GaussLegendreRule.hpp

// This file is part of Bembel, the higher order C++ boundary element library.
//
// Copyright (C) 2022 see <http://www.bembel.eu>
//
// It was written as part of a cooperation of J. Doelz, H. Harbrecht, S. Kurz,
// M. Multerer, S. Schoeps, and F. Wolf at Technische Universitaet Darmstadt,
// Universitaet Basel, and Universita della Svizzera italiana, Lugano. This
// source code is subject to the GNU General Public License version 3 and
// provided WITHOUT ANY WARRANTY, see <http://www.bembel.eu> for further
// information.
 
#ifndef BEMBEL_SRC_QUADRATURE_GAUSSLEGENDRERULE_HPP_
#define BEMBEL_SRC_QUADRATURE_GAUSSLEGENDRERULE_HPP_
 
namespace Bembel {
template <unsigned int Order>
struct GaussLegendreRule {
  GaussLegendreRule() {
    Eigen::MatrixXd A(Order, Order);
    A.setZero();
    for (auto i = 1; i < Order; ++i) A(i, i - 1) = i / (sqrt(4. * i * i - 1));
    Eigen::SelfAdjointEigenSolver<Eigen::MatrixXd> es;
    es.compute(A);
    xi_.resize(Order);
    w_.resize(Order);
    for (auto i = 0; i < Order; ++i) {
      xi_[i] = 0.5 * (es.eigenvalues()(i) + 1);
      w_[i] = es.eigenvectors()(0, i) * es.eigenvectors()(0, i);
    }
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<0> {
  GaussLegendreRule() {
    assert(false && "Order of GaussLegendre Quadrature has to be at least 1");
  }
};
template <>
struct GaussLegendreRule<1> {
  GaussLegendreRule() {
    xi_ = {0.5};
    w_ = {1.0};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<2> {
  GaussLegendreRule() {
    xi_ = {0.2113248654051871177454256097490212721762,
           0.7886751345948128822545743902509787278238};
    w_ = {0.5, 0.5};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<3> {
  GaussLegendreRule() {
    xi_ = {0.1127016653792583114820734600217600389167, 0.5,
           0.8872983346207416885179265399782399610833};
 
    w_ = {.2777777777777777777777777777777777777778,
          .4444444444444444444444444444444444444444,
          .2777777777777777777777777777777777777778};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<4> {
  GaussLegendreRule() {
    xi_ = {0.06943184420297371238802675555359524745214,
           0.3300094782075718675986671204483776563997,
           0.6699905217924281324013328795516223436003,
           0.9305681557970262876119732444464047525479};
 
    w_ = {.1739274225687269286865319746109997036177,
          .3260725774312730713134680253890002963823,
          .3260725774312730713134680253890002963823,
          .1739274225687269286865319746109997036177};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<5> {
  GaussLegendreRule() {
    xi_ = {0.04691007703066800360118656085030351743717,
           0.2307653449471584544818427896498955975164, 0.5,
           0.7692346550528415455181572103501044024836,
           0.9530899229693319963988134391496964825628};
 
    w_ = {.1184634425280945437571320203599586813216,
          .2393143352496832340206457574178190964561,
          .2844444444444444444444444444444444444444,
          .2393143352496832340206457574178190964561,
          .1184634425280945437571320203599586813216};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<6> {
  GaussLegendreRule() {
    xi_ = {0.03376524289842398609384922275300269543262,
           0.1693953067668677431693002024900473264968,
           0.3806904069584015456847491391596440322907,
           0.6193095930415984543152508608403559677093,
           0.8306046932331322568306997975099526735032,
           0.9662347571015760139061507772469973045674};
 
    w_ = {.08566224618958517252014807108636644676341,
          .1803807865240693037849167569188580558308,
          .2339569672863455236949351719947754974058,
          .2339569672863455236949351719947754974058,
          .1803807865240693037849167569188580558308,
          .08566224618958517252014807108636644676341};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<7> {
  GaussLegendreRule() {
    xi_ = {0.02544604382862073773690515797607436879961,
           0.1292344072003027800680676133596057964629,
           0.2970774243113014165466967939615192683263,
           0.5,
           0.7029225756886985834533032060384807316737,
           0.8707655927996972199319323866403942035371,
           0.9745539561713792622630948420239256312004};
 
    w_ = {0.06474248308443484663530571633954100916429,
          0.1398526957446383339507338857118897912435,
          0.1909150252525594724751848877444875669392,
          0.2089795918367346938775510204081632653061,
          0.1909150252525594724751848877444875669392,
          0.1398526957446383339507338857118897912435,
          0.06474248308443484663530571633954100916429};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<8> {
  GaussLegendreRule() {
    xi_ = {0.01985507175123188415821956571526350478588,
           0.1016667612931866302042230317620847815814,
           0.237233795041835507091130475405376825479,
           0.4082826787521750975302619288199080096666,
           0.5917173212478249024697380711800919903334,
           0.762766204958164492908869524594623174521,
           0.8983332387068133697957769682379152184186,
           0.9801449282487681158417804342847364952141};
 
    w_ = {0.0506142681451881295762656771549810950577,
          0.1111905172266872352721779972131204422151,
          0.1568533229389436436689811009933006566302,
          0.1813418916891809914825752246385978060971,
          0.1813418916891809914825752246385978060971,
          0.1568533229389436436689811009933006566302,
          0.1111905172266872352721779972131204422151,
          0.0506142681451881295762656771549810950577};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<9> {
  GaussLegendreRule() {
    xi_ = {0.0159198802461869550822118985481635649753,
           0.08198444633668210285028510596513256172795,
           0.1933142836497048013456489803292629076071,
           0.337873288298095535480730992678331695714,
           0.5,
           0.662126711701904464519269007321668304286,
           0.8066857163502951986543510196707370923929,
           0.9180155536633178971497148940348674382721,
           0.9840801197538130449177881014518364350247};
 
    w_ = {.04063719418078720598594607905526182533783,
          .09032408034742870202923601562145640475717,
          .1303053482014677311593714347093164248859,
          .1561735385200014200343152032922218327994,
          .1651196775006298815822625346434870244394,
          .1561735385200014200343152032922218327994,
          .1303053482014677311593714347093164248859,
          .09032408034742870202923601562145640475717,
          .04063719418078720598594607905526182533783};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<10> {
  GaussLegendreRule() {
    xi_ = {0.01304673574141413996101799395777397328587,
           0.06746831665550774463395165578825347573623,
           0.1602952158504877968828363174425632121154,
           0.2833023029353764046003670284171079189,
           0.4255628305091843945575869994351400076912,
           0.5744371694908156054424130005648599923088,
           0.7166976970646235953996329715828920811,
           0.8397047841495122031171636825574367878846,
           0.9325316833444922553660483442117465242638,
           0.9869532642585858600389820060422260267141};
 
    w_ = {0.03333567215434406879678440494666589642893,
          0.07472567457529029657288816982884866620128,
          0.1095431812579910219977674671140815962294,
          0.1346333596549981775456134607847346764299,
          0.1477621123573764350869464973256691647105,
          0.1477621123573764350869464973256691647105,
          0.1346333596549981775456134607847346764299,
          0.1095431812579910219977674671140815962294,
          0.07472567457529029657288816982884866620128,
          0.03333567215434406879678440494666589642893};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<11> {
  GaussLegendreRule() {
    xi_ = {0.01088567092697150359803099943857130461429,
           0.05646870011595235046242111534803636668416,
           0.1349239972129753379532918739844232709752,
           0.2404519353965940920371371652706952227599,
           0.3652284220238275138342340072995692376602,
           0.5,
           0.6347715779761724861657659927004307623398,
           0.7595480646034059079628628347293047772401,
           0.8650760027870246620467081260155767290248,
           0.9435312998840476495375788846519636333158,
           0.9891143290730284964019690005614286953857};
 
    w_ = {0.02783428355808683324137686022127428936426,
          0.06279018473245231231734714961197005009881,
          0.09314510546386712571304882071582794584564,
          0.1165968822959952399592618524215875697159,
          0.1314022722551233310903444349452545976864,
          0.136462543388950315357241764168171094578,
          0.1314022722551233310903444349452545976864,
          0.1165968822959952399592618524215875697159,
          0.09314510546386712571304882071582794584564,
          0.06279018473245231231734714961197005009881,
          0.02783428355808683324137686022127428936426};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<12> {
  GaussLegendreRule() {
    xi_ = {0.009219682876640374654725454925359588519922,
           0.0479413718147625716607670669404519037312,
           0.1150486629028476564815530833935909620075,
           0.2063410228566912763516487905297328598155,
           0.3160842505009099031236542316781412193718,
           0.4373832957442655422637793152680734350083,
           0.5626167042557344577362206847319265649917,
           0.6839157494990900968763457683218587806282,
           0.7936589771433087236483512094702671401845,
           0.8849513370971523435184469166064090379925,
           0.9520586281852374283392329330595480962688,
           0.9907803171233596253452745450746404114801};
 
    w_ = {0.02358766819325591359730798074250853015851,
          0.05346966299765921548012735909699811210729,
          0.08003916427167311316732626477167953593601,
          0.1015837133615329608745322279048991882533,
          0.1167462682691774043804249494624390281297,
          0.1245735229067013925002812180214756054152,
          0.1245735229067013925002812180214756054152,
          0.1167462682691774043804249494624390281297,
          0.1015837133615329608745322279048991882533,
          0.08003916427167311316732626477167953593601,
          0.05346966299765921548012735909699811210729,
          0.02358766819325591359730798074250853015851};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<13> {
  GaussLegendreRule() {
    xi_ = {0.007908472640705925263585275596445194467505,
           0.04120080038851101739672608174964024380476,
           0.09921095463334504360289675520857005484719,
           0.1788253302798298896780076965022421749642,
           0.275753624481776573561043573936180066099,
           0.3847708420224326029672359394510055823942,
           0.5,
           0.6152291579775673970327640605489944176058,
           0.724246375518223426438956426063819933901,
           0.8211746697201701103219923034977578250358,
           0.9007890453666549563971032447914299451528,
           0.9587991996114889826032739182503597561952,
           0.9920915273592940747364147244035548055325};
 
    w_ = {0.02024200238265793976001079610049303002099,
          0.04606074991886422395721088797689856046184,
          0.06943675510989361923180088843443573381093,
          0.08907299038097286914002334599804899775641,
          0.1039080237684442511562616096530263816933,
          0.1131415901314486192060450930198883092174,
          0.1162757766154369550972947576344179740783,
          0.1131415901314486192060450930198883092174,
          0.1039080237684442511562616096530263816933,
          0.08907299038097286914002334599804899775641,
          0.06943675510989361923180088843443573381093,
          0.04606074991886422395721088797689856046184,
          0.02024200238265793976001079610049303002099};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<14> {
  GaussLegendreRule() {
    xi_ = {0.006858095651593830579201366647973599161954,
           0.03578255816821324133180443031106286776148,
           0.08639934246511750340510262867480251948015,
           0.1563535475941572649259900984903329312308,
           0.2423756818209229540173546407244056688456,
           0.3404438155360551197821640879157622665829,
           0.445972525646328168966877674890082626194,
           0.554027474353671831033122325109917373806,
           0.6595561844639448802178359120842377334171,
           0.7576243181790770459826453592755943311544,
           0.8436464524058427350740099015096670687692,
           0.9136006575348824965948973713251974805199,
           0.9642174418317867586681955696889371322385,
           0.993141904348406169420798633352026400838};
 
    w_ = {0.01755973016587593151591643806909589030985,
          0.04007904357988010490281663853142715479185,
          0.06075928534395159234470740453623831297833,
          0.07860158357909676728480096931192107830283,
          0.09276919873896890687085829506257851812446,
          0.1025992318606478019829620328306090278552,
          0.1076319267315788950979382216581300176375,
          0.1076319267315788950979382216581300176375,
          0.1025992318606478019829620328306090278552,
          0.09276919873896890687085829506257851812446,
          0.07860158357909676728480096931192107830283,
          0.06075928534395159234470740453623831297833,
          0.04007904357988010490281663853142715479185,
          0.01755973016587593151591643806909589030985};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<15> {
  GaussLegendreRule() {
    xi_ = {0.006003740989757285755217140706693709426514,
           0.031363303799647047846120526144895264378,
           0.07589670829478639189967583961289157431687,
           0.1377911343199149762919069726930309951846,
           0.2145139136957305762313866313730446793808,
           0.3029243264612183150513963145094772658186,
           0.3994029530012827388496858483027018960936,
           0.5,
           0.6005970469987172611503141516972981039064,
           0.6970756735387816849486036854905227341814,
           0.7854860863042694237686133686269553206192,
           0.8622088656800850237080930273069690048154,
           0.9241032917052136081003241603871084256831,
           0.968636696200352952153879473855104735622,
           0.9939962590102427142447828592933062905735};
 
    w_ = {0.01537662099805863417731419678860220886087,
          0.03518302374405406235463370822533366923335,
          0.05357961023358596750593477334293465170777,
          0.06978533896307715722390239725551416126043,
          0.08313460290849696677660043024060440556545,
          0.09308050000778110551340028093321141225311,
          0.09921574266355578822805916322191966240935,
          0.1012891209627806364403100999837596574193,
          0.09921574266355578822805916322191966240935,
          0.09308050000778110551340028093321141225311,
          0.08313460290849696677660043024060440556545,
          0.06978533896307715722390239725551416126043,
          0.05357961023358596750593477334293465170777,
          0.03518302374405406235463370822533366923335,
          0.01537662099805863417731419678860220886087};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<16> {
  GaussLegendreRule() {
    xi_ = {0.005299532504175033701922913274833686286863,
           0.02771248846338371196100579223269582745443,
           0.06718439880608412805976605114380343380633,
           0.1222977958224984830524494025762788658231,
           0.1910618777986781257766641179756044905041,
           0.27099161117138630682879027850821121323,
           0.359198224610370543384769749269751946757,
           0.4524937450811812799073403322875209684348,
           0.5475062549188187200926596677124790315652,
           0.640801775389629456615230250730248053243,
           0.72900838882861369317120972149178878677,
           0.8089381222013218742233358820243955094959,
           0.8777022041775015169475505974237211341769,
           0.9328156011939158719402339488561965661937,
           0.9722875115366162880389942077673041725456,
           0.9947004674958249662980770867251663137131};
 
    w_ = {0.01357622970587704742589028622800905175613,
          0.03112676196932394643142191849718884713749,
          0.04757925584124639240496255380112311317763,
          0.06231448562776693602623814109600821007244,
          0.07479799440828836604075086527373927448525,
          0.08457825969750126909465603951517998110582,
          0.09130170752246179443338183398460996969178,
          0.09472530522753424814269836160414155257345,
          0.09472530522753424814269836160414155257345,
          0.09130170752246179443338183398460996969178,
          0.08457825969750126909465603951517998110582,
          0.07479799440828836604075086527373927448525,
          0.06231448562776693602623814109600821007244,
          0.04757925584124639240496255380112311317763,
          0.03112676196932394643142191849718884713749,
          0.01357622970587704742589028622800905175613};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<17> {
  GaussLegendreRule() {
    xi_ = {0.004712262342791332162282990029667361746105,
           0.02466223911561611938864152105209848927831,
           0.05988042313650704893852215275592215368829,
           0.109242998051599296537384972239761974888,
           0.1711644203916546170748488916784988324261,
           0.2436547314567615160568767156852240627085,
           0.3243841182730618423514072414523269974797,
           0.4107579092520760720746612531729672212623,
           0.5,
           0.5892420907479239279253387468270327787377,
           0.6756158817269381576485927585476730025203,
           0.7563452685432384839431232843147759372915,
           0.8288355796083453829251511083215011675739,
           0.890757001948400703462615027760238025112,
           0.9401195768634929510614778472440778463117,
           0.9753377608843838806113584789479015107217,
           0.9952877376572086678377170099703326382539};
 
    w_ = {0.01207415143427396598005501314378266234585,
          0.02772976468699360056472008267912233025642,
          0.04251807415858959044176768509553103692525,
          0.05594192359670198554739419281317796336792,
          0.06756818423426273664315999085117509868606,
          0.07702288053840514404071579740097930597024,
          0.0840020510782250222549853318941615775106,
          0.08828135268349632316263549505659861957546,
          0.08972323517810326272913282213094281072439,
          0.08828135268349632316263549505659861957546,
          0.0840020510782250222549853318941615775106,
          0.07702288053840514404071579740097930597024,
          0.06756818423426273664315999085117509868606,
          0.05594192359670198554739419281317796336792,
          0.04251807415858959044176768509553103692525,
          0.02772976468699360056472008267912233025642,
          0.01207415143427396598005501314378266234585};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<18> {
  GaussLegendreRule() {
    xi_ = {0.004217415789534526634991997646924614873711,
           0.02208802521430112240940205353511184501358,
           0.05369876675122213039696970443642724229605,
           0.09814752051373844215879127249270460144835,
           0.1541564784698233960625544593555758052739,
           0.2201145844630262326960642257373354315362,
           0.2941244192685786769820341030834741814605,
           0.3740568871542472452055135725610443849186,
           0.4576124934791323493788690735321080941333,
           0.5423875065208676506211309264678919058667,
           0.6259431128457527547944864274389556150814,
           0.7058755807314213230179658969165258185395,
           0.7798854155369737673039357742626645684638,
           0.8458435215301766039374455406444241947261,
           0.9018524794862615578412087275072953985516,
           0.9463012332487778696030302955635727577039,
           0.9779119747856988775905979464648881549864,
           0.9957825842104654733650080023530753851263};
 
    w_ = {0.01080800676324165515667135513322623469384,
          0.02485727444748489822666747310131932084043,
          0.03821286512744452826456483880831826280266,
          0.0504710220531435827814069924624173035314,
          0.06127760335573923009225956340010077761408,
          0.07032145733532532560236565187597361404775,
          0.07734233756313262246270900191818738609661,
          0.08213824187291636149302688823296379520617,
          0.08457119148157179592032823506749330516705,
          0.08457119148157179592032823506749330516705,
          0.08213824187291636149302688823296379520617,
          0.07734233756313262246270900191818738609661,
          0.07032145733532532560236565187597361404775,
          0.06127760335573923009225956340010077761408,
          0.0504710220531435827814069924624173035314,
          0.03821286512744452826456483880831826280266,
          0.02485727444748489822666747310131932084043,
          0.01080800676324165515667135513322623469384};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<19> {
  GaussLegendreRule() {
    xi_ = {0.003796578078207798405491164873369753205342,
           0.01989592393258498457361057965617423669245,
           0.0484220481925910491786695357338437560953,
           0.08864267173142858751053875664364304911273,
           0.1395169113323853106914520695881091851714,
           0.1997273476691594882651809175268803600658,
           0.2677146293120195271413664259479488160119,
           0.3417179500181850840049413355750775410539,
           0.4198206771798873120659519421296282252476,
           0.5,
           0.5801793228201126879340480578703717747524,
           0.6582820499818149159950586644249224589461,
           0.7322853706879804728586335740520511839881,
           0.8002726523308405117348190824731196399342,
           0.8604830886676146893085479304118908148286,
           0.9113573282685714124894612433563569508873,
           0.9515779518074089508213304642661562439047,
           0.9801040760674150154263894203438257633075,
           0.9962034219217922015945088351266302467947};
 
    w_ = {0.009730894114863238518156020732219217876453,
          0.02240711338284980016641907870099710597588,
          0.03452227136882061329035412900300652248092,
          0.04574501081122499973223104706191982633046,
          0.05578332277366699735801195084088299874067,
          0.06437698126966811383775789242843855852792,
          0.07130335108680330588787305472095148623783,
          0.0763830210329298333894277004488314992305,
          0.07948442169697717382497821973252360083939,
          0.08052722492439184798958181266045836751995,
          0.07948442169697717382497821973252360083939,
          0.0763830210329298333894277004488314992305,
          0.07130335108680330588787305472095148623783,
          0.06437698126966811383775789242843855852792,
          0.05578332277366699735801195084088299874067,
          0.04574501081122499973223104706191982633046,
          0.03452227136882061329035412900300652248092,
          0.02240711338284980016641907870099710597588,
          0.009730894114863238518156020732219217876453};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
template <>
struct GaussLegendreRule<20> {
  GaussLegendreRule() {
    xi_ = {0.003435700407452537606938805764339860888676,
           0.01801403636104310436616693440136138904397,
           0.04388278587433704706612377939835094347541,
           0.08044151408889058830273546914923965733519,
           0.1268340467699246036928474648221792048446,
           0.1819731596367424872735816518868570316283,
           0.2445664990245864509978179745223745007873,
           0.3131469556422902196637259114875363813021,
           0.3861070744291774609597519023157126876285,
           0.4617367394332513331226797953005808944976,
           0.5382632605667486668773202046994191055024,
           0.6138929255708225390402480976842873123715,
           0.6868530443577097803362740885124636186979,
           0.7554335009754135490021820254776254992127,
           0.8180268403632575127264183481131429683717,
           0.8731659532300753963071525351778207951554,
           0.9195584859111094116972645308507603426648,
           0.9561172141256629529338762206016490565246,
           0.981985963638956895633833065598638610956,
           0.9965642995925474623930611942356601391113};
 
    w_ = {0.008807003569576059155930981175926408181072,
          0.02030071490019347066551997613746605493955,
          0.0313360241670545317847532675935208031758,
          0.04163837078835237436237907161102310305009,
          0.05096505990862021751837506774017493808335,
          0.05909726598075920865618868885569114350252,
          0.06584431922458831344924724987408156745806,
          0.07104805465919102566464916253358246651726,
          0.07458649323630187339391436850098471834634,
          0.07637669356536292534904216597754879674597,
          0.07637669356536292534904216597754879674597,
          0.07458649323630187339391436850098471834634,
          0.07104805465919102566464916253358246651726,
          0.06584431922458831344924724987408156745806,
          0.05909726598075920865618868885569114350252,
          0.05096505990862021751837506774017493808335,
          0.04163837078835237436237907161102310305009,
          0.0313360241670545317847532675935208031758,
          0.02030071490019347066551997613746605493955,
          0.008807003569576059155930981175926408181072};
  }
  std::vector<double> xi_;
  std::vector<double> w_;
};
 
}  // namespace Bembel
#endif  // BEMBEL_SRC_QUADRATURE_GAUSSLEGENDRERULE_HPP_