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

\biref operator to calculate and assemble Cmat More...

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

Inherits FlatPrismElementForcesAndSurcesCore::UserDataOperator.

Collaboration diagram for BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalAssemCmat:
[legend]

Public Member Functions

 OpRVEBCsPeriodicCalAssemCmat (const string field_name, const string lagrang_field_name, Mat aij, RVEBC_Data_Periodic &data, CommonFunctionsPeriodic &common_functions_periodic, CommonDataPeriodic &common_data_periodic, 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)
 

Public Attributes

Mat Aij
 
RVEBC_Data_PeriodicdAta
 
CommonFunctionsPeriodic commonFunctionsPeriodic
 
CommonDataPeriodiccommonDataPeriodic
 
bool hoGeometry
 
MatrixDouble Kmat
 
MatrixDouble transKmat
 
ublas::vector< int > rrvec
 
ublas::vector< int > ccvec
 
MatrixDouble NTN
 

Detailed Description

\biref operator to calculate and assemble Cmat

Definition at line 338 of file BCs_RVELagrange_Periodic.hpp.

Constructor & Destructor Documentation

◆ OpRVEBCsPeriodicCalAssemCmat()

BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalAssemCmat::OpRVEBCsPeriodicCalAssemCmat ( const string  field_name,
const string  lagrang_field_name,
Mat  aij,
RVEBC_Data_Periodic data,
CommonFunctionsPeriodic common_functions_periodic,
CommonDataPeriodic common_data_periodic,
bool  ho_geometry = false 
)

Definition at line 345 of file BCs_RVELagrange_Periodic.hpp.

353 :
354
356 lagrang_field_name, field_name, UserDataOperator::OPROWCOL
357 ),Aij(aij),
358 dAta(data),
359 commonFunctionsPeriodic(common_functions_periodic),
360 commonDataPeriodic(common_data_periodic),
361 hoGeometry(ho_geometry){
362 sYmm = false; //This will make sure to loop over all intities (e.g. for order=2 it will make doWork to loop 16 time)
363 }
ForcesAndSourcesCore::UserDataOperator UserDataOperator
constexpr auto field_name

Member Function Documentation

◆ doWork()

PetscErrorCode BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalAssemCmat::doWork ( int  row_side,
int  col_side,
EntityType  row_type,
EntityType  col_type,
DataForcesAndSourcesCore::EntData row_data,
DataForcesAndSourcesCore::EntData col_data 
)

Definition at line 369 of file BCs_RVELagrange_Periodic.hpp.

