v0.14.0
Public Member Functions | Public Attributes | List of all members
EshelbianPlasticity::OpCalculateRotationAndSpatialGradient Struct Reference

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

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

Public Member Functions

 OpCalculateRotationAndSpatialGradient (boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
 
MoFEMErrorCode doWork (int side, EntityType type, EntData &data)
 

Public Attributes

boost::shared_ptr< DataAtIntegrationPts > dataAtPts
 data at integration pts More...
 

Detailed Description

Definition at line 488 of file EshelbianPlasticity.hpp.

Constructor & Destructor Documentation

◆ OpCalculateRotationAndSpatialGradient()

EshelbianPlasticity::OpCalculateRotationAndSpatialGradient::OpCalculateRotationAndSpatialGradient ( boost::shared_ptr< DataAtIntegrationPts > &  data_ptr)
inline

Definition at line 491 of file EshelbianPlasticity.hpp.

493  : VolUserDataOperator(NOSPACE, OPSPACE), dataAtPts(data_ptr) {}

Member Function Documentation

◆ doWork()

MoFEMErrorCode EshelbianPlasticity::OpCalculateRotationAndSpatialGradient::doWork ( int  side,
EntityType  type,
EntData &  data 
)

Definition at line 257 of file EshelbianOperators.cpp.

259  {
260  MoFEMFunctionBegin;
261 
262  FTensor::Index<'L', size_symm> L;
263  auto t_L = symm_L_tensor();
264 
265  int nb_integration_pts = getGaussPts().size2();
266  FTensor::Index<'i', 3> i;
267  FTensor::Index<'j', 3> j;
268  FTensor::Index<'k', 3> k;
269  FTensor::Index<'l', 3> l;
270  FTensor::Index<'m', 3> m;
271  FTensor::Index<'l', 3> n;
272 
273  dataAtPts->hAtPts.resize(9, nb_integration_pts, false);
274  dataAtPts->hdOmegaAtPts.resize(9 * 3, nb_integration_pts, false);
275  dataAtPts->hdLogStretchAtPts.resize(9 * 6, nb_integration_pts, false);
276  dataAtPts->leviKirchhoffAtPts.resize(3, nb_integration_pts, false);
277  dataAtPts->leviKirchhoffPAtPts.resize(9 * 3, nb_integration_pts, false);
278  dataAtPts->leviKirchhoffOmegaAtPts.resize(9, nb_integration_pts, false);
279 
280  dataAtPts->adjointPdstretchAtPts.resize(9, nb_integration_pts, false);
281  dataAtPts->adjointPdUAtPts.resize(size_symm, nb_integration_pts, false);
282  dataAtPts->adjointPdUdPAtPts.resize(9 * size_symm, nb_integration_pts, false);
283  dataAtPts->adjointPdUdOmegaAtPts.resize(3 * size_symm, nb_integration_pts,
284  false);
285  dataAtPts->rotMatAtPts.resize(9, nb_integration_pts, false);
286  dataAtPts->diffRotMatAtPts.resize(27, nb_integration_pts, false);
287  dataAtPts->stretchTensorAtPts.resize(6, nb_integration_pts, false);
288  dataAtPts->diffStretchTensorAtPts.resize(36, nb_integration_pts, false);
289  dataAtPts->eigenVals.resize(3, nb_integration_pts, false);
290  dataAtPts->eigenVecs.resize(9, nb_integration_pts, false);
291  dataAtPts->nbUniq.resize(nb_integration_pts, false);
292 
293  dataAtPts->logStretchTotalTensorAtPts.resize(6, nb_integration_pts, false);
294 
295  auto t_h = getFTensor2FromMat<3, 3>(dataAtPts->hAtPts);
296  auto t_h_domega = getFTensor3FromMat<3, 3, 3>(dataAtPts->hdOmegaAtPts);
297  auto t_h_dlog_u =
298  getFTensor3FromMat<3, 3, size_symm>(dataAtPts->hdLogStretchAtPts);
299  auto t_levi_kirchoff = getFTensor1FromMat<3>(dataAtPts->leviKirchhoffAtPts);
300  auto t_levi_kirchoff_dP =
301  getFTensor3FromMat<3, 3, 3>(dataAtPts->leviKirchhoffPAtPts);
302  auto t_levi_kirchoff_domega =
303  getFTensor2FromMat<3, 3>(dataAtPts->leviKirchhoffOmegaAtPts);
304  auto t_approx_P_adjont_dstretch =
305  getFTensor2FromMat<3, 3>(dataAtPts->adjointPdstretchAtPts);
306  auto t_approx_P_adjont_log_du =
307  getFTensor1FromMat<size_symm>(dataAtPts->adjointPdUAtPts);
308  auto t_approx_P_adjont_log_du_dP =
309  getFTensor3FromMat<3, 3, size_symm>(dataAtPts->adjointPdUdPAtPts);
310  auto t_approx_P_adjont_log_du_domega =
311  getFTensor2FromMat<3, size_symm>(dataAtPts->adjointPdUdOmegaAtPts);
312  auto t_omega = getFTensor1FromMat<3>(dataAtPts->rotAxisAtPts);
313  auto t_R = getFTensor2FromMat<3, 3>(dataAtPts->rotMatAtPts);
314  auto t_diff_R = getFTensor3FromMat<3, 3, 3>(dataAtPts->diffRotMatAtPts);
315  auto t_log_u =
316  getFTensor2SymmetricFromMat<3>(dataAtPts->logStretchTensorAtPts);
317  auto t_u = getFTensor2SymmetricFromMat<3>(dataAtPts->stretchTensorAtPts);
318  auto t_approx_P = getFTensor2FromMat<3, 3>(dataAtPts->approxPAtPts);
319  auto t_diff_u =
320  getFTensor4DdgFromMat<3, 3, 1>(dataAtPts->diffStretchTensorAtPts);
321  auto t_eigen_vals = getFTensor1FromMat<3>(dataAtPts->eigenVals);
322  auto t_eigen_vecs = getFTensor2FromMat<3, 3>(dataAtPts->eigenVecs);
323 
324  auto &nbUniq = dataAtPts->nbUniq;
325  auto t_grad_h1 = getFTensor2FromMat<3, 3>(dataAtPts->wGradH1AtPts);
326  auto t_log_stretch_total =
327  getFTensor2SymmetricFromMat<3>(dataAtPts->logStretchTotalTensorAtPts);
328  auto t_kd = FTensor::Kronecker_Delta<int>();
329 
330  for (int gg = 0; gg != nb_integration_pts; ++gg) {
331 
332  FTensor::Tensor1<double, 3> eig;
333  FTensor::Tensor2<double, 3, 3> eigen_vec;
334  FTensor::Tensor3<double, 3, 3, size_symm> t_Ldiff_u;
335  FTensor::Tensor2<double, 3, 3> t_h1;
336  FTensor::Tensor2<double, 3, 3> t_Ru;
337  FTensor::Tensor2<double, 3, 3> t_approx_P_intermidiate;
338 
339  t_h1(i, j) = t_grad_h1(i, j) + t_kd(i, j);
340 
341  auto calulate_rotation = [&]() {
342  auto t0_diff =
343  get_diff_rotation_form_vector(t_omega, EshelbianCore::rotSelector);
344  auto t0 = get_rotation_form_vector(t_omega, EshelbianCore::rotSelector);
345  t_diff_R(i, j, k) = t0_diff(i, j, k);
346  t_R(i, j) = t0(i, j);
347  };
348 
349  auto calulate_streach = [&]() {
350  eigen_vec(i, j) = t_log_u(i, j);
351  CHKERR computeEigenValuesSymmetric(eigen_vec, eig);
352  // rare case when two eigen values are equal
353  nbUniq[gg] = get_uniq_nb<3>(&eig(0));
354  if (nbUniq[gg] < 3) {
355  sort_eigen_vals<3>(eig, eigen_vec);
356  }
357  t_eigen_vals(i) = eig(i);
358  t_eigen_vecs(i, j) = eigen_vec(i, j);
359  auto t_u_tmp =
360  EigenMatrix::getMat(t_eigen_vals, t_eigen_vecs, EshelbianCore::f);
361  t_u(i, j) = t_u_tmp(i, j);
362  auto t_diff_u_tmp =
363  EigenMatrix::getDiffMat(t_eigen_vals, t_eigen_vecs, EshelbianCore::f,
364  EshelbianCore::d_f, nbUniq[gg]);
365  t_diff_u(i, j, k, l) = t_diff_u_tmp(i, j, k, l);
366  t_Ldiff_u(i, j, L) = t_diff_u(i, j, m, n) * t_L(m, n, L);
367  };
368 
369  calulate_rotation();
370  calulate_streach();
371 
372  switch (EshelbianCore::gradApperoximator) {
373  case LARGE_ROT:
374 
375  // rotation can be large, moderate or small
376  t_Ru(i, m) = t_R(i, l) * t_u(l, m);
377  t_h(i, j) = t_Ru(i, m) * t_h1(m, j);
378  t_h_domega(i, j, k) = (t_diff_R(i, l, k) * t_u(l, m)) * t_h1(m, j);
379  t_h_dlog_u(i, j, L) = (t_R(i, l) * t_Ldiff_u(l, m, L)) * t_h1(m, j);
380 
381  // conservation of angular momentum at current configuration
382  t_approx_P_intermidiate(i, m) = t_approx_P(i, j) * t_h1(m, j);
383  t_approx_P_adjont_dstretch(l, m) =
384  t_approx_P_intermidiate(i, m) * t_R(i, l);
385 
386  t_levi_kirchoff(k) =
387  levi_civita(l, m, k) * (t_approx_P_adjont_dstretch(l, m));
388  t_levi_kirchoff_dP(i, j, k) =
389  (levi_civita(l, m, k) * t_R(i, l)) * t_h1(m, j);
390  t_levi_kirchoff_domega(k, n) =
391  levi_civita(l, m, k) *
392  (t_approx_P_intermidiate(i, m) * t_diff_R(i, l, n));
393 
394  t_approx_P_adjont_log_du(L) =
395  t_Ldiff_u(l, m, L) * t_approx_P_adjont_dstretch(l, m);
396  t_approx_P_adjont_log_du_dP(i, j, L) = t_h_dlog_u(i, j, L);
397  t_approx_P_adjont_log_du_domega(n, L) =
398  t_Ldiff_u(l, m, L) *
399  (t_approx_P_intermidiate(i, m) * t_diff_R(i, l, n));
400 
401  break;
402  case MODERATE_ROT:
403 
404  {
405 
406  // assume small current rotation
407  FTensor::Tensor2<double, 3, 3> t_Omega;
408  t_Omega(i, j) = FTensor::levi_civita<double>(i, j, k) * t_omega(k);
409  t_Ru(i, m) = t_Omega(i, m) + t_u(i, m);
410  t_h(i, j) = t_Ru(i, m) * t_h1(m, j);
411  t_h_domega(i, j, k) = FTensor::levi_civita<double>(i, m, k) * t_h1(m, j);
412  t_h_dlog_u(i, j, L) = t_Ldiff_u(i, m, L) * t_h1(m, j);
413 
414  // conservation of angular momentum at intermediate configuration
415  t_approx_P_intermidiate(i, m) = t_approx_P(i, j) * t_h1(m, j);
416  t_approx_P_adjont_dstretch(i, m) = t_approx_P_intermidiate(i, m);
417 
418  t_levi_kirchoff(k) =
419  levi_civita(i, m, k) * (t_approx_P_adjont_dstretch(i, m));
420  t_levi_kirchoff_dP(i, j, k) = levi_civita(i, m, k) * t_h1(m, j);
421  t_levi_kirchoff_domega(k, n) = 0;
422 
423  t_approx_P_adjont_log_du(L) =
424  t_Ldiff_u(i, m, L) * t_approx_P_adjont_dstretch(i, m);
425  t_approx_P_adjont_log_du_dP(i, j, L) = t_h_dlog_u(i, j, L);
426  t_approx_P_adjont_log_du_domega(n, L) = 0;
427 
428  }
429 
430  break;
431 
432  case SMALL_ROT:
433 
434  {
435 
436  // assume small current rotation
437  FTensor::Tensor2<double, 3, 3> t_Omega;
438  t_Omega(i, j) = FTensor::levi_civita<double>(i, j, k) * t_omega(k);
439  t_h(i, j) = t_Omega(i, j) + t_u(i, j);
440  t_h_domega(i, j, k) = FTensor::levi_civita<double>(i, j, k);
441  t_h_dlog_u(i, j, L) = t_Ldiff_u(i, j, L);
442 
443  // conservation of angular momentum at reference condition
444  t_levi_kirchoff(k) = levi_civita(i, j, k) * t_approx_P(i, j);
445  t_levi_kirchoff_dP(i, j, k) = levi_civita(i, j, k);
446  t_levi_kirchoff_domega(k, l) = 0;
447 
448  t_approx_P_adjont_dstretch(i, j) = t_approx_P(i, j);
449  t_approx_P_adjont_log_du(L) =
450  t_Ldiff_u(i, j, L) * t_approx_P_adjont_dstretch(i, j);
451  t_approx_P_adjont_log_du_dP(i, j, L) = t_h_dlog_u(i, j, L);
452  t_approx_P_adjont_log_du_domega(m, L) = 0;
453 
454  }
455 
456  break;
457  }
458 
459  FTensor::Tensor2<double, 3, 3> t_C_h1;
460  t_C_h1(i, j) = t_h1(k, i) * t_h1(k, j);
461  eigen_vec(i, j) = t_C_h1(i, j);
462  CHKERR computeEigenValuesSymmetric(eigen_vec, eig);
463  switch (EshelbianCore::stretchSelector) {
464  case StretchSelector::LINEAR:
465  break;
466  case StretchSelector::LOG:
467  for (int ii = 0; ii != 3; ++ii) {
468  eig(ii) = std::max(eig(ii),
469  2 * EshelbianCore::inv_f(EshelbianCore::f(
470  std::numeric_limits<double>::min_exponent)));
471  }
472  break;
473  }
474 
475  auto t_log_u2_h1_tmp =
476  EigenMatrix::getMat(eig, eigen_vec, EshelbianCore::inv_f);
477 
478  switch (EshelbianCore::gradApperoximator) {
479  case LARGE_ROT:
480  case MODERATE_ROT:
481  t_log_stretch_total(i, j) = t_log_u2_h1_tmp(i, j) / 2 + t_log_u(i, j);
482  break;
483  case SMALL_ROT:
484  t_log_stretch_total(i, j) = t_log_u(i, j);
485  break;
486  };
487 
488  ++t_h;
489  ++t_h_domega;
490  ++t_h_dlog_u;
491  ++t_levi_kirchoff;
492  ++t_levi_kirchoff_dP;
493  ++t_levi_kirchoff_domega;
494  ++t_approx_P_adjont_dstretch;
495  ++t_approx_P_adjont_log_du;
496  ++t_approx_P_adjont_log_du_dP;
497  ++t_approx_P_adjont_log_du_domega;
498  ++t_approx_P;
499  ++t_R;
500  ++t_diff_R;
501  ++t_log_u;
502  ++t_u;
503  ++t_diff_u;
504  ++t_eigen_vals;
505  ++t_eigen_vecs;
506  ++t_omega;
507  ++t_grad_h1;
508  ++t_log_stretch_total;
509  }
510 
511  MoFEMFunctionReturn(0);
512 }

Member Data Documentation

◆ dataAtPts

boost::shared_ptr<DataAtIntegrationPts> EshelbianPlasticity::OpCalculateRotationAndSpatialGradient::dataAtPts

data at integration pts

Definition at line 490 of file EshelbianPlasticity.hpp.

The documentation for this struct was generated from the following files:
NOSPACE
@ NOSPACE
Definition: definitions.h:83
EigenMatrix::getMat
FTensor::Tensor2_symmetric< double, 3 > getMat(Val< double, 3 > &t_val, Vec< double, 3 > &t_vec, Fun< double > f)
Get the Mat object.
Definition: MatrixFunction.cpp:53
EshelbianPlasticity::get_rotation_form_vector
auto get_rotation_form_vector(FTensor::Tensor1< T, 3 > &t_omega, RotSelector rotSelector=LARGE_ROT)
Definition: EshelbianOperators.cpp:83
EshelbianPlasticity::LINEAR
@ LINEAR
Definition: EshelbianPlasticity.hpp:21
EshelbianPlasticity::size_symm
constexpr static auto size_symm
Definition: EshelbianOperators.cpp:47
FTensor::Tensor1< double, 3 >
EshelbianPlasticity::OpCalculateRotationAndSpatialGradient::dataAtPts
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
Definition: EshelbianPlasticity.hpp:490
FTensor::Kronecker_Delta
Kronecker Delta class.
Definition: Kronecker_Delta.hpp:15
EshelbianPlasticity::EshelbianCore::d_f
static boost::function< double(const double)> d_f
Definition: EshelbianPlasticity.hpp:870
FTensor::levi_civita
constexpr std::enable_if<(Dim0<=2 &&Dim1<=2), Tensor2_Expr< Levi_Civita< T >, T, Dim0, Dim1, i, j > >::type levi_civita(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
levi_civita functions to make for easy adhoc use
Definition: Levi_Civita.hpp:617
FTensor::Tensor2< double, 3, 3 >
EshelbianPlasticity::EshelbianCore::gradApperoximator
static enum RotSelector gradApperoximator
Definition: EshelbianPlasticity.hpp:865
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:535
FTensor::Tensor3
Definition: Tensor3_value.hpp:12
EshelbianPlasticity::EshelbianCore::stretchSelector
static enum StretchSelector stretchSelector
Definition: EshelbianPlasticity.hpp:866
EshelbianPlasticity::EshelbianCore::rotSelector
static enum RotSelector rotSelector
Definition: EshelbianPlasticity.hpp:864
EshelbianPlasticity::EshelbianCore::f
static boost::function< double(const double)> f
Definition: EshelbianPlasticity.hpp:869
EshelbianPlasticity::VolUserDataOperator
VolumeElementForcesAndSourcesCore::UserDataOperator VolUserDataOperator
Definition: EshelbianPlasticity.hpp:28
MoFEM::L
VectorDouble L
Definition: Projection10NodeCoordsOnField.cpp:124
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
t_kd
constexpr auto t_kd
Definition: free_surface.cpp:137
EshelbianPlasticity::symm_L_tensor
auto symm_L_tensor()
Definition: EshelbianOperators.cpp:60
FTensor::Index
Definition: Index.hpp:23
convert.n
n
Definition: convert.py:82
EshelbianPlasticity::EshelbianCore::inv_f
static boost::function< double(const double)> inv_f
Definition: EshelbianPlasticity.hpp:872
EshelbianPlasticity::SMALL_ROT
@ SMALL_ROT
Definition: EshelbianPlasticity.hpp:20
EshelbianPlasticity::LARGE_ROT
@ LARGE_ROT
Definition: EshelbianPlasticity.hpp:20
EshelbianPlasticity::LOG
@ LOG
Definition: EshelbianPlasticity.hpp:21
j
FTensor::Index< 'j', 3 > j
Definition: matrix_function.cpp:19
EshelbianPlasticity::get_diff_rotation_form_vector
auto get_diff_rotation_form_vector(FTensor::Tensor1< T, 3 > &t_omega, RotSelector rotSelector=LARGE_ROT)
Definition: EshelbianOperators.cpp:119
m
FTensor::Index< 'm', 3 > m
Definition: shallow_wave.cpp:80
EshelbianPlasticity::MODERATE_ROT
@ MODERATE_ROT
Definition: EshelbianPlasticity.hpp:20
k
FTensor::Index< 'k', 3 > k
Definition: matrix_function.cpp:20
EigenMatrix::getDiffMat
FTensor::Ddg< double, 3, 3 > getDiffMat(Val< double, 3 > &t_val, Vec< double, 3 > &t_vec, Fun< double > f, Fun< double > d_f, const int nb)
Get the Diff Mat object.
Definition: MatrixFunction.cpp:64
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
l
FTensor::Index< 'l', 3 > l
Definition: matrix_function.cpp:21
MoFEM::computeEigenValuesSymmetric
MoFEMErrorCode computeEigenValuesSymmetric(const MatrixDouble &mat, VectorDouble &eig, MatrixDouble &eigen_vec)
compute eigenvalues of a symmetric matrix using lapack dsyev
Definition: Templates.hpp:1423
Generated by Doxygen 1.8.17 and hosted at University of Glasgow