EshelbianPlasticity::OpTopologicalObjectivePythonImpl Struct Reference

Inheritance diagram for EshelbianPlasticity::OpTopologicalObjectivePythonImpl:

Collaboration diagram for EshelbianPlasticity::OpTopologicalObjectivePythonImpl:

Public Types
using	OP = ForcesAndSourcesCore::UserDataOperator

Public Member Functions
	OpTopologicalObjectivePythonImpl (boost::shared_ptr< DataAtIntegrationPts > data_at_pts_ptr, boost::shared_ptr< TopologicalData > topo_p, boost::shared_ptr< ObjectiveFunctionData > python_ptr, const ObjectiveModelType eval_energy_model=PYTHON_MODEL)
MoFEMErrorCode	doWork (int side, EntityType type, EntData &data)

Private Attributes
ObjectiveModelType	evalEnergyModel = PYTHON_MODEL
boost::shared_ptr< DataAtIntegrationPts >	dataAtPts
boost::shared_ptr< TopologicalData >	topoData
boost::shared_ptr< ObjectiveFunctionData >	pythonPtr

Detailed Description

Definition at line 333 of file EshelbianTopologicalDerivativeOperators.cpp.

Member Typedef Documentation

OP

using EshelbianPlasticity::OpTopologicalObjectivePythonImpl::OP = ForcesAndSourcesCore::UserDataOperator

Definition at line 335 of file EshelbianTopologicalDerivativeOperators.cpp.

Constructor & Destructor Documentation

OpTopologicalObjectivePythonImpl()

EshelbianPlasticity::OpTopologicalObjectivePythonImpl::OpTopologicalObjectivePythonImpl ( boost::shared_ptr< DataAtIntegrationPts >  data_at_pts_ptr, boost::shared_ptr< TopologicalData >  topo_p, boost::shared_ptr< ObjectiveFunctionData >  python_ptr, const ObjectiveModelType  eval_energy_model = PYTHON_MODEL  )

inline

Definition at line 337 of file EshelbianTopologicalDerivativeOperators.cpp.

      : OP(NOSPACE, OP::OPSPACE), dataAtPts(data_at_pts_ptr), topoData(topo_p),
        pythonPtr(python_ptr), evalEnergyModel(eval_energy_model) {}

Member Function Documentation

doWork()

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

Definition at line 354 of file EshelbianTopologicalDerivativeOperators.cpp.

                                                                       {
  MoFEMFunctionBegin;
 
#ifndef NDEBUG
  if (!dataAtPts)
    SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
            "DataAtIntegrationPts pointer is null");
  if (!topoData)
    SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
            "Topological data pointer is null");
#endif // NDEBUG
 
  const int nb_gauss_pts = getGaussPts().size2();
  if (!nb_gauss_pts)
    MoFEMFunctionReturnHot(0);
 
  MatrixSizeHelper<GetFTensor1FromMatType<3, -1, DL>, DL>::size(
      *topoData->objConstantRotation, nb_gauss_pts);
  topoData->objConstantRotation->clear();
  const auto &features = dataAtPts->physicsPtr->getFeatures();
  const auto variable_compliance_mask =
      PhysicalEquations::Features{}
          .set(PhysicalEquations::NON_HOMOGENEOUS_MATERIAL)
          .set(PhysicalEquations::NO_STRETCH_NONLINEAR);
 
  auto stress_full_ptr = boost::make_shared<MatrixDouble>();
  auto get_stress_full =
      MatrixSizeHelper<GetFTensor2FromMatType<SPACE_DIM, SPACE_DIM, -1, DL>,
                       DL>::size(*stress_full_ptr, nb_gauss_pts);
  stress_full_ptr->clear();
  auto strain_full_ptr = boost::make_shared<MatrixDouble>();
  auto get_strain_full =
      MatrixSizeHelper<GetFTensor2FromMatType<SPACE_DIM, SPACE_DIM, -1, DL>,
                       DL>::size(*strain_full_ptr, nb_gauss_pts);
  strain_full_ptr->clear();
 
  auto t_stress = get_stress_full();
  auto t_strain = get_strain_full();
 
  auto t_biot = dataAtPts->getFTensorAdjointPdstretch(nb_gauss_pts);
  auto t_u = dataAtPts->getFTensorStretch(nb_gauss_pts);
 
  auto next = [&]() {
    ++t_stress;
    ++t_strain;
    ++t_biot;
    ++t_u;
  };
 
  FTENSOR_INDEXES(SPACE_DIM, i, j, k, l);
 
  for (auto gg = 0; gg != nb_gauss_pts; ++gg) {
    // we have to handle all variants, that will render how the Jacobian
    // gradient is evaluated.
    t_stress(i, j) = t_biot(i, j);
    t_strain(i, j) = t_u(i, j);
    next();
  }
 
  auto evaluate_python_objective = [&]() {
    MoFEMFunctionBegin;
    if (!pythonPtr)
      SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
              "ObjectiveFunctionData pointer is null");
 
    auto &coords = OP::getCoordsAtGaussPts();
    CHKERR pythonPtr->evalInteriorObjectiveFunction(
        coords, dataAtPts->getSmallWL2AtPts(), stress_full_ptr, strain_full_ptr,
        topoData->getObjAtPts(), false);
    CHKERR pythonPtr->evalInteriorObjectiveGradientStrain(
        coords, dataAtPts->getSmallWL2AtPts(), stress_full_ptr, strain_full_ptr,
        topoData->getObjDStrainAtPts(), false);
    CHKERR pythonPtr->evalInteriorObjectiveGradientU(
        coords, dataAtPts->getSmallWL2AtPts(), stress_full_ptr, strain_full_ptr,
        topoData->getObjDDisplacementAtPts(), false);
    CHKERR pythonPtr->evalInteriorObjectiveGradientStress(
        coords, dataAtPts->getSmallWL2AtPts(), stress_full_ptr, strain_full_ptr,
        topoData->getObjDStressAtPts(), false);
    MoFEMFunctionReturn(0);
  };
 
  auto evaluate_energy_of_hencky_model = [&]() {
    MoFEMFunctionBegin;
 
    MatrixSizeHelper<GetFTensor1FromMatType<SPACE_DIM, -1, DL>, DL>::size(
        *topoData->getObjDDisplacementAtPts(), nb_gauss_pts);
    topoData->getObjDDisplacementAtPts()->clear();
    MatrixSizeHelper<GetFTensor2FromMatType<SPACE_DIM, SPACE_DIM, -1, DL>,
                     DL>::size(*topoData->getObjDStressAtPts(), nb_gauss_pts);
    topoData->getObjDStressAtPts()->clear();
    MatrixSizeHelper<GetFTensor1FromMatType<SPACE_DIM, -1, DL>, DL>::size(
        *topoData->getObjDRotationAtPts(), nb_gauss_pts);
 
    auto eval_evergy = [&](auto &&t_D) {
      auto get_obj =
          MatrixSizeHelper<GetFTensor1FromMatType<1, -1, DL>, DL>::size(
              *topoData->getObjAtPts(), nb_gauss_pts);
      auto get_dstrain_obj =
          MatrixSizeHelper<GetFTensor2FromMatType<SPACE_DIM, SPACE_DIM, -1, DL>,
                           DL>::size(*topoData->getObjDStrainAtPts(),
                                     nb_gauss_pts);
      auto t_obj = get_obj();
      auto t_dstrain_obj = get_dstrain_obj();
      auto t_log_u = dataAtPts->getFTensorLogStretchTotal(nb_gauss_pts);
      for (int gg = 0; gg != nb_gauss_pts; ++gg) {
        t_obj(0) = 0.5 * (t_log_u(i, j) * t_D(i, j, k, l) * t_log_u(k, l));
        t_dstrain_obj(i, j) = t_D(i, j, k, l) * t_log_u(k, l);
        ++t_log_u;
        ++t_obj;
        ++t_dstrain_obj;
        ++t_D;
      }
    };
 
    if (!features.test(PhysicalEquations::NON_HOMOGENEOUS_MATERIAL)) {
      eval_evergy(
          getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, 0>(dataAtPts->matD));
    } else {
      eval_evergy(getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM>(dataAtPts->matD));
    }
 
    topoData->getObjDRotationAtPts()->clear();
    MoFEMFunctionReturn(0);
  };
 
  auto evaluate_energy_of_hencky_model_nostreach = [&]() {
    MoFEMFunctionBegin;
 
    MatrixSizeHelper<GetFTensor1FromMatType<SPACE_DIM, -1, DL>, DL>::size(
        *topoData->getObjDDisplacementAtPts(), nb_gauss_pts);
    topoData->getObjDDisplacementAtPts()->clear();
    MatrixSizeHelper<GetFTensor2FromMatType<SPACE_DIM, SPACE_DIM, -1, DL>,
                     DL>::size(*topoData->getObjDStrainAtPts(), nb_gauss_pts);
    topoData->getObjDStrainAtPts()->clear();
    MatrixSizeHelper<GetFTensor1FromMatType<SPACE_DIM, -1, DL>, DL>::size(
        *topoData->getObjDRotationAtPts(), nb_gauss_pts);
    topoData->getObjDRotationAtPts()->clear();
 
    auto eval_evergy = [&](auto &&t_inv_D) {
      auto get_obj =
          MatrixSizeHelper<GetFTensor1FromMatType<1, -1, DL>, DL>::size(
              *topoData->getObjAtPts(), nb_gauss_pts);
      auto get_dstress_obj =
          MatrixSizeHelper<GetFTensor2FromMatType<SPACE_DIM, SPACE_DIM, -1, DL>,
                           DL>::size(*topoData->getObjDStressAtPts(),
                                     nb_gauss_pts);
      auto t_obj = get_obj();
      auto t_dstress_obj = get_dstress_obj();
      auto t_stress = dataAtPts->getFTensorAdjointPdstretch(nb_gauss_pts);
      for (int gg = 0; gg != nb_gauss_pts; ++gg) {
        t_obj(0) =
            0.5 * (t_stress(i, j) * t_inv_D(i, j, k, l) * t_stress(k, l));
        t_dstress_obj(i, j) = t_inv_D(i, j, k, l) * t_stress(k, l);
        ++t_stress;
        ++t_obj;
        ++t_dstress_obj;
        ++t_inv_D;
      }
    };
 
    if ((features & variable_compliance_mask).none()) {
      eval_evergy(
          getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, 0>(dataAtPts->matInvD));
    } else {
      eval_evergy(
          getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM>(dataAtPts->matInvD));
    }
    MoFEMFunctionReturn(0);
  };
 
  auto conversion_of_biot_stress = [&]() {
    MoFEMFunctionBegin;
    // Pull Python's Biot-stress covector back to Piola and both rotation
    // factors.
    MatrixSizeHelper<GetFTensor1FromMatType<3, -1, DL>, DL>::size(
        *topoData->getObjDRotationAtPts(), nb_gauss_pts);
    topoData->getObjDRotationAtPts()->clear();
    auto objective = topoData->getFTensorObjDStress(nb_gauss_pts);
    auto R = dataAtPts->getFTensorRotMat(nb_gauss_pts);
    auto P = dataAtPts->getFTensorApproxP(nb_gauss_pts);
    auto gradient = dataAtPts->getFTensorSmallWGradH1(nb_gauss_pts);
    FTENSOR_INDEXES(3, i, k, l, m, o);
    constexpr auto identity = FTensor::Kronecker_Delta<double>();
    for (int gg = 0; gg != nb_gauss_pts; ++gg) {
      FTensor::Tensor2<double, 3, 3> covector, h1;
      covector(l, o) = objective(l, o);
      h1(o, k) = identity(o, k);
      if (EshelbianCore::gradApproximator != SMALL_ROT) {
        if (EshelbianCore::gradApproximator != NO_H1_CONFIGURATION)
          h1(o, k) += gradient(o, k);
        objective(i, k) = R(i, l) * (covector(l, o) * h1(o, k));
        for (int mode = 0; mode != 2; ++mode) {
          FTensor::Tensor3<double, 3, 3, 3> derivative;
          if (EshelbianCore::rotSelector == LARGE_ROT) {
            using namespace LogarithmicDirectionalDerivatives;
            const auto read = [&](MatrixPtr values) {
              Vector result;
              auto view = GetFTensor1FromMatImpl<3, -1, DL, MatrixDouble>::get(
                  *values, gg, 0);
              result(i) = view(i);
              return result;
            };
            const SplitRotation rotation(read(dataAtPts->constantRotation),
                                         read(dataAtPts->higherRotation));
            derivative(i, l, m) = rotation.firstTensor(mode)(i, l, m);
          } else {
            derivative(i, l, m) = levi_civita(i, l, m);
          }
          auto rotation_covector =
              topoData->getFTensorObjDRotation(nb_gauss_pts, mode, gg);
          rotation_covector(m) =
              covector(l, o) * (derivative(i, l, m) * P(i, k)) * h1(o, k);
        }
      }
      ++objective;
      ++R;
      ++P;
      ++gradient;
    }
    MoFEMFunctionReturn(0);
  };
 
  auto conversion_of_stretch = [&]() {
    MoFEMFunctionBegin;
    // Python sees the physical stretch tensor. The stretch field stores log
    // stretch, so convert dJ/dU to dJ/dlogU once here.
    auto t_obj_dstretch = topoData->getFTensorObjDStrain(nb_gauss_pts);
    auto t_diff_stretch = dataAtPts->getFTensorDiffStretch(nb_gauss_pts);
 
    FTENSOR_INDEX(SPACE_DIM, i);
    FTENSOR_INDEX(SPACE_DIM, j);
    FTENSOR_INDEX(SPACE_DIM, k);
    FTENSOR_INDEX(SPACE_DIM, l);
 
    for (auto gg = 0; gg != nb_gauss_pts; ++gg) {
      FTensor::Tensor2<double, SPACE_DIM, SPACE_DIM> t_dJ_dstretch;
      t_dJ_dstretch(i, j) = t_obj_dstretch(i, j);
 
      t_obj_dstretch(k, l) = t_dJ_dstretch(i, j) * t_diff_stretch(i, j, k, l);
 
      ++t_obj_dstretch;
      ++t_diff_stretch;
    }
 
    MoFEMFunctionReturn(0);
  };
 
  auto conversion_of_stretch_to_stress_for_no_stretch = [&](auto t_inv_D) {
    MoFEMFunctionBegin;
    // In no-stretch mode, log stretch is computed from the stress field. Fold
    // the objective stretch derivative into the stress derivative buffer.
    auto t_obj_dstress = topoData->getFTensorObjDStress(nb_gauss_pts);
    auto t_obj_dstretch = topoData->getFTensorObjDStrain(nb_gauss_pts);
 
    FTENSOR_INDEX(SPACE_DIM, i);
    FTENSOR_INDEX(SPACE_DIM, j);
    FTENSOR_INDEX(SPACE_DIM, k);
    FTENSOR_INDEX(SPACE_DIM, l);
 
    for (auto gg = 0; gg != nb_gauss_pts; ++gg) {
      FTensor::Tensor2_symmetric<double, SPACE_DIM> t_dstretch_dstress;
      t_dstretch_dstress(i, j) =
          ((t_obj_dstretch(k, l) || t_obj_dstretch(l, k)) / 2.) *
          t_inv_D(k, l, i, j);
 
      t_obj_dstress(i, j) += t_dstretch_dstress(i, j);
 
      ++t_obj_dstress;
      ++t_obj_dstretch;
      ++t_inv_D;
    }
 
    MoFEMFunctionReturn(0);
  };
 
  switch (evalEnergyModel) {
  case PYTHON_MODEL:
    CHKERR evaluate_python_objective();
    CHKERR conversion_of_biot_stress();
    if ((features & PhysicalEquations::noStretchMask).any()) {
      if ((features & variable_compliance_mask).none()) {
        CHKERR conversion_of_stretch_to_stress_for_no_stretch(
            getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, 0>(dataAtPts->matInvD));
      } else {
        CHKERR conversion_of_stretch_to_stress_for_no_stretch(
            getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM>(dataAtPts->matInvD));
      }
    } else {
      CHKERR conversion_of_stretch();
    }
    break;
  case HENCKY_MODEL:
    if ((features & PhysicalEquations::noStretchMask).any()) {
      CHKERR evaluate_energy_of_hencky_model_nostreach();
      CHKERR conversion_of_biot_stress();
    } else {
      CHKERR evaluate_energy_of_hencky_model();
    }
    break;
  default:
    SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
            "Objective model type not handled");
  }
 
  MoFEMFunctionReturn(0);
}

Member Data Documentation

dataAtPts

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

private

Definition at line 349 of file EshelbianTopologicalDerivativeOperators.cpp.

evalEnergyModel

ObjectiveModelType EshelbianPlasticity::OpTopologicalObjectivePythonImpl::evalEnergyModel = PYTHON_MODEL

private

Definition at line 348 of file EshelbianTopologicalDerivativeOperators.cpp.

pythonPtr

boost::shared_ptr<ObjectiveFunctionData> EshelbianPlasticity::OpTopologicalObjectivePythonImpl::pythonPtr

private

Definition at line 351 of file EshelbianTopologicalDerivativeOperators.cpp.

topoData

boost::shared_ptr<TopologicalData> EshelbianPlasticity::OpTopologicalObjectivePythonImpl::topoData

private

Definition at line 350 of file EshelbianTopologicalDerivativeOperators.cpp.

---

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

• users_modules/eshelbian_plasticity/src/impl/EshelbianTopologicalDerivativeOperators.cpp