Assemble P *. More...

#include "tutorials/cor-7/src/ElasticityMixedFormulation.hpp"

Inheritance diagram for OpAssembleP:

Collaboration diagram for OpAssembleP:

Public Member Functions
	OpAssembleP (DataAtIntegrationPts &common_data, BlockData &data)

PetscErrorCode	doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)

Public Attributes
DataAtIntegrationPts &	commonData

MatrixDouble	locP

MatrixDouble	translocP

BlockData &	dAta

Detailed Description

Assemble P *.

*

\[ {\bf{P}} = - \int\limits_\Omega {{\bf{N}}_p^T\frac{1}{\lambda }{{\bf{N}}_p}d\Omega } \]

Examples: mofem/tutorials/cor-7/elasticity_mixed_formulation.cpp, and mofem/tutorials/cor-7/src/ElasticityMixedFormulation.hpp.

Definition at line 113 of file ElasticityMixedFormulation.hpp.

Constructor & Destructor Documentation

◆ OpAssembleP()

OpAssembleP::OpAssembleP	(	DataAtIntegrationPts &	common_data,
		BlockData &	data
	)

inline

Definition at line 121 of file ElasticityMixedFormulation.hpp.

      : VolumeElementForcesAndSourcesCore::UserDataOperator(
            "P", "P", UserDataOperator::OPROWCOL),
        commonData(common_data), dAta(data) {
    sYmm = true;
  }

Member Function Documentation

◆ doWork()

PetscErrorCode OpAssembleP::doWork	(	int	row_side,
		int	col_side,
		EntityType	row_type,
		EntityType	col_type,
		EntitiesFieldData::EntData &	row_data,
		EntitiesFieldData::EntData &	col_data
	)

inline

Examples: mofem/tutorials/cor-7/src/ElasticityMixedFormulation.hpp.

Definition at line 128 of file ElasticityMixedFormulation.hpp.

                                                            {
 
    MoFEMFunctionBegin;
    const int row_nb_dofs = row_data.getIndices().size();
    if (!row_nb_dofs)
      MoFEMFunctionReturnHot(0);
    const int col_nb_dofs = col_data.getIndices().size();
    if (!col_nb_dofs)
      MoFEMFunctionReturnHot(0);
 
    if (dAta.tEts.find(getNumeredEntFiniteElementPtr()->getEnt()) ==
        dAta.tEts.end()) {
      MoFEMFunctionReturnHot(0);
    }
    CHKERR commonData.getBlockData(dAta);
    // Set size can clear local tangent matrix
    locP.resize(row_nb_dofs, col_nb_dofs, false);
    locP.clear();
 
    const int row_nb_gauss_pts = row_data.getN().size1();
    if (!row_nb_gauss_pts)
      MoFEMFunctionReturnHot(0);
    const int row_nb_base_functions = row_data.getN().size2();
    auto row_base_functions = row_data.getFTensor0N();
 
    // get data
    FTensor::Index<'i', 3> i;
    FTensor::Index<'j', 3> j;
    const double lambda = commonData.lAmbda;
 
    double coefficient = commonData.pOisson == 0.5 ? 0. : 1 / lambda;
 
    // integration
    if (coefficient != 0.) {
      for (int gg = 0; gg != row_nb_gauss_pts; gg++) {
 
        // Get volume and integration weight
        double w = getVolume() * getGaussPts()(3, gg);
 
        // INTEGRATION
        int row_bb = 0;
        for (; row_bb != row_nb_dofs; row_bb++) {
          auto col_base_functions = col_data.getFTensor0N(gg, 0);
          for (int col_bb = 0; col_bb != col_nb_dofs; col_bb++) {
 
            locP(row_bb, col_bb) -=
                w * row_base_functions * col_base_functions * coefficient;
 
            ++col_base_functions;
          }
          ++row_base_functions;
        }
        for (; row_bb != row_nb_base_functions; row_bb++) {
          ++row_base_functions;
        }
      }
    }
    
 
    CHKERR MatSetValues(
        getFEMethod()->ksp_B, row_nb_dofs, &*row_data.getIndices().begin(),
        col_nb_dofs, &*col_data.getIndices().begin(), &locP(0, 0), ADD_VALUES);
 
    // is symmetric
    if (row_side != col_side || row_type != col_type) {
      translocP.resize(col_nb_dofs, row_nb_dofs, false);
      noalias(translocP) = trans(locP);
      CHKERR MatSetValues(getFEMethod()->ksp_B, col_nb_dofs,
                          &*col_data.getIndices().begin(), row_nb_dofs,
                          &*row_data.getIndices().begin(), &translocP(0, 0),
                          ADD_VALUES);
    }
 
    MoFEMFunctionReturn(0);
  }