v0.14.0
gel_constitutive_equation_test.cpp
Go to the documentation of this file.
1 /** \file gel_constitutive_equation_test.cpp
2  \brief Atom test verifying implementation of gel constitutive equations
3  \ingroup gel
4 */
5 
6 /* This file is part of MoFEM.
7  * MoFEM is free software: you can redistribute it and/or modify it under
8  * the terms of the GNU Lesser General Public License as published by the
9  * Free Software Foundation, either version 3 of the License, or (at your
10  * option) any later version.
11  *
12  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
13  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
14  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
15  * License for more details.
16  *
17  * You should have received a copy of the GNU Lesser General Public
18  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
19 
20 #include <MoFEM.hpp>
21 using namespace MoFEM;
22 #include <PostProcOnRefMesh.hpp>
23 
24 #include <boost/numeric/ublas/vector_proxy.hpp>
25 #include <boost/numeric/ublas/matrix.hpp>
26 #include <boost/numeric/ublas/matrix_proxy.hpp>
27 #include <boost/numeric/ublas/vector.hpp>
28 #include <adolc/adolc.h>
29 #include <Gels.hpp>
30 using namespace GelModule;
31 
32 #include <boost/iostreams/tee.hpp>
33 #include <boost/iostreams/stream.hpp>
34 #include <fstream>
35 #include <iostream>
36 namespace bio = boost::iostreams;
37 using bio::tee_device;
38 using bio::stream;
39 
40 
41 
42 static char help[] = "...\n\n";
43 
44 int main(int argc, char *argv[]) {
45 
46  MoFEM::Core::Initialize(&argc,&argv,(char *)0,help);
47 
48  try {
49 
50  moab::Core mb_instance;
51  moab::Interface& moab = mb_instance;
52  int rank;
53  MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
54  MoFEM::Core core(moab);
55  MoFEM::Interface& m_field = core;
56 
57  // Create gel instance
58  Gel gel(m_field);
59 
60  map<int,Gel::BlockMaterialData> material_data;
61 
62  // Set material parameters
63  material_data[0].gAlpha = 1;
64  material_data[0].vAlpha = 0.3;
65  material_data[0].gBeta = 1;
66  material_data[0].vBeta = 0.3;
67  material_data[0].gBetaHat = 1;
68  material_data[0].vBetaHat = 0.3;
69  material_data[0].oMega = 1;
70  material_data[0].vIscosity = 1;
71  material_data[0].pErmeability = 2.;
72 
73  Gel::ConstitutiveEquation<double> ce(material_data);
74 
75  // Set state variables
76  MatrixDouble3by3 &F = ce.F;
77  F.resize(3,3);
78  F.clear();
79  F(0,0) = 0.01; F(0,1) = 0.02; F(0,2) = 0.03;
80  F(1,0) = 0.04; F(1,1) = 0.05; F(1,2) = 0.06;
81  F(2,0) = 0.07; F(2,1) = 0.08; F(2,2) = 0.09;
82 
83  MatrixDouble3by3 &F_dot = ce.FDot;
84  F_dot.resize(3,3);
85  F_dot.clear();
86  F_dot(0,0) = 0.01; F_dot(0,1) = 0.02; F_dot(0,2) = 0.03;
87  F_dot(1,0) = 0.04; F_dot(1,1) = 0.05; F_dot(1,2) = 0.06;
88  F_dot(2,0) = 0.07; F_dot(2,1) = 0.08; F_dot(2,2) = 0.09;
89 
90  MatrixDouble3by3 &strainHat = ce.strainHat;
91  strainHat.resize(3,3);
92  strainHat.clear();
93  strainHat(0,0) = 0.01; strainHat(0,1) = 0.02; strainHat(0,2) = 0.03;
94  strainHat(1,0) = 0.04; strainHat(1,1) = 0.05; strainHat(1,2) = 0.06;
95  strainHat(2,0) = 0.07; strainHat(2,1) = 0.08; strainHat(2,2) = 0.09;
96 
97  MatrixDouble3by3 &strainHatDot = ce.strainHatDot;
98  strainHatDot.resize(3,3);
99  strainHatDot.clear();
100  strainHatDot(0,0) = 0.01; strainHatDot(0,1) = 0.02; strainHatDot(0,2) = 0.03;
101  strainHatDot(1,0) = 0.04; strainHatDot(1,1) = 0.05; strainHatDot(1,2) = 0.06;
102  strainHatDot(2,0) = 0.07; strainHatDot(2,1) = 0.08; strainHatDot(2,2) = 0.09;
103 
104  ce.mU = 1;
105  VectorDouble3 &gradient_mu = ce.gradientMu;
106  gradient_mu.resize(3);
107  gradient_mu.clear();
108  gradient_mu[0] = 1;
109  gradient_mu[1] = 0;
110  gradient_mu[2] = 0;
111 
112  // Creeate tee like output, printing results on screen and simultaneously
113  // writing results to file
114  typedef tee_device<ostream, ofstream> TeeDevice;
115  typedef stream<TeeDevice> TeeStream;
116  ofstream ofs("gel_constitutive_model_test.txt");
117  TeeDevice my_tee(cout,ofs);
118  TeeStream my_split(my_tee);
119 
120  // Do calculations
121 
122  ierr = ce.calculateCauchyDefromationTensor(); CHKERRQ(ierr);
123  my_split << "C\n" << ce.C << endl << endl;
124 
125  ierr = ce.calculateStrainTotal(); CHKERRQ(ierr);
126  my_split << "strainTotal\n" << ce.strainTotal << endl << endl;
127 
128  ierr = ce.calculateTraceStrainTotalDot(); CHKERRQ(ierr);
129  my_split << "traceStrianTotalDot\n" << ce.traceStrainTotalDot << endl << endl;
130 
131  ierr = ce.calculateSolventConcentrationDot(); CHKERRQ(ierr);
132  my_split << "solventConcentrationDot\n" << ce.solventConcentrationDot << endl << endl;
133 
134  ierr = ce.calculateStressAlpha(); CHKERRQ(ierr);
135  my_split << "stressAlpha\n" << ce.stressAlpha << endl << endl;
136 
137  ierr = ce.calculateStressBeta(); CHKERRQ(ierr);
138  my_split << "stressBeta\n" << ce.stressBeta << endl << endl;
139 
140  ierr = ce.calculateStressBetaHat(); CHKERRQ(ierr);
141  my_split << "stressBeta\n" << ce.stressBetaHat << endl << endl;
142 
143  ierr = ce.calculateStressTotal(); CHKERRQ(ierr);
144  my_split << "stressTotal\n" << ce.stressTotal << endl << endl;
145 
146  ierr = ce.calculateSolventFlux(); CHKERRQ(ierr);
147  my_split << "solventFlux\n" << ce.solventFlux << endl << endl;
148 
149  ierr = ce.calculateStrainHatFlux(); CHKERRQ(ierr);
150  my_split << "strainHatDot\n" << ce.strainHatDot << endl << endl;
151 
152  ierr = ce.calculateResidualStrainHat(); CHKERRQ(ierr);
153  my_split << "residualStrainHat\n" << ce.residualStrainHat << endl << endl;
154 
155  }
156  CATCH_ERRORS;
157 
159 
160  return 0;
161 }
TeeStream
stream< TeeDevice > TeeStream
Definition: forces_and_sources_testing_contact_prism_element.cpp:10
GelModule::Gel::ConstitutiveEquation::calculateCauchyDefromationTensor
virtual PetscErrorCode calculateCauchyDefromationTensor()
Definition: Gels.hpp:156
GelModule::Gel::ConstitutiveEquation::solventConcentrationDot
TYPE solventConcentrationDot
Volume rate change.
Definition: Gels.hpp:153
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition: Core.hpp:82
GelModule::Gel::ConstitutiveEquation::solventFlux
ublas::vector< TYPE, ublas::bounded_array< TYPE, 9 > > solventFlux
Solvent flux.
Definition: Gels.hpp:152
MoFEM::Types::VectorDouble3
VectorBoundedArray< double, 3 > VectorDouble3
Definition: Types.hpp:92
help
static char help[]
Definition: gel_constitutive_equation_test.cpp:42
GelModule::Gel::ConstitutiveEquation::calculateStrainHatFlux
virtual PetscErrorCode calculateStrainHatFlux()
Calculate rate of strain hat.
Definition: Gels.hpp:255
GelModule
Definition: Gels.hpp:27
MoFEM.hpp
GelModule::Gel::ConstitutiveEquation::gradientMu
ublas::vector< TYPE, ublas::bounded_array< TYPE, 3 > > gradientMu
Gradient of solvent concentration.
Definition: Gels.hpp:131
MoFEM::CoreTmp< 0 >::Finalize
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:112
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::calculateStressTotal
virtual PetscErrorCode calculateStressTotal()
Definition: Gels.hpp:292
GelModule::Gel
Implementation of Gel constitutive model.
Definition: Gels.hpp:74
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
TeeDevice
tee_device< std::ostream, std::ofstream > TeeDevice
Definition: forces_and_sources_testing_contact_prism_element.cpp:9
MoFEM::Interface
DeprecatedCoreInterface Interface
Definition: Interface.hpp:2010
GelModule::Gel::ConstitutiveEquation::calculateStressAlpha
virtual PetscErrorCode calculateStressAlpha()
Calculate stress in spring alpha.
Definition: Gels.hpp:200
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::calculateStrainTotal
virtual PetscErrorCode calculateStrainTotal()
Calculate total strain.
Definition: Gels.hpp:166
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
main
int main(int argc, char *argv[])
Definition: gel_constitutive_equation_test.cpp:44
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
GelModule::Gel::ConstitutiveEquation::calculateStressBetaHat
virtual PetscErrorCode calculateStressBetaHat()
Calculate stress due to concentration of solvent molecules.
Definition: Gels.hpp:280
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::calculateSolventFlux
virtual PetscErrorCode calculateSolventFlux()
Calculate flux.
Definition: Gels.hpp:319
PostProcOnRefMesh.hpp
Post-process fields on refined mesh.
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::calculateStressBeta
virtual PetscErrorCode calculateStressBeta()
Calculate stress in spring beta.
Definition: Gels.hpp:226
Gels.hpp
Implementation of Gel finite element.
MoFEM::CoreTmp< 0 >::Initialize
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:72
CATCH_ERRORS
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:385
MoFEM::Core
CoreTmp< 0 > Core
Definition: Core.hpp:1148
GelModule::Gel::ConstitutiveEquation::calculateTraceStrainTotalDot
virtual PetscErrorCode calculateTraceStrainTotalDot()
Definition: Gels.hpp:181
MoFEM::Exceptions::ierr
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76
GelModule::Gel::ConstitutiveEquation::mU
TYPE mU
Solvent concentration.
Definition: Gels.hpp:130
GelModule::Gel::ConstitutiveEquation::calculateSolventConcentrationDot
virtual PetscErrorCode calculateSolventConcentrationDot()
Calculate solvent concentration rate.
Definition: Gels.hpp:332
MoFEM::Types::MatrixDouble3by3
MatrixBoundedArray< double, 9 > MatrixDouble3by3
Definition: Types.hpp:105
GelModule::Gel::ConstitutiveEquation
Constitutive model functions.
Definition: Gels.hpp:114
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
GelModule::Gel::ConstitutiveEquation::calculateResidualStrainHat
virtual PetscErrorCode calculateResidualStrainHat()
Definition: Gels.hpp:301
F
@ F
Definition: free_surface.cpp:394