v0.9.2
Public Member Functions | Public Attributes | List of all members
MixTransport::UnsaturatedFlowElement::OpIntegrateFluxes Struct Reference

#include <users_modules/basic_finite_elements/mix_transport/src/UnsaturatedFlow.hpp>

Inheritance diagram for MixTransport::UnsaturatedFlowElement::OpIntegrateFluxes:
[legend]
Collaboration diagram for MixTransport::UnsaturatedFlowElement::OpIntegrateFluxes:
[legend]

Public Member Functions

 OpIntegrateFluxes (UnsaturatedFlowElement &ctx, const std::string fluxes_name)
 Constructor. More...
 
MoFEMErrorCode doWork (int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
 Integrate boundary condition. More...
 

Public Attributes

UnsaturatedFlowElementcTx
 
FTensor::Index< 'i', 3 > i
 

Detailed Description

Definition at line 901 of file UnsaturatedFlow.hpp.

Constructor & Destructor Documentation

◆ OpIntegrateFluxes()

MixTransport::UnsaturatedFlowElement::OpIntegrateFluxes::OpIntegrateFluxes ( UnsaturatedFlowElement ctx,
const std::string  fluxes_name 
)

Constructor.

Examples
UnsaturatedFlow.hpp.

Definition at line 908 of file UnsaturatedFlow.hpp.

911  fluxes_name, UserDataOperator::OPROW),
912  cTx(ctx) {}
ForcesAndSourcesCore::UserDataOperator UserDataOperator

Member Function Documentation

◆ doWork()

MoFEMErrorCode MixTransport::UnsaturatedFlowElement::OpIntegrateFluxes::doWork ( int  side,
EntityType  type,
DataForcesAndSourcesCore::EntData data 
)

Integrate boundary condition.

Parameters
sidelocal index of entity
typetype of entity
datadata on entity
Returns
error code
Examples
UnsaturatedFlow.hpp.

Definition at line 923 of file UnsaturatedFlow.hpp.

924  {
926  int nb_dofs = data.getFieldData().size();
927  if (nb_dofs == 0)
929  // Get base function
930  auto t_n_hdiv = data.getFTensor1N<3>();
931  // get normal of face
932  auto t_normal = getFTensor1Normal();
933  // Integration weight
934  auto t_w = getFTensor0IntegrationWeight();
935  double flux_on_entity = 0;
936  int nb_gauss_pts = data.getN().size1();
937  for (int gg = 0; gg < nb_gauss_pts; gg++) {
938  auto t_data = data.getFTensor0FieldData();
939  for (int rr = 0; rr != nb_dofs; rr++) {
940  flux_on_entity -= (0.5 * t_data * t_w) * (t_n_hdiv(i) * t_normal(i));
941  ++t_n_hdiv;
942  ++t_data;
943  }
944  ++t_w;
945  }
946  CHKERR VecSetValue(cTx.ghostFlux, 0, flux_on_entity, ADD_VALUES);
948  }
Vec ghostFlux
Ghost Vector of integrated fluxes.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

Member Data Documentation

◆ cTx

UnsaturatedFlowElement& MixTransport::UnsaturatedFlowElement::OpIntegrateFluxes::cTx
Examples
UnsaturatedFlow.hpp.

Definition at line 903 of file UnsaturatedFlow.hpp.

◆ i

FTensor::Index<'i', 3> MixTransport::UnsaturatedFlowElement::OpIntegrateFluxes::i
Examples
UnsaturatedFlow.hpp.

Definition at line 914 of file UnsaturatedFlow.hpp.


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