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