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

Assemble residual. More...

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

Inherits UserDataOperator.

Collaboration diagram for MinimalSurfaceEquation::MinimalSurfaceElement::OpAssembleResidaul:
[legend]

Public Member Functions

 OpAssembleResidaul (const string field_name, CommonData &common_data)
 
PetscErrorCode doWork (int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
 

Public Attributes

CommonDatacommonData
 
VectorDouble rEsidual
 

Detailed Description

Assemble residual.

Examples
minimal_surface_area.cpp.

Definition at line 296 of file MinimalSurfaceElement.hpp.

Constructor & Destructor Documentation

◆ OpAssembleResidaul()

MinimalSurfaceEquation::MinimalSurfaceElement::OpAssembleResidaul::OpAssembleResidaul ( const string  field_name,
CommonData common_data 
)

Definition at line 299 of file MinimalSurfaceElement.hpp.

301  field_name, UserDataOperator::OPROW),
302  commonData(common_data) {}
ForcesAndSourcesCore::UserDataOperator UserDataOperator

Member Function Documentation

◆ doWork()

PetscErrorCode MinimalSurfaceEquation::MinimalSurfaceElement::OpAssembleResidaul::doWork ( int  side,
EntityType  type,
DataForcesAndSourcesCore::EntData &  data 
)

Definition at line 304 of file MinimalSurfaceElement.hpp.

305  {
307 
308  int nb_dofs = data.getIndices().size();
309  if (nb_dofs == 0) {
311  }
312  rEsidual.resize(nb_dofs, false);
313  rEsidual.clear();
314  int nb_gauss_pts = data.getN().size1();
315  for (int gg = 0; gg != nb_gauss_pts; gg++) {
316  double val = getGaussPts()(2, gg) * getArea();
317  ublas::matrix_row<MatrixDouble> an_by_grad_u(commonData.aNbyGradU, gg);
318  noalias(rEsidual) += val * prod(data.getDiffN(gg), an_by_grad_u);
319  }
320  CHKERR VecSetValues(getFEMethod()->snes_f, data.getIndices().size(),
321  &*data.getIndices().data().begin(), &rEsidual[0],
322  ADD_VALUES);
324  }
#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 VecSetValues(Vec V, const DataForcesAndSourcesCore::EntData &data, const double *ptr, InsertMode iora)
Assemble PETSc vector.
#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

◆ commonData

CommonData& MinimalSurfaceEquation::MinimalSurfaceElement::OpAssembleResidaul::commonData

Definition at line 298 of file MinimalSurfaceElement.hpp.

◆ rEsidual

VectorDouble MinimalSurfaceEquation::MinimalSurfaceElement::OpAssembleResidaul::rEsidual

Definition at line 303 of file MinimalSurfaceElement.hpp.


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