tPlot timeseries of kinematic properties - 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
       ---
 (DIR) commit c07cff509228069d18a5d0b0f93a744e503b0167
 (DIR) parent 6c22ac62c66defbacc06fbcee65028f9790e5337
 (HTM) Author: Anders Damsgaard <anders@adamsgaard.dk>
       Date:   Tue, 26 Nov 2019 11:23:24 +0100
       
       Plot timeseries of kinematic properties
       
       Diffstat:
         M python/supraglacial-plots.py        |      59 +++++++++++++++++++++++++++++++
       
       1 file changed, 59 insertions(+), 0 deletions(-)
       ---
 (DIR) diff --git a/python/supraglacial-plots.py b/python/supraglacial-plots.py
       t@@ -16,6 +16,7 @@ dp_dz_vals = np.empty(len(dpdz_values)*len(slope_angle_values))
        slope_angles = np.empty_like(dp_dz_vals)
        fluxes = np.empty_like(dp_dz_vals)
        
       +'''
        j = 0
        for dpdz in dpdz_values:
            outfiles = ''
       t@@ -84,3 +85,61 @@ plt.xlabel('Slope [$^\circ$]')
        plt.ylabel('Specific sediment flux [m$^2$/s]')
        plt.savefig('supraglacial_flux.png')
        plt.savefig('supraglacial_flux.pdf')
       +'''
       +
       +# time series
       +
       +for dpdz in dpdz_values:
       +    for slope_angle in slope_angle_values:
       +                print('### ' + sim.id())
       +        sim = sphere.sim("supraglacial-slope{}-dpdz{}".format(slope_angle, dpdz), fluid=True)
       +        sim.readlast()
       +
       +        N_time = sim.status() - 1
       +        timesteps = np.linspace(1, sim.time_current[0], N_time)
       +        porosity = np.empty(N_time)
       +        velocity = np.empty(N_time)
       +        displacement = np.empty(N_time)
       +        flux = np.empty(N_time)
       +        z = np.zeros(20)
       +        v_x_space_avg = np.empty_like(z)
       +        xsum_space_avg = np.empty_like(z)
       +
       +        for it in np.arange(N_time):
       +            sim.readstep(it)
       +
       +            dz = np.max(sim.x[:,2])/len(z)
       +            for i in range(len(z)):
       +                z[i] = i*dz + 0.5*dz
       +                I = np.nonzero((sim.x[:,2] >= i*dz) & (sim.x[:,2] < (i+1)*dz))
       +                xsum_space_avg[i] = np.mean(sim.xyzsum[I,0])
       +
       +            porosity[it] = np.mean(sim.phi)
       +            velocity[it] = np.mean(sim.vel[:,0])
       +            displacement[it] = np.mean(sim.xyzsum[:,0])
       +            flux[it] = np.trapz(xsum_space_avg/sim.time_current, dx=dz)
       +
       +        ax1 = plt.subplot(4,1,1)
       +        plt.plot(timesteps, porosity, '+-')
       +        plt.ylabel('Porosity [-]')
       +        plt.setp(ax1.get_xticklabels(), visible=False)
       +
       +        ax2 = plt.subplot(4,1,2)
       +        plt.plot(timesteps, velocity, '+-')
       +        plt.ylabel('Avg. velocity [m/s]')
       +        plt.setp(ax2.get_xticklabels(), visible=False)
       +
       +        ax3 = plt.subplot(4,1,3)
       +        plt.plot(timesteps, displacement, '+-')
       +        plt.ylabel('Cumulative displacement [m]')
       +        plt.setp(ax3.get_xticklabels(), visible=False)
       +
       +        ax1 = plt.subplot(4,1,4)
       +        plt.plot(timesteps, flux, '+-')
       +        plt.ylabel('Cumulative flux [m$^2$/s]')
       +        plt.xlabel('Time [s]')
       +
       +        plt.savefig(sim.id() + '-timeseries.png')
       +        plt.savefig(sim.id() + '-timeseries.pdf')
       +        plt.clf()
       +        plt.close()