Commit 35b5f9e6 authored by Bernd Flemisch's avatar Bernd Flemisch
Browse files

[el2p, multidomain] fix residual calculation, generalize multidomain

convergence writer

- The residual calculation for el2p was wrong, probably never used. It now
  forwards the call to the GridOperator functionality.

- The multidomain convergence writer assumed that sub-domain solution
  vectors can be default-contructed. This is not the case if a PDELab
  model is used for the subdomain.

Reviewed by Dennis.



git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/trunk@14200 2fb0f335-1f38-0410-981e-8018bf24f1b0
parent 3bb14e46
......@@ -629,6 +629,8 @@ private:
{ return problem_().vertexMapper(); }
const ElementMapper &elementMapper_() const
{ return problem_().elementMapper(); }
const GridOperator &gridOperator() const
{ return *gridOperator_;}
Problem *problemPtr_;
......
......@@ -225,46 +225,7 @@ public:
Scalar globalResidual(SolutionVector &residual,
const SolutionVector &u)
{
SolutionVector tmp(curSol());
curSol() = u;
Scalar res = globalResidual(residual);
curSol() = tmp;
return res;
}
/*!
* \brief Compute the global residual for the current solution
* vector.
*
* \param residual Stores the result
*/
Scalar globalResidual(SolutionVector &residual)
{
residual = 0;
ElementIterator eIt = gridView_().template begin<0>();
const ElementIterator eEndIt = gridView_().template end<0>();
for (; eIt != eEndIt; ++eIt) {
localResidual().eval(*eIt);
if (isBox)
{
#if DUNE_VERSION_NEWER(DUNE_COMMON, 2, 4)
for (int i = 0; i < eIt->subEntities(dim); ++i)
#else
for (int i = 0; i < eIt->template count<dim>(); ++i)
#endif
{
int globalI = vertexMapper().map(*eIt, i, dim);
residual[globalI] += localResidual().residual(i);
}
}
else
{
int globalI = elementMapper().map(*eIt);
residual[globalI] = localResidual().residual(0);
}
}
jacAsm_.gridOperator().residual(u, residual);
// calculate the square norm of the residual
Scalar result2 = residual.two_norm2();
......@@ -283,6 +244,17 @@ public:
return std::sqrt(result2);
}
/*!
* \brief Compute the global residual for the current solution
* vector.
*
* \param residual Stores the result
*/
Scalar globalResidual(SolutionVector &residual)
{
return globalResidual(residual, curSol());
}
/*!
* \brief Compute the integral over the domain of the storage
* terms of all conservation quantities.
......@@ -711,13 +683,13 @@ public:
def[eqIdx] = writer.allocateManagedBuffer(numDofs);
}
for (unsigned int vIdxGlobal = 0; vIdxGlobal < u.size(); vIdxGlobal++)
for (unsigned int vIdxGlobal = 0; vIdxGlobal < u.base().size(); vIdxGlobal++)
{
for (int eqIdx = 0; eqIdx < numEq; ++eqIdx)
for (int eqIdx = 0; eqIdx < numEq; ++eqIdx)
{
(*x[eqIdx])[vIdxGlobal] = u[vIdxGlobal][eqIdx];
(*delta[eqIdx])[vIdxGlobal] = - deltaU[vIdxGlobal][eqIdx];
(*def[eqIdx])[vIdxGlobal] = residual[vIdxGlobal][eqIdx];
(*x[eqIdx])[vIdxGlobal] = u.base()[vIdxGlobal][eqIdx];
(*delta[eqIdx])[vIdxGlobal] = - deltaU.base()[vIdxGlobal][eqIdx];
(*def[eqIdx])[vIdxGlobal] = residual.base()[vIdxGlobal][eqIdx];
}
}
......
......@@ -116,21 +116,13 @@ struct MultiDomainConvergenceWriter
void writeFields(const SolutionVector &uLastIter,
const SolutionVector &deltaU)
{
SolutionVector1 uLastIter1;
SolutionVector2 uLastIter2;
SolutionVector1 deltaU1;
SolutionVector2 deltaU2;
uLastIter1.resize(ctl_.method().model().sdModel1().numDofs());
uLastIter2.resize(ctl_.method().model().sdModel2().numDofs());
deltaU1.resize(ctl_.method().model().sdModel1().numDofs());
deltaU2.resize(ctl_.method().model().sdModel2().numDofs());
typedef Dumux::SplitAndMerge<TypeTag> Common;
Common::splitSolVector(uLastIter, uLastIter1, uLastIter2);
Common::splitSolVector(deltaU, deltaU1, deltaU2);
SolutionVector1 uLastIter1(ctl_.method().model().sdModel1().curSol());
SolutionVector2 uLastIter2(ctl_.method().model().sdModel2().curSol());
SolutionVector1 deltaU1(uLastIter1);
SolutionVector2 deltaU2(uLastIter2);
SplitAndMerge<TypeTag>::splitSolVector(uLastIter, uLastIter1, uLastIter2);
SplitAndMerge<TypeTag>::splitSolVector(deltaU, deltaU1, deltaU2);
std::cout << "\n writing convergence file of current Newton iteration \n";
ctl_.method().model().sdModel1().addConvergenceVtkFields(*vtkMultiWriter1_, uLastIter1, deltaU1);
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment