Skip to content

GitLab

Commit 879cea41 authored by Praetorius, Simon's avatar Praetorius, Simon
Browse files

Merge branch 'example/stokes' into 'master' 

Cleanup of examples for stokes equation

See merge request !17
parents 2494d6df b3062bde
Hide whitespace changes
Inline Side-by-side
examples/init/ellipt.json.2d deleted 100644 → 0
View file @ 2494d6df
{
"dimension of world": 2,
"elliptMesh": {
"macro file name": "./macro/macro.stand.2d",
"global refinements": 5
},
"ellipt": {
"mesh": "elliptMesh",
"names": "u",
"solver" : {
"name": "cg",
"max iteration": 1000,
"absolute tolerance": 1e-6,
"info": 1,
"left precon": "diag"
},
"output[0]": {
"filename": "ellipt.2d",
"name": "u",
"output directory": "output"
}
}
}
examples/init/ellipt.json.3d deleted 100644 → 0
View file @ 2494d6df
{
"dimension of world": 3,
"elliptMesh": {
"macro file name": "./macro/BG1_withBoundary.3d",
"global refinements": 0,
"min corner": "0,0,0",
"max corner": "1,1,1",
"num cells": "2,2,2",
"dimension": "1,1,1"
},
"ellipt": {
"mesh": "elliptMesh",
"names": "u",
"solver" : {
"name": "cg",
"max iteration": 1000,
"absolute tolerance": 1e-6,
"info": 1,
"left precon": "diag"
},
"output": {
"filename": "ellipt.3d",
"output directory": "output"
}
}
}
examples/init/heat.json.2d deleted 100644 → 0
View file @ 2494d6df
{
"dimension of world": 2,
"heatMesh": {
"macro file name": "./macro/macro.stand.2d",
"global refinements": 4,
"min corner": "0,0",
"max corner": "1,1",
"num cells": "2,2",
"dimension": "2,2"
},
"heat": {
"mesh": "heatMesh",
"names": "u",
"solver" : {
"name": "cg",
"max iteration": 1000,
"absolute tolerance": 1e-6,
"break if tolerance not reached": 1,
"info": 1,
"left precon": "diag"
},
"output": {
"filename": "heat.2d",
"output directory": "output",
"ParaView format": 1,
"ParaView mode": 1
}
},
"adapt": {
"timestep": 0.1,
"start time": 0.0,
"end time": 1.0
}
}
examples/init/heat.json.3d deleted 100644 → 0
View file @ 2494d6df
{
"dimension of world": 3,
"heatMesh": {
"macro file name": "./macro/macro.stand.3d",
"global refinements": 2
},
"heat": {
"mesh": "heatMesh",
"names": "u",
"solver" : {
"name": "cg",
"max iteration": 1000,
"absolute tolerance": 1e-6,
"info": 1,
"left precon": "diag"
},
"output": {
"filename": "heat.3d",
"output directory": "output"
}
},
"adapt": {
"timestep": 0.1,
"start time": 0.0,
"end time": 0.2
}
}
examples/init/navier_stokes.json.2d deleted 100644 → 0
View file @ 2494d6df
{
"dimension of world": 2,
"stokesMesh": {
"macro file name": "./macro/macro.small2.2d",
"global refinements": 4,
"min corner": "0,0",
"max corner": "1,1",
"num cells": "2,2"
},
"stokes": {
"mesh": "stokesMesh",
"names": ["u0", "u1", "p"],
"solver" : {
"name": "direct",
"ell": 3,
"max iteration": 5000,
"relative tolerance": 1e-6,
"info": 1
},
"output": {
"filename": "navier_stokes.2d",
"output directory": "output",
"ParaView format": 1
},
"viscosity": 0.1,
"density": 1.0,
"boundary velocity": 1.0
},
"adapt": {
"timestep": 1.0e-3,
"start time": 0.0,
"end time": 1.0
}
}
examples/init/pfc.json.2d deleted 100644 → 0
View file @ 2494d6df
{
"dimension of world": 2,
"pfcMesh": {
"macro file name": "./macro/macro.square.2d",
"global refinements": 8,
"min corner": "0,0",
"max corner": "1,1",
"num cells": "2,2",
"dimension": "50,50"
},
"pfc": {
"mesh": "pfcMesh",
"solver" : {
"name": "fgmres",
"max iteration": 1000,
"relative tolerance": 1e-6,
"info": 10,
"ell": 3,
"print cycle": 10,
"right precon": "pfc"
},
"output": {
"filename": "pfc.2d",
"output directory": "output"
}
},
"precon_pfc": {
"M": {
"solver": "cg",
"max iteration": 1000,
"relative tolerance": 1e-3,
"left precon": "diag",
"info": 0
},
"MpL": {
"solver": "cg",
"max iteration": 1000,
"relative tolerance": 1e-3,
"left precon": "diag",
"info": 0
},
"MpL2": {
"solver": "cg",
"max iteration": 1000,
"relative tolerance": 1e-3,
"left precon": "diag",
"info": 0
}
},
"adapt": {
"timestep": 0.1,
"start time": 0.0,
"end time": 10.0
}
}
examples/init/stokes.dat.2d
View file @ 879cea41
dimension of world: 2
stokesMesh->global refinements: 0
stokesMesh->dimension: 1.0 1.0
stokesMesh->num cells: 4 4
stokes->mesh: stokesMesh
stokes->solver->name: direct
stokes->solver->name: bicgstag
stokes->solver->ell: 3
stokes->solver->max iteration: 1000
stokes->solver->relative tolerance: 1e-8
stokes->solver->restart: 50
stokes->solver->precon->name: ilu
stokes->solver->info: 2
stokes->output[0]->filename: stokes_u.2d
......
examples/init/stokes.json.2d deleted 100644 → 0
View file @ 2494d6df
{
"dimension of world": 2,
"stokesMesh": {
"macro file name": "./macro/macro.stand.2d",
"global refinements": 2
},
"stokes": {
"mesh": "stokesMesh",
"solver" : {
"name": "minres",
"max iteration": 2000,
"reduction": 1e-8,
"restart": 50,
"info": 2
},
"output[0]": {
"filename": "stokes_u.2d",
"name": "u",
"subsampling": 1,
"output directory": "output",
"ParaView format": 1,
"ParaView mode": 1
},
"output[1]": {
"filename": "stokes_p.2d",
"name": "p",
"output directory": "output",
"ParaView format": 1,
"ParaView mode": 1
}
}
}
examples/init/test.json.2d deleted 100644 → 0
View file @ 2494d6df
{
"dimension of world": 2,
"testMesh": {
"macro file name": "macro/macro.stand.2d",
"global refinements": 10
},
"test": {
"mesh": "testMesh",
"dim": 2,
"components": 1,
"polynomial degree[0]": 1,
"solver": "cg",
"solver" : {
"max iteration": 1000,
"tolerance": 1e-8,
"info": 10,
"left precon": "diag"
},
"output": {
"filename": "output/test.2d",
"ParaView format": 1,
"ParaView mode": 1
}
}
}
examples/init/vecellipt.json.2d deleted 100644 → 0
View file @ 2494d6df
{
"dimension of world": 2,
"elliptMesh": {
"macro file name": "./macro/macro.stand.2d",
"global refinements": 5
},
"ellipt": {
"mesh": "elliptMesh",
"names": "u",
"solver" : {
"name": "direct",
"max iteration": 1000,
"absolute tolerance": 1e-6,
"info": 1,
"left precon": "diag"
},
"output[0]": {
"filename": "vecellipt0.2d",
"name": "u0",
"output directory": "output"
},
"output[1]": {
"filename": "vecellipt1.2d",
"name": "u1",
"output directory": "output"
}
}
}
examples/navier_stokes.cc
View file @ 879cea41
#ifdef HAVE_CONFIG_H
#include "config.h"
#endif
#include <iostream>
#include <ctime>
#include <cmath>
#include <amdis/AMDiS.hpp>
#include <amdis/common/FieldMatVec.hpp>
#include <amdis/AdaptInstationary.hpp>
#include <amdis/LocalOperators.hpp>
#include <amdis/ProblemStat.hpp>
......@@ -16,23 +7,17 @@
using namespace AMDiS;
struct NavierStokesBasis
{
using Grid = Dune::YaspGrid<GRIDDIM>;
using GlobalBasis = typename TaylorHoodBasis<Grid::LeafGridView>::GlobalBasis;
};
using StokesProblem = ProblemStat<NavierStokesBasis>;
using StokesProblemInstat = ProblemInstat<NavierStokesBasis>;
using Grid = Dune::YaspGrid<2>;
using Basis = TaylorHoodBasis<Grid::LeafGridView>;
int main(int argc, char** argv)
{
AMDiS::init(argc, argv);
StokesProblem prob("stokes");
ProblemStat<Basis> prob("stokes");
prob.initialize(INIT_ALL);
StokesProblemInstat probInstat("stokes", prob);
ProblemInstat<Basis> probInstat("stokes", prob);
probInstat.initialize(INIT_UH_OLD);
double viscosity = 1.0;
......@@ -82,31 +67,21 @@ int main(int argc, char** argv)
trans(gradientAtQP(prob.solution(_v))) * prob.solution(_v));
prob.addVectorOperator(opNonlin3, _v);
// define boundary regions
auto left = [](auto const& x) { return x[0] < 1.e-8; };
auto not_left = [](auto const& x) { return x[0] > 1.0 - 1.e-8 || x[1] < 1.e-8 || x[1] > 1.0 - 1.e-8; };
// define boundary values
auto parabolic_y = [](auto const& x) -> Dune::FieldVector<double,2>
auto parabolic_y = [](auto const& x)
{
return {0.0, x[1]*(1.0 - x[1])};
return FieldVector<double,2>{0.0, x[1]*(1.0 - x[1])};
};
auto zero = [](auto const& x) -> Dune::FieldVector<double,2>
{
return {0.0, 0.0};
};
FieldVector<double,2> zero(0);
// set boundary conditions for velocity
prob.addDirichletBC(left, _v, _v, parabolic_y);
prob.addDirichletBC(not_left, _v, _v, zero);
// set point constraint for pressure
prob.addDirichletBC([](auto const& x) { return x[0] < 1.e-8 && x[1] < 1.e-8; }, _p, _p, 0.0);
prob.boundaryManager()->setBoxBoundary({1,2,2,2});
prob.addDirichletBC(1, _v, _v, parabolic_y);
prob.addDirichletBC(2, _v, _v, zero);
// set initial conditions
prob.solution(_v).interpolate(parabolic_y);
prob.solution(_p).interpolate(0.0);
AdaptInstationary adapt("adapt", prob, adaptInfo, probInstat, adaptInfo);
adapt.adapt();
......
examples/stokes0.cc
View file @ 879cea41
#ifdef HAVE_CONFIG_H
#include "config.h"
#endif
#include <iostream>
#include <ctime>
#include <cmath>
#include <amdis/AMDiS.hpp>
#include <amdis/LocalOperators.hpp>
#include <amdis/ProblemStat.hpp>
#ifdef DOW
#undef DOW
#endif
#define DOW WORLDDIM
using namespace AMDiS;
// 3 components: velocity with polynomial degree 2 and pressure with polynomial degree 1
using Grid = Dune::YaspGrid<GRIDDIM>;
using StokesParam = TaylorHoodBasis<Grid::LeafGridView>;
using Grid = Dune::YaspGrid<2>;
using StokesParam = TaylorHoodBasis<Grid::LeafGridView>;
using StokesProblem = ProblemStat<StokesParam>;
int main(int argc, char** argv)
......@@ -37,8 +22,8 @@ int main(int argc, char** argv)
auto _v = Dune::Indices::_0;
auto _p = Dune::Indices::_1;
// <viscosity*grad(u_i), grad(v_i)>
for (std::size_t i = 0; i < DOW; ++i) {
for (std::size_t i = 0; i < WORLDDIM; ++i) {
// <viscosity*grad(u_i), grad(v_i)>
auto opL = makeOperator(tag::gradtest_gradtrial{}, viscosity);
prob.addMatrixOperator(opL, treepath(_v,i), treepath(_v,i));
......@@ -51,42 +36,23 @@ int main(int argc, char** argv)
prob.addMatrixOperator(opDiv, _p, treepath(_v,i));
}
auto opZero = makeOperator(tag::test_trial{}, 0.0);
prob.addMatrixOperator(opZero, _p, _p);
// define boundary regions
auto left = [](auto const& x) { return x[0] < 1.e-8; };
auto not_left = [](auto const& x) { return x[0] > 1.0 - 1.e-8 || x[1] < 1.e-8 || x[1] > 1.0 - 1.e-8; };
// define boundary values
auto parabolic_y = [](auto const& x) -> Dune::FieldVector<double,DOW>
auto parabolic_y = [](auto const& x)
{
return {0.0, x[1]*(1.0 - x[1])};
return FieldVector<double,2>{0.0, x[1]*(1.0 - x[1])};
};
auto zero = [](auto const& x) -> Dune::FieldVector<double,DOW>
{
return {0.0, 0.0};
};
FieldVector<double,2> zero(0);
// set boundary conditions for velocity
prob.addDirichletBC(left, _v, _v, parabolic_y);
prob.addDirichletBC(not_left, _v, _v, zero);
prob.addDirichletBC([](auto const& x) { return x[0] < 1.e-8 && x[1] < 1.e-8; }, _p, _p, 0.0);
prob.boundaryManager()->setBoxBoundary({1,2,2,2});
prob.addDirichletBC(1, _v, _v, parabolic_y);
prob.addDirichletBC(2, _v, _v, zero);
AdaptInfo adaptInfo("adapt");
// assemble and solve system
prob.assemble(adaptInfo);
#ifdef DEBUG_MTL
// write matrix to file
mtl::io::matrix_market_ostream out("matrix_stokes0.mtx");
out << prob.systemMatrix().matrix();
std::cout << prob.systemMatrix().matrix() << '\n';
#endif
prob.solve(adaptInfo);
// output solution
......
examples/stokes1.cc
View file @ 879cea41
#ifdef HAVE_CONFIG_H
#include "config.h"
#endif
#include <iostream>
#include <ctime>
#include <cmath>
#include <amdis/AMDiS.hpp>
#include <amdis/LocalOperators.hpp>
#include <amdis/ProblemStat.hpp>
#ifdef DOW
#undef DOW
#endif
#define DOW WORLDDIM
using namespace AMDiS;
// 3 components: velocity with polynomial degree 2 and pressure with polynomial degree 1
using Grid = Dune::YaspGrid<GRIDDIM>;
using StokesParam = TaylorHoodBasis<Grid::LeafGridView>;
using Grid = Dune::YaspGrid<2>;
using StokesParam = TaylorHoodBasis<Grid::LeafGridView>;
using StokesProblem = ProblemStat<StokesParam>;
int main(int argc, char** argv)
......@@ -37,8 +22,8 @@ int main(int argc, char** argv)
auto _v = Dune::Indices::_0;
auto _p = Dune::Indices::_1;
// <viscosity*grad(u_i), grad(v_i)>
for (std::size_t i = 0; i < DOW; ++i) {
for (std::size_t i = 0; i < WORLDDIM; ++i) {
// <viscosity*grad(u_i), grad(v_i)>
auto opL = makeOperator(tag::gradtest_gradtrial{}, viscosity);
prob.addMatrixOperator(opL, treepath(_v,i), treepath(_v,i));
}
......@@ -51,43 +36,23 @@ int main(int argc, char** argv)
auto opDiv = makeOperator(tag::test_divtrialvec{}, 1.0);
prob.addMatrixOperator(opDiv, _p, _v);
auto opZero = makeOperator(tag::test_trial{}, 0.0);
prob.addMatrixOperator(opZero, _p, _p);
// define boundary regions
auto left = [](auto const& x) { return x[0] < 1.e-8; };
auto not_left = [](auto const& x) { return x[0] > 1.0 - 1.e-8 || x[1] < 1.e-8 || x[1] > 1.0 - 1.e-8; };
// define boundary values
auto parabolic_y = [](auto const& x) -> Dune::FieldVector<double,DOW>
auto parabolic_y = [](auto const& x)
{
return {0.0, x[1]*(1.0 - x[1])};
return FieldVector<double,2>{0.0, x[1]*(1.0 - x[1])};
};
auto zero = [](auto const& x) -> Dune::FieldVector<double,DOW>
{
return {0.0, 0.0};
};
FieldVector<double,2> zero(0);
// set boundary conditions for velocity
prob.addDirichletBC(left, _v, _v, parabolic_y);
prob.addDirichletBC(not_left, _v, _v, zero);
prob.addDirichletBC([](auto const& x) { return x[0] < 1.e-8 && x[1] < 1.e-8; }, _p, _p, 0.0);
prob.boundaryManager()->setBoxBoundary({1,2,2,2});
prob.addDirichletBC(1, _v, _v, parabolic_y);
prob.addDirichletBC(2, _v, _v, zero);
AdaptInfo adaptInfo("adapt");
// assemble and solve system
prob.assemble(adaptInfo);
#ifdef DEBUG_MTL
// write matrix to file
mtl::io::matrix_market_ostream out("matrix_stokes1.mtx");
out << prob.systemMatrix().matrix();
std::cout << prob.systemMatrix().matrix() << '\n';
#endif
prob.solve(adaptInfo);
// output solution
......
examples/stokes3.cc
View file @ 879cea41
#ifdef HAVE_CONFIG_H
#include "config.h"
#endif
#include <iostream>
#include <ctime>
#include <cmath>
#include <amdis/AMDiS.hpp>
#include <amdis/LocalOperators.hpp>
#include <amdis/ProblemStat.hpp>
......@@ -13,8 +5,8 @@