v0.13.1
Public Member Functions | Private Attributes | List of all members
PoissonOps::OpLhsUU Struct Reference

#include <users_modules/softmech/poisson_problem/src/PoissonOperators.hpp>

Inheritance diagram for PoissonOps::OpLhsUU:
[legend]
Collaboration diagram for PoissonOps::OpLhsUU:
[legend]

Public Member Functions

 OpLhsUU (std::string domain_field1, std::string domain_field2)
 
MoFEMErrorCode doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
 

Private Attributes

MatrixDouble locLhs
 
MatrixDouble transLocLhs
 

Detailed Description

Definition at line 23 of file PoissonOperators.hpp.

Constructor & Destructor Documentation

◆ OpLhsUU()

PoissonOps::OpLhsUU::OpLhsUU ( std::string  domain_field1,
std::string  domain_field2 
)

Definition at line 25 of file PoissonOperators.hpp.

26 : OpVolEle(domain_field1, domain_field2, OpVolEle::OPROWCOL) {
27 sYmm = true;
28 }
FaceElementForcesAndSourcesCore::UserDataOperator OpVolEle
@ OPROWCOL
operator doWork is executed on FE rows &columns

Member Function Documentation

◆ doWork()

MoFEMErrorCode PoissonOps::OpLhsUU::doWork ( int  row_side,
int  col_side,
EntityType  row_type,
EntityType  col_type,
EntData row_data,
EntData col_data 
)

Definition at line 30 of file PoissonOperators.hpp.

32 {
34
35 const int nb_row_dofs = row_data.getIndices().size();
36 const int nb_col_dofs = col_data.getIndices().size();
37
38 if (nb_row_dofs && nb_col_dofs) {
39 locLhs.resize(nb_row_dofs, nb_col_dofs, false);
40 locLhs.clear();
41 const int nb_integration_pts = getGaussPts().size2();
42 auto t_w = getFTensor0IntegrationWeight();
43
44 auto t_row_diff_base = row_data.getFTensor1DiffN<2>();
45
46 const double vol = getMeasure();
47 for (int gg = 0; gg != nb_integration_pts; ++gg) {
48 const double a = t_w * vol;
49 for (int rr = 0; rr != nb_row_dofs; ++rr) {
50 auto t_col_diff_base = col_data.getFTensor1DiffN<2>(gg, 0);
51 for (int cc = 0; cc != nb_col_dofs; ++cc) {
52 locLhs(rr, cc) += a * t_row_diff_base(i) * t_col_diff_base(i);
53 ++t_col_diff_base;
54 } // endFor cc
55 ++t_row_diff_base;
56 } // endFor rr
57 ++t_w;
58 } // endFor gg
59 CHKERR MatSetValues(getFEMethod()->ksp_B, row_data, col_data,
60 &locLhs(0, 0), ADD_VALUES);
61 if (row_side != col_side || row_type != col_type) {
62 transLocLhs.resize(nb_col_dofs, nb_row_dofs, false);
63 noalias(transLocLhs) = trans(locLhs);
64 CHKERR MatSetValues(getFEMethod()->ksp_B, col_data, row_data,
65 &transLocLhs(0, 0), ADD_VALUES);
66 } // (row_side != col_side || row_type != col_type)
67 } // endIf (nb_row_dofs && nb_col_dofs)
68
70 }
constexpr double a
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#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
FTensor::Index< 'i', 2 > i
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1DiffN(const FieldApproximationBase base)
Get derivatives of base functions.
const VectorInt & getIndices() const
Get global indices of dofs on entity.

Member Data Documentation

◆ locLhs

MatrixDouble PoissonOps::OpLhsUU::locLhs
private

Definition at line 73 of file PoissonOperators.hpp.

◆ transLocLhs

MatrixDouble PoissonOps::OpLhsUU::transLocLhs
private

Definition at line 73 of file PoissonOperators.hpp.


The documentation for this struct was generated from the following file: