[i448e] in devtools, provide program that prints code with Gauss-Hermite...

[i448e] in devtools, provide program that prints code with Gauss-Hermite quadrature nodes and weights ()

Merging branch 'i448e'  into 'main'.

See merge request !1331

devtools/tabulate/GaussHermiteConstants.cpp

+/* Nodes and weight for Gauss-Hermite quadrature.
+   Computed using the program GaussHermite_NodesWeights.cpp
+   by John Burkardt and Joachim Wuttke [1], based on algorithms by
+   Roger Martin and James Wilkinson [2] and Sylvan Elhay, Jaroslav Kautsky [3]
+
+   [1] https://jugit.fz-juelich.de/mlz/bornagain/-/tree/main/devtools/tabulate/GaussHermite_NodesWeights.cpp
+   [2] The Implicit QL Algorithm, Numerische Mathematik, Volume 12, Number 5, December 1968, pages 377-383.
+   [3] Algorithm 655: IQPACK, FORTRAN Subroutines for the Weights of Interpolatory Quadrature, ACM Transactions on Mathematical Software, Volume 13, Number 4, December 1987, pages 399-415.
+*/
+static const int maxOrderGaussHermite = 25;
+static const double GaussHermiteNodesWeights[650] = {
+               /* order 1 */ 0, 1, 
+
+/* order 2 */ -0.707106781186547, 0.5, 0.707106781186547, 0.5, 
+
+/* order 3 */ -1.22474487139159, 0.166666666666667, 0, 0.666666666666667, 1.22474487139159, 0.166666666666667, 
+
+/* order 4 */ -1.65068012388578, 0.0458758547680685, -0.52464762327529, 0.454124145231931, 0.52464762327529, 0.454124145231931, 1.65068012388578, 0.0458758547680685, 
+
+/* order 5 */ -2.02018287045609, 0.0112574113277207, -0.958572464613819, 0.222075922005613, 0, 0.533333333333334, 0.958572464613819, 0.222075922005613, 2.02018287045609, 0.0112574113277207, 
+
+/* order 6 */ -2.35060497367449, 0.00255578440205624, -1.3358490740137, 0.0886157460419146, -0.436077411927616, 0.408828469556029, 0.436077411927616, 0.408828469556029, 1.3358490740137, 0.0886157460419146, 2.35060497367449, 0.00255578440205624, 
+
+/* order 7 */ -2.65196135683523, 0.000548268855972217, -1.67355162876747, 0.0307571239675865, -0.816287882858965, 0.240123178605013, 0, 0.457142857142858, 0.816287882858965, 0.240123178605013, 1.67355162876747, 0.0307571239675865, 2.65196135683523, 0.000548268855972217, 
+
+/* order 8 */ -2.93063742025724, 0.000112614538375368, -1.98165675669584, 0.00963522012078825, -1.15719371244678, 0.117239907661759, -0.381186990207322, 0.373012257679077, 0.381186990207322, 0.373012257679077, 1.15719371244678, 0.117239907661759, 1.98165675669584, 0.00963522012078825, 2.93063742025724, 0.000112614538375368, 
+
+/* order 9 */ -3.19099320178153, 2.23458440077466e-05, -2.26658058453184, 0.00278914132123177, -1.46855328921667, 0.049916406765218, -0.723551018752838, 0.24409750289494, 0, 0.406349206349207, 0.723551018752838, 0.24409750289494, 1.46855328921667, 0.049916406765218, 2.26658058453184, 0.00278914132123177, 3.19099320178153, 2.23458440077466e-05, 
+
+/* order 10 */ -3.43615911883774, 4.31065263071829e-06, -2.53273167423279, 0.000758070934312214, -1.75668364929988, 0.0191115805007703, -1.03661082978951, 0.135483702980268, -0.342901327223705, 0.344642334932019, 0.342901327223705, 0.344642334932019, 1.03661082978951, 0.135483702980268, 1.75668364929988, 0.0191115805007703, 2.53273167423279, 0.000758070934312214, 3.43615911883774, 4.31065263071829e-06, 
+
+/* order 11 */ -3.66847084655958, 8.12184979021493e-07, -2.78329009978165, 0.000195671930271223, -2.02594801582576, 0.00672028523553728, -1.32655708449493, 0.0661387460710579, -0.6568095668821, 0.24224029987397, 0, 0.36940836940837, 0.6568095668821, 0.24224029987397, 1.32655708449493, 0.0661387460710579, 2.02594801582576, 0.00672028523553728, 2.78329009978165, 0.000195671930271223, 3.66847084655958, 8.12184979021493e-07, 
+
+/* order 12 */ -3.88972489786978, 1.49992716763717e-07, -3.02063702512089, 4.83718492259062e-05, -2.27950708050106, 0.00220338068753321, -1.5976826351526, 0.0291166879123642, -0.947788391240164, 0.14696704804533, -0.314240376254359, 0.32166436151283, 0.314240376254359, 0.32166436151283, 0.947788391240164, 0.14696704804533, 1.5976826351526, 0.0291166879123642, 2.27950708050106, 0.00220338068753321, 3.02063702512089, 4.83718492259062e-05, 3.88972489786978, 1.49992716763717e-07, 
+
+/* order 13 */ -4.10133759617864, 2.72262764280589e-08, -3.24660897837241, 1.15265965273339e-05, -2.51973568567824, 0.000681236350442929, -1.85310765160151, 0.0117705605059965, -1.22005503659075, 0.0791689558604499, -0.60576387917106, 0.237871522964136, 0, 0.340992340992341, 0.60576387917106, 0.237871522964136, 1.22005503659075, 0.0791689558604499, 1.85310765160151, 0.0117705605059965, 2.51973568567824, 0.000681236350442929, 3.24660897837241, 1.15265965273339e-05, 4.10133759617864, 2.72262764280589e-08, 
+
+/* order 14 */ -4.30444857047363, 4.86816125774837e-09, -3.46265693360227, 2.66099134406762e-06, -2.7484707249854, 0.000200339553760744, -2.09518325850772, 0.00442891910694742, -1.47668273114114, 0.0386501088242534, -0.878713787329399, 0.154083339842513, -0.291745510672562, 0.30263462681302, 0.291745510672562, 0.30263462681302, 0.878713787329399, 0.154083339842513, 1.47668273114114, 0.0386501088242534, 2.09518325850772, 0.00442891910694742, 2.7484707249854, 0.000200339553760744, 3.46265693360227, 2.66099134406762e-06, 4.30444857047363, 4.86816125774837e-09, 
+
+/* order 15 */ -4.49999070730939, 8.58964989963332e-10, -3.66995037340445, 5.97541959792063e-07, -2.9671669279056, 5.64214640518902e-05, -2.32573248617386, 0.00156735750354995, -1.71999257518649, 0.0173657744921376, -1.13611558521092, 0.0894177953998445, -0.565069583255576, 0.232462293609732, 0, 0.318259518259518, 0.565069583255576, 0.232462293609732, 1.13611558521092, 0.0894177953998445, 1.71999257518649, 0.0173657744921376, 2.32573248617386, 0.00156735750354995, 2.9671669279056, 5.64214640518902e-05, 3.66995037340445, 5.97541959792063e-07, 4.49999070730939, 8.58964989963332e-10, 
+
+/* order 16 */ -4.68873893930582, 1.49781472316184e-10, -3.86944790486012, 1.30947321628682e-07, -3.17699916197996, 1.53000321624873e-05, -2.54620215784748, 0.000525984926573909, -1.95178799091625, 0.00726693760118475, -1.38025853919888, 0.0472847523540141, -0.822951449144656, 0.15833837275095, -0.273481046138152, 0.286568521238013, 0.273481046138152, 0.286568521238013, 0.822951449144656, 0.15833837275095, 1.38025853919888, 0.0472847523540141, 1.95178799091625, 0.00726693760118475, 2.54620215784748, 0.000525984926573909, 3.17699916197996, 1.53000321624873e-05, 3.86944790486012, 1.30947321628682e-07, 4.68873893930582, 1.49781472316184e-10, 
+
+/* order 17 */ -4.8713451936744, 2.58431491937492e-11, -4.06194667587548, 2.80801611793057e-08, -3.37893209114149, 4.01267944797983e-06, -2.75776291570389, 0.000168491431551338, -2.17350282666662, 0.00285894606228462, -1.61292431422123, 0.023086657025711, -1.06764872574345, 0.0974063711627213, -0.531633001342655, 0.226706308468977, 0, 0.299538370126606, 0.531633001342655, 0.226706308468977, 1.06764872574345, 0.0974063711627213, 1.61292431422123, 0.023086657025711, 2.17350282666662, 0.00285894606228462, 2.75776291570389, 0.000168491431551338, 3.37893209114149, 4.01267944797983e-06, 4.06194667587548, 2.80801611793057e-08, 4.8713451936744, 2.58431491937492e-11, 
+
+/* order 18 */ -5.04836400887447, 4.4165887693587e-12, -4.24811787356813, 5.90548847883656e-09, -3.57376906848627, 1.02155239763698e-06, -2.96137750553161, 5.17989614411622e-05, -2.38629908916669, 0.00106548479629165, -1.83553160426163, 0.0105165177519413, -1.30092085838962, 0.0548966324802226, -0.776682919267412, 0.160685303893513, -0.258267750519097, 0.272783234654288, 0.258267750519097, 0.272783234654288, 0.776682919267412, 0.160685303893513, 1.30092085838962, 0.0548966324802226, 1.83553160426163, 0.0105165177519413, 2.38629908916669, 0.00106548479629165, 2.96137750553161, 5.17989614411622e-05, 3.57376906848627, 1.02155239763698e-06, 4.24811787356813, 5.90548847883656e-09, 5.04836400887447, 4.4165887693587e-12, 
+
+/* order 19 */ -5.22027169053748, 7.48283005405713e-13, -4.42853280660378, 1.22037084844749e-09, -3.76218735196402, 2.5322200320929e-07, -3.1578488183476, 1.53511459546667e-05, -2.59113378979454, 0.000378502109414269, -2.04923170985062, 0.00450723542034204, -1.52417061939353, 0.0286666910301185, -1.01036838713431, 0.103603657276144, -0.503520163423888, 0.220941712199144, 0, 0.283773192751521, 0.503520163423888, 0.220941712199144, 1.01036838713431, 0.103603657276144, 1.52417061939353, 0.0286666910301185, 2.04923170985062, 0.00450723542034204, 2.59113378979454, 0.000378502109414269, 3.1578488183476, 1.53511459546667e-05, 3.76218735196402, 2.5322200320929e-07, 4.42853280660378, 1.22037084844749e-09, 5.22027169053748, 7.48283005405713e-13, 
+
+/* order 20 */ -5.38748089001124, 1.25780067243789e-13, -4.60368244955074, 2.48206236231517e-10, -3.94476404011562, 6.12749025998293e-08, -3.34785456738321, 4.40212109023087e-06, -2.78880605842813, 0.00012882627996193, -2.25497400208927, 0.0018301031310805, -1.73853771211659, 0.013997837447101, -1.23407621539532, 0.0615063720639769, -0.737473728545395, 0.161739333984, -0.245340708300901, 0.260793063449555, 0.245340708300901, 0.260793063449555, 0.737473728545395, 0.161739333984, 1.23407621539532, 0.0615063720639769, 1.73853771211659, 0.013997837447101, 2.25497400208927, 0.0018301031310805, 2.78880605842813, 0.00012882627996193, 3.34785456738321, 4.40212109023087e-06, 3.94476404011562, 6.12749025998293e-08, 4.60368244955074, 2.48206236231517e-10, 5.38748089001124, 1.25780067243789e-13, 
+
+/* order 21 */ -5.55035187326468, 2.09899121956567e-14, -4.77399234341122, 4.97536860412163e-11, -4.12199554749184, 1.4506612844931e-08, -3.53197287713768, 1.22535483614826e-06, -2.9799912077046, 4.21923474255168e-05, -2.45355212451284, 0.000708047795481538, -1.94496294918625, 0.00643969705140878, -1.44893425065073, 0.0339527297865429, -0.961499634418369, 0.108392285626419, -0.479450707079108, 0.21533371569506, 0, 0.270260183572876, 0.479450707079108, 0.21533371569506, 0.961499634418369, 0.108392285626419, 1.44893425065073, 0.0339527297865429, 1.94496294918625, 0.00643969705140878, 2.45355212451284, 0.000708047795481538, 2.9799912077046, 4.21923474255168e-05, 3.53197287713768, 1.22535483614826e-06, 4.12199554749184, 1.4506612844931e-08, 4.77399234341122, 4.97536860412163e-11, 5.55035187326468, 2.09899121956567e-14, 
+
+/* order 22 */ -5.70920135320527, 3.4794606478771e-15, -4.93983413106017, 9.84137898234611e-12, -4.29431248059316, 3.36651415945817e-09, -3.7107015328778, 3.31985374981405e-07, -3.16526590920214, 1.33459771268095e-05, -2.64563744105817, 0.000262283303255964, -2.14423359279854, 0.00280876104757721, -1.65587437328642, 0.0175690728808057, -1.17671395848124, 0.0671963114288902, -0.703686097170007, 0.161906293413675, -0.234179139930991, 0.250243596586935, 0.234179139930991, 0.250243596586935, 0.703686097170007, 0.161906293413675, 1.17671395848124, 0.0671963114288902, 1.65587437328642, 0.0175690728808057, 2.14423359279854, 0.00280876104757721, 2.64563744105817, 0.000262283303255964, 3.16526590920214, 1.33459771268095e-05, 3.7107015328778, 3.31985374981405e-07, 4.29431248059316, 3.36651415945817e-09, 4.93983413106017, 9.84137898234611e-12, 5.70920135320527, 3.4794606478771e-15, 
+
+/* order 23 */ -5.86430949898457, 5.73238316780196e-16, -5.10153461047668, 1.92293531156779e-12, -4.46209117374, 7.67088886239997e-10, -3.8844727081061, 8.77506248386196e-08, -3.34512715994122, 4.08997724499214e-06, -2.83180378712616, 9.34081860903131e-05, -2.33701621147445, 0.00116762863749786, -1.85567703767137, 0.00857967839146562, -1.38403958568249, 0.0388671837034811, -0.919151465442564, 0.112073382602621, -0.458538350068105, 0.209959669577543, 0, 0.258509740808839, 0.458538350068105, 0.209959669577543, 0.919151465442564, 0.112073382602621, 1.38403958568249, 0.0388671837034811, 1.85567703767137, 0.00857967839146562, 2.33701621147445, 0.00116762863749786, 2.83180378712616, 9.34081860903131e-05, 3.34512715994122, 4.08997724499214e-06, 3.8844727081061, 8.77506248386196e-08, 4.46209117374, 7.67088886239997e-10, 5.10153461047668, 1.92293531156779e-12, 5.86430949898457, 5.73238316780196e-16, 
+
+/* order 24 */ -6.01592556142574, 9.39019368904176e-17, -5.25938292766804, 3.71497415276242e-13, -4.62566275642379, 1.71866492796486e-10, -4.05366440244815, 2.26746167348067e-08, -3.52000681303453, 1.21765974544258e-06, -3.01254613756556, 3.20950056527461e-05, -2.52388101701143, 0.000464718718779398, -2.0490035736617, 0.00397660892918131, -1.58425001096169, 0.0211263444089675, -1.12676081761125, 0.0720693640171782, -0.674171107037212, 0.161459512867, -0.224414547472515, 0.24087011554664, 0.224414547472515, 0.24087011554664, 0.674171107037212, 0.161459512867, 1.12676081761125, 0.0720693640171782, 1.58425001096169, 0.0211263444089675, 2.0490035736617, 0.00397660892918131, 2.52388101701143, 0.000464718718779398, 3.01254613756556, 3.20950056527461e-05, 3.52000681303453, 1.21765974544258e-06, 4.05366440244815, 2.26746167348067e-08, 4.62566275642379, 1.71866492796486e-10, 5.25938292766804, 3.71497415276242e-13, 6.01592556142574, 9.39019368904176e-17, 
+
+/* order 25 */ -6.16427243405245, 1.53003899799871e-17, -5.41363635528004, 7.10210303700393e-14, -4.78532036735222, 3.79115000047724e-11, -4.21860944438656, 5.73802386889943e-09, -3.69028287699836, 3.5301525602455e-07, -3.1882949244251, 1.06721949052025e-05, -2.70532023717303, 0.000177766906926527, -2.23642013026728, 0.0017578504052638, -1.77800112433715, 0.0108567559914623, -1.32728070207308, 0.0433799701676449, -0.881982756213821, 0.114880924303952, -0.440147298645308, 0.204851025650341, 0, 0.248169351176485, 0.440147298645308, 0.204851025650341, 0.881982756213821, 0.114880924303952, 1.32728070207308, 0.0433799701676449, 1.77800112433715, 0.0108567559914623, 2.23642013026728, 0.0017578504052638, 2.70532023717303, 0.000177766906926527, 3.1882949244251, 1.06721949052025e-05, 3.69028287699836, 3.5301525602455e-07, 4.21860944438656, 5.73802386889943e-09, 4.78532036735222, 3.79115000047724e-11, 5.41363635528004, 7.10210303700393e-14, 6.16427243405245, 1.53003899799871e-17, 
+
+};

devtools/tabulate/GaussHermite_Demo.cpp

+// Demonstrates usage of the Gauss-Hermite nodes and weights
+// generated by our program GaussHermite_NodesWeights.cpp.
+//
+// Joachim Wuttke, 2023, public domain
+
+#include "GaussHermiteConstants.cpp"
+#include <cmath>
+#include <iostream>
+#include <iomanip>
+
+// With f(x)=cos(x), analytical integration yields exp(-1/4).
+// See https://www.chebfun.org/examples/quad/HermiteQuad.html, with different normalization.
+double f(double x)
+{
+    return cos(x);
+}
+
+int main(int argc, char* argv[])
+{
+    std::cout << std::setprecision(15);
+    for (int n=0; n<12; ++n) {
+
+        double sum = 0;
+        for (int i = 0; i<n; ++i) {
+            int k = n * (n-1) + 2*i;
+            sum += GaussHermiteNodesWeights[k+1] * f(GaussHermiteNodesWeights[k]);
+        }
+        std::cout << sum << "\n";
+    }
+}