Simulation of viscoelastic flows

Submitted by ROAC on Mon, 11/07/2011 - 00:27

Hi,

I want to implement a simulation of viscoelastic flow that consists of a polymer melt with suspensed microscale particles . Does CFDEM have the module for viscoelastic flow?

AC

Forums:

cgoniva's picture

cgoniva | Mon, 11/07/2011 - 11:18

Hi AC!

That sounds like a very interesting application.

Currently we have no such solver, but that should not be a big deal.

I'd suggest to take a viscoeslastic solver coming with the ext version and add the CFDEM functionality following the example of cfdemSolverPiso_shared.

Cheers,
Chris

ROAC | Fri, 12/02/2011 - 19:38

Thank you Chris. I am thinking ...would it be possible to simulate viscoelastic flow just using the SPH method (without coupling CFD)?

If SPH works well, should I be able to put solid particles into the flow to simulate migration of particles?...

AC

cgoniva's picture

cgoniva | Mon, 12/05/2011 - 10:47

Hi!

to be honest, I am not an SPH expert, thus I cannot give you any advise on SPH with viscoelastic flow.

Cheers,
Chris

heliana60 | Fri, 04/15/2016 - 15:46

Hello Christoph,

I want to use the resolved method using viscoelastic fluids. Do you think you could help me out solving some errors I have found so far :S

I am not really good at programming, but so far what I have done is to create a new solver based on the cfdemSolverIB and introduced a new divergence of the stress (depending on viscoelastic models), which is based on models already implemented in OpenFOAM (viscoelasticFluidFoam Laws). In any case, the flow is solved correctly, my problem comes when calculating according to the immersed boundary conditions method the force acting on the particles, which is according to Shirgaonkar (2009) and the one implemented for these IB method f=-pI+div(NewtonianStress). ShirgaonkarIB calls the IBDragPerV in the forceSubModels (volumetric force) which is get, but I want to include also the divergence of the stress that comes from the non-Newtonian part of the fluid. So I have something like:

IBViscoDragPerV_ = rhoField()*(fvc::div(tau_/rho_)+nuField()*fvc::laplacian(U)-fvc::grad(p));
return IBViscoDragPerV_;

the fvc::div(tau_) part is so far the problem for me :(

tau_ is created in the viscoelastic models I use, and it comes from the viscoelasticLaw. Actually it is the non-newtonian stress (i should rename it because it is misleading).
I am really lost at introducing tau_ in the IBViscoDragPerV_ equation, since tau is created at each time step I just want to read it and calculate fvc::div(tau_/rho_) .... When I try to run the model it crashes after a few seconds showing something like:

Step Atoms KinEng rke Volume
21 2 0.00083081303 0 4
[3] #0 Foam::error::printStack(Foam::Ostream&) at ??:?
[3] #1 Foam::sigFpe::sigHandler(int) at ??:?
[3] #2 ? in "/lib/x86_64-linux-gnu/libc.so.6"
[3] #3 LAMMPS_NS::Neighbor::check_distance() at ??:?
[3] #4 LAMMPS_NS::Verlet::run(int) at ??:?
[3] #5 LAMMPS_NS::Run::command(int, char**) at ??:?
[3] #6 void LAMMPS_NS::Input::command_creator(LAMMPS_NS::LAMMPS*, int, char**) at ??:?
[3] #7 LAMMPS_NS::Input::execute_command() at ??:?
[3] #8 LAMMPS_NS::Input::one(char const*) at ??:?
[3] #9 Foam::twoWayMPI::couple(int) const at ??:?
[3] #10 Foam::cfdemCloudIB::evolve(Foam::GeometricField&) at ??:?
[3] #11 ? at ??:?
[3] #12 __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
[3] #13 ? at ??:?
[mp-oldroyd:05316] *** Process received signal ***
[mp-oldroyd:05316] Signal: Floating point exception (8)
[mp-oldroyd:05316] Signal code: (-6)

Which is super scary. I know obviously I am doing something wrong with the code :) , so if you could give me an advice/hint how to include this in the IBViscoDragPerV_ force I would really appreciate it :) I think I am going mad :(

Have a nice weekend!

Heliana

heliana60 | Wed, 05/04/2016 - 13:41

Hello there,

First of all sorry for the last entry which was so inconclusive. I did not put the effort in the programming part as I should have had.

Now I have. I "created" a solver using the cfdemSolverIB based on the viscoelasticFluidFoam solver included in openFoam-1.6-ext. It solves viscoelastic flows depending on the constitutive equation you need. There are different models already in openFoam-1.6-ext that can be used, like Oldroyd-B model. I am using this new cfdemSolverIBVisco to solve the viscoelastic flow and to calculate the viscous forces coming from the non-Newtonian fluid I use to immersed bodies.

I created a new force instead of the IBDragPerV_ already implemented because my flow is non-Newtonian. So my new force will be:

IBViscoDragPerV_ = rhoField()*(divStress-fvc::grad(p));

Where divStress is the divergence of the stress that sums the divergence of the non-Newtonian stress and the divergence of the Newtonian stress. (I created a new subModel where I called this force just so I do not modify the already existing ShirgaonkarIB model. )

So far so good, the cfdemParticle compiles without errors and the new solver as well. The problem comes when I run the case that is in the tutorial (1 or 2 particles) the solver runs good with all contCumulative_0 in the order of 10^-8 until the simulation kind of explotes and Co number goes to 10^29 ..... When I use 1 particle this happens just before touching the bottom of the box. When I use 2 particles (like in the twoSpheresGlowinskiMPI tutorial) the system crashes a bit further from the bottom. Force increases there rapidly (more due to the grad(P)) at this point before crashing.

I am checking which the cause could be and as soon as I know I will share here. If someone has done something similar or could give a hint/advice it would be really nice!

Have a nice week,

Heliana