v0.9.0
Public Member Functions | Public Attributes | List of all members
CellEngineering::OpVirtualPotentialRho Struct Reference

Post-process tractions. More...

#include <users_modules/cell_engineering/src/CellForces.hpp>

Inheritance diagram for CellEngineering::OpVirtualPotentialRho:
[legend]
Collaboration diagram for CellEngineering::OpVirtualPotentialRho:
[legend]

Public Member Functions

 OpVirtualPotentialRho (std::string field_name, CommonData &common_data)
 
MoFEMErrorCode doWork (int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
 

Public Attributes

CommonDatacommonData
 

Detailed Description

Post-process tractions.

Todo:
Should be done with Prisms, no need to have three kinds of elements, tets for solid, prisms for constrains and triangles to do post-processing. We can do post-pocessing with prisms.
Examples
cell_forces.cpp.

Definition at line 460 of file CellForces.hpp.

Constructor & Destructor Documentation

◆ OpVirtualPotentialRho()

CellEngineering::OpVirtualPotentialRho::OpVirtualPotentialRho ( std::string  field_name,
CommonData common_data 
)

Definition at line 465 of file CellForces.hpp.

467  field_name, UserDataOperator::OPROW),
468  commonData(common_data) {}
ForcesAndSourcesCore::UserDataOperator UserDataOperator

Member Function Documentation

◆ doWork()

MoFEMErrorCode CellEngineering::OpVirtualPotentialRho::doWork ( int  side,
EntityType  type,
DataForcesAndSourcesCore::EntData data 
)

Real work is done here, i.e. assemble tractions

Definition at line 473 of file CellForces.hpp.

474  {
475  //
476  //
478 
479  if (data.getFieldData().empty())
481 
482  const int nb_gauss_pts = data.getN().size1();
483  const int nb_dofs = data.getFieldData().size();
484 
485  VectorDouble *vec;
486  MatrixDouble *mat;
487  MatrixDouble *grad;
488 
489  if (rowFieldName == "RHO") {
493  }
494 
495  if (type == MBVERTEX) {
496  mat->resize(nb_gauss_pts, 2, false);
497  mat->clear();
498  vec->resize(nb_gauss_pts, false);
499  vec->clear();
500  grad->resize(0, 0, false);
501  }
502 
503  // calculate forces at integration points
504  for (unsigned int gg = 0; gg < nb_gauss_pts; gg++) {
505  for (int rr = 0; rr != nb_dofs; rr++) {
506  const double diff_base_x = data.getDiffN()(gg, 2 * rr + 0);
507  const double diff_base_y = data.getDiffN()(gg, 2 * rr + 1);
508  const double val = data.getFieldData()[rr];
509  (*mat)(gg, 0) += diff_base_x * val;
510  (*mat)(gg, 1) += diff_base_y * val;
511  const double base = data.getN()(gg, rr);
512  (*vec)[gg] += base * val;
513  }
514  }
515 
517  }
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:74
#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
VectorDouble forcePotentialAtGaussPoints
Definition: CellForces.hpp:34
MatrixDouble gradForceAtGaussPtrs
Definition: CellForces.hpp:35
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:72

Member Data Documentation

◆ commonData

CommonData& CellEngineering::OpVirtualPotentialRho::commonData

Definition at line 463 of file CellForces.hpp.


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