tdiffusivity-results.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
       ---
       tdiffusivity-results.py (2024B)
       ---
            1 #!/usr/bin/env python
            2 import matplotlib
            3 matplotlib.use('Agg')
            4 matplotlib.rcParams.update({'font.size': 18, 'font.family': 'serif'})
            5 import os
            6 import shutil
            7 
            8 import sphere
            9 import numpy
           10 import matplotlib.pyplot as plt
           11 #import diffusivitycalc
           12 
           13 c_phi = 1.0
           14 c_grad_p = 1.0
           15 #sigma0_list = numpy.array([5.0e3, 10.0e3, 20.0e3, 40.0e3, 80.0e3, 160.0e3])
           16 sigma0_list = numpy.array([10.0e3, 20.0e3, 40.0e3, 80.0e3, 160.0e3])
           17 #alpha = numpy.empty_like(sigma0_list)
           18 #phi_bar = numpy.empty_like(sigma0_list)
           19 load = numpy.array([])
           20 alpha = numpy.array([])
           21 phi_bar = numpy.array([])
           22 
           23 #dc = diffusivitycalc.DiffusivityCalc()
           24 
           25 i = 0
           26 for sigma0 in sigma0_list:
           27 
           28     sid = 'cons-sigma0=' + str(sigma0) + '-c_phi=' + \
           29                      str(c_phi) + '-c_grad_p=' + str(c_grad_p)
           30     if os.path.isfile('../output/' + sid + '.status.dat'):
           31         sim = sphere.sim(sid, fluid=True)
           32 
           33         #sim.visualize('walls')
           34         sim.plotLoadCurve()
           35         load = numpy.append(load, sigma0)
           36         alpha = numpy.append(alpha, sim.c_v)
           37         phi_bar = numpy.append(phi_bar, sim.phi_bar)
           38         #sim.writeVTKall()
           39 
           40     #else:
           41         #print(sid + ' not found')
           42 
           43     i += 1
           44 
           45 fig, ax1 = plt.subplots()
           46 load /= 1000.0
           47 
           48 #plt.ticklabel_format(style='sci', axis='y', scilimits=(0,0))
           49 ax1.plot(load, alpha, 'o-k')
           50 ax1.set_xlabel('Normal stress $\\sigma_0$ [kPa]')
           51 ax1.set_ylabel('Hydraulic diffusivity $\\alpha$ [m$^2$s$^{-1}$]')
           52 #ax1.ticklabel_format(style='plain', axis='y')
           53 #ax1.grid()
           54 
           55 ax2 = ax1.twinx()
           56 #color = 'black'
           57 #ax2.plot(load, phi_bar, 'o--' + color)
           58 ax2.plot(load, phi_bar, 'o--', color='black')
           59 ax2.set_ylabel('Mean porosity $\\bar{\\phi}$ [-]')
           60 #ax2.set_ylabel('Mean porosity $\\bar{\\phi}$ [-]', color=color)
           61 ax2.get_yaxis().get_major_formatter().set_useOffset(False)
           62 #for tl in ax2.get_yticklabels():
           63     #tl.set_color(color)
           64 
           65 filename = 'diffusivity-sigma0-vs-alpha.pdf'
           66 plt.tight_layout()
           67 plt.savefig(filename)
           68 shutil.copyfile(filename, '/home/adc/articles/own/2-org/' + filename)
           69 #print(os.getcwd() + '/' + filename)
           70 print(filename)