v0.10.0
Public Member Functions | Public Attributes | List of all members
MoFEM::OpSetCovariantPiolaTransformOnFace Struct Reference

transform Hcurl base fluxes from reference element to physical triangle More...

#include <src/finite_elements/DataOperators.hpp>

Inheritance diagram for MoFEM::OpSetCovariantPiolaTransformOnFace:
[legend]
Collaboration diagram for MoFEM::OpSetCovariantPiolaTransformOnFace:
[legend]

Public Member Functions

 OpSetCovariantPiolaTransformOnFace (const VectorDouble &normal, const MatrixDouble &normals_at_pts, const VectorDouble &tangent0, const MatrixDouble &tangent0_at_pts, const VectorDouble &tangent1, const MatrixDouble &tangent1_at_pts)
 
MoFEMErrorCode doWork (int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
 Operator for linear form, usually to calculate values on right hand side. More...
 
- Public Member Functions inherited from MoFEM::DataOperator
 DataOperator (const bool symm=true)
 
virtual ~DataOperator ()=default
 
virtual MoFEMErrorCode doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
 Operator for bi-linear form, usually to calculate values on left hand side. More...
 
virtual MoFEMErrorCode opLhs (DataForcesAndSourcesCore &row_data, DataForcesAndSourcesCore &col_data)
 
virtual MoFEMErrorCode opRhs (DataForcesAndSourcesCore &data, const bool error_if_no_base=false)
 
bool getSymm () const
 Get if operator uses symmetry of DOFs or not. More...
 
void setSymm ()
 set if operator is executed taking in account symmetry More...
 
void unSetSymm ()
 unset if operator is executed for non symmetric problem More...
 

Public Attributes

const VectorDoublenOrmal
 
const MatrixDoublenormalsAtGaussPts
 
const VectorDoubletAngent0
 
const MatrixDoubletangent0AtGaussPt
 
const VectorDoubletAngent1
 
const MatrixDoubletangent1AtGaussPt
 
- Public Attributes inherited from MoFEM::DataOperator
bool sYmm
 If true assume that matrix is symmetric structure. More...
 
std::array< bool, MBMAXTYPE > doEntities
 If true operator is executed for entity. More...
 
booldoVertices
 \deprectaed If false skip vertices More...
 
booldoEdges
 \deprectaed If false skip edges More...
 
booldoQuads
 \deprectaed More...
 
booldoTris
 \deprectaed More...
 
booldoTets
 \deprectaed More...
 
booldoPrisms
 \deprectaed More...
 

Detailed Description

transform Hcurl base fluxes from reference element to physical triangle

Definition at line 544 of file DataOperators.hpp.

Constructor & Destructor Documentation

◆ OpSetCovariantPiolaTransformOnFace()

MoFEM::OpSetCovariantPiolaTransformOnFace::OpSetCovariantPiolaTransformOnFace ( const VectorDouble normal,
const MatrixDouble normals_at_pts,
const VectorDouble tangent0,
const MatrixDouble tangent0_at_pts,
const VectorDouble tangent1,
const MatrixDouble tangent1_at_pts 
)

Definition at line 553 of file DataOperators.hpp.

559  : nOrmal(normal), normalsAtGaussPts(normals_at_pts), tAngent0(tangent0),
560  tangent0AtGaussPt(tangent0_at_pts), tAngent1(tangent1),
561  tangent1AtGaussPt(tangent1_at_pts) {}
FTensor::Tensor1< double, 2 > normal(FTensor::Tensor1< T1, 3 > &t_coords, FTensor::Tensor1< T2, 2 > &t_disp)
Definition: ContactOps.hpp:137

Member Function Documentation

◆ doWork()

MoFEMErrorCode MoFEM::OpSetCovariantPiolaTransformOnFace::doWork ( int  side,
EntityType  type,
DataForcesAndSourcesCore::EntData data 
)
virtual

Operator for linear form, usually to calculate values on right hand side.

Reimplemented from MoFEM::DataOperator.

Definition at line 988 of file DataOperators.cpp.

989  {
991 
992  if (type != MBEDGE && type != MBTRI)
994 
998 
1000  &tAngent0[0], &tAngent1[0], &nOrmal[0],
1001 
1002  &tAngent0[1], &tAngent1[1], &nOrmal[1],
1003 
1004  &tAngent0[2], &tAngent1[2], &nOrmal[2],
1005 
1006  3);
1007  double det;
1009  CHKERR determinantTensor3by3(t_m, det);
1010  CHKERR invertTensor3by3(t_m, det, t_inv_m);
1011 
1012  for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; ++b) {
1013 
1014  FieldApproximationBase base = static_cast<FieldApproximationBase>(b);
1015 
1016  auto &baseN = data.getN(base);
1017  auto &diffBaseN = data.getDiffN(base);
1018 
1019  int nb_dofs = baseN.size2() / 3;
1020  int nb_gauss_pts = baseN.size1();
1021 
1022  MatrixDouble piola_n(baseN.size1(), baseN.size2());
1023  MatrixDouble diff_piola_n(diffBaseN.size1(), diffBaseN.size2());
1024 
1025  if (nb_dofs > 0 && nb_gauss_pts > 0) {
1026 
1028  &baseN(0, HVEC0), &baseN(0, HVEC1), &baseN(0, HVEC2));
1030  &diffBaseN(0, HVEC0_0), &diffBaseN(0, HVEC0_1),
1031  &diffBaseN(0, HVEC1_0), &diffBaseN(0, HVEC1_1),
1032  &diffBaseN(0, HVEC2_0), &diffBaseN(0, HVEC2_1));
1033  FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 3> t_transformed_h_curl(
1034  &piola_n(0, HVEC0), &piola_n(0, HVEC1), &piola_n(0, HVEC2));
1036  t_transformed_diff_h_curl(
1037  &diff_piola_n(0, HVEC0_0), &diff_piola_n(0, HVEC0_1),
1038  &diff_piola_n(0, HVEC1_0), &diff_piola_n(0, HVEC1_1),
1039  &diff_piola_n(0, HVEC2_0), &diff_piola_n(0, HVEC2_1));
1040 
1041  int cc = 0;
1042  if (normalsAtGaussPts.size1() == (unsigned int)nb_gauss_pts) {
1043  // HO geometry is set, so jacobian is different at each gauss point
1045  &tangent0AtGaussPt(0, 0), &tangent1AtGaussPt(0, 0),
1046  &normalsAtGaussPts(0, 0), &tangent0AtGaussPt(0, 1),
1047  &tangent1AtGaussPt(0, 1), &normalsAtGaussPts(0, 1),
1048  &tangent0AtGaussPt(0, 2), &tangent1AtGaussPt(0, 2),
1049  &normalsAtGaussPts(0, 2), 3);
1050  for (int gg = 0; gg < nb_gauss_pts; ++gg) {
1051  CHKERR determinantTensor3by3(t_m_at_pts, det);
1052  CHKERR invertTensor3by3(t_m_at_pts, det, t_inv_m);
1053  for (int ll = 0; ll != nb_dofs; ll++) {
1054  t_transformed_h_curl(i) = t_inv_m(j, i) * t_h_curl(j);
1055  t_transformed_diff_h_curl(i, k) =
1056  t_inv_m(j, i) * t_diff_h_curl(j, k);
1057  ++t_h_curl;
1058  ++t_transformed_h_curl;
1059  ++t_diff_h_curl;
1060  ++t_transformed_diff_h_curl;
1061  ++cc;
1062  }
1063  ++t_m_at_pts;
1064  }
1065  } else {
1066  for (int gg = 0; gg < nb_gauss_pts; ++gg) {
1067  for (int ll = 0; ll != nb_dofs; ll++) {
1068  t_transformed_h_curl(i) = t_inv_m(j, i) * t_h_curl(j);
1069  t_transformed_diff_h_curl(i, k) =
1070  t_inv_m(j, i) * t_diff_h_curl(j, k);
1071  ++t_h_curl;
1072  ++t_transformed_h_curl;
1073  ++t_diff_h_curl;
1074  ++t_transformed_diff_h_curl;
1075  ++cc;
1076  }
1077  }
1078  }
1079  if (cc != nb_gauss_pts * nb_dofs)
1080  SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSIBLE_CASE, "Data inconsistency");
1081 
1082  baseN.data().swap(piola_n.data());
1083  diffBaseN.data().swap(diff_piola_n.data());
1084  }
1085  }
1086 
1088 }
MoFEMErrorCode determinantTensor3by3(T1 &t, T2 &det)
Calculate determinant 3 by 3.
Definition: Templates.hpp:571
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:76
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
MoFEMErrorCode invertTensor3by3(ublas::matrix< T, L, A > &jac_data, ublas::vector< T, A > &det_data, ublas::matrix< T, L, A > &inv_jac_data)
Calculate inverse of tensor rank 2 at integration points.
Definition: Templates.hpp:552
static Index< 'i', 3 > i
FieldApproximationBase
approximation base
Definition: definitions.h:150
static Index< 'j', 3 > j
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:152
#define CHKERR
Inline error check.
Definition: definitions.h:604
static Index< 'k', 3 > k
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415

Member Data Documentation

◆ nOrmal

const VectorDouble& MoFEM::OpSetCovariantPiolaTransformOnFace::nOrmal

Definition at line 546 of file DataOperators.hpp.

◆ normalsAtGaussPts

const MatrixDouble& MoFEM::OpSetCovariantPiolaTransformOnFace::normalsAtGaussPts

Definition at line 547 of file DataOperators.hpp.

◆ tAngent0

const VectorDouble& MoFEM::OpSetCovariantPiolaTransformOnFace::tAngent0

Definition at line 548 of file DataOperators.hpp.

◆ tangent0AtGaussPt

const MatrixDouble& MoFEM::OpSetCovariantPiolaTransformOnFace::tangent0AtGaussPt

Definition at line 549 of file DataOperators.hpp.

◆ tAngent1

const VectorDouble& MoFEM::OpSetCovariantPiolaTransformOnFace::tAngent1

Definition at line 550 of file DataOperators.hpp.

◆ tangent1AtGaussPt

const MatrixDouble& MoFEM::OpSetCovariantPiolaTransformOnFace::tangent1AtGaussPt

Definition at line 551 of file DataOperators.hpp.


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