Skip to content
Snippets Groups Projects
Commit 210485c4 authored by Kilian Weishaupt's avatar Kilian Weishaupt Committed by Timo Koch
Browse files

[staggeredGrid][test] Adapt donea test to new vtk interface

parent f40c4cca
No related branches found
No related tags found
2 merge requests!617[WIP] Next,!370Feature/staggered grid
...@@ -34,9 +34,6 @@ ...@@ -34,9 +34,6 @@
#include <dumux/material/components/constant.hh> #include <dumux/material/components/constant.hh>
#include <dumux/linear/amgbackend.hh>
namespace Dumux namespace Dumux
{ {
template <class TypeTag> template <class TypeTag>
...@@ -291,51 +288,44 @@ public: ...@@ -291,51 +288,44 @@ public:
} }
/*! /*!
* \brief Append all quantities of interest which can be derived * \brief Adds additional VTK output data to the VTKWriter. Function is called by the output module on every write.
* from the solution of the current time step to the VTK
* writer.
*/ */
void addOutputVtkFields() template<class VtkOutputModule>
void addVtkOutputFields(VtkOutputModule& outputModule) const
{ {
//Here we calculate the analytical solution auto& pressureExact = outputModule.createScalarField("pressureExact", 0);
const auto numElements = this->gridView().size(0); auto& velocityExact = outputModule.createVectorField("velocityExact", 0);
auto& pExact = *(this->resultWriter().allocateManagedBuffer(numElements));
auto& velocityExact = *(this->resultWriter()).template allocateManagedBuffer<double, dimWorld>(numElements);
using DofTypeIndices = typename GET_PROP(TypeTag, DofTypeIndices);
typename DofTypeIndices::FaceIdx faceIdx;
const auto someElement = *(elements(this->gridView()).begin());
auto someFvGeometry = localView(this->model().globalFvGeometry()); auto& scalarFaceVelocityExact = outputModule.createFaceScalarField("scalarFaceVelocityExact");
someFvGeometry.bindElement(someElement); auto& vectorFaceVelocityExact = outputModule.createFaceVectorField("vectorFaceVelocityExact");
const auto someScv = *(scvs(someFvGeometry).begin());
Scalar time = std::max(this->timeManager().time() + this->timeManager().timeStepSize(), 1e-10);
for (const auto& element : elements(this->gridView())) for (const auto& element : elements(this->gridView()))
{ {
auto fvGeometry = localView(this->model().globalFvGeometry()); auto fvGeometry = localView(this->model().globalFvGeometry());
fvGeometry.bindElement(element); fvGeometry.bindElement(element);
for (auto&& scv : scvs(fvGeometry)) for (auto&& scv : scvs(fvGeometry))
{ {
auto dofIdxGlobal = scv.dofIndex(); auto ccDofIdx = scv.dofIndex();
const auto& globalPos = scv.dofPosition(); auto ccDofPosition = scv.dofPosition();
auto analyticalSolutionAtCc = dirichletAtPos(ccDofPosition);
pExact[dofIdxGlobal] = dirichletAtPos(globalPos)[pressureIdx];
GlobalPosition velocityVector(0.0); GlobalPosition velocityVector(0.0);
for (auto&& scvf : scvfs(fvGeometry)) for (auto&& scvf : scvfs(fvGeometry))
{ {
auto faceDofIdx = scvf.dofIndex();
auto faceDofPosition = scvf.center();
auto dirIdx = scvf.directionIndex(); auto dirIdx = scvf.directionIndex();
velocityVector[dirIdx] += 0.5*dirichletAtPos(globalPos)[faceIdx][dirIdx]; auto analyticalSolutionAtFace = dirichletAtPos(faceDofPosition);
scalarFaceVelocityExact[faceDofIdx] = analyticalSolutionAtFace[faceIdx][dirIdx];
GlobalPosition tmp(0.0);
tmp[dirIdx] = analyticalSolutionAtFace[faceIdx][dirIdx];
vectorFaceVelocityExact[faceDofIdx] = std::move(tmp);
} }
velocityExact[dofIdxGlobal] = velocityVector; pressureExact[ccDofIdx] = analyticalSolutionAtCc[pressureIdx];
velocityExact[ccDofIdx] = analyticalSolutionAtCc[faceIdx];
} }
} }
this->resultWriter().attachDofData(pExact, "p_exact", false);
this->resultWriter().attachDofData(velocityExact, "velocity_exact", false, dim);
} }
/*! /*!
......
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