7#ifndef __SURFACE_SLIDING_CONSTRAINS_HPP__
8#define __SURFACE_SLIDING_CONSTRAINS_HPP__
25 boost::shared_ptr<VectorDouble> &active_variables_ptr)
30 EntitiesFieldData::EntData &data) {
32 if (type == MBVERTEX) {
33 for (
unsigned int dd = 0; dd != data.getFieldData().size(); ++dd)
40 template <
int SizeLambda>
46 boost::shared_ptr<VectorDouble> &active_variables_ptr)
51 EntitiesFieldData::EntData &data) {
53 if (type == MBVERTEX) {
54 for (
unsigned int dd = 0; dd != data.getFieldData().size(); ++dd)
61 template <
int SizeLambda,
int SizePositions>
67 boost::shared_ptr<VectorDouble> &results_ptr)
73 EntitiesFieldData::EntData &data) {
77 VectorInt &indices = data.getIndices();
81 else if (indices.size() == SizeLambda)
83 else if (indices.size() == SizePositions)
87 "Element %s: Data inconsistency nb of indices %d",
90 CHKERR VecSetOption(
getSNESf(), VEC_IGNORE_NEGATIVE_INDICES, PETSC_TRUE);
96 template <
int SizeLambda,
int SizePositions>
102 const std::string field_name_col,
103 boost::shared_ptr<MatrixDouble> &jacobian_ptr)
112 EntitiesFieldData::EntData &row_data,
113 EntitiesFieldData::EntData &col_data) {
115 if (row_type != MBVERTEX)
117 if (col_type != MBVERTEX)
119 VectorInt &row_indices = row_data.getIndices();
120 VectorInt &col_indices = col_data.getIndices();
121 if (row_indices.empty() || col_indices.empty())
125 if (row_indices.size() == SizeLambda)
127 else if (row_indices.size() == SizePositions)
128 shift_row = SizeLambda;
131 "Data inconsistency");
134 if (col_indices.size() == SizeLambda)
136 else if (col_indices.size() == SizePositions)
137 shift_col = SizeLambda;
140 "Data inconsistency");
142 MatrixDouble jac(row_indices.size(), col_indices.size());
143 for (
unsigned int rr = 0; rr != row_indices.size(); ++rr) {
144 for (
unsigned int cc = 0; cc != col_indices.size(); ++cc) {
145 jac(rr, cc) = (*jacobianPtr)(shift_row + rr, shift_col + cc);
146 if (jac(rr, cc) != jac(rr, cc))
148 "Jacobian assemble not a number jac(rr,cc) = %3.4f",
152 CHKERR MatSetValues(
getSNESB(), row_indices.size(), &row_indices[0],
153 col_indices.size(), &col_indices[0], &jac(0, 0),
338 virtual MoFEMErrorCode
370 if (
B != PETSC_NULL) {
374 if (
F != PETSC_NULL) {
394 ierr = PetscOptionsBegin(PETSC_COMM_WORLD,
"",
395 "Get surface sliding constrains element scaling",
398 CHKERR PetscOptionsScalar(
"-surface_sliding_alpha",
"scaling parameter",
"",
400 ierr = PetscOptionsEnd();
412 CHKERRABORT(PETSC_COMM_WORLD,
ierr);
430 boost::shared_ptr<VectorDouble> &active_variables_ptr,
431 boost::shared_ptr<VectorDouble> &results_ptr,
432 boost::shared_ptr<MatrixDouble> &jacobian_ptr,
443 EntitiesFieldData::EntData &data) {
445 if (type != MBVERTEX)
448 VectorInt &indices = data.getIndices();
452 ublas::vector<adouble> lambda_dofs(3);
453 for (
int dd = 0; dd != 3; ++dd) {
454 lambda_dofs[dd] <<= (*activeVariablesPtr)[dd];
456 ublas::vector<adouble> position_dofs(9);
457 for (
int dd = 0; dd != 9; ++dd) {
458 position_dofs[dd] <<= (*activeVariablesPtr)[3 + dd];
472 int nb_gauss_pts = data.getN().size1();
473 int nb_base_functions = data.getN().size2();
474 auto t_base1 = data.getFTensor0N();
475 auto t_base2 = data.getFTensor0N();
476 auto t_diff_base = data.getFTensor1DiffN<2>();
484 ublas::vector<adouble> c_vec(3);
485 ublas::vector<adouble> f_vec(9);
491 for (
int gg = 0; gg != nb_gauss_pts; ++gg) {
494 &position_dofs[0], &position_dofs[1], &position_dofs[2], 3);
498 t_position_ksi(
i) = 0;
499 t_position_eta(
i) = 0;
502 for (
int bb = 0; bb != nb_base_functions; ++bb) {
503 t_position(
i) += t_base1 * t_position_dofs(
i);
504 t_position_ksi(
i) += t_diff_base(
N0) * t_position_dofs(
i);
505 t_position_eta(
i) += t_diff_base(
N1) * t_position_dofs(
i);
506 lambda += t_base1 * t_lambda_dof;
513 t_delta(
i) = t_position(
i) - t_coord_ref(
i);
521 adouble area = sqrt(t_normal(
i) * t_normal(
i));
523 for (
int bb = 0; bb != nb_base_functions; ++bb) {
524 if (indices[bb] != -1) {
525 val = w * eo * t_base2;
526 t_c += val * t_normal(
i) * t_delta(
i);
528 t_f(
i) += val * t_normal(
i);
538 for (
int rr = 0; rr != 3; ++rr)
541 for (
int rr = 0; rr != 9; ++rr)
542 f_vec(rr) >>= (*resultsPtr)[3 + rr];
547 double *jac_ptr[3 + 9];
548 for (
int rr = 0; rr != 3 + 9; ++rr)
556 "ADOL-C function evaluation with error");
564 const std::string lagrange_multipliers_field_name,
565 const std::string material_field_name,
566 const double *alpha =
nullptr) {
569 if (alpha !=
nullptr) {
573 boost::shared_ptr<VectorDouble> active_variables_ptr(
574 new VectorDouble(3 + 9));
575 boost::shared_ptr<VectorDouble> results_ptr(
new VectorDouble(3 + 9));
576 boost::shared_ptr<MatrixDouble> jacobian_ptr(
577 new MatrixDouble(3 + 9, 3 + 9));
581 lagrange_multipliers_field_name, active_variables_ptr));
583 material_field_name, active_variables_ptr));
585 tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
595 lagrange_multipliers_field_name, active_variables_ptr));
597 material_field_name, active_variables_ptr));
599 tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
602 lagrange_multipliers_field_name, material_field_name, jacobian_ptr));
604 material_field_name, lagrange_multipliers_field_name, jacobian_ptr));
606 material_field_name, material_field_name, jacobian_ptr));
613 const std::string lagrange_multipliers_field_name,
614 const std::string material_field_name) {
617 boost::shared_ptr<VectorDouble> active_variables_ptr(
618 new VectorDouble(3 + 9));
619 boost::shared_ptr<VectorDouble> results_ptr(
new VectorDouble(3 + 9));
620 boost::shared_ptr<MatrixDouble> jacobian_ptr(
621 new MatrixDouble(3 + 9, 3 + 9));
626 lagrange_multipliers_field_name, active_variables_ptr));
628 material_field_name, active_variables_ptr));
630 tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
633 lagrange_multipliers_field_name, material_field_name, jacobian_ptr));
643 static MoFEMErrorCode
createTag(moab::Interface &moab, Tag &th0, Tag &th1,
644 Tag &th2, Tag &th3) {
646 double def_val[] = {0, 0, 0};
647 CHKERR moab.tag_get_handle(
"EDGE_BASE0", 3, MB_TYPE_DOUBLE, th0,
648 MB_TAG_CREAT | MB_TAG_SPARSE, def_val);
649 CHKERR moab.tag_get_handle(
"EDGE_BASE1", 3, MB_TYPE_DOUBLE, th1,
650 MB_TAG_CREAT | MB_TAG_SPARSE, def_val);
651 int def_numb_val[] = {-1};
652 CHKERR moab.tag_get_handle(
"PATCH_NUMBER", 1, MB_TYPE_INTEGER, th2,
653 MB_TAG_CREAT | MB_TAG_SPARSE, def_numb_val);
654 int def_orientation_val[] = {1};
655 CHKERR moab.tag_get_handle(
"PATCH_ORIENTATION", 1, MB_TYPE_INTEGER, th3,
656 MB_TAG_CREAT | MB_TAG_SPARSE,
657 def_orientation_val);
666 auto get_edges = [&](
const Range &ents) {
668 CHKERR moab.get_adjacencies(ents, 1,
false, edges,
669 moab::Interface::UNION);
673 auto get_face_adj = [&, get_edges](
const Range &faces) {
675 CHKERR moab.get_adjacencies(subtract(get_edges(faces), edges), 2,
false,
676 adj_faces, moab::Interface::UNION);
677 return intersect(adj_faces, tris);
680 auto get_patch = [&, get_face_adj](
const EntityHandle face) {
682 patch_ents.insert(face);
685 nb0 = patch_ents.size();
686 patch_ents.merge(get_face_adj(patch_ents));
687 }
while (nb0 != patch_ents.size());
691 auto get_patches = [&]() {
692 std::vector<Range> patches;
693 while (!tris.empty()) {
694 patches.push_back(get_patch(tris[0]));
695 tris = subtract(tris, patches.back());
700 Tag th0, th1, th2, th3;
703 auto patches = get_patches();
705 for (
auto patch : patches) {
708 std::vector<int> tags_vals(patch.size(), pp);
709 CHKERR moab.tag_set_data(th2, patch, &*tags_vals.begin());
717 moab::Interface &moab,
Range edges,
Range tris,
718 bool number_pathes =
true,
719 boost::shared_ptr<SurfaceSlidingConstrains::DriverElementOrientation>
720 surface_orientation =
nullptr,
723 Tag th0, th1, th2, th3;
728 for (
auto edge : edges) {
730 CHKERR moab.get_adjacencies(&edge, 1, 2,
false, adj_faces);
731 adj_faces = intersect(adj_faces, tris);
732 if (adj_faces.size() != 2) {
734 "Should be 2 faces adjacent to edge but is %d",
737 VectorDouble3
v[2] = {VectorDouble3(3), VectorDouble3(3)};
738 auto get_tensor_from_vec = [](VectorDouble3 &
v) {
740 &
v[0], &
v[1], &
v[2]);
747 std::array<double, 3> areas;
748 auto calculate_normals = [&, get_tensor_from_vec]() {
750 for (
auto face : adj_faces) {
751 double &x = (
v[ff])[0];
752 double &y = (
v[ff])[1];
753 double &z = (
v[ff])[2];
754 moab::Util::normal(&moab, face, x, y, z);
755 auto t_n = get_tensor_from_vec(
v[ff]);
756 areas[ff] = sqrt(t_n(
i) * t_n(
i));
760 CHKERR moab.tag_get_data(th3, &face, 1, &orientation);
761 if (orientation == -1) {
764 if (surface_orientation && m_field_ptr) {
765 CHKERR surface_orientation->getElementOrientation(*m_field_ptr,
767 int eo = surface_orientation->elementOrientation;
777 auto get_patch_number = [&]() {
778 std::vector<int>
p = {0, 0};
779 CHKERR moab.tag_get_data(th2, adj_faces, &*
p.begin());
783 std::array<EntityHandle, 2> faces_handles = {adj_faces[0],
785 auto order_normals = [&, get_patch_number]() {
786 auto p = get_patch_number();
789 faces_handles[0] = adj_faces[1];
790 faces_handles[1] = adj_faces[0];
796 auto t_n0 = get_tensor_from_vec(
v[0]);
797 auto t_n1 = get_tensor_from_vec(
v[1]);
800 auto get_base_for_coplanar_faces = [&]() {
803 std::vector<EntityHandle> face_conn;
804 CHKERR moab.get_connectivity(&faces_handles[1], 1, face_conn,
true);
805 std::vector<EntityHandle> edge_conn;
806 CHKERR moab.get_connectivity(&edge, 1, edge_conn,
true);
807 std::sort(face_conn.begin(), face_conn.end());
808 std::sort(edge_conn.begin(), edge_conn.end());
811 if (face_conn[
n] != edge_conn[
n])
813 VectorDouble3 coords_face(3);
814 CHKERR moab.get_coords(&face_conn[
n], 1,
815 &*coords_face.data().begin());
816 VectorDouble6 coords_edge(6);
817 CHKERR moab.get_coords(&edge_conn[0], 2,
818 &*coords_edge.data().begin());
820 coords_edge[0], coords_edge[1], coords_edge[2]);
822 coords_edge[3], coords_edge[4], coords_edge[5]);
823 auto t_face_n = get_tensor_from_vec(coords_face);
824 t_face_n(
i) -= t_edge_n0(
i);
826 t_delta(
i) = t_edge_n1(
i) - t_edge_n0(
i);
827 t_delta(
i) /= sqrt(t_delta(
i) * t_delta(
i));
829 if (t_n1(
i) * t_face_n(
i) < 0)
835 auto get_base_for_not_planar_faces = [&]() {
838 t_n1(
i) /= sqrt(t_n1(
i) * t_n1(
i));
843 constexpr double tol = 1e-6;
844 if ((t_cross(
i) * t_cross(
i)) <
tol)
845 CHKERR get_base_for_coplanar_faces();
847 CHKERR get_base_for_not_planar_faces();
849 VectorDouble3 &v0 =
v[0];
850 VectorDouble3 &v1 =
v[1];
851 CHKERR moab.tag_set_data(th0, &edge, 1, &v0[0]);
852 CHKERR moab.tag_set_data(th1, &edge, 1, &v1[0]);
857 static MoFEMErrorCode
saveEdges(moab::Interface &moab, std::string name,
863 CHKERR moab.create_meshset(MESHSET_SET, meshset);
864 CHKERR moab.add_entities(meshset, edges);
865 if (faces !=
nullptr) {
866 CHKERR moab.add_entities(meshset, *faces);
868 CHKERR moab.write_file(name.c_str(),
"VTK",
"", &meshset, 1, ths, 4);
869 CHKERR moab.delete_entities(&meshset, 1);
891 if (
B != PETSC_NULL) {
895 if (
F != PETSC_NULL) {
914 ierr = PetscOptionsBegin(PETSC_COMM_WORLD,
"",
915 "Get edge sliding constrains element scaling",
918 CHKERR PetscOptionsScalar(
"-edge_sliding_alpha",
"scaling parameter",
"",
920 ierr = PetscOptionsEnd();
930 CHKERRABORT(PETSC_COMM_WORLD,
ierr);
945 boost::shared_ptr<VectorDouble> &active_variables_ptr,
946 boost::shared_ptr<VectorDouble> &results_ptr,
947 boost::shared_ptr<MatrixDouble> &jacobian_ptr,
948 bool evaluate_jacobian,
const double &alpha)
956 EntitiesFieldData::EntData &data) {
958 if (type != MBVERTEX)
966 Tag th0, th1, th2, th3;
976 VectorInt &indices = data.getIndices();
980 ublas::vector<adouble> lambda_dofs(4);
981 for (
int dd = 0; dd != 4; ++dd) {
982 lambda_dofs[dd] <<= (*activeVariablesPtr)[dd];
984 ublas::vector<adouble> position_dofs(6);
985 for (
int dd = 0; dd != 6; ++dd) {
986 position_dofs[dd] <<= (*activeVariablesPtr)[4 + dd];
990 &position_dofs[0], &position_dofs[1], &position_dofs[2]);
992 &position_dofs[3], &position_dofs[4], &position_dofs[5]);
995 t_tangent(
i) = t_node1(
i) - t_node0(
i);
996 adouble l = sqrt(t_tangent(
i) * t_tangent(
i));
1000 t_dot0 = t_edge_base0(
i) * t_tangent(
i);
1001 t_dot1 = t_edge_base1(
i) * t_tangent(
i);
1004 t_base0(
i) = t_edge_base0(
i) - t_dot0 * t_tangent(
i);
1005 t_base1(
i) = t_edge_base1(
i) - t_dot1 * t_tangent(
i);
1006 t_base0(
i) /= sqrt(t_base0(
i) * t_base0(
i));
1007 t_base1(
i) /= sqrt(t_base1(
i) * t_base1(
i));
1009 auto t_base_fun1 = data.getFTensor0N();
1010 auto t_base_fun2 = data.getFTensor0N();
1016 ublas::vector<adouble> c_vec(4);
1017 ublas::vector<adouble> f_vec(6);
1021 int nb_gauss_pts = data.getN().size1();
1022 int nb_base_functions = data.getN().size2();
1023 for (
int gg = 0; gg != nb_gauss_pts; ++gg) {
1026 &position_dofs[0], &position_dofs[1], &position_dofs[2], 3);
1028 &lambda_dofs[1], 2);
1032 for (
int bb = 0; bb != nb_base_functions; ++bb) {
1033 t_position(
i) += t_base_fun1 * t_position_dofs(
i);
1034 t_lambda(
j) += t_base_fun1 * t_lambda_dof(
j);
1040 t_delta(
i) = t_position(
i) - t_coord_ref(
i);
1048 for (
int bb = 0; bb != nb_base_functions; ++bb) {
1049 if (indices[2 * bb] != -1) {
1050 val = w * t_base_fun2;
1051 t_c(
N0) += val * dot0;
1052 t_c(
N1) += val * dot1;
1053 val1 = val * t_lambda(
N0);
1054 val2 = val * t_lambda(
N1);
1055 t_f(
i) += val1 * t_base0(
i) + val2 * t_base1(
i);
1065 for (
int rr = 0; rr != 4; ++rr) {
1066 c_vec[rr] >>= (*resultsPtr)[rr];
1068 for (
int rr = 0; rr != 6; ++rr) {
1069 f_vec(rr) >>= (*resultsPtr)[4 + rr];
1075 double *jac_ptr[4 + 6];
1076 for (
int rr = 0; rr != 4 + 6; ++rr) {
1077 jac_ptr[rr] = &(*jacobianPtr)(rr, 0);
1084 "ADOL-C function evaluation with error");
1093 const std::string lagrange_multipliers_field_name,
1094 const std::string material_field_name) {
1099 material_field_name);
1104 const std::string lagrange_multipliers_field_name,
1105 const std::string material_field_name,
1106 const double *alpha =
nullptr) {
1113 boost::shared_ptr<VectorDouble> active_variables_ptr(
1114 new VectorDouble(4 + 6));
1115 boost::shared_ptr<VectorDouble> results_ptr(
new VectorDouble(4 + 6));
1116 boost::shared_ptr<MatrixDouble> jacobian_ptr(
1117 new MatrixDouble(4 + 6, 4 + 6));
1121 lagrange_multipliers_field_name, active_variables_ptr));
1123 material_field_name, active_variables_ptr));
1125 tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
1126 jacobian_ptr,
false,
aLpha));
1135 lagrange_multipliers_field_name, active_variables_ptr));
1137 material_field_name, active_variables_ptr));
1139 tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
1140 jacobian_ptr,
true,
aLpha));
1142 lagrange_multipliers_field_name, material_field_name, jacobian_ptr));
1144 material_field_name, lagrange_multipliers_field_name, jacobian_ptr));
1146 material_field_name, material_field_name, jacobian_ptr));
1153 const std::string lagrange_multipliers_field_name,
1154 const std::string material_field_name) {
1160 material_field_name);
1166 const std::string lagrange_multipliers_field_name,
1167 const std::string material_field_name) {
1170 boost::shared_ptr<VectorDouble> active_variables_ptr(
1171 new VectorDouble(4 + 6));
1172 boost::shared_ptr<VectorDouble> results_ptr(
new VectorDouble(4 + 6));
1173 boost::shared_ptr<MatrixDouble> jacobian_ptr(
1174 new MatrixDouble(4 + 6, 4 + 6));
1179 lagrange_multipliers_field_name, active_variables_ptr));
1181 material_field_name, active_variables_ptr));
1183 tag, lagrange_multipliers_field_name, active_variables_ptr, results_ptr,
1184 jacobian_ptr,
true,
aLpha));
1186 lagrange_multipliers_field_name, material_field_name, jacobian_ptr));
ForcesAndSourcesCore::UserDataOperator UserDataOperator
static PetscErrorCode ierr
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_OPERATION_UNSUCCESSFUL
@ MOFEM_DATA_INCONSISTENCY
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
FTensor::Index< 'n', SPACE_DIM > n
constexpr double lambda
surface tension
FTensor::Index< 'i', SPACE_DIM > i
const double v
phase velocity of light in medium (cm/ns)
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
auto cross(const Tensor1_Expr< A, T, 3, i > &a, const Tensor1_Expr< B, U, 3, j > &b, const Index< k, 3 > &)
implementation of Data Operators for Forces and Sources
constexpr auto field_name
static MoFEMErrorCode numberSurfaces(moab::Interface &moab, Range edges, Range tris)
static MoFEMErrorCode createTag(moab::Interface &moab, Tag &th0, Tag &th1, Tag &th2, Tag &th3)
static MoFEMErrorCode setTags(moab::Interface &moab, Range edges, Range tris, bool number_pathes=true, boost::shared_ptr< SurfaceSlidingConstrains::DriverElementOrientation > surface_orientation=nullptr, MoFEM::Interface *m_field_ptr=nullptr)
static MoFEMErrorCode saveEdges(moab::Interface &moab, std::string name, Range edges, Range *faces=nullptr)
MoFEMErrorCode preProcess()
function is run at the beginning of loop
MyEdgeFE(MoFEM::Interface &m_field)
boost::shared_ptr< VectorDouble > activeVariablesPtr
boost::shared_ptr< MatrixDouble > jacobianPtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
boost::shared_ptr< VectorDouble > resultsPtr
OpJacobian(int tag, const std::string field_name, boost::shared_ptr< VectorDouble > &active_variables_ptr, boost::shared_ptr< VectorDouble > &results_ptr, boost::shared_ptr< MatrixDouble > &jacobian_ptr, bool evaluate_jacobian, const double &alpha)
MoFEMErrorCode setOperators(int tag, Range edges, Range faces, const std::string lagrange_multipliers_field_name, const std::string material_field_name)
boost::shared_ptr< MyEdgeFE > feRhsPtr
EdgeSlidingConstrains(MoFEM::Interface &m_field)
boost::shared_ptr< MyEdgeFE > feLhsPtr
MyEdgeFE & getLoopFeRhs()
MyEdgeFE & getLoopFeLhs()
MoFEM::Interface & mField
MoFEMErrorCode setOperators(int tag, const std::string lagrange_multipliers_field_name, const std::string material_field_name, const double *alpha=nullptr)
MoFEMErrorCode getOptions()
MoFEMErrorCode setOperatorsConstrainOnly(int tag, Range edges, Range faces, const std::string lagrange_multipliers_field_name, const std::string material_field_name)
MoFEMErrorCode setOperatorsConstrainOnly(int tag, const std::string lagrange_multipliers_field_name, const std::string material_field_name)
OpAssembleLhs(const std::string field_name_row, const std::string field_name_col, boost::shared_ptr< MatrixDouble > &jacobian_ptr)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
boost::shared_ptr< MatrixDouble > jacobianPtr
OpAssembleRhs(const std::string field_name, boost::shared_ptr< VectorDouble > &results_ptr)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
boost::shared_ptr< VectorDouble > resultsPtr
boost::shared_ptr< VectorDouble > activeVariablesPtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
OpGetActiveDofsLambda(const std::string field_name, boost::shared_ptr< VectorDouble > &active_variables_ptr)
OpGetActiveDofsPositions(const std::string field_name, boost::shared_ptr< VectorDouble > &active_variables_ptr)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
boost::shared_ptr< VectorDouble > activeVariablesPtr
Implementation of surface sliding constrains.
virtual ~GenericSliding()=default
virtual moab::Interface & get_moab()=0
bool sYmm
If true assume that matrix is symmetric structure.
Deprecated interface functions.
default operator for EDGE element
const EdgeElementForcesAndSourcesCore * getEdgeFE()
get pointer to this finite element
default operator for TRI element
FaceElementForcesAndSourcesCore * getFaceFE()
return pointer to Generic Triangle Finite Element object
EntityHandle getFEEntityHandle() const
Return finite element entity handle.
auto getFTensor1CoordsAtGaussPts()
Get coordinates at integration points assuming linear geometry.
friend class ForcesAndSourcesCore
std::string getFEName() const
Get name of the element.
@ OPCOL
operator doWork function is executed on FE columns
@ OPROW
operator doWork function is executed on FE rows
@ OPROWCOL
operator doWork is executed on FE rows &columns
const FEMethod * getFEMethod() const
Return raw pointer to Finite Element Method object.
MatrixDouble & getGaussPts()
matrix of integration (Gauss) points for Volume Element
virtual MoFEMErrorCode preProcess()
function is run at the beginning of loop
boost::ptr_deque< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.
Mat & snes_B
preconditioner of jacobian matrix
Class implemented by user to detect face orientation.
virtual MoFEMErrorCode getElementOrientation(MoFEM::Interface &m_field, const EntityHandle face)
virtual MoFEMErrorCode getElementOrientation(MoFEM::Interface &m_field, const FEMethod *fe_method_ptr)
MoFEMErrorCode preProcess()
function is run at the beginning of loop
MyTriangleFE(MoFEM::Interface &m_field)
OpJacobian(int tag, const std::string field_name, boost::shared_ptr< VectorDouble > &active_variables_ptr, boost::shared_ptr< VectorDouble > &results_ptr, boost::shared_ptr< MatrixDouble > &jacobian_ptr, DriverElementOrientation &orientation, bool evaluate_jacobian, const double &alpha)
boost::shared_ptr< VectorDouble > activeVariablesPtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
boost::shared_ptr< VectorDouble > resultsPtr
DriverElementOrientation & oRientation
boost::shared_ptr< MatrixDouble > jacobianPtr
Shape preserving constrains, i.e. nodes sliding on body surface.
virtual ~SurfaceSlidingConstrains()=default
boost::shared_ptr< MyTriangleFE > feLhsPtr
MyTriangleFE & getLoopFeLhs()
SurfaceSlidingConstrains(MoFEM::Interface &m_field, DriverElementOrientation &orientation)
MoFEMErrorCode getOptions()
MoFEMErrorCode setOperators(int tag, const std::string lagrange_multipliers_field_name, const std::string material_field_name, const double *alpha=nullptr)
MyTriangleFE & getLoopFeRhs()
MoFEMErrorCode setOperatorsConstrainOnly(int tag, const std::string lagrange_multipliers_field_name, const std::string material_field_name)
MoFEM::Interface & mField
boost::shared_ptr< MyTriangleFE > feRhsPtr
DriverElementOrientation & crackFrontOrientation