v0.14.0
Public Member Functions | Public Attributes | List of all members
GelModule::Gel::ConstitutiveEquation< TYPE > Struct Template Reference

Constitutive model functions. More...

#include <users_modules/gels/src/Gels.hpp>

Inheritance diagram for GelModule::Gel::ConstitutiveEquation< TYPE >:
[legend]
Collaboration diagram for GelModule::Gel::ConstitutiveEquation< TYPE >:
[legend]

Public Member Functions

 ConstitutiveEquation (map< int, BlockMaterialData > &data)
 
virtual ~ConstitutiveEquation ()
 
virtual PetscErrorCode calculateCauchyDefromationTensor ()
 
virtual PetscErrorCode calculateStrainTotal ()
 Calculate total strain. More...
 
virtual PetscErrorCode calculateTraceStrainTotalDot ()
 
virtual PetscErrorCode calculateStressAlpha ()
 Calculate stress in spring alpha. More...
 
virtual PetscErrorCode calculateStressBeta ()
 Calculate stress in spring beta. More...
 
virtual PetscErrorCode calculateStrainHatFlux ()
 Calculate rate of strain hat. More...
 
virtual PetscErrorCode calculateStressBetaHat ()
 Calculate stress due to concentration of solvent molecules. More...
 
virtual PetscErrorCode calculateStressTotal ()
 
virtual PetscErrorCode calculateResidualStrainHat ()
 
virtual PetscErrorCode calculateSolventFlux ()
 Calculate flux. More...
 
virtual PetscErrorCode calculateSolventConcentrationDot ()
 Calculate solvent concentration rate. More...
 

Public Attributes

int iD
 
map< int, BlockMaterialData > & dAta
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > F
 Gradient of deformation. More...
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > FDot
 Rate of gradient of deformation. More...
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > strainHat
 Internal variable, strain in dashpot beta. More...
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > strainHatDot
 Internal variable, strain in dashpot beta. More...
 
TYPE mU
 Solvent concentration. More...
 
ublas::vector< TYPE, ublas::bounded_array< TYPE, 3 > > gradientMu
 Gradient of solvent concentration. More...
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > C
 Cauchy deformation. More...
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > gradientU
 Gradient of displacements. More...
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > strainTotal
 Total strain applied at integration point. More...
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > stressAlpha
 Stress generated by spring alpha. More...
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > stressBeta
 Stress generated by spring beta. More...
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > strainHatFlux
 Rate of dashpot (beta) strain. More...
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > stressBetaHat
 Stress as result of volume change due to solvent concentration. More...
 
TYPE traceStrainTotal
 
TYPE traceStrainHat
 
TYPE traceStressBeta
 
TYPE traceStrainTotalDot
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > stressTotal
 Total stress. More...
 
ublas::vector< TYPE, ublas::bounded_array< TYPE, 9 > > solventFlux
 Solvent flux. More...
 
TYPE solventConcentrationDot
 Volume rate change. More...
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > residualStrainHat
 Residual for calculation epsilon hat. More...
 

Detailed Description

template<typename TYPE>
struct GelModule::Gel::ConstitutiveEquation< TYPE >

Constitutive model functions.

Gel model

Definition at line 114 of file Gels.hpp.

Constructor & Destructor Documentation

◆ ConstitutiveEquation()

template<typename TYPE >
GelModule::Gel::ConstitutiveEquation< TYPE >::ConstitutiveEquation ( map< int, BlockMaterialData > &  data)
inline

Definition at line 119 of file Gels.hpp.

119  :
120  dAta(data) {
121  }

◆ ~ConstitutiveEquation()

template<typename TYPE >
virtual GelModule::Gel::ConstitutiveEquation< TYPE >::~ConstitutiveEquation ( )
inlinevirtual

Definition at line 122 of file Gels.hpp.

122 {}

Member Function Documentation

◆ calculateCauchyDefromationTensor()

template<typename TYPE >
virtual PetscErrorCode GelModule::Gel::ConstitutiveEquation< TYPE >::calculateCauchyDefromationTensor ( )
inlinevirtual

Definition at line 156 of file Gels.hpp.

156  {
157  PetscFunctionBegin;
158  C.resize(3,3);
159  noalias(C) = prod(trans(F),F);
160  PetscFunctionReturn(0);
161  }

◆ calculateResidualStrainHat()

template<typename TYPE >
virtual PetscErrorCode GelModule::Gel::ConstitutiveEquation< TYPE >::calculateResidualStrainHat ( )
inlinevirtual

Definition at line 301 of file Gels.hpp.

301  {
302  PetscFunctionBegin;
303  residualStrainHat.resize(3,3,false);
305  PetscFunctionReturn(0);
306  }

◆ calculateSolventConcentrationDot()

template<typename TYPE >
virtual PetscErrorCode GelModule::Gel::ConstitutiveEquation< TYPE >::calculateSolventConcentrationDot ( )
inlinevirtual

Calculate solvent concentration rate.

FIXME: For simplicity as first approximation set volume change as trace of gradient total strain

Definition at line 332 of file Gels.hpp.

332  {
333  PetscFunctionBegin;
334  solventConcentrationDot = traceStrainTotalDot/*/dAta[iD].oMega*/;
335  PetscFunctionReturn(0);
336  }

◆ calculateSolventFlux()

template<typename TYPE >
virtual PetscErrorCode GelModule::Gel::ConstitutiveEquation< TYPE >::calculateSolventFlux ( )
inlinevirtual

Calculate flux.

\[ J_k = -\left( \frac{\kappa}{\eta\Omega^2} \right) \mu_{,i} \]

dAta[iD].oMega

Definition at line 319 of file Gels.hpp.

319  {
320  PetscFunctionBegin;
321  double a = dAta[iD].pErmeability/(dAta[iD].vIscosity*dAta[iD].oMega/**dAta[iD].oMega*/);
323  PetscFunctionReturn(0);
324  }

◆ calculateStrainHatFlux()

template<typename TYPE >
virtual PetscErrorCode GelModule::Gel::ConstitutiveEquation< TYPE >::calculateStrainHatFlux ( )
inlinevirtual

Calculate rate of strain hat.

\[ \frac{\partial \hat{\varepsilon}_{ij}}{\partial t} = \frac{1}{2\hat{G}^\beta} \left( \sigma_{ij}^\beta - \frac{\hat{v}^\beta}{1+\hat{v}^\beta}\sigma_{kk}^\beta\delta_{ij} \right) \]

Definition at line 255 of file Gels.hpp.

255  {
256  PetscFunctionBegin;
257  traceStressBeta = 0;
258  for(int ii = 0;ii<3;ii++) {
259  traceStressBeta = stressBeta(ii,ii);
260  }
261  strainHatFlux.resize(3,3,false);
262  double a = (1.0/(2.0*dAta[iD].gBetaHat));
263  noalias(strainHatFlux) = a*stressBeta;
264  double b = a*(dAta[iD].vBetaHat/(1.0+dAta[iD].vBetaHat));
265  for(int ii = 0;ii<3;ii++) {
266  strainHatFlux(ii,ii) -= b*traceStressBeta;
267  }
268  PetscFunctionReturn(0);
269  }

◆ calculateStrainTotal()

template<typename TYPE >
virtual PetscErrorCode GelModule::Gel::ConstitutiveEquation< TYPE >::calculateStrainTotal ( )
inlinevirtual

Calculate total strain.

Definition at line 166 of file Gels.hpp.

166  {
167  PetscFunctionBegin;
168  gradientU.resize(3,3,false);
169  noalias(gradientU) = F;
170  for(int ii = 0;ii<3;ii++) {
171  gradientU(ii,ii) -= 1;
172  }
173  strainTotal.resize(3,3,false);
174  noalias(strainTotal) = gradientU + trans(gradientU);
175  strainTotal *= 0.5;
176  PetscFunctionReturn(0);
177  }

◆ calculateStressAlpha()

template<typename TYPE >
virtual PetscErrorCode GelModule::Gel::ConstitutiveEquation< TYPE >::calculateStressAlpha ( )
inlinevirtual

Calculate stress in spring alpha.

\[ \sigma^\alpha_{ij} = 2G^\alpha\left(\varepsilon_{ij} + \frac{v^\alpha}{1-2v^\alpha}\varepsilon_{kk}\delta_{ij}\right) \]

Note: In general implementation for large strain this function should calculate Piola-Kirchhoff Stress I.

Definition at line 200 of file Gels.hpp.

200  {
201  PetscFunctionBegin;
202  traceStrainTotal = 0;
203  for(int ii = 0;ii<3;ii++) {
204  traceStrainTotal += strainTotal(ii,ii);
205  }
206  stressAlpha.resize(3,3,false);
207  double a = 2.0*dAta[iD].gAlpha;
208  double b = a*(dAta[iD].vAlpha/(1.0-2.0*dAta[iD].vAlpha));
209  noalias(stressAlpha) = a*strainTotal;
210  for(int ii = 0; ii<3; ii++){
211  stressAlpha(ii,ii) += b*traceStrainTotal;
212  }
213  PetscFunctionReturn(0);
214  }

◆ calculateStressBeta()

template<typename TYPE >
virtual PetscErrorCode GelModule::Gel::ConstitutiveEquation< TYPE >::calculateStressBeta ( )
inlinevirtual

Calculate stress in spring beta.

\[ \sigma^\beta_{ij} = 2G^\beta\left[ (\varepsilon_{ij}-\hat{\varepsilon}_{ij}) + \frac{v^\beta}{1-2v^\beta}(\varepsilon_{kk}-\hat{\varepsilon}_{kk})\delta_{ij} \right] \]

Definition at line 226 of file Gels.hpp.

226  {
227  PetscFunctionBegin;
228  traceStrainTotal = 0;
229  traceStrainHat = 0;
230  for(int ii = 0;ii<3;ii++) {
231  traceStrainHat += strainHat(ii,ii);
232  traceStrainTotal += strainTotal(ii,ii);
233  }
234  stressBeta.resize(3,3,false);
235  double a = 2.0*dAta[iD].gBeta;
236  double b = a*(dAta[iD].vBeta/(1.0-2.0*dAta[iD].vBeta));
237  noalias(stressBeta) = a*(strainTotal-strainHat);
238  for(int ii = 0;ii<3;ii++) {
240  }
241  PetscFunctionReturn(0);
242  }

◆ calculateStressBetaHat()

template<typename TYPE >
virtual PetscErrorCode GelModule::Gel::ConstitutiveEquation< TYPE >::calculateStressBetaHat ( )
inlinevirtual

Calculate stress due to concentration of solvent molecules.

\[ \hat{\sigma}_{ij}^\beta = \frac{\Delta\mu}{\Omega}\delta_{ij} \]

note that \(\Delta\mu = \mu-\mu_0\). In current version of code value of \(\Delta\mu\) is approximated.

Definition at line 280 of file Gels.hpp.

280  {
281  PetscFunctionBegin;
282  stressBetaHat.resize(3,3,false);
283  stressBetaHat.clear();
284  for(int ii=0; ii<3; ii++){
285  stressBetaHat(ii,ii) = -(mU-dAta[iD].mU0)/dAta[iD].oMega;
286  }
287  PetscFunctionReturn(0);
288  }

◆ calculateStressTotal()

template<typename TYPE >
virtual PetscErrorCode GelModule::Gel::ConstitutiveEquation< TYPE >::calculateStressTotal ( )
inlinevirtual

Definition at line 292 of file Gels.hpp.

292  {
293  PetscFunctionBegin;
294  stressTotal.resize(3,3,false);
295  noalias(stressTotal) = stressAlpha;
296  noalias(stressTotal) += stressBeta;
297  noalias(stressTotal) += stressBetaHat;
298  PetscFunctionReturn(0);
299  }

◆ calculateTraceStrainTotalDot()

template<typename TYPE >
virtual PetscErrorCode GelModule::Gel::ConstitutiveEquation< TYPE >::calculateTraceStrainTotalDot ( )
inlinevirtual

\berief Calculate trace of rate of gradient of deformation

Definition at line 181 of file Gels.hpp.

181  {
182  PetscFunctionBegin;
184  for(int ii = 0;ii<3;ii++) {
185  traceStrainTotalDot += FDot(ii,ii);
186  }
187  PetscFunctionReturn(0);
188  }

Member Data Documentation

◆ C

template<typename TYPE >
ublas::matrix<TYPE,ublas::row_major,ublas::bounded_array<TYPE,9> > GelModule::Gel::ConstitutiveEquation< TYPE >::C

Cauchy deformation.

Definition at line 135 of file Gels.hpp.

◆ dAta

template<typename TYPE >
map<int,BlockMaterialData>& GelModule::Gel::ConstitutiveEquation< TYPE >::dAta

Definition at line 117 of file Gels.hpp.

◆ F

template<typename TYPE >
ublas::matrix<TYPE,ublas::row_major,ublas::bounded_array<TYPE,9> > GelModule::Gel::ConstitutiveEquation< TYPE >::F

Gradient of deformation.

Definition at line 126 of file Gels.hpp.

◆ FDot

template<typename TYPE >
ublas::matrix<TYPE,ublas::row_major,ublas::bounded_array<TYPE,9> > GelModule::Gel::ConstitutiveEquation< TYPE >::FDot

Rate of gradient of deformation.

Definition at line 127 of file Gels.hpp.

◆ gradientMu

template<typename TYPE >
ublas::vector<TYPE,ublas::bounded_array<TYPE,3> > GelModule::Gel::ConstitutiveEquation< TYPE >::gradientMu

Gradient of solvent concentration.

Definition at line 131 of file Gels.hpp.

◆ gradientU

template<typename TYPE >
ublas::matrix<TYPE,ublas::row_major,ublas::bounded_array<TYPE,9> > GelModule::Gel::ConstitutiveEquation< TYPE >::gradientU

Gradient of displacements.

Definition at line 136 of file Gels.hpp.

◆ iD

template<typename TYPE >
int GelModule::Gel::ConstitutiveEquation< TYPE >::iD

Definition at line 116 of file Gels.hpp.

◆ mU

template<typename TYPE >
TYPE GelModule::Gel::ConstitutiveEquation< TYPE >::mU

Solvent concentration.

Definition at line 130 of file Gels.hpp.

◆ residualStrainHat

template<typename TYPE >
ublas::matrix<TYPE,ublas::row_major,ublas::bounded_array<TYPE,9> > GelModule::Gel::ConstitutiveEquation< TYPE >::residualStrainHat

Residual for calculation epsilon hat.

Definition at line 154 of file Gels.hpp.

◆ solventConcentrationDot

template<typename TYPE >
TYPE GelModule::Gel::ConstitutiveEquation< TYPE >::solventConcentrationDot

Volume rate change.

Definition at line 153 of file Gels.hpp.

◆ solventFlux

template<typename TYPE >
ublas::vector<TYPE,ublas::bounded_array<TYPE,9> > GelModule::Gel::ConstitutiveEquation< TYPE >::solventFlux

Solvent flux.

Definition at line 152 of file Gels.hpp.

◆ strainHat

template<typename TYPE >
ublas::matrix<TYPE,ublas::row_major,ublas::bounded_array<TYPE,9> > GelModule::Gel::ConstitutiveEquation< TYPE >::strainHat

Internal variable, strain in dashpot beta.

Definition at line 128 of file Gels.hpp.

◆ strainHatDot

template<typename TYPE >
ublas::matrix<TYPE,ublas::row_major,ublas::bounded_array<TYPE,9> > GelModule::Gel::ConstitutiveEquation< TYPE >::strainHatDot

Internal variable, strain in dashpot beta.

Definition at line 129 of file Gels.hpp.

◆ strainHatFlux

template<typename TYPE >
ublas::matrix<TYPE,ublas::row_major,ublas::bounded_array<TYPE,9> > GelModule::Gel::ConstitutiveEquation< TYPE >::strainHatFlux

Rate of dashpot (beta) strain.

Definition at line 141 of file Gels.hpp.

◆ strainTotal

template<typename TYPE >
ublas::matrix<TYPE,ublas::row_major,ublas::bounded_array<TYPE,9> > GelModule::Gel::ConstitutiveEquation< TYPE >::strainTotal

Total strain applied at integration point.

Definition at line 137 of file Gels.hpp.

◆ stressAlpha

template<typename TYPE >
ublas::matrix<TYPE,ublas::row_major,ublas::bounded_array<TYPE,9> > GelModule::Gel::ConstitutiveEquation< TYPE >::stressAlpha

Stress generated by spring alpha.

Definition at line 139 of file Gels.hpp.

◆ stressBeta

template<typename TYPE >
ublas::matrix<TYPE,ublas::row_major,ublas::bounded_array<TYPE,9> > GelModule::Gel::ConstitutiveEquation< TYPE >::stressBeta

Stress generated by spring beta.

Definition at line 140 of file Gels.hpp.

◆ stressBetaHat

template<typename TYPE >
ublas::matrix<TYPE,ublas::row_major,ublas::bounded_array<TYPE,9> > GelModule::Gel::ConstitutiveEquation< TYPE >::stressBetaHat

Stress as result of volume change due to solvent concentration.

Definition at line 142 of file Gels.hpp.

◆ stressTotal

template<typename TYPE >
ublas::matrix<TYPE,ublas::row_major,ublas::bounded_array<TYPE,9> > GelModule::Gel::ConstitutiveEquation< TYPE >::stressTotal

Total stress.

Definition at line 151 of file Gels.hpp.

◆ traceStrainHat

template<typename TYPE >
TYPE GelModule::Gel::ConstitutiveEquation< TYPE >::traceStrainHat

Definition at line 145 of file Gels.hpp.

◆ traceStrainTotal

template<typename TYPE >
TYPE GelModule::Gel::ConstitutiveEquation< TYPE >::traceStrainTotal

Definition at line 144 of file Gels.hpp.

◆ traceStrainTotalDot

template<typename TYPE >
TYPE GelModule::Gel::ConstitutiveEquation< TYPE >::traceStrainTotalDot

Definition at line 147 of file Gels.hpp.

◆ traceStressBeta

template<typename TYPE >
TYPE GelModule::Gel::ConstitutiveEquation< TYPE >::traceStressBeta

Definition at line 146 of file Gels.hpp.


The documentation for this struct was generated from the following file:
GelModule::Gel::ConstitutiveEquation::dAta
map< int, BlockMaterialData > & dAta
Definition: Gels.hpp:117
GelModule::Gel::ConstitutiveEquation::solventConcentrationDot
TYPE solventConcentrationDot
Volume rate change.
Definition: Gels.hpp:153
GelModule::Gel::ConstitutiveEquation::solventFlux
ublas::vector< TYPE, ublas::bounded_array< TYPE, 9 > > solventFlux
Solvent flux.
Definition: Gels.hpp:152
GelModule::Gel::ConstitutiveEquation::iD
int iD
Definition: Gels.hpp:116
GelModule::Gel::ConstitutiveEquation::gradientU
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > gradientU
Gradient of displacements.
Definition: Gels.hpp:136
GelModule::Gel::ConstitutiveEquation::traceStrainTotal
TYPE traceStrainTotal
Definition: Gels.hpp:144
GelModule::Gel::ConstitutiveEquation::gradientMu
ublas::vector< TYPE, ublas::bounded_array< TYPE, 3 > > gradientMu
Gradient of solvent concentration.
Definition: Gels.hpp:131
GelModule::Gel::ConstitutiveEquation::stressTotal
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > stressTotal
Total stress.
Definition: Gels.hpp:151
GelModule::Gel::ConstitutiveEquation::traceStressBeta
TYPE traceStressBeta
Definition: Gels.hpp:146
GelModule::Gel::ConstitutiveEquation::C
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > C
Cauchy deformation.
Definition: Gels.hpp:135
GelModule::Gel::ConstitutiveEquation::stressAlpha
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > stressAlpha
Stress generated by spring alpha.
Definition: Gels.hpp:139
a
constexpr double a
Definition: approx_sphere.cpp:30
GelModule::Gel::ConstitutiveEquation::residualStrainHat
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > residualStrainHat
Residual for calculation epsilon hat.
Definition: Gels.hpp:154
GelModule::Gel::ConstitutiveEquation::traceStrainHat
TYPE traceStrainHat
Definition: Gels.hpp:145
GelModule::Gel::ConstitutiveEquation::strainTotal
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > strainTotal
Total strain applied at integration point.
Definition: Gels.hpp:137
GelModule::Gel::ConstitutiveEquation::strainHatDot
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > strainHatDot
Internal variable, strain in dashpot beta.
Definition: Gels.hpp:129
GelModule::Gel::ConstitutiveEquation::mU
TYPE mU
Solvent concentration.
Definition: Gels.hpp:130
GelModule::Gel::ConstitutiveEquation::strainHatFlux
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > strainHatFlux
Rate of dashpot (beta) strain.
Definition: Gels.hpp:141
GelModule::Gel::ConstitutiveEquation::stressBetaHat
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > stressBetaHat
Stress as result of volume change due to solvent concentration.
Definition: Gels.hpp:142
GelModule::Gel::ConstitutiveEquation::FDot
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > FDot
Rate of gradient of deformation.
Definition: Gels.hpp:127
GelModule::Gel::ConstitutiveEquation::strainHat
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > strainHat
Internal variable, strain in dashpot beta.
Definition: Gels.hpp:128
GelModule::Gel::ConstitutiveEquation::stressBeta
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > stressBeta
Stress generated by spring beta.
Definition: Gels.hpp:140
GelModule::Gel::ConstitutiveEquation::F
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > F
Gradient of deformation.
Definition: Gels.hpp:126
GelModule::Gel::ConstitutiveEquation::traceStrainTotalDot
TYPE traceStrainTotalDot
Definition: Gels.hpp:147