v0.13.1
Public Member Functions | Public Attributes | Protected Member Functions | List of all members
FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX Struct Reference

#include <users_modules/fracture_mechanics/src/CrackFrontElement.hpp>

Inherits VolumeElementForcesAndSourcesCore::UserDataOperator.

Collaboration diagram for FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX:
[legend]

Public Member Functions

 OpAleLhsWithDensitySingularElement_dx_dX (const std::string row_field, const std::string col_field, boost::shared_ptr< HookeElement::DataAtIntegrationPts > &data_at_pts, boost::shared_ptr< VectorDouble > rho_at_gauss_pts, boost::shared_ptr< MatrixDouble > rho_grad_at_gauss_pts, const double rho_n, const double rho_0, boost::shared_ptr< MatrixDouble > singular_displacement)
 

Public Attributes

boost::shared_ptr< VectorDouble > rhoAtGaussPtsPtr
 
boost::shared_ptr< MatrixDouble > rhoGradAtGaussPtsPtr
 
boost::shared_ptr< MatrixDouble > singularDisplacement
 
const double rhoN
 
const double rHo0
 
MatrixDouble K
 
MatrixDouble transK
 
VectorDouble nF
 
boost::shared_ptr< HookeElement::DataAtIntegrationPts > dataAtPts
 
VectorInt rowIndices
 
VectorInt colIndices
 
int nbRows
 number of dofs on rows More...
 
int nbCols
 number if dof on column More...
 
int nbIntegrationPts
 number of integration points More...
 
bool isDiag
 true if this block is on diagonal More...
 

Protected Member Functions

MoFEMErrorCode doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
 
