v0.14.0
Loading...
Searching...
No Matches
Public Member Functions | Public Attributes | List of all members
BCs_RVELagrange_Disp::OpRVEBCsLhs Struct Reference

\biref operator to calculate the LHS for the RVE bounary conditions More...

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

Inheritance diagram for BCs_RVELagrange_Disp::OpRVEBCsLhs:
[legend]
Collaboration diagram for BCs_RVELagrange_Disp::OpRVEBCsLhs:
[legend]

Public Member Functions

 OpRVEBCsLhs (const string field_name, const string lagrang_field_name, Mat aij, RVEBC_Data &data, bool ho_geometry=false)
 
PetscErrorCode doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
 calculate thermal convection term in the lhs of equations More...
 

Public Attributes

RVEBC_DatadAta
 
bool hoGeometry
 
Mat Aij
 
MatrixDouble K
 
MatrixDouble transK
 

Detailed Description

\biref operator to calculate the LHS for the RVE bounary conditions

Definition at line 101 of file BCs_RVELagrange_Disp.hpp.

Constructor & Destructor Documentation

◆ OpRVEBCsLhs()

BCs_RVELagrange_Disp::OpRVEBCsLhs::OpRVEBCsLhs ( const string  field_name,
const string  lagrang_field_name,
Mat  aij,
RVEBC_Data data,
bool  ho_geometry = false 
)
inline

Definition at line 106 of file BCs_RVELagrange_Disp.hpp.

109 :
110 FaceElementForcesAndSourcesCore::UserDataOperator(
111 lagrang_field_name, field_name, UserDataOperator::OPROWCOL
112 ),
113 dAta(data),hoGeometry(ho_geometry),Aij(aij) {
114 // This will make sure to loop over all entities
115 // (e.g. for order=2 it will make doWork to loop 16 time)
116 sYmm = false;
117 }
constexpr auto field_name

Member Function Documentation

◆ doWork()

PetscErrorCode BCs_RVELagrange_Disp::OpRVEBCsLhs::doWork ( int  row_side,
int  col_side,
EntityType  row_type,
EntityType  col_type,
DataForcesAndSourcesCore::EntData &  row_data,
DataForcesAndSourcesCore::EntData &  col_data 
)
inline

calculate thermal convection term in the lhs of equations

C = intS N^T N dS

Definition at line 125 of file BCs_RVELagrange_Disp.hpp.

130 {
131 PetscFunctionBegin;
132
133 try {
134 int nb_row = row_data.getIndices().size();
135 int nb_col = col_data.getIndices().size();
136
137 if(nb_row == 0) PetscFunctionReturn(0);
138 if(nb_col == 0) PetscFunctionReturn(0);
139
140 auto weak_ptr_dof =
141 getNumeredEntFiniteElementPtr()->getRowDofsByPetscGlobalDofIdx(
142 row_data.getIndices()[0]);
143 const FENumeredDofEntity *dof_ptr;
144 if(auto ptr = weak_ptr_dof.lock())
145 dof_ptr = ptr.get();
146 else
147 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Dof not found");
148
149 int rank = dof_ptr->getNbOfCoeffs();
150
151 K.resize(nb_row/rank,nb_col/rank);
152 K.clear();
153
154 for(unsigned int gg = 0;gg<row_data.getN().size1();gg++) {
155 double area;
156 if(hoGeometry) {
157 area = norm_2(getNormalsAtGaussPts(gg))*0.5;
158 } else {
159 area = getArea();
160 }
161 double val = getGaussPts()(2,gg)*area;
162 noalias(K) += val*outer_prod(
163 row_data.getN(gg,nb_row/rank),col_data.getN(gg,nb_col/rank)
164 );
165 }
166
167 VectorInt row_indices,col_indices;
168 row_indices.resize(nb_row/rank);
169 col_indices.resize(nb_col/rank);
170
171 for(int rr = 0;rr < rank; rr++) {
172 unsigned int nb_rows;
173 unsigned int nb_cols;
174 int *rows;
175 int *cols;
176 if(rank > 1) {
177
178 ublas::noalias(row_indices) = ublas::vector_slice<VectorInt >
179 (row_data.getIndices(), ublas::slice(rr, rank, row_data.getIndices().size()/rank));
180
181 ublas::noalias(col_indices) = ublas::vector_slice<VectorInt >
182 (col_data.getIndices(), ublas::slice(rr, rank, col_data.getIndices().size()/rank));
183
184 nb_rows = row_indices.size();
185 nb_cols = col_indices.size();
186 rows = &*row_indices.data().begin();
187 cols = &*col_indices.data().begin();
188
189 } else {
190
191 nb_rows = row_data.getIndices().size();
192 nb_cols = col_data.getIndices().size();
193 rows = &*row_data.getIndices().data().begin();
194 cols = &*col_data.getIndices().data().begin();
195
196 }
197
198 // Matrix C
199 ierr = MatSetValues(Aij,nb_rows,rows,nb_cols,cols,&K(0,0),ADD_VALUES); CHKERRQ(ierr);
200
201 // Matrix CT
202 transK.resize(nb_col/rank,nb_row/rank);
203 noalias(transK) = trans(K);
204 ierr = MatSetValues(Aij,nb_cols,cols,nb_rows,rows,&transK(0,0),ADD_VALUES); CHKERRQ(ierr);
205
206 }
207
208 } catch (const std::exception& ex) {
209 ostringstream ss;
210 ss << "throw in method: " << ex.what() << endl;
211 SETERRQ(PETSC_COMM_SELF,1,ss.str().c_str());
212 }
213
214 PetscFunctionReturn(0);
215 }
static PetscErrorCode ierr
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
UBlasVector< int > VectorInt
Definition: Types.hpp:67
MoFEMErrorCode MatSetValues(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.

Member Data Documentation

◆ Aij

Mat BCs_RVELagrange_Disp::OpRVEBCsLhs::Aij

Definition at line 105 of file BCs_RVELagrange_Disp.hpp.

◆ dAta

RVEBC_Data& BCs_RVELagrange_Disp::OpRVEBCsLhs::dAta

Definition at line 103 of file BCs_RVELagrange_Disp.hpp.

◆ hoGeometry

bool BCs_RVELagrange_Disp::OpRVEBCsLhs::hoGeometry

Definition at line 104 of file BCs_RVELagrange_Disp.hpp.

◆ K

MatrixDouble BCs_RVELagrange_Disp::OpRVEBCsLhs::K

Definition at line 119 of file BCs_RVELagrange_Disp.hpp.

◆ transK

MatrixDouble BCs_RVELagrange_Disp::OpRVEBCsLhs::transK

Definition at line 119 of file BCs_RVELagrange_Disp.hpp.


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