Intrepid
Intrepid_Polylib.hpp
Go to the documentation of this file.
1/*
2// @HEADER
3// ************************************************************************
4//
5// Intrepid Package
6// Copyright (2007) Sandia Corporation
7//
8// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
9// license for use of this work by or on behalf of the U.S. Government.
10//
11// Redistribution and use in source and binary forms, with or without
12// modification, are permitted provided that the following conditions are
13// met:
14//
15// 1. Redistributions of source code must retain the above copyright
16// notice, this list of conditions and the following disclaimer.
17//
18// 2. Redistributions in binary form must reproduce the above copyright
19// notice, this list of conditions and the following disclaimer in the
20// documentation and/or other materials provided with the distribution.
21//
22// 3. Neither the name of the Corporation nor the names of the
23// contributors may be used to endorse or promote products derived from
24// this software without specific prior written permission.
25//
26// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37//
38// Questions? Contact Pavel Bochev (pbboche@sandia.gov)
39// Denis Ridzal (dridzal@sandia.gov), or
40// Kara Peterson (kjpeter@sandia.gov)
41//
42// ************************************************************************
43// @HEADER
44*/
45
47//
48// File: Intrepid_Polylib.hpp
49//
50// For more information, please see: http://www.nektar.info
51//
52// The MIT License
53//
54// Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
55// Department of Aeronautics, Imperial College London (UK), and Scientific
56// Computing and Imaging Institute, University of Utah (USA).
57//
58// License for the specific language governing rights and limitations under
59// Permission is hereby granted, free of charge, to any person obtaining a
60// copy of this software and associated documentation files (the "Software"),
61// to deal in the Software without restriction, including without limitation
62// the rights to use, copy, modify, merge, publish, distribute, sublicense,
63// and/or sell copies of the Software, and to permit persons to whom the
64// Software is furnished to do so, subject to the following conditions:
65//
66// The above copyright notice and this permission notice shall be included
67// in all copies or substantial portions of the Software.
68//
69// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
70// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
71// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
72// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
73// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
74// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
75// DEALINGS IN THE SOFTWARE.
76//
77// Description:
78// This file is redistributed with the Intrepid package. It should be used
79// in accordance with the above MIT license, at the request of the authors.
80// This file is NOT covered by the usual Intrepid/Trilinos LGPL license.
81//
82// Origin: Nektar++ library, http://www.nektar.info, downloaded on
83// March 10, 2009.
84//
86
87
95#ifndef INTREPID_POLYLIB_HPP
96#define INTREPID_POLYLIB_HPP
97
98#include "Intrepid_ConfigDefs.hpp"
99#include "Intrepid_Types.hpp"
100#include "Teuchos_Assert.hpp"
101
102namespace Intrepid {
103
197 enum EIntrepidPLPoly {
198 PL_GAUSS=0,
199 PL_GAUSS_RADAU_LEFT,
200 PL_GAUSS_RADAU_RIGHT,
201 PL_GAUSS_LOBATTO,
202 PL_MAX
203 };
204
205 inline EIntrepidPLPoly & operator++(EIntrepidPLPoly &type) {
206 return type = static_cast<EIntrepidPLPoly>(type+1);
207 }
208
209 inline EIntrepidPLPoly operator++(EIntrepidPLPoly &type, int) {
210 EIntrepidPLPoly oldval = type;
211 ++type;
212 return oldval;
213 }
214
215
224
225 public:
226
227 /* Points and weights */
228
236 template<class Scalar>
237 static void zwgj (Scalar *z, Scalar *w, const int np, const Scalar alpha, const Scalar beta);
238
239
247 template<class Scalar>
248 static void zwgrjm (Scalar *z, Scalar *w, const int np, const Scalar alpha, const Scalar beta);
249
250
258 template<class Scalar>
259 static void zwgrjp (Scalar *z, Scalar *w, const int np, const Scalar alpha, const Scalar beta);
260
261
269 template<class Scalar>
270 static void zwglj (Scalar *z, Scalar *w, const int np, const Scalar alpha, const Scalar beta);
271
272
273
274 /* Derivative operators */
275
284 template<class Scalar>
285 static void Dgj (Scalar *D, const Scalar *z, const int np, const Scalar alpha, const Scalar beta);
286
287
296 template<class Scalar>
297 static void Dgrjm (Scalar *D, const Scalar *z, const int np, const Scalar alpha, const Scalar beta);
298
299
308 template<class Scalar>
309 static void Dgrjp (Scalar *D, const Scalar *z, const int np, const Scalar alpha, const Scalar beta);
310
311
320 template<class Scalar>
321 static void Dglj (Scalar *D, const Scalar *z, const int np, const Scalar alpha, const Scalar beta);
322
323
324
325 /* Lagrangian interpolants */
326
346 template<class Scalar>
347 static Scalar hgj (const int i, const Scalar z, const Scalar *zgj,
348 const int np, const Scalar alpha, const Scalar beta);
349
350
370 template<class Scalar>
371 static Scalar hgrjm (const int i, const Scalar z, const Scalar *zgrj,
372 const int np, const Scalar alpha, const Scalar beta);
373
374
394 template<class Scalar>
395 static Scalar hgrjp (const int i, const Scalar z, const Scalar *zgrj,
396 const int np, const Scalar alpha, const Scalar beta);
397
398
418 template<class Scalar>
419 static Scalar hglj (const int i, const Scalar z, const Scalar *zglj,
420 const int np, const Scalar alpha, const Scalar beta);
421
422
423
424 /* Interpolation operators */
425
436 template<class Scalar>
437 static void Imgj (Scalar *im, const Scalar *zgj, const Scalar *zm, const int nz,
438 const int mz, const Scalar alpha, const Scalar beta);
439
440
451 template<class Scalar>
452 static void Imgrjm(Scalar *im, const Scalar *zgrj, const Scalar *zm, const int nz,
453 const int mz, const Scalar alpha, const Scalar beta);
454
455
466 template<class Scalar>
467 static void Imgrjp(Scalar *im, const Scalar *zgrj, const Scalar *zm, const int nz,
468 const int mz, const Scalar alpha, const Scalar beta);
469
470
481 template<class Scalar>
482 static void Imglj (Scalar *im, const Scalar *zglj, const Scalar *zm, const int nz,
483 const int mz, const Scalar alpha, const Scalar beta);
484
485
486 /* Polynomial functions */
487
527 template<class Scalar>
528 static void jacobfd (const int np, const Scalar *z, Scalar *poly_in, Scalar *polyd,
529 const int n, const Scalar alpha, const Scalar beta);
530
531
545 template<class Scalar>
546 static void jacobd (const int np, const Scalar *z, Scalar *polyd, const int n,
547 const Scalar alpha, const Scalar beta);
548
549
550
551 /* Helper functions. */
552
559 template<class Scalar>
560 static void Jacobz (const int n, Scalar *z, const Scalar alpha, const Scalar beta);
561
562
585 template<class Scalar>
586 static void JacZeros (const int n, Scalar *a, const Scalar alpha, const Scalar beta);
587
588
612 template<class Scalar>
613 static void TriQL (const int n, Scalar *d, Scalar *e);
614
615
625 template<class Scalar>
626 static Scalar gammaF (const Scalar x);
627
628
629 }; // class IntrepidPolylib
630
631} // end of Intrepid namespace
632
633// include templated definitions
635
636#endif
Definition file for a set of functions providing orthogonal polynomial polynomial calculus and interp...
Contains definitions of custom data types in Intrepid.
Providing orthogonal polynomial calculus and interpolation, created by Spencer Sherwin,...
static void Jacobz(const int n, Scalar *z, const Scalar alpha, const Scalar beta)
Calculate the n zeros, z, of the Jacobi polynomial, i.e. .
static void zwgrjp(Scalar *z, Scalar *w, const int np, const Scalar alpha, const Scalar beta)
Gauss-Radau-Jacobi zeros and weights with end point at z=1.
static void JacZeros(const int n, Scalar *a, const Scalar alpha, const Scalar beta)
Zero determination through the eigenvalues of a tridiagonal matrix from the three term recursion rela...
static void Imgj(Scalar *im, const Scalar *zgj, const Scalar *zm, const int nz, const int mz, const Scalar alpha, const Scalar beta)
Interpolation Operator from Gauss-Jacobi points to an arbitrary distribution at points zm.
static Scalar hgrjp(const int i, const Scalar z, const Scalar *zgrj, const int np, const Scalar alpha, const Scalar beta)
Compute the value of the i th Lagrangian interpolant through the np Gauss-Radau-Jacobi points zgrj at...
static void jacobfd(const int np, const Scalar *z, Scalar *poly_in, Scalar *polyd, const int n, const Scalar alpha, const Scalar beta)
Routine to calculate Jacobi polynomials, , and their first derivative, .
static void zwgrjm(Scalar *z, Scalar *w, const int np, const Scalar alpha, const Scalar beta)
Gauss-Radau-Jacobi zeros and weights with end point at z=-1.
static Scalar hgrjm(const int i, const Scalar z, const Scalar *zgrj, const int np, const Scalar alpha, const Scalar beta)
Compute the value of the i th Lagrangian interpolant through the np Gauss-Radau-Jacobi points zgrj at...
static void Imgrjp(Scalar *im, const Scalar *zgrj, const Scalar *zm, const int nz, const int mz, const Scalar alpha, const Scalar beta)
Interpolation Operator from Gauss-Radau-Jacobi points (including z=1) to an arbitrary distrubtion at ...
static void Dgrjp(Scalar *D, const Scalar *z, const int np, const Scalar alpha, const Scalar beta)
Compute the Derivative Matrix associated with the Gauss-Radau-Jacobi zeros with a zero at z=1.
static void Imglj(Scalar *im, const Scalar *zglj, const Scalar *zm, const int nz, const int mz, const Scalar alpha, const Scalar beta)
Interpolation Operator from Gauss-Lobatto-Jacobi points to an arbitrary distrubtion at points zm.
static void zwglj(Scalar *z, Scalar *w, const int np, const Scalar alpha, const Scalar beta)
Gauss-Lobatto-Jacobi zeros and weights with end point at z=-1,1.
static void Dgj(Scalar *D, const Scalar *z, const int np, const Scalar alpha, const Scalar beta)
Compute the Derivative Matrix and its transpose associated with the Gauss-Jacobi zeros.
static void TriQL(const int n, Scalar *d, Scalar *e)
QL algorithm for symmetric tridiagonal matrix.
static void Dgrjm(Scalar *D, const Scalar *z, const int np, const Scalar alpha, const Scalar beta)
Compute the Derivative Matrix and its transpose associated with the Gauss-Radau-Jacobi zeros with a z...
static void Dglj(Scalar *D, const Scalar *z, const int np, const Scalar alpha, const Scalar beta)
Compute the Derivative Matrix associated with the Gauss-Lobatto-Jacobi zeros.
static void Imgrjm(Scalar *im, const Scalar *zgrj, const Scalar *zm, const int nz, const int mz, const Scalar alpha, const Scalar beta)
Interpolation Operator from Gauss-Radau-Jacobi points (including z=-1) to an arbitrary distrubtion at...
static Scalar hglj(const int i, const Scalar z, const Scalar *zglj, const int np, const Scalar alpha, const Scalar beta)
Compute the value of the i th Lagrangian interpolant through the np Gauss-Lobatto-Jacobi points zglj ...
static void jacobd(const int np, const Scalar *z, Scalar *polyd, const int n, const Scalar alpha, const Scalar beta)
Calculate the derivative of Jacobi polynomials.
static Scalar gammaF(const Scalar x)
Calculate the Gamma function , , for integer values x and halves.
static Scalar hgj(const int i, const Scalar z, const Scalar *zgj, const int np, const Scalar alpha, const Scalar beta)
Compute the value of the i th Lagrangian interpolant through the np Gauss-Jacobi points zgj at the ar...
static void zwgj(Scalar *z, Scalar *w, const int np, const Scalar alpha, const Scalar beta)
Gauss-Jacobi zeros and weights.