OGS
MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars > Class Template Reference

Detailed Description

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
class MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >

Uses a material model provided by MFront (via MFront's generic interface and the MGIS library).

Definition at line 247 of file MFrontGeneric.h.

#include <MFrontGeneric.h>

Collaboration diagram for MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >:
[legend]

Public Types

using KelvinVector
using KelvinMatrix
using InternalVariable

Public Member Functions

 MFrontGeneric (mgis::behaviour::Behaviour &&behaviour, std::vector< ParameterLib::Parameter< double > const * > &&material_properties, std::map< std::string, ParameterLib::Parameter< double > const * > &&state_variables_initial_properties, std::optional< ParameterLib::CoordinateSystem > const &local_coordinate_system)
std::unique_ptr< typename MechanicsBase< DisplacementDim >::MaterialStateVariables > createMaterialStateVariables () const
void initializeInternalStateVariables (double const t, ParameterLib::SpatialPosition const &x, typename MechanicsBase< DisplacementDim >::MaterialStateVariables &material_state_variables) const
std::optional< std::tuple< OGSMFrontThermodynamicForcesData, std::unique_ptr< typename MechanicsBase< DisplacementDim >::MaterialStateVariables >, OGSMFrontTangentOperatorData > > integrateStress (MaterialPropertyLib::VariableArray const &variable_array_prev, MaterialPropertyLib::VariableArray const &variable_array, double const t, ParameterLib::SpatialPosition const &x, double const dt, typename MechanicsBase< DisplacementDim >::MaterialStateVariables const &material_state_variables) const
std::vector< InternalVariablegetInternalVariables () const
template<typename ForcesGradsCombinations = typename ForcesGradsCombinations<Gradients, TDynForces, ExtStateVars>::type>
OGSMFrontTangentOperatorBlocksView< DisplacementDim, ForcesGradsCombinationscreateTangentOperatorBlocksView () const
OGSMFrontThermodynamicForcesView< DisplacementDim, TDynForces > createThermodynamicForcesView () const
double getBulkModulus (double const, ParameterLib::SpatialPosition const &, KelvinMatrix const *const C) const
double computeFreeEnergyDensity (double const, ParameterLib::SpatialPosition const &, double const, KelvinVector const &, KelvinVector const &, typename MechanicsBase< DisplacementDim >::MaterialStateVariables const &) const

Private Types

using GradientsAndExtStateVars

Private Attributes

mgis::behaviour::Behaviour _behaviour
int const equivalent_plastic_strain_offset_
std::vector< ParameterLib::Parameter< double > const * > _material_properties
std::map< std::string, ParameterLib::Parameter< double > const * > _state_variables_initial_properties
ParameterLib::CoordinateSystem const *const _local_coordinate_system

Member Typedef Documentation

◆ GradientsAndExtStateVars

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
using MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::GradientsAndExtStateVars
private
Initial value:
boost::mp11::mp_append<Gradients, ExtStateVars>

Definition at line 253 of file MFrontGeneric.h.

◆ InternalVariable

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
using MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::InternalVariable
Initial value:

Definition at line 263 of file MFrontGeneric.h.

◆ KelvinMatrix

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
using MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::KelvinMatrix
Initial value:
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), kelvin_vector_dimensions(DisplacementDim), Eigen::RowMajor > KelvinMatrixType

Definition at line 261 of file MFrontGeneric.h.

◆ KelvinVector

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
using MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::KelvinVector
Initial value:
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType

Definition at line 259 of file MFrontGeneric.h.

Constructor & Destructor Documentation

◆ MFrontGeneric()

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::MFrontGeneric ( mgis::behaviour::Behaviour && behaviour,
std::vector< ParameterLib::Parameter< double > const * > && material_properties,
std::map< std::string, ParameterLib::Parameter< double > const * > && state_variables_initial_properties,
std::optional< ParameterLib::CoordinateSystem > const & local_coordinate_system )
inline

Definition at line 266 of file MFrontGeneric.h.

281 ? &local_coordinate_system.value()
282 : nullptr)
283 {
284 {
285 auto check_gradient = [&grads = _behaviour.gradients,
286 hyp = _behaviour.hypothesis,
287 i = 0]<typename Grad>(Grad) mutable
288 {
289 // TODO allow reordering of gradients and thermodynamic forces?
290 if (grads[i].name != Grad::name)
291 {
292 OGS_FATAL(
293 "OGS expects the {}th gradient to be {} but MFront "
294 "provides {}.",
295 i, Grad::name, grads[i].name);
296 }
297
298 if (grads[i].type != Grad::type)
299 {
300 OGS_FATAL(
301 "The behaviour's {}th driver ({}) must be of type {}.",
303 }
304
307 {
308 OGS_FATAL(
309 "The behaviour's {}th driver's ({}) size in OGS is {} "
310 "but {} in MFront.",
311 i, grads[i].name,
314 }
315
316 i++;
317 };
318
319 if (_behaviour.gradients.size() !=
321 OGS_FATAL(
322 "The behaviour must have exactly {} gradients as input.",
324
326 }
327
328 {
329 auto check_tdyn_force = [&tdfs = _behaviour.thermodynamic_forces,
330 hyp = _behaviour.hypothesis,
331 i = 0]<typename TDF>(TDF) mutable
332 {
333 if (tdfs[i].name != TDF::name)
334 {
335 OGS_FATAL(
336 "OGS expects the {}th thermodynamic force to be {} but "
337 "MFront provides {}.",
338 i, TDF::name, tdfs[i].name);
339 }
340
341 if (tdfs[i].type != TDF::type)
342 {
343 OGS_FATAL(
344 "The behaviour's {}th thermodynamic force ({}) must be "
345 "of type {}.",
347 }
348
351 {
352 OGS_FATAL(
353 "The behaviour's {}th thermodynamic force ({}) must "
354 "have size {} instead of {}.",
357 }
358
359 i++;
360 };
361
362 if (_behaviour.thermodynamic_forces.size() !=
364 OGS_FATAL(
365 "The behaviour must compute exactly {} thermodynamic "
366 "forces.",
368
370 }
371
372 auto const hypothesis = _behaviour.hypothesis;
373
374 static_assert(
376 "Temperature is the only allowed external state variable.");
377
378 if (!_behaviour.esvs.empty())
379 {
380 if (_behaviour.esvs[0].name != "Temperature")
381 {
382 OGS_FATAL(
383 "Only temperature is supported as external state "
384 "variable.");
385 }
386
388 hypothesis) != 1)
389 OGS_FATAL(
390 "Temperature must be a scalar instead of having {:d} "
391 "components.",
393 _behaviour.thermodynamic_forces[0], hypothesis));
394 }
395
396 if (_behaviour.mps.size() != _material_properties.size())
397 {
398 ERR("There are {:d} material properties in the loaded behaviour:",
399 _behaviour.mps.size());
400 for (auto const& mp : _behaviour.mps)
401 {
402 ERR("\t{:s}", mp.name);
403 }
404 OGS_FATAL("But the number of passed material properties is {:d}.",
405 _material_properties.size());
406 }
407 }
#define OGS_FATAL(...)
Definition Error.h:26
void ERR(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:48
std::vector< ParameterLib::Parameter< double > const * > _material_properties
std::map< std::string, ParameterLib::Parameter< double > const * > _state_variables_initial_properties
ParameterLib::CoordinateSystem const *const _local_coordinate_system
const char * varTypeToString(int v)
int getEquivalentPlasticStrainOffset(mgis::behaviour::Behaviour const &b)

Member Function Documentation

◆ computeFreeEnergyDensity()

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
double MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::computeFreeEnergyDensity ( double const ,
ParameterLib::SpatialPosition const & ,
double const ,
KelvinVector const & ,
KelvinVector const & ,
typename MechanicsBase< DisplacementDim >::MaterialStateVariables const &  ) const
inline

Definition at line 716 of file MFrontGeneric.h.

724 {
725 // TODO implement
727 }

◆ createMaterialStateVariables()

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
std::unique_ptr< typename MechanicsBase< DisplacementDim >::MaterialStateVariables > MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::createMaterialStateVariables ( ) const
inline

◆ createTangentOperatorBlocksView()

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
template<typename ForcesGradsCombinations = typename ForcesGradsCombinations<Gradients, TDynForces, ExtStateVars>::type>
OGSMFrontTangentOperatorBlocksView< DisplacementDim, ForcesGradsCombinations > MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::createTangentOperatorBlocksView ( ) const
inline

Definition at line 686 of file MFrontGeneric.h.

◆ createThermodynamicForcesView()

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
OGSMFrontThermodynamicForcesView< DisplacementDim, TDynForces > MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::createThermodynamicForcesView ( ) const
inline

Definition at line 694 of file MFrontGeneric.h.

◆ getBulkModulus()

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
double MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::getBulkModulus ( double const ,
ParameterLib::SpatialPosition const & ,
KelvinMatrix const *const C ) const
inline

Definition at line 699 of file MFrontGeneric.h.

702 {
703 if (C == nullptr)
704 {
705 OGS_FATAL(
706 "MFront::getBulkModulus() requires the tangent stiffness C "
707 "input "
708 "argument to be valid.");
709 }
713 return 1. / 9. * identity2.transpose() * *C * identity2;
714 }
constexpr int kelvin_vector_dimensions(int const displacement_dim)
Kelvin vector dimensions for given displacement dimension.

◆ getInternalVariables()

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
std::vector< InternalVariable > MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::getInternalVariables ( ) const
inline

Definition at line 625 of file MFrontGeneric.h.

626 {
628
629 for (auto const& iv : _behaviour.isvs)
630 {
631 auto const name = iv.name;
633 _behaviour.isvs, name, _behaviour.hypothesis);
634 auto const size =
636
637 // TODO (naumov): For orthotropic materials the internal variables
638 // should be rotated to the global coordinate system before output.
639 // MFront stores the variables in local coordinate system.
640 // The `size` variable could be used to find out the type of
641 // variable.
643 name, static_cast<int>(size),
644 [offset, size](
645 typename MechanicsBase<
648 {
650 DisplacementDim> const*>(&state) != nullptr);
651 auto const& internal_state_variables =
653 DisplacementDim> const&>(state)
654 ._behaviour_data.s1.internal_state_variables;
655
656 cache.resize(size);
658 size,
659 begin(cache));
660 return cache;
661 },
662 [offset, size](typename MechanicsBase<
665 {
667 DisplacementDim> const*>(&state) != nullptr);
669 static_cast<
671 state)
672 ._behaviour_data.s1.internal_state_variables;
673
674 return {internal_state_variables.data() + offset, size};
675 }};
677 }
678
679 return internal_variables;
680 }
typename MechanicsBase< DisplacementDim >::InternalVariable InternalVariable

◆ initializeInternalStateVariables()

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
void MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::initializeInternalStateVariables ( double const t,
ParameterLib::SpatialPosition const & x,
typename MechanicsBase< DisplacementDim >::MaterialStateVariables & material_state_variables ) const
inline

Definition at line 417 of file MFrontGeneric.h.

422 {
425
426 auto& state =
429
430 auto const& ivs = getInternalVariables();
431
433 {
434 // find corresponding internal variable
435 auto const& iv = BaseLib::findElementOrError(
436 ivs,
437 [name = name](InternalVariable const& iv)
438 { return iv.name == name; },
439 [name = name]()
440 { OGS_FATAL("Internal variable `{:s}' not found.", name); });
441
442 // evaluate parameter
443 std::vector<double> values = (*parameter)(t, x);
444
445 // copy parameter data into iv
446 auto const values_span = iv.reference(state);
447 assert(values.size() == values_span.size());
449 }
450
451 auto const& s1 = state._behaviour_data.s1.internal_state_variables;
452 auto& s0 = state._behaviour_data.s0.internal_state_variables;
454 }
std::vector< InternalVariable > getInternalVariables() const
ranges::range_reference_t< Range > findElementOrError(Range &range, std::predicate< ranges::range_reference_t< Range > > auto &&predicate, std::invocable auto error_callback)
Definition Algorithm.h:81

◆ integrateStress()

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
std::optional< std::tuple< OGSMFrontThermodynamicForcesData, std::unique_ptr< typename MechanicsBase< DisplacementDim >::MaterialStateVariables >, OGSMFrontTangentOperatorData > > MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::integrateStress ( MaterialPropertyLib::VariableArray const & variable_array_prev,
MaterialPropertyLib::VariableArray const & variable_array,
double const t,
ParameterLib::SpatialPosition const & x,
double const dt,
typename MechanicsBase< DisplacementDim >::MaterialStateVariables const & material_state_variables ) const
inline

Definition at line 460 of file MFrontGeneric.h.

468 {
469 using namespace MathLib::KelvinVector;
470
471 assert(
474 // New state, copy of current one, packed in unique_ptr for return.
475 auto state = std::make_unique<
479 auto& behaviour_data = state->_behaviour_data;
480
481 behaviour_data.dt = dt;
482 behaviour_data.rdt = 1.0;
483 behaviour_data.K[0] =
484 4.0; // if K[0] is greater than 3.5, the consistent
485 // tangent operator must be computed.
486 if (behaviour_data.s0.b.btype ==
488 {
489 behaviour_data.K[1] = 1.0; // The second Piola-Kirchoff stress is
490 // used for the stress measure
491 behaviour_data.K[2] =
492 1.0; // The derivative of the second
493 // Piola-Kirchoff stress with respect to
494 // the Green-Lagrange strain is returned.
495 }
496
497 // evaluate parameters at (t, x)
498 {
499 {
500 auto out = behaviour_data.s0.material_properties.begin();
501 for (auto* param : _material_properties)
502 {
503 auto const& vals = (*param)(t - dt, x);
504 out = std::copy(vals.begin(), vals.end(), out);
505 }
506 assert(out == behaviour_data.s0.material_properties.end());
507 }
508
509 {
510 auto out = behaviour_data.s1.material_properties.begin();
511 for (auto* param : _material_properties)
512 {
513 auto const& vals = (*param)(t, x);
514 out = std::copy(vals.begin(), vals.end(), out);
515 }
516 assert(out == behaviour_data.s1.material_properties.end());
517 }
518 }
519
520 // TODO unify with gradient handling? Make gradient and external state
521 // var input both optional?
522 if (!behaviour_data.s0.external_state_variables.empty())
523 {
524 // assuming that there is only temperature
525 // NOTE the temperature can be NaN, e.g., if OGS's process does not
526 // have a temperature defined
527 behaviour_data.s0.external_state_variables[0] =
528 variable_array_prev.temperature;
529 }
530
531 if (!behaviour_data.s1.external_state_variables.empty())
532 {
533 // assuming that there is only temperature
534 // NOTE the temperature can be NaN, e.g., if OGS's process does not
535 // have a temperature defined
536 behaviour_data.s1.external_state_variables[0] =
537 variable_array.temperature;
538 }
539
540 // rotation tensor
543 {
545 _local_coordinate_system->transformation<DisplacementDim>(x));
546 }
547
550 variable_array_prev, Q, behaviour_data.s0.gradients.data()});
553 variable_array, Q, behaviour_data.s1.gradients.data()});
554
555 // previous and current state of thermodynamic forces are both set from
556 // variable_array_prev. TODO optimization potential compute Q * grad_ogs
557 // only once for both cases.
561 behaviour_data.s0.thermodynamic_forces.data()});
565 behaviour_data.s1.thermodynamic_forces.data()});
566
569 if (status != 1)
570 {
572 "MFront: integration failed with status " +
573 std::to_string(status) + ".");
574 }
575
577 tdyn_forces_data.data.resize( // TODO data stored on heap but size is
578 // compile-time constant
579 behaviour_data.s1.thermodynamic_forces.size());
580
583 &in_data = behaviour_data.s1.thermodynamic_forces,
584 &Q]<typename TDF>(TDF tdf)
585 {
586 OGSMFrontThermodynamicForcesView<DisplacementDim, TDynForces>
587 view;
588
589 if constexpr (TDF::type ==
590 mgis::behaviour::Variable::Type::STENSOR)
591 {
592 if (Q)
593 {
594 view.block(tdf, out_data) =
595 *Q * eigenSwap45View(view.block(tdf, in_data));
596 }
597 else
598 {
599 view.block(tdf, out_data) =
600 eigenSwap45View(view.block(tdf, in_data));
601 }
602 }
603 else if constexpr (TDF::type ==
605 {
606 view.block(tdf, out_data) = view.block(tdf, in_data);
607 }
608 else
609 {
610 OGS_FATAL("Not yet implemented.");
611 }
612 });
613
614 return std::make_optional(
623 }

Member Data Documentation

◆ _behaviour

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
mgis::behaviour::Behaviour MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::_behaviour
private

Definition at line 730 of file MFrontGeneric.h.

◆ _local_coordinate_system

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
ParameterLib::CoordinateSystem const* const MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::_local_coordinate_system
private

Definition at line 735 of file MFrontGeneric.h.

◆ _material_properties

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
std::vector<ParameterLib::Parameter<double> const*> MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::_material_properties
private

Definition at line 732 of file MFrontGeneric.h.

◆ _state_variables_initial_properties

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
std::map<std::string, ParameterLib::Parameter<double> const*> MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::_state_variables_initial_properties
private

Definition at line 734 of file MFrontGeneric.h.

◆ equivalent_plastic_strain_offset_

template<int DisplacementDim, typename Gradients, typename TDynForces, typename ExtStateVars>
int const MaterialLib::Solids::MFront::MFrontGeneric< DisplacementDim, Gradients, TDynForces, ExtStateVars >::equivalent_plastic_strain_offset_
private

Definition at line 731 of file MFrontGeneric.h.


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