76 std::vector<long double> points =
77 Polynomials::jacobi_polynomial_roots<long double>(n, 0, 0);
79 for (
unsigned int i = 0; i < (points.size() + 1) / 2; ++i)
85 const long double pp =
88 const long double x = -1. + 2. * points[i];
89 const double w = 1. / ((1. - x * x) * pp * pp);
106 long double result = n - 1;
107 for (
int i = n - 2; i > 1; --i)
121 std::vector<long double>
126 const unsigned int q = x.size();
127 std::vector<long double>
w(q);
129 const long double factor =
130 std::pow(2., alpha + beta + 1) *
gamma(alpha + q) *
gamma(beta + q) /
131 ((q - 1) *
gamma(q) *
gamma(alpha + beta + q + 1));
132 for (
unsigned int i = 0; i < q; ++i)
134 const long double s =
136 w[i] = factor / (s * s);
139 w[q - 1] *= (alpha + 1);
163 std::vector<double> q_points(n);
171 q_points[1] = 2. / 3.;
175 q_points[1] = (6. - std::sqrt(6)) * 0.1;
176 q_points[2] = (6. + std::sqrt(6)) * 0.1;
180 q_points[0] = 0.000000000000000000;
181 q_points[1] = 0.212340538239152943;
182 q_points[2] = 0.590533135559265343;
183 q_points[3] = 0.911412040487296071;
186 q_points[0] = 0.000000000000000000;
187 q_points[1] = 0.139759864343780571;
188 q_points[2] = 0.416409567631083166;
189 q_points[3] = 0.723156986361876197;
190 q_points[4] = 0.942895803885482331;
193 q_points[0] = 0.000000000000000000;
194 q_points[1] = 0.098535085798826416;
195 q_points[2] = 0.304535726646363913;
196 q_points[3] = 0.562025189752613841;
197 q_points[4] = 0.801986582126391845;
198 q_points[5] = 0.960190142948531222;
201 q_points[0] = 0.000000000000000000;
202 q_points[1] = 0.073054328680258851;
203 q_points[2] = 0.230766137969945495;
204 q_points[3] = 0.441328481228449865;
205 q_points[4] = 0.663015309718845702;
206 q_points[5] = 0.851921400331515644;
207 q_points[6] = 0.970683572840215114;
210 q_points[0] = 0.000000000000000000;
211 q_points[1] = 0.056262560536922135;
212 q_points[2] = 0.180240691736892389;
213 q_points[3] = 0.352624717113169672;
214 q_points[4] = 0.547153626330555420;
215 q_points[5] = 0.734210177215410598;
216 q_points[6] = 0.885320946839095790;
217 q_points[7] = 0.977520613561287499;
228 const ::QGaussRadau<1>::EndPoint end_point)
233 case ::QGaussRadau<1>::EndPoint::left:
235 case ::QGaussRadau<1>::EndPoint::right:
237 std::vector<double> points(n);
238 for (
unsigned int i = 0; i < n; ++i)
240 points[n - i - 1] = 1. - left_points[i];
248 "This constructor can only be called with either "
249 "QGaussRadau::left or QGaussRadau::right as second argument."));
265 std::vector<double> weights(n);
276 weights[0] = 1. / 9.;
277 weights[1] = (16. + std::sqrt(6)) / 36.;
278 weights[2] = (16. - std::sqrt(6)) / 36.;
281 weights[0] = 0.062500000000000000;
282 weights[1] = 0.328844319980059696;
283 weights[2] = 0.388193468843171852;
284 weights[3] = 0.220462211176768369;
287 weights[0] = 0.040000000000000001;
288 weights[1] = 0.223103901083570894;
289 weights[2] = 0.311826522975741427;
290 weights[3] = 0.281356015149462124;
291 weights[4] = 0.143713560791225797;
294 weights[0] = 0.027777777777777776;
295 weights[1] = 0.159820376610255471;
296 weights[2] = 0.242693594234484888;
297 weights[3] = 0.260463391594787597;
298 weights[4] = 0.208450667155953895;
299 weights[5] = 0.100794192626740456;
302 weights[0] = 0.020408163265306121;
303 weights[1] = 0.119613744612656100;
304 weights[2] = 0.190474936822115581;
305 weights[3] = 0.223554914507283209;
306 weights[4] = 0.212351889502977870;
307 weights[5] = 0.159102115733650767;
308 weights[6] = 0.074494235556010341;
311 weights[0] = 0.015625000000000000;
312 weights[1] = 0.092679077401489660;
313 weights[2] = 0.152065310323392683;
314 weights[3] = 0.188258772694559262;
315 weights[4] = 0.195786083726246729;
316 weights[5] = 0.173507397817250691;
317 weights[6] = 0.124823950664932445;
318 weights[7] = 0.057254407372128648;
330 const ::QGaussRadau<1>::EndPoint end_point)
335 case ::QGaussRadau<1>::EndPoint::left:
337 case ::QGaussRadau<1>::EndPoint::right:
339 std::vector<double> weights(n);
340 for (
unsigned int i = 0; i < n; ++i)
342 weights[n - i - 1] = left_weights[i];
349 "This constructor can only be called with either "
350 "QGaussRadau::EndPoint::left or "
351 "QGaussRadau::EndPoint::right as second argument."));
362 , end_point(end_point)
365 std::vector<double> p =
367 std::vector<double>
w =
370 for (
unsigned int i = 0; i < this->size(); ++i)
373 this->weights[i] =
w[i];
386 std::vector<long double> points =
387 Polynomials::jacobi_polynomial_roots<long double>(n - 2, 1, 1);
388 points.insert(points.begin(), 0);
389 points.push_back(1.);
390 std::vector<long double>
w =
394 for (
unsigned int i = 0; i < points.size(); ++i)
397 this->weights[i] = 0.5 *
w[i];
417 static const double xpts[] = {0.0, 1.0};
418 static const double wts[] = {0.5, 0.5};
420 for (
unsigned int i = 0; i < this->
size(); ++i)
433 static const double xpts[] = {0.0, 0.5, 1.0};
434 static const double wts[] = {1. / 6., 2. / 3., 1. / 6.};
436 for (
unsigned int i = 0; i < this->
size(); ++i)
449 static const double xpts[] = {0.0, .25, .5, .75, 1.0};
450 static const double wts[] = {
451 7. / 90., 32. / 90., 12. / 90., 32. / 90., 7. / 90.};
453 for (
unsigned int i = 0; i < this->
size(); ++i)
466 static const double xpts[] = {
467 0.0, 1. / 6., 1. / 3., .5, 2. / 3., 5. / 6., 1.0};
468 static const double wts[] = {41. / 840.,
476 for (
unsigned int i = 0; i < this->
size(); ++i)
491 for (
unsigned int i = 0; i < this->
size(); ++i)
498 this->
weights[i] = revert ?
w[n - 1 - i] :
w[i];
506 std::vector<double> q_points(n);
511 q_points[0] = 0.3333333333333333;
515 q_points[0] = 0.1120088061669761;
516 q_points[1] = 0.6022769081187381;
520 q_points[0] = 0.06389079308732544;
521 q_points[1] = 0.3689970637156184;
522 q_points[2] = 0.766880303938942;
526 q_points[0] = 0.04144848019938324;
527 q_points[1] = 0.2452749143206022;
528 q_points[2] = 0.5561654535602751;
529 q_points[3] = 0.848982394532986;
533 q_points[0] = 0.02913447215197205;
534 q_points[1] = 0.1739772133208974;
535 q_points[2] = 0.4117025202849029;
536 q_points[3] = 0.6773141745828183;
537 q_points[4] = 0.89477136103101;
541 q_points[0] = 0.02163400584411693;
542 q_points[1] = 0.1295833911549506;
543 q_points[2] = 0.3140204499147661;
544 q_points[3] = 0.5386572173517997;
545 q_points[4] = 0.7569153373774084;
546 q_points[5] = 0.922668851372116;
551 q_points[0] = 0.0167193554082585;
552 q_points[1] = 0.100185677915675;
553 q_points[2] = 0.2462942462079286;
554 q_points[3] = 0.4334634932570557;
555 q_points[4] = 0.6323509880476823;
556 q_points[5] = 0.81111862674023;
557 q_points[6] = 0.940848166743287;
561 q_points[0] = 0.01332024416089244;
562 q_points[1] = 0.07975042901389491;
563 q_points[2] = 0.1978710293261864;
564 q_points[3] = 0.354153994351925;
565 q_points[4] = 0.5294585752348643;
566 q_points[5] = 0.7018145299391673;
567 q_points[6] = 0.849379320441094;
568 q_points[7] = 0.953326450056343;
572 q_points[0] = 0.01086933608417545;
573 q_points[1] = 0.06498366633800794;
574 q_points[2] = 0.1622293980238825;
575 q_points[3] = 0.2937499039716641;
576 q_points[4] = 0.4466318819056009;
577 q_points[5] = 0.6054816627755208;
578 q_points[6] = 0.7541101371585467;
579 q_points[7] = 0.877265828834263;
580 q_points[8] = 0.96225055941096;
584 q_points[0] = 0.00904263096219963;
585 q_points[1] = 0.05397126622250072;
586 q_points[2] = 0.1353118246392511;
587 q_points[3] = 0.2470524162871565;
588 q_points[4] = 0.3802125396092744;
589 q_points[5] = 0.5237923179723384;
590 q_points[6] = 0.6657752055148032;
591 q_points[7] = 0.7941904160147613;
592 q_points[8] = 0.898161091216429;
593 q_points[9] = 0.9688479887196;
598 q_points[0] = 0.007643941174637681;
599 q_points[1] = 0.04554182825657903;
600 q_points[2] = 0.1145222974551244;
601 q_points[3] = 0.2103785812270227;
602 q_points[4] = 0.3266955532217897;
603 q_points[5] = 0.4554532469286375;
604 q_points[6] = 0.5876483563573721;
605 q_points[7] = 0.7139638500230458;
606 q_points[8] = 0.825453217777127;
607 q_points[9] = 0.914193921640008;
608 q_points[10] = 0.973860256264123;
612 q_points[0] = 0.006548722279080035;
613 q_points[1] = 0.03894680956045022;
614 q_points[2] = 0.0981502631060046;
615 q_points[3] = 0.1811385815906331;
616 q_points[4] = 0.2832200676673157;
617 q_points[5] = 0.398434435164983;
618 q_points[6] = 0.5199526267791299;
619 q_points[7] = 0.6405109167754819;
620 q_points[8] = 0.7528650118926111;
621 q_points[9] = 0.850240024421055;
622 q_points[10] = 0.926749682988251;
623 q_points[11] = 0.977756129778486;
639 std::vector<double> quadrature_weights(n);
644 quadrature_weights[0] = -1.0;
647 quadrature_weights[0] = -0.7185393190303845;
648 quadrature_weights[1] = -0.2814606809696154;
652 quadrature_weights[0] = -0.5134045522323634;
653 quadrature_weights[1] = -0.3919800412014877;
654 quadrature_weights[2] = -0.0946154065661483;
658 quadrature_weights[0] = -0.3834640681451353;
659 quadrature_weights[1] = -0.3868753177747627;
660 quadrature_weights[2] = -0.1904351269501432;
661 quadrature_weights[3] = -0.03922548712995894;
665 quadrature_weights[0] = -0.2978934717828955;
666 quadrature_weights[1] = -0.3497762265132236;
667 quadrature_weights[2] = -0.234488290044052;
668 quadrature_weights[3] = -0.0989304595166356;
669 quadrature_weights[4] = -0.01891155214319462;
673 quadrature_weights[0] = -0.2387636625785478;
674 quadrature_weights[1] = -0.3082865732739458;
675 quadrature_weights[2] = -0.2453174265632108;
676 quadrature_weights[3] = -0.1420087565664786;
677 quadrature_weights[4] = -0.05545462232488041;
678 quadrature_weights[5] = -0.01016895869293513;
683 quadrature_weights[0] = -0.1961693894252476;
684 quadrature_weights[1] = -0.2703026442472726;
685 quadrature_weights[2] = -0.239681873007687;
686 quadrature_weights[3] = -0.1657757748104267;
687 quadrature_weights[4] = -0.0889432271377365;
688 quadrature_weights[5] = -0.03319430435645653;
689 quadrature_weights[6] = -0.005932787015162054;
693 quadrature_weights[0] = -0.164416604728002;
694 quadrature_weights[1] = -0.2375256100233057;
695 quadrature_weights[2] = -0.2268419844319134;
696 quadrature_weights[3] = -0.1757540790060772;
697 quadrature_weights[4] = -0.1129240302467932;
698 quadrature_weights[5] = -0.05787221071771947;
699 quadrature_weights[6] = -0.02097907374214317;
700 quadrature_weights[7] = -0.003686407104036044;
704 quadrature_weights[0] = -0.1400684387481339;
705 quadrature_weights[1] = -0.2097722052010308;
706 quadrature_weights[2] = -0.211427149896601;
707 quadrature_weights[3] = -0.1771562339380667;
708 quadrature_weights[4] = -0.1277992280331758;
709 quadrature_weights[5] = -0.07847890261203835;
710 quadrature_weights[6] = -0.0390225049841783;
711 quadrature_weights[7] = -0.01386729555074604;
712 quadrature_weights[8] = -0.002408041036090773;
716 quadrature_weights[0] = -0.12095513195457;
717 quadrature_weights[1] = -0.1863635425640733;
718 quadrature_weights[2] = -0.1956608732777627;
719 quadrature_weights[3] = -0.1735771421828997;
720 quadrature_weights[4] = -0.135695672995467;
721 quadrature_weights[5] = -0.0936467585378491;
722 quadrature_weights[6] = -0.05578772735275126;
723 quadrature_weights[7] = -0.02715981089692378;
724 quadrature_weights[8] = -0.00951518260454442;
725 quadrature_weights[9] = -0.001638157633217673;
730 quadrature_weights[0] = -0.1056522560990997;
731 quadrature_weights[1] = -0.1665716806006314;
732 quadrature_weights[2] = -0.1805632182877528;
733 quadrature_weights[3] = -0.1672787367737502;
734 quadrature_weights[4] = -0.1386970574017174;
735 quadrature_weights[5] = -0.1038334333650771;
736 quadrature_weights[6] = -0.06953669788988512;
737 quadrature_weights[7] = -0.04054160079499477;
738 quadrature_weights[8] = -0.01943540249522013;
739 quadrature_weights[9] = -0.006737429326043388;
740 quadrature_weights[10] = -0.001152486965101561;
744 quadrature_weights[0] = -0.09319269144393;
745 quadrature_weights[1] = -0.1497518275763289;
746 quadrature_weights[2] = -0.166557454364573;
747 quadrature_weights[3] = -0.1596335594369941;
748 quadrature_weights[4] = -0.1384248318647479;
749 quadrature_weights[5] = -0.1100165706360573;
750 quadrature_weights[6] = -0.07996182177673273;
751 quadrature_weights[7] = -0.0524069547809709;
752 quadrature_weights[8] = -0.03007108900074863;
753 quadrature_weights[9] = -0.01424924540252916;
754 quadrature_weights[10] = -0.004899924710875609;
755 quadrature_weights[11] = -0.000834029009809656;
763 return quadrature_weights;
771 const bool factor_out_singularity)
773 ((origin[0] == 0) || (origin[0] == 1)) ? (alpha == 1 ? n : 2 * n) : 4 * n)
774 , fraction(((origin[0] == 0) || (origin[0] == 1.)) ? 1. : origin[0])
799 unsigned int ns_offset = (
fraction == 1) ? n : 2 * n;
801 for (
unsigned int i = 0, j = ns_offset; i < n; ++i, ++j)
809 if ((alpha != 1) || (
fraction != 1))
829 if (factor_out_singularity ==
true)
830 for (
unsigned int i = 0; i <
size(); ++i)
835 "The singularity cannot be on a Gauss point of the same order!"));
838 Assert(denominator != 0.0,
840 "The quadrature formula you are using does not allow to "
841 "factor out the singularity, which is zero at one point."));
842 this->
weights[i] /= denominator;
851 const double eps = 1
e-8;
852 bool on_edge =
false;
853 for (
unsigned int d = 0;
d < 2; ++
d)
854 on_edge = on_edge || (std::abs(singularity[
d]) <
eps ||
855 std::abs(singularity[
d] - 1.0) <
eps);
856 const bool on_vertex =
858 std::abs((singularity -
Point<2>(.5, .5)).norm_square() - .5) <
eps;
870 const bool factor_out_singularity)
879 std::vector<QGaussOneOverR<2>> quads;
880 std::vector<Point<2>> origins;
883 quads.emplace_back(n, 3, factor_out_singularity);
884 quads.emplace_back(n, 2, factor_out_singularity);
885 quads.emplace_back(n, 1, factor_out_singularity);
886 quads.emplace_back(n, 0, factor_out_singularity);
888 origins.emplace_back(0., 0.);
889 origins.emplace_back(singularity[0], 0.);
890 origins.emplace_back(0., singularity[1]);
891 origins.push_back(singularity);
896 unsigned int q_id = 0;
899 for (
unsigned int box = 0; box < 4; ++box)
902 dist =
Point<2>(std::abs(dist[0]), std::abs(dist[1]));
903 double area = dist[0] * dist[1];
905 for (
unsigned int q = 0; q < quads[box].size(); ++q, ++q_id)
907 const Point<2> &qp = quads[box].point(q);
909 origins[box] +
Point<2>(dist[0] * qp[0], dist[1] * qp[1]);
910 this->
weights[q_id] = quads[box].weight(q) * area;
918 const unsigned int vertex_index,
919 const bool factor_out_singularity)
955 std::vector<double> &ws = this->
weights;
958 for (
unsigned int q = 0; q < gauss.
size(); ++q)
962 ps[q][1] = gp[0] * std::tan(pi4 * gp[1]);
963 ws[q] = gauss.
weight(q) * pi4 / std::cos(pi4 * gp[1]);
964 if (factor_out_singularity)
968 ws[gauss.
size() + q] = ws[q];
969 ps[gauss.
size() + q][0] = ps[q][1];
970 ps[gauss.
size() + q][1] = ps[q][0];
975 switch (vertex_index)
992 double R00 = std::cos(theta), R01 = -std::sin(theta);
993 double R10 = std::sin(theta), R11 = std::cos(theta);
995 if (vertex_index != 0)
996 for (
unsigned int q = 0; q <
size(); ++q)
998 double x = ps[q][0] - .5, y = ps[q][1] - .5;
1000 ps[q][0] = R00 * x + R01 * y + .5;
1001 ps[q][1] = R10 * x + R11 * y + .5;
1010 std::vector<unsigned int> permutation(quad.
size());
1011 for (
unsigned int i = 0; i < quad.
size(); ++i)
1014 std::sort(permutation.begin(),
1016 [
this](
const unsigned int x,
const unsigned int y) {
1017 return this->compare_weights(x, y);
1027 for (
unsigned int i = 0; i < quad.
size(); ++i)
1031 if (permutation[i] != i)
1041 return (this->weights[a] < this->weights[
b]);
1059 , end_point(end_point)
1136 const double eta_bar = singularity[0] * 2. - 1.;
1137 const double eta_star = eta_bar * eta_bar - 1.;
1141 std::vector<double> weights_dummy(
weights.size());
1142 unsigned int cont = 0;
1143 const double tol = 1
e-10;
1161 weights.resize(weights_dummy.size() - 1);
1169 if (std::abs(eta_star) <= tol)
1172 std::pow((eta_bar * eta_star + std::abs(eta_star)), 1.0 / 3.0) +
1173 std::pow((eta_bar * eta_star - std::abs(eta_star)), 1.0 / 3.0) +
1178 gamma_bar = (eta_bar * eta_star + std::abs(eta_star)) /
1179 std::abs(eta_bar * eta_star + std::abs(eta_star)) *
1180 std::pow(std::abs(eta_bar * eta_star + std::abs(eta_star)),
1182 (eta_bar * eta_star - std::abs(eta_star)) /
1183 std::abs(eta_bar * eta_star - std::abs(eta_star)) *
1184 std::pow(std::abs(eta_bar * eta_star - std::abs(eta_star)),
1191 double eta = (Utilities::fixed_power<3>(
gamma - gamma_bar) +
1192 gamma_bar * (gamma_bar * gamma_bar + 3)) /
1193 (1 + 3 * gamma_bar * gamma_bar);
1195 double J = 3 * ((
gamma - gamma_bar) * (
gamma - gamma_bar)) /
1196 (1 + 3 * gamma_bar * gamma_bar);
1213 std::vector<double> points(n);
1215 for (
unsigned short i = 0; i < n; ++i)
1222 (1. +
double(2 * i + 1) /
double(2 * (n - 1) + 2))));
1235 std::vector<double> weights(n);
1237 for (
unsigned short i = 0; i < n; ++i)
1253 Assert(n > 0,
ExcMessage(
"Need at least one point for the quadrature rule"));
1257 for (
unsigned int i = 0; i < this->
size(); ++i)
1278 const unsigned int n,
1279 const ::QGaussRadauChebyshev<1>::EndPoint end_point)
1281 std::vector<double> points(n);
1283 for (
unsigned short i = 0; i < n; ++i)
1289 case ::QGaussRadauChebyshev<1>::EndPoint::left:
1295 (1 + 2 *
double(i) / (2 *
double(n - 1) + 1.))));
1299 case ::QGaussRadauChebyshev<1>::EndPoint::right:
1303 (1. - std::cos(
numbers::PI * (2 *
double(n - 1 - i) /
1304 (2 *
double(n - 1) + 1.))));
1312 "This constructor can only be called with either "
1313 "QGaussRadauChebyshev::EndPoint::left or "
1314 "QGaussRadauChebyshev::EndPoint:right as second argument."));
1325 const unsigned int n,
1326 const ::QGaussRadauChebyshev<1>::EndPoint end_point)
1328 std::vector<double> weights(n);
1330 for (
unsigned short i = 0; i < n; ++i)
1337 else if (end_point ==
1353 , end_point(end_point)
1355 Assert(n > 0,
ExcMessage(
"Need at least one point for quadrature rules."));
1356 std::vector<double> points =
1361 for (
unsigned int i = 0; i < this->
size(); ++i)
1364 this->weights[i] =
weights[i];
1375 , end_point(end_point)
1388 std::vector<double> points(n);
1390 for (
unsigned short i = 0; i < n; ++i)
1396 (1. + std::cos(
numbers::PI * (1 +
double(i) /
double(n - 1))));
1405 std::vector<double> weights(n);
1407 for (
unsigned short i = 0; i < n; ++i)
1411 if (i == 0 || i == (n - 1))
1428 "Need at least two points for Gauss-Lobatto quadrature rule"));
1429 std::vector<double> p =
1431 std::vector<double>
w =
1434 for (
unsigned int i = 0; i < this->
size(); ++i)
1452 std::vector<Point<dim>> qpoints;
1453 std::vector<double> weights;
1455 for (
unsigned int i = 0; i < quad.
size(); ++i)
1462 for (
int d = 0;
d < dim; ++
d)
1467 this->weights.push_back(quad.
weight(i));
1475 template <
int spacedim>
1481 ExcMessage(
"Invalid combination of dim and spacedim ."));
1483 for (
unsigned int d = 0;
d < dim; ++
d)
1487 const double J = std::abs(B.determinant());
1493 std::vector<Point<spacedim>> qp(this->size());
1494 std::vector<double>
w(this->size());
1496 for (
unsigned int i = 0; i < this->size(); ++i)
1500 w[i] = J * this->weight(i);
1509 template <
int spacedim>
1512 const std::vector<std::array<
Point<spacedim>, dim + 1>> &simplices)
const
1514 Assert(!(dim == 1 && spacedim == 1),
1515 ExcMessage(
"This function is not supposed to work in 1D-1d case."));
1517 ExcMessage(
"Invalid combination of dim and spacedim ."));
1519 std::vector<Point<spacedim>> qp;
1520 std::vector<double> ws;
1521 for (
const auto &
simplex : simplices)
1523 const auto rule = this->compute_affine_transformation(
simplex);
1525 rule.get_points().end(),
1526 std::back_inserter(qp),
1529 rule.get_weights().end(),
1530 std::back_inserter(ws),
1531 [&](
const double w) { return w; });
1545 for (
unsigned int i = 0; i < base.
size(); ++i)
1547 const auto &q = base.
point(i);
1548 const auto w = base.
weight(i);
1550 const auto xhat = q[0];
1551 const auto yhat = q[1];
1555 const double st = std::sin(t);
1556 const double ct = std::cos(t);
1557 const double r = xhat / (st + ct);
1559 const double J = pi * xhat / (2 * (std::sin(pi * yhat) + 1));
1582 for (
unsigned int i = 0; i < base.
size(); ++i)
1584 const auto &q = base.
point(i);
1585 const auto w = base.
weight(i);
1587 const auto xhat = q[0];
1588 const auto yhat = q[1];
1590 const double x = std::pow(xhat, beta) * (1 - yhat);
1591 const double y = std::pow(xhat, beta) * yhat;
1593 const double J = beta * std::pow(xhat, 2. * beta - 1.);
1613 "The split point should be inside the unit reference cell."));
1615 std::array<Point<dim>, dim + 1>
vertices;
1623 for (
unsigned int start = 0; start < (dim > 2 ? 2 : 1); ++start)
1625 for (
unsigned int i = 0; i < dim; ++i)
1632 quad.get_points().begin(),
1633 quad.get_points().end());
1634 this->weights.insert(this->weights.end(),
1635 quad.get_weights().begin(),
1636 quad.get_weights().end());
1648 if (dim == 0 || dim == 1)
1650 const ::QGauss<dim> quad(n_points_1D);
1653 this->
weights = quad.get_weights();
1657 if (n_points_1D == 1)
1659 const double p = 1.0 / 3.0;
1661 this->
weights.emplace_back(0.5);
1663 else if (n_points_1D == 2)
1669 const double Q12 = 1.0 / 2.0;
1671 0.1550510257216822);
1673 0.6449489742783178);
1675 0.1550510257216822);
1677 0.6449489742783178);
1679 this->
weights.emplace_back(0.31804138174397717 * Q12);
1680 this->
weights.emplace_back(0.18195861825602283 * Q12);
1681 this->
weights.emplace_back(0.31804138174397717 * Q12);
1682 this->
weights.emplace_back(0.18195861825602283 * Q12);
1684 else if (n_points_1D == 3)
1687 const double p0 = 2.0 / 7.0 - std::sqrt(15.0) / 21.0;
1688 const double p1 = 2.0 / 7.0 + std::sqrt(15.0) / 21.0;
1689 const double p2 = 3.0 / 7.0 - 2.0 * std::sqrt(15.0) / 21.0;
1690 const double p3 = 3.0 / 7.0 + 2.0 * std::sqrt(15.0) / 21.0;
1699 const double q12 = 0.5;
1700 const double w0 = 9.0 / 40.0;
1701 const double w1 = 31.0 / 240.0 - std::sqrt(15.0) / 1200.0;
1702 const double w2 = 31.0 / 240.0 + std::sqrt(15.0) / 1200.0;
1703 this->
weights.emplace_back(q12 * w0);
1704 this->
weights.emplace_back(q12 * w1);
1705 this->
weights.emplace_back(q12 * w1);
1706 this->
weights.emplace_back(q12 * w1);
1707 this->
weights.emplace_back(q12 * w2);
1708 this->
weights.emplace_back(q12 * w2);
1709 this->
weights.emplace_back(q12 * w2);
1711 else if (n_points_1D == 4)
1719 if (n_points_1D == 1)
1721 const double Q14 = 1.0 / 4.0;
1722 const double Q16 = 1.0 / 6.0;
1725 this->
weights.emplace_back(Q16);
1732 else if (n_points_1D == 2)
1734 const double Q16 = 1.0 / 6.0;
1735 this->
weights.emplace_back(0.1223220027573451 * Q16);
1736 this->
weights.emplace_back(0.1280664127107469 * Q16);
1737 this->
weights.emplace_back(0.1325680271444452 * Q16);
1738 this->
weights.emplace_back(0.1406244096604032 * Q16);
1739 this->
weights.emplace_back(0.2244151669175574 * Q16);
1740 this->
weights.emplace_back(0.2520039808095023 * Q16);
1744 0.01271836631368145);
1747 0.3621268299455338);
1749 0.01140332944455717,
1750 0.3586207204668839);
1753 0.6384932999617267);
1756 0.1308471689520965);
1759 0.1403728057942107);
1763 else if (n_points_1D == 3)
1768 else if (n_points_1D == 4)
1778 "QGaussSimplex is currently only implemented for "
1779 "n_points_1D = 1, 2, 3, and 4 while you are asking for "
1786 template <std::
size_t b_dim>
1787 std::vector<std::array<double, b_dim>>
1788 all_permutations(
const std::array<double, b_dim> &b_point)
1790 std::vector<std::array<double, b_dim>> output;
1795 std::array<double, b_dim> temp = b_point;
1796 std::sort(temp.begin(), temp.end());
1799 output.push_back(temp);
1801 while (std::next_permutation(temp.begin(), temp.end()));
1811 const unsigned int n_points_1D,
1812 const bool use_odd_order)
1824 std::array<double, dim + 1> centroid;
1825 std::fill(centroid.begin(), centroid.end(), 1.0 / (dim + 1.0));
1826 std::vector<std::vector<std::array<double, dim + 1>>> b_point_permutations;
1827 std::vector<double> b_weights;
1837 auto process_point_1 = [&](
const double a,
const double w) {
1838 const double b = 1.0 - dim * a;
1839 std::array<double, dim + 1> b_point;
1840 std::fill(b_point.begin(), b_point.begin() + dim, a);
1843 b_weights.push_back(
w);
1844 b_point_permutations.push_back(all_permutations(b_point));
1849 auto process_point_2 = [&](
const double a,
const double w) {
1851 const double b = (1.0 - 2.0 * a) / 2.0;
1852 std::array<double, dim + 1> b_point;
1853 std::fill(b_point.begin(), b_point.begin() + dim - 1, a);
1854 b_point[dim - 1] =
b;
1857 b_weights.push_back(
w);
1858 b_point_permutations.push_back(all_permutations(b_point));
1864 auto process_point_3 = [&](
const double a,
const double b,
const double w) {
1865 const double c = 1.0 - (dim - 1.0) * a -
b;
1866 std::array<double, dim + 1> b_point;
1867 std::fill(b_point.begin(), b_point.begin() + dim - 1, a);
1868 b_point[dim - 1] =
b;
1871 b_weights.push_back(
w);
1872 b_point_permutations.push_back(all_permutations(b_point));
1875 switch (n_points_1D)
1884 b_point_permutations.push_back({centroid});
1885 b_weights.push_back(1.0000000000000000e+00);
1890 process_point_1(1.6666666666666669e-01,
1891 3.3333333333333331e-01);
1898 b_point_permutations.push_back({centroid});
1899 b_weights.push_back(1.0000000000000000e+00);
1904 process_point_1(1.3819660112501050e-01,
1905 2.5000000000000000e-01);
1917 process_point_1(9.1576213509770743e-02, 1.0995174365532187e-01);
1918 process_point_1(4.4594849091596489e-01, 2.2338158967801147e-01);
1924 process_point_1(3.2816330251638171e-01,
1925 1.3621784253708741e-01);
1926 process_point_1(1.0804724989842859e-01,
1927 1.1378215746291261e-01);
1946 b_point_permutations.push_back({centroid});
1947 b_weights.push_back(2.2500000000000001e-01);
1948 process_point_1(1.0128650732345634e-01,
1949 1.2593918054482714e-01);
1950 process_point_1(4.7014206410511511e-01,
1951 1.3239415278850619e-01);
1956 process_point_1(6.3089014491502227e-02,
1957 5.0844906370206819e-02);
1958 process_point_1(2.4928674517091043e-01,
1959 1.1678627572637937e-01);
1960 process_point_3(5.3145049844816938e-02,
1961 3.1035245103378439e-01,
1962 8.2851075618373571e-02);
1969 process_point_1(3.1088591926330061e-01,
1970 1.1268792571801590e-01);
1971 process_point_1(9.2735250310891248e-02,
1972 7.3493043116361956e-02);
1973 process_point_2(4.5503704125649642e-02,
1974 4.2546020777081472e-02);
1979 process_point_1(4.0673958534611372e-02,
1980 1.0077211055320640e-02);
1981 process_point_1(3.2233789014227548e-01,
1982 5.5357181543654717e-02);
1983 process_point_1(2.1460287125915201e-01,
1984 3.9922750258167487e-02);
1985 process_point_3(6.3661001875017442e-02,
1986 6.0300566479164919e-01,
1987 4.8214285714285710e-02);
2001 process_point_1(3.3730648554587850e-02,
2002 1.6545050110792131e-02);
2003 process_point_1(4.7430969250471822e-01,
2004 7.7086646185986069e-02);
2005 process_point_1(2.4157738259540357e-01,
2006 1.2794417123015558e-01);
2007 process_point_3(4.7036644652595216e-02,
2008 1.9868331479735168e-01,
2009 5.5878732903199779e-02);
2014 b_point_permutations.push_back({centroid});
2015 b_weights.push_back(1.4431560767778717e-01);
2016 process_point_1(5.0547228317030957e-02,
2017 3.2458497623198079e-02);
2018 process_point_1(4.5929258829272313e-01,
2019 9.5091634267284619e-02);
2020 process_point_1(1.7056930775176021e-01,
2021 1.0321737053471824e-01);
2022 process_point_3(8.3947774099575878e-03,
2023 2.6311282963463811e-01,
2024 2.7230314174434993e-02);
2031 b_point_permutations.push_back({centroid});
2032 b_weights.push_back(9.5485289464130846e-02);
2033 process_point_1(3.1570114977820279e-01,
2034 4.2329581209967028e-02);
2035 process_point_2(5.0489822598396350e-02,
2036 3.1896927832857580e-02);
2037 process_point_3(1.8883383102600099e-01,
2038 5.7517163758699996e-01,
2039 3.7207130728334620e-02);
2040 process_point_3(2.1265472541483140e-02,
2041 8.1083024109854862e-01,
2042 8.1107708299033420e-03);
2047 process_point_1(1.0795272496221089e-01,
2048 2.6426650908408830e-02);
2049 process_point_1(1.8510948778258660e-01,
2050 5.2031747563738531e-02);
2051 process_point_1(4.2316543684767283e-02,
2052 7.5252561535401989e-03);
2053 process_point_1(3.1418170912403898e-01,
2054 4.1763782856934897e-02);
2055 process_point_2(4.3559132858383021e-01,
2056 3.6280930261308818e-02);
2057 process_point_3(2.1433930127130570e-02,
2058 7.1746406342630831e-01,
2059 7.1569028908444327e-03);
2060 process_point_3(2.0413933387602909e-01,
2061 5.8379737830214440e-01,
2062 1.5453486150960340e-02);
2076 b_point_permutations.push_back({centroid});
2077 b_weights.push_back(9.7135796282798836e-02);
2078 process_point_1(4.4729513394452691e-02,
2079 2.5577675658698031e-02);
2080 process_point_1(4.8968251919873762e-01,
2081 3.1334700227139071e-02);
2082 process_point_1(4.3708959149293664e-01,
2083 7.7827541004774278e-02);
2084 process_point_1(1.8820353561903275e-01,
2085 7.9647738927210249e-02);
2086 process_point_3(3.6838412054736258e-02,
2087 2.2196298916076568e-01,
2088 4.3283539377289376e-02);
2093 b_point_permutations.push_back({centroid});
2094 b_weights.push_back(8.1743329146285973e-02);
2095 process_point_1(3.2055373216943517e-02,
2096 1.3352968813149567e-02);
2097 process_point_1(1.4216110105656438e-01,
2098 4.5957963604744731e-02);
2099 process_point_3(2.8367665339938453e-02,
2100 1.6370173373718250e-01,
2101 2.5297757707288385e-02);
2102 process_point_3(2.9619889488729734e-02,
2103 3.6914678182781102e-01,
2104 3.4184648162959429e-02);
2105 process_point_3(1.4813288578382056e-01,
2106 3.2181299528883545e-01,
2107 6.3904906396424044e-02);
2114 b_point_permutations.push_back({centroid});
2115 b_weights.push_back(5.8010548912480253e-02);
2116 process_point_1(6.1981697552226933e-10,
2117 6.4319281759256394e-05);
2118 process_point_1(1.6077453539526160e-01,
2119 2.3173338462425461e-02);
2120 process_point_1(3.2227652182142102e-01,
2121 2.9562912335429289e-02);
2122 process_point_1(4.5108918345413578e-02,
2123 8.0639799796161822e-03);
2124 process_point_2(1.1229654600437609e-01,
2125 3.8134080103702457e-02);
2126 process_point_3(4.5887144875245922e-01,
2127 2.5545792330413102e-03,
2128 8.3844221982985519e-03);
2129 process_point_3(3.3775870685338598e-02,
2130 7.1835032644207453e-01,
2131 1.0234559352745330e-02);
2132 process_point_3(1.8364136980992790e-01,
2133 3.4415910578175279e-02,
2134 2.0524915967988139e-02);
2139 b_point_permutations.push_back({centroid});
2140 b_weights.push_back(4.7399773556020743e-02);
2141 process_point_1(3.1225006869518868e-01,
2142 2.6937059992268701e-02);
2143 process_point_1(1.1430965385734609e-01,
2144 9.8691597167933822e-03);
2145 process_point_3(4.1043073921896539e-01,
2146 1.6548602561961109e-01,
2147 1.1393881220195230e-02);
2148 process_point_3(6.1380088247906528e-03,
2149 9.4298876734520487e-01,
2150 3.6194434433925362e-04);
2151 process_point_3(1.2105018114558939e-01,
2152 4.7719037990428043e-01,
2153 2.5739731980456069e-02);
2154 process_point_3(3.2779468216442620e-02,
2155 5.9425626948000698e-01,
2156 1.0135871679755789e-02);
2157 process_point_3(3.2485281564823047e-02,
2158 8.0117728465834437e-01,
2159 6.5761472770359038e-03);
2160 process_point_3(1.7497934218393901e-01,
2161 6.2807184547536599e-01,
2162 1.2907035798861989e-02);
2175 b_point_permutations.push_back({centroid});
2176 b_weights.push_back(8.5761179732224219e-02);
2177 process_point_1(2.8485417614371900e-02, 1.0431870512894697e-02);
2178 process_point_1(4.9589190096589092e-01, 1.6606273054585369e-02);
2179 process_point_1(1.0263548271224643e-01, 3.8630759237019321e-02);
2180 process_point_1(4.3846592676435220e-01, 6.7316154079468296e-02);
2181 process_point_1(2.1021995670317828e-01, 7.0515684111716576e-02);
2182 process_point_3(7.3254276860644785e-03,
2183 1.4932478865208237e-01,
2184 1.0290289572953278e-02);
2185 process_point_3(4.6010500165429957e-02,
2186 2.8958112563770588e-01,
2187 4.0332476640500554e-02);
2192 process_point_1(2.4646363436335583e-02, 7.9316425099736389e-03);
2193 process_point_1(4.8820375094554153e-01, 2.4266838081452032e-02);
2194 process_point_1(1.0925782765935427e-01, 2.8486052068877544e-02);
2195 process_point_1(4.4011164865859309e-01, 4.9918334928060942e-02);
2196 process_point_1(2.7146250701492608e-01, 6.2541213195902765e-02);
2197 process_point_3(2.1382490256170616e-02,
2198 1.2727971723358933e-01,
2199 1.5083677576511438e-02);
2200 process_point_3(2.3034156355267121e-02,
2201 2.9165567973834094e-01,
2202 2.1783585038607559e-02);
2203 process_point_3(1.1629601967792658e-01,
2204 2.5545422863851736e-01,
2205 4.3227363659414209e-02);
2214 b_point_permutations.push_back({centroid});
2215 b_weights.push_back(6.7960036586831640e-02);
2216 process_point_1(2.1509681108843159e-02, 6.0523371035391717e-03);
2217 process_point_1(4.8907694645253935e-01, 2.3994401928894731e-02);
2218 process_point_1(4.2694141425980042e-01, 5.5601967530453329e-02);
2219 process_point_1(2.2137228629183292e-01, 5.8278485119199981e-02);
2220 process_point_3(5.1263891023823893e-03,
2221 2.7251581777342970e-01,
2222 9.5906810035432631e-03);
2223 process_point_3(2.4370186901093827e-02,
2224 1.1092204280346341e-01,
2225 1.4965401105165668e-02);
2226 process_point_3(8.7895483032197297e-02,
2227 1.6359740106785048e-01,
2228 2.4179039811593819e-02);
2229 process_point_3(6.8012243554206653e-02,
2230 3.0844176089211778e-01,
2231 3.4641276140848373e-02);
2236 process_point_1(1.9390961248701044e-02, 4.9234036024000819e-03);
2237 process_point_1(6.1799883090872587e-02, 1.4433699669776668e-02);
2238 process_point_1(4.8896391036217862e-01, 2.1883581369428889e-02);
2239 process_point_1(4.1764471934045394e-01, 3.2788353544125348e-02);
2240 process_point_1(1.7720553241254344e-01, 4.2162588736993016e-02);
2241 process_point_1(2.7347752830883865e-01, 5.1774104507291585e-02);
2242 process_point_3(1.2683309328720416e-03,
2243 1.1897449769695684e-01,
2244 5.0102288385006719e-03);
2245 process_point_3(1.4646950055654417e-02,
2246 2.9837288213625779e-01,
2247 1.4436308113533840e-02);
2248 process_point_3(5.7124757403647919e-02,
2249 1.7226668782135557e-01,
2250 2.4665753212563674e-02);
2251 process_point_3(9.2916249356971847e-02,
2252 3.3686145979634496e-01,
2253 3.8571510787060684e-02);
2261 for (
unsigned int permutation_n = 0; permutation_n < b_weights.size();
2264 for (
const std::array<double, dim + 1> &b_point :
2265 b_point_permutations[permutation_n])
2267 const double volume = (dim == 2 ? 1.0 / 2.0 : 1.0 / 6.0);
2268 this->
weights.emplace_back(volume * b_weights[permutation_n]);
2270 for (
int d = 0;
d < dim; ++
d)
2271 c_point[
d] = b_point[
d];
2293 const unsigned int n_copies)
2303 const unsigned int n_copies)
2309 setup_qiterated_1D(base_quad, n_copies);
2314 const auto n_refinements =
2315 static_cast<unsigned int>(std::round(std::log2(n_copies)));
2316 Assert((1u << n_refinements) == n_copies,
2317 ExcMessage(
"The number of copies must be a power of 2."));
2330 std::vector<Point<dim>> points;
2331 std::vector<double> weights;
2335 for (
unsigned int qp = 0; qp < base_quad.
size(); ++qp)
2338 weights.push_back(fe_values.
JxW(qp));
2363 for (
unsigned int i = 0; i < quad_line.
size(); ++i)
2364 for (
unsigned int j = 0; j < quad_tri.
size(); ++j)
2367 quad_tri.
point(j)[1],
2368 quad_line.
point(i)[0]);
2385 if (n_points_1D == 1)
2387 const double Q14 = 1.0 / 4.0;
2388 const double Q43 = 4.0 / 3.0;
2391 this->
weights.emplace_back(Q43);
2393 else if (n_points_1D == 2)
2396 this->
quadrature_points.emplace_back(-0.26318405556971, -0.26318405556971, 0.54415184401122);
2397 this->
quadrature_points.emplace_back(-0.50661630334979, -0.50661630334979, 0.12251482265544);
2398 this->
quadrature_points.emplace_back(-0.26318405556971, +0.26318405556971, 0.54415184401122);
2399 this->
quadrature_points.emplace_back(-0.50661630334979, +0.50661630334979, 0.12251482265544);
2400 this->
quadrature_points.emplace_back(+0.26318405556971, -0.26318405556971, 0.54415184401122);
2401 this->
quadrature_points.emplace_back(+0.50661630334979, -0.50661630334979, 0.12251482265544);
2402 this->
quadrature_points.emplace_back(+0.26318405556971, +0.26318405556971, 0.54415184401122);
2403 this->
quadrature_points.emplace_back(+0.50661630334979, +0.50661630334979, 0.12251482265544);
2406 this->
weights.emplace_back(0.10078588207983);
2407 this->
weights.emplace_back(0.23254745125351);
2408 this->
weights.emplace_back(0.10078588207983);
2409 this->
weights.emplace_back(0.23254745125351);
2410 this->
weights.emplace_back(0.10078588207983);
2411 this->
weights.emplace_back(0.23254745125351);
2412 this->
weights.emplace_back(0.10078588207983);
2413 this->
weights.emplace_back(0.23254745125351);
2519 const std::vector<std::array<
Point<2>, 1 + 1>> &simplices)
const;
2523 const std::vector<std::array<
Point<3>, 1 + 1>> &simplices)
const;
2527 const std::vector<std::array<
Point<2>, 2 + 1>> &simplices)
const;
2531 const std::vector<std::array<
Point<3>, 2 + 1>> &simplices)
const;
2535 const std::vector<std::array<
Point<3>, 3 + 1>> &simplices)
const;
const Point< spacedim > & quadrature_point(const unsigned int q_point) const
double JxW(const unsigned int q_point) const
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, level_dof_access >> &cell)
Abstract base class for mapping classes.
QDuffy(const Quadrature< 1 > &radial_quadrature, const Quadrature< 1 > &angular_quadrature, const double beta=1.0)
QGaussChebyshev(const unsigned int n)
Generate a formula with n quadrature points.
QGaussLobattoChebyshev(const unsigned int n)
Generate a formula with n quadrature points.
QGaussLobatto(const unsigned int n)
QGaussLogR(const unsigned int n, const Point< dim > &x0=Point< dim >(), const double alpha=1, const bool factor_out_singular_weight=false)
QGaussLog(const unsigned int n, const bool revert=false)
static std::vector< double > get_quadrature_points(const unsigned int n)
static std::vector< double > get_quadrature_weights(const unsigned int n)
QGaussOneOverR(const unsigned int n, const Point< dim > &singularity, const bool factor_out_singular_weight=false)
static unsigned int quad_size(const Point< dim > &singularity, const unsigned int n)
QGaussPyramid(const unsigned int n_points_1D)
QGaussRadauChebyshev(const unsigned int n, const EndPoint end_point=QGaussRadauChebyshev::EndPoint::left)
Generate a formula with n quadrature points.
QGaussRadau(const unsigned int n, const EndPoint end_point=QGaussRadau::EndPoint::left)
QGaussSimplex(const unsigned int n_points_1D)
QGaussWedge(const unsigned int n_points_1D)
QGauss(const unsigned int n)
QIteratedSimplex(const Quadrature< dim > &base_quadrature, const unsigned int n_copies)
Quadrature< spacedim > compute_affine_transformation(const std::array< Point< spacedim >, dim+1 > &vertices) const
QSimplex(const Quadrature< dim > &quad)
Quadrature< spacedim > mapped_quadrature(const std::vector< std::array< Point< spacedim >, dim+1 >> &simplices) const
QSorted(const Quadrature< dim > &quad)
bool compare_weights(const unsigned int a, const unsigned int b) const
QSplit(const QSimplex< dim > &base, const Point< dim > &split_point)
QTelles(const Quadrature< 1 > &base_quad, const Point< dim > &singularity)
QTrianglePolar(const Quadrature< 1 > &radial_quadrature, const Quadrature< 1 > &angular_quadrature)
QWitherdenVincentSimplex(const unsigned int n_points_1D, const bool use_odd_order=true)
std::vector< Point< dim > > quadrature_points
Quadrature & operator=(const Quadrature< dim > &)
bool is_tensor_product_flag
double weight(const unsigned int i) const
const Point< dim > & point(const unsigned int i) const
std::vector< double > weights
unsigned int size() const
void refine_global(const unsigned int times=1)
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
IteratorRange< active_cell_iterator > active_cell_iterators() const
static ::ExceptionBase & ExcInternalError()
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotImplemented()
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
@ update_JxW_values
Transformed quadrature weights.
@ update_quadrature_points
Transformed quadrature points.
void simplex(Triangulation< dim, dim > &tria, const std::vector< Point< dim >> &vertices)
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
double norm(const FEValuesBase< dim > &fe, const ArrayView< const std::vector< Tensor< 1, dim >>> &Du)
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
void quadrature_points(const Triangulation< dim, spacedim > &triangulation, const Quadrature< dim > &quadrature, const std::vector< std::vector< BoundingBox< spacedim >>> &global_bounding_boxes, ParticleHandler< dim, spacedim > &particle_handler, const Mapping< dim, spacedim > &mapping=(ReferenceCells::get_hypercube< dim >() .template get_default_linear_mapping< dim, spacedim >()), const std::vector< std::vector< double >> &properties={})
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
Tensor< 2, dim, Number > w(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
Number jacobi_polynomial_value(const unsigned int degree, const int alpha, const int beta, const Number x)
std::string to_string(const number value, const unsigned int digits=numbers::invalid_unsigned_int)
std::vector< double > get_quadrature_weights(const unsigned int n)
std::vector< double > get_quadrature_points(const unsigned int n)
std::vector< double > get_quadrature_weights(const unsigned int n)
std::vector< double > get_quadrature_points(const unsigned int n)
std::vector< long double > compute_quadrature_weights(const std::vector< long double > &x, const int alpha, const int beta)
long double gamma(const unsigned int n)
std::vector< double > get_quadrature_points(const unsigned int n, const ::QGaussRadauChebyshev< 1 >::EndPoint end_point)
std::vector< double > get_quadrature_weights(const unsigned int n, const ::QGaussRadauChebyshev< 1 >::EndPoint end_point)
std::vector< double > get_quadrature_points(const unsigned int n, const ::QGaussRadau< 1 >::EndPoint end_point)
std::vector< double > get_quadrature_weights(const unsigned int n, const ::QGaussRadau< 1 >::EndPoint end_point)
std::vector< double > get_left_quadrature_weights(const unsigned int n)
std::vector< double > get_left_quadrature_points(const unsigned int n)
void split_point(const Point< dim1+dim2 > &source, Point< dim1 > &p1, Point< dim2 > &p2)
static constexpr double PI_2
static constexpr double PI
void transform(const InputIterator &begin_in, const InputIterator &end_in, OutputIterator out, const Function &function, const unsigned int grainsize)
static Point< dim > unit_cell_vertex(const unsigned int vertex)
const ::Triangulation< dim, spacedim > & tria