MoFEMErrorCode iNtegrate (DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
 
MoFEMErrorCode aSsemble (DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
 

Detailed Description

Definition at line 710 of file CrackFrontElement.hpp.

Constructor & Destructor Documentation

◆ OpAleLhsWithDensitySingularElement_dx_dX()

FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::OpAleLhsWithDensitySingularElement_dx_dX ( const std::string  row_field,
const std::string  col_field,
boost::shared_ptr< HookeElement::DataAtIntegrationPts > &  data_at_pts,
boost::shared_ptr< VectorDouble >  rho_at_gauss_pts,
boost::shared_ptr< MatrixDouble >  rho_grad_at_gauss_pts,
const double  rho_n,
const double  rho_0,
boost::shared_ptr< MatrixDouble >  singular_displacement 
)

Definition at line 1335 of file CrackFrontElement.cpp.

1343 : UserDataOperator(row_field, col_field, OPROWCOL, false),
1344 dataAtPts(data_at_pts), rhoAtGaussPtsPtr(rho_at_gauss_pts),
1345 rhoGradAtGaussPtsPtr(rho_grad_at_gauss_pts),
1346 singularDisplacement(singular_displacement), rhoN(rho_n), rHo0(rho_0) {}
ForcesAndSourcesCore::UserDataOperator UserDataOperator
boost::shared_ptr< HookeElement::DataAtIntegrationPts > dataAtPts

Member Function Documentation

◆ aSsemble()

MoFEMErrorCode FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::aSsemble ( DataForcesAndSourcesCore::EntData row_data,
DataForcesAndSourcesCore::EntData col_data 
)
protected

Definition at line 1458 of file CrackFrontElement.cpp.

1460 {
1462
1463 const int *row_indices = &*row_data.getIndices().data().begin();
1464 const int *col_indices = &*col_data.getIndices().data().begin();
1465
1466 auto &data = *dataAtPts;
1467 Mat B = getFEMethod()->ksp_B != PETSC_NULL ? getFEMethod()->ksp_B
1468 : getFEMethod()->snes_B;
1469 // assemble local matrix
1470 CHKERR MatSetValues(B, nbRows, row_indices, nbCols, col_indices,
1471 &*K.data().begin(), ADD_VALUES);
1472
1473 if (!isDiag && sYmm) {
1474 // if not diagonal term and since global matrix is symmetric assemble
1475 // transpose term.
1476 transK.resize(K.size2(), K.size1(), false);
1477 noalias(transK) = trans(K);
1478 CHKERR MatSetValues(B, nbCols, col_indices, nbRows, row_indices,
1479 &*transK.data().begin(), ADD_VALUES);
1480 }
1482}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535

◆ doWork()

MoFEMErrorCode FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::doWork ( int  row_side,
int  col_side,
EntityType  row_type,
EntityType  col_type,
DataForcesAndSourcesCore::EntData row_data,
DataForcesAndSourcesCore::EntData col_data 
)
protected

Definition at line 1348 of file CrackFrontElement.cpp.

1351 {
1352
1354 nbRows = row_data.getIndices().size();
1355 if (!nbRows)
1357
1358 nbCols = col_data.getIndices().size();
1359 if (!nbCols)
1361
1362 K.resize(nbRows, nbCols, false);
1363 K.clear();
1364
1365 nbIntegrationPts = getGaussPts().size2();
1366 if (row_side == col_side && row_type == col_type) {
1367 isDiag = true;
1368 } else {
1369 isDiag = false;
1370 }
1371
1372 CHKERR iNtegrate(row_data, col_data);
1373 CHKERR aSsemble(row_data, col_data);
1374
1376}
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
MoFEMErrorCode iNtegrate(DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
MoFEMErrorCode aSsemble(DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)

◆ iNtegrate()

MoFEMErrorCode FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::iNtegrate ( DataForcesAndSourcesCore::EntData row_data,
DataForcesAndSourcesCore::EntData col_data 
)
protected

Definition at line 1378 of file CrackFrontElement.cpp.

1380 {
1382 auto get_tensor2 = [](MatrixDouble &m, const int r, const int c) {
1384 &m(r + 0, c + 0), &m(r + 0, c + 1), &m(r + 0, c + 2), &m(r + 1, c + 0),
1385 &m(r + 1, c + 1), &m(r + 1, c + 2), &m(r + 2, c + 0), &m(r + 2, c + 1),
1386 &m(r + 2, c + 2));
1387 };
1388
1393
1394 double vol = getVolume();
1395
1396 auto t_w = getFTensor0IntegrationWeight();
1397
1398 auto t_row_diff_base = row_data.getFTensor1DiffN<3>();
1399 const int row_nb_base_fun = row_data.getN().size2();
1400
1402 auto t_grad_rho = getFTensor1FromMat<3>(*rhoGradAtGaussPtsPtr);
1403
1404 auto t_cauchy_stress =
1405 getFTensor2SymmetricFromMat<3>(*(dataAtPts->cauchyStressMat));
1406 auto t_invH = getFTensor2FromMat<3, 3>(*dataAtPts->invHMat);
1407 auto &det_H = *dataAtPts->detHVec;
1408
1409 auto t_initial_singular_displacement =
1410 getFTensor1FromMat<3>(*singularDisplacement);
1411
1412 for (int gg = 0; gg != nbIntegrationPts; ++gg) {
1413
1414 double a = t_w * vol * det_H[gg];
1415 const double stress_dho_coef = (rhoN / rho);
1416
1417 int rr = 0;
1418 for (; rr != nbRows / 3; ++rr) {
1419
1420 auto t_m = get_tensor2(K, 3 * rr, 0);
1421
1422 FTensor::Tensor1<double, 3> t_row_diff_base_pulled;
1423 t_row_diff_base_pulled(i) = t_row_diff_base(j) * t_invH(j, i);
1424
1425 FTensor::Tensor1<double, 3> t_row_stress;
1426 t_row_stress(i) = a * t_row_diff_base_pulled(j) * t_cauchy_stress(i, j);
1428 t_a(i, k) = t_row_stress(i) * stress_dho_coef * t_grad_rho(k);
1429
1430 auto t_col_diff_base = col_data.getFTensor1DiffN<3>(gg, 0);
1431 auto t_col_base = col_data.getFTensor0N(gg, 0);
1432 for (int cc = 0; cc != nbCols / 3; ++cc) {
1433
1434 t_m(i, k) +=
1435 t_a(i, k) * t_col_diff_base(l) * t_initial_singular_displacement(l);
1436
1437 ++t_col_base;
1438 ++t_col_diff_base;
1439 ++t_m;
1440 }
1441 ++t_row_diff_base;
1442 }
1443
1444 for (; rr != row_nb_base_fun; ++rr)
1445 ++t_row_diff_base;
1446
1447 ++t_initial_singular_displacement;
1448 ++t_w;
1449 ++t_cauchy_stress;
1450 ++t_invH;
1451 ++rho;
1452 ++t_grad_rho;
1453 }
1454
1456}
constexpr double a
FTensor::Index< 'm', SPACE_DIM > m
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:77
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
Definition: Templates.hpp:135
const double r
rate factor
double rho
Definition: plastic.cpp:98

Member Data Documentation

◆ colIndices

VectorInt FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::colIndices

Definition at line 728 of file CrackFrontElement.hpp.

◆ dataAtPts

boost::shared_ptr<HookeElement::DataAtIntegrationPts> FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::dataAtPts

Definition at line 725 of file CrackFrontElement.hpp.

◆ isDiag

bool FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::isDiag

true if this block is on diagonal

Definition at line 733 of file CrackFrontElement.hpp.

◆ K

MatrixDouble FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::K

Definition at line 721 of file CrackFrontElement.hpp.

◆ nbCols

int FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::nbCols

number if dof on column

Definition at line 731 of file CrackFrontElement.hpp.

◆ nbIntegrationPts

int FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::nbIntegrationPts

number of integration points

Definition at line 732 of file CrackFrontElement.hpp.

◆ nbRows

int FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::nbRows

number of dofs on rows

Definition at line 730 of file CrackFrontElement.hpp.

◆ nF

VectorDouble FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::nF

Definition at line 723 of file CrackFrontElement.hpp.

◆ rHo0

const double FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::rHo0

Definition at line 718 of file CrackFrontElement.hpp.

◆ rhoAtGaussPtsPtr

boost::shared_ptr<VectorDouble> FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::rhoAtGaussPtsPtr

Definition at line 713 of file CrackFrontElement.hpp.

◆ rhoGradAtGaussPtsPtr

boost::shared_ptr<MatrixDouble> FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::rhoGradAtGaussPtsPtr

Definition at line 714 of file CrackFrontElement.hpp.

◆ rhoN

const double FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::rhoN

Definition at line 717 of file CrackFrontElement.hpp.

◆ rowIndices

VectorInt FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::rowIndices

Definition at line 727 of file CrackFrontElement.hpp.

◆ singularDisplacement

boost::shared_ptr<MatrixDouble> FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::singularDisplacement

Definition at line 715 of file CrackFrontElement.hpp.

◆ transK

MatrixDouble FractureMechanics::OpAleLhsWithDensitySingularElement_dx_dX::transK

Definition at line 722 of file CrackFrontElement.hpp.


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