v0.10.0
Public Member Functions | Public Attributes | List of all members
EshelbianPlasticity::OpSpatialPhysical_du_du Struct Reference

#include <users_modules/eshelbian_plasticty/src/EshelbianPlasticity.hpp>

Inheritance diagram for EshelbianPlasticity::OpSpatialPhysical_du_du:
[legend]
Collaboration diagram for EshelbianPlasticity::OpSpatialPhysical_du_du:
[legend]

Public Member Functions

 OpSpatialPhysical_du_du (const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const double alpha)
 
MoFEMErrorCode integrate (EntData &row_data, EntData &col_data)
 
- Public Member Functions inherited from EshelbianPlasticity::OpAssembleVolume
 OpAssembleVolume (const std::string &field, boost::shared_ptr< DataAtIntegrationPts > data_ptr, const char type)
 
 OpAssembleVolume (const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > data_ptr, const char type, const bool assemble_symmetry)
 
- Public Member Functions inherited from EshelbianPlasticity::OpAssembleBasic< VolUserDataOperator >
 OpAssembleBasic (const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > data_ptr, const char type)
 
 OpAssembleBasic (const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > data_ptr, const char type, const bool assemble_symmetry)
 
virtual MoFEMErrorCode integrate (EntData &data)
 
virtual MoFEMErrorCode integrate (int row_side, EntityType row_type, EntData &data)
 
virtual MoFEMErrorCode assemble (EntData &data)
 
virtual MoFEMErrorCode assemble (int row_side, EntityType row_type, EntData &data)
 
virtual MoFEMErrorCode assemble (int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
 
MoFEMErrorCode doWork (int side, EntityType type, EntData &data)
 
MoFEMErrorCode doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
 

Public Attributes

const double alpha_u
 
- Public Attributes inherited from EshelbianPlasticity::OpAssembleBasic< VolUserDataOperator >
const bool assembleSymmetry
 
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
 data at integration pts More...
 
VectorDouble nF
 local right hand side vector More...
 
MatrixDouble K
 local tangent matrix More...
 
MatrixDouble transposeK
 

Detailed Description

Definition at line 959 of file EshelbianPlasticity.hpp.

Constructor & Destructor Documentation

◆ OpSpatialPhysical_du_du()

EshelbianPlasticity::OpSpatialPhysical_du_du::OpSpatialPhysical_du_du ( const std::string &  row_field,
const std::string &  col_field,
boost::shared_ptr< DataAtIntegrationPts > &  data_ptr,
const double  alpha 
)

Definition at line 961 of file EshelbianPlasticity.hpp.

965  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false),
966  alpha_u(alpha) {
967  sYmm = false;
968  }
EshelbianPlasticity::OpSpatialPhysical_du_du::alpha_u
const double alpha_u
Definition: EshelbianPlasticity.hpp:960
EshelbianPlasticity::OpAssembleVolume::OpAssembleVolume
OpAssembleVolume(const std::string &field, boost::shared_ptr< DataAtIntegrationPts > data_ptr, const char type)
Definition: EshelbianPlasticity.hpp:768

Member Function Documentation

◆ integrate()

MoFEMErrorCode EshelbianPlasticity::OpSpatialPhysical_du_du::integrate ( EntData &  row_data,
EntData &  col_data 
)
virtual

Reimplemented from EshelbianPlasticity::OpAssembleBasic< VolUserDataOperator >.

Definition at line 1280 of file EshelbianOperators.cpp.

1281  {
1282  MoFEMFunctionBegin;
1283 
1284  int nb_integration_pts = row_data.getN().size1();
1285  int row_nb_dofs = row_data.getIndices().size();
1286  int col_nb_dofs = col_data.getIndices().size();
1287  auto get_ftensor4 = [](MatrixDouble &m, const int r, const int c) {
1288  return FTensor::Ddg<FTensor::PackPtr<double *, 6>, 3, 3>(
1289 
1290  &m(r + 0, c + 0), &m(r + 0, c + 1), &m(r + 0, c + 2), &m(r + 0, c + 3),
1291  &m(r + 0, c + 4), &m(r + 0, c + 5),
1292 
1293  &m(r + 1, c + 0), &m(r + 1, c + 1), &m(r + 1, c + 2), &m(r + 1, c + 3),
1294  &m(r + 1, c + 4), &m(r + 1, c + 5),
1295 
1296  &m(r + 2, c + 0), &m(r + 2, c + 1), &m(r + 2, c + 2), &m(r + 2, c + 3),
1297  &m(r + 2, c + 4), &m(r + 2, c + 5),
1298 
1299  &m(r + 3, c + 0), &m(r + 3, c + 1), &m(r + 3, c + 2), &m(r + 3, c + 3),
1300  &m(r + 3, c + 4), &m(r + 3, c + 5),
1301 
1302  &m(r + 4, c + 0), &m(r + 4, c + 1), &m(r + 4, c + 2), &m(r + 4, c + 3),
1303  &m(r + 4, c + 4), &m(r + 4, c + 5),
1304 
1305  &m(r + 5, c + 0), &m(r + 5, c + 1), &m(r + 5, c + 2), &m(r + 5, c + 3),
1306  &m(r + 5, c + 4), &m(r + 5, c + 5)
1307 
1308  );
1309  };
1310  FTensor::Index<'i', 3> i;
1311  FTensor::Index<'j', 3> j;
1312  FTensor::Index<'k', 3> k;
1313  FTensor::Index<'l', 3> l;
1314  FTensor::Index<'m', 3> m;
1315  FTensor::Index<'n', 3> n;
1316 
1317  FTensor::Number<0> N0;
1318  FTensor::Number<1> N1;
1319  FTensor::Number<2> N2;
1320 
1321  auto v = getVolume();
1322  auto t_w = getFTensor0IntegrationWeight();
1323  auto t_P_dh0 = getFTensor3FromMat(dataAtPts->P_dh0);
1324  auto t_P_dh1 = getFTensor3FromMat(dataAtPts->P_dh1);
1325  auto t_P_dh2 = getFTensor3FromMat(dataAtPts->P_dh2);
1326  auto t_R = getFTensor2FromMat<3, 3>(dataAtPts->rotMatAtPts);
1327  auto t_approx_P = getFTensor2FromMat<3, 3>(dataAtPts->approxPAtPts);
1328  auto t_P = getFTensor2FromMat<3, 3>(dataAtPts->PAtPts);
1329  auto t_dot_log_u =
1330  getFTensor2SymmetricFromMat<3>(dataAtPts->logStreachDotTensorAtPts);
1331  auto t_u = getFTensor2SymmetricFromMat<3>(dataAtPts->streachTensorAtPts);
1332  auto t_diff_u = FTensor::Tensor4<FTensor::PackPtr<double *, 1>, 3, 3, 3, 3>(
1333  &dataAtPts->diffStreachTensorAtPts(0, 0),
1334  dataAtPts->diffStreachTensorAtPts.size2());
1335  std::array<FTensor::Tensor4<FTensor::PackPtr<double *, 1>, 3, 3, 3, 3>, 6>
1336  t_diff2_u = {
1337 
1338  FTensor::Tensor4<FTensor::PackPtr<double *, 1>, 3, 3, 3, 3>(
1339  &dataAtPts->expLogUHessian[0](0, 0),
1340  dataAtPts->expLogUHessian[0].size2()),
1341  FTensor::Tensor4<FTensor::PackPtr<double *, 1>, 3, 3, 3, 3>(
1342  &dataAtPts->expLogUHessian[1](0, 0),
1343  dataAtPts->expLogUHessian[1].size2()),
1344  FTensor::Tensor4<FTensor::PackPtr<double *, 1>, 3, 3, 3, 3>(
1345  &dataAtPts->expLogUHessian[2](0, 0),
1346  dataAtPts->expLogUHessian[2].size2()),
1347  FTensor::Tensor4<FTensor::PackPtr<double *, 1>, 3, 3, 3, 3>(
1348  &dataAtPts->expLogUHessian[3](0, 0),
1349  dataAtPts->expLogUHessian[3].size2()),
1350  FTensor::Tensor4<FTensor::PackPtr<double *, 1>, 3, 3, 3, 3>(
1351  &dataAtPts->expLogUHessian[4](0, 0),
1352  dataAtPts->expLogUHessian[4].size2()),
1353  FTensor::Tensor4<FTensor::PackPtr<double *, 1>, 3, 3, 3, 3>(
1354  &dataAtPts->expLogUHessian[5](0, 0),
1355  dataAtPts->expLogUHessian[5].size2())
1356 
1357  };
1358 
1359  constexpr auto t_kd = FTensor::Kronecker_Delta<int>();
1360  FTensor::Tensor4<double, 3, 3, 3, 3> t_one4;
1361  t_one4(i, j, k, l) = t_kd(j, l) * t_kd(i, k);
1362 
1363  FTensor::Tensor4<double, 3, 3, 3, 3> t_one4_symmetric;
1364  t_one4_symmetric(i, j, k, l) = 0;
1365  for (auto ii : {0, 1, 2})
1366  for (auto jj : {0, 1, 2})
1367  for (auto kk : {0, 1, 2})
1368  for (auto ll : {0, 1, 2}) {
1369 
1370  if (ll != kk)
1371  t_one4_symmetric(ii, jj, kk, ll) =
1372  t_one4(ii, jj, kk, ll) + t_one4(ii, jj, ll, kk);
1373  else
1374  t_one4_symmetric(ii, jj, kk, ll) = t_one4(ii, jj, kk, ll);
1375  }
1376 
1377  int row_nb_base_functions = row_data.getN().size2();
1378  auto t_row_base_fun = row_data.getFTensor0N();
1379 
1380  FTensor::Tensor2_symmetric<double, 3> t_one;
1381  t_one(i, j) = 0;
1382  for (auto ii : {0, 1, 2})
1383  t_one(ii, ii) = 1;
1384 
1385  const double ts_a = getTSa();
1386  for (int gg = 0; gg != nb_integration_pts; ++gg) {
1387  double a = v * t_w;
1388 
1389  FTensor::Tensor4<double, 3, 3, 3, 3> t_P_dh;
1390  t_P_dh(i, j, N0, k) = t_P_dh0(i, j, k);
1391  t_P_dh(i, j, N1, k) = t_P_dh1(i, j, k);
1392  t_P_dh(i, j, N2, k) = t_P_dh2(i, j, k);
1393 
1394  FTensor::Tensor2<double, 3, 3> t_RTP;
1395  t_RTP(i, j) = t_R(k, i) * t_approx_P(k, j);
1396  FTensor::Tensor2<double, 3, 3> t_deltaP;
1397  t_deltaP(i, j) = t_RTP(i, j) - t_P(i, j);
1398 
1399  FTensor::Tensor2<double, 3, 3> t_dot_u;
1400  t_dot_u(i, j) = t_u(i, k) * t_dot_log_u(k, j);
1401 
1402  FTensor::Tensor4<double, 3, 3, 3, 3> t_stress_diff;
1403  t_stress_diff(i, j, k, l) = 0;
1404  for (auto ii : {0, 1, 2})
1405  for (auto jj : {0, 1, 2})
1406  for (auto kk : {0, 1, 2})
1407  for (auto ll : {0, 1, 2})
1408  for (auto mm : {0, 1, 2})
1409  for (auto nn : {0, 1, 2})
1410  t_stress_diff(ii, jj, mm, nn) -=
1411  t_P_dh(ii, jj, kk, ll) * t_diff_u(kk, ll, mm, nn);
1412 
1413  // FTensor::Tensor4<double, 3, 3, 3, 3> t_dot_u_diff;
1414  // t_dot_u_diff(i, j, k, l) = 0;
1415  // for (auto ii : {0, 1, 2})
1416  // for (auto jj : {0, 1, 2})
1417  // for (auto kk : {0, 1, 2})
1418  // for (auto ll : {0, 1, 2})
1419  // for (auto mm : {0, 1, 2})
1420  // t_dot_u_diff(ii, jj, kk, ll) +=
1421  // t_diff_u(ii, mm, kk, ll) * t_dot_log_u(mm, jj) +
1422  // ts_a * t_u(ii, mm) * t_one4_symmetric(mm, jj, kk, ll);
1423 
1424  // auto t_viscous_stress = calculate_viscous_stress(t_dot_u);
1425  // t_viscous_stress(i, j) = t_dot_u(i, j);
1426  // t_viscous_stress(i, j) *= alpha_u * t_dot_log_u(i, j);
1427 
1428  // auto t_viscous_stress_diff =
1429  // calculate_diff_viscous_stress(t_dot_u, t_dot_u_diff);
1430  // t_viscous_stress_diff(i, j, k, l) = t_dot_u_diff(i, j, k, l);
1431  // t_viscous_stress_diff(i, j, k, l) *= alpha_u * t_one4(i, j, k, l);
1432 
1433  FTensor::Tensor2_symmetric<double, 3> t_viscous_stress;
1434  t_viscous_stress(i, j) = alpha_u * t_dot_log_u(i, j);
1435  FTensor::Tensor4<double, 3, 3, 3, 3> t_viscous_stress_diff;
1436  t_viscous_stress_diff(i, j, k, l) = t_one4_symmetric(i, j, k, l);
1437  t_viscous_stress_diff(i, j, k, l) *= (ts_a * alpha_u);
1438 
1439  int rr = 0;
1440  for (; rr != row_nb_dofs / 6; ++rr) {
1441 
1442  auto t_col_base_fun = col_data.getFTensor0N(gg, 0);
1443  auto t_m = get_ftensor4(K, 6 * rr, 0);
1444 
1445  for (int cc = 0; cc != col_nb_dofs / 6; ++cc) {
1446 
1447  const double b = a * t_row_base_fun * t_col_base_fun;
1448 
1449  for (int ii : {0, 1, 2})
1450  for (int jj = 0; jj <= ii; ++jj)
1451  for (int kk : {0, 1, 2})
1452  for (int ll = 0; ll <= kk; ++ll)
1453  for (int mm : {0, 1, 2})
1454  for (int nn : {0, 1, 2})
1455  t_m(ii, jj, kk, ll) +=
1456  b * t_diff_u(mm, nn, ii, jj) *
1457  (t_stress_diff(mm, nn, kk, ll) -
1458  t_viscous_stress_diff(mm, nn, kk, ll));
1459 
1460  for (int ii : {0, 1, 2})
1461  for (int jj = 0; jj <= ii; ++jj) {
1462  int ss = 0;
1463  for (int kk : {0, 1, 2})
1464  for (int ll = 0; ll <= kk; ++ll, ++ss)
1465  for (int mm : {0, 1, 2})
1466  for (int nn : {0, 1, 2})
1467  t_m(ii, jj, kk, ll) +=
1468  b * (t_diff2_u[ss])(mm, nn, ii, jj) *
1469  (t_deltaP(mm, nn) - t_viscous_stress(mm, nn));
1470  }
1471 
1472  ++t_m;
1473  ++t_col_base_fun;
1474  }
1475 
1476  ++t_row_base_fun;
1477  }
1478  for (; rr != row_nb_base_functions; ++rr)
1479  ++t_row_base_fun;
1480 
1481  ++t_w;
1482  ++t_P_dh0;
1483  ++t_P_dh1;
1484  ++t_P_dh2;
1485  ++t_R;
1486  ++t_approx_P;
1487  ++t_P;
1488  ++t_dot_log_u;
1489  ++t_u;
1490  ++t_diff_u;
1491  for (auto ll : {0, 1, 2, 3, 4, 5})
1492  ++(t_diff2_u[ll]);
1493  }
1494  MoFEMFunctionReturn(0);
1495 }
MoFEM::Types::MatrixDouble
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:76
FTensor::Kronecker_Delta
Kronecker Delta class.
Definition: Kronecker_Delta.hpp:15
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
EshelbianPlasticity::OpAssembleBasic< VolUserDataOperator >::dataAtPts
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
Definition: EshelbianPlasticity.hpp:638
l
static Index< 'l', 3 > l
Definition: BasicFeTools.hpp:56
n
static Index< 'n', 3 > n
Definition: BasicFeTools.hpp:58
EshelbianPlasticity::OpSpatialPhysical_du_du::alpha_u
const double alpha_u
Definition: EshelbianPlasticity.hpp:960
MoFEM::DataForcesAndSourcesCore::EntData::getIndices
const VectorInt & getIndices() const
Get global indices of dofs on entity.
Definition: DataStructures.hpp:1178
m
static Index< 'm', 3 > m
Definition: BasicFeTools.hpp:57
EshelbianPlasticity::getFTensor3FromMat
FTensor::Tensor3< FTensor::PackPtr< double *, 1 >, 3, 3, 3 > getFTensor3FromMat(MatrixDouble &m)
Definition: EshelbianPlasticity.hpp:24
i
static Index< 'i', 3 > i
Definition: BasicFeTools.hpp:53
N2
static Number< 2 > N2
Definition: BasicFeTools.hpp:64
j
static Index< 'j', 3 > j
Definition: BasicFeTools.hpp:54
ReactionDiffusionEquation::r
const double r
rate factor
Definition: reaction_diffusion.cpp:33
N1
static Number< 1 > N1
Definition: BasicFeTools.hpp:63
FTensor::Tensor2_symmetric
Definition: Tensor2_symmetric_value.hpp:13
k
static Index< 'k', 3 > k
Definition: BasicFeTools.hpp:55
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
MoFEM::DataForcesAndSourcesCore::EntData::getN
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
Definition: DataStructures.hpp:1288
N0
static Number< 0 > N0
Definition: BasicFeTools.hpp:62
MoFEM::DataForcesAndSourcesCore::EntData::getFTensor0N
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.
Definition: DataStructures.hpp:1459

Member Data Documentation

◆ alpha_u

const double EshelbianPlasticity::OpSpatialPhysical_du_du::alpha_u

Definition at line 960 of file EshelbianPlasticity.hpp.

The documentation for this struct was generated from the following files:
Generated on Sun Sep 27 2020 23:35:04 for MoFEM by Doxygen 1.8.15 and hosted at University of Glasgow