GopherProxy

	tadd script to test deformation depth with/without rotation - sphere - GPU-base…
	git clone git://src.adamsgaard.dk/sphere
	Log
	Files
	Refs
	LICENSE
	---
	commit 29ef77af0b602c4b4402f6c86b0673e266ebfe39
	parent 9c3b2327798617bfd18112c2cf55ad3ceb26fa11
	Author: Anders Damsgaard <[email protected]>
	Date: Thu, 14 Apr 2022 16:04:53 +0200

	add script to test deformation depth with/without rotation

	Diffstat:
	A python/deformationdepth.py \| 196 +++++++++++++++++++++++++++++…

	1 file changed, 196 insertions(+), 0 deletions(-)
	---
	diff --git a/python/deformationdepth.py b/python/deformationdepth.py
	t@@ -0,0 +1,196 @@
	+#!/usr/bin/env python
	+
	+# Import sphere functionality
	+import sphere
	+
	+### EXPERIMENT SETUP ###
	+initialization = True
	+consolidation = True
	+shearing_norot = True
	+shearing = False
	+rendering = False
	+plots = True
	+
	+# Number of particles
	+np = 1e4
	+
	+# Common simulation id
	+sim_id = "defdepth"
	+
	+# Deviatoric stress [Pa]
	+Nlist = [100e3, 10e3, 20e3, 50e3]
	+
	+### INITIALIZATION ###
	+
	+# New class
	+init = sphere.sim(np = np, nd = 3, nw = 0, sid = sim_id + "-init")
	+
	+# Save radii
	+init.generateRadii(mean = 0.02)
	+
	+# Use default params
	+init.defaultParams(gamma_n = 100.0, mu_s = 0.6, mu_d = 0.6)
	+
	+# Add gravity
	+init.g[2] = -9.81
	+
	+# Periodic x and y boundaries
	+init.periodicBoundariesXY()
	+
	+# Initialize positions in random grid (also sets world size)
	+hcells = np**(1.0/3.0)
	+init.initRandomGridPos(gridnum = [hcells, hcells, 1e9])
	+
	+# Set duration of simulation
	+init.initTemporal(total = 5.0)
	+
	+if (initialization == True):
	+
	+ # Run sphere
	+ init.run(dry = True)
	+ init.run()
	+
	+ if (plots == True):
	+ # Make a graph of energies
	+ init.visualize('energy')
	+
	+ init.writeVTKall()
	+
	+ if (rendering == True):
	+ # Render images with raytracer
	+ init.render(method = "angvel", max_val = 0.3, verbose = False)
	+
	+
	+
	+# For each normal stress, consolidate and subsequently shear the material
	+for N in Nlist:
	+
	+ ### CONSOLIDATION ###
	+
	+ # New class
	+ cons = sphere.sim(np = init.np, nw = 1, sid = sim_id +
	+ "-cons-N{}".format(N))
	+
	+ # Read last output file of initialization step
	+ lastf = status(sim_id + "-init")
	+ cons.readbin("../output/" + sim_id + "-init.output{:0=5}.bin".format(lastf…
	+
	+ # Periodic x and y boundaries
	+ cons.periodicBoundariesXY()
	+
	+ # Setup consolidation experiment
	+ cons.consolidate(normal_stress = N, periodic = init.periodic)
	+ cons.adaptiveGrid()
	+
	+
	+ # Set duration of simulation
	+ cons.initTemporal(total = 1.5)
	+
	+ """
	+ cons.w_m[0] *= 0.001
	+ cons.mu_s[0] = 0.0
	+ cons.mu_d[0] = 0.0
	+ cons.gamma_wn[0] = 1e4
	+ cons.gamma_wt[0] = 1e4
	+ cons.contactmodel[0] = 1
	+ """
	+
	+ if (consolidation == True):
	+
	+ # Run sphere
	+ cons.run(dry = True) # show values, don't run
	+ cons.run() # run
	+
	+ if (plots == True):
	+ # Make a graph of energies
	+ cons.visualize('energy')
	+ cons.visualize('walls')
	+
	+ cons.writeVTKall()
	+
	+ if (rendering == True):
	+ # Render images with raytracer
	+ cons.render(method = "pres", max_val = 2.0*N, verbose = False)
	+
	+
	+ ### SHEARING WITHOUT ROTATION ###
	+
	+ # New class
	+ shear_nrt = sphere.sim(np = cons.np, nw = cons.nw, sid = sim_id +
	+ "-shear_nrt-N{}".format(N))
	+
	+ # Read last output file of initialization step
	+ lastf = status(sim_id + "-cons-N{}".format(N))
	+ shear_nrt.readbin("../output/" + sim_id +
	+ "-cons-N{}.output{:0=5}.bin".format(N, lastf),
	+ verbose = False)
	+
	+ # Periodic x and y boundaries
	+ shear_nrt.periodicBoundariesXY()
	+
	+ # Setup shear_nrt experiment
	+ shear_nrt.shear(shear_strain_rate = 0.05, periodic = init.periodic)
	+ shear_nrt.adaptiveGrid()
	+
	+ # Set duration of simulation
	+ shear_nrt.initTemporal(total = 20.0)
	+
	+ # Disable rotation for regular grains
	+ shear_nrt.fixvel[numpy.nonzero(shear_nrt.fixvel == 0.0)] = -10.0
	+
	+ if (shearing == True):
	+
	+ # Run sphere
	+ shear_nrt.run(dry = True)
	+ shear_nrt.run()
	+
	+ if (plots == True):
	+ # Make a graph of energies
	+ shear_nrt.visualize('energy')
	+ shear_nrt.visualize('shear')
	+
	+ shear_nrt.writeVTKall()
	+
	+ if (rendering == True):
	+ # Render images with raytracer
	+ shear_nrt.render(method = "pres", max_val = 2.0*N, verbose = False)
	+
	+
	+ ### SHEARING ###
	+
	+ # New class
	+ shear = sphere.sim(np = cons.np, nw = cons.nw, sid = sim_id +
	+ "-shear-N{}".format(N))
	+
	+ # Read last output file of initialization step
	+ lastf = status(sim_id + "-cons-N{}".format(N))
	+ shear.readbin("../output/" + sim_id +
	+ "-cons-N{}.output{:0=5}.bin".format(N, lastf),
	+ verbose = False)
	+
	+ # Periodic x and y boundaries
	+ shear.periodicBoundariesXY()
	+
	+ # Setup shear experiment
	+ shear.shear(shear_strain_rate = 0.05, periodic = init.periodic)
	+ shear.adaptiveGrid()
	+
	+ # Set duration of simulation
	+ shear.initTemporal(total = 20.0)
	+
	+ if (shearing == True):
	+
	+ # Run sphere
	+ shear.run(dry = True)
	+ shear.run()
	+
	+ if (plots == True):
	+ # Make a graph of energies
	+ shear.visualize('energy')
	+ shear.visualize('shear')
	+
	+ shear.writeVTKall()
	+
	+ if (rendering == True):
	+ # Render images with raytracer
	+ shear.render(method = "pres", max_val = 2.0*N, verbose = False)