v0.9.1
Public Member Functions | Public Attributes | List of all members
MoFEM::OpGetCoordsAndNormalsOnPrism Struct Reference

calculate normals at Gauss points of triangle element More...

#include <src/finite_elements/DataOperators.hpp>

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

Public Member Functions

 OpGetCoordsAndNormalsOnPrism (MatrixDouble &coords_at_gaussptf3, MatrixDouble &normals_at_gaussptf3, MatrixDouble &tangent1_at_gaussptf3, MatrixDouble &tangent2_at_gaussptf3, MatrixDouble &coords_at_gaussptf4, MatrixDouble &normals_at_gaussptf4, MatrixDouble &tangent1_at_gaussptf4, MatrixDouble &tangent2_at_gaussptf4)
 
MoFEMErrorCode doWork (int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
 Operator for linear form, usually to calculate values on right hand side. More...
 
MoFEMErrorCode calculateNormals ()
 
- 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

MatrixDoublecOords_at_GaussPtF3
 
MatrixDoublenOrmals_at_GaussPtF3
 
MatrixDoubletAngent1_at_GaussPtF3
 
MatrixDoubletAngent2_at_GaussPtF3
 
MatrixDoublecOords_at_GaussPtF4
 
MatrixDoublenOrmals_at_GaussPtF4
 
MatrixDoubletAngent1_at_GaussPtF4
 
MatrixDoubletAngent2_at_GaussPtF4
 
MatrixDouble sPin
 
- 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

calculate normals at Gauss points of triangle element

Definition at line 471 of file DataOperators.hpp.

Constructor & Destructor Documentation

◆ OpGetCoordsAndNormalsOnPrism()

MoFEM::OpGetCoordsAndNormalsOnPrism::OpGetCoordsAndNormalsOnPrism ( MatrixDouble coords_at_gaussptf3,
MatrixDouble normals_at_gaussptf3,
MatrixDouble tangent1_at_gaussptf3,
MatrixDouble tangent2_at_gaussptf3,
MatrixDouble coords_at_gaussptf4,
MatrixDouble normals_at_gaussptf4,
MatrixDouble tangent1_at_gaussptf4,
MatrixDouble tangent2_at_gaussptf4 
)

Definition at line 482 of file DataOperators.hpp.

487  : cOords_at_GaussPtF3(coords_at_gaussptf3),
488  nOrmals_at_GaussPtF3(normals_at_gaussptf3),
489  tAngent1_at_GaussPtF3(tangent1_at_gaussptf3),
490  tAngent2_at_GaussPtF3(tangent2_at_gaussptf3),
491  cOords_at_GaussPtF4(coords_at_gaussptf4),
492  nOrmals_at_GaussPtF4(normals_at_gaussptf4),
493  tAngent1_at_GaussPtF4(tangent1_at_gaussptf4),
494  tAngent2_at_GaussPtF4(tangent2_at_gaussptf4) {}

Member Function Documentation

◆ calculateNormals()

MoFEMErrorCode MoFEM::OpGetCoordsAndNormalsOnPrism::calculateNormals ( )

Definition at line 871 of file DataOperators.cpp.

871  {
873 
874  sPin.resize(3, 3);
875  sPin.clear();
876  nOrmals_at_GaussPtF3.resize(tAngent1_at_GaussPtF3.size1(), 3, false);
877  for (unsigned int gg = 0; gg < tAngent1_at_GaussPtF3.size1(); ++gg) {
878  ierr = Spin(&*sPin.data().begin(), &tAngent1_at_GaussPtF3(gg, 0));
879  CHKERRG(ierr);
880  cblas_dgemv(CblasRowMajor, CblasNoTrans, 3, 3, 1., &*sPin.data().begin(), 3,
881  &tAngent2_at_GaussPtF3(gg, 0), 1, 0.,
882  &nOrmals_at_GaussPtF3(gg, 0), 1);
883  }
884  sPin.clear();
885  nOrmals_at_GaussPtF4.resize(tAngent1_at_GaussPtF4.size1(), 3, false);
886  for (unsigned int gg = 0; gg < tAngent1_at_GaussPtF4.size1(); ++gg) {
887  ierr = Spin(&*sPin.data().begin(), &tAngent1_at_GaussPtF4(gg, 0));
888  CHKERRG(ierr);
889  cblas_dgemv(CblasRowMajor, CblasNoTrans, 3, 3, 1., &*sPin.data().begin(), 3,
890  &tAngent2_at_GaussPtF4(gg, 0), 1, 0.,
891  &nOrmals_at_GaussPtF4(gg, 0), 1);
892  }
894 }
void cblas_dgemv(const enum CBLAS_ORDER order, const enum CBLAS_TRANSPOSE TransA, const int M, const int N, const double alpha, const double *A, const int lda, const double *X, const int incX, const double beta, double *Y, const int incY)
Definition: cblas_dgemv.c:11
PetscErrorCode Spin(double *spinOmega, double *vecOmega)
calculate spin matrix from vector
Definition: fem_tools.c:558
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:550
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ doWork()

MoFEMErrorCode MoFEM::OpGetCoordsAndNormalsOnPrism::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 787 of file DataOperators.cpp.

788  {
790 
791  if (data.getFieldData().size() == 0)
793  const int valid_edges3[] = {1, 1, 1, 0, 0, 0, 0, 0, 0};
794  const int valid_faces3[] = {0, 0, 0, 1, 0, 0, 0, 0, 0};
795  const int valid_edges4[] = {0, 0, 0, 0, 0, 0, 1, 1, 1};
796  const int valid_faces4[] = {0, 0, 0, 0, 1, 0, 0, 0, 0};
797 
798  if (type == MBEDGE) {
799  if (!valid_edges3[side] || valid_edges4[side])
801  } else if (type == MBTRI) {
802  if (!valid_faces3[side] || valid_faces4[side])
804  }
805 
806  switch (type) {
807  case MBVERTEX: {
808  for (unsigned int gg = 0; gg < data.getN().size1(); ++gg) {
809  for (int dd = 0; dd < 3; dd++) {
810  cOords_at_GaussPtF3(gg, dd) =
811  cblas_ddot(3, &data.getN(gg)[0], 1, &data.getFieldData()[dd], 3);
813  3, &data.getDiffN()(gg, 0), 2, &data.getFieldData()[dd], 3);
815  3, &data.getDiffN()(gg, 1), 2, &data.getFieldData()[dd], 3);
817  3, &data.getN(gg)[0], 1, &data.getFieldData()[9 + dd], 3);
819  3, &data.getDiffN()(gg, 6 + 0), 2, &data.getFieldData()[9 + dd], 3);
821  3, &data.getDiffN()(gg, 6 + 1), 2, &data.getFieldData()[9 + dd], 3);
822  }
823  }
824  } break;
825  case MBEDGE:
826  case MBTRI: {
827  if (2 * data.getN().size2() != data.getDiffN().size2()) {
828  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
829  }
830  unsigned int nb_dofs = data.getFieldData().size();
831  if (nb_dofs % 3 != 0) {
832  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
833  }
834  if (nb_dofs > 3 * data.getN().size2()) {
835  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
836  "data inconsistency, side %d type %d", side, type);
837  }
838  for (unsigned int gg = 0; gg < data.getN().size1(); ++gg) {
839  for (int dd = 0; dd < 3; dd++) {
840  if ((type == MBTRI && valid_faces3[side]) ||
841  (type == MBEDGE && valid_edges3[side])) {
843  nb_dofs / 3, &data.getN(gg)[0], 1, &data.getFieldData()[dd], 3);
844  tAngent1_at_GaussPtF3(gg, dd) +=
845  cblas_ddot(nb_dofs / 3, &data.getDiffN()(gg, 0), 2,
846  &data.getFieldData()[dd], 3);
847  tAngent2_at_GaussPtF3(gg, dd) +=
848  cblas_ddot(nb_dofs / 3, &data.getDiffN()(gg, 1), 2,
849  &data.getFieldData()[dd], 3);
850  } else if ((type == MBTRI && valid_faces4[side]) ||
851  (type == MBEDGE && valid_edges4[side])) {
853  nb_dofs / 3, &data.getN(gg)[0], 1, &data.getFieldData()[dd], 3);
854  tAngent1_at_GaussPtF4(gg, dd) +=
855  cblas_ddot(nb_dofs / 3, &data.getDiffN()(gg, 0), 2,
856  &data.getFieldData()[dd], 3);
857  tAngent2_at_GaussPtF4(gg, dd) +=
858  cblas_ddot(nb_dofs / 3, &data.getDiffN()(gg, 1), 2,
859  &data.getFieldData()[dd], 3);
860  }
861  }
862  }
863  } break;
864  default:
865  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
866  }
867 
869 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
double cblas_ddot(const int N, const double *X, const int incX, const double *Y, const int incY)
Definition: cblas_ddot.c:12
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

Member Data Documentation

◆ cOords_at_GaussPtF3

MatrixDouble& MoFEM::OpGetCoordsAndNormalsOnPrism::cOords_at_GaussPtF3

Definition at line 473 of file DataOperators.hpp.

◆ cOords_at_GaussPtF4

MatrixDouble& MoFEM::OpGetCoordsAndNormalsOnPrism::cOords_at_GaussPtF4

Definition at line 477 of file DataOperators.hpp.

◆ nOrmals_at_GaussPtF3

MatrixDouble& MoFEM::OpGetCoordsAndNormalsOnPrism::nOrmals_at_GaussPtF3

Definition at line 474 of file DataOperators.hpp.

◆ nOrmals_at_GaussPtF4

MatrixDouble& MoFEM::OpGetCoordsAndNormalsOnPrism::nOrmals_at_GaussPtF4

Definition at line 478 of file DataOperators.hpp.

◆ sPin

MatrixDouble MoFEM::OpGetCoordsAndNormalsOnPrism::sPin

Definition at line 496 of file DataOperators.hpp.

◆ tAngent1_at_GaussPtF3

MatrixDouble& MoFEM::OpGetCoordsAndNormalsOnPrism::tAngent1_at_GaussPtF3

Definition at line 475 of file DataOperators.hpp.

◆ tAngent1_at_GaussPtF4

MatrixDouble& MoFEM::OpGetCoordsAndNormalsOnPrism::tAngent1_at_GaussPtF4

Definition at line 479 of file DataOperators.hpp.

◆ tAngent2_at_GaussPtF3

MatrixDouble& MoFEM::OpGetCoordsAndNormalsOnPrism::tAngent2_at_GaussPtF3

Definition at line 476 of file DataOperators.hpp.

◆ tAngent2_at_GaussPtF4

MatrixDouble& MoFEM::OpGetCoordsAndNormalsOnPrism::tAngent2_at_GaussPtF4

Definition at line 480 of file DataOperators.hpp.


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