Hi Dennis,
Thanks a lot, I got through elemFluxVarsCache.bind runtime error!
However, a new problem has popped up đ
When I get to computing (or recovering) the advective flux with
fluxVars.advectiveFlux(nPhaseIdx, upwindTermN), the execution brings me to
caluclating Dracy's flux in CCTpfaDarcysLaw::flux(...) since
- !advFluxIsCached_[phaseIdx] is true in
PorousMediumFluxVariables::advectiveFlux(..) and
- (!problem.couplingManager().isOnInteriorBoundary(element, scvf)) is true in
CCTpfaFacetCouplingDarcysLawImpl::flux(..)
So I am ending up in CCTpfaDarcysLaw::flux(...) and the outsideVolVars =
elemVolVars[scvf.outsideScvIdx()] cannot be computed because scvf is a boundary
face, resulting in the error message "Assertion `it != volVarIndices_.end() &&
"Could not find the current volume variables for volVarIdx!"' failed."
Looking forward for your feedback.
Thanks,
Nikolai
________________________________
From: Dumux <dumux-boun...@listserv.uni-stuttgart.de> on behalf of Dennis
Gläser <dennis.glae...@iws.uni-stuttgart.de>
Sent: Tuesday, August 21, 2018 11:39:34 AM
To: dumux@listserv.uni-stuttgart.de
Subject: Re: [DuMuX] Check global conservation error
Hi Nikolai,
I think your implementation is almost correct. If you want to balance masses,
you have to add the density to your upwind terms though, i.e.:
auto upwindTermN = [](const auto& volVars) { return
volVars.mobility(nPhaseIdx)*volVars.density(nPhaseIdx); };
...
Also, you have to call bind() instead of bindElement() on the geometry and the
volume variables, i.e.:
auto fvGeometry = localView(this->fvGridGeometry());
fvGeometry.bind(element);
auto elemVolVars = localView(gridVariables.curGridVolVars());
elemVolVars.bind(element, fvGeometry, curSol);
This has the following background:
If you want to do computations on the element without assembling fluxes, you do
not have to prepare all variables within the entire stencil. For these purposes
you call bindElement(). This only prepares primary/secondary variables on the
element. This is enough e.g. to evaluate the mass inside the domain after a
time step. But, if you want to compute fluxes, you need to prepare the entire
stencil. This stencil includes neighboring elements as well for cell-centered
schemes. In this case you prepare you local views accordingly by calling bind()
instead of bindElement().
I hope this makes it work!
Best wishes,
Dennis
On 21.08.2018 10:48, Nikolai Andrianov wrote:
Dear DuMuX experts,
Please advise how can I implement a global conservation error check, that is,
evaluate the residual ( mass(t^{n+1}) - mass(t^n) ) / dt + sum of fluxes over
the domain boundaries.
So far I can evaluate the mass in the domain by looping through the elements as
in void getVolumes of
https://git.iws.uni-stuttgart.de/andrian/rate-sens-nofrac/blob/master/src/rate_problem.hh.
I tried to use the following to access the fluxes in rate_problem.hh (inspired
from porousmediumflow/velocityoutput.hh and
porousmediumflow/nonequilibrium/gridvariables.hh):
// the upwind term to be used for the volume flux evaluation
auto upwindTermN = [](const auto& volVars) { return
volVars.mobility(nPhaseIdx); };
auto upwindTermW = [](const auto& volVars) { return
volVars.mobility(wPhaseIdx); };
auto elemFluxVarsCache =
localView(gridVariables.gridFluxVarsCache());
elemFluxVarsCache.bind(element, fvGeometry, elemVolVars);
for (auto&& scvf : scvfs(fvGeometry)) {
if (scvf.boundary()) {
//auto bcTypes = problemBoundaryTypes(element, scvf);
//if (bcTypes.hasOnlyDirichlet())
{
FluxVariables fluxVars;
fluxVars.init(*this, element, fvGeometry, elemVolVars,
scvf, elemFluxVarsCache);
Scalar extrusionFactor =
problem_.extrusionFactor(element, fvGeometry.scv(scvf.insideScvIdx()),
elementSolution(element, sol_, fvGridGeometry_));
FO += fluxVars.advectiveFlux(nPhaseIdx, upwindTermN);//
/ extrusionFactor;
FW += fluxVars.advectiveFlux(wPhaseIdx, upwindTermW);//
/ extrusionFactor;
}
}
}
However, I get compilation error at elemFluxVarsCache.bind(element, fvGeometry,
elemVolVars) saying that Assertion `it != indices.end() && "Could not find the
scv/scvf! Make sure to properly bind this class!"' failed.
Your feedback is greatly appreciated.
Many thanks,
Nikolai
PS: please disregard my previous mail with the same subject..
_______________________________________________
Dumux mailing list
Dumux@listserv.uni-stuttgart.de<mailto:Dumux@listserv.uni-stuttgart.de>
https://listserv.uni-stuttgart.de/mailman/listinfo/dumux
_______________________________________________
Dumux mailing list
Dumux@listserv.uni-stuttgart.de
https://listserv.uni-stuttgart.de/mailman/listinfo/dumux
