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
223 class IntrepidPolylib {
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.
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.