UEqnAddPorosity.H
Go to the documentation of this file.
1// Including porosity effects in UEqn following:
2// Jensen, B., Jacobsen, N. G., & Christensen, E. D. (2014).
3// Investigations on the porous media equations and resistance
4// coefficients for coastal structures. Coastal Engineering, 84, 56-72.
5
7{
8 const volScalarField& porosity = tporosity.cref();
9
10 const word porosityModel("JensenEtAl2014");
11 const dictionary& dict =
12 porosityProperties.subDict(porosityModel + "Coeffs");
13 const dimensionedScalar alpha(dimless/dimArea, dict.get<scalar>("alpha"));
14 const dimensionedScalar beta(dimless/dimLength, dict.get<scalar>("beta"));
15 const dimensionedScalar d50(dimless, dict.get<scalar>("d50"));
16 const dimensionedScalar KC(dimless, dict.get<scalar>("KC"));
17
18 // Generating Darcy-Forchheimer coefficient: F = rho*U*(a + b*|U|)
19 // Shoud it be mu or muEff in the equation below?
20 {
21 // Darcy term
22 volScalarField DarcyForchheimerCoeff
23 (
24 alpha*sqr(1 - porosity)*mixture.mu()/sqr(porosity)/sqr(d50)
25 );
26
27 // Adding Forchheimer term
28 DarcyForchheimerCoeff += rho*mag(U)
29 *beta*(1 + pos(KC)*7.5/KC)*(1 - porosity)/sqr(porosity)/d50;
30
31 // Adding Darcy-Forchheimer term as implicit source term
32 UEqn += fvm::Sp(DarcyForchheimerCoeff, U);
33 }
34
35 {
36 // Generating added mass force coefficient
37 const dimensionedScalar gamma_p(dimless, dict.get<scalar>("gamma_p"));
38 const volScalarField Cm(gamma_p*(1 - porosity));
39
40 UEqn += Cm*fvm::ddt(rho, U);
41 UEqn += Cm*MRF.DDt(rho, U);
42 }
43
44 // Dividing both matrix entries and source term by porosity to compensate
45 // for the fact that the FVM cell volume averages use division by cell
46 // volume V whereas only the cell pore volume, porosity*V, is accessible.
47 UEqn *= scalar(1)/porosity;
48}
IOMRFZoneList & MRF
U
Definition: pEqn.H:72
fvVectorMatrix & UEqn
Definition: UEqn.H:13
IOdictionary porosityProperties(IOobject("porosityProperties", runTime.constant(), runTime, IOobject::READ_IF_PRESENT, IOobject::NO_WRITE))
const bool porosityEnabled(porosityProperties.getOrDefault< bool >("porosityEnabled", false))
tmp< volScalarField > tporosity
volScalarField & alpha
dictionary dict
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)
Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n"<< endl;autoPtr< temperaturePhaseChangeTwoPhaseMixture > mixture
Definition: createFields.H:39