Skip to content
Snippets Groups Projects
Commit 0468f826 authored by Elias Pipping's avatar Elias Pipping
Browse files

[Cleanup] Rename elasticity parameters

parent d1101c68
No related branches found
No related tags found
No related merge requests found
...@@ -110,8 +110,10 @@ int main(int argc, char *argv[]) { ...@@ -110,8 +110,10 @@ int main(int argc, char *argv[]) {
using ScalarVector = Dune::BlockVector<Dune::FieldVector<double, 1>>; using ScalarVector = Dune::BlockVector<Dune::FieldVector<double, 1>>;
using Nonlinearity = GlobalNonlinearity<Matrix, Vector>; using Nonlinearity = GlobalNonlinearity<Matrix, Vector>;
auto const E = parset.get<double>("body.E"); auto const youngModulus = parset.get<double>("body.youngModulus");
auto const nu = parset.get<double>("body.nu"); auto const poissonRatio = parset.get<double>("body.poissonRatio");
auto const shearViscosity = parset.get<double>("body.shearViscosity");
auto const bulkViscosity = parset.get<double>("body.bulkViscosity");
auto const density = parset.get<double>("body.density"); auto const density = parset.get<double>("body.density");
double const gravity = 9.81; double const gravity = 9.81;
...@@ -228,9 +230,8 @@ int main(int argc, char *argv[]) { ...@@ -228,9 +230,8 @@ int main(int argc, char *argv[]) {
MyAssembler myAssembler(leafView); MyAssembler myAssembler(leafView);
Matrix A, C, M; Matrix A, C, M;
myAssembler.assembleElasticity(E, nu, A); myAssembler.assembleElasticity(youngModulus, poissonRatio, A);
myAssembler.assembleViscosity(parset.get<double>("body.shearViscosity"), myAssembler.assembleViscosity(shearViscosity, bulkViscosity, C);
parset.get<double>("body.bulkViscosity"), C);
myAssembler.assembleMass(density, M); myAssembler.assembleMass(density, M);
EnergyNorm<Matrix, Vector> const ANorm(A), MNorm(M); EnergyNorm<Matrix, Vector> const ANorm(A), MNorm(M);
// Q: Does it make sense to weigh them in this manner? // Q: Does it make sense to weigh them in this manner?
...@@ -490,7 +491,7 @@ int main(int argc, char *argv[]) { ...@@ -490,7 +491,7 @@ int main(int argc, char *argv[]) {
myAssembler.vertexBasis, u); myAssembler.vertexBasis, u);
ScalarVector vonMisesStress; ScalarVector vonMisesStress;
VonMisesStressAssembler<Grid, P0Basis::LocalFiniteElement> VonMisesStressAssembler<Grid, P0Basis::LocalFiniteElement>
localStressAssembler(E, nu, gridDisplacement); localStressAssembler(youngModulus, poissonRatio, gridDisplacement);
p0Assembler.assembleFunctional(localStressAssembler, vonMisesStress); p0Assembler.assembleFunctional(localStressAssembler, vonMisesStress);
writeVtk(myAssembler.vertexBasis, u, alpha, p0Basis, vonMisesStress, writeVtk(myAssembler.vertexBasis, u, alpha, p0Basis, vonMisesStress,
......
...@@ -7,8 +7,8 @@ writeVTK = false ...@@ -7,8 +7,8 @@ writeVTK = false
finalTime = 15 finalTime = 15
[body] [body]
E = 5e7 youngModulus = 5e7
nu = 0.3 # The closer we get to 0.5, the more wiggly everything gets poissonRatio = 0.3 # The closer we get to 0.5, the more wiggly everything gets
density = 5000 density = 5000
height = 1 height = 1
width = 5 width = 5
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment