1681 {
1683
1684 EntityHandle fe_ent = OP::getFEEntityHandle();
1686 if (bc.faces.find(fe_ent) != bc.faces.end()) {
1689 auto &v_analytical_expr = std::get<1>(analytical_data);
1690
1692 int nb_integration_pts = OP::getGaussPts().size2();
1693 auto t_w = OP::getFTensor0IntegrationWeight();
1694 int nb_base_functions = data.
getN().size2();
1696
1697#ifndef NDEBUG
1698 if (!this->sourceVec) {
1700 "Source vector for OpTauStabilizationOpAnalyticalDispBc is not "
1701 "set");
1702 }
1703 if (data.
getN().size1() != nb_integration_pts) {
1705 "Number of integration points in data should be %d but is %d",
1706 nb_integration_pts, (
int)data.
getN().size1());
1707 }
1708 if (nb_base_functions < nb_dofs /
SPACE_DIM) {
1710 "Number of base functions in data should be at least %d but is "
1711 "%d",
1712 nb_dofs /
SPACE_DIM, nb_base_functions);
1713 }
1714#endif
1715
1716 auto t_disp_val =
1718 *this->sourceVec, nb_integration_pts)();
1720
1722
1723 auto area = getMeasure();
1724 auto t_coords = OP::getFTensor1CoordsAtGaussPts();
1725 for (int gg = 0; gg != nb_integration_pts; ++gg) {
1726 auto tau_scale =
1727 area * t_w * OP::betaCoeff(t_coords(0), t_coords(1), t_coords(2));
1730 t_bc_residual(ii) =
1731 bc.flags[ii] ? t_disp_val(ii) - t_bc_disp(ii) : 0.;
1732 auto t_nf = getFTensor1FromPtr<3, 3>(OP::locF.data().data());
1733 int bb = 0;
1734 for (; bb != nb_dofs /
SPACE_DIM; ++bb) {
1736 (tau_scale * t_row_base_fun) * t_bc_residual(
i);
1737 ++t_nf;
1738 ++t_row_base_fun;
1739 }
1740 for (; bb != nb_base_functions; ++bb)
1741 ++t_row_base_fun;
1742
1743 ++t_w;
1744 ++t_coords;
1745 ++t_disp_val;
1746 ++t_bc_disp;
1747 }
1748 }
1749 }
1750
1752}
#define FTENSOR_INDEX(DIM, I)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_DATA_INCONSISTENCY
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
FTensor::Index< 'i', SPACE_DIM > i
std::tuple< std::string, MatrixDouble > getAnalyticalExpr(OP_PTR op_ptr, MatrixDouble &analytical_expr, const std::string block_name)
DataLayoutTraits< DataLayout::GaussByCoeffs > DL
auto getFTensor1FromMat(M &data, int rr=0, int cc=0)
Get tensor rank 1 (vector) form data matrix.
decltype(GetFTensor1FromMatImpl< Tensor_Dim, S, DL, M >::get(std::declval< M & >(), 0, 0)) GetFTensor1FromMatType
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
const VectorInt & getIndices() const
Get global indices of degrees of freedom on entity.