speccon example code: speccon1d_vrc_stone_column_luetal2010.pyΒΆ
# speccon1d_vrc example (if viewing this in docs, plots are at bottom of page)
# Stone column with vertical and radial drainage in both the column and soil
# Comapre with Lu et al. (2010)
# The orignal solution of Lu et al. (2010)
# is implemented separately in
# geotecha.consolidation.luetal2010.
# Lu, Meng-Meng, Kang-He Xie, and Biao Guo. 2010. 'Consolidation
# Theory for a Composite Foundation Considering Radial and Vertical
# Flows within the Column and the Variation of Soil Permeability
# within the Disturbed Soil Zone'. Canadian Geotechnical
# Journal 47 (2): 207-17. doi:10.1139/T09-086.
# This file should be run with python. It will not work if run with the
# speccon1d_vr.exe script program.
from __future__ import division, print_function
import numpy as np
from geotecha.speccon.speccon1d_vrc import Speccon1dVRC
import matplotlib.pyplot as plt
#Expected values
#t = time values
#avp = average excess pore pressure
#z = depth values
#por = excess pore pressure at time t and depth z.
#settle = settlement
z = np.array(
[ 0. , 0.52631579, 1.05263158, 1.57894737,
2.10526316, 2.63157895, 3.15789474, 3.68421053,
4.21052632, 4.73684211, 5.26315789, 5.78947368,
6.31578947, 6.84210526, 7.36842105, 7.89473684,
8.42105263, 8.94736842, 9.47368421, 10. ])
t = np.array(
[ 0.01 , 0.02682696, 0.07196857, 0.19306977,
0.51794747, 1.38949549, 3.72759372, 10. ])
por = np.array(
[[ 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
0.00000000e+00, 0.00000000e+00],
[ 9.17485805e-01, 7.01719236e-01, 4.21755393e-01,
1.94786515e-01, 6.54056592e-02, 1.33355230e-02,
3.89712328e-04, 3.16237477e-08],
[ 9.66168183e-01, 8.98486152e-01, 6.80015128e-01,
3.60540868e-01, 1.28467752e-01, 2.65620291e-02,
7.76761900e-04, 6.30314731e-08],
[ 9.69825825e-01, 9.20901559e-01, 7.83513686e-01,
4.82747535e-01, 1.87257571e-01, 3.95727628e-02,
1.15850419e-03, 9.40086297e-08],
[ 9.73001762e-01, 9.29271016e-01, 8.19657575e-01,
5.64302374e-01, 2.40538758e-01, 5.22652992e-02,
1.53233099e-03, 1.24343612e-07],
[ 9.75806088e-01, 9.36464754e-01, 8.38815166e-01,
6.17506711e-01, 2.87833953e-01, 6.45432601e-02,
1.89568828e-03, 1.53829201e-07],
[ 9.78277402e-01, 9.42814343e-01, 8.54129319e-01,
6.54605021e-01, 3.29297198e-01, 7.63175703e-02,
2.24609368e-03, 1.82263981e-07],
[ 9.80449715e-01, 9.48405034e-01, 8.67548553e-01,
6.83416142e-01, 3.65476307e-01, 8.75072151e-02,
2.58115346e-03, 2.09453713e-07],
[ 9.82352927e-01, 9.53310968e-01, 8.79367998e-01,
7.07654855e-01, 3.97066262e-01, 9.80395238e-02,
2.89857882e-03, 2.35212664e-07],
[ 9.84013230e-01, 9.57597206e-01, 8.89734934e-01,
7.28756552e-01, 4.24725876e-01, 1.07850060e-01,
3.19620163e-03, 2.59364874e-07],
[ 9.85453477e-01, 9.61320595e-01, 8.98773593e-01,
7.47245201e-01, 4.48982558e-01, 1.16882229e-01,
3.47198913e-03, 2.81745359e-07],
[ 9.86693490e-01, 9.64530519e-01, 9.06592203e-01,
7.63358521e-01, 4.70211358e-01, 1.25086718e-01,
3.72405786e-03, 3.02201237e-07],
[ 9.87750336e-01, 9.67269553e-01, 9.13284407e-01,
7.77252517e-01, 4.88657500e-01, 1.32420888e-01,
3.95068649e-03, 3.20592775e-07],
[ 9.88638561e-01, 9.69574028e-01, 9.18930380e-01,
7.89052617e-01, 5.04473596e-01, 1.38848192e-01,
4.15032750e-03, 3.36794338e-07],
[ 9.89370390e-01, 9.71474512e-01, 9.23597781e-01,
7.98864200e-01, 5.17753731e-01, 1.44337675e-01,
4.32161779e-03, 3.50695255e-07],
[ 9.89955896e-01, 9.72996216e-01, 9.27342552e-01,
8.06774875e-01, 5.28557651e-01, 1.48863591e-01,
4.46338787e-03, 3.62200568e-07],
[ 9.90403137e-01, 9.74159324e-01, 9.30209586e-01,
8.12855360e-01, 5.36925201e-01, 1.52405122e-01,
4.57466988e-03, 3.71231683e-07],
[ 9.90718270e-01, 9.74979268e-01, 9.32233259e-01,
8.17160070e-01, 5.42883792e-01, 1.54946203e-01,
4.65470412e-03, 3.77726910e-07],
[ 9.90905631e-01, 9.75466923e-01, 9.33437836e-01,
8.19727561e-01, 5.46451805e-01, 1.56475410e-01,
4.70294426e-03, 3.81641879e-07],
[ 9.90967799e-01, 9.75628759e-01, 9.33837762e-01,
8.20580842e-01, 5.47639964e-01, 1.56985904e-01,
4.71906097e-03, 3.82949847e-07]])
pors = np.array(
[[ 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
0.00000000e+00, 0.00000000e+00],
[ 9.27633737e-01, 7.09244763e-01, 4.25937878e-01,
1.96408780e-01, 6.57859407e-02, 1.33918199e-02,
3.91305874e-04, 3.17530525e-08],
[ 9.75739340e-01, 9.07326067e-01, 6.86393712e-01,
3.63445060e-01, 1.29204350e-01, 2.66740512e-02,
7.79938099e-04, 6.32891995e-08],
[ 9.78392063e-01, 9.28999691e-01, 7.90242233e-01,
4.86429783e-01, 1.88307407e-01, 3.97393868e-02,
1.16324133e-03, 9.43930170e-08],
[ 9.80677548e-01, 9.36565657e-01, 8.25988279e-01,
5.68305980e-01, 2.41847694e-01, 5.24848888e-02,
1.53859670e-03, 1.24852035e-07],
[ 9.82694742e-01, 9.43043805e-01, 8.44624631e-01,
6.21533374e-01, 2.89345525e-01, 6.48137183e-02,
1.90343973e-03, 1.54458186e-07],
[ 9.84471652e-01, 9.48759787e-01, 8.59451193e-01,
6.58509924e-01, 3.30959720e-01, 7.66364016e-02,
2.25527791e-03, 1.83009232e-07],
[ 9.86032950e-01, 9.53791006e-01, 8.72433990e-01,
6.87150555e-01, 3.67247069e-01, 8.78715909e-02,
2.59170769e-03, 2.10310138e-07],
[ 9.87400313e-01, 9.58204661e-01, 8.83865707e-01,
7.11214114e-01, 3.98912768e-01, 9.84463472e-02,
2.91043094e-03, 2.36174414e-07],
[ 9.88592726e-01, 9.62059704e-01, 8.93890009e-01,
7.32151484e-01, 4.26624937e-01, 1.08296025e-01,
3.20927064e-03, 2.60425378e-07],
[ 9.89626745e-01, 9.65407626e-01, 9.02627867e-01,
7.50490874e-01, 4.50918321e-01, 1.17363873e-01,
3.48618575e-03, 2.82897374e-07],
[ 9.90516726e-01, 9.68293148e-01, 9.10184627e-01,
7.66470876e-01, 4.72173178e-01, 1.25600462e-01,
3.73928509e-03, 3.03436893e-07],
[ 9.91275027e-01, 9.70754815e-01, 9.16651423e-01,
7.80247573e-01, 4.90638166e-01, 1.32963073e-01,
3.96684030e-03, 3.21903631e-07],
[ 9.91912175e-01, 9.72825511e-01, 9.22106262e-01,
7.91946311e-01, 5.06468070e-01, 1.39415101e-01,
4.16729755e-03, 3.38171440e-07],
[ 9.92437017e-01, 9.74532898e-01, 9.26614960e-01,
8.01672362e-01, 5.19758346e-01, 1.44925551e-01,
4.33928815e-03, 3.52129196e-07],
[ 9.92856841e-01, 9.75899786e-01, 9.30231927e-01,
8.09513227e-01, 5.30569632e-01, 1.49468650e-01,
4.48163784e-03, 3.63681552e-07],
[ 9.93177475e-01, 9.76944435e-01, 9.33000820e-01,
8.15539526e-01, 5.38942373e-01, 1.53023563e-01,
4.59337481e-03, 3.72749594e-07],
[ 9.93403372e-01, 9.77680801e-01, 9.34955063e-01,
8.19805595e-01, 5.44904396e-01, 1.55574210e-01,
4.67373627e-03, 3.79271379e-07],
[ 9.93537667e-01, 9.78118721e-01, 9.36118251e-01,
8.22349926e-01, 5.48474361e-01, 1.57109160e-01,
4.72217362e-03, 3.83202356e-07],
[ 9.93582226e-01, 9.78264048e-01, 9.36504425e-01,
8.23195491e-01, 5.49663152e-01, 1.57621569e-01,
4.73835623e-03, 3.84515672e-07]])
porc = np.array(
[[ 0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
0.00000000e+00, 0.00000000e+00, 0.00000000e+00,
0.00000000e+00, 0.00000000e+00],
[ 1.05651237e-01, 9.96770796e-02, 8.71566174e-02,
6.50053734e-02, 3.49831418e-02, 8.83176806e-03,
2.62228672e-04, 2.12793636e-08],
[ 2.00475571e-01, 1.91293001e-01, 1.69728466e-01,
1.28205500e-01, 6.95399456e-02, 1.76002592e-02,
5.22665940e-04, 4.24133677e-08],
[ 2.84526766e-01, 2.73050992e-01, 2.45229852e-01,
1.88167662e-01, 1.03270757e-01, 2.62428385e-02,
7.79532644e-04, 6.32576454e-08],
[ 3.58938854e-01, 3.45699726e-01, 3.13201228e-01,
2.44013857e-01, 1.35823923e-01, 3.46981320e-02,
1.03107403e-03, 8.36698092e-08],
[ 4.24713829e-01, 4.10140609e-01, 3.74057906e-01,
2.95373688e-01, 1.66908194e-01, 4.29066031e-02,
1.27557173e-03, 1.03510423e-07],
[ 4.82737407e-01, 4.67178819e-01, 4.28379409e-01,
3.42212830e-01, 1.96295431e-01, 5.08110698e-02,
1.51135554e-03, 1.22643956e-07],
[ 5.33790958e-01, 5.17527283e-01, 4.76713556e-01,
3.84663146e-01, 2.23815396e-01, 5.83571532e-02,
1.73681479e-03, 1.40939707e-07],
[ 5.78562021e-01, 5.61815510e-01, 5.19551308e-01,
4.22914138e-01, 2.49345817e-01, 6.54936514e-02,
1.95040936e-03, 1.58272696e-07],
[ 6.17653556e-01, 6.00597404e-01, 5.57328872e-01,
4.57162026e-01, 2.72801017e-01, 7.21728422e-02,
2.15068022e-03, 1.74524523e-07],
[ 6.51592052e-01, 6.34358143e-01, 5.90431721e-01,
4.87591411e-01, 2.94121507e-01, 7.83507177e-02,
2.33625938e-03, 1.89584170e-07],
[ 6.80834613e-01, 6.63520194e-01, 6.19198268e-01,
5.14370101e-01, 3.13265738e-01, 8.39871626e-02,
2.50587924e-03, 2.03348764e-07],
[ 7.05775110e-01, 6.88448572e-01, 6.43923104e-01,
5.37648014e-01, 3.30204229e-01, 8.90460831e-02,
2.65838120e-03, 2.15724281e-07],
[ 7.26749478e-01, 7.09455365e-01, 6.64859820e-01,
5.57557059e-01, 3.44915673e-01, 9.34955001e-02,
2.79272365e-03, 2.26626182e-07],
[ 7.44040239e-01, 7.26803618e-01, 6.82223449e-01,
5.74211227e-01, 3.57384471e-01, 9.73076133e-02,
2.90798900e-03, 2.35979996e-07],
[ 7.57880303e-01, 7.40710602e-01, 6.96192539e-01,
5.87706707e-01, 3.67599173e-01, 1.00458846e-01,
3.00339000e-03, 2.43721828e-07],
[ 7.68456103e-01, 7.51350513e-01, 7.06910875e-01,
5.98122006e-01, 3.75551464e-01, 1.02929876e-01,
3.07827505e-03, 2.49798792e-07],
[ 7.75910102e-01, 7.58856641e-01, 7.14488869e-01,
6.05518064e-01, 3.81235486e-01, 1.04705656e-01,
3.13213271e-03, 2.54169378e-07],
[ 7.80342711e-01, 7.63323027e-01, 7.19004632e-01,
6.09938335e-01, 3.84647364e-01, 1.05775430e-01,
3.16459515e-03, 2.56803729e-07],
[ 7.81813642e-01, 7.64805630e-01, 7.20504729e-01,
6.11408866e-01, 3.85784895e-01, 1.06132739e-01,
3.17544065e-03, 2.57683850e-07]])
avp = np.array(
[[ 9.61190784e-01, 9.20180000e-01, 8.32864585e-01,
6.60717386e-01, 3.83698642e-01, 1.00835010e-01,
3.00428723e-03, 2.43793444e-07]])
settle = np.array(
[[ 0.01122695, 0.02309079, 0.04834989, 0.09814961, 0.17828718,
0.26011559, 0.28841662, 0.28928564]])
########################################################
##luetal2010 input to generate expected data
# import numpy as np
# from geotecha.consolidation.luetal2010 import luetal2010
# rc=0.1
# re=0.9
# H=10
# rs=0.15
# ks=0.5
# kv=1.0
# kvc=2000
# kh=1.0
# khc=2
# mvs=0.1
# mvc=0.0005
# gamw=10
# utop=1
# ubot=1
# nterms=100
# z = np.linspace(0,10,20)# np.array([0, 3.0, 6.0, 9.0, 10.0])
# t = np.logspace(np.log10(0.01), np.log10(10), 8)#np.array([0.1, 0.8, 2.0])
#
# por, pors, porc, avp, settle = luetal2010(
# z,t, rc, re, H, rs, ks,
# kv, kvc, kh, khc,
# mvs, mvc, gamw, utop, ubot, nterms)
# #Note for avp and settle you'll need to add an extra [] around the existing [].
# for v in ['z', 't', 'por', 'pors', 'porc', 'avp', 'settle']:
# print('{} = np.{}\n'.format(v, repr(eval(v))).replace("array([","array(\n ["))
#########################################################
reader = ("""\
H = 10
drn = 1
#re=0.9, rc=0.1, rs=0.15, kh/ks=2, n=9, s=1.5, kap=2.0
#mu=1.87284145852, eta = 2/re**2/mu=1.3183901879341955
mvs = 0.1
mvc = 0.0005
mvref = mvs*9**2/(9**2 - 1 + mvs/mvc)
mv = PolyLine([0,1], [1,1])
khref = 1.0
etref = 1.3183901879341955
n=9
dTh = khref * etref/mvref/10 #khref/mvref/gamw * etref
kh = khc=PolyLine([0,1], [1,1])
et = PolyLine([0,1], [1,1])
khcref = 2.0
dThc = 8*khcref / (mvref * 10) /0.1**2
khc = PolyLine([0,1], [1,1])
kvref = 1.0
dTv = kvref/mvref/10/ H**2
kv = PolyLine([0,1], [1,1])
kvcref = 2000
kvc = PolyLine([0, 1], [1, 1])
dTvc=kvcref/H**2/mvref/10
neig = 40
surcharge_vs_depth = [PolyLine([0,1], [1,1])]
surcharge_vs_time = [PolyLine([0,0,10], [0,1,1])]
ppress_z = np.%s
tvals = np.%s
avg_ppress_z_pairs = [[0,1]]
settlement_z_pairs = [[0,1]]
""" % (repr(z/10), repr(t)))
a = Speccon1dVRC(reader)
a.make_all()
# custom plots
title = ("Lu et al (2010) Stone column, load is uniform with depth")
fig = plt.figure(figsize=(13,6))
fig.suptitle(title)
#z vs u
ax1 = fig.add_subplot("141")
ax1.set_xlabel('Excess pore pressure \nin soil, kPa')
ax1.set_ylabel('Depth')
ax1.invert_yaxis()
ax1.plot(pors, z,
ls=".", color='Blue', marker="+", ms=5,
label='expected')
ax1.plot(a.pors, z,
ls='-', color='red', marker='o', ms=5, markerfacecolor='None',
markeredgecolor='red',
label='calculated')
# avp vs t
ax2 = fig.add_subplot("143")
ax2.set_xlabel('Time')
ax2.set_ylabel('Average excess pore press \n(soil and column), kPa')
ax2.set_xscale('log')
ax2.set_xlim((0.01, 10))
ax2.plot(t, avp[0],
ls=".", color='Blue', marker="+", ms=5,
label='expected')
ax2.plot(t, a.avp[0],
ls='-', color='red', marker='o', ms=5, markerfacecolor='None',
markeredgecolor='red',
label='calculated')
# settlement vs t
ax3 = fig.add_subplot("144")
ax3.set_xlabel('Time')
ax3.set_ylabel('Settlement')
ax3.invert_yaxis()
ax3.set_xscale('log')
ax3.set_xlim((0.01, 10))
ax3.plot(t, settle[0],
ls=".", color='Blue', marker="+", ms=5,
label='expected')
ax3.plot(t, a.set[0],
ls='-', color='red', marker='o', ms=5, markerfacecolor='None',
markeredgecolor='red',
label='calculated')
leg = ax3.legend()
leg.draggable()
ax4 = fig.add_subplot("142")
ax4.set_xlabel('Excess pore pressure \nin column, kPa')
ax4.set_ylabel('Depth')
ax4.invert_yaxis()
ax4.plot(porc, z,
ls=".", color='Blue', marker="+", ms=5,
label='expected')
ax4.plot(a.porc, z,
ls='-', color='red', marker='o', ms=5, markerfacecolor='None',
markeredgecolor='red',
label='calculated')
fig.subplots_adjust(top=0.90, bottom=0.15, left=0.1, right=0.94, wspace=0.4)
#fig.tight_layout()
plt.show()
(Source code, png, hires.png, pdf)