v0.15.0
Loading...
Searching...
No Matches
FreeSurfaceOps::OpLhsH_dH< I > Struct Template Reference

Lhs for H dH. More...

#include "tutorials/vec-5/src/FreeSurfaceOps.hpp"

Inheritance diagram for FreeSurfaceOps::OpLhsH_dH< I >:
[legend]
Collaboration diagram for FreeSurfaceOps::OpLhsH_dH< I >:
[legend]

Public Member Functions

 OpLhsH_dH (const std::string field_name, boost::shared_ptr< MatrixDouble > u_ptr, boost::shared_ptr< VectorDouble > h_ptr, boost::shared_ptr< MatrixDouble > grad_g_ptr)
 
MoFEMErrorCode iNtegrate (EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
 

Private Attributes

boost::shared_ptr< MatrixDouble > uPtr
 
boost::shared_ptr< VectorDouble > hPtr
 
boost::shared_ptr< MatrixDouble > gradGPtr
 

Detailed Description

template<bool I>
struct FreeSurfaceOps::OpLhsH_dH< I >

Lhs for H dH.

Examples
free_surface.cpp.

Definition at line 973 of file FreeSurfaceOps.hpp.

Constructor & Destructor Documentation

◆ OpLhsH_dH()

template<bool I>
FreeSurfaceOps::OpLhsH_dH< I >::OpLhsH_dH ( const std::string field_name,
boost::shared_ptr< MatrixDouble > u_ptr,
boost::shared_ptr< VectorDouble > h_ptr,
boost::shared_ptr< MatrixDouble > grad_g_ptr )
inline

Definition at line 975 of file FreeSurfaceOps.hpp.

979 AssemblyDomainEleOp::OPROWCOL),
980 uPtr(u_ptr), hPtr(h_ptr), gradGPtr(grad_g_ptr) {
981 sYmm = false;
982 }
constexpr auto field_name
boost::shared_ptr< MatrixDouble > uPtr
boost::shared_ptr< MatrixDouble > gradGPtr
boost::shared_ptr< VectorDouble > hPtr
FormsIntegrators< DomainEleOp >::Assembly< A >::OpBase AssemblyDomainEleOp

Member Function Documentation

◆ iNtegrate()

template<bool I>
MoFEMErrorCode FreeSurfaceOps::OpLhsH_dH< I >::iNtegrate ( EntitiesFieldData::EntData & row_data,
EntitiesFieldData::EntData & col_data )
inline

Definition at line 984 of file FreeSurfaceOps.hpp.

985 {
987
988 const double vol = getMeasure();
989 auto t_w = getFTensor0IntegrationWeight();
990 auto t_coords = getFTensor1CoordsAtGaussPts();
991 auto t_row_base = row_data.getFTensor0N();
992 auto t_row_diff_base = row_data.getFTensor1DiffN<SPACE_DIM>();
993
994#ifndef NDEBUG
995 if (row_data.getDiffN().size1() != row_data.getN().size1())
996 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "wrong size 1");
997 if (row_data.getDiffN().size2() != row_data.getN().size2() * SPACE_DIM) {
998 MOFEM_LOG("SELF", Sev::error)
999 << "Side " << rowSide << " " << CN::EntityTypeName(rowType);
1000 MOFEM_LOG("SELF", Sev::error) << row_data.getN();
1001 MOFEM_LOG("SELF", Sev::error) << row_data.getDiffN();
1002 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "wrong size 2");
1003 }
1004
1005 if (col_data.getDiffN().size1() != col_data.getN().size1())
1006 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "wrong size 1");
1007 if (col_data.getDiffN().size2() != col_data.getN().size2() * SPACE_DIM) {
1008 MOFEM_LOG("SELF", Sev::error)
1009 << "Side " << rowSide << " " << CN::EntityTypeName(rowType);
1010 MOFEM_LOG("SELF", Sev::error) << col_data.getN();
1011 MOFEM_LOG("SELF", Sev::error) << col_data.getDiffN();
1012 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "wrong size 2");
1013 }
1014#endif
1015
1016 if constexpr (I) {
1017
1018 auto t_h = getFTensor0FromVec(*hPtr);
1019 auto t_grad_g = getFTensor1FromMat<SPACE_DIM>(*gradGPtr);
1020
1021 for (int gg = 0; gg != nbIntegrationPts; gg++) {
1022
1023 const double r = t_coords(0);
1024 const double alpha = t_w * vol * cylindrical(r);
1025
1026 int rr = 0;
1027 for (; rr != nbRows; ++rr) {
1028
1029 auto t_col_base = col_data.getFTensor0N(gg, 0);
1030 auto t_col_diff_base = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
1031
1032 for (int cc = 0; cc != nbCols; ++cc) {
1033
1034 locMat(rr, cc) += (t_row_base * t_col_base * alpha);
1035
1036 ++t_col_base;
1037 ++t_col_diff_base;
1038 }
1039
1040 ++t_row_base;
1041 ++t_row_diff_base;
1042 }
1043
1044 for (; rr < nbRowBaseFunctions; ++rr) {
1045 ++t_row_base;
1046 ++t_row_diff_base;
1047 }
1048
1049 ++t_h;
1050 ++t_grad_g;
1051 ++t_w;
1052 ++t_coords;
1053 }
1054
1055 } else {
1056
1057 auto t_h = getFTensor0FromVec(*hPtr);
1058 auto t_grad_g = getFTensor1FromMat<SPACE_DIM>(*gradGPtr);
1059 auto t_u = getFTensor1FromMat<U_FIELD_DIM>(*uPtr);
1060
1061 auto ts_a = getTSa();
1062
1063 for (int gg = 0; gg != nbIntegrationPts; gg++) {
1064
1065 const double r = t_coords(0);
1066 const double alpha = t_w * vol * cylindrical(r);
1067
1068 auto m_dh = get_M_dh(t_h) * alpha;
1069 int rr = 0;
1070 for (; rr != nbRows; ++rr) {
1071
1072 auto t_col_base = col_data.getFTensor0N(gg, 0);
1073 auto t_col_diff_base = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
1074
1075 for (int cc = 0; cc != nbCols; ++cc) {
1076
1077 locMat(rr, cc) += (t_row_base * t_col_base * alpha) * ts_a;
1078 locMat(rr, cc) +=
1079 (t_row_base * alpha) * (t_col_diff_base(i) * t_u(i));
1080 locMat(rr, cc) +=
1081 (t_row_diff_base(i) * t_grad_g(i)) * (t_col_base * m_dh);
1082
1083 ++t_col_base;
1084 ++t_col_diff_base;
1085 }
1086
1087 ++t_row_base;
1088 ++t_row_diff_base;
1089 }
1090
1091 for (; rr < nbRowBaseFunctions; ++rr) {
1092 ++t_row_base;
1093 ++t_row_diff_base;
1094 }
1095
1096 ++t_u;
1097 ++t_h;
1098 ++t_grad_g;
1099 ++t_w;
1100 ++t_coords;
1101 }
1102 }
1103
1105 }
constexpr int SPACE_DIM
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_DATA_INCONSISTENCY
Definition definitions.h:31
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
auto cylindrical
auto get_M_dh
#define MOFEM_LOG(channel, severity)
Log.
FTensor::Index< 'i', SPACE_DIM > i
int r
Definition sdf.py:8
constexpr IntegrationType I

Member Data Documentation

◆ gradGPtr

template<bool I>
boost::shared_ptr<MatrixDouble> FreeSurfaceOps::OpLhsH_dH< I >::gradGPtr
private

Definition at line 1110 of file FreeSurfaceOps.hpp.

◆ hPtr

template<bool I>
boost::shared_ptr<VectorDouble> FreeSurfaceOps::OpLhsH_dH< I >::hPtr
private

Definition at line 1109 of file FreeSurfaceOps.hpp.

◆ uPtr

template<bool I>
boost::shared_ptr<MatrixDouble> FreeSurfaceOps::OpLhsH_dH< I >::uPtr
private

Definition at line 1108 of file FreeSurfaceOps.hpp.


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