MarmotKelvinChain.h
Go to the documentation of this file.
1 /* ---------------------------------------------------------------------
2  * _
3  * _ __ ___ __ _ _ __ _ __ ___ ___ | |_
4  * | '_ ` _ \ / _` | '__| '_ ` _ \ / _ \| __|
5  * | | | | | | (_| | | | | | | | | (_) | |_
6  * |_| |_| |_|\__,_|_| |_| |_| |_|\___/ \__|
7  *
8  * Unit of Strength of Materials and Structural Analysis
9  * University of Innsbruck,
10  * 2020 - today
11  *
12  * festigkeitslehre@uibk.ac.at
13  *
14  * Alexander Dummer alexander.dummer@uibk.ac.at
15  *
16  * This file is part of the MAteRialMOdellingToolbox (marmot).
17  *
18  * This library is free software; you can redistribute it and/or
19  * modify it under the terms of the GNU Lesser General Public
20  * License as published by the Free Software Foundation; either
21  * version 2.1 of the License, or (at your option) any later version.
22  *
23  * The full text of the license can be found in the file LICENSE.md at
24  * the top level directory of marmot.
25  * ---------------------------------------------------------------------
26  */
27 
28 #pragma once
30 #include "Marmot/MarmotTypedefs.h"
31 #include "autodiff/forward/real.hpp"
32 #include <functional>
33 
34 namespace Marmot::Materials {
35 
36  namespace KelvinChain {
37 
38  typedef Eigen::VectorXd Properties;
39  typedef Eigen::Map< Properties > mapProperties;
40 
41  typedef Eigen::Matrix< double, 6, Eigen::Dynamic > StateVarMatrix;
42  typedef Eigen::Map< StateVarMatrix > mapStateVarMatrix;
43 
44  template < int N >
45  struct Factorial {
46  enum { value = N * Factorial< N - 1 >::value };
47  };
48 
49  template <>
50  struct Factorial< 0 > {
51  enum { value = 1 };
52  };
53 
54  template < int k >
55  double evaluatePostWidderFormula( std::function< autodiff::Real< k, double >( autodiff::Real< k, double > ) > phi,
56  double tau )
57  {
58 
59  autodiff::Real< k, double > tau_( tau * k );
60 
61  double val = -pow( -tau * k, k ) / double( Factorial< k - 1 >::value );
62  val *= autodiff::derivatives( phi, autodiff::along( 1. ), autodiff::at( tau_ ) )[k];
63  return val;
64  }
65 
66  template < int k >
67  double approximateZerothCompliance( std::function< autodiff::Real< k, double >( autodiff::Real< k, double > ) > phi,
68  double tauMin,
69  double spacing = 10. )
70  {
72  double val_ = -pow( -k, k ) * pow( tau, k - 1 ) / double( Factorial< k - 1 >::value );
73  autodiff::Real< k, double > tau_( tau * k );
74  val_ *= autodiff::derivatives( phi, autodiff::along( 1. ), autodiff::at( tau_ ) )[k];
75  return val_;
76  };
77 
78  double
80  { 1e-14, tauMin / sqrt( spacing ) },
81  100,
83  return val;
84  }
85 
86  template < int k >
87  Properties computeElasticModuli( std::function< autodiff::Real< k, double >( autodiff::Real< k, double > ) > phi,
88  Properties retardationTimes,
89  bool gaussQuadrature = false )
90  {
91  Properties elasticModuli( retardationTimes.size() );
92  double spacing = retardationTimes( 1 ) / retardationTimes( 0 );
93 
94  for ( int i = 0; i < retardationTimes.size(); i++ ) {
95  double tau = retardationTimes( i );
96  if ( !gaussQuadrature ) {
97  elasticModuli( i ) = 1. / ( log( spacing ) * evaluatePostWidderFormula< k >( phi, tau ) );
98  }
99  else {
100  elasticModuli( i ) = 1. /
101  ( log( spacing ) / 2. *
102  ( evaluatePostWidderFormula< k >( phi, tau * pow( spacing, -sqrt( 3. ) / 6. ) ) +
103  evaluatePostWidderFormula< k >( phi, tau * pow( spacing, sqrt( 3. ) / 6. ) ) ) );
104  }
105  }
106 
107  return elasticModuli;
108  }
109 
110  Properties generateRetardationTimes( int n, double min, double spacing );
111 
112  void updateStateVarMatrix( const double dT,
113  Properties elasticModuli,
114  Properties retardationTimes,
115  Eigen::Ref< StateVarMatrix > stateVars,
116  const Marmot::Vector6d& dStress,
117  const Marmot::Matrix6d& unitComplianceMatrix );
118 
119  void evaluateKelvinChain( const double dT,
120  Properties elasticModuli,
121  Properties retardationTimes,
122  StateVarMatrix stateVars,
123  double& uniaxialCompliance,
124  Marmot::Vector6d& dStrain,
125  const double factor );
126 
127  void computeLambdaAndBeta( double dT, double tau, double& lambda, double& beta );
128 
129  } // namespace KelvinChain
130 } // namespace Marmot::Materials
Marmot::Materials::KelvinChain::evaluatePostWidderFormula
double evaluatePostWidderFormula(std::function< autodiff::Real< k, double >(autodiff::Real< k, double >) > phi, double tau)
Definition: MarmotKelvinChain.h:55
Marmot::Materials::KelvinChain::evaluateKelvinChain
void evaluateKelvinChain(const double dT, Properties elasticModuli, Properties retardationTimes, StateVarMatrix stateVars, double &uniaxialCompliance, Marmot::Vector6d &dStrain, const double factor)
Definition: MarmotKelvinChain.cpp:18
Marmot::Materials::KelvinChain::computeLambdaAndBeta
void computeLambdaAndBeta(double dT, double tau, double &lambda, double &beta)
Definition: MarmotKelvinChain.cpp:57
Marmot::NumericalAlgorithms::Integration::simpson
@ simpson
Definition: MarmotNumericalIntegration.h:37
MarmotTypedefs.h
Marmot::Materials::KelvinChain::mapStateVarMatrix
Eigen::Map< StateVarMatrix > mapStateVarMatrix
Definition: MarmotKelvinChain.h:42
Marmot::Matrix6d
Eigen::Matrix< double, 6, 6 > Matrix6d
Definition: MarmotTypedefs.h:35
Marmot::Materials
Definition: MarmotKelvinChain.h:34
MarmotNumericalIntegration.h
Marmot::Materials::KelvinChain::mapProperties
Eigen::Map< Properties > mapProperties
Definition: MarmotKelvinChain.h:39
Marmot::Materials::KelvinChain::approximateZerothCompliance
double approximateZerothCompliance(std::function< autodiff::Real< k, double >(autodiff::Real< k, double >) > phi, double tauMin, double spacing=10.)
Definition: MarmotKelvinChain.h:67
Marmot::Materials::KelvinChain::updateStateVarMatrix
void updateStateVarMatrix(const double dT, Properties elasticModuli, Properties retardationTimes, Eigen::Ref< StateVarMatrix > stateVars, const Marmot::Vector6d &dStress, const Marmot::Matrix6d &unitComplianceMatrix)
Marmot::NumericalAlgorithms::Integration::integrateScalarFunction
double integrateScalarFunction(scalar_to_scalar_function_type f, const std::tuple< double, double > integrationLimits, const int n, const integrationRule intRule)
Definition: MarmotNumericalIntegration.cpp:8
Marmot::Materials::KelvinChain::computeElasticModuli
Properties computeElasticModuli(std::function< autodiff::Real< k, double >(autodiff::Real< k, double >) > phi, Properties retardationTimes, bool gaussQuadrature=false)
Definition: MarmotKelvinChain.h:87
Marmot::Materials::KelvinChain::generateRetardationTimes
Properties generateRetardationTimes(int n, double min, double spacing)
Definition: MarmotKelvinChain.cpp:10
Marmot::Materials::KelvinChain::StateVarMatrix
Eigen::Matrix< double, 6, Eigen::Dynamic > StateVarMatrix
Definition: MarmotKelvinChain.h:41
Marmot::Materials::KelvinChain::Properties
Eigen::VectorXd Properties
Definition: MarmotKelvinChain.h:38
Marmot::Vector6d
Eigen::Matrix< double, 6, 1 > Vector6d
Definition: MarmotTypedefs.h:43
Marmot::NumericalAlgorithms::Integration::scalar_to_scalar_function_type
std::function< double(const double x) > scalar_to_scalar_function_type
Definition: MarmotNumericalIntegration.h:35
Marmot::Materials::KelvinChain::Factorial::value
@ value
Definition: MarmotKelvinChain.h:46
Marmot::Materials::KelvinChain::Factorial
Definition: MarmotKelvinChain.h:45
Marmot::FiniteElement::Spatial1D::Bar2::N
NSized N(double xi)
Definition: MarmotFiniteElement1D.cpp:15