tshear-test.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
       ---
       tshear-test.py (3592B)
       ---
            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      = True
           11 plots          = True
           12 
           13 # Number of particles
           14 np = 1e4
           15 
           16 # Common simulation id
           17 sim_id = "shear-test"
           18 
           19 # Deviatoric stress [Pa]
           20 Nlist = [80e3]
           21 
           22 ### INITIALIZATION ###
           23 
           24 # New class
           25 init = sphere.sim(np = np, nd = 3, nw = 0, sid = sim_id + "-init")
           26 
           27 # Save radii
           28 init.generateRadii(mean = 0.02)
           29 
           30 # Use default params
           31 init.defaultParams(gamma_n = 100.0, mu_s = 0.6, mu_d = 0.6)
           32 
           33 # Add gravity
           34 init.g[2] = -9.81
           35 
           36 # Periodic x and y boundaries
           37 init.periodicBoundariesXY()
           38 
           39 # Initialize positions in random grid (also sets world size)
           40 hcells = np**(1.0/3.0)
           41 init.initRandomGridPos(gridnum = [hcells, hcells, 1e9])
           42 
           43 # Set duration of simulation
           44 init.initTemporal(total = 5.0)
           45 
           46 if (initialization == True):
           47 
           48     # Run sphere
           49     init.run(dry = True)
           50     init.run()
           51 
           52     if (plots == True):
           53         # Make a graph of energies
           54         init.visualize('energy')
           55 
           56     init.writeVTKall()
           57 
           58     if (rendering == True):
           59         # Render images with raytracer
           60         init.render(method = "angvel", max_val = 0.3, verbose = False)
           61 
           62 
           63 
           64 # For each normal stress, consolidate and subsequently shear the material
           65 for N in Nlist:
           66 
           67     ### CONSOLIDATION ###
           68 
           69     # New class
           70     cons = sphere.sim(np = init.np, nw = 1, sid = sim_id +
           71                       "-cons-N{}".format(N))
           72 
           73     # Read last output file of initialization step
           74     lastf = status(sim_id + "-init")
           75     cons.readbin("../output/" + sim_id + "-init.output{:0=5}.bin".format(lastf), verbose=False)
           76 
           77     # Periodic x and y boundaries
           78     cons.periodicBoundariesXY()
           79 
           80     # Setup consolidation experiment
           81     cons.consolidate(normal_stress = N, periodic = init.periodic)
           82     cons.adaptiveGrid()
           83 
           84 
           85     # Set duration of simulation
           86     cons.initTemporal(total = 1.5)
           87 
           88     """
           89     cons.w_m[0] *= 0.001
           90     cons.mu_s[0] = 0.0
           91     cons.mu_d[0] = 0.0
           92     cons.gamma_wn[0] = 1e4
           93     cons.gamma_wt[0] = 1e4
           94     cons.contactmodel[0] = 1
           95     """
           96 
           97     if (consolidation == True):
           98 
           99         # Run sphere
          100         cons.run(dry = True) # show values, don't run
          101         cons.run() # run
          102 
          103         if (plots == True):
          104             # Make a graph of energies
          105             cons.visualize('energy')
          106             cons.visualize('walls')
          107 
          108         cons.writeVTKall()
          109 
          110         if (rendering == True):
          111             # Render images with raytracer
          112             cons.render(method = "pres", max_val = 2.0*N, verbose = False)
          113 
          114 
          115     ### SHEARING ###
          116 
          117     # New class
          118     shear = sphere.sim(np = cons.np, nw = cons.nw, sid = sim_id +
          119                        "-shear-N{}".format(N))
          120 
          121     # Read last output file of initialization step
          122     lastf = status(sim_id + "-cons-N{}".format(N))
          123     shear.readbin("../output/" + sim_id +
          124                   "-cons-N{}.output{:0=5}.bin".format(N, lastf),
          125                   verbose = False)
          126 
          127     # Periodic x and y boundaries
          128     shear.periodicBoundariesXY()
          129 
          130     # Setup shear experiment
          131     shear.shear(shear_strain_rate = 0.05, periodic = init.periodic)
          132     shear.adaptiveGrid()
          133 
          134     # Set duration of simulation
          135     shear.initTemporal(total = 20.0)
          136 
          137     if (shearing == True):
          138 
          139         # Run sphere
          140         shear.run(dry = True)
          141         shear.run()
          142 
          143         if (plots == True):
          144             # Make a graph of energies
          145             shear.visualize('energy')
          146             shear.visualize('shear')
          147 
          148         shear.writeVTKall()
          149 
          150         if (rendering == True):
          151             # Render images with raytracer
          152             shear.render(method = "pres", max_val = 2.0*N, verbose = False)