Change coupling execution order (CFD first, then DEM)

Nucleophobe's picture
Submitted by Nucleophobe on Fri, 03/04/2016 - 06:01

Hi,

From what I can tell, in all of the CFDEMcoupling solvers, the execution starts with a call to the DEM solver, after which the flow field is solved using CFD.

For instance, in cfdemSolverIB, particleCloud.evolve() is called at the beginning of the time loop, which then runs the function dataExchangeM().couple(). This function then executes LIGGGHTS, which, during the first timestep, initializes the number of particles in the simulation, etc., and ALSO advances the DEM solver the specified number of DEM timesteps (e.g., 10).

I am working on a custom solver similar to 'cfdemSolverIB' and I would like to reverse the order of the solvers (CFD first, THEN DEM). The algorithm would look something like this:
1) Call LIGGGHTS to initialize the particleCloud variables (number of particles, particle positions, and particle velocities) BUT don't execute LIGGGHTs for any timesteps.
2) Calculate the forces on the particles using the initial conditions of the flow field.
3) Pass the particle forces to LIGGGHTS and run for the appropriate number of DEM timesteps (e.g., 10).
4) Advance the CFD solver one time step.
5) Perform the pressure and velocity correction steps to couple the CFD and DEM.
6) Advance in time.

Is there any built-in way to accomplish this (in particular, step (1)), or do I need to alter the data exchange code(s) significantly?

I can post more details if needed; I'm just trying to get a general idea of how I should proceed.

Thanks,
-Nuc

Forums:

Nucleophobe's picture

Nucleophobe | Tue, 03/08/2016 - 18:41

Hi,

It's been a few days, so I just thought that I would report back that I have figured out a solution for now.

LAMMPS/LIGGGHTS can accept an input of run 0 to initialize simulations without taking any timesteps:
http://lammps.sandia.gov/doc/run.html

So, to accomplish the goals I outlined above, I modified the data exchange model twoWayMPI so that on the first coupling call, run 0 is performed. Then, on subsequent calls, the default behavior is restored (i.e., run X number of timesteps, according to the input files).

I also modified cfdemCloudIB.C so that dataExchangeM().couple() is called at the beginning and end of the coupling operation on the first call, but only at the end of the coupling operation for subsequent calls. This isn't the most elegant solution, but it works.

Hopefully this helps someone else.
-Nuc

j-kerbl's picture

j-kerbl | Thu, 03/17/2016 - 16:48

Hi Nuc,

thanks for sharing!

Cheers
Josef

NTT1508 | Wed, 03/30/2016 - 02:25

Hi Nuc,

Thank for your contribution. Do you think this will result in more accuracy? I think it does not much effect on result because the initial step is normally very small compared to the whole computation time in many cases.

Regards,

Nucleophobe's picture

Nucleophobe | Thu, 03/31/2016 - 22:29

NTT1508:

As originally implemented, there are no fluid drag forces during the first time step. Furthermore, the buoyancy force is also not considered (!).

So, in some cases (e.g., when considering gravity), a non-physically large acceleration can take place during the first N number of DEM timesteps before CFD coupling occurs.

I noticed this problem when examining the forces on a falling sphere. During the first set of DEM timesteps, the full weight of the particle is the only force applied, so the particle may accelerate very quickly. Then, during the CFD coupling step, the fluid drag forces that are calculated are non-physically large. The oscillation in the forces is eventually damped out, but the accuracy of the velocity vs. time curve is definitely effected.

You may want to examine the forces vs. time and the velocity vs. time to determine whether you need to worry about this problem.
-Nuc

Nucleophobe's picture

Nucleophobe | Thu, 03/31/2016 - 22:57

*...but the accuracy of the velocity vs. time curve is definitely affected.

The root of the problem here is that the CFDEM coupling occurs in a partitioned or "staggered" manner. That is, the solution of the fluid motion and the solution of the solid motion take place separately, which means there is always some error in the forces and in the particle motion.

There is some discussion on loose vs. tight coupling of on the Fluid--Structure Interaction Wikipedia page:
https://en.wikipedia.org/wiki/Fluid%E2%80%93structure_interaction

Anyway, the solution is to implement sub iterations to converge the forces and displacements to the correct solution. You can find more details here in this thesis:
https://etda.libraries.psu.edu/paper/11068/5747

"Fixed-point fluid-structure interaction solvers with dynamic relaxation" by Küttler and Wall (Computational Mechanics, 2008) is another good resource.
-Nuc

NTT1508 | Mon, 04/04/2016 - 01:17

Hi Nuc,

Yes the gap could be significant when the number of steps for exchange between DEM-CFD is large. In that case, a dummy run (0) should work well.
Thank for your clarification.

Regards,