Nektar++
Polylib.h
Go to the documentation of this file.
1 #ifndef H_PLYLIB
2 
5 #include <complex>
6 
7 /*
8  * LIBRARY ROUTINES FOR POLYNOMIAL CALCULUS AND INTERPOLATION
9  */
10 
11 /** \brief The namespace associated with the the Polylib library
12  * (\ref pagePolylib "Polylib introduction")
13  */
14 namespace Polylib
15 {
16 
17 /**
18  \page pagePolylib The Polylib library
19  \section sectionPolyLib Routines For Orthogonal Polynomial Calculus and
20  Interpolation
21 
22  Spencer Sherwin,
23  Aeronautics, Imperial College London
24 
25  Based on codes by Einar Ronquist and Ron Henderson
26 
27  Abbreviations
28  - z - Set of collocation/quadrature points
29  - w - Set of quadrature weights
30  - D - Derivative matrix
31  - h - Lagrange Interpolant
32  - I - Interpolation matrix
33  - g - Gauss
34  - k - Kronrod
35  - gr - Gauss-Radau
36  - gl - Gauss-Lobatto
37  - j - Jacobi
38  - m - point at minus 1 in Radau rules
39  - p - point at plus 1 in Radau rules
40 
41  -----------------------------------------------------------------------\n
42  MAIN ROUTINES\n
43  -----------------------------------------------------------------------\n
44 
45  Points and Weights:
46 
47  - zwgj Compute Gauss-Jacobi points and weights
48  - zwgrjm Compute Gauss-Radau-Jacobi points and weights (z=-1)
49  - zwgrjp Compute Gauss-Radau-Jacobi points and weights (z= 1)
50  - zwglj Compute Gauss-Lobatto-Jacobi points and weights
51  - zwgk Compute Gauss-Kronrod-Jacobi points and weights
52  - zwrk Compute Radau-Kronrod points and weights
53  - zwlk Compute Lobatto-Kronrod points and weights
54  - JacZeros Compute Gauss-Jacobi points and weights
55 
56  Derivative Matrices:
57 
58  - Dgj Compute Gauss-Jacobi derivative matrix
59  - Dgrjm Compute Gauss-Radau-Jacobi derivative matrix (z=-1)
60  - Dgrjp Compute Gauss-Radau-Jacobi derivative matrix (z= 1)
61  - Dglj Compute Gauss-Lobatto-Jacobi derivative matrix
62 
63  Lagrange Interpolants:
64 
65  - hgj Compute Gauss-Jacobi Lagrange interpolants
66  - hgrjm Compute Gauss-Radau-Jacobi Lagrange interpolants (z=-1)
67  - hgrjp Compute Gauss-Radau-Jacobi Lagrange interpolants (z= 1)
68  - hglj Compute Gauss-Lobatto-Jacobi Lagrange interpolants
69 
70  Interpolation Operators:
71 
72  - Imgj Compute interpolation operator gj->m
73  - Imgrjm Compute interpolation operator grj->m (z=-1)
74  - Imgrjp Compute interpolation operator grj->m (z= 1)
75  - Imglj Compute interpolation operator glj->m
76 
77  Polynomial Evaluation:
78 
79  - jacobfd Returns value and derivative of Jacobi poly. at point z
80  - jacobd Returns derivative of Jacobi poly. at point z (valid at z=-1,1)
81 
82  -----------------------------------------------------------------------\n
83  LOCAL ROUTINES\n
84  -----------------------------------------------------------------------\n
85 
86  - jacobz Returns Jacobi polynomial zeros
87  - gammaf Gamma function for integer values and halves
88  - RecCoeff Calculates the recurrence coefficients for orthogonal poly
89  - TriQL QL algorithm for symmetrix tridiagonal matrix
90  - JKMatrix Generates the Jacobi-Kronrod matrix
91 
92  ------------------------------------------------------------------------\n
93 
94  Useful references:
95 
96  - [1] Gabor Szego: Orthogonal Polynomials, American Mathematical Society,
97  Providence, Rhode Island, 1939.
98  - [2] Abramowitz \& Stegun: Handbook of Mathematical Functions,
99  Dover, New York, 1972.
100  - [3] Canuto, Hussaini, Quarteroni \& Zang: Spectral Methods in Fluid
101  Dynamics, Springer-Verlag, 1988.
102  - [4] Ghizzetti \& Ossicini: Quadrature Formulae, Academic Press, 1970.
103  - [5] Karniadakis \& Sherwin: Spectral/hp element methods for CFD, 1999
104 
105 
106  NOTES
107 
108  -# Legendre polynomial \f$ \alpha = \beta = 0 \f$
109  -# Chebychev polynomial \f$ \alpha = \beta = -0.5 \f$
110  -# All routines are double precision.
111  -# All array subscripts start from zero, i.e. vector[0..N-1]
112 */
113 
114 /*-----------------------------------------------------------------------
115  M A I N R O U T I N E S
116  -----------------------------------------------------------------------*/
117 
118 /* Points and weights */
119 LIB_UTILITIES_EXPORT void zwgj(double *, double *, const int, const double,
120  const double);
121 LIB_UTILITIES_EXPORT void zwgrjm(double *, double *, const int, const double,
122  const double);
123 LIB_UTILITIES_EXPORT void zwgrjp(double *, double *, const int, const double,
124  const double);
125 LIB_UTILITIES_EXPORT void zwglj(double *, double *, const int, const double,
126  const double);
127 LIB_UTILITIES_EXPORT void zwgk(double *, double *, const int, const double,
128  const double);
129 LIB_UTILITIES_EXPORT void zwrk(double *, double *, const int, const double,
130  const double);
131 LIB_UTILITIES_EXPORT void zwlk(double *, double *, const int, const double,
132  const double);
133 LIB_UTILITIES_EXPORT void JacZeros(const int, double *, double *, const double,
134  const double);
135 
136 /* Derivative operators */
137 LIB_UTILITIES_EXPORT void Dgj(double *, const double *, const int, const double,
138  const double);
139 LIB_UTILITIES_EXPORT void Dgrjm(double *, const double *, const int,
140  const double, const double);
141 LIB_UTILITIES_EXPORT void Dgrjp(double *, const double *, const int,
142  const double, const double);
143 LIB_UTILITIES_EXPORT void Dglj(double *, const double *, const int,
144  const double, const double);
145 
146 /* Lagrangian interpolants */
147 LIB_UTILITIES_EXPORT double hgj(const int, const double, const double *,
148  const int, const double, const double);
149 LIB_UTILITIES_EXPORT double hgrjm(const int, const double, const double *,
150  const int, const double, const double);
151 LIB_UTILITIES_EXPORT double hgrjp(const int, const double, const double *,
152  const int, const double, const double);
153 LIB_UTILITIES_EXPORT double hglj(const int, const double, const double *,
154  const int, const double, const double);
155 
156 /* Interpolation operators */
157 LIB_UTILITIES_EXPORT void Imgj(double *, const double *, const double *,
158  const int, const int, const double,
159  const double);
160 LIB_UTILITIES_EXPORT void Imgrjm(double *, const double *, const double *,
161  const int, const int, const double,
162  const double);
163 LIB_UTILITIES_EXPORT void Imgrjp(double *, const double *, const double *,
164  const int, const int, const double,
165  const double);
166 LIB_UTILITIES_EXPORT void Imglj(double *, const double *, const double *,
167  const int, const int, const double,
168  const double);
169 
170 /* Polynomial functions */
171 LIB_UTILITIES_EXPORT void jacobfd(const int, const double *, double *, double *,
172  const int, const double, const double);
173 LIB_UTILITIES_EXPORT void jacobd(const int, const double *, double *, const int,
174  const double, const double);
175 
176 LIB_UTILITIES_EXPORT std::complex<Nektar::NekDouble> ImagBesselComp(
177  const int, std::complex<Nektar::NekDouble>);
178 
179 /* Gamma function routines */
180 LIB_UTILITIES_EXPORT double gammaF(const double);
181 LIB_UTILITIES_EXPORT double gammaFracGammaF(const int, const double, const int,
182  const double);
183 
184 } // namespace Polylib
185 
186 #define H_PLYLIB
187 #endif /* END OF POLYLIB.H DECLARATIONS */
#define LIB_UTILITIES_EXPORT
The namespace associated with the the Polylib library (Polylib introduction)
Definition: Polylib.cpp:18
double gammaF(const double)
Calculate the Gamma function , , for integer values and halves.
Definition: Polylib.cpp:1216
void Dgj(double *D, const double *z, const int np, const double alpha, const double beta)
Compute the Derivative Matrix and its transpose associated with the Gauss-Jacobi zeros.
Definition: Polylib.cpp:558
void zwgrjm(double *z, double *w, const int np, const double alpha, const double beta)
Gauss-Radau-Jacobi zeros and weights with end point at z=-1.
Definition: Polylib.cpp:127
void JacZeros(const int n, double *a, double *b, const double alpha, const double beta)
Zero and Weight determination through the eigenvalues and eigenvectors of a tridiagonal matrix from t...
Definition: Polylib.cpp:1406
void zwgrjp(double *z, double *w, const int np, const double alpha, const double beta)
Gauss-Radau-Jacobi zeros and weights with end point at z=1.
Definition: Polylib.cpp:168
void Dglj(double *D, const double *z, const int np, const double alpha, const double beta)
Compute the Derivative Matrix associated with the Gauss-Lobatto-Jacobi zeros.
Definition: Polylib.cpp:718
void Imgrjp(double *im, const double *zgrj, const double *zm, const int nz, const int mz, const double alpha, const double beta)
Interpolation Operator from Gauss-Radau-Jacobi points (including z=1) to an arbitrary distrubtion at ...
Definition: Polylib.cpp:985
void zwglj(double *z, double *w, const int np, const double alpha, const double beta)
Gauss-Lobatto-Jacobi zeros and weights with end point at z=-1,1.
Definition: Polylib.cpp:207
void zwgj(double *z, double *w, const int np, const double alpha, const double beta)
Gauss-Jacobi zeros and weights.
Definition: Polylib.cpp:99
void Imgj(double *im, const double *zgj, const double *zm, const int nz, const int mz, const double alpha, const double beta)
Interpolation Operator from Gauss-Jacobi points to an arbitrary distribution at points zm.
Definition: Polylib.cpp:923
double hgrjm(const int i, const double z, const double *zgrj, const int np, const double alpha, const double beta)
Compute the value of the i th Lagrangian interpolant through the np Gauss-Radau-Jacobi points zgrj at...
Definition: Polylib.cpp:824
double gammaFracGammaF(const int, const double, const int, const double)
Calculate fraction of two Gamma functions, , for integer values and halves.
Definition: Polylib.cpp:1265
double hgj(const int i, const double z, const double *zgj, const int np, const double alpha, const double beta)
Compute the value of the i th Lagrangian interpolant through the np Gauss-Jacobi points zgj at the ar...
Definition: Polylib.cpp:789
void Imgrjm(double *im, const double *zgrj, const double *zm, const int nz, const int mz, const double alpha, const double beta)
Interpolation Operator from Gauss-Radau-Jacobi points (including z=-1) to an arbitrary distrubtion at...
Definition: Polylib.cpp:954
void zwgk(double *z, double *w, const int npt, const double alpha, const double beta)
Gauss-Kronrod-Jacobi zeros and weights.
Definition: Polylib.cpp:256
void Dgrjm(double *D, const double *z, const int np, const double alpha, const double beta)
Compute the Derivative Matrix and its transpose associated with the Gauss-Radau-Jacobi zeros with a z...
Definition: Polylib.cpp:605
std::complex< Nektar::NekDouble > ImagBesselComp(int n, std::complex< Nektar::NekDouble > y)
Calcualte the bessel function of the first kind with complex double input y. Taken from Numerical Rec...
Definition: Polylib.cpp:1731
void Dgrjp(double *D, const double *z, const int np, const double alpha, const double beta)
Compute the Derivative Matrix associated with the Gauss-Radau-Jacobi zeros with a zero at z=1.
Definition: Polylib.cpp:661
double hgrjp(const int i, const double z, const double *zgrj, const int np, const double alpha, const double beta)
Compute the value of the i th Lagrangian interpolant through the np Gauss-Radau-Jacobi points zgrj at...
Definition: Polylib.cpp:860
void zwrk(double *z, double *w, const int npt, const double alpha, const double beta)
Gauss-Radau-Kronrod-Jacobi zeros and weights.
Definition: Polylib.cpp:328
double hglj(const int i, const double z, const double *zglj, const int np, const double alpha, const double beta)
Compute the value of the i th Lagrangian interpolant through the np Gauss-Lobatto-Jacobi points zgrj ...
Definition: Polylib.cpp:895
void Imglj(double *im, const double *zglj, const double *zm, const int nz, const int mz, const double alpha, const double beta)
Interpolation Operator from Gauss-Lobatto-Jacobi points to an arbitrary distrubtion at points zm.
Definition: Polylib.cpp:1016
void jacobd(const int np, const double *z, double *polyd, const int n, const double alpha, const double beta)
Calculate the derivative of Jacobi polynomials.
Definition: Polylib.cpp:1187
void jacobfd(const int np, const double *z, double *poly_in, double *polyd, const int n, const double alpha, const double beta)
Routine to calculate Jacobi polynomials, , and their first derivative, .
Definition: Polylib.cpp:1074
void zwlk(double *z, double *w, const int npt, const double alpha, const double beta)
Gauss-Lobatto-Kronrod-Jacobi zeros and weights.
Definition: Polylib.cpp:434