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 664 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 1186 of file CrackFrontElement.cpp.

1194 : UserDataOperator(row_field, col_field, OPROWCOL, true),
1195 dataAtPts(data_at_pts), rhoAtGaussPtsPtr(rho_at_gauss_pts),
1196 rhoGradAtGaussPtsPtr(rho_grad_at_gauss_pts),
1197 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 1309 of file CrackFrontElement.cpp.

1311 {
1313
1314 const int *row_indices = &*row_data.getIndices().data().begin();
1315 const int *col_indices = &*col_data.getIndices().data().begin();
1316
1317 auto &data = *dataAtPts;
1318 Mat B = getFEMethod()->ksp_B != PETSC_NULL ? getFEMethod()->ksp_B
1319 : getFEMethod()->snes_B;
1320 // assemble local matrix
1321 CHKERR MatSetValues(B, nbRows, row_indices, nbCols, col_indices,
1322 &*K.data().begin(), ADD_VALUES);
1323
1324 if (!isDiag && sYmm) {
1325 // if not diagonal term and since global matrix is symmetric assemble
1326 // transpose term.
1327 transK.resize(K.size2(), K.size1(), false);
1328 noalias(transK) = trans(K);
1329 CHKERR MatSetValues(B, nbCols, col_indices, nbRows, row_indices,
1330 &*transK.data().begin(), ADD_VALUES);
1331 }
1333}
#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 1199 of file CrackFrontElement.cpp.

1202 {
1203
1205 nbRows = row_data.getIndices().size();
1206 if (!nbRows)
1208
1209 nbCols = col_data.getIndices().size();
1210 if (!nbCols)
1212
1213 K.resize(nbRows, nbCols, false);
1214 K.clear();
1215
1216 nbIntegrationPts = getGaussPts().size2();
1217 if (row_side == col_side && row_type == col_type) {
1218 isDiag = true;
1219 } else {
1220 isDiag = false;
1221 }
1222
1223 CHKERR iNtegrate(row_data, col_data);
1224 CHKERR aSsemble(row_data, col_data);
1225
1227}
#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 1229 of file CrackFrontElement.cpp.

1231 {
1233 auto get_tensor2 = [](MatrixDouble &m, const int r, const int c) {
1235 &m(r + 0, c + 0), &m(r + 0, c + 1), &m(r + 0, c + 2), &m(r + 1, c + 0),
1236 &m(r + 1, c + 1), &m(r + 1, c + 2), &m(r + 2, c + 0), &m(r + 2, c + 1),
1237 &m(r + 2, c + 2));
1238 };
1239
1244
1245 double vol = getVolume();
1246
1247 auto t_w = getFTensor0IntegrationWeight();
1248
1249 auto t_row_diff_base = row_data.getFTensor1DiffN<3>();
1250 const int row_nb_base_fun = row_data.getN().size2();
1251
1253 auto t_grad_rho = getFTensor1FromMat<3>(*rhoGradAtGaussPtsPtr);
1254
1255 auto t_eshelby_stress =
1256 getFTensor2FromMat<3, 3>(*dataAtPts->eshelbyStressMat);
1257 auto t_invH = getFTensor2FromMat<3, 3>(*dataAtPts->invHMat);
1258 auto &det_H = *dataAtPts->detHVec;
1259
1260 auto t_initial_singular_displacement =
1261 getFTensor1FromMat<3>(*singularDisplacement);
1262
1263 for (int gg = 0; gg != nbIntegrationPts; ++gg) {
1264
1265 double a = t_w * vol * det_H[gg];
1266 const double stress_dho_coef = (rhoN / rho);
1267
1268 int rr = 0;
1269 for (; rr != nbRows / 3; ++rr) {
1270
1271 auto t_m = get_tensor2(K, 3 * rr, 0);
1272
1273 FTensor::Tensor1<double, 3> t_row_diff_base_pulled;
1274 t_row_diff_base_pulled(i) = t_row_diff_base(j) * t_invH(j, i);
1275
1276 FTensor::Tensor1<double, 3> t_row_stress;
1277 t_row_stress(i) = a * t_row_diff_base_pulled(j) * t_eshelby_stress(i, j);
1279 t_a(i, k) = t_row_stress(i) * stress_dho_coef * t_grad_rho(k);
1280
1281 auto t_col_diff_base = col_data.getFTensor1DiffN<3>(gg, 0);
1282 auto t_col_base = col_data.getFTensor0N(gg, 0);
1283 for (int cc = 0; cc != nbCols / 3; ++cc) {
1284
1285 t_m(i, k) +=
1286 t_a(i, k) * t_col_diff_base(l) * t_initial_singular_displacement(l);
1287
1288 ++t_col_base;
1289 ++t_col_diff_base;
1290 ++t_m;
1291 }
1292 ++t_row_diff_base;
1293 }
1294
1295 for (; rr != row_nb_base_fun; ++rr)
1296 ++t_row_diff_base;
1297
1298 ++t_initial_singular_displacement;
1299 ++t_w;
1300 ++t_eshelby_stress;
1301 ++t_invH;
1302 ++rho;
1303 ++t_grad_rho;
1304 }
1305
1307}
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 682 of file CrackFrontElement.hpp.

◆ dataAtPts

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

Definition at line 679 of file CrackFrontElement.hpp.

◆ isDiag

bool FractureMechanics::OpAleLhsWithDensitySingularElement_dX_dX::isDiag

true if this block is on diagonal

Definition at line 687 of file CrackFrontElement.hpp.

◆ K

MatrixDouble FractureMechanics::OpAleLhsWithDensitySingularElement_dX_dX::K

Definition at line 675 of file CrackFrontElement.hpp.

◆ nbCols

int FractureMechanics::OpAleLhsWithDensitySingularElement_dX_dX::nbCols

number if dof on column

Definition at line 685 of file CrackFrontElement.hpp.

◆ nbIntegrationPts

int FractureMechanics::OpAleLhsWithDensitySingularElement_dX_dX::nbIntegrationPts

number of integration points

Definition at line 686 of file CrackFrontElement.hpp.

◆ nbRows

int FractureMechanics::OpAleLhsWithDensitySingularElement_dX_dX::nbRows

number of dofs on rows

Definition at line 684 of file CrackFrontElement.hpp.

◆ nF

VectorDouble FractureMechanics::OpAleLhsWithDensitySingularElement_dX_dX::nF

Definition at line 677 of file CrackFrontElement.hpp.

◆ rHo0

const double FractureMechanics::OpAleLhsWithDensitySingularElement_dX_dX::rHo0

Definition at line 672 of file CrackFrontElement.hpp.

◆ rhoAtGaussPtsPtr

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

Definition at line 667 of file CrackFrontElement.hpp.

◆ rhoGradAtGaussPtsPtr

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

Definition at line 668 of file CrackFrontElement.hpp.

◆ rhoN

const double FractureMechanics::OpAleLhsWithDensitySingularElement_dX_dX::rhoN

Definition at line 671 of file CrackFrontElement.hpp.

◆ rowIndices

VectorInt FractureMechanics::OpAleLhsWithDensitySingularElement_dX_dX::rowIndices

Definition at line 681 of file CrackFrontElement.hpp.

◆ singularDisplacement

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

Definition at line 669 of file CrackFrontElement.hpp.

◆ transK

MatrixDouble FractureMechanics::OpAleLhsWithDensitySingularElement_dX_dX::transK

Definition at line 676 of file CrackFrontElement.hpp.


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