tadded function to read a specified time - 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 3f5dbc246a357c9339173ebd483f807b2e27b1eb
(DIR) parent 7e337442d9b68f38a76a8f2479de8e216fd19c9c
(HTM) Author: Anders Damsgaard <anders.damsgaard@geo.au.dk>
Date: Fri, 27 Feb 2015 11:42:32 +0100
added function to read a specified time
Diffstat:
M python/halfshear-darcy-combined.py | 29 +++++++++++++++++------------
M python/sphere.py | 29 +++++++++++++++++++++++++++++
2 files changed, 46 insertions(+), 12 deletions(-)
---
(DIR) diff --git a/python/halfshear-darcy-combined.py b/python/halfshear-darcy-combined.py
t@@ -21,6 +21,9 @@ calculateforcechains = False
legend_alpha=0.5
linewidth=0.5
+t_DEM_to_t_real = 5.787e-5
+
+
###################
#### DATA READ ####
t@@ -128,11 +131,13 @@ horizontalalignment='left'
fontweight='bold'
bbox={'facecolor':'white', 'alpha':1.0, 'pad':3}
+t = t/t_DEM_to_t_real / (60.*60.*24.)
+
fig = plt.figure(figsize=[3.5,8])
## ax1: N, tau, ax2: p_f
ax1 = plt.subplot(5, 1, 1)
-lns0 = ax1.plot(t, sigma_def/1000., '-k', label="$\\sigma_0$",
+lns0 = ax1.plot(t, sigma_def/1000., '-k', label="$N$",
linewidth=linewidth)
#lns1 = ax1.plot(t, sigma_eff/1000., '-k', label="$\\sigma'$")
#lns2 = ax1.plot(t, tau_def/1000., '-r', label="$\\tau$")
t@@ -145,10 +150,10 @@ ax2color = 'blue'
#lns4 = ax2.plot(t, p_f_top/1000.0 + 80.0, '-',
#color=ax2color,
#label='$p_\\text{f}^\\text{forcing}$')
-lns5 = ax2.plot(t, p_f_bar/1000.0 + 80.0, '--',
+lns5 = ax2.plot(t, p_f_bar/1000.0, '--',
color=ax2color,
- label='$\\bar{p}_\\text{f}$', linewidth=linewidth)
-ax2.set_ylabel('Mean fluid pressure [kPa]')
+ label='$\\Delta\\bar{p}_\\text{f}$', linewidth=linewidth)
+ax2.set_ylabel('Mean change in fluid pressure [kPa]')
ax2.yaxis.label.set_color(ax2color)
for tl in ax2.get_yticklabels():
tl.set_color(ax2color)
t@@ -161,7 +166,7 @@ labs = [l.get_label() for l in lns]
ax2.legend(lns, labs, loc='upper right', ncol=3,
fancybox=True, framealpha=legend_alpha)
ax1.set_ylim([-30, 200])
-ax2.set_ylim(ax1.get_ylim())
+#ax2.set_ylim(ax1.get_ylim())
ax1.text(bbox_x, bbox_y, 'a',
horizontalalignment=horizontalalignment,
t@@ -219,7 +224,7 @@ ax7.set_ylim([1.0e1, 2.0e4])
ax8 = ax7.twinx()
ax8color='green'
ax8.plot(t, coordinationnumber, color=ax8color, linewidth=linewidth)
-ax8.set_ylabel('Coordination number [-]')
+ax8.set_ylabel('Contacts per particle [-]')
ax8.yaxis.label.set_color(ax8color)
for tl in ax8.get_yticklabels():
tl.set_color(ax8color)
t@@ -249,19 +254,19 @@ im9 = ax9.pcolormesh(t, zpos_c, poros,
rasterized=True)
ax9.set_ylim([zpos_c[0], sim.w_x[0]])
ax9.set_ylabel('Vertical position [m]')
-#cb9 = plt.colorbar(im9, orientation='horizontal')#, pad=0.20)
+
cbaxes = fig.add_axes([0.32, 0.1, 0.4, 0.01]) # x,y,w,h
-#cb9 = plt.colorbar(im9, orientation='horizontal', shrink=0.7)
+
#ax9.add_patch(matplotlib.patches.Rectangle(
- #(3.0, 0.06), # x,y
+ #(3.0, 0.04), # x,y
#15., # dx
#.15, # dy
#fill=True,
#linewidth=1,
#facecolor='white'))
ax9.add_patch(matplotlib.patches.Rectangle(
- (3.0, 0.04), # x,y
- 15., # dx
+ (0.6, 0.04), # x,y
+ 3., # dx
.15, # dy
fill=True,
linewidth=1,
t@@ -298,7 +303,7 @@ plt.setp(ax5.get_xticklabels(), visible=False)
plt.setp(ax7.get_xticklabels(), visible=False)
#plt.setp(ax8.get_xticklabels(), visible=False)
-ax9.set_xlabel('Time [s]')
+ax9.set_xlabel('Time [d]')
fig.tight_layout()
plt.subplots_adjust(hspace=0.05)
(DIR) diff --git a/python/sphere.py b/python/sphere.py
t@@ -2177,6 +2177,35 @@ class sim:
fn = "../output/{0}.output{1:0=5}.bin".format(self.sid, lastfile)
self.readbin(fn, verbose)
+ def readTime(self, time, verbose=True):
+ '''
+ Read the output file most closely corresponding to the time given as an
+ argument.
+
+ :param time: The desired current time [s]
+ :type time: float
+
+ See also :func:`readbin()`, :func:`readfirst()`, :func:`readsecond`, and
+ :func:`readstep`.
+ '''
+
+ self.readfirst(verbose=False)
+ t_first = self.currentTime()
+ n_first = self.time_step_count[0]
+
+ self.readlast(verbose=False)
+ t_last = self.currentTime()
+ n_last = self.time_step_count[0]
+
+ if time < t_first | time > t_last:
+ raise Exception('Error: The specified time {} s is outside the ' +
+ 'range of output files [{}; {}] s.'.format(time, \
+ t_first, t_last))
+
+ dt_dn = (t_last - t_first)/(n_last - n_first)
+ step = int((time - t_first)/dt_dn) + n_first
+ sim.readstep(step, verbose=verbose)
+
def generateRadii(self, psd = 'logn',
mean = 440e-6,
variance = 8.8e-9,