View Issue Details

IDProjectCategoryView StatusLast Update
0002985OpenFOAM[All Projects] Bugpublic2018-08-08 15:44
ReporterTZirwesAssigned Tohenry 
PrioritynormalSeveritycrashReproducibilityalways
Status closedResolutionsuspended 
PlatformGNU/LinuxOSOpenSuSEOS Version13.1
Product Version5.x 
Fixed in Version 
Summary0002985: Problems with DPMFoam/MPPICFoam
DescriptionSetup:
The setup is a simple cylindrical pipe (created with blockMesh) with an inlet (U = 0.01 m/s) and an outlet. No particles are inserted during the simulation.

Running the case with pimpleFoam works fine. Using MPPICFoam or DPMFoam however shows four problems:

1) MPPICFoam+DPMFoam: Maximum velocity escalates if gravity has a component in x or y direction

If the case is started with pimpleFoam, the maximum velocity stays at about 0.01 m/s.
If the case is started with DPMFoam or MPPICFoam, the maximum velocity increases to very high values. After ten time steps, the velocity is already at about 0.1 m/s (while still increasing). See the output of the function object during the run, for example after about 10 time steps:

volFieldValue minmaxvalue write:
    max() of U.gas = (0.0116891 0.0048455 0.105317)

This velocity is way too high since the inlet and initial velocity is |U|=0.01 m/s.

The velocity seems to increase (and decrease) first the outlet near the wall, and the pressure adjusts correspondingly at the inlet (see the pictures in the attached archive).

We tried to change:
* boundary conditions at inlet/outlet according to all available MPPIC and DPM tutorials (see 0/U.gas and 0/p) => no change
* the PIMPLE parameters (outerCorrections, innerCorrections...) => no change
* the discretization => no change

The problem can be resolved by changing the gravitational accelaration so that it only has a component in the z direction. Then, the maximum velocity stays constant at 0.01 m/s.
The problem can also be resolved by removing the term " + rAUcf*(g & mesh.Sf())" from UcEqn.H, thereby removing the gravitation from the continuous phase:
https://github.com/OpenFOAM/OpenFOAM-5.x/blob/master/applications/solvers/lagrangian/DPMFoam/UcEqn.H#L17

How to reproduce:

Start the case with either MPPICFoam or DPMFoam in serial. Obeserve how the maximum velocity increases with every time step for the first few time steps.
Then go to constant/g and change the value of g to value (0 0 -9.81) (the exact value does not matter as long as the component in x and y direction is approximately zero). Then the maximum velocity stays at about 0.01 m/s.


2) MPPICFoam + DPMFoam: Case cannot be started if particles are set to inactive together with manualInection:

In constant/kinematicCloudProperties, if "active" is set to false , the solver will not start because it reports a missing "massFlowRate" entry in the manualInjection model.

How to reproduce:

With either DPMFoam or MPPICFoam, start the case with "active true;" (working) or "active false;" (not working) in constant/kinematicCloudProperties.

3) MPPICFoam: GAMG solver for kinematicCloud:alpha fails without particles

If no particles are present (as in this case), using MPPICFoam slows down if the solver for kinematicCloud:alpha in system/fvSolution is changed from PBiCGStab to GAMG.
With PBiCGStab, everything works fine with and without particles.
With GAMG, everything works fine with at least one particle.
But without particles, GAMG becomes very slow taking 1000 iterations each time step for the solution of kinematicCloud:alpha:

GAMG: Solving for kinematicCloud:alpha, Initial residual = 4.58125e-05, Final residual = 1.11705e-06, No Iterations 1000

How to reproduce:
Using MPPICFoam, change the solver in line 72 of system/fvSolution to GAMG and see the difference between no particles (using constant/kinematicCloudPositions) and one particle (using constant/kinematicCloudPositions_oneParcel).


4) The case will not start, if coordinates in manualInjection are not valid anymore (probably for all lagrangian solver)

This problem is relevant for rotating/moving meshes. When the simulation is started using manualInjection, particles are added into the domain at SOI.
When the simulation is restarted at a later time, no particles are added because SOI is exeeded. But if the mesh is moving and positions in the constant/kinematicCloudPositions file do not lie in the current domain anymore, the solver will not start, forcing the user to remove the particle positions file.
The error looks like this:
--> FOAM FATAL ERROR:
Cannot find parcel injection cell. Parcel position = (100 0 0.5)

How to reproduce:

Start the case with either MPPICFoam or DPMFoam using constant/kinematicCloudPositions_oneParcles.
After some time, write out the results and end the simulation. Change the position in kinematicCloudPositions_oneParcle to (100 0 0.5) (or any other location that is not in the domain).
Restart the simulation. Although no particles should be added (because SOI is exeeded), the simulation will not start due to invalid particle locations (this corresponds to a change in mesh topology over time where the initially valid locations in the positions file are not valid anymore).
Additional InformationMost recent version of 5.x (5.x-7f7d351b741b from Jun 6, 2018)
All runs performed in serial
Before running the case, execute blockMesh
TagsLagrangian, MPPICFoam

Activities

TZirwes

2018-06-18 22:53

reporter  

BugReport.tar.gz (160,526 bytes)

will

2018-06-25 10:27

manager   ~0009821

1) The outlet boundary condition you have applied is not physical. You have a hydrostatic pressure variation in the bulk but a uniform pressure on the boundary. That difference is driving the flow. It doesn't happen in pimpleFoam because pimpleFoam doesn't add gravitational acceleration.

2) I agree that this is illogical. What's happening is that the "active" switch is disabling all reading, which defaults transient to false (despite the entry saying "true" in the same file). This changes the behaviour of the injection model, requiring a flow rate to be read, rather than a total mass.

I think that if the injection model is read, then the transient flag should also be read. If not, not. I would tend towards enabling reading in all cases and limiting the effect of the "active" flag. This could be done by removing the "if (active_)" clause around "read" in cloudSolution.C. I suspect other cases will break as a result of this change, though, so care will be required to correctly update all the examples.

3) This is a limitation of the GAMG solver

4) Agreed. This is a problem. More significantly, there doesn't seem to be any update of the injection models triggered as a result of mesh motion. There is a hook which updates as a result of topological change, but this probably needs to be applied in motion cases, too, as most models store cell indices and similar which will change as the mesh moves around the injection position.

TZirwes

2018-07-02 08:46

reporter   ~0009842

Dear Will,
thank you for your answer. Regarding 1):

I now tried all available boundary conditions for the pressure outlet in addition to the ones tested in the original post:

totalPressure and uniformTotalPressure:
maximum velocity increases just like in the uploaded case.
 
uniformDensityHydrostaticPressure and phaseHydrostaticPressure:
the maximum velocity still increases indefinitely, although slower than with the other BCs (just a remark: the keyword "rho" from these BCs requires a scalar, not a word as indicated by the documentation in https://github.com/OpenFOAM/OpenFOAM-5.x/blob/019c14b8cf0a051e8e2cc2ec5d8b2ab566fd3cd8/src/finiteVolume/fields/fvPatchFields/derived/uniformDensityHydrostaticPressure/uniformDensityHydrostaticPressureFvPatchScalarField.H and https://github.com/OpenFOAM/OpenFOAM-5.x/blob/019c14b8cf0a051e8e2cc2ec5d8b2ab566fd3cd8/src/finiteVolume/fields/fvPatchFields/derived/phaseHydrostaticPressure/phaseHydrostaticPressureFvPatchScalarField.H)

prghPressure and prghTotalPressure:
Cannot be used because an "hRef" object is not found even if the entry is provided in the outlet BC in 0/p.

prghTotalHydrostaticPressure:
ph_rgh field cannot be found.

Does this mean DPMFoam and MPPICFoam solvers can only be used if gravity is normal to all outlet patches?

henry

2018-08-08 15:43

manager   ~0009924

It is possible that this issue can be resolved by the specification of suitable boundary conditions which would take a bit of time to test and would need user support funding. If this proves unreliable the solvers would need to be changed to solve for p_rgh to handle strong buoyancy more accurately and stably and this would require development or maintenance funding:

https://openfoam.org/news/funding-2018/
https://openfoam.org/maintenance/

henry

2018-08-08 15:44

manager   ~0009925

Requires funding.

Issue History

Date Modified Username Field Change
2018-06-18 22:53 TZirwes New Issue
2018-06-18 22:53 TZirwes File Added: BugReport.tar.gz
2018-06-18 22:53 TZirwes Tag Attached: Lagrangian
2018-06-18 22:53 TZirwes Tag Attached: MPPICFoam
2018-06-25 10:27 will Note Added: 0009821
2018-07-02 08:46 TZirwes Note Added: 0009842
2018-08-08 15:43 henry Note Added: 0009924
2018-08-08 15:44 henry Assigned To => henry
2018-08-08 15:44 henry Status new => closed
2018-08-08 15:44 henry Resolution open => suspended
2018-08-08 15:44 henry Note Added: 0009925