pub fn gauss_legendre_quadrature<F>(f: F, n: usize, (a, b): (f64, f64)) -> f64
where
F: Fn(f64) -> f64,
{
(b - a) / 2_f64
* unit_gauss_legendre_quadrature(
|x| f(x * (b - a) / 2_f64 + (a + b) / 2_f64),
n,
)
}
pub fn gauss_legendre_quadrature_cached<F>(
f: F,
(a, b): (f64, f64),
weights: &[f64],
roots: &[f64],
) -> f64
where
F: Fn(f64) -> f64,
{
(b - a) / 2_f64
* unit_gauss_legendre_quadrature_cached(
|x| f(x * (b - a) / 2_f64 + (a + b) / 2_f64),
weights,
roots,
)
}
fn unit_gauss_legendre_quadrature<F>(f: F, n: usize) -> f64
where
F: Fn(f64) -> f64,
{
let (weights, roots) = gauss_legendre_table(n);
unit_gauss_legendre_quadrature_cached(f, &weights, &roots)
}
fn unit_gauss_legendre_quadrature_cached<F>(
f: F,
weights: &[f64],
roots: &[f64],
) -> f64
where
F: Fn(f64) -> f64,
{
weights
.iter()
.zip(roots.iter())
.map(|(weight, &root)| weight * f(root))
.sum::<f64>()
}
pub fn gauss_legendre_table(n: usize) -> (Vec<f64>, Vec<f64>) {
let ref_root: &[f64] = match n {
2 => &LEGENDRE_ROOT_2[..],
3 => &LEGENDRE_ROOT_3[..],
4 => &LEGENDRE_ROOT_4[..],
5 => &LEGENDRE_ROOT_5[..],
6 => &LEGENDRE_ROOT_6[..],
7 => &LEGENDRE_ROOT_7[..],
8 => &LEGENDRE_ROOT_8[..],
9 => &LEGENDRE_ROOT_9[..],
10 => &LEGENDRE_ROOT_10[..],
11 => &LEGENDRE_ROOT_11[..],
12 => &LEGENDRE_ROOT_12[..],
13 => &LEGENDRE_ROOT_13[..],
14 => &LEGENDRE_ROOT_14[..],
15 => &LEGENDRE_ROOT_15[..],
16 => &LEGENDRE_ROOT_16[..],
17 => &LEGENDRE_ROOT_17[..],
18 => &LEGENDRE_ROOT_18[..],
19 => &LEGENDRE_ROOT_19[..],
20 => &LEGENDRE_ROOT_20[..],
21 => &LEGENDRE_ROOT_21[..],
22 => &LEGENDRE_ROOT_22[..],
23 => &LEGENDRE_ROOT_23[..],
24 => &LEGENDRE_ROOT_24[..],
25 => &LEGENDRE_ROOT_25[..],
26 => &LEGENDRE_ROOT_26[..],
27 => &LEGENDRE_ROOT_27[..],
28 => &LEGENDRE_ROOT_28[..],
29 => &LEGENDRE_ROOT_29[..],
30 => &LEGENDRE_ROOT_30[..],
_ => panic!("Legendre quadrature is limited up to n = 30"),
};
let ref_weight: &[f64] = match n {
2 => &LEGENDRE_WEIGHT_2[..],
3 => &LEGENDRE_WEIGHT_3[..],
4 => &LEGENDRE_WEIGHT_4[..],
5 => &LEGENDRE_WEIGHT_5[..],
6 => &LEGENDRE_WEIGHT_6[..],
7 => &LEGENDRE_WEIGHT_7[..],
8 => &LEGENDRE_WEIGHT_8[..],
9 => &LEGENDRE_WEIGHT_9[..],
10 => &LEGENDRE_WEIGHT_10[..],
11 => &LEGENDRE_WEIGHT_11[..],
12 => &LEGENDRE_WEIGHT_12[..],
13 => &LEGENDRE_WEIGHT_13[..],
14 => &LEGENDRE_WEIGHT_14[..],
15 => &LEGENDRE_WEIGHT_15[..],
16 => &LEGENDRE_WEIGHT_16[..],
17 => &LEGENDRE_WEIGHT_17[..],
18 => &LEGENDRE_WEIGHT_18[..],
19 => &LEGENDRE_WEIGHT_19[..],
20 => &LEGENDRE_WEIGHT_20[..],
21 => &LEGENDRE_WEIGHT_21[..],
22 => &LEGENDRE_WEIGHT_22[..],
23 => &LEGENDRE_WEIGHT_23[..],
24 => &LEGENDRE_WEIGHT_24[..],
25 => &LEGENDRE_WEIGHT_25[..],
26 => &LEGENDRE_WEIGHT_26[..],
27 => &LEGENDRE_WEIGHT_27[..],
28 => &LEGENDRE_WEIGHT_28[..],
29 => &LEGENDRE_WEIGHT_29[..],
30 => &LEGENDRE_WEIGHT_30[..],
_ => panic!("Legendre quadrature is limited up to n = 30"),
};
let mut result_root = vec![0_f64; n];
let mut result_weight = vec![0_f64; n];
result_root[..ref_root.len()].clone_from_slice(ref_root);
result_weight[..ref_weight.len()].clone_from_slice(ref_weight);
match n % 2 {
0 => {
for i in ref_root.len()..n {
result_root[i] = -ref_root[n - i - 1];
result_weight[i] = ref_weight[n - i - 1];
}
}
1 => {
for i in ref_root.len()..n {
result_root[i] = -ref_root[n - i];
result_weight[i] = ref_weight[n - i];
}
}
_ => unreachable!(),
}
(result_weight, result_root)
}
const LEGENDRE_ROOT_2: [f64; 1] = [0.577_350_269_189_626];
const LEGENDRE_ROOT_3: [f64; 2] = [0_f64, 0.774_596_669_241_483];
const LEGENDRE_ROOT_4: [f64; 2] =
[0.339_981_043_584_856, 0.861_136_311_594_053];
const LEGENDRE_ROOT_5: [f64; 3] =
[0_f64, 0.538_469_310_105_683, 0.906_179_845_938_664];
const LEGENDRE_ROOT_6: [f64; 3] = [
0.238_619_186_083_197,
0.661_209_386_466_265,
0.932_469_514_203_152,
];
const LEGENDRE_ROOT_7: [f64; 4] = [
0_f64,
0.405_845_151_377_397,
0.741_531_185_599_394,
0.949_107_912_342_759,
];
const LEGENDRE_ROOT_8: [f64; 4] = [
0.183_434_642_495_650,
0.525_532_409_916_329,
0.796_666_477_413_627,
0.960_289_856_497_536,
];
const LEGENDRE_ROOT_9: [f64; 5] = [
0_f64,
0.324_253_423_403_809,
0.613_371_432_700_590,
0.836_031_107_326_636,
0.968_160_239_507_626,
];
const LEGENDRE_ROOT_10: [f64; 5] = [
0.148_874_338_981_631,
0.433_395_394_129_247,
0.679_409_568_299_024,
0.865_063_366_688_985,
0.973_906_528_517_172,
];
const LEGENDRE_ROOT_11: [f64; 6] = [
0_f64,
0.269_543_155_952_345,
0.519_096_129_110_681,
0.730_152_005_574_049,
0.887_062_599_768_095,
0.978_228_658_146_057,
];
const LEGENDRE_ROOT_12: [f64; 6] = [
0.125_333_408_511_469,
0.367_831_498_918_180,
0.587_317_954_286_617,
0.769_902_674_194_305,
0.904_117_256_370_475,
0.981_560_634_246_719,
];
const LEGENDRE_ROOT_13: [f64; 7] = [
0_f64,
0.230_458_315_955_135,
0.448_492_751_036_447,
0.642_349_339_440_340,
0.801_578_090_733_310,
0.917_598_399_222_978,
0.984_183_054_718_588,
];
const LEGENDRE_ROOT_14: [f64; 7] = [
0.108_054_948_707_344,
0.319_112_368_927_890,
0.515_248_636_358_154,
0.687_292_904_811_685,
0.827_201_315_069_765,
0.928_434_883_663_574,
0.986_283_808_696_812,
];
const LEGENDRE_ROOT_15: [f64; 8] = [
0_f64,
0.201_194_093_997_435,
0.394_151_347_077_563,
0.570_972_172_608_539,
0.724_417_731_360_170,
0.848_206_583_410_427,
0.937_273_392_400_706,
0.987_992_518_020_485,
];
const LEGENDRE_ROOT_16: [f64; 8] = [
0.095_012_509_837_637,
0.281_603_550_779_259,
0.458_016_777_657_227,
0.617_876_244_402_644,
0.755_404_408_355_003,
0.865_631_202_387_832,
0.944_575_023_073_233,
0.989_400_934_991_650,
];
const LEGENDRE_ROOT_17: [f64; 9] = [
0_f64,
0.178_484_181_495_847_85,
0.351_231_763_453_876_3,
0.512_690_537_086_476_9,
0.657_671_159_216_690_7,
0.781_514_003_896_801_4,
0.880_239_153_726_985_9,
0.950_675_521_768_767_8,
0.990_575_475_314_417_4,
];
const LEGENDRE_ROOT_18: [f64; 9] = [
0.084_775_013_041_735_3,
0.251_886_225_691_505_5,
0.411_751_161_462_842_63,
0.559_770_831_073_947_5,
0.691_687_043_060_353_2,
0.803_704_958_972_523_1,
0.892_602_466_497_555_7,
0.955_823_949_571_397_7,
0.991_565_168_420_930_9,
];
const LEGENDRE_ROOT_19: [f64; 10] = [
0_f64,
0.160_358_645_640_225_37,
0.316_564_099_963_629_83,
0.464_570_741_375_960_94,
0.600_545_304_661_681,
0.720_966_177_335_229_4,
0.822_714_656_537_142_8,
0.903_155_903_614_817_9,
0.960_208_152_134_83,
0.992_406_843_843_584_4,
];
const LEGENDRE_ROOT_20: [f64; 10] = [
0.076_526_521_133_497_34,
0.227_785_851_141_645_07,
0.373_706_088_715_419_55,
0.510_867_001_950_827_1,
0.636_053_680_726_515,
0.746_331_906_460_150_8,
0.839_116_971_822_218_8,
0.912_234_428_251_326,
0.963_971_927_277_913_8,
0.993_128_599_185_094_9,
];
const LEGENDRE_ROOT_21: [f64; 11] = [
0_f64,
0.145_561_854_160_895_1,
0.288_021_316_802_401_1,
0.424_342_120_207_438_8,
0.551_618_835_887_219_8,
0.667_138_804_197_412_3,
0.768_439_963_475_677_9,
0.853_363_364_583_317_3,
0.920_099_334_150_400_8,
0.967_226_838_566_306_3,
0.993_752_170_620_389_5,
];
const LEGENDRE_ROOT_22: [f64; 11] = [
0.069_739_273_319_722_23,
0.207_860_426_688_221_27,
0.341_935_820_892_084_24,
0.469_355_837_986_757,
0.587_640_403_506_911_6,
0.694_487_263_186_682_7,
0.787_816_805_979_208_1,
0.865_812_577_720_300_2,
0.926_956_772_187_174,
0.970_060_497_835_428_7,
0.994_294_585_482_399_2,
];
const LEGENDRE_ROOT_23: [f64; 12] = [
0_f64,
0.133_256_824_298_466_1,
0.264_135_680_970_344_93,
0.390_301_038_030_290_8,
0.509_501_477_846_007_5,
0.619_609_875_763_646_1,
0.718_661_363_131_950_2,
0.804_888_401_618_839_9,
0.876_752_358_270_441_6,
0.932_971_086_826_016_1,
0.972_542_471_218_115_2,
0.994_769_334_997_552_2,
];
const LEGENDRE_ROOT_24: [f64; 12] = [
0.064_056_892_862_605_63,
0.191_118_867_473_616_3,
0.315_042_679_696_163_4,
0.433_793_507_626_045_1,
0.545_421_471_388_839_6,
0.648_093_651_936_975_5,
0.740_124_191_578_554_4,
0.820_001_985_973_903,
0.886_415_527_004_401_1,
0.938_274_552_002_732_8,
0.974_728_555_971_309_5,
0.995_187_219_997_021_3,
];
const LEGENDRE_ROOT_25: [f64; 13] = [
0_f64,
0.122_864_692_610_710_4,
0.243_866_883_720_988_44,
0.361_172_305_809_387_86,
0.473_002_731_445_715,
0.577_662_930_241_223,
0.673_566_368_473_468_4,
0.759_259_263_037_357_6,
0.833_442_628_760_834,
0.894_991_997_878_275_3,
0.942_974_571_228_974_3,
0.976_663_921_459_517_5,
0.995_556_969_790_498_1,
];
const LEGENDRE_ROOT_26: [f64; 13] = [
0.059_230_093_429_313_21,
0.176_858_820_356_890_18,
0.292_004_839_485_956_9,
0.403_051_755_123_486_3,
0.508_440_714_824_505_7,
0.606_692_293_017_618_1,
0.696_427_260_419_957_3,
0.776_385_948_820_678_9,
0.845_445_942_788_498,
0.902_637_861_984_307_1,
0.947_159_066_661_714_2,
0.978_385_445_956_471,
0.995_885_701_145_616_9,
];
const LEGENDRE_ROOT_27: [f64; 14] = [
0_f64,
0.113_972_585_609_529_97,
0.226_459_365_439_536_85,
0.335_993_903_638_508_9,
0.441_148_251_750_026_87,
0.540_551_564_579_456_9,
0.632_907_971_946_495_2,
0.717_013_473_739_423_7,
0.791_771_639_070_508_2,
0.856_207_908_018_294_5,
0.909_482_320_677_491_1,
0.950_900_557_814_705_1,
0.979_923_475_961_501_2,
0.996_179_262_888_988_6,
];
const LEGENDRE_ROOT_28: [f64; 14] = [
0.055_079_289_884_034_27,
0.164_569_282_133_380_76,
0.272_061_627_635_178_1,
0.376_251_516_089_078_7,
0.475_874_224_955_118_3,
0.569_720_471_811_401_7,
0.656_651_094_038_865,
0.735_610_878_013_631_8,
0.805_641_370_917_179_1,
0.865_892_522_574_395_1,
0.915_633_026_392_132_1,
0.954_259_280_628_938_2,
0.981_303_165_370_872_7,
0.996_442_497_573_954_4,
];
const LEGENDRE_ROOT_29: [f64; 15] = [
0_f64,
0.106_278_230_132_679_23,
0.211_352_286_166_001_07,
0.314_031_637_867_639_93,
0.413_152_888_174_008_64,
0.507_592_955_124_227_6,
0.596_281_797_138_227_8,
0.678_214_537_602_686_5,
0.752_462_851_734_477_1,
0.818_185_487_615_252_4,
0.874_637_804_920_102_8,
0.921_180_232_953_058_7,
0.957_285_595_778_087_7,
0.982_545_505_261_413_2,
0.996_679_442_260_596_6,
];
const LEGENDRE_ROOT_30: [f64; 15] = [
0.051_471_842_555_317_7,
0.153_869_913_608_583_54,
0.254_636_926_167_889_85,
0.352_704_725_530_878_1,
0.447_033_769_538_089_15,
0.536_624_148_142_019_9,
0.620_526_182_989_242_9,
0.697_850_494_793_315_8,
0.767_777_432_104_826_2,
0.829_565_762_382_768_4,
0.882_560_535_792_052_7,
0.926_200_047_429_274_3,
0.960_021_864_968_307_5,
0.983_668_123_279_747_2,
0.996_893_484_074_649_5,
];
const LEGENDRE_WEIGHT_2: [f64; 1] = [1_f64];
const LEGENDRE_WEIGHT_3: [f64; 2] =
[0.888_888_888_888_889, 0.555_555_555_555_556];
const LEGENDRE_WEIGHT_4: [f64; 2] =
[0.652_145_154_862_546, 0.347_854_845_137_454];
const LEGENDRE_WEIGHT_5: [f64; 3] = [
0.568_888_888_888_889,
0.478_628_670_499_366,
0.236_926_885_056_189,
];
const LEGENDRE_WEIGHT_6: [f64; 3] = [
0.467_913_934_572_691,
0.360_761_573_048_139,
0.171_324_492_379_170,
];
const LEGENDRE_WEIGHT_7: [f64; 4] = [
0.417_959_183_673_469,
0.381_830_050_505_119,
0.279_705_391_489_277,
0.129_484_966_168_870,
];
const LEGENDRE_WEIGHT_8: [f64; 4] = [
0.362_683_783_378_362,
0.313_706_645_877_887,
0.222_381_034_453_374,
0.101_228_536_290_376,
];
const LEGENDRE_WEIGHT_9: [f64; 5] = [
0.330_239_355_001_260,
0.312_347_077_040_003,
0.260_610_696_402_935,
0.180_648_160_694_857,
0.081_274_388_361_574,
];
const LEGENDRE_WEIGHT_10: [f64; 5] = [
0.295_524_224_714_753,
0.269_266_719_309_996,
0.219_086_362_515_982,
0.149_451_349_150_581,
0.066_671_344_308_688,
];
const LEGENDRE_WEIGHT_11: [f64; 6] = [
0.272_925_086_777_901,
0.262_804_544_510_247,
0.233_193_764_591_990,
0.186_290_210_927_734,
0.125_580_369_464_905,
0.055_668_567_116_174,
];
const LEGENDRE_WEIGHT_12: [f64; 6] = [
0.249_147_045_813_403,
0.233_492_536_538_355,
0.203_167_426_723_066,
0.160_078_328_543_346,
0.106_939_325_995_318,
0.047_175_336_386_512,
];
const LEGENDRE_WEIGHT_13: [f64; 7] = [
0.232_551_553_230_874,
0.226_283_180_262_897,
0.207_816_047_536_889,
0.178_145_980_761_946,
0.138_873_510_219_787,
0.092_121_499_837_728,
0.040_484_004_765_316,
];
const LEGENDRE_WEIGHT_14: [f64; 7] = [
0.215_263_853_463_158,
0.205_198_463_721_290,
0.185_538_397_477_938,
0.157_203_167_158_194,
0.121_518_570_687_903,
0.080_158_087_159_760,
0.035_119_460_331_752,
];
const LEGENDRE_WEIGHT_15: [f64; 8] = [
0.202_578_241_925_561,
0.198_431_485_327_111,
0.186_161_000_015_562,
0.166_269_205_816_994,
0.139_570_677_926_154,
0.107_159_220_467_172,
0.070_366_047_488_108,
0.030_753_241_996_117,
];
const LEGENDRE_WEIGHT_16: [f64; 8] = [
0.189_450_610_455_069,
0.182_603_415_044_924,
0.169_156_519_395_003,
0.149_595_988_816_577,
0.124_628_971_255_534,
0.095_158_511_682_493,
0.062_253_523_938_648,
0.027_152_459_411_754,
];
const LEGENDRE_WEIGHT_17: [f64; 9] = [
0.179_446_470_356_206_53,
0.176_562_705_366_992_65,
0.168_004_102_156_450_04,
0.154_045_761_076_810_28,
0.135_136_368_468_525_47,
0.111_883_847_193_403_97,
0.085_036_148_317_179_18,
0.055_459_529_373_987_2,
0.024_148_302_868_547_93,
];
const LEGENDRE_WEIGHT_18: [f64; 9] = [
0.169_142_382_963_143_6,
0.164_276_483_745_832_72,
0.154_684_675_126_265_24,
0.140_642_914_670_650_65,
0.122_555_206_711_478_46,
0.100_942_044_106_287_17,
0.076_425_730_254_889_06,
0.049_714_548_894_969_8,
0.021_616_013_526_483_31,
];
const LEGENDRE_WEIGHT_19: [f64; 10] = [
0.161_054_449_848_783_7,
0.158_968_843_393_954_34,
0.152_766_042_065_859_67,
0.142_606_702_173_606_6,
0.128_753_962_539_336_21,
0.111_566_645_547_333_99,
0.091_490_021_622_45,
0.069_044_542_737_641_23,
0.044_814_226_765_699_6,
0.019_461_788_229_726_477,
];
const LEGENDRE_WEIGHT_20: [f64; 10] = [
0.152_753_387_130_725_85,
0.149_172_986_472_603_74,
0.142_096_109_318_382_04,
0.131_688_638_449_176_63,
0.118_194_531_961_518_41,
0.101_930_119_817_240_44,
0.083_276_741_576_704_75,
0.062_672_048_334_109_07,
0.040_601_429_800_386_94,
0.017_614_007_139_152_118,
];
const LEGENDRE_WEIGHT_21: [f64; 11] = [
0.146_081_133_649_690_41,
0.144_524_403_989_970_05,
0.139_887_394_791_073_15,
0.132_268_938_633_337_46,
0.121_831_416_053_728_53,
0.108_797_299_167_148_38,
0.093_444_423_456_033_86,
0.076_100_113_628_379_3,
0.057_134_425_426_857_21,
0.036_953_789_770_852_494,
0.016_017_228_257_774_335,
];
const LEGENDRE_WEIGHT_22: [f64; 11] = [
0.139_251_872_855_631_98,
0.136_541_498_346_015_17,
0.131_173_504_787_062_38,
0.123_252_376_810_512_42,
0.112_932_296_080_539_22,
0.100_414_144_442_880_96,
0.085_941_606_217_067_73,
0.069_796_468_424_520_49,
0.052_293_335_152_683_286,
0.033_774_901_584_814_16,
0.014_627_995_298_272_201,
];
const LEGENDRE_WEIGHT_23: [f64; 12] = [
0.133_654_572_186_106_18,
0.132_462_039_404_696_6,
0.128_905_722_188_082_15,
0.123_049_084_306_729_53,
0.114_996_640_222_411_36,
0.104_892_091_464_541_41,
0.092_915_766_060_035_15,
0.079_281_411_776_718_95,
0.064_232_421_408_525_85,
0.048_037_671_731_084_67,
0.030_988_005_856_979_444,
0.013_411_859_487_141_772,
];
const LEGENDRE_WEIGHT_24: [f64; 12] = [
0.127_938_195_346_752_16,
0.125_837_456_346_828_3,
0.121_670_472_927_803_39,
0.115_505_668_053_725_6,
0.107_444_270_115_965_63,
0.097_618_652_104_113_89,
0.086_190_161_531_953_28,
0.073_346_481_411_080_31,
0.059_298_584_915_436_78,
0.044_277_438_817_419_81,
0.028_531_388_628_933_663,
0.012_341_229_799_987_2,
];
const LEGENDRE_WEIGHT_25: [f64; 13] = [
0.123_176_053_726_715_45,
0.122_242_442_990_310_04,
0.119_455_763_535_784_77,
0.114_858_259_145_711_64,
0.108_519_624_474_263_65,
0.100_535_949_067_050_64,
0.091_028_261_982_963_65,
0.080_140_700_335_001_02,
0.068_038_333_812_356_91,
0.054_904_695_975_835_192,
0.040_939_156_701_306_313,
0.026_354_986_615_032_137,
0.011_393_798_501_026_288,
];
const LEGENDRE_WEIGHT_26: [f64; 13] = [
0.118_321_415_279_262_28,
0.116_660_443_485_296_58,
0.113_361_816_546_319_67,
0.108_471_840_528_576_59,
0.102_059_161_094_425_42,
0.094_213_800_355_914_15,
0.085_045_894_313_485_23,
0.074_684_149_765_659_75,
0.063_274_046_329_574_84,
0.050_975_825_297_147_81,
0.037_962_383_294_362_76,
0.024_417_851_092_631_91,
0.010_551_372_617_343_01,
];
const LEGENDRE_WEIGHT_27: [f64; 14] = [
0.114_220_867_378_957,
0.113_476_346_108_965_15,
0.111_252_488_356_845_19,
0.107_578_285_788_533_19,
0.102_501_637_817_745_8,
0.096_088_727_370_028_5,
0.088_423_158_543_756_95,
0.079_604_867_773_057_77,
0.069_748_823_766_245_6,
0.058_983_536_859_833_6,
0.047_449_412_520_615_06,
0.035_297_053_757_419_71,
0.022_686_231_596_180_62,
0.009_798_996_051_294_36,
];
const LEGENDRE_WEIGHT_28: [f64; 14] = [
0.110_047_013_016_475_2,
0.108_711_192_258_294_13,
0.106_055_765_922_846_42,
0.102_112_967_578_060_77,
0.096_930_657_997_929_92,
0.090_571_744_393_032_84,
0.083_113_417_228_901_21,
0.074_646_214_234_568_78,
0.065_272_923_966_999_6,
0.055_107_345_675_716_75,
0.044_272_934_759_004_23,
0.032_901_427_782_304_38,
0.021_132_112_592_771_26,
0.009_124_282_593_094_52,
];
const LEGENDRE_WEIGHT_29: [f64; 15] = [
0.106_479_381_718_314_24,
0.105_876_155_097_320_94,
0.104_073_310_077_729_38,
0.101_091_273_759_914_97,
0.096_963_834_094_408_6,
0.091_737_757_139_258_76,
0.085_472_257_366_172_53,
0.078_238_327_135_763_79,
0.070_117_933_255_051_28,
0.061_203_090_657_079_139,
0.051_594_826_902_497_92,
0.041_402_062_518_682_836,
0.030_740_492_202_093_623,
0.019_732_085_056_122_706,
0.008_516_903_878_746_41,
];
const LEGENDRE_WEIGHT_30: [f64; 15] = [
0.102_852_652_893_558_84,
0.101_762_389_748_405_5,
0.099_593_420_586_795_27,
0.096_368_737_174_644_26,
0.092_122_522_237_786_13,
0.086_899_787_201_082_98,
0.080_755_895_229_420_22,
0.073_755_974_737_705_2,
0.065_974_229_882_180_5,
0.057_493_156_217_619_066,
0.048_402_672_830_594_05,
0.038_799_192_569_627_05,
0.028_784_707_883_323_37,
0.018_466_468_311_090_96,
0.007_968_192_496_166_606,
];
#[cfg(test)]
mod tests {
use super::*;
#[test]
fn gauss_legendre_quadrature_3rd_order_poly() {
let range = (-1.0, 2.0);
let f = |x: f64| {
2_f64.mul_add(x.powi(3), x * x) - 3_f64.mul_add(x, -4.0)
};
let ans = 18.0;
for n in 2..=30 {
let q = gauss_legendre_quadrature(f, n, range);
match n {
11 => assert::close(q, ans, 1e-9),
12 => assert::close(q, ans, 1e-3),
_ => assert::close(q, ans, 1e-13),
}
}
}
#[test]
fn gauss_legendre_quadrature_cached_3rd_order_poly() {
let range = (-1.0, 2.0);
let f = |x: f64| {
2_f64.mul_add(x.powi(3), x * x) - 3_f64.mul_add(x, -4.0)
};
let ans = 18.0;
for n in 2..=30 {
let (weights, roots) = gauss_legendre_table(n);
let q =
gauss_legendre_quadrature_cached(f, range, &weights, &roots);
match n {
11 => assert::close(q, ans, 1e-9),
12 => assert::close(q, ans, 1e-3),
_ => assert::close(q, ans, 1e-13),
}
}
}
}