tscripts to interpret permeability and diffusivity experiments - 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 96eeffa6fcc7b60bcc548030d18412780338bfe5
 (DIR) parent 98dd28151787ad19ca169cf0c8f76bc378d3108e
 (HTM) Author: Anders Damsgaard <anders.damsgaard@geo.au.dk>
       Date:   Wed, 27 Aug 2014 09:50:25 +0200
       
       scripts to interpret permeability and diffusivity experiments
       
       Diffstat:
         M python/consolidation-curve.py       |       2 +-
         A python/diffusivity-results.py       |      44 +++++++++++++++++++++++++++++++
         M python/permeability-results.py      |      38 +++++++++++++++++--------------
       
       3 files changed, 66 insertions(+), 18 deletions(-)
       ---
 (DIR) diff --git a/python/consolidation-curve.py b/python/consolidation-curve.py
       t@@ -74,8 +74,8 @@ plt.ylabel('Thickness change [m]')
        for c in range(len(c_grad_p_list)):
            plt.semilogx(t[c], H[c], 'o-', label='$c$ = %.2f' % (c_grad_p_list[c]))
        plt.grid()
       -plt.legend(loc=0)
        
       +plt.legend(loc=0)
        plt.tight_layout()
        filename = 'cons-curves.pdf'
        plt.savefig(filename)
 (DIR) diff --git a/python/diffusivity-results.py b/python/diffusivity-results.py
       t@@ -0,0 +1,44 @@
       +#!/usr/bin/env python
       +import matplotlib
       +matplotlib.use('Agg')
       +matplotlib.rcParams.update({'font.size': 18, 'font.family': 'serif'})
       +import os
       +
       +import sphere
       +import numpy
       +import matplotlib.pyplot as plt
       +import diffusivitycalc
       +
       +
       +c_phi = 1.0
       +c_grad_p = 1.0
       +sigma0_list = numpy.array([5.0e3, 10.0e3, 20.0e3, 40.0e3, 80.0e3, 160.0e3])
       +alpha = numpy.empty(len(sigma0_list))
       +
       +dc = diffusivitycalc.DiffusivityCalc()
       +
       +i = 0
       +for sigma0 in sigma0_list:
       +
       +    sim = sphere.sim('cons-sigma0=' + str(sigma0) + '-c_phi=' + \
       +                     str(c_phi) + '-c_grad_p=' + str(c_grad_p), fluid=True)
       +
       +    sim.visualize('walls')
       +    sim.plotLoadCurve()
       +    #sim.writeVTKall()
       +
       +    i += 1
       +
       +
       +plt.xlabel('Normal stress $\\sigma_0$ [kPa]')
       +plt.ylabel('Hydraulic diffusivity $\\alpha$ [m$^2$s$^{-1}$]')
       +#plt.ticklabel_format(style='sci', axis='y', scilimits=(0,0))
       +sigma0 /= 1000.0
       +plt.plot(sigma0, alpha, 'o-k')
       +plt.grid()
       +
       +plt.tight_layout()
       +filename = 'diffusivity-sigma0-vs-alpha.pdf'
       +plt.savefig(filename)
       +print(os.getcwd() + '/' + filename)
       +plt.savefig(filename)
 (DIR) diff --git a/python/permeability-results.py b/python/permeability-results.py
       t@@ -27,22 +27,23 @@ for c_grad_p in cvals:
                sids.append('permeability-dp=' + str(dp) + '-c_phi=' + \
                        str(c_phi) + '-c_grad_p=' + str(c_grad_p))
        
       -    K[c] = numpy.empty(len(sids))
       -    dpdz[c] = numpy.empty_like(K)
       -    Q[c] = numpy.empty_like(K)
       -    phi_bar[c] = numpy.empty_like(K)
       +    K[c] = numpy.zeros(len(sids))
       +    dpdz[c] = numpy.zeros_like(K[c])
       +    Q[c] = numpy.zeros_like(K[c])
       +    phi_bar[c] = numpy.zeros_like(K[c])
            i = 0
        
            for sid in sids:
       -        pc = PermeabilityCalc(sid, plot_evolution=False, print_results=False,
       -                verbose=False)
       -        K[c][i] = pc.conductivity()
       -        pc.findPressureGradient()
       -        pc.findCrossSectionalFlux()
       -        dpdz[c][i] = pc.dPdL[2]
       -        Q[c][i] = pc.Q[2]
       -        pc.findMeanPorosity()
       -        phi_bar[c][i] = pc.phi_bar
       +        if os.path.isfile('../output/' + sid + '.status.dat'):
       +            pc = PermeabilityCalc(sid, plot_evolution=False, print_results=False,
       +                    verbose=False)
       +            K[c][i] = pc.conductivity()
       +            pc.findPressureGradient()
       +            pc.findCrossSectionalFlux()
       +            dpdz[c][i] = pc.dPdL[2]
       +            Q[c][i] = pc.Q[2]
       +            pc.findMeanPorosity()
       +            phi_bar[c][i] = pc.phi_bar
        
                i += 1
        
       t@@ -58,9 +59,12 @@ fig = plt.figure()
        plt.xlabel('Pressure gradient $\\Delta p/\\Delta z$ [kPa m$^{-1}$]')
        plt.ylabel('Hydraulic conductivity $K$ [ms$^{-1}$]')
        plt.ticklabel_format(style='sci', axis='y', scilimits=(0,0))
       -for c in range(len(c_vals)):
       -    dpdz /= 1000.0
       -    plt.plot(dpdz[c], K[c], 'o-', label='$c$ = %.2f' % (c_vals[c]))
       +for c in range(len(cvals)):
       +    dpdz[c] /= 1000.0
       +    #plt.plot(dpdz[c], K[c], 'o-', label='$c$ = %.2f' % (cvals[c]))
       +    #plt.semilogx(dpdz[c], K[c], 'o-', label='$c$ = %.2f' % (cvals[c]))
       +    #plt.semilogy(dpdz[c], K[c], 'o-', label='$c$ = %.2f' % (cvals[c]))
       +    plt.loglog(dpdz[c], K[c], 'o-', label='$c$ = %.2f' % (cvals[c]))
        plt.grid()
        
        #plt.subplot(3,1,2)
       t@@ -75,8 +79,8 @@ plt.grid()
        #plt.plot(dpdz, phi_bar, '+')
        #plt.grid()
        
       +plt.legend()
        plt.tight_layout()
        filename = 'permeability-dpdz-vs-K-vs-c.pdf'
       -plt.savefig(filename)
        print(os.getcwd() + '/' + filename)
        plt.savefig(filename)