v0.8.23
Public Member Functions | Public Attributes | List of all members
MinimalSurfaceEquation::MinimalSurfaceElement::OpAssmebleBcLhs Struct Reference

Integrate vector on lhs,. More...

#include <users_modules/minimal_surface_equation/src/MinimalSurfaceElement.hpp>

Inherits UserDataOperator.

Public Member Functions

 OpAssmebleBcLhs (const string field_name, Mat m_a)
 
PetscErrorCode doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
 

Public Attributes

Mat mA
 
MatrixDouble K
 
MatrixDouble transK
 

Detailed Description

Integrate vector on lhs,.

\[ \mathbf{A} = \int_L \mathbf{N}^\textrm{T}\mathbf{N} \textrm{d}L \]

Examples
minimal_surface_area.cpp.

Definition at line 140 of file MinimalSurfaceElement.hpp.

Constructor & Destructor Documentation

◆ OpAssmebleBcLhs()

MinimalSurfaceEquation::MinimalSurfaceElement::OpAssmebleBcLhs::OpAssmebleBcLhs ( const string  field_name,
Mat  m_a 
)

Definition at line 143 of file MinimalSurfaceElement.hpp.

145  field_name, UserDataOperator::OPROWCOL),
146  mA(m_a) {}
ForcesAndSourcesCore::UserDataOperator UserDataOperator

Member Function Documentation

◆ doWork()

PetscErrorCode MinimalSurfaceEquation::MinimalSurfaceElement::OpAssmebleBcLhs::doWork ( int  row_side,
int  col_side,
EntityType  row_type,
EntityType  col_type,
DataForcesAndSourcesCore::EntData &  row_data,
DataForcesAndSourcesCore::EntData &  col_data 
)

Definition at line 148 of file MinimalSurfaceElement.hpp.

151  {
153 
154  int row_nb_dofs = row_data.getIndices().size();
155  if (row_nb_dofs == 0) {
157  }
158  int col_nb_dofs = col_data.getIndices().size();
159  if (col_nb_dofs == 0) {
161  }
162  K.resize(row_nb_dofs, col_nb_dofs, false);
163  K.clear();
164  int nb_gauss_pts = row_data.getN().size1();
165  for (int gg = 0; gg != nb_gauss_pts; gg++) {
166  double val = getLength() * getGaussPts()(1, gg);
167  VectorAdaptor row_n = row_data.getN(gg);
168  VectorAdaptor col_n = col_data.getN(gg);
169  noalias(K) += val * outer_prod(row_n, col_n);
170  }
171  CHKERR MatSetValues(mA, row_nb_dofs,
172  &*row_data.getIndices().data().begin(), col_nb_dofs,
173  &*col_data.getIndices().data().begin(),
174  &*K.data().begin(), ADD_VALUES);
175  // Matrix is symmetric, assemble off diagonal part
176  if (row_type != col_type || row_side != col_side) {
177  transK.resize(col_nb_dofs, row_nb_dofs, false);
178  noalias(transK) = trans(K);
179  CHKERR MatSetValues(mA, col_nb_dofs,
180  &*col_data.getIndices().data().begin(), row_nb_dofs,
181  &*row_data.getIndices().data().begin(),
182  &*transK.data().begin(), ADD_VALUES);
183  }
184 
186  }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
MoFEMErrorCode MatSetValues(Mat M, const DataForcesAndSourcesCore::EntData &row_data, const DataForcesAndSourcesCore::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition: Types.hpp:101
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

Member Data Documentation

◆ K

MatrixDouble MinimalSurfaceEquation::MinimalSurfaceElement::OpAssmebleBcLhs::K

Definition at line 147 of file MinimalSurfaceElement.hpp.

◆ mA

Mat MinimalSurfaceEquation::MinimalSurfaceElement::OpAssmebleBcLhs::mA

Definition at line 142 of file MinimalSurfaceElement.hpp.

◆ transK

MatrixDouble MinimalSurfaceEquation::MinimalSurfaceElement::OpAssmebleBcLhs::transK

Definition at line 147 of file MinimalSurfaceElement.hpp.


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