Commit a4c463df authored by Thomas Witkowski's avatar Thomas Witkowski
Browse files

And some more work on demos and tutorial.

parent 4b9e369e
......@@ -28,12 +28,9 @@
namespace AMDiS {
/// Contains boundary constants
/// Flag to denote interior boundaryies.
typedef enum {
INTERIOR = 0, /**< interior => no boundary (b = 0) */
DIRICHLET = 1, /**< dirichlet boundary (b > 0) */
NEUMANN = -1, /**< neumann boundary (-100 < b < 0) */
ROBIN = -100 /**< robin boundary (b <= -100) */
INTERIOR = 0
} BoundaryConstants;
/// Type specifier for the different boundary types
......
......@@ -491,12 +491,16 @@ namespace AMDiS {
void MacroInfo::dirichletBoundary()
{
// === Traverse all elements and set either interior boundaries, if the ===
// === element has a neighbour on this bound. Otherwise set some arbitrary ===
// === Dirichlet boundary conditions. ===
for (int i = 0; i < static_cast<int>(mel.size()); i++)
for (int k = 0; k < mesh->getGeo(NEIGH); k++)
if (mel[i]->neighbour[k])
mel[i]->boundary[k] = INTERIOR;
else
mel[i]->boundary[k] = DIRICHLET;
mel[i]->boundary[k] = 1;
}
......@@ -521,13 +525,17 @@ namespace AMDiS {
++melIt, i++) {
for (int j = 0; j < mesh->getGeo(NEIGH); j++) {
if ((*melIt)->getBoundary(j) != INTERIOR) {
if ((*melIt)->getBoundary(j) >= DIRICHLET) {
if ((*melIt)->getBoundary(j) > 0) {
// Here we have some Dirichlet boundary conditions.
int j1 = mel_vertex[i][(j + 1) % 3];
int j2 = mel_vertex[i][(j + 2) % 3];
bound[j1] = std::max(bound[j1], (*melIt)->getBoundary(j));
bound[j2] = std::max(bound[j2], (*melIt)->getBoundary(j));
} else if ((*melIt)->getBoundary(j) <= NEUMANN) {
} else {
// Otherwise we have some Neumann boundary conditions.
int j1 = mel_vertex[i][(j + 1) % 3];
int j2 = mel_vertex[i][(j + 2) % 3];
......
......@@ -26,16 +26,7 @@ ball->adapt->max iteration: 4
ball->adapt->info: 8
ball->output->filename: ball
ball->output->ParaView format: 1
ball->output->TecPlot format: 0
ball->output->TecPlot ext: .tec
ball->output->AMDiS format: 0
ball->output->AMDiS mesh ext: .mesh
ball->output->AMDiS data ext: .dat
ball->output->append index: 0
ball->output->index length: 6
ball->output->index decimals: 3
......
......@@ -16,10 +16,6 @@ ball->estimator: 0
ball->marker->strategy: 0 % 0: no adaption 1: GR 2: MS 3: ES 4:GERS
ball->output->filename: output/ball
ball->output->ParaView format: 1
ball->output->AMDiS format: 0
ball->output->AMDiS mesh ext: .mesh
ball->output->AMDiS data ext: .dat
......@@ -28,16 +28,7 @@ bunny->adapt->max iteration: 0
bunny->adapt->info: 8
bunny->output->filename: bunny_fixed
bunny->output->ParaView format: 1
bunny->output->TecPlot format: 0
bunny->output->TecPlot ext: .tec
bunny->output->AMDiS format: 0
bunny->output->AMDiS mesh ext: .mesh
bunny->output->AMDiS data ext: .dat
bunny->output->append index: 0
bunny->output->index length: 6
bunny->output->index decimals: 3
......
......@@ -25,9 +25,4 @@ neumann->adapt->max iteration: 100
neumann->adapt->refine bisections: 2
neumann->output->filename: output/neumann
neumann->output->ParaView format: 1
neumann->output->AMDiS format: 0
neumann->output->AMDiS mesh ext: .mesh
neumann->output->AMDiS data ext: .dat
......@@ -27,16 +27,7 @@ parametric->adapt->max iteration: 8
parametric->adapt->info: 8
parametric->output->filename: parametric
parametric->output->ParaView format: 1
parametric->output->TecPlot format: 0
parametric->output->TecPlot ext: .tec
parametric->output->AMDiS format: 0
parametric->output->AMDiS mesh ext: .mesh
parametric->output->AMDiS data ext: .dat
parametric->output->append index: 0
parametric->output->index length: 6
parametric->output->index decimals: 3
......
dimension of world: 1
periodicMesh->macro file name: ./macro/periodic.macro.1d
periodicMesh->global refinements: 0
periodicMesh->periodic file: ./init/periodic.per.1d
periodicMesh->preserve coarse dofs: 0
periodicMesh->global refinements: 0
periodic->mesh: periodicMesh
periodic->dim: 1
......@@ -32,16 +29,7 @@ periodic->adapt->refine bisections: 2
periodic->adapt->coarsen allowed: 0
periodic->output->filename: periodic
periodic->output->ParaView format: 1
periodic->output->TecPlot format: 0
periodic->output->TecPlot ext: .tec
periodic->output->AMDiS format: 0
periodic->output->AMDiS mesh ext: .mesh
periodic->output->AMDiS data ext: .dat
periodic->output->append index: 0
periodic->output->index length: 6
periodic->output->index decimals: 3
......
dimension of world: 2
periodicMesh->macro file name: ./macro/periodic.macro.2d
periodicMesh->global refinements: 0
periodicMesh->periodic file: ./init/periodic.per.2d
periodicMesh->preserve coarse dofs: 0
periodicMesh->global refinements: 0
periodic->mesh: periodicMesh
periodic->dim: 2
......@@ -32,17 +29,7 @@ periodic->adapt->refine bisections: 2
periodic->adapt->coarsen allowed: 0
periodic->output->filename: output/periodic
periodic->output->ParaView format: 1
periodic->output->TecPlot format: 0
periodic->output->TecPlot ext: .tec
periodic->output->AMDiS format: 0
periodic->output->AMDiS mesh ext: .mesh
periodic->output->AMDiS data ext: .dat
periodic->output->append index: 0
periodic->output->index length: 6
periodic->output->index decimals: 3
......
dimension of world: 3
periodicMesh->macro file name: ./macro/periodic.macro.3d
periodicMesh->global refinements: 0
periodicMesh->periodic file: ./init/periodic.per.3d
periodicMesh->preserve coarse dofs: 0
periodicMesh->global refinements: 0
periodic->mesh: periodicMesh
periodic->dim: 3
......@@ -32,16 +29,7 @@ periodic->adapt->refine bisections: 3
periodic->adapt->coarsen allowed: 0
periodic->output->filename: periodic
periodic->output->ParaView format: 1
periodic->output->TecPlot format: 0
periodic->output->TecPlot ext: .tec
periodic->output->AMDiS format: 0
periodic->output->AMDiS mesh ext: .mesh
periodic->output->AMDiS data ext: .dat
periodic->output->append index: 0
periodic->output->index length: 6
periodic->output->index decimals: 3
......
......@@ -17,8 +17,5 @@ sphere->marker->strategy: 0
sphere->output->filename: output/sphere
sphere->output->ParaView format: 1
sphere->output->AMDiS format: 0
sphere->output->AMDiS mesh ext: .mesh
sphere->output->AMDiS data ext: .dat
......@@ -16,7 +16,4 @@ torus->marker->strategy: 0
torus->output->filename: output/torus
torus->output->ParaView format: 1
torus->output->AMDiS format: 0
torus->output->AMDiS mesh ext: .mesh
torus->output->AMDiS data ext: .dat
......@@ -5,7 +5,7 @@ vecheatMesh->global refinements: 0
vecheat->space->polynomial degree[0]: 1
vecheat->space->polynomial degree[1]: 1
vecheat->space->dim: 1
vecheat->space->mesh: vecheatMesh
vecheat->space->components: 2
......@@ -32,22 +32,12 @@ vecheat->theta: 1.0
vecheat->adapt[0]->tolerance: 0.01
vecheat->adapt[1]->tolerance: 0.01
vecheat->adapt->timestep: 0.01
vecheat->adapt->max iteration: 2
vecheat->adapt->refine bisections: 1
vecheat->adapt->coarsen bisections: 1
vecheat->adapt[0]->time tolerance: 0.01
vecheat->adapt[1]->time tolerance: 0.01
vecheat->adapt[0]->rel space error: 0.5
vecheat->adapt[0]->rel time error: 0.5
vecheat->adapt[0]->info: 8
vecheat->adapt[0]->coarsen allowed: 1 % 0|1
vecheat->adapt[1]->rel space error: 0.5
vecheat->adapt[1]->rel time error: 0.5
vecheat->adapt[1]->info: 8
vecheat->adapt[1]->coarsen allowed: 1 % 0|1
vecheat->adapt->timestep: 0.01
vecheat->adapt->start time: 0.0
vecheat->adapt->end time: 1.0
vecheat->adapt->strategy: 1 % 0=explicit, 1=implicit
vecheat->adapt->max iteration: 1
......@@ -58,8 +48,6 @@ vecheat->initial->marker[1]->strategy: 2 % 0=none, 1=GR, 2=MS, 3=ES, 4=
vecheat->initial->adapt->max iteration: 10
vecheat->initial->adapt->info: 8
vecheat->space->dim: 1
vecheat->space->marker[0]->strategy: 3 % 0=none, 1=GR, 2=MS, 3=ES, 4=GERS
vecheat->space->marker[0]->ESTheta: 0.9
vecheat->space->marker[0]->ESThetaC: 0.05
......@@ -71,31 +59,15 @@ vecheat->space->marker[1]->ESThetaC: 0.05
vecheat->space->marker[1]->info: 8
vecheat->space->output[0]->filename: vecheat0_
vecheat->space->output[0]->ParaView format: 0
vecheat->space->output[0]->TecPlot format: 0
vecheat->space->output[0]->TecPlot ext: .tec
vecheat->space->output[0]->AMDiS format: 0
vecheat->space->output[0]->AMDiS mesh ext: .mesh
vecheat->space->output[0]->AMDiS data ext: .dat
vecheat->space->output[0]->ParaView format: 1
vecheat->space->output[0]->ParaView animation: 1
vecheat->space->output[0]->append index: 1
vecheat->space->output[0]->index length: 6
vecheat->space->output[0]->index decimals: 3
vecheat->space->output[1]->filename: vecheat1_
vecheat->space->output[1]->ParaView format: 0
vecheat->space->output[1]->TecPlot format: 0
vecheat->space->output[1]->TecPlot ext: .tec
vecheat->space->output[1]->AMDiS format: 0
vecheat->space->output[1]->AMDiS mesh ext: .mesh
vecheat->space->output[1]->AMDiS data ext: .dat
vecheat->space->output[1]->ParaView format: 1
vecheat->space->output[1]->ParaView animation: 1
vecheat->space->output[1]->append index: 1
vecheat->space->output[1]->index length: 6
vecheat->space->output[1]->index decimals: 3
......
......@@ -5,7 +5,7 @@ vecheatMesh->global refinements: 5
vecheat->space->polynomial degree[0]: 1
vecheat->space->polynomial degree[1]: 1
vecheat->space->dim: 2
vecheat->space->mesh: vecheatMesh
vecheat->space->components: 2
......@@ -34,19 +34,13 @@ vecheat->adapt->coarsen bisections: 2
vecheat->adapt[0]->tolerance: 0.001
vecheat->adapt[0]->tolerance: 0.001
vecheat->adapt[0]->time tolerance: 0.1
vecheat->adapt[1]->time tolerance: 0.1
vecheat->adapt->timestep: 0.01
vecheat->adapt->start time: 0.0
vecheat->adapt->end time: 1.0
vecheat->adapt->max iteration: 2
vecheat->adapt[0]->info: 8
vecheat->adapt[0]->coarsen allowed: 1 % 0|1
vecheat->adapt[1]->info: 8
vecheat->adapt[1]->coarsen allowed: 1 % 0|1
vecheat->adapt->strategy: 1 % 0=explicit, 1=implicit
vecheat->adapt->max iteration: 1
......@@ -57,8 +51,6 @@ vecheat->initial->marker[1]->strategy: 2 % 0=none, 1=GR, 2=MS, 3=ES, 4=
vecheat->initial->adapt->max iteration: 10
vecheat->initial->adapt->info: 8
vecheat->space->dim: 2
vecheat->space->marker[0]->strategy: 3 % 0=none, 1=GR, 2=MS, 3=ES, 4=GERS
vecheat->space->marker[0]->ESTheta: 0.9
vecheat->space->marker[0]->ESThetaC: 0.05
......@@ -69,11 +61,20 @@ vecheat->space->marker[1]->ESTheta: 0.9
vecheat->space->marker[1]->ESThetaC: 0.05
vecheat->space->marker[1]->info: 8
vecheat->space->output->filename: vecheat_
vecheat->space->output->ParaView format: 1
vecheat->space->output->append index: 1
vecheat->space->output->index length: 6
vecheat->space->output->index decimals: 3
vecheat->space->output[0]->filename: vecheat0_
vecheat->space->output[0]->ParaView format: 1
vecheat->space->output[0]->ParaView animation: 1
vecheat->space->output[0]->append index: 1
vecheat->space->output[0]->index length: 6
vecheat->space->output[0]->index decimals: 3
vecheat->space->output[1]->filename: vecheat1_
vecheat->space->output[1]->ParaView format: 1
vecheat->space->output[1]->ParaView animation: 1
vecheat->space->output[1]->append index: 1
vecheat->space->output[1]->index length: 6
vecheat->space->output[1]->index decimals: 3
WAIT: 0
......
dimension of world: 3
vecheatMesh->macro file name: ./macro/macro.stand.3d
vecheatMesh->global refinements: 0
vecheatMesh->global refinements: 3
vecheat->space->polynomial degree[0]: 1
vecheat->space->polynomial degree[1]: 1
vecheat->space->dim: 3
vecheat->space->mesh: vecheatMesh
vecheat->space->components: 2
vecheat->space->solver: cg
vecheat->space->solver: umfpack
vecheat->space->solver->max iteration: 1000
vecheat->space->solver->tolerance: 1.e-8
vecheat->space->solver->info: 8
vecheat->space->solver->left precon: diag
vecheat->space->solver->info: 2
vecheat->space->solver->left precon: no
vecheat->space->solver->right precon: no
vecheat->space->estimator[0]: residual
......@@ -29,25 +29,18 @@ vecheat->space->estimator[1]->C3: 1.0
vecheat->theta: 1.0
vecheat->adapt->refine bisections: 3
vecheat->adapt->coarsen bisections: 3
vecheat->adapt->refine bisections: 2
vecheat->adapt->coarsen bisections: 2
vecheat->adapt[0]->tolerance: 0.01
vecheat->adapt[1]->tolerance: 0.01
vecheat->adapt[0]->tolerance: 0.001
vecheat->adapt[0]->tolerance: 0.001
vecheat->adapt[0]->time tolerance: 0.1
vecheat->adapt[1]->time tolerance: 0.1
vecheat->adapt->timestep: 0.01
vecheat->adapt->start time: 0.0
vecheat->adapt->end time: 1.0
vecheat->adapt->max iteration: 2
vecheat->adapt[0]->rel space error: 0.5
vecheat->adapt[0]->rel time error: 0.5
vecheat->adapt[0]->info: 8
vecheat->adapt[0]->coarsen allowed: 1 % 0|1
vecheat->adapt[1]->rel space error: 0.5
vecheat->adapt[1]->rel time error: 0.5
vecheat->adapt[1]->info: 8
vecheat->adapt[1]->coarsen allowed: 1 % 0|1
vecheat->adapt->strategy: 1 % 0=explicit, 1=implicit
vecheat->adapt->max iteration: 1
......@@ -58,8 +51,6 @@ vecheat->initial->marker[1]->strategy: 2 % 0=none, 1=GR, 2=MS, 3=ES, 4=
vecheat->initial->adapt->max iteration: 10
vecheat->initial->adapt->info: 8
vecheat->space->dim: 3
vecheat->space->marker[0]->strategy: 3 % 0=none, 1=GR, 2=MS, 3=ES, 4=GERS
vecheat->space->marker[0]->ESTheta: 0.9
vecheat->space->marker[0]->ESThetaC: 0.05
......@@ -71,35 +62,20 @@ vecheat->space->marker[1]->ESThetaC: 0.05
vecheat->space->marker[1]->info: 8
vecheat->space->output[0]->filename: vecheat0_
vecheat->space->output[0]->ParaView format: 0
vecheat->space->output[0]->TecPlot format: 0
vecheat->space->output[0]->TecPlot ext: .tec
vecheat->space->output[0]->AMDiS format: 0
vecheat->space->output[0]->AMDiS mesh ext: .mesh
vecheat->space->output[0]->AMDiS data ext: .dat
vecheat->space->output[0]->ParaView format: 1
vecheat->space->output[0]->ParaView animation: 1
vecheat->space->output[0]->append index: 1
vecheat->space->output[0]->index length: 6
vecheat->space->output[0]->index decimals: 3
vecheat->space->output[1]->filename: vecheat1_
vecheat->space->output[1]->ParaView format: 0
vecheat->space->output[1]->TecPlot format: 0
vecheat->space->output[1]->TecPlot ext: .tec
vecheat->space->output[1]->AMDiS format: 0
vecheat->space->output[1]->AMDiS mesh ext: .mesh
vecheat->space->output[1]->AMDiS data ext: .dat
vecheat->space->output[1]->ParaView format: 1
vecheat->space->output[1]->ParaView animation: 1
vecheat->space->output[1]->append index: 1
vecheat->space->output[1]->index length: 6
vecheat->space->output[1]->index decimals: 3
WAIT: 0
......
......@@ -174,13 +174,6 @@ int main(int argc, char** argv)
adaptInfoInitial);
// ===== create boundary functions =====
G *boundaryFct = new G;
boundaryFct->setTimePtr(heat.getTime());
heat.setExactSolution(boundaryFct);
heatSpace.addDirichletBC(1, boundaryFct);
// ===== create rhs functions =====
int degree = heatSpace.getFeSpace()->getBasisFcts()->getDegree();
F *rhsFct = new F(degree);
......@@ -213,6 +206,13 @@ int main(int argc, char** argv)
heatSpace.addVectorOperator(F);
// ===== create boundary functions =====
G *boundaryFct = new G;
boundaryFct->setTimePtr(heat.getTime());
heat.setExactSolution(boundaryFct);
heatSpace.addDirichletBC(1, boundaryFct);
// ===== start adaption loop =====
int errorCode = adaptInstat.adapt();
......
......@@ -30,7 +30,7 @@ public:
/// Implementation of AbstractFunction::operator().
double operator()(const WorldVector<double>& x) const
{
int dim = x.getSize();
int dim = Global::getGeo(WORLD);
double r2 = x * x;
double ux = sin(M_PI * (*timePtr)) * exp(-10.0 * r2);
double ut = M_PI * cos(M_PI * (*timePtr)) * exp(-10.0 * r2);
......@@ -68,9 +68,8 @@ public:
/// set the time in all needed functions!
void setTime(AdaptInfo *adaptInfo)
{
ProblemInstat::setTime(adaptInfo);
rhsTime = adaptInfo->getTime() - (1-theta) * adaptInfo->getTimestep();
boundaryTime = adaptInfo->getTime();
tau1 = 1.0 / adaptInfo->getTimestep();
}
// ===== initial problem methods =====================================
......@@ -79,7 +78,7 @@ public:
void solveInitialProblem(AdaptInfo *adaptInfo)
{
int size = problemStat->getNumComponents();
boundaryTime = rhsTime = 0.0;
rhsTime = 0.0;
for (int i = 0; i < size; i++) {
problemStat->getMesh(i)->dofCompress();
problemStat->getSolution()->getDOFVector(i)->interpol(boundaryFct);
......@@ -91,7 +90,6 @@ public:
{
int size = problemStat->getNumComponents();
double errMax, errSum;
boundaryTime = 0.0;
for (int i = 0; i < size; i++) {
errSum = Error<double>::L2Err(*boundaryFct,
......@@ -157,24 +155,12 @@ public:
return &theta1;
}
/// Returns pointer to \ref tau1
double *getTau1Ptr()
{
return &tau1;
}
/// Returns pointer to \ref rhsTime.
double *getRHSTimePtr()
double *getRhsTimePtr()
{
return &rhsTime;
}
/// Returns pointer to \ref theta1.
double *getBoundaryTimePtr()
{
return &boundaryTime;
}
/// Returns \ref boundaryFct;
AbstractFunction<double, WorldVector<double> > *getBoundaryFct()
{
......@@ -188,15 +174,9 @@ protected:
/// theta - 1
double theta1;
/// 1.0 / timestep
double tau1;
/// time for right hand side functions.
double rhsTime;
/// time for boundary functions.
double boundaryTime;
/// Pointer to boundary function. Needed for initial problem.
AbstractFunction<double, WorldVector<double> > *boundaryFct;
};
......@@ -237,22 +217,14 @@ int main(int argc, char** argv)
vecheat,
adaptInfoInitial);
double fac = *(vecheat.getThetaPtr()) != 0.0 ? 1.0 : 1.0e-3;
int nRefine = 0;
int dim = vecheatSpace.getMesh(0)->getDim();
GET_PARAMETER(0, vecheatSpace.getMesh(0)->getName()
+ "->global refinements", "%d", &nRefine);
if (*(vecheat.getThetaPtr()) == 0.5)
*(adaptInfo.getTimestepPtr()) *= fac * pow(2.0, static_cast<double>(-nRefine) / dim);
else
*(adaptInfo.getTimestepPtr()) *= fac * pow(2.0, -nRefine);
// ===== create rhs functions =====
F *rhsFct0 = new F(vecheatSpace.getFeSpace(0)->getBasisFcts()->getDegree());
rhsFct0->setTimePtr(vecheat.getRHSTimePtr());
rhsFct0->setTimePtr(vecheat.getRhsTimePtr());
F *rhsFct1 = new F(vecheatSpace.getFeSpace(1)->getBasisFcts()->getDegree());
rhsFct1->setTimePtr(vecheat.getRHSTimePtr());
rhsFct1->setTimePtr(vecheat.getRhsTimePtr());
// ===== create operators =====
double one = 1.0;
......@@ -282,17 +254,17 @@ int main(int argc, char** argv)
C00.addZeroOrderTerm(new Simple_ZOT);
C00.setUhOld(vecheat.getOldSolution()->getDOFVector(0));
vecheatSpace.addMatrixOperator(C00, 0, 0,
vecheat.getTau1Ptr(), vecheat.getTau1Ptr());
vecheatSpace.addVectorOperator(C00, 0, vecheat.getTau1Ptr(), vecheat.getTau1Ptr());
vecheat.getInvTau(), vecheat.getInvTau());
vecheatSpace.addVectorOperator(C00, 0, vecheat.getInvTau(), vecheat.getInvTau());
Operator C11(vecheatSpace.getFeSpace(1), vecheatSpace.getFeSpace(1));
C11.addZeroOrderTerm(new Simple_ZOT);
C11.setUhOld(vecheat.getOldSolution()->getDOFVector(1));
vecheatSpace.addMatrixOperator(C11, 1, 1,
vecheat.getTau1Ptr(), vecheat.getTau1Ptr());
vecheat.getInvTau(), vecheat.getInvTau());