Dear members of the DuMux community,
I am contacting you because I would like to get some pointers on how to
implement a source term [1] valid in the whole domain, [2] which would depend
on the solution vector.
>From what I understood, I need to use the "PointSource" function rather than
>the source() function.
After looking at "dumux\test\porousmediumflow\richardsnc\implicit", I added the
following:
1) in the main.cc:
before "timeLoop->start();":
problem->computePointSourceMap()
2) in the problem.hh:
using PointSource = GetPropType<TypeTag, Properties::PointSource>;
using Problem = GetPropType<TypeTag, Properties::Problem>;
[...]
void addPointSources(std::vector<PointSource>& pointSources) const
{
// loop over all vertices of the 1D domain to get one source per
subcontrol volume.
for (const auto& vertex : Dune::vertices( this->fvGridGeometry().gridView() ) )
{
auto globalPos = vertex.geometry().center();//1D space
pointSources.emplace_back(globalPos,
[this](const Problem &problem,
const Element &element,
const FVElementGeometry &fvGeometry,
const ElementVolumeVariables &elemVolVars,
const SubControlVolume &scv)
{
return bioChemicalReaction(elemVolVars,scv);
}
);
}
}
PrimaryVariables biochemicalReactions( const ElementVolumeVariables
&elemVolVars,
const SubControlVolume &scv)
{
const auto& volVars = elemVolVars[scv];
double volume = scv.volume(); //to have m3 scv.
const auto massOrMoleDensity = [](const auto& volVars, const int compIdx, const
bool isFluid)
{
return (useMoles ? volVars.molarDensity(compIdx) : volVars.density(compIdx) );
};
const auto massOrMoleFraction= [](const auto& volVars, const int phaseIdx,
const int compIdx)
{
return (useMoles ? volVars.moleFraction(phaseIdx, compIdx) :
volVars.massFraction(phaseIdx, compIdx) )
};
//mol comp/m3 scv
double C_S = massOrMoleDensity(volVars, h2OIdx) //mol phase/ m3 phase
* massOrMoleFraction(volVars,h2OIdx, mucilIdx) //mol comp/ mol phase
* volVars.saturation(h2OIdx) //m3 phase / m3 pores
* volVars.porosity(); // m3 pores/m3 scv
double rate = - 0.5; //1/s
double dwater_dt = 0.;
// reaction rate in
mol/s for all components
// of the subcontrol volume
double dC_S_dt = rate * C_S * volume; // example of a simple reaction
return PrimaryVariables({dwater_dt,dC_S_dt});
}
3) in properties.hh:
template<class TypeTag>
struct PointSource<TypeTag, TTag::RichardsNCTT> { using type =
SolDependentPointSource<TypeTag>; };
####
For this simplified toy model, [1] there are no flux BCs, [2] the (molar)
density of the phase is always equal to that of pure water,
[3] there is a simple equation for the reaction to make the comparison with a
reference analytical solution easier,
(the biochemical reactions implemented later will be more complex).
Therefore the only change in the concentration of the component comes from the
reaction = [d C_S/dt = -0.5 * C_S].
I have the issue that the solver yields at each time step:
[ C_S_end = C_S_init/(1 - rate ) ]
a) Do I need to overload the "addSourceDerivatives()" function and/or change
the solver settings to get the correct result?
b) If yes, do you have one or several examples using the
"addSourceDerivatives()"?
c) Are there other functions which need to be overloaded/changed (like
"pointSourceDerivative()")?
d) If I understood correctly, using "addPointSources()" or "pointSource()"
yield the same results? cf. comments in the function
"dumux\common\fvproblem.hh:scvPointSources()"
NB: I use version 3.0 of DuMux
Many thanks for your help!
With best regards,
Mona Giraud
------------------------------------------------------------------------------------------------
------------------------------------------------------------------------------------------------
Forschungszentrum Juelich GmbH
52425 Juelich
Sitz der Gesellschaft: Juelich
Eingetragen im Handelsregister des Amtsgerichts Dueren Nr. HR B 3498
Vorsitzender des Aufsichtsrats: MinDir Stefan Müller
Geschaeftsfuehrung: Prof. Dr. Astrid Lambrecht (Vorsitzende),
Karsten Beneke (stellv. Vorsitzender), Dr. Ir. Pieter Jansens, Prof. Dr. Frauke
Melchior
------------------------------------------------------------------------------------------------
------------------------------------------------------------------------------------------------
_______________________________________________
DuMux mailing list
[email protected]
https://listserv.uni-stuttgart.de/mailman/listinfo/dumux
