v0.13.1
Public Member Functions | Public Attributes | List of all members
BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalCTLam Struct Reference

#include <users_modules/homogenisation/src/BCs_RVELagrange_Periodic.hpp>

Inherits FlatPrismElementForcesAndSurcesCore::UserDataOperator.

Collaboration diagram for BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalCTLam:
[legend]

Public Member Functions

 OpRVEBCsPeriodicCalCTLam (const string field_name, RVEBC_Data_Periodic &data, Vec f, CommonDataPeriodic &common_data_periodic, CommonFunctionsPeriodic &common_functions_periodic, bool ho_geometry=false)
 
PetscErrorCode doWork (int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
 

Public Attributes

Vec F
 
RVEBC_Data_PeriodicdAta
 
bool hoGeometry
 
CommonDataPeriodiccommonDataPeriodic
 
CommonFunctionsPeriodiccommonFunctionsPeriodic
 
ublas::vector< ublas::vector< int > > col_ind
 
VectorDouble fe1
 
VectorDouble fe2
 
MatrixDouble N_mat
 

Detailed Description

Definition at line 981 of file BCs_RVELagrange_Periodic.hpp.

Constructor & Destructor Documentation

◆ OpRVEBCsPeriodicCalCTLam()

BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalCTLam::OpRVEBCsPeriodicCalCTLam ( const string  field_name,
RVEBC_Data_Periodic data,
Vec  f,
CommonDataPeriodic common_data_periodic,
CommonFunctionsPeriodic common_functions_periodic,
bool  ho_geometry = false 
)

Definition at line 989 of file BCs_RVELagrange_Periodic.hpp.

996 :
998 field_name,UserDataOperator::OPCOL
999 ),
1000 F(f),
1001 dAta(data),
1002 hoGeometry(ho_geometry),
1003 commonDataPeriodic(common_data_periodic),
1004 commonFunctionsPeriodic(common_functions_periodic) {
1005 }
ForcesAndSourcesCore::UserDataOperator UserDataOperator
constexpr auto field_name

Member Function Documentation

◆ doWork()

PetscErrorCode BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalCTLam::doWork ( int  side,
EntityType  type,
DataForcesAndSourcesCore::EntData data 
)

Definition at line 1013 of file BCs_RVELagrange_Periodic.hpp.

