tjp-long-shear-soft.py - sphere - GPU-based 3D discrete element method algorithm with optional fluid coupling
 (HTM) git clone git://src.adamsgaard.dk/sphere
 (DIR) Log
 (DIR) Files
 (DIR) Refs
 (DIR) LICENSE
       ---
       tjp-long-shear-soft.py (3635B)
       ---
            1 #!/usr/bin/env python
            2 
            3 # Import sphere functionality
            4 import sphere
            5 
            6 ### EXPERIMENT SETUP ###
            7 initialization = True
            8 consolidation = True
            9 shearing = True
           10 rendering = False
           11 plots = True
           12 
           13 # CUDA device to use
           14 device = 0
           15 
           16 # Number of particles
           17 np = 1e4
           18 
           19 # Common simulation id
           20 sim_id = "jp-long-shear-soft"
           21 
           22 # Deviatoric stress [Pa]
           23 devslist = [100e3]
           24 
           25 ### INITIALIZATION ###
           26 
           27 # New class
           28 init = sphere.sim(np=np, nd=3, nw=0, sid=sim_id + "-init")
           29 
           30 # Save radii
           31 init.generateRadii(mean=0.01)
           32 
           33 # Use default params
           34 init.defaultParams(k_n=1.16e7, k_t=1.16e7, gamma_n=100.0, mu_s=0.6, mu_d=0.6)
           35 
           36 # Add gravity
           37 init.g[2] = -9.81
           38 
           39 # Periodic x and y boundaries
           40 init.periodicBoundariesXY()
           41 
           42 # Initialize positions in random grid (also sets world size)
           43 hcells = np**(1.0/3.0)
           44 init.initRandomGridPos(gridnum=[hcells, hcells, 1e9])
           45 
           46 # Set duration of simulation
           47 init.initTemporal(total=10.0, epsilon=0.07)
           48 
           49 if (initialization):
           50 
           51     # Run sphere
           52     init.run(dry=True)
           53     init.run(device=device)
           54 
           55     if (plots):
           56         # Make a graph of energies
           57         init.visualize('energy')
           58 
           59     init.writeVTKall()
           60 
           61     if (rendering):
           62         # Render images with raytracer
           63         init.render(method="angvel", max_val=0.3, verbose=False)
           64 
           65 
           66 # For each normal stress, consolidate and subsequently shear the material
           67 for devs in devslist:
           68 
           69     ### CONSOLIDATION ###
           70 
           71     # New class
           72     cons = sphere.sim(
           73         np=init.np,
           74         nw=1,
           75         sid=sim_id +
           76         "-cons-devs{}".format(devs))
           77 
           78     # Read last output file of initialization step
           79     lastf = sphere.status(sim_id + "-init")
           80     cons.readbin(
           81         "../output/" +
           82         sim_id +
           83         "-init.output{:0=5}.bin".format(lastf),
           84         verbose=False)
           85 
           86     # Periodic x and y boundaries
           87     cons.periodicBoundariesXY()
           88 
           89     # Setup consolidation experiment
           90     cons.consolidate(normal_stress=devs)
           91 
           92     # Set duration of simulation
           93     cons.initTemporal(total=3.0, epsilon=0.07)
           94 
           95     """
           96     cons.w_m[0] *= 0.001
           97     cons.mu_s[0] = 0.0
           98     cons.mu_d[0] = 0.0
           99     cons.gamma_wn[0] = 1e4
          100     cons.gamma_wt[0] = 1e4
          101     cons.contactmodel[0] = 1
          102     """
          103 
          104     if (consolidation):
          105 
          106         # Run sphere
          107         cons.run(dry=True)  # show values, don't run
          108         cons.run(device=device)  # run
          109 
          110         if (plots):
          111             # Make a graph of energies
          112             cons.visualize('energy')
          113             cons.visualize('walls')
          114 
          115         cons.writeVTKall()
          116 
          117         if (rendering):
          118             # Render images with raytracer
          119             cons.render(method="pres", max_val=2.0*devs, verbose=False)
          120 
          121     ### SHEARING ###
          122 
          123     # New class
          124     shear = sphere.sim(
          125         np=cons.np,
          126         nw=cons.nw,
          127         sid=sim_id +
          128         "-shear-devs{}".format(devs))
          129 
          130     # Read last output file of initialization step
          131     lastf = sphere.status(sim_id + "-cons-devs{}".format(devs))
          132     shear.readbin(
          133         "../output/" + sim_id + "-cons-devs{}.output{:0=5}.bin".format(devs,
          134                                                                        lastf), verbose=False)
          135 
          136     # Periodic x and y boundaries
          137     shear.periodicBoundariesXY()
          138 
          139     # Setup shear experiment
          140     shear.shear(shear_strain_rate=0.05)
          141 
          142     # Set duration of simulation
          143     shear.initTemporal(total=240.0, epsilon=0.07)
          144 
          145     if (shearing):
          146 
          147         # Run sphere
          148         shear.run(dry=True)
          149         shear.run(device=device)
          150 
          151         if (plots):
          152             # Make a graph of energies
          153             shear.visualize('energy')
          154             shear.visualize('shear')
          155 
          156         shear.writeVTKall()
          157 
          158         if (rendering):
          159             # Render images with raytracer
          160             shear.render(method="pres", max_val=2.0*devs, verbose=False)