[mikeP]

Hello,
I have a momentum transport equation which is not pressure-driven, in addition to the continuity equation. Instead it has gravity and drag terms as the driving force.

I can include it into OpenFOAM as:
div(phi,U) == gravity + drag

However, I can not include the continuity equation with this one. The solver does not compile when I try
fvm::div(phi) or fvc::div(phi)

I know that in simple algorithm, the continuity equation is included in the pressure equation. But i don't have to solve for pressure, it is not anywhere in the equations. How can I do it?

[mikeP]

I need to include the equation shown in attachment, which is a mass conservation equation, where "a" is phase volume fraction in a multiphase flow.

I defined "a" as a volScalarField.

When I do at as:

fvScalarMatrix alphaEqn
            (
                fvm::ddt(alpha)
                + fvm::div(alpha,U)
            );
            solve(alphaEqn);

I get:
error: no matching function for call to ‘div(Foam::volScalarField&, Foam::volVectorField&)’

I suppose, I can not use two volume fields as argument in fvm::div

Then I switch to

fvScalarMatrix alphaEqn
            (
                fvm::ddt(alpha)
                + fvc::div(phi)
            );
            solve(alphaEqn);

where phi is defined as linearInterpolate(alpha*U) & mesh.Sf()
This time solver converges, but when I solve my case, it gives floating point exception if my internal field value for "a" is different than the voundary condition value.

[nimasam]

try this :

fvScalarMatrix alphaEqn
(
fvm::ddt(alpha)
+ fvc::div(phi,alpha)
);
solve(alphaEqn);

phi here is U & mesh.Sf()

[mikeP]

I modified it as
linearInterpolate(alpha*U) & mesh.Sf()
in order to get phase flux. I am not sure if this makes sense.