Commit ad3473ed authored by Jonathan Lambrechts's avatar Jonathan Lambrechts
Browse files

remove state in all .py files

parent 8668533f
Pipeline #8535 failed with stages
in 1 minute and 33 seconds
......@@ -74,14 +74,12 @@ p = scontact.ParticleProblem(2,friction_enabled = True,rotation_enabled=True)
p2 = scontact.ParticleProblem(2,friction_enabled = True,rotation_enabled=True)
p.read_vtk(outputdir, 0)
p2.read_vtk(outputdir, 0)
state = p.state()
state2 = p2.state()
p2.remove_particles_flag((state2.x[:,1] + p2.r()[:,0] <1.15))
p2.remove_particles_flag((state2.x[:,1] + p2.r()[:,0] >1.1))
p2.remove_particles_flag((state2.x[:,0] + p2.r()[:,0] <0.15))
p2.remove_particles_flag((state2.x[:,0] + p2.r()[:,0] >-0.15))
state2.x[:,1] -= 0.6
p2.set_state(state2)
x2 = p2.position()
p2.remove_particles_flag((x2[:,1] + p2.r()[:,0] <1.15))
p2.remove_particles_flag((x2[:,1] + p2.r()[:,0] >1.1))
p2.remove_particles_flag((x2[:,0] + p2.r()[:,0] <0.15))
p2.remove_particles_flag((x2[:,0] + p2.r()[:,0] >-0.15))
p2.position()[:,1] -= 0.6
p.clear_boundaries()
p.load_msh_boundaries("mesh.msh", ["Inner"], material = "Steel")
p.load_msh_boundaries("mesh.msh", ["Right","Left","RightUp","RightDown","LeftUp","LeftDown"], material = "Plexi")
......@@ -103,13 +101,13 @@ while t < tEnd :
p.load_msh_boundaries("mesh.msh", ["Inner"],material="Steel")
p.load_msh_boundaries("mesh.msh", ["Right","Left","RightUp","RightDown","LeftUp","LeftDown"],material="Plexi")
#Adding new particles
if (max(p.state().x[:,1])+max(p.r()) <0.5) and t > 3 and t < 25:
p.add_particles(p2.state().x,p2.r(),p2.mass(),v=p2.state().v,tag="Sand",forced=p2.forced_flag(),contact_forces=p2.contact_forces())
if (max(p.position()[:,1])+max(p.r()) <0.5) and t > 3 and t < 25:
p.add_particles(p2.position(),p2.r(),p2.mass(),v=p2.velocity(),tag="Sand",forced=p2.forced_flag(),contact_forces=p2.contact_forces())
#Solving the contacts
F = np.zeros(2)
time_integration.iterate(None,p,dt,min_nsub=1,contact_tol=1e-7,external_particles_forces=g*p.mass(),after_sub_iter= lambda nsub : accumulate(F,nsub))
#Removing particles outside the hopper
p.remove_particles_flag( (p.state().x[:,1] >-0.6))
p.remove_particles_flag( (p.position()[:,1] >-0.6))
#Computing mean velocity and writing to file
with open(filename,"a") as file1:
file1.write(str(F[0])+";"+str(F[1])+";"
......
......@@ -82,13 +82,12 @@ else :
p2 = scontact.ParticleProblem(3,friction_enabled = True,rotation_enabled=True)
p.read_vtk(outputdir, 0)
p2.read_vtk(outputdir, 0)
p2.remove_particles_flag((p2.state().x[:,1] + p2.r()[:,0] <1.25))
p2.remove_particles_flag((p2.state().x[:,1] + p2.r()[:,0] >1.1))
p2.remove_particles_flag((p2.state().x[:,0] + p2.r()[:,0] <0.18))
p2.remove_particles_flag((p2.state().x[:,0] + p2.r()[:,0] >-0.18))
state2 = p2.state()
state2.x[:,1] -= 0.5
p2.set_state(state2)
x2 = p2.position()
p2.remove_particles_flag((x2[:,1] + p2.r()[:,0] <1.25))
p2.remove_particles_flag((x2[:,1] + p2.r()[:,0] >1.1))
p2.remove_particles_flag((x2[:,0] + p2.r()[:,0] <0.18))
p2.remove_particles_flag((x2[:,0] + p2.r()[:,0] >-0.18))
p2.position()[:,1] -= 0.5
p.clear_boundaries()
p.load_msh_boundaries("mesh" + ".msh", ["Inner","Right","Left","Front","Rear","LockDown"], material = "Steel")
......@@ -105,13 +104,13 @@ while t < tEnd :
p.clear_boundaries()
p.load_msh_boundaries("mesh" + ".msh", ["Inner","Right","Left","Front","Rear"],material="Steel")
#Adding new particles
if (max(state.x[:,1])+max(p.r()) <0.6) and t > 10 and t < 25:
p.add_particles(p2.state().x,p2.r(),p2.mass(),v=p2.state().v,tag="Sand",forced=p2.forced_flag(),contact_forces=p2.contact_forces())
if (max(p.position()[:,1])+max(p.r()) <0.6) and t > 10 and t < 25:
p.add_particles(p2.position(),p2.r(),p2.mass(),v=p2.velocity(),tag="Sand",forced=p2.forced_flag(),contact_forces=p2.contact_forces())
#Solving the contacts
F = np.zeros(3)
time_integration.iterate(None,p,dt,min_nsub=1,contact_tol=1e-6,external_particles_forces=g*p.mass(),after_sub_iter= lambda subdt : accumulate(F))
#Removing particles outside the hopper
p.remove_particles_flag( (p.state().x[:,1] >-0.6))
p.remove_particles_flag( (p.position()[:,1] >-0.6))
#Computing mean velocity and writing to file
with open(filename,"a") as file1:
file1.write(str(F[0])+";"+str(F[1])+";"+str(F[2])+";"
......
......@@ -104,7 +104,7 @@ fluid.set_open_boundary("Top",pressure=0,porosity=1,concentration=1)
fluid.set_open_boundary("Bottom",velocity=[0,outerBndV],porosity=1,concentration=1)
#fluid.set_strong_boundary("Injection",1,outerBndV)
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
#set initial_condition
c = np.zeros_like(fluid.coordinates()[:,1])
......@@ -115,7 +115,7 @@ fluid.export_vtk(outputdir,0,0)
ii = 0
tic = time.time()
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
while t < tEnd :
#Adaptation of the mesh. Args are minimal mesh radius, maximal mesh radius and number of elements
#Fluid solver
......@@ -133,7 +133,7 @@ while t < tEnd :
for i in range(nsub) :
p.iterate(dt/nsub, forces)
t += dt
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
#Output files writting
if ii %outf == 0 :
ioutput = int(ii/outf) + 1
......
......@@ -88,9 +88,9 @@ fluid.set_open_boundary("Bottom",velocity=[0,outerBndV])
ii = 0
t = 0
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
fluid.export_vtk(outputdir,0,0)
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
tic = time.clock()
while t < tEnd :
#Fluid solver
......
......@@ -89,7 +89,7 @@ fluid.set_wall_boundary("Lateral")
fluid.set_open_boundary("Top",pressure=0, concentration=1)
fluid.set_open_boundary("Injection",velocity=[0, outerBndV], concentration=1)
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
#set initial_condition
print(fluid.solution().shape, fluid.coordinates().shape)
......@@ -100,7 +100,7 @@ fluid.export_vtk(outputdir,0,0)
ii = 0
tic = time.time()
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
while t < tEnd :
#Adaptation of the mesh. Args are minimal mesh radius, maximal mesh radius and number of elements
#Fluid solver
......@@ -118,7 +118,7 @@ while t < tEnd :
for i in range(nsub) :
p.iterate(dt/nsub, forces)
t += dt
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
#Output files writting
if ii %outf == 0 :
ioutput = int(ii/outf) + 1
......
......@@ -75,7 +75,7 @@ fluid.import_vtk(outputrel+"/fluid_%05d.vtu"%iReload)
rel = 1
fluid.set_particles(p.mass(), p.volume(), p.state(),rel)
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity(),rel)
ii = (iReload-1)*outf+1
t=ii*dt
tic = time.clock()
......@@ -99,7 +99,7 @@ while t < tEnd :
#Contact solver
for i in range(nsub) :
p.iterate(dt/nsub, forces)
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
#Output files writting
if ii %outf1 == 0 :
ioutput = int(ii1/outf1) + 1
......@@ -109,4 +109,4 @@ while t < tEnd :
t += dt
ii += 1
ii1 += 1
print("%i : %.2g/%.2g (cpu %.6g)" % (ii, t, tEnd, time.time() - tic))
\ No newline at end of file
print("%i : %.2g/%.2g (cpu %.6g)" % (ii, t, tEnd, time.time() - tic))
......@@ -154,7 +154,7 @@ fluid.set_strong_boundary("Lateral",0,0)
#fluid.set_strong_boundary("Lateral",1,0)
if not no_particle :
fluid.set_particles(p.mass(), p.volume(), p.state(), p.contact_forces(),init=True)
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity(), p.contact_forces(),init=True)
#set initial_condition
#c = np.zeros_like(fluid.coordinates()[:,1])
......@@ -204,7 +204,7 @@ while t < tEnd :
p.iterate(dt/nsub, f0)
p.position()[:,:] = pp0
p.contact_forces()[:,:] = cp0
fluid.set_particles(p.mass(), p.volume(), p.state(),p.contact_forces())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity(),p.contact_forces())
fluid.solution()[:,:] = sf0
#fluid.concentration_dg()[:] = af0
......@@ -227,7 +227,7 @@ while t < tEnd :
for i in range(nsub) :
#tol = 5e-9
p.iterate(dt/nsub, (alpha*f0+(1-alpha)*f1))
fluid.set_particles(p.mass(), p.volume(), p.state(),p.contact_forces())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity(),p.contact_forces())
t += dt
flux = fluid.volume_flux("Injection")
......
......@@ -152,7 +152,7 @@ fluid.set_open_boundary("Injection",pressure=inBndP,porosity=1)
#fluid.set_strong_boundary("Lateral",0,0)
if not no_particle :
fluid.set_particles(p.mass(), p.volume(), p.state(),0*p.velocity())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity(),0*p.velocity())
#set initial_condition
#c = np.zeros_like(fluid.coordinates()[:,1])
......@@ -166,7 +166,7 @@ fluid.solution()[:,2] = (H-fluid.coordinates()[:,1])*(-g)*rho0
ii = 0
tic = time.time()
if not no_particle :
fluid.set_particles(p.mass(), p.volume(), p.state(),0*p.velocity())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity(),0*p.velocity())
fluid.export_vtk(outputdir,0,0)
ioutput = 0
......@@ -195,7 +195,7 @@ while t < tEnd :
vn = p.velocity().copy()
p.iterate(dt/nsub, forces)
fc = (p.velocity()-vn)*p.mass()/(dt/nsub) - forces
fluid.set_particles(p.mass(), p.volume(), p.state(),fc)
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity(),fc)
t += dt
flux = fluid.volume_flux("Injection")
injector.iterate((Q+flux*w)*rho0*dt)
......
......@@ -151,7 +151,7 @@ fluid.set_open_boundary("Injection",pressure=inBndP,porosity=1,concentration=1)
#fluid.set_strong_boundary("Lateral",0,0)
if not no_particle :
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
#set initial_condition
c = np.zeros_like(fluid.coordinates()[:,1])
......@@ -166,7 +166,7 @@ ii = 0
tic = time.time()
if not no_particle :
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
fluid.export_vtk(outputdir,0,0)
ioutput = 0
......@@ -193,7 +193,7 @@ while t < tEnd :
#Contact solver
for i in range(nsub) :
p.iterate(dt/nsub, forces)
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
t += dt
flux = fluid.volume_flux("Injection")
injector.iterate((Q+flux*w)*rho0*dt)
......
......@@ -145,7 +145,7 @@ fluid.set_strong_boundary("Top",2,0)
#fluid.set_strong_boundary("Lateral",0,0)
if not no_particle :
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
#set initial_condition
c = np.zeros_like(fluid.coordinates()[:,1])
......@@ -158,7 +158,7 @@ ii = 0
tic = time.time()
if not no_particle :
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
while t < tEnd :
#Adaptation of the mesh. Args are minimal mesh radius, maximal mesh radius and number of elements
......@@ -177,7 +177,7 @@ while t < tEnd :
#Contact solver
for i in range(nsub) :
p.iterate(dt/nsub, forces)
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
t += dt
flux = fluid.volume_flux("Injection")
injector.iterate((Q+flux*w)*rho0*dt)
......
......@@ -93,7 +93,7 @@ fluid.load_msh("mesh.msh")
fluid.set_strong_boundary("Outer",0,lambda x : (x[:,1]/rout)*-v )
fluid.set_strong_boundary("Outer",1,lambda x : (-x[:,0]/rout)*-v)
fluid.set_wall_boundary("Top",0)
fluid.set_particles(p.mass(), p.volume(), p.state(),p.contact_forces(),init=True)
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity(),p.contact_forces(),init=True)
fluid.solution()[:,2] = g[1]*rho*(fluid.coordinates()[:,1] - rout/2)
fluid.export_vtk(outputdir, 0.,0)
#
......
......@@ -74,7 +74,7 @@ fluid.set_strong_boundary("Box",1,0)
fluid.set_open_boundary("Left",velocity=[0.001, 0])
fluid.set_open_boundary("Right",pressure=0)
fluid.set_wall_boundary("Box")
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
fluid.export_vtk(outputdir,0,0)
tic = time.time()
......@@ -87,7 +87,7 @@ while t < tEnd :
print("NSUB", nsub,"VMAX",vmax, "VMAX * dt", vmax * dt, "r", min(p.r()))
for i in range(nsub) :
p.iterate(dt/nsub, forces)
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
t += dt
if ii %outf == 0 :
ioutput = int(ii/outf) + 1
......
......@@ -97,7 +97,7 @@ fluid = fluid.FluidProblem(2,g,nu*rho,rho)
fluid.load_msh("mesh.msh")
fluid.set_wall_boundary("Top", pressure=0)
fluid.set_wall_boundary("Box")
fluid.set_particles(p.mass(), p.volume(), p.state(), p.contact_forces(), init=True)
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity(), p.contact_forces(), init=True)
fluid.export_vtk(outputdir,0,0)
tic = time.time()
......
......@@ -105,7 +105,7 @@ fluid.set_open_boundary("BottomOut",velocity=[0, outerBndV])
fluid.set_open_boundary("TopOut",pressure=0)
fluid.set_open_boundary("Top",pressure=0)
# Set location of the particles in the mesh and compute the porosity in each computation cell
fluid.set_particles(p.mass(), p.volume(), p.state(), p.contact_forces(), init=True)
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity(), p.contact_forces(), init=True)
fluid.export_vtk(outputdir,0,0)
tic = time.time()
......
......@@ -105,7 +105,7 @@ fluid.set_wall_boundary("X")
fluid.set_wall_boundary("Z")
fluid.set_wall_boundary("Bottom")
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
fluid.export_vtk(outputdir,0,0)
tic = time.time()
......@@ -129,7 +129,7 @@ while t < tEnd :
tol = 1e-6
p.iterate(dt/nsub, forces, tol)
#Localisation of the grains in the fluid
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
t += dt
#Output files writting
if ii %outf == 0 :
......
......@@ -72,12 +72,12 @@ fluid.set_strong_boundary("Bottom",2,0)
#fluid.set_weak_boundary("Top","Inflow")
#fluid.set_weak_boundary("XVel","Wall")
#fluid.set_particles(p.mass(), p.volume(), p.state())
#fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
fluid.export_vtk(outputdir,0,0)
tic = time.clock()
fluid.export_vtk(outputdir,0,0)
#fluid.set_particles(p.mass(), p.volume(), p.state())
#fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
while t < tEnd :
fluid.implicit_euler(dt)
#forces = fluid.compute_node_force(dt)
......@@ -87,7 +87,7 @@ while t < tEnd :
#print("NSUB", nsub,"VMAX",vmax, "VMAX * dt", vmax * dt, "r", min(p.r()))
#for i in range(nsub) :
# p.iterate(dt/nsub, forces,1e-7)
#fluid.set_particles(p.mass(), p.volume(), p.state())
#fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
t += dt
if ii %outf == 0 :
ioutput = int(ii/outf) + 1
......
......@@ -94,12 +94,12 @@ fluid.set_wall_boundary("XVel")
L=.33/2
#fluid.solution()[:,:] = np.zeros_like(fluid.solution()[:,:])#0*1/L**4 *(fluid.coordinates()[:,0]-L)*(fluid.coordinates()[:,0]+L)*(fluid.coordinates()[:,2]-L)*(fluid.coordinates()[:,2]+L)*np.sin(np.pi*1/2*dt)
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
fluid.export_vtk(outputdir,0,0)
tic = time.time()
fluid.export_vtk(outputdir,0,0)
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
while t < tEnd :
fluid.implicit_euler(dt)
forces = fluid.compute_node_force(dt)
......@@ -110,7 +110,7 @@ while t < tEnd :
print("NSUB", nsub,"VMAX",vmax, "VMAX * dt", vmax * dt, "r", min(p.r()))
for i in range(nsub) :
p.iterate(dt/nsub, forces,1e-7)
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
t += dt
if ii %outf == 0 :
ioutput = int(ii/outf) + 1
......
......@@ -92,11 +92,11 @@ fluid.set_wall_boundary("XVel")
L=.33/2
#fluid.solution()[:,:] = np.zeros_like(fluid.solution()[:,:])#0*1/L**4 *(fluid.coordinates()[:,0]-L)*(fluid.coordinates()[:,0]+L)*(fluid.coordinates()[:,2]-L)*(fluid.coordinates()[:,2]+L)*np.sin(np.pi*1/2*dt)
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
fluid.export_vtk(outputdir,0,0)
tic = time.time()
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
while t < tEnd :
fluid.implicit_euler(dt)
forces = fluid.compute_node_force(dt)
......@@ -107,7 +107,7 @@ while t < tEnd :
print("NSUB", nsub,"VMAX",vmax, "VMAX * dt", vmax * dt, "r", min(p.r()))
for i in range(nsub) :
p.iterate(dt/nsub, forces,1e-7)
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
t += dt
if ii %outf == 0 :
ioutput = int(ii/outf) + 1
......
......@@ -96,11 +96,11 @@ fluid.set_wall_boundary("XVel")
L=.33/2
#fluid.solution()[:,:] = np.zeros_like(fluid.solution()[:,:])#0*1/L**4 *(fluid.coordinates()[:,0]-L)*(fluid.coordinates()[:,0]+L)*(fluid.coordinates()[:,2]-L)*(fluid.coordinates()[:,2]+L)*np.sin(np.pi*1/2*dt)
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
fluid.export_vtk(outputdir,0,0)
tic = time.time()
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
while t < tEnd :
fluid.implicit_euler(dt)
forces = fluid.compute_node_force(dt)
......@@ -111,7 +111,7 @@ while t < tEnd :
print("NSUB", nsub,"VMAX",vmax, "VMAX * dt", vmax * dt, "r", min(p.r()))
for i in range(nsub) :
p.iterate(dt/nsub, forces,1e-7)
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
t += dt
if ii %outf == 0 :
ioutput = int(ii/outf) + 1
......
......@@ -105,13 +105,13 @@ fluid.set_wall_boundary("InnerRight")
fluid.set_wall_boundary("Top",pressure=0)
fluid.set_wall_boundary("Bottom")
fluid.set_wall_boundary("Vertical")
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
fluid.export_vtk(outputdir,0,0)
tic = time.clock()
#computation loop
while t < tEnd :
fluid.set_particles(p.mass(), p.volume(), p.state())
fluid.set_particles(p.mass(), p.volume(), p.position(), p.velocity())
fluid.implicit_euler(dt)
forces = fluid.compute_node_force(dt)
vn = p.velocity() + forces * dt / p.mass()
......
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