374 {
375 PetscFunctionBegin;
376
377 try {
378// cout<<"Hello from OpRVEBCsPeriodicLhsAssemCmat "<<endl;
379 if(row_data.getIndices().size()==0) PetscFunctionReturn(0);
380 if(col_data.getIndices().size()==0) PetscFunctionReturn(0);
381 if(row_type == MBEDGE && row_side >= 3) PetscFunctionReturn(0); //ignore second triangel edges (cananical number 6 to 8)
382 if(col_type == MBEDGE && col_side >= 3) PetscFunctionReturn(0); //ignore second triangel edges (cananical number 6 to 8)
383 if(row_type == MBTRI && row_side == 4) PetscFunctionReturn(0); //ignore second triangel face (cananical number 4)
384 if(col_type == MBTRI && col_side == 4) PetscFunctionReturn(0); //ignore second triangel face (cananical number 4)
385
386
387 //Cmat calculation
388
389 int rank = col_data.getFieldDofs()[0]->getNbOfCoeffs();
390
391 MatrixDouble N_mat_row;
392 MatrixDouble N_mat_col;
393
394
395 //1st face
396 for(unsigned int gg = 0;gg<col_data.getN().size1();gg++) {
397 double area;
398 VectorDouble normal_f3 = getNormalsAtGaussPtsF3(gg);
399 area = cblas_dnrm2(3,&normal_f3[0],1)*0.5;
400 double val = getGaussPts()(2,gg)*area;
401
402 if(row_type == MBVERTEX) {
403 ierr = commonFunctionsPeriodic.shapeMatNew(rank, gg, commonDataPeriodic.RowN[0][row_type][row_side], N_mat_row, 2); CHKERRQ(ierr);
404// cerr << "N_mat_row "<< N_mat_row << endl;
405 }
406 else {
407 ierr = commonFunctionsPeriodic.shapeMatNew(rank, gg, commonDataPeriodic.RowN[0][row_type][row_side], N_mat_row); CHKERRQ(ierr);}
408
409 if(col_type == MBVERTEX) {
410 ierr = commonFunctionsPeriodic.shapeMatNew(rank, gg, commonDataPeriodic.ColN[0][col_type][col_side], N_mat_col, 2); CHKERRQ(ierr);
411// cerr << "N_mat_row "<< N_mat_row << endl;
412 }
413 else {
414 ierr = commonFunctionsPeriodic.shapeMatNew(rank, gg, commonDataPeriodic.ColN[0][col_type][col_side], N_mat_col); CHKERRQ(ierr);}
415 if(gg==0){
416// cerr << "row_type "<< row_type <<" row_side "<<row_side << " "<< endl;
417// cerr << "col_type "<< col_type <<" col_side "<< col_side<< " "<< endl;
418 NTN=prod(trans(N_mat_row),N_mat_col); //we don't need to define size of NTN
419 Kmat=-1.0*val*NTN; //we don't need to defien size of K
420 }else{
421 NTN=prod(trans(N_mat_row),N_mat_col);
422 Kmat+=-1.0*val*NTN;
423 }
424 }
425 rrvec=commonDataPeriodic.RowInd[0][row_type][row_side];
426 ccvec=commonDataPeriodic.ColInd[0][col_type][col_side];
427 ierr = MatSetValues(Aij,rrvec.size(),&rrvec[0],ccvec.size(),&ccvec[0],&Kmat(0,0),ADD_VALUES); CHKERRQ(ierr);
428 transKmat = trans(Kmat);
429 ierr = MatSetValues(Aij,ccvec.size(),&ccvec[0],rrvec.size(),&rrvec[0],&transKmat(0,0),ADD_VALUES); CHKERRQ(ierr);
430
431// cerr << "rrvec "<< rrvec << endl;
432// cerr << "ccvec "<< ccvec << endl;
433// cerr << "Kmat "<< Kmat << endl;
434
435
436 //2nd face
437 for(unsigned int gg = 0;gg<col_data.getN().size1();gg++) {
438 double area;
439 VectorDouble normal_f3 = getNormalsAtGaussPtsF3(gg);
440 area = cblas_dnrm2(3,&normal_f3[0],1)*0.5;
441 double val = getGaussPts()(2,gg)*area;
442
443 if(row_type == MBVERTEX) {
444 ierr = commonFunctionsPeriodic.shapeMatNew(rank, gg, commonDataPeriodic.RowN[1][row_type][row_side], N_mat_row, 2); CHKERRQ(ierr);
445// cerr << "N_mat_row "<< N_mat_row << endl;
446 }
447 else {
448 ierr = commonFunctionsPeriodic.shapeMatNew(rank, gg, commonDataPeriodic.RowN[1][row_type][row_side], N_mat_row); CHKERRQ(ierr);}
449
450 if(col_type == MBVERTEX) {
451 ierr = commonFunctionsPeriodic.shapeMatNew(rank, gg, commonDataPeriodic.ColN[1][col_type][col_side], N_mat_col, 2); CHKERRQ(ierr);
452// cerr << "N_mat_row "<< N_mat_row << endl;
453 }
454 else {
455 ierr = commonFunctionsPeriodic.shapeMatNew(rank, gg, commonDataPeriodic.ColN[1][col_type][col_side], N_mat_col); CHKERRQ(ierr);}
456 if(gg==0){
457// cerr << "row_type "<< row_type <<" row_side "<<row_side << " "<< endl;
458// cerr << "col_type "<< col_type <<" col_side "<< col_side<< " "<< endl;
459 NTN=prod(trans(N_mat_row),N_mat_col); //we don't need to define size of NTN
460 Kmat=val*NTN; //we don't need to defien size of K
461 }else{
462 NTN=prod(trans(N_mat_row),N_mat_col);
463 Kmat+=val*NTN;
464 }
465 }
466 rrvec=commonDataPeriodic.RowInd[1][row_type][row_side];
467 ccvec=commonDataPeriodic.ColInd[1][col_type][col_side];
468 ierr = MatSetValues(Aij,rrvec.size(),&rrvec[0],ccvec.size(),&ccvec[0],&Kmat(0,0),ADD_VALUES); CHKERRQ(ierr);
469 transKmat = trans(Kmat);
470 ierr = MatSetValues(Aij,ccvec.size(),&ccvec[0],rrvec.size(),&rrvec[0],&transKmat(0,0),ADD_VALUES); CHKERRQ(ierr);
471
472
473// cerr << "\n\n\nrrvec "<< rrvec << endl;
474// cerr << "ccvec "<< ccvec << endl;
475// cerr << "Kmat "<< Kmat << endl;
476
477
478 } catch (const std::exception& ex) {
479 ostringstream ss;
480 ss << "throw in method: " << ex.what() << endl;
481 SETERRQ(PETSC_COMM_SELF,1,ss.str().c_str());
482 }
483// string wait;
484// cin>>wait;
485 PetscFunctionReturn(0);
486 }
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:88
UBlasVector< double > VectorDouble
Definition: Types.hpp:79
ublas::vector< map< EntityType, map< int, ublas::vector< int > > > > RowInd
ublas::vector< map< EntityType, map< int, MatrixDouble > > > ColN
ublas::vector< map< EntityType, map< int, MatrixDouble > > > RowN
ublas::vector< map< EntityType, map< int, ublas::vector< int > > > > ColInd
PetscErrorCode shapeMatNew(int rank, unsigned int gg, MatrixDouble &RowN, MatrixDouble &N_mat, int div=1)

Member Data Documentation

◆ Aij

Mat BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalAssemCmat::Aij

Definition at line 339 of file BCs_RVELagrange_Periodic.hpp.

◆ ccvec

ublas::vector<int> BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalAssemCmat::ccvec

Definition at line 366 of file BCs_RVELagrange_Periodic.hpp.

◆ commonDataPeriodic

CommonDataPeriodic& BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalAssemCmat::commonDataPeriodic

Definition at line 342 of file BCs_RVELagrange_Periodic.hpp.

◆ commonFunctionsPeriodic

CommonFunctionsPeriodic BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalAssemCmat::commonFunctionsPeriodic

Definition at line 341 of file BCs_RVELagrange_Periodic.hpp.

◆ dAta

RVEBC_Data_Periodic& BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalAssemCmat::dAta

Definition at line 340 of file BCs_RVELagrange_Periodic.hpp.

◆ hoGeometry

bool BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalAssemCmat::hoGeometry

Definition at line 343 of file BCs_RVELagrange_Periodic.hpp.

◆ Kmat

MatrixDouble BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalAssemCmat::Kmat

Definition at line 365 of file BCs_RVELagrange_Periodic.hpp.

◆ NTN

MatrixDouble BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalAssemCmat::NTN

Definition at line 367 of file BCs_RVELagrange_Periodic.hpp.

◆ rrvec

ublas::vector<int> BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalAssemCmat::rrvec

Definition at line 366 of file BCs_RVELagrange_Periodic.hpp.

◆ transKmat

MatrixDouble BCs_RVELagrange_Periodic::OpRVEBCsPeriodicCalAssemCmat::transKmat

Definition at line 365 of file BCs_RVELagrange_Periodic.hpp.


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