1013 {
1014 PetscFunctionBegin;
1015// cout<<" OpRVEBCsPeriodicCalCTLam "<<endl;
1016 if(data.getIndices().size()==0) PetscFunctionReturn(0);
1017 int rank = data.getFieldDofs()[0]->getNbOfCoeffs();
1018 if(type == MBVERTEX) {
1019 col_ind.resize(2);
1020 int nb=data.getIndices().size()/2; //for one face for nodes only (total indieces are for both faces)
1021 col_ind[0].resize(nb); //1st face of prism -- nodes
1022 col_ind[1].resize(nb); //2nd face of prism -- nodes
1023 for(int ii=0; ii<nb; ii++){
1024 col_ind[0][ii]=data.getIndices()[ii];
1025 col_ind[1][ii]=data.getIndices()[ii+nb];
1026 }
1027 }
1028
1029
1030// cout<<"data.getN().size1() "<<data.getN().size1()<<endl;
1031 for(unsigned int gg = 0;gg<data.getN().size1();gg++) {
1032// cout<<"gg "<<gg<<endl;
1033 double area;
1034 VectorDouble normal_f3 = getNormalsAtGaussPtsF3(gg);
1035 area = cblas_dnrm2(3,&normal_f3[0],1)*0.5;
1036 double val = getGaussPts()(2,gg)*area;
1037
1038 if(type == MBVERTEX) {
1039 ierr = commonFunctionsPeriodic.shapeMat(rank, gg, data, N_mat, 2); CHKERRQ(ierr);}
1040 else {
1041 ierr = commonFunctionsPeriodic.shapeMat(rank, gg, data, N_mat); CHKERRQ(ierr);}
1042
1043
1044 if(gg==0){
1045 if(type == MBVERTEX) {
1046 fe1=-val*prod(trans(N_mat), commonDataPeriodic.LagMulAtGaussPts[gg]);
1047 fe2= val*prod(trans(N_mat), commonDataPeriodic.LagMulAtGaussPts[gg]);}
1048 else if(type == MBEDGE && side < 3){
1049 fe1=-val*prod(trans(N_mat), commonDataPeriodic.LagMulAtGaussPts[gg]);}
1050 else if(type == MBEDGE && side >= 6){
1051 fe2= val*prod(trans(N_mat), commonDataPeriodic.LagMulAtGaussPts[gg]);}
1052 else if(type == MBTRI && side == 3){
1053 fe1=-val*prod(trans(N_mat), commonDataPeriodic.LagMulAtGaussPts[gg]);}
1054 else if(type == MBTRI && side == 4){
1055 fe2= val*prod(trans(N_mat), commonDataPeriodic.LagMulAtGaussPts[gg]);}
1056 }else{
1057 if(type == MBVERTEX) {
1058 fe1+=-val*prod(trans(N_mat), commonDataPeriodic.LagMulAtGaussPts[gg]);
1059 fe2+= val*prod(trans(N_mat), commonDataPeriodic.LagMulAtGaussPts[gg]);}
1060 else if(type == MBEDGE && side < 3){
1061 fe1+=-val*prod(trans(N_mat), commonDataPeriodic.LagMulAtGaussPts[gg]);}
1062 else if(type == MBEDGE && side >= 6){
1063 fe2+= val*prod(trans(N_mat), commonDataPeriodic.LagMulAtGaussPts[gg]);}
1064 else if(type == MBTRI && side == 3){
1065 fe1+=-val*prod(trans(N_mat), commonDataPeriodic.LagMulAtGaussPts[gg]);}
1066 else if(type == MBTRI && side == 4){
1067 fe2+= val*prod(trans(N_mat), commonDataPeriodic.LagMulAtGaussPts[gg]);}
1068 }
1069
1070 }
1071// cout<<" Before assembley of F"<<endl;
1072
1073 if(type == MBVERTEX) {
1074 ierr = VecSetValues(F,col_ind[0].size(),&col_ind[0][0],&fe1[0],ADD_VALUES); CHKERRQ(ierr);
1075 ierr = VecSetValues(F,col_ind[1].size(),&col_ind[1][0],&fe2[0],ADD_VALUES); CHKERRQ(ierr);}
1076 else if(type == MBEDGE && side < 3){
1077 ierr = VecSetValues(F,data.getIndices().size(),&data.getIndices()[0],&fe1[0],ADD_VALUES); CHKERRQ(ierr);}
1078 else if(type == MBEDGE && side >= 6){
1079 ierr = VecSetValues(F,data.getIndices().size(),&data.getIndices()[0],&fe2[0],ADD_VALUES); CHKERRQ(ierr);}
1080 else if(type == MBTRI && side == 3){
1081 ierr = VecSetValues(F,data.getIndices().size(),&data.getIndices()[0],&fe1[0],ADD_VALUES); CHKERRQ(ierr);}
1082 else if(type == MBTRI && side == 4){
1083 ierr = VecSetValues(F,data.getIndices().size(),&data.getIndices()[0],&fe2[0],ADD_VALUES); CHKERRQ(ierr);}
1084 PetscFunctionReturn(0);
1085 }
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
UBlasVector< double > VectorDouble
Definition: Types.hpp:79
PetscErrorCode shapeMat(int rank, unsigned int gg, DataForcesAndSourcesCore::EntData &col_data, MatrixDouble &N_mat, int div=1)

Member Data Documentation

◆ col_ind

ublas::vector<ublas::vector<int> > BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalCTLam::col_ind

Definition at line 1008 of file BCs_RVELagrange_Periodic.hpp.

◆ commonDataPeriodic

CommonDataPeriodic& BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalCTLam::commonDataPeriodic

Definition at line 986 of file BCs_RVELagrange_Periodic.hpp.

◆ commonFunctionsPeriodic

CommonFunctionsPeriodic& BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalCTLam::commonFunctionsPeriodic

Definition at line 987 of file BCs_RVELagrange_Periodic.hpp.

◆ dAta

RVEBC_Data_Periodic& BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalCTLam::dAta

Definition at line 984 of file BCs_RVELagrange_Periodic.hpp.

◆ F

Vec BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalCTLam::F

Definition at line 983 of file BCs_RVELagrange_Periodic.hpp.

◆ fe1

VectorDouble BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalCTLam::fe1

Definition at line 1009 of file BCs_RVELagrange_Periodic.hpp.

◆ fe2

VectorDouble BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalCTLam::fe2

Definition at line 1009 of file BCs_RVELagrange_Periodic.hpp.

◆ hoGeometry

bool BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalCTLam::hoGeometry

Definition at line 985 of file BCs_RVELagrange_Periodic.hpp.

◆ N_mat

MatrixDouble BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalCTLam::N_mat

Definition at line 1011 of file BCs_RVELagrange_Periodic.hpp.


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