1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155

//! # gauss-quad
//!
//! **gauss-quad** is a Gaussian quadrature library for numerical integration.
//!
//! ## Quadrature rules
//! **gauss-quad** implements the following quadrature rules:
//! * Gauss-Legendre
//! * Gauss-Jacobi
//! * Gauss-Laguerre
//! * Gauss-Hermite
//! * Midpoint
//! * Simpson
//!
//! ## Using **gauss-quad**
//! To use any of the quadrature rules in your project, first initialize the rule with
//! a specified degree and then you can use it for integration, e.g.:
//! ```
//! use gauss_quad::GaussLegendre;
//! // This macro is used in these docs to compare floats.
//! // The assertion succeeds if the two sides are within floating point error,
//! // or an optional epsilon.
//! use approx::assert_abs_diff_eq;
//!
//! // initialize the quadrature rule
//! let degree = 10;
//! let quad = GaussLegendre::init(degree);
//!
//! // use the rule to integrate a function
//! let left_bound = 0.0;
//! let right_bound = 1.0;
//! let integral = quad.integrate(left_bound, right_bound, |x| x * x);
//! assert_abs_diff_eq!(integral, 1.0 / 3.0);
//! ```
//!
//! ## Setting up a quadrature rule
//! Using a quadrature rule takes two steps:
//! 1. Initialization
//! 2. Integration
//!
//! First, rules must be initialized using some specific input parameters.
//!
//! Then, you can integrate functions using those rules:
//! ```
//! # use gauss_quad::*;
//! # let degree = 5;
//! # let alpha = 1.2;
//! # let beta = 1.2;
//! # let a = 0.0;
//! # let b = 1.0;
//! # let c = -10.;
//! # let d = 100.;
//! let gauss_legendre = GaussLegendre::init(degree);
//! // Integrate on the domain [a, b]
//! let x_cubed = gauss_legendre.integrate(a, b, |x| x * x * x);
//!
//! let gauss_jacobi = GaussJacobi::init(degree, alpha, beta);
//! // Integrate on the domain [c, d]
//! let double_x = gauss_jacobi.integrate(c, d, |x| 2.0 * x);
//!
//! let gauss_laguerre = GaussLaguerre::init(degree, alpha);
//! // no explicit domain, Gauss-Laguerre integration is done on the domain [0, ∞).
//! let piecewise = gauss_laguerre.integrate(|x| if x > 0.0 && x < 2.0 { x } else { 0.0 });
//!
//! let gauss_hermite = GaussHermite::init(degree);
//! // again, no explicit domain since integration is done over the domain (-∞, ∞).
//! let constant = gauss_hermite.integrate(|x| if x > -1.0 && x < 1.0 { 2.0 } else { 1.0 });
//! ```
//!
//! ## Specific quadrature rules
//! Different rules may take different parameters.
//!
//! For example, the `GaussLaguerre` rule requires both a `degree` and an `alpha`
//! parameter.
//!
//! `GaussLaguerre` is also defined as an improper integral over the domain [0, ∞).
//! This means no domain bounds are needed in the `integrate` call.
//! ```
//! # use gauss_quad::GaussLaguerre;
//! // initialize the quadrature rule
//! let degree = 10;
//! let alpha = 0.5;
//! let quad = GaussLaguerre::init(degree, alpha);
//!
//! // use the rule to integrate a function
//! let integral = quad.integrate(|x| x * x);
//! ```
//!
//! ## Panics and errors
//! Quadrature rules are only defined for a certain set of input values.
//! For example, every rule is only defined for degrees where `degree > 1`.
//! ```should_panic
//! # use gauss_quad::GaussLaguerre;
//! let degree = 1;
//! let quad = GaussLaguerre::init(degree, 0.1); // panics!
//! ```
//!
//! Specific rules may have other requirements.
//! `GaussJacobi` for example, requires alpha and beta parameters larger than -1.0.
//! ```should_panic
//! # use gauss_quad::GaussJacobi;
//! let degree = 10;
//! let alpha = 0.1;
//! let beta = -1.1;
//!
//! let quad = GaussJacobi::init(degree, alpha, beta); // panics!
//! ```
//! Make sure to read the specific quadrature rule's documentation before using it.
//!
//! Error handling is very simple: bad input values will cause the program to panic
//! and abort with a short error message.
//!
//! ## Passing functions to quadrature rules
//! The `integrate` method expects functions of the form `Fn(f64) -> f64`, i.e. functions of
//! one parameter.
//!
//! ```
//! # use gauss_quad::GaussLegendre;
//!
//! // initialize the quadrature rule
//! let degree = 10;
//! let quad = GaussLegendre::init(degree);
//!
//! // use the rule to integrate a function
//! let left_bound = 0.0;
//! let right_bound = 1.0;
//! let integral = quad.integrate(left_bound, right_bound, |x| x * x);
//! ```

use nalgebra::{Dynamic, Matrix, VecStorage};
pub type DMatrixf64 = Matrix<f64, Dynamic, Dynamic, VecStorage<f64, Dynamic, Dynamic>>;
#[doc(inline)]
pub use core::f64::consts::PI;

mod gamma;
#[cfg(test)]
mod gaussian_quadrature;
pub mod hermite;
pub mod jacobi;
pub mod laguerre;
pub mod legendre;
pub mod midpoint;
pub mod simpson;

#[doc(inline)]
pub use hermite::GaussHermite;
#[doc(inline)]
pub use jacobi::GaussJacobi;
#[doc(inline)]
pub use laguerre::GaussLaguerre;
#[doc(inline)]
pub use legendre::GaussLegendre;
#[doc(inline)]
pub use midpoint::Midpoint;
#[doc(inline)]
pub use simpson::Simpson;