[Diro7]

Hello!

I'm working on a custom solver derived from buoyantBoussinesqSimpleFoam (OF6).
I need a simple sink for energy (linear in temperature), so until now what I did was something like:

fvScalarMatrix TEqn
(
    fvm::div(phi, T)
  - fvm::laplacian(alphaEff, T)
==
    radiation->ST(rhoCpRef, T)
  - heatExchanger*(HXCoeff/rhoCp)*(fvm::Sp(1, T) - TExt)
  + fvOptions(T)
);

where heatExchanger is a flag-like nondimensional field which is 1 in the cells which are part of the heat exchanging region and 0 elsewhere (in the case I'm using as test, it is actually 1 everywhere). HXCoeff is a sort of volumetric heat transfer coefficient (W m-3 K-1), TExt is the heat sink temperature and rhoCp is read from a custom dictionary (since the original rhoCpRef is needed by the dummy radiation model and is initialised to 1.0 in the original solver). All these parameters in bold are of dimensionedScalar type.

Everything runs ok so far and the results are exactly as expected.

The problems come when I try to get rid of this hard-coded solution and use the semiImplicitSource fvOptions source, which in principle should allow to model the same sink at run-time.
I removed the sink term from T equation, and added a fvOptions file which contains:

heatExchanger
{
    type    scalarSemiImplicitSource;

    selectionMode   all; 
    volumeMode      specific;

    injectionRateSuSp
    {
        T   (Su Sp);
    }
}

where

Su = TExt * HXCoeff / rhoCp
Sp = - HXCoeff / rhoCp

The solver diverges showing a totally different behaviour with respect the other implementation.
I tried everything, such as switching from specific to absolute mode, changing signs and other stuff, but it seems a can't replicate the behaviour I desire.

Can anyone tell me what I'm doing wrong? I have no clue, except (maybe) for having noticed that TEqn has dimensions [0 3 -1 1 0 0 0] and I don't understand how it is related to the way I have to implement a scalarSemiImplicitSource for the T equation.

Thank you in advance for the help!

Andrea


EDIT:

I figured out that the dimensions of TEqn are ok. They are basically K/s (the dimensions of the "physical" equation) multiplied by volume. This feature seems common to all kinds of equations and I guess it being related to the finite volume discretisation (they are the dimensions of the discretised/volume-integrated equation, not the physical one). This also explains why the addSup member function from SemiImplicitSource.C reads

template<class Type>
void Foam::fv::SemiImplicitSource<Type>::addSup
(
    fvMatrix<Type>& eqn,
    const label fieldi
)
{
    [...]

    const GeometricField<Type, fvPatchField, volMesh>& psi = eqn.psi();

    typename GeometricField<Type, fvPatchField, volMesh>::Internal Su
    (
        IOobject
        (
            name_ + fieldNames_[fieldi] + "Su",
            mesh_.time().timeName(),
            mesh_,
            IOobject::NO_READ,
            IOobject::NO_WRITE
        ),
        mesh_,
        dimensioned<Type>
        (
            "zero",
            eqn.dimensions()/dimVolume,
            Zero
        ),
        false
    );

    [...]
 }

where equation dimensions are divided by volume to retrieve the correct physical dimensions before adding the SuSp term to the original equation.

Anyway, this explains the "strange" dimensions but doesn't solve my problem.
Indeed, this also confirms that the fvOptions (Su Sp) coefficients for the T equation are supposed to have dimensions K/s and 1/s respectively, which means the I'm not messing up with the coefficients.

The specific/absolute switch, if I got it right, basically determines if the coefficients are good as they are provided or if they need to be normalised by the cellZone volume before being applied to the equation.
So I really don't get what I'm doing wrong, since I only have to set two simple coefficients...

[Eubo]

Hi!
Did you figure out why it was diverging?
I am running a similar case using semiImplicitSource and buoyantBoussinesqSimpleFoam. I want to confirm if the method of calculating Su Sp and explanations that you have provided are correct.

[Diro7]

I don't exactly remember the issue, but my calculations were right.
The issue, if I'm correct, was related with how semiImplicitSource switches between implicit and explicit treatment of the linear term of the source based on the sign of the coefficient.
If you are experiencing similar issues I can give you more detail (I'm pretty sure there was a thread specifically on this matter but I can't find it right now).