1416 {
1418
1421 if (bc.faces.find(fe_ent) != bc.faces.end()) {
1424 auto &v_analytical_expr = std::get<1>(analytical_data);
1425
1427 int nb_integration_pts = OP::getGaussPts().size2();
1428 auto t_w = OP::getFTensor0IntegrationWeight();
1429 int nb_base_functions = data.
getN().size2();
1431
1432#ifndef NDEBUG
1433 if (!this->sourceVec) {
1435 "Source vector for OpTauStabilizationOpAnalyticalDispBc is not "
1436 "set");
1437 }
1438 if (data.
getN().size1() != nb_integration_pts) {
1440 "Number of integration points in data should be %d but is %d",
1441 nb_integration_pts, (
int)data.
getN().size1());
1442 }
1443 if (nb_base_functions < nb_dofs /
SPACE_DIM) {
1445 "Number of base functions in data should be at least %d but is "
1446 "%d",
1447 nb_dofs /
SPACE_DIM, nb_base_functions);
1448 }
1449#endif
1450
1451 auto t_disp_val =
1453 *this->sourceVec, nb_integration_pts)();
1455
1456 auto area = getMeasure();
1457 auto t_coords = OP::getFTensor1CoordsAtGaussPts();
1458 for (int gg = 0; gg != nb_integration_pts; ++gg) {
1459 auto tau_scale =
1460 area * t_w * OP::betaCoeff(t_coords(0), t_coords(1), t_coords(2));
1461 auto t_nf = getFTensor1FromPtr<3, 3>(OP::locF.data().data());
1462 int bb = 0;
1463 for (; bb != nb_dofs /
SPACE_DIM; ++bb) {
1464 for (
auto ii = 0; ii !=
SPACE_DIM; ++ii) {
1465 if (bc.flags[ii]) {
1466 t_nf(ii) +=
1467 (tau_scale * t_row_base_fun) * (t_disp_val(ii) - t_bc_disp(ii));
1468 }
1469 }
1470 ++t_nf;
1471 ++t_row_base_fun;
1472 }
1473 for (; bb != nb_base_functions; ++bb)
1474 ++t_row_base_fun;
1475
1476 ++t_w;
1477 ++t_coords;
1478 ++t_disp_val;
1479 ++t_bc_disp;
1480 }
1481 }
1482 }
1483
1485}
#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()
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.