Commit c77308a6 authored by Nathan Coppin's avatar Nathan Coppin
Browse files

allowing 0 particles

parent 16668f34
Pipeline #7496 passed with stage
in 1 minute and 59 seconds
......@@ -270,7 +270,7 @@ class ParticleProblem :
tol -- optional argument defining the interpenetration tolerance to stop the NLGS iterations of the NSCD
"""
self._lib.particleProblemIterate.restype = c_int
return self._lib.particleProblemIterate(self._p, c_double(np.min(self.r())), c_double(dt), c_double(tol), c_int(-1),c_int(force_motion),_np2c(forces))
return self._lib.particleProblemIterate(self._p, c_double(np.min(self.r()) if self.r().size != 0 else tol), c_double(dt), c_double(tol), c_int(-1),c_int(force_motion),_np2c(forces))
def get_contacts(self,contact_type) :
ctype = {"particle_particle":0,"particle_disk":1,"particle_segment":2,"particle_triangle":3}[contact_type]
......
......@@ -16,14 +16,16 @@ def _advance_particles(particles, f, dt, min_nsub,contact_tol,iteration=0,after_
# Compute free velocities of the grains
if f is None:
f = np.zeros_like(particles.velocity())
vn = particles.velocity() + f*dt / particles.mass()
# Estimation of the solid sub-time steps
if particles.r() is not None:
if particles.r() is not None and particles.r().size != 0:
vn = particles.velocity() + f*dt / particles.mass()
print(vn,particles.r())
vmax = np.max(np.linalg.norm(vn,axis=1))
nsub = max(min_nsub, int(np.ceil((vmax * dt * 8)/min(particles.r()))))
else:
nsub = 1
print("NSUB", nsub,"VMAX",vmax, "VMAX * dt", vmax * dt, "r", min(particles.r()) if particles.r() is not None else 0 )
vmax = 0
print("NSUB", nsub,"VMAX",vmax, "VMAX * dt", vmax * dt, "r", min(particles.r()) if (particles.r() is not None and particles.r().size !=0) else 0 )
#If time step was split too many times, ignore the convergence and move the grains
if iteration == max_split:
for i in range(nsub) :
......@@ -119,7 +121,7 @@ def iterate(fluid, particles, dt, min_nsub=1, contact_tol=1e-8, external_particl
if particles is not None and external_particles_forces is not None:
if external_particles_forces.shape!=particles.velocity().shape:
raise ValueError("external_particles_forces must have shape (number of particles,dimension!")
if (particles is None or particles.r() is None) and fluid is not None:
if (particles is None or particles.r() is None or particles.r().size == 0) and fluid is not None:
fluid.implicit_euler(dt)
# For two fluids flows, advance the concentration field with the fluid velocity.
# The number of sub-time steps for the advection is automatically computed.
......@@ -127,7 +129,7 @@ def iterate(fluid, particles, dt, min_nsub=1, contact_tol=1e-8, external_particl
fluid.advance_concentration(dt)
if particles is not None and fluid is None:
_advance_particles(particles,external_particles_forces,dt,min_nsub,contact_tol,after_sub_iter=after_sub_iter,max_split=max_split)
if fluid is not None and particles is not None and particles.r() is not None:
if fluid is not None and particles is not None and particles.r() is not None and particles.r().size !=0:
if fixed_grains:
fluid.implicit_euler(dt)
fluid.set_particles(particles.mass(), particles.volume(), particles.position(), particles.velocity(),-fluid.compute_node_force(dt))
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment