\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:

Collaboration diagram for BCs_RVELagrange_Disp::OpRVEBCsLhs:

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_Data &	dAta

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.

      :
     FaceElementForcesAndSourcesCore::UserDataOperator(
       lagrang_field_name, field_name, UserDataOperator::OPROWCOL
     ),
     dAta(data),hoGeometry(ho_geometry),Aij(aij) {
       // This will make sure to loop over all entities
       // (e.g. for order=2 it will make doWork to loop 16 time)
       sYmm = false;
     }

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.

       {
       PetscFunctionBegin;
  
       try {
         int nb_row = row_data.getIndices().size();
         int nb_col = col_data.getIndices().size();
  
         if(nb_row == 0) PetscFunctionReturn(0);
         if(nb_col == 0) PetscFunctionReturn(0);
  
         auto weak_ptr_dof =
             getNumeredEntFiniteElementPtr()->getRowDofsByPetscGlobalDofIdx(
                 row_data.getIndices()[0]);
         const FENumeredDofEntity *dof_ptr;
         if(auto ptr = weak_ptr_dof.lock()) 
           dof_ptr = ptr.get();
         else
           SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Dof not found");
  
         int rank = dof_ptr->getNbOfCoeffs();
  
         K.resize(nb_row/rank,nb_col/rank);
         K.clear();
  
         for(unsigned int gg = 0;gg<row_data.getN().size1();gg++) {
           double area;
           if(hoGeometry) {
             area = norm_2(getNormalsAtGaussPts(gg))*0.5;
           }   else {
             area = getArea();
           }
           double val = getGaussPts()(2,gg)*area;
           noalias(K) += val*outer_prod(
             row_data.getN(gg,nb_row/rank),col_data.getN(gg,nb_col/rank)
           );
         }
  
         VectorInt row_indices,col_indices;
         row_indices.resize(nb_row/rank);
         col_indices.resize(nb_col/rank);
  
         for(int rr = 0;rr < rank; rr++) {
           unsigned int nb_rows;
           unsigned int nb_cols;
           int *rows;
           int *cols;
           if(rank > 1) {
  
             ublas::noalias(row_indices) = ublas::vector_slice<VectorInt >
             (row_data.getIndices(), ublas::slice(rr, rank, row_data.getIndices().size()/rank));
  
             ublas::noalias(col_indices) = ublas::vector_slice<VectorInt >
             (col_data.getIndices(), ublas::slice(rr, rank, col_data.getIndices().size()/rank));
  
             nb_rows = row_indices.size();
             nb_cols = col_indices.size();
             rows = &*row_indices.data().begin();
             cols = &*col_indices.data().begin();
  
           } else {
  
             nb_rows = row_data.getIndices().size();
             nb_cols = col_data.getIndices().size();
             rows = &*row_data.getIndices().data().begin();
             cols = &*col_data.getIndices().data().begin();
  
           }
  
           // Matrix C
           ierr = MatSetValues(Aij,nb_rows,rows,nb_cols,cols,&K(0,0),ADD_VALUES); CHKERRQ(ierr);
  
           // Matrix CT
           transK.resize(nb_col/rank,nb_row/rank);
           noalias(transK) = trans(K);
           ierr = MatSetValues(Aij,nb_cols,cols,nb_rows,rows,&transK(0,0),ADD_VALUES); CHKERRQ(ierr);
  
         }
  
       } catch (const std::exception& ex) {
         ostringstream ss;
         ss << "throw in method: " << ex.what() << endl;
         SETERRQ(PETSC_COMM_SELF,1,ss.str().c_str());
       }
  
       PetscFunctionReturn(0);
     }

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