pPorousFluidEqn.H 1.11 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
porousRho = porousThermo.rho();
porousRho = max(porousRho, rhoMin);
porousRho = min(porousRho, rhoMax);
porousRho.relax();

volScalarField rAUPorous(1.0/porousUEqn().A());

porousU = rAUPorous*porousUEqn().H();
porousUEqn.clear();

bool closedVolume = false;

porousPhi =
    fvc::interpolate(porousRho)
    *(fvc::interpolate(porousU) & porousMesh.Sf());

for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{

    fvScalarMatrix pEqn
    (
        fvm::laplacian(porousRho*rAUPorous, porousP) ==  fvc::div(porousPhi)
    );

    pEqn.setReference(pRefCell, pRefValue);

    pEqn.solve();

    if (nonOrth == nNonOrthCorr)
    {
        porousPhi -= pEqn.flux();

    }
}

porousP.relax();

porousU -= rAUPorous*fvc::grad(porousP);
porousU.correctBoundaryConditions();

if (closedVolume)
{
    porousP += (initialMass - fvc::domainIntegrate(porousPsi*porousP))
        /fvc::domainIntegrate(porousPsi);
}

porousRho = porousThermo.rho();
porousRho = max(porousRho, rhoMin);
porousRho = min(porousRho, rhoMax);
porousRho.relax();

Info<< "rho max/min : "
    << max(porousRho).value() << " "
    << min(porousRho).value() << endl;