v0.8.23
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, const bool do_vertices=true, const bool do_edges=true, const bool do_quads=true, const bool do_tris=true, const bool do_tets=true, const bool do_prisms=true)
 
virtual ~DataOperator ()
 
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, bool symm=true)
 
virtual MoFEMErrorCode opLhs (DataForcesAndSourcesCore &row_data, DataForcesAndSourcesCore &col_data)
 
virtual MoFEMErrorCode opRhs (DataForcesAndSourcesCore &data, const bool do_vertices, const bool do_edges, const bool do_quads, const bool do_tris, const bool do_tets, const bool do_prisms, const bool error_if_no_base=true)
 
virtual MoFEMErrorCode opRhs (DataForcesAndSourcesCore &data, const bool error_if_no_base=true)
 
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...
 
bool doVertices
 If false skip vertices. More...
 
bool doEdges
 If false skip edges. More...
 
bool doQuads
 
bool doTris
 
bool doTets
 
bool doPrisms
 

Detailed Description

transform Hcurl base fluxes from reference element to physical triangle

Definition at line 527 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 536 of file DataOperators.hpp.

542  : nOrmal(normal), normalsAtGaussPts(normals_at_pts), tAngent0(tangent0),
543  tangent0AtGaussPt(tangent0_at_pts), tAngent1(tangent1),
544  tangent1AtGaussPt(tangent1_at_pts) {}

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 1408 of file DataOperators.cpp.

1409  {
1411 
1412  if (type != MBEDGE && type != MBTRI)
1414 
1418 
1420  &tAngent0[0], &tAngent1[0], &nOrmal[0],
1421 
1422  &tAngent0[1], &tAngent1[1], &nOrmal[1],
1423 
1424  &tAngent0[2], &tAngent1[2], &nOrmal[2],
1425 
1426  3);
1427  double det;
1429  CHKERR determinantTensor3by3(t_m, det);
1430  CHKERR invertTensor3by3(t_m, det, t_inv_m);
1431 
1432  for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; ++b) {
1433 
1434  FieldApproximationBase base = static_cast<FieldApproximationBase>(b);
1435 
1436  auto &baseN = data.getN(base);
1437  auto &diffBaseN = data.getDiffN(base);
1438 
1439  int nb_dofs = baseN.size2() / 3;
1440  int nb_gauss_pts = baseN.size1();
1441 
1442  MatrixDouble piola_n(baseN.size1(), baseN.size2());
1443  MatrixDouble diff_piola_n(diffBaseN.size1(), diffBaseN.size2());
1444 
1445  if (nb_dofs > 0 && nb_gauss_pts > 0) {
1446 
1448  &baseN(0, HVEC0), &baseN(0, HVEC1), &baseN(0, HVEC2));
1450  &diffBaseN(0, HVEC0_0), &diffBaseN(0, HVEC0_1),
1451  &diffBaseN(0, HVEC1_0), &diffBaseN(0, HVEC1_1),
1452  &diffBaseN(0, HVEC2_0), &diffBaseN(0, HVEC2_1));
1453  FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 3> t_transformed_h_curl(
1454  &piola_n(0, HVEC0), &piola_n(0, HVEC1), &piola_n(0, HVEC2));
1456  t_transformed_diff_h_curl(
1457  &diff_piola_n(0, HVEC0_0), &diff_piola_n(0, HVEC0_1),
1458  &diff_piola_n(0, HVEC1_0), &diff_piola_n(0, HVEC1_1),
1459  &diff_piola_n(0, HVEC2_0), &diff_piola_n(0, HVEC2_1));
1460 
1461  int cc = 0;
1462  if (normalsAtGaussPts.size1() == (unsigned int)nb_gauss_pts) {
1463  // HO geometry is set, so jacobian is different at each gauss point
1465  &tangent0AtGaussPt(0, 0), &tangent1AtGaussPt(0, 0),
1466  &normalsAtGaussPts(0, 0), &tangent0AtGaussPt(0, 1),
1467  &tangent1AtGaussPt(0, 1), &normalsAtGaussPts(0, 1),
1468  &tangent0AtGaussPt(0, 2), &tangent1AtGaussPt(0, 2),
1469  &normalsAtGaussPts(0, 2), 3);
1470  for (int gg = 0; gg < nb_gauss_pts; ++gg) {
1471  CHKERR determinantTensor3by3(t_m_at_pts, det);
1472  CHKERR invertTensor3by3(t_m_at_pts, det, t_inv_m);
1473  for (int ll = 0; ll != nb_dofs; ll++) {
1474  t_transformed_h_curl(i) = t_inv_m(j, i) * t_h_curl(j);
1475  t_transformed_diff_h_curl(i, k) =
1476  t_inv_m(j, i) * t_diff_h_curl(j, k);
1477  ++t_h_curl;
1478  ++t_transformed_h_curl;
1479  ++t_diff_h_curl;
1480  ++t_transformed_diff_h_curl;
1481  ++cc;
1482  }
1483  ++t_m_at_pts;
1484  }
1485  } else {
1486  for (int gg = 0; gg < nb_gauss_pts; ++gg) {
1487  for (int ll = 0; ll != nb_dofs; ll++) {
1488  t_transformed_h_curl(i) = t_inv_m(j, i) * t_h_curl(j);
1489  t_transformed_diff_h_curl(i, k) =
1490  t_inv_m(j, i) * t_diff_h_curl(j, k);
1491  ++t_h_curl;
1492  ++t_transformed_h_curl;
1493  ++t_diff_h_curl;
1494  ++t_transformed_diff_h_curl;
1495  ++cc;
1496  }
1497  }
1498  }
1499  if (cc != nb_gauss_pts * nb_dofs)
1500  SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSIBLE_CASE, "Data inconsistency");
1501 
1502  baseN.data().swap(piola_n.data());
1503  diffBaseN.data().swap(diff_piola_n.data());
1504  }
1505  }
1506 
1508 }
MoFEMErrorCode determinantTensor3by3(T1 &t, T2 &det)
Calculate determinant 3 by 3.
Definition: Templates.hpp:418
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:77
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
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:399
FieldApproximationBase
approximation base
Definition: definitions.h:143
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:145
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

Member Data Documentation

◆ nOrmal

const VectorDouble& MoFEM::OpSetCovariantPiolaTransformOnFace::nOrmal

Definition at line 529 of file DataOperators.hpp.

◆ normalsAtGaussPts

const MatrixDouble& MoFEM::OpSetCovariantPiolaTransformOnFace::normalsAtGaussPts

Definition at line 530 of file DataOperators.hpp.

◆ tAngent0

const VectorDouble& MoFEM::OpSetCovariantPiolaTransformOnFace::tAngent0

Definition at line 531 of file DataOperators.hpp.

◆ tangent0AtGaussPt

const MatrixDouble& MoFEM::OpSetCovariantPiolaTransformOnFace::tangent0AtGaussPt

Definition at line 532 of file DataOperators.hpp.

◆ tAngent1

const VectorDouble& MoFEM::OpSetCovariantPiolaTransformOnFace::tAngent1

Definition at line 533 of file DataOperators.hpp.

◆ tangent1AtGaussPt

const MatrixDouble& MoFEM::OpSetCovariantPiolaTransformOnFace::tangent1AtGaussPt

Definition at line 534 of file DataOperators.hpp.


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