[mco]

Hi everyone,


I've consulted the programming guide for OpenFoam and I'm working on translating an equation for a decay constant to put into a solver. This is a solver for cyclic heat transfer with constant wall temperature.

Suppose I have a field that is a normalized temperature field.

The equation to find the decay constant, , is:





here, is the velocity in the direction of flow, is the normal to the wall, is the thermal diffusivity of the fluid, and is the length of the domain.


The numerator represents the heat flux leaving the control volume over the wall surface, and the denominator represents the net heat flux at the inlet boundary condition.

My solution to this is the following code in a custom solver:

// Calculate decay constant lambda

volScalarField alphaEff("alphaEff",turbulence->nu()/Pr + turbulence->nut()/Prt);

vector flowDir = (Ubar.value()/mag(Ubar.value()));
label inletID = mesh.boundaryMesh().findPatchID("inlet");
label wallsID = mesh.boundaryMesh().findPatchID("walls");

if(runTime.value() >= 20)
{
    // Calculate decay constant lambda Numerator
    dimensionedScalar lambdaNum
    (
        "lambdaNum",
        dimensionSet(0, 3, -1, 1, 0, 0, 0),
        gSum(alphaEff.boundaryField()[wallsID]*theta.boundaryField()[wallsID].snGrad()*mesh.boundary()[wallsID].magSf())
    );

     // Calculate decay constant lambda Denomenator
    dimensionedScalar lambdaDen
    (
        "lambdaDen",
        dimensionSet(0, 3, -1, 1, 0, 0, 0),
        gSum((theta.boundaryField()[inletID]*(flowDir & U.boundaryField()[inletID])*mesh.boundary()[inletID].magSf()) + 
        (alphaEff.boundaryField()[inletID]*(flowDir & fvc::grad(theta)().boundaryField()[inletID])*mesh.boundary()[inletID].magSf()))
    );
        
    lambda = -(1/Lc)*Foam::log(1 - (lambdaNum/lambdaDen));
}

Info << "Lambda: " << lambda.value() << endl;

All the constants are defined in the createFields.H file. The solver waits until time steps to let the cyclic pressure field stabilize before calculating

I suspect that it might be something to do with the second term in lambdaDen. Specifically how I take the gradient of the field value at the inlet.


Any help with this would be greatly appreciated.


Best Regards,

Michael