v0.14.0
Loading...
Searching...
No Matches
JacobiPolynomial.cpp
Go to the documentation of this file.
1/** \file legendrepolynomial.cpp
2 * \brief implementation of multi-grid solver for p- adaptivity
3 */
4
5
6
7namespace MoFEM {
8
10JacobiPolynomialCtx::query_interface(boost::typeindex::type_index type_index,
11 UnknownInterface **iface) const {
12 *iface = const_cast<JacobiPolynomialCtx *>(this);
13 return 0;
14}
15
17JacobiPolynomial::query_interface(boost::typeindex::type_index type_index,
18 UnknownInterface **iface) const {
20 *iface = const_cast<JacobiPolynomial *>(this);
22}
23
24template <class TYPE>
26 TYPE *ctx) {
28 ctx->baseFunPtr->resize(pts_x.size2(), ctx->P + 1, false);
29 ctx->baseDiffFunPtr->resize(pts_x.size2(), ctx->dIm * (ctx->P + 1), false);
30 if (pts_x.size1() != pts_t.size1() || pts_x.size2() != pts_t.size2()) {
31 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
32 "Inconsistent size of arguments");
33 }
34 double *l = NULL;
35 double *diff_l = NULL;
36 for (unsigned int gg = 0; gg < pts_x.size2(); gg++) {
37 if (ctx->baseFunPtr)
38 l = &((*ctx->baseFunPtr)(gg, 0));
39 if (ctx->baseDiffFunPtr)
40 diff_l = &((*ctx->baseDiffFunPtr)(gg, 0));
41 ierr = (ctx->basePolynomialsType1)(ctx->P, ctx->aLpha, pts_x(0, gg),
42 pts_t(0, gg), ctx->diffX, ctx->diffT, l,
43 diff_l, ctx->dIm);
45 }
47}
48
51 boost::shared_ptr<BaseFunctionCtx> ctx_ptr) {
53 auto ctx = ctx_ptr->getInterface<JacobiPolynomialCtx>();
54 CHKERR get_value(pts_x, pts_t, ctx);
56}
57
59 boost::typeindex::type_index type_index, UnknownInterface **iface) const {
60 *iface = const_cast<IntegratedJacobiPolynomialCtx *>(this);
61 return 0;
62}
63
65 boost::typeindex::type_index type_index, UnknownInterface **iface) const {
66 *iface = const_cast<IntegratedJacobiPolynomial *>(this);
67 return 0;
68}
69
71 MatrixDouble &pts_x, MatrixDouble &pts_t,
72 boost::shared_ptr<BaseFunctionCtx> ctx_ptr) {
74 auto ctx = ctx_ptr->getInterface<IntegratedJacobiPolynomialCtx>();
75 CHKERR get_value(pts_x, pts_t, ctx);
77}
78
79} // namespace MoFEM
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:483
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
FTensor::Index< 'l', 3 > l
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:77
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
static MoFEMErrorCode get_value(MatrixDouble &pts_x, MatrixDouble &pts_t, TYPE *ctx)
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
MoFEMErrorCode getValue(MatrixDouble &pts_x, MatrixDouble &pts_t, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Class used to give arguments to Legendre base functions.
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Calculating Legendre base functions.
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
MoFEMErrorCode getValue(MatrixDouble &pts_x, MatrixDouble &pts_t, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
base class for all interface classes