[mk1]

Hello everyone,

I'm trying to simulate a system which consists of three (immiscible) fluids
using OpenFOAM 4.0 and the multiphaseEulerFoam solver.
The phases are air, "particles" and water. The air and particles are dispersed in the water.

The geometry is a rectangle (2D-problem), which contains 2 walls and an inlet and an outlet. The region near the inlet is filled solely by air, while
the region near the outlet is occupied by both particles and water, with phase fractions 0.5 and 0.5, respectively.

What I'd like to do is as follows:
A flow is induced by setting a pressure difference between the inlet and the outlet (pres. diff. ~ 50kPa). The particle phase should not flow through the outlet at all (semi-permeable wall). Hence, only the water phase exits through the outlet and only air should flow in through the inlet to conserve the volume. Intuitively, the phase fraction of the particles should reach 1 at the outlet at some point, and water can no longer exit (the outlet becomes clogged).

For some reason though, my velocity fields of all the phases near the phase interface start to blow up after around 0.05 seconds. I've tried
various drag force models to no avail, which makes me think my BCs might be incorrect.

BCs:

p_rgh:

FoamFile
{
version 2.0;
format ascii;
class volScalarField;
object p_rgh;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions [1 -1 -2 0 0 0 0];

internalField uniform 1e5;

boundaryField
{
inlet
{
type fixedValue;
value $internalField;
}

outlet
{
type prghPressure;
p uniform 5.0e4;
value uniform 5.0e4;
}

walls
{
type fixedFluxPressure;
value $internalField;
}

frontAndBackPlanes
{
type empty;
}
}

and U.air and U.water:

FoamFile
{
version 2.0;
format ascii;
class volVectorField;
object U.air;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions [0 1 -1 0 0 0 0];

internalField uniform (0 0 0);

boundaryField
{
inlet
{
type pressureInletOutletVelocity;
value $internalField;
}

outlet
{
type pressureInletOutletVelocity;
value $internalField;
}

walls
{
type noSlip;
}

frontAndBackPlanes
{
type empty;
}
}

and U.particles:

FoamFile
{
version 2.0;
format ascii;
class volVectorField;
object U.particles;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions [0 1 -1 0 0 0 0];

internalField uniform (0 0 0);

boundaryField
{
inlet
{
type zeroGradient;
}

outlet
{
type fixedValue;
value uniform (0 0 0);
}

walls
{
type fixedValue;
value uniform (0 0 0);
}

frontAndBackPlanes
{
type empty;
}
}

Are these BCs correct for what I'm trying to accomplish? Any other ideas
on why the simulation blows up? I'm attaching a minimal case for testing.

Any help is greatly appreciated!