25#ifndef DUMUX_ELASTICPROBLEM_HH
26#define DUMUX_ELASTICPROBLEM_HH
28#include <dune/common/fmatrix.hh>
29#include <dune/grid/yaspgrid.hh>
39template <
class TypeTag>
48template<
class TypeTag>
49struct Grid<TypeTag, TTag::TestElastic> {
using type = Dune::YaspGrid<2>; };
51template<
class TypeTag>
54template<
class TypeTag>
65template<
class TypeTag>
77 using ElementVolumeVariables =
typename GridVariables::GridVolumeVariables::LocalView;
78 using ElementFluxVariablesCache =
typename GridVariables::GridFluxVariablesCache::LocalView;
81 using FVElementGeometry =
typename GridGeometry::LocalView;
82 using SubControlVolume =
typename GridGeometry::SubControlVolume;
83 using SubControlVolumeFace =
typename GridGeometry::SubControlVolumeFace;
86 using Element =
typename GridView::template Codim<0>::Entity;
87 using GlobalPosition =
typename Element::Geometry::GlobalCoordinate;
89 static constexpr Scalar pi = M_PI;
90 static constexpr int dim = GridView::dimension;
91 static constexpr int dimWorld = GridView::dimensionworld;
92 using GradU = Dune::FieldMatrix<Scalar, dim, dimWorld>;
104 {
return PrimaryVariables(0.0); }
108 {
return PrimaryVariables(0.0); }
118 BoundaryTypes values;
119 values.setAllDirichlet();
121 if (globalPos[0] < eps_ && globalPos[1] > eps_ && globalPos[1] < 1.0 - eps_)
122 values.setNeumann(Indices::momentum(0));
123 if (globalPos[1] > 1.0 - eps_ && globalPos[0] > eps_ && globalPos[0] < 1.0 - eps_)
124 values.setNeumann(Indices::momentum(1));
131 const FVElementGeometry& fvGeometry,
132 const ElementVolumeVariables& elemVolvars,
133 const ElementFluxVariablesCache& elemFluxVarsCache,
134 const SubControlVolumeFace& scvf)
const
139 for (
int i = 0; i < dim; ++i)
140 for (
int j = 0; j < dimWorld; ++j)
141 epsilon[i][j] = 0.5*(gradU[i][j] + gradU[j][i]);
143 const auto& lameParams = this->
spatialParams().lameParamsAtPos(scvf.ipGlobal());
145 const auto traceEpsilon =
trace(epsilon);
146 for (
int i = 0; i < dim; ++i)
148 sigma[i][i] = lameParams.lambda()*traceEpsilon;
149 for (
int j = 0; j < dimWorld; ++j)
150 sigma[i][j] += 2.0*lameParams.mu()*epsilon[i][j];
154 sigma.mv(scvf.unitOuterNormal(), values);
162 NumEqVector
source(
const Element& element,
163 const FVElementGeometry& fvGeometry,
164 const ElementVolumeVariables& elemVolVars,
165 const SubControlVolume& scv)
const
170 const auto ipGlobal = scv.center();
171 const auto x = ipGlobal[0];
172 const auto y = ipGlobal[1];
175 const auto& lameParams = this->
spatialParams().lameParamsAtPos(scv.center());
176 const auto lambda = lameParams.lambda();
177 const auto mu = lameParams.mu();
180 const Scalar pi_2 = 2.0*pi;
181 const Scalar pi_2_square = pi_2*pi_2;
182 const Scalar cos_2pix = cos(pi_2*x);
183 const Scalar sin_2pix = sin(pi_2*x);
184 const Scalar cos_2piy = cos(pi_2*y);
185 const Scalar sin_2piy = sin(pi_2*y);
187 const Scalar dE11_dx = -2.0*sin_2piy;
188 const Scalar dE22_dx = pi_2_square*cos_2pix*cos_2piy;
189 const Scalar dE11_dy = pi_2*(1.0-2.0*x)*cos_2piy;
190 const Scalar dE22_dy = -1.0*pi_2_square*sin_2pix*sin_2piy;
191 const Scalar dE12_dy = 0.5*pi_2_square*(cos_2pix*cos_2piy - (x-x*x)*sin_2piy);
192 const Scalar dE21_dx = 0.5*((1.0-2*x)*pi_2*cos_2piy - pi_2_square*sin_2pix*sin_2piy);
195 PrimaryVariables divSigma(0.0);
196 divSigma[Indices::momentum(0)] = lambda*(dE11_dx + dE22_dx) + 2*mu*(dE11_dx + dE12_dy);
197 divSigma[Indices::momentum(1)] = lambda*(dE11_dy + dE22_dy) + 2*mu*(dE21_dx + dE22_dy);
208 const auto x = globalPos[0];
209 const auto y = globalPos[1];
211 PrimaryVariables exact(0.0);
212 exact[Indices::momentum(0)] = (x-x*x)*sin(2*pi*y);
213 exact[Indices::momentum(1)] = sin(2*pi*x)*sin(2*pi*y);
225 const auto x = globalPos[0];
226 const auto y = globalPos[1];
228 static constexpr int xIdx = Indices::momentum(0);
229 static constexpr int yIdx = Indices::momentum(1);
231 GradU exactGrad(0.0);
232 exactGrad[xIdx][xIdx] = (1-2*x)*sin(2*pi*y);
233 exactGrad[xIdx][yIdx] = (x - x*x)*2*pi*cos(2*pi*y);
234 exactGrad[yIdx][xIdx] = 2*pi*cos(2*pi*x)*sin(2*pi*y);
235 exactGrad[yIdx][yIdx] = 2*pi*sin(2*pi*x)*cos(2*pi*y);
240 static constexpr Scalar eps_ = 3e-6;
Defines a type tag and some properties for models using the box scheme.
Dune::DenseMatrix< MatrixType >::field_type trace(const Dune::DenseMatrix< MatrixType > &M)
Trace of a dense matrix.
Definition: math.hh:760
make the local view function available whenever we use the grid geometry
Definition: adapt.hh:29
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:149
const GridGeometry & gridGeometry() const
The finite volume grid geometry.
Definition: common/fvproblem.hh:588
The DUNE grid type.
Definition: common/properties.hh:57
Property to specify the type of a problem which has to be solved.
Definition: common/properties.hh:69
The type of the spatial parameters object.
Definition: common/properties.hh:221
Base class for all geomechanical problems.
Definition: geomechanics/fvproblem.hh:68
Base class for all fully implicit porous media problems.
Definition: dumux/porousmediumflow/problem.hh:39
SpatialParams & spatialParams()
Returns the spatial parameters object.
Definition: dumux/porousmediumflow/problem.hh:146
Problem definition for the deformation of an elastic body.
Definition: test/geomechanics/elastic/problem.hh:67
NumEqVector neumann(const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolvars, const ElementFluxVariablesCache &elemFluxVarsCache, const SubControlVolumeFace &scvf) const
Evaluates the boundary conditions for a Neumann boundary segment.
Definition: test/geomechanics/elastic/problem.hh:130
PrimaryVariables dirichletAtPos(const GlobalPosition &globalPos) const
Evaluates the boundary conditions for a Dirichlet boundary segment.
Definition: test/geomechanics/elastic/problem.hh:107
NumEqVector source(const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const SubControlVolume &scv) const
Evaluates the source term for all phases within a given sub-control volume.
Definition: test/geomechanics/elastic/problem.hh:162
PrimaryVariables exactSolution(const GlobalPosition &globalPos) const
Evaluates the exact displacement to this problem at a given position.
Definition: test/geomechanics/elastic/problem.hh:204
BoundaryTypes boundaryTypesAtPos(const GlobalPosition &globalPos) const
Specifies which kind of boundary condition should be used for which equation on a given boundary segm...
Definition: test/geomechanics/elastic/problem.hh:116
static constexpr Scalar temperature()
Returns the temperature in the domain.
Definition: test/geomechanics/elastic/problem.hh:99
GradU exactGradient(const GlobalPosition &globalPos) const
Evaluates the exact displacement gradient to this problem at a given position.
Definition: test/geomechanics/elastic/problem.hh:220
PrimaryVariables initialAtPos(const GlobalPosition &globalPos) const
Evaluates the initial value for a control volume.
Definition: test/geomechanics/elastic/problem.hh:103
ElasticProblem(std::shared_ptr< const GridGeometry > gridGeometry)
Definition: test/geomechanics/elastic/problem.hh:95
Definition: test/geomechanics/elastic/problem.hh:45
std::tuple< Elastic, BoxModel > InheritsFrom
Definition: test/geomechanics/elastic/problem.hh:45
Dune::YaspGrid< 2 > type
Definition: test/geomechanics/elastic/problem.hh:49
Definition of the spatial parameters for the linear elasticity problem.
Definition: geomechanics/elastic/spatialparams.hh:41
Base class for all geomechanical problems.
Defines a type tag and some properties for the elastic geomechanical model.
Definition of the spatial parameters for the MaxwellStefan problem.