14#include <boost/math/quadrature/gauss_kronrod.hpp>
15using namespace boost::math::quadrature;
51 GAUSS>::OpSource<1, 3>;
53 GAUSS>::OpSource<1, 1>;
62constexpr double omega = 7.292 * 1e-5;
63constexpr double g = 9.80616;
64constexpr double mu = 1e4;
65constexpr double h0 = 1e4;
69constexpr double phi_0 = M_PI / 7;
71constexpr double phi_2 = M_PI / 4;
85 boost::shared_ptr<MatrixDouble> u_dot_ptr,
86 boost::shared_ptr<MatrixDouble> grad_u_ptr,
87 boost::shared_ptr<MatrixDouble> grad_h_ptr)
95 const double vol = getMeasure();
96 auto t_w = getFTensor0IntegrationWeight();
99 auto t_dot_u = getFTensor1FromMat<3>(*
uDotPtr);
100 auto t_u = getFTensor1FromMat<3>(*
uPtr);
101 auto t_grad_u = getFTensor2FromMat<3, 3>(*
uGradPtr);
102 auto t_grad_h = getFTensor1FromMat<3>(*
hGradPtr);
103 auto t_coords = getFTensor1CoordsAtGaussPts();
104 auto t_normal = getFTensor1NormalsAtGaussPts();
108 const double alpha = t_w * vol;
109 auto t_nf = getFTensor1FromArray<3, 3>(
locF);
114 const auto a = std::sqrt(t_coords(
i) * t_coords(
i));
115 const auto sin_fi = t_coords(2) /
a;
116 const auto f = 2 *
omega * sin_fi;
119 t_r(
i) = t_normal(
i);
124 t_P(
i,
j) = t_r(
i) * t_r(
j);
128 t_A(
i,
m) = levi_civita(
i,
j,
m) * t_r(
j);
130 t_rhs(
m) = t_Q(
m,
i) * (t_dot_u(
i) + t_grad_u(
i,
j) * t_u(
j) +
131 f * t_A(
i,
j) * t_u(
j) +
g * t_grad_h(
i));
132 t_rhs_grad(
m,
j) = t_Q(
m,
i) * (
mu * t_grad_u(
i,
j));
134 t_rhs(
m) += t_P(
m,
j) * t_u(
j);
137 for (; rr !=
nbRows / 3; ++rr) {
138 t_nf(
i) += alpha * t_row_base * t_rhs(
i);
139 t_nf(
i) += alpha * t_row_diff_base(
j) * t_rhs_grad(
i,
j);
162 boost::shared_ptr<MatrixDouble>
uPtr;
170 OpULhs_dU(
const std::string field_name_row,
const std::string field_name_col,
171 boost::shared_ptr<MatrixDouble> u_ptr,
172 boost::shared_ptr<MatrixDouble> grad_u_ptr)
183 const double vol = getMeasure();
184 auto t_w = getFTensor0IntegrationWeight();
187 auto t_coords = getFTensor1CoordsAtGaussPts();
188 auto t_normal = getFTensor1NormalsAtGaussPts();
190 auto t_u = getFTensor1FromMat<3>(*
uPtr);
191 auto t_grad_u = getFTensor2FromMat<3, 3>(*
uGradPtr);
193 auto get_t_mat = [&](
const int rr) {
202 const auto ts_a = getFEMethod()->ts_a;
206 const auto a = std::sqrt(t_coords(
i) * t_coords(
i));
207 const auto sin_fi = t_coords(2) /
a;
208 const auto f = 2 *
omega * sin_fi;
211 t_r(
i) = t_normal(
i);
216 t_P(
i,
j) = t_r(
i) * t_r(
j);
220 t_A(
i,
m) = levi_civita(
i,
j,
m) * t_r(
j);
224 t_Q(
m,
i) * (ts_a *
t_kd(
i,
j) + t_grad_u(
i,
j) + f * t_A(
i,
j)) +
227 const double alpha = t_w * vol;
230 for (; rr !=
nbRows / 3; rr++) {
234 auto t_mat = get_t_mat(3 * rr);
236 for (
int cc = 0; cc !=
nbCols / 3; cc++) {
237 t_mat(
i,
j) += (alpha * t_row_base * t_col_base) * t_rhs_du(
i,
j);
238 t_mat(
i,
j) += (alpha *
mu) * t_Q(
i,
j) *
239 (t_row_diff_base(
m) * t_col_diff_base(
m));
263 boost::shared_ptr<MatrixDouble>
uPtr;
269 OpULhs_dH(
const std::string field_name_row,
const std::string field_name_col)
280 const double vol = getMeasure();
282 auto t_w = getFTensor0IntegrationWeight();
286 auto t_normal = getFTensor1NormalsAtGaussPts();
288 auto get_t_vec = [&](
const int rr) {
301 t_r(
i) = t_normal(
i);
306 t_P(
i,
j) = t_r(
i) * t_r(
j);
309 const double alpha = t_w * vol;
312 for (; rr !=
nbRows / 3; rr++) {
313 auto t_vec = get_t_vec(3 * rr);
315 const double a = alpha *
g * t_row_base;
317 for (
int cc = 0; cc !=
nbCols; cc++) {
318 t_vec(
i) +=
a * (t_Q(
i,
m) * t_col_diff_base(
m));
416 PetscBool is_restart = PETSC_FALSE;
420 auto restart_vector = [&]() {
423 auto dm =
simple->getDM();
425 <<
"reading vector in binary from vector.dat ...";
427 PetscViewerBinaryOpen(PETSC_COMM_WORLD,
"vector.dat", FILE_MODE_READ,
441 const double e_n = exp(-4 / pow(
phi_1 -
phi_0, 2));
448 auto get_phi = [&](
const double x,
const double y,
const double z) {
450 const double r = std::sqrt(t_r(
i) * t_r(
i));
454 auto init_u_phi = [&](
const double phi) {
462 auto init_u = [&](
const double x,
const double y,
const double z) {
464 const double u_phi = init_u_phi(get_phi(x, y, z));
468 t_u(
i) = ((t_A(
i,
j) * t_r(
j)) * u_phi);
473 auto init_h = [&](
const double x,
const double y,
const double z) {
474 const double a = std::sqrt(x * x + y * y + z * z);
476 auto integral = [&](
const double fi) {
477 const double u_phi = init_u_phi(fi);
478 const auto f = 2 *
omega * sin(fi);
479 return a * u_phi * (
f + (tan(fi) /
a) * u_phi);
482 auto montain = [&](
const double lambda,
const double fi) {
490 const double fi = get_phi(x, y, z);
491 const double lambda = atan2(x, y);
495 h1 = gauss_kronrod<double, 32>::integrate(
496 integral,
phi_0, fi, 0, std::numeric_limits<float>::epsilon());
498 return h0 + (montain(
lambda, fi) - (h1 /
g));
501 auto set_domain_general = [&](
auto &pipeline) {
502 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
508 auto set_domain_rhs = [&](
auto &pipeline) {
513 auto set_domain_lhs = [&](
auto &pipeline) {
515 new OpMassUU(
"U",
"U", [](
double,
double,
double) {
return 1; }));
517 new OpMassHH(
"H",
"H", [](
double,
double,
double) {
return 1; }));
520 auto post_proc = [&]() {
523 auto dm =
simple->getDM();
525 auto post_proc_fe = boost::make_shared<PostProcEle>(
mField);
527 auto det_ptr = boost::make_shared<VectorDouble>();
528 auto jac_ptr = boost::make_shared<MatrixDouble>();
529 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
531 post_proc_fe->getOpPtrVector().push_back(
533 post_proc_fe->getOpPtrVector().push_back(
535 post_proc_fe->getOpPtrVector().push_back(
537 post_proc_fe->getOpPtrVector().push_back(
540 auto u_ptr = boost::make_shared<MatrixDouble>();
541 auto h_ptr = boost::make_shared<VectorDouble>();
542 auto pos_ptr = boost::make_shared<MatrixDouble>();
544 post_proc_fe->getOpPtrVector().push_back(
546 post_proc_fe->getOpPtrVector().push_back(
548 post_proc_fe->getOpPtrVector().push_back(
553 post_proc_fe->getOpPtrVector().push_back(
555 new OpPPMap(post_proc_fe->getPostProcMesh(),
556 post_proc_fe->getMapGaussPts(),
560 {{
"U", u_ptr}, {
"HO_POSITIONS", pos_ptr}},
569 CHKERR post_proc_fe->writeFile(
"out_init.h5m");
574 auto solve_init = [&]() {
579 auto solver = pipeline_mng->createKSP();
580 CHKERR KSPSetFromOptions(solver);
582 CHKERR KSPGetPC(solver, &pc);
583 PetscBool is_pcfs = PETSC_FALSE;
584 PetscObjectTypeCompare((PetscObject)pc, PCFIELDSPLIT, &is_pcfs);
585 if (is_pcfs == PETSC_TRUE) {
587 auto name_prb =
simple->getProblemName();
590 name_prb,
ROW,
"U", 0, 3, is_u);
591 CHKERR PCFieldSplitSetIS(pc, PETSC_NULL, is_u);
592 CHKERR PCFieldSplitSetType(pc, PC_COMPOSITE_ADDITIVE);
597 auto dm =
simple->getDM();
602 CHKERR VecGhostUpdateBegin(
D, INSERT_VALUES, SCATTER_FORWARD);
603 CHKERR VecGhostUpdateEnd(
D, INSERT_VALUES, SCATTER_FORWARD);
617 set_domain_general(pipeline_mng->getOpDomainRhsPipeline());
618 set_domain_rhs(pipeline_mng->getOpDomainRhsPipeline());
619 set_domain_general(pipeline_mng->getOpDomainLhsPipeline());
620 set_domain_lhs(pipeline_mng->getOpDomainLhsPipeline());
626 pipeline_mng->getOpDomainRhsPipeline().clear();
627 pipeline_mng->getOpDomainLhsPipeline().clear();
640 auto set_domain_general = [&](
auto &pipeline) {
641 auto det_ptr = boost::make_shared<VectorDouble>();
642 auto jac_ptr = boost::make_shared<MatrixDouble>();
643 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
652 auto set_domain_rhs = [&](
auto &pipeline) {
653 auto dot_u_ptr = boost::make_shared<MatrixDouble>();
654 auto u_ptr = boost::make_shared<MatrixDouble>();
655 auto grad_u_ptr = boost::make_shared<MatrixDouble>();
656 auto div_u_ptr = boost::make_shared<VectorDouble>();
657 auto dot_h_ptr = boost::make_shared<VectorDouble>();
658 auto grad_h_ptr = boost::make_shared<MatrixDouble>();
671 "H", dot_h_ptr, [](
double,
double,
double) {
return 1.; }));
673 "H", div_u_ptr, [](
double,
double,
double) {
return h0; }));
674 pipeline.push_back(
new OpConvectiveH(
"H", u_ptr, grad_h_ptr));
676 new OpURhs(
"U", u_ptr, dot_u_ptr, grad_u_ptr, grad_h_ptr));
679 auto set_domain_lhs = [&](
auto &pipeline) {
680 auto u_ptr = boost::make_shared<MatrixDouble>();
681 auto grad_u_ptr = boost::make_shared<MatrixDouble>();
682 auto grad_h_ptr = boost::make_shared<MatrixDouble>();
689 pipeline.push_back(
new OpMassHH(
"H",
"H", [&](
double,
double,
double) {
693 "H",
"U", [](
const double,
const double,
const double) {
return h0; },
699 pipeline.push_back(
new OpULhs_dU(
"U",
"U", u_ptr, grad_u_ptr));
700 pipeline.push_back(
new OpULhs_dH(
"U",
"H"));
711 set_domain_general(pipeline_mng->getOpDomainRhsPipeline());
712 set_domain_general(pipeline_mng->getOpDomainLhsPipeline());
714 set_domain_rhs(pipeline_mng->getOpDomainRhsPipeline());
715 set_domain_lhs(pipeline_mng->getOpDomainLhsPipeline());
739 "out_step_" + boost::lexical_cast<std::string>(
ts_step) +
".h5m");
741 <<
"writing vector in binary to vector.dat ...";
743 PetscViewerBinaryOpen(PETSC_COMM_WORLD,
"vector.dat", FILE_MODE_WRITE,
745 VecView(
ts_u, viewer);
746 PetscViewerDestroy(&viewer);
753 boost::shared_ptr<PostProcEle>
postProc;
762 auto dm =
simple->getDM();
764 auto set_initial_step = [&](
auto ts) {
769 CHKERR TSSetStepNumber(ts, step);
773 auto set_fieldsplit_preconditioner_ksp = [&](
auto ksp) {
776 CHKERR KSPGetPC(ksp, &pc);
777 PetscBool is_pcfs = PETSC_FALSE;
778 PetscObjectTypeCompare((PetscObject)pc, PCFIELDSPLIT, &is_pcfs);
779 if (is_pcfs == PETSC_TRUE) {
781 auto name_prb =
simple->getProblemName();
784 name_prb,
ROW,
"U", 0, 3, is_u);
785 CHKERR PCFieldSplitSetIS(pc, PETSC_NULL, is_u);
791 auto get_fe_post_proc = [&]() {
792 auto post_proc_fe = boost::make_shared<PostProcEle>(
mField);
794 auto det_ptr = boost::make_shared<VectorDouble>();
795 auto jac_ptr = boost::make_shared<MatrixDouble>();
796 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
798 post_proc_fe->getOpPtrVector().push_back(
800 post_proc_fe->getOpPtrVector().push_back(
802 post_proc_fe->getOpPtrVector().push_back(
804 post_proc_fe->getOpPtrVector().push_back(
806 post_proc_fe->getOpPtrVector().push_back(
809 auto u_ptr = boost::make_shared<MatrixDouble>();
810 auto h_ptr = boost::make_shared<VectorDouble>();
811 auto pos_ptr = boost::make_shared<MatrixDouble>();
813 auto grad_u_ptr = boost::make_shared<MatrixDouble>();
814 auto grad_h_ptr = boost::make_shared<MatrixDouble>();
816 post_proc_fe->getOpPtrVector().push_back(
818 post_proc_fe->getOpPtrVector().push_back(
820 post_proc_fe->getOpPtrVector().push_back(
823 post_proc_fe->getOpPtrVector().push_back(
825 post_proc_fe->getOpPtrVector().push_back(
830 post_proc_fe->getOpPtrVector().push_back(
833 post_proc_fe->getPostProcMesh(), post_proc_fe->getMapGaussPts(),
837 {{
"U", u_ptr}, {
"HO_POSITIONS", pos_ptr}, {
"GRAD_H", grad_h_ptr}},
839 {{
"GRAD_U", grad_u_ptr}}, {}
848 auto set_fieldsplit_preconditioner_ts = [&](
auto solver) {
851 CHKERR TSGetSNES(solver, &snes);
853 CHKERR SNESGetKSP(snes, &ksp);
854 CHKERR set_fieldsplit_preconditioner_ksp(ksp);
859 ts = pipeline_mng->createTSIM();
861 boost::shared_ptr<FEMethod> null_fe;
862 auto monitor_ptr = boost::make_shared<Monitor>(dm, get_fe_post_proc());
864 null_fe, monitor_ptr);
869 CHKERR TSSetExactFinalTime(ts, TS_EXACTFINALTIME_MATCHSTEP);
874 CHKERR TSSetSolution(ts, T);
875 CHKERR TSSetFromOptions(ts);
876 CHKERR set_fieldsplit_preconditioner_ts(ts);
878 CHKERR set_initial_step(ts);
880 CHKERR TSGetTime(ts, &ftime);
902int main(
int argc,
char *argv[]) {
909 auto core_log = logging::core::get();
917 DMType dm_name =
"DMMOFEM";
922 moab::Core mb_instance;
923 moab::Interface &moab = mb_instance;
void simple(double P1[], double P2[], double P3[], double c[], const int N)
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
Tensor1< T, Tensor_Dim > normalize()
#define CATCH_ERRORS
Catch errors.
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, SPACE_DIM > OpMass
auto init_h
Initialisation function.
FormsIntegrators< DomainEleOp >::Assembly< A >::BiLinearForm< I >::OpMixScalarTimesDiv< SPACE_DIM, COORD_TYPE > OpMixScalarTimesDiv
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
PetscErrorCode DMoFEMLoopFiniteElements(DM dm, const char fe_name[], MoFEM::FEMethod *method, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr())
Executes FEMethod for finite elements in DM.
auto createDMVector(DM dm)
Get smart vector from DM.
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
#define MOFEM_LOG(channel, severity)
Log.
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
constexpr std::enable_if<(Dim0<=2 &&Dim1<=2), Tensor2_Expr< Levi_Civita< T >, T, Dim0, Dim1, i, j > >::type levi_civita(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
levi_civita functions to make for easy adhoc use
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
implementation of Data Operators for Forces and Sources
PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
Set Monitor To TS solver.
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
SmartPetscObj< Vec > vectorDuplicate(Vec vec)
Create duplicate vector of smart vector.
OpCalculateScalarFieldValuesFromPetscVecImpl< PetscData::CTX_SET_X_T > OpCalculateScalarFieldValuesDot
OpCalculateHOJacForFaceImpl< 3 > OpCalculateHOJacForFaceEmbeddedIn3DSpace
const double u0
inital vale on blocksets
OpPostProcMapInMoab< SPACE_DIM, SPACE_DIM > OpPPMap
constexpr auto field_name
static constexpr int approx_order
OpBaseImpl< PETSC, EdgeEleOp > OpBase
static char help[]
[Check]
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, 1 > OpMassHH
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, 3 > OpSourceU
constexpr double beta_montain
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, 3 > OpMassUU
FTensor::Index< 'l', 3 > l
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpConvectiveTermLhsDu< 1, 1, 3 > OpConvectiveH_dU
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpConvectiveTermRhs< 1, 1, 3 > OpConvectiveH
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMixScalarTimesDiv< 3 > OpBaseDivU
OpBaseTimesDotH OpBaseTimesDivU
FTensor::Index< 'j', 3 > j
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpConvectiveTermLhsDy< 1, 1, 3 > OpConvectiveH_dGradH
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesScalar< 1 > OpBaseTimesDotH
FTensor::Index< 'i', 3 > i
FTensor::Index< 'm', 3 > m
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpSource< 1, 1 > OpSourceH
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMixVectorTimesGrad< 1, 3, 3 > OpBaseGradH
constexpr double alpha_montain
MoFEMErrorCode boundaryCondition()
MoFEMErrorCode assembleSystem()
MoFEMErrorCode readMesh()
MoFEMErrorCode checkResults()
MoFEMErrorCode solveSystem()
MoFEMErrorCode createCommonData()
Example(MoFEM::Interface &m_field)
MoFEMErrorCode runProblem()
MoFEM::Interface & mField
boost::shared_ptr< FEMethod > domianRhsFEPtr
MoFEMErrorCode setupProblem()
MoFEMErrorCode outputResults()
boost::shared_ptr< FEMethod > domianLhsFEPtr
Simple interface for fast problem set-up.
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1DiffN(const FieldApproximationBase base)
Get derivatives of base functions.
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.
structure for User Loop Methods on finite elements
default operator for TRI element
Section manager is used to create indexes and sections.
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
VectorDouble locF
local entity vector
int nbRows
number of dofs on rows
int nbIntegrationPts
number of integration points
MatrixDouble locMat
local entity block matrix
int nbCols
number if dof on column
int nbRowBaseFunctions
number or row base functions
Calculate field values (template specialization) for tensor field rank 1, i.e. vector field.
Calculate HO coordinates at gauss points.
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
Get value at integration points for scalar field.
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
Approximate field values for given petsc vector.
Get values at integration pts for tensor filed rank 1, i.e. vector field.
Calculate normals at Gauss points of triangle element.
Post post-proc data at points from hash maps.
Modify integration weights on face to take in account higher-order geometry.
PipelineManager interface.
Simple interface for fast problem set-up.
intrusive_ptr for managing petsc objects
PetscInt ts_step
time step number
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
[Push operators to pipeline]
Monitor(SmartPetscObj< DM > dm, boost::shared_ptr< PostProcEle > post_proc)
MoFEMErrorCode postProcess()
function is run at the end of loop
boost::shared_ptr< PostProcEle > postProc
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
Integrate grad-grad operator.
OpULhs_dH(const std::string field_name_row, const std::string field_name_col)
boost::shared_ptr< MatrixDouble > uPtr
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
Integrate grad-grad operator.
boost::shared_ptr< MatrixDouble > uGradPtr
OpULhs_dU(const std::string field_name_row, const std::string field_name_col, boost::shared_ptr< MatrixDouble > u_ptr, boost::shared_ptr< MatrixDouble > grad_u_ptr)
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data)
Class dedicated to integrate operator.
boost::shared_ptr< MatrixDouble > uDotPtr
boost::shared_ptr< MatrixDouble > uGradPtr
boost::shared_ptr< MatrixDouble > hGradPtr
boost::shared_ptr< MatrixDouble > uPtr
OpURhs(const std::string field_name, boost::shared_ptr< MatrixDouble > u_ptr, boost::shared_ptr< MatrixDouble > u_dot_ptr, boost::shared_ptr< MatrixDouble > grad_u_ptr, boost::shared_ptr< MatrixDouble > grad_h_ptr)