← Back to team overview

yade-users team mailing list archive

Re: [Question #683336]: Yade-OpenFoam-coupling delete the particles outside the fluid cells

 

Question #683336 on Yade changed:
https://answers.launchpad.net/yade/+question/683336

Deepak proposed the following answer:
and also in the PyRunner for deletePar, set the iterPeriod=1.

On Wed, Sep 4, 2019 at 11:19 AM Deepak Kn <deepak.kn1990@xxxxxxxxx>
wrote:

> Hello,
>
> There could be several reasons for this :
>
> Have you checked for the stability of the system ? Do particles 'fly' out
> of the domain within a few timesteps ? Have you made sure the courant
> number/time step in the fluid side is an acceptable value (<0.5) ?
> If all these are okay, then I suspect it is with the deletePar function as
> it may not be in sync with coupling.  One way to ensure this is to add an
> if condition to  check whether coupling occurs in the present time step..
>
> def deletePar(self) :
>      if (O.iter%fluidCoupling.dataExchangeInterval==0) :
>                ...
> try it and let me know..
>
>
>
>
>
>
>
>
> On Tue, Sep 3, 2019 at 10:33 AM Anqi H <
> question683336@xxxxxxxxxxxxxxxxxxxxx> wrote:
>
>> Question #683336 on Yade changed:
>> https://answers.launchpad.net/yade/+question/683336
>>
>>     Status: Answered => Open
>>
>> Anqi H is still having a problem:
>> Hi Deepak,
>>
>> Thank you so much for your help so far. I've modified my script based on
>> your comment to check if a body does not have any cohesive bonds in its
>> interactions, if true then delete this sphere from fluidcoupling IDs. In
>> the first openFoam time step, some spheres were deleted as expected,
>> however, after the deletions and in the same foam timestep, it seems the
>> program hangs after 500 yade iterations and in the fluidCoupling step (I
>> have put the particle deletion pyrunner before calling fluidCoupling). I
>> have copied my python script below. I am sorry it's probably difficult
>> to run this script as it requires two other txt files that stored the
>> sphere locations. I was wondering if it was due to some obvious errors
>> in my script that I haven't noticed?
>>
>>
>> from __future__ import print_function
>> import sys
>> from yadeimport import *
>> from yade.utils import *
>> from yade import ymport
>>
>> initMPI()                           #Initialize the mpi environment,
>> always required.
>> fluidCoupling = yade.FoamCoupling();     #Initialize the engine
>> fluidCoupling.getRank();            #part of Initialization.
>>
>>
>> #example of spheres in shear flow : two-way point force coupling
>> class simulation():
>>
>>   def __init__(self):
>>
>>     O.periodic = True;
>>
>>     #proppant properties
>>     FrictAng_p = 0.9
>>     Density_p = 2650
>>     Young_p = 100e6
>>     TensileStr_p=3000
>>     Cohesion_p=3000
>>
>>     Young = 57e8 #apparent modulus
>>     FrictAng = 0.5
>>     Density = 2650
>>     Poisson = 0.28
>>     Cohesion = 38e6 # pa
>>     TensileStr = 38e6 # pa
>>
>>     rock =
>> JCFpmMat(type=1,young=Young,frictionAngle=FrictAng,density=Density,poisson=Poisson,tensileStrength=TensileStr,cohesion=Cohesion,label='rock')
>>     proppant =
>> JCFpmMat(type=2,young=Young_p,frictionAngle=FrictAng_p,density=Density_p,tensileStrength=TensileStr_p,cohesion=Cohesion_p,label='proppant')
>>
>> O.materials.append(JCFpmMat(young=Young_p,frictionAngle=0,density=0,label='wallmat'))
>>
>>     O.materials.append(proppant)
>>     O.materials.append(rock)
>>
>>     rock_assembly =
>>
>> O.bodies.append(ymport.textExt('new_rock.txt','x_y_z_r',shift=Vector3(0,0,0),scale=1,material=rock))
>>
>>     for b in O.bodies:
>>       if b.state.pos[1]>0.017 or b.state.pos[1]<0.0135:
>>         O.bodies.erase(b.id)
>>       #if (b.state.pos[1]>0.0135 and b.state.pos[1]<0.015) :
>>        # b.groupMask=2
>>       else:
>>         b.dynamic = False
>>         b.state.vel=(0, 0, 0)
>>
>>     proppant_assembly =
>> O.bodies.append(ymport.textExt('new_prop.txt','x_y_z_r',shift=Vector3(0,0,0),scale=1,material=proppant,color=(1.00,0.67,0.50)))
>>     print ('length of bodies proppant '+str(len(O.bodies)))
>>
>>     for b in proppant_assembly:
>>       O.bodies[b].groupMask=2
>>       if O.bodies[b].state.pos[0] < 0:
>>         print("found it ")
>>         O.bodies.erase(b)
>>
>>
>>     sphereIDs = [b.id for b in O.bodies if (type(b.shape)==Sphere and
>> b.material.label=='proppant')]
>>
>>     #coupling engine settings
>>     fluidCoupling.setNumParticles(len(sphereIDs))
>>     fluidCoupling.setIdList(sphereIDs)
>>     fluidCoupling.isGaussianInterp=False;  #use pimpleFoamYade for
>> gaussianInterp
>>
>>     # Integrator
>>     newton=NewtonIntegrator(damping=0.0, gravity = (0.0 ,0.0, 0.0))
>>      # add small damping in case of stability issues.. ~ 0.1 max, also
>> note : If gravity is needed, set it in constant/g dir.
>>
>>
>>     O.engines=[
>>   ForceResetter(),
>>   InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=1.2,
>> label="bols"),Bo1_Facet_Aabb()], allowBiggerThanPeriod=True),
>>   InteractionLoop(
>>     [Ig2_Sphere_Sphere_ScGeom(interactionDetectionFactor=1.2,
>> label="ig2s"),Ig2_Facet_Sphere_ScGeom()],
>>     [Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(cohesiveTresholdIteration=1)],
>>
>> [Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(recordCracks=True,label='lawFunctor')]
>>   ),
>>   GlobalStiffnessTimeStepper(timestepSafetyCoefficient=0.7, label = "ts"),
>>         PyRunner(command='sim.deletePar()',iterPeriod=30,
>> label='checkPar'),
>>         fluidCoupling, #to be called after timestepper
>>         PyRunner(command='sim.printMessage()', iterPeriod= 1000,
>> label='outputMessage'),
>>   newton,
>>         #PyRunner(command='sim.deletePar()',iterPeriod=50,
>> label='checkPar'),
>>         VTKRecorder(fileName='yadep/3d-vtk-',mask =
>> 2,recorders=['spheres','colors'],iterPeriod=1000)
>>     ]
>>
>>
>>   def printMessage(self):
>>
>>      print("********************************YADE-ITER = " + str(O.iter)
>> +" **********************************")
>>      if O.iter == 4000:
>>          maxVel = 0.05
>>          for b in O.bodies:
>>              if type(b.shape)==Sphere:
>>                  bodyVel = abs(b.state.vel.norm())
>>                  if bodyVel > maxVel:
>>                      raise ValueError("Body velocity exceeds imposed
>> shear velocity by ", abs(bodyVel-maxVel))
>>
>>   def deletePar(self):
>>     ids = fluidCoupling.getIdList()
>>     for b in ids:
>>       temp=0
>>       for i in O.bodies[b].intrs():
>>         if i.phys.isCohesive==True:
>>           temp+=1
>>       if temp==0:
>>         fluidCoupling.eraseId(b)
>>         print('deleted '+str(b))
>>     print("******YADE-ITER = " + str(O.iter) +" **********")
>>
>>     #if O.bodies[b].state.nbInitBonds==0 or
>> (O.bodies[b].state.nbInitBonds==O.bodies[b].state.nbBrokenBonds):
>>
>>
>>   def steprun(self):
>>     O.step()
>>     bols.aabbEnlargeFactor = 1
>>     ig2s.interactionDetectionFactor = 1
>>     print('-------------------reset aabbEnlargeFactor------------------')
>>
>>   def irun(self,num):
>>       O.run(num,1)
>>
>>
>> if __name__=="__main__":
>>     sim = simulation()
>>     sim.steprun()
>>     sim.irun(5000) #run 5000 iteration and wait
>>     # print("body id = ", O.bodies[34].id)
>>     fluidCoupling.killMPI()
>>
>> import builtins
>> builtins.sim=sim
>>
>> --
>> You received this question notification because your team yade-users is
>> an answer contact for Yade.
>>
>> _______________________________________________
>> Mailing list: https://launchpad.net/~yade-users
>> Post to     : yade-users@xxxxxxxxxxxxxxxxxxx
>> Unsubscribe : https://launchpad.net/~yade-users
>> More help   : https://help.launchpad.net/ListHelp
>>
>

-- 
You received this question notification because your team yade-users is
an answer contact for Yade.