14.4.1. Laterally-Loaded Pile Foundation¶
The original model can be found here.
The Python code is converted by Pavan Chigullapally from University of Auckland, Auckland (pchi893@aucklanduni.ac.nz), and shown below, which can be downloaded
here
.
1# -*- coding: utf-8 -*-
2"""
3Created on Thu Jan 10 18:24:47 2019
4
5@author: pchi893
6"""
7
8
9
10
11##########################################################
12# #
13# Procedure to compute ultimate lateral resistance, p_u, #
14# and displacement at 50% of lateral capacity, y50, for #
15# p-y springs representing cohesionless soil. #
16# Converted to openseespy by: Pavan Chigullapally #
17# University of Auckland #
18# #
19# Created by: Hyung-suk Shin #
20# University of Washington #
21# Modified by: Chris McGann #
22# Pedro Arduino #
23# Peter Mackenzie-Helnwein #
24# University of Washington #
25# #
26###########################################################
27
28# references
29# American Petroleum Institute (API) (1987). Recommended Practice for Planning, Designing and
30# Constructing Fixed Offshore Platforms. API Recommended Practice 2A(RP-2A), Washington D.C,
31# 17th edition.
32#
33# Brinch Hansen, J. (1961). "The ultimate resistance of rigid piles against transversal forces."
34# Bulletin No. 12, Geoteknisk Institute, Copenhagen, 59.
35#
36# Boulanger, R. W., Kutter, B. L., Brandenberg, S. J., Singh, P., and Chang, D. (2003). Pile
37# Foundations in liquefied and laterally spreading ground during earthquakes: Centrifuge experiments
38# and analyses. Center for Geotechnical Modeling, University of California at Davis, Davis, CA.
39# Rep. UCD/CGM-03/01.
40#
41# Reese, L.C. and Van Impe, W.F. (2001), Single Piles and Pile Groups Under Lateral Loading.
42# A.A. Balkema, Rotterdam, Netherlands.
43
44import math
45
46def get_pyParam ( pyDepth, gamma, phiDegree, b, pEleLength, puSwitch, kSwitch, gwtSwitch):
47
48 #----------------------------------------------------------
49 # define ultimate lateral resistance, pult
50 #----------------------------------------------------------
51
52 # pult is defined per API recommendations (Reese and Van Impe, 2001 or API, 1987) for puSwitch = 1
53 # OR per the method of Brinch Hansen (1961) for puSwitch = 2
54
55 pi = 3.14159265358979
56 phi = phiDegree * (pi/180)
57 zbRatio = pyDepth / b
58
59 #-------API recommended method-------
60
61 if puSwitch == 1:
62
63 # obtain loading-type coefficient A for given depth-to-diameter ratio zb
64 # ---> values are obtained from a figure and are therefore approximate
65 zb = []
66 dataNum = 41
67 for i in range(dataNum):
68 b1 = i * 0.125
69 zb.append(b1)
70 As = [2.8460, 2.7105, 2.6242, 2.5257, 2.4271, 2.3409, 2.2546, 2.1437, 2.0575, 1.9589, 1.8973, 1.8111, 1.7372, 1.6632, 1.5893, 1.5277, 1.4415, 1.3799, 1.3368, 1.2690, 1.2074, 1.1581,
71 1.1211, 1.0780, 1.0349, 1.0164, 0.9979, 0.9733, 0.9610, 0.9487, 0.9363, 0.9117, 0.8994, 0.8994, 0.8871, 0.8871, 0.8809, 0.8809, 0.8809, 0.8809, 0.8809]
72
73 # linear interpolation to define A for intermediate values of depth:diameter ratio
74 for i in range(dataNum):
75 if zbRatio >= 5.0:
76 A = 0.88
77 elif zb[i] <= zbRatio and zbRatio <= zb[i+1]:
78 A = (As[i+1] - As[i])/(zb[i+1] - zb[i]) * (zbRatio-zb[i]) + As[i]
79
80 # define common terms
81 alpha = phi / 2
82 beta = pi / 4 + phi / 2
83 K0 = 0.4
84
85 tan_1 = math.tan(pi / 4 - phi / 2)
86 Ka = math.pow(tan_1 , 2)
87
88 # terms for Equation (3.44), Reese and Van Impe (2001)
89 tan_2 = math.tan(phi)
90 tan_3 = math.tan(beta - phi)
91 sin_1 = math.sin(beta)
92 cos_1 = math.cos(alpha)
93 c1 = K0 * tan_2 * sin_1 / (tan_3*cos_1)
94
95 tan_4 = math.tan(beta)
96 tan_5 = math.tan(alpha)
97 c2 = (tan_4/tan_3)*tan_4 * tan_5
98
99 c3 = K0 * tan_4 * (tan_2 * sin_1 - tan_5)
100
101 c4 = tan_4 / tan_3 - Ka
102
103 # terms for Equation (3.45), Reese and Van Impe (2001)
104 pow_1 = math.pow(tan_4,8)
105 pow_2 = math.pow(tan_4,4)
106 c5 = Ka * (pow_1-1)
107 c6 = K0 * tan_2 * pow_2
108
109 # Equation (3.44), Reese and Van Impe (2001)
110 pst = gamma * pyDepth * (pyDepth * (c1 + c2 + c3) + b * c4)
111
112 # Equation (3.45), Reese and Van Impe (2001)
113 psd = b * gamma * pyDepth * (c5 + c6)
114
115 # pult is the lesser of pst and psd. At surface, an arbitrary value is defined
116 if pst <=psd:
117 if pyDepth == 0:
118 pu = 0.01
119
120 else:
121 pu = A * pst
122
123 else:
124 pu = A * psd
125
126 # PySimple1 material formulated with pult as a force, not force/length, multiply by trib. length
127 pult = pu * pEleLength
128
129 #-------Brinch Hansen method-------
130 elif puSwitch == 2:
131 # pressure at ground surface
132 cos_2 = math.cos(phi)
133
134 tan_6 = math.tan(pi/4+phi/2)
135
136 sin_2 = math.sin(phi)
137 sin_3 = math.sin(pi/4 + phi/2)
138
139 exp_1 = math.exp((pi/2+phi)*tan_2)
140 exp_2 = math.exp(-(pi/2-phi) * tan_2)
141
142 Kqo = exp_1 * cos_2 * tan_6 - exp_2 * cos_2 * tan_1
143 Kco = (1/tan_2) * (exp_1 * cos_2 * tan_6 - 1)
144
145 # pressure at great depth
146 exp_3 = math.exp(pi * tan_2)
147 pow_3 = math.pow(tan_2,4)
148 pow_4 = math.pow(tan_6,2)
149 dcinf = 1.58 + 4.09 * (pow_3)
150 Nc = (1/tan_2)*(exp_3)*(pow_4 - 1)
151 Ko = 1 - sin_2
152 Kcinf = Nc * dcinf
153 Kqinf = Kcinf * Ko * tan_2
154
155 # pressure at an arbitrary depth
156 aq = (Kqo/(Kqinf - Kqo))*(Ko*sin_2/sin_3)
157 KqD = (Kqo + Kqinf * aq * zbRatio)/(1 + aq * zbRatio)
158
159 # ultimate lateral resistance
160 if pyDepth == 0:
161 pu = 0.01
162 else:
163 pu = gamma * pyDepth * KqD * b
164
165 # PySimple1 material formulated with pult as a force, not force/length, multiply by trib. length
166 pult = pu * pEleLength
167
168 #----------------------------------------------------------
169 # define displacement at 50% lateral capacity, y50
170 #----------------------------------------------------------
171
172 # values of y50 depend of the coefficent of subgrade reaction, k, which can be defined in several ways.
173 # for gwtSwitch = 1, k reflects soil above the groundwater table
174 # for gwtSwitch = 2, k reflects soil below the groundwater table
175 # a linear variation of k with depth is defined for kSwitch = 1 after API (1987)
176 # a parabolic variation of k with depth is defined for kSwitch = 2 after Boulanger et al. (2003)
177
178 # API (1987) recommended subgrade modulus for given friction angle, values obtained from figure (approximate)
179
180 ph = [28.8, 29.5, 30.0, 31.0, 32.0, 33.0, 34.0, 35.0, 36.0, 37.0, 38.0, 39.0, 40.0]
181
182 # subgrade modulus above the water table
183 if gwtSwitch == 1:
184 k = [10, 23, 45, 61, 80, 100, 120, 140, 160, 182, 215, 250, 275]
185
186 else:
187 k = [10, 20, 33, 42, 50, 60, 70, 85, 95, 107, 122, 141, 155]
188
189 dataNum = 13
190 for i in range(dataNum):
191 if ph[i] <= phiDegree and phiDegree <= ph[i+1]:
192 khat = (k[i+1]-k[i])/(ph[i+1]-ph[i])*(phiDegree - ph[i]) + k[i]
193
194 # change units from (lb/in^3) to (kN/m^3)
195 k_SIunits = khat * 271.45
196
197 # define parabolic distribution of k with depth if desired (i.e. lin_par switch == 2)
198 sigV = pyDepth * gamma
199
200 if sigV == 0:
201 sigV = 0.01
202
203 if kSwitch == 2:
204 # Equation (5-16), Boulanger et al. (2003)
205 cSigma = math.pow(50 / sigV , 0.5)
206 # Equation (5-15), Boulanger et al. (2003)
207 k_SIunits = cSigma * k_SIunits
208
209 # define y50 based on pult and subgrade modulus k
210
211 # based on API (1987) recommendations, p-y curves are described using tanh functions.
212 # tcl does not have the atanh function, so must define this specifically
213
214 # i.e. atanh(x) = 1/2*ln((1+x)/(1-x)), |x| < 1
215
216 # when half of full resistance has been mobilized, p(y50)/pult = 0.5
217 x = 0.5
218 log_1 = math.log((1+x)/(1-x))
219 atanh_value = 0.5 * log_1
220
221 # need to be careful at ground surface (don't want to divide by zero)
222 if pyDepth == 0.0:
223 pyDepth = 0.01
224
225 y50 = 0.5 * (pu/ A)/(k_SIunits * pyDepth) * atanh_value
226 # return pult and y50 parameters
227 outResult = []
228 outResult.append(pult)
229 outResult.append(y50)
230
231 return outResult
232
233#########################################################################################################################################################################
234
235#########################################################################################################################################################################
236
237###########################################################
238# #
239# Procedure to compute ultimate tip resistance, qult, and #
240# displacement at 50% mobilization of qult, z50, for #
241# use in q-z curves for cohesionless soil. #
242# Converted to openseespy by: Pavan Chigullapally #
243# University of Auckland #
244# Created by: Chris McGann #
245# Pedro Arduino #
246# University of Washington #
247# #
248###########################################################
249
250# references
251# Meyerhof G.G. (1976). "Bearing capacity and settlement of pile foundations."
252# J. Geotech. Eng. Div., ASCE, 102(3), 195-228.
253#
254# Vijayvergiya, V.N. (1977). "Load-movement characteristics of piles."
255# Proc., Ports 77 Conf., ASCE, New York.
256#
257# Kulhawy, F.H. ad Mayne, P.W. (1990). Manual on Estimating Soil Properties for
258# Foundation Design. Electrical Power Research Institute. EPRI EL-6800,
259# Project 1493-6 Final Report.
260
261def get_qzParam (phiDegree, b, sigV, G):
262
263 # define required constants; pi, atmospheric pressure (kPa), pa, and coeff. of lat earth pressure, Ko
264 pi = 3.14159265358979
265 pa = 101
266 sin_4 = math.sin(phiDegree * (pi/180))
267 Ko = 1 - sin_4
268
269 # ultimate tip pressure can be computed by qult = Nq*sigV after Meyerhof (1976)
270 # where Nq is a bearing capacity factor, phi is friction angle, and sigV is eff. overburden
271 # stress at the pile tip.
272 phi = phiDegree * (pi/180)
273
274 # rigidity index
275 tan_7 = math.tan(phi)
276 Ir = G/(sigV * tan_7)
277 # bearing capacity factor
278 tan_8 = math.tan(pi/4+phi/2)
279 sin_5 = math.sin(phi)
280 pow_4 = math.pow(tan_8,2)
281 pow_5 = math.pow(Ir,(4*sin_5)/(3*(1+sin_5)))
282 exp_4 = math.exp(pi/2-phi)
283
284 Nq = (1+2*Ko)*(1/(3-sin_5))*exp_4*(pow_4)*(pow_5)
285 # tip resistance
286 qu = Nq * sigV
287 # QzSimple1 material formulated with qult as force, not stress, multiply by area of pile tip
288 pow_6 = math.pow(b, 2)
289 qult = qu * pi*pow_6/4
290
291 # the q-z curve of Vijayvergiya (1977) has the form, q(z) = qult*(z/zc)^(1/3)
292 # where zc is critical tip deflection given as ranging from 3-9% of the
293 # pile diameter at the tip.
294
295 # assume zc is 5% of pile diameter
296 zc = 0.05 * b
297
298 # based on Vijayvergiya (1977) curve, z50 = 0.125*zc
299 z50 = 0.125 * zc
300
301 # return values of qult and z50 for use in q-z material
302 outResult = []
303 outResult.append(qult)
304 outResult.append(z50)
305
306 return outResult
307
308#########################################################################################################################################################################
309
310#########################################################################################################################################################################
311##########################################################
312# #
313# Procedure to compute ultimate resistance, tult, and #
314# displacement at 50% mobilization of tult, z50, for #
315# use in t-z curves for cohesionless soil. #
316# Converted to openseespy by: Pavan Chigullapally #
317# University of Auckland #
318# Created by: Chris McGann #
319# University of Washington #
320# #
321###########################################################
322
323def get_tzParam ( phi, b, sigV, pEleLength):
324
325# references
326# Mosher, R.L. (1984). "Load transfer criteria for numerical analysis of
327# axial loaded piles in sand." U.S. Army Engineering and Waterways
328# Experimental Station, Automatic Data Processing Center, Vicksburg, Miss.
329#
330# Kulhawy, F.H. (1991). "Drilled shaft foundations." Foundation engineering
331# handbook, 2nd Ed., Chap 14, H.-Y. Fang ed., Van Nostrand Reinhold, New York
332
333 pi = 3.14159265358979
334
335 # Compute tult based on tult = Ko*sigV*pi*dia*tan(delta), where
336 # Ko is coeff. of lateral earth pressure at rest,
337 # taken as Ko = 0.4
338 # delta is interface friction between soil and pile,
339 # taken as delta = 0.8*phi to be representative of a
340 # smooth precast concrete pile after Kulhawy (1991)
341
342 delta = 0.8 * phi * pi/180
343
344 # if z = 0 (ground surface) need to specify a small non-zero value of sigV
345
346 if sigV == 0.0:
347 sigV = 0.01
348
349 tan_9 = math.tan(delta)
350 tu = 0.4 * sigV * pi * b * tan_9
351
352 # TzSimple1 material formulated with tult as force, not stress, multiply by tributary length of pile
353 tult = tu * pEleLength
354
355 # Mosher (1984) provides recommended initial tangents based on friction angle
356 # values are in units of psf/in
357 kf = [6000, 10000, 10000, 14000, 14000, 18000]
358 fric = [28, 31, 32, 34, 35, 38]
359
360 dataNum = len(fric)
361
362
363 # determine kf for input value of phi, linear interpolation for intermediate values
364 if phi < fric[0]:
365 k = kf[0]
366 elif phi > fric[5]:
367 k = kf[5]
368 else:
369 for i in range(dataNum):
370 if fric[i] <= phi and phi <= fric[i+1]:
371 k = ((kf[i+1] - kf[i])/(fric[i+1] - fric[i])) * (phi - fric[i]) + kf[i]
372
373
374 # need to convert kf to units of kN/m^3
375 kSIunits = k * 1.885
376
377 # based on a t-z curve of the shape recommended by Mosher (1984), z50 = tult/kf
378 z50 = tult / kSIunits
379
380 # return values of tult and z50 for use in t-z material
381 outResult = []
382 outResult.append(tult)
383 outResult.append(z50)
384
385 return outResult
386
387
388#########################################################################################################################################################################
389
390#########################################################################################################################################################################
391
392###########################################################
393# #
394# Static pushover of a single pile, modeled as a beam on #
395# a nonlinear Winkler foundation. Lateral soil response #
396# is described by p-y springs. Vertical soil response #
397# described by t-z and q-z springs. #
398# Converted to openseespy by: Pavan Chigullapally #
399# University of Auckland #
400# Created by: Chris McGann #
401# HyungSuk Shin #
402# Pedro Arduino #
403# Peter Mackenzie-Helnwein #
404# --University of Washington-- #
405# #
406# ---> Basic units are kN and meters #
407# #
408###########################################################
409
410
411import openseespy.opensees as op
412
413op.wipe()
414
415#########################################################################################################################################################################
416
417#########################################################################################################################################################################
418
419# all the units are in SI units N and mm
420
421#----------------------------------------------------------
422# pile geometry and mesh
423#----------------------------------------------------------
424
425# length of pile head (above ground surface) (m)
426L1 = 1.0
427# length of embedded pile (below ground surface) (m)
428L2 = 20.0
429# pile diameter
430diameter = 1.0
431
432# number of pile elements
433nElePile = 84
434# pile element length
435eleSize = (L1+L2)/nElePile
436
437# number of total pile nodes
438nNodePile = 1 + nElePile
439
440#----------------------------------------------------------
441# create spring nodes
442#----------------------------------------------------------
443# spring nodes created with 3 dim, 3 dof
444op.model('basic', '-ndm', 3, '-ndf', 3)
445
446# counter to determine number of embedded nodes
447count = 0
448
449# create spring nodes
450
451#1 to 85 are spring nodes
452
453pile_nodes = dict()
454
455for i in range(nNodePile):
456 zCoord = eleSize * i
457 if zCoord <= L2:
458 op.node(i+1, 0.0, 0.0, zCoord)
459 op.node(i+101, 0.0, 0.0, zCoord)
460 pile_nodes[i+1] = (0.0, 0.0, zCoord)
461 pile_nodes[i+101] = (0.0, 0.0, zCoord)
462 count = count + 1
463
464print("Finished creating all spring nodes...")
465
466# number of embedded nodes
467nNodeEmbed = count
468
469# spring node fixities
470for i in range(nNodeEmbed):
471 op.fix(i+1, 1, 1, 1)
472 op.fix(i+101, 0, 1, 1)
473
474print("Finished creating all spring node fixities...")
475
476#----------------------------------------------------------
477# soil properties
478#----------------------------------------------------------
479
480# soil unit weight (kN/m^3)
481gamma = 17.0
482# soil internal friction angle (degrees)
483phi = 36.0
484# soil shear modulus at pile tip (kPa)
485Gsoil = 150000.0
486
487# select pult definition method for p-y curves
488# API (default) --> 1
489# Brinch Hansen --> 2
490puSwitch = 1
491
492# variation in coefficent of subgrade reaction with depth for p-y curves
493# API linear variation (default) --> 1
494# modified API parabolic variation --> 2
495kSwitch = 1
496
497# effect of ground water on subgrade reaction modulus for p-y curves
498# above gwt --> 1
499# below gwt --> 2
500gwtSwitch = 1
501
502#----------------------------------------------------------
503# create spring material objects
504#----------------------------------------------------------
505
506# p-y spring material
507
508for i in range(1 , nNodeEmbed+1):
509 # depth of current py node
510 pyDepth = L2 - eleSize * (i-1)
511 # procedure to define pult and y50
512 pyParam = get_pyParam(pyDepth, gamma, phi, diameter, eleSize, puSwitch, kSwitch, gwtSwitch)
513 pult = pyParam [0]
514 y50 = pyParam [1]
515 op.uniaxialMaterial('PySimple1', i, 2, pult, y50, 0.0)
516
517
518# t-z spring material
519for i in range(2, nNodeEmbed+1):
520 # depth of current tz node
521 pyDepth = eleSize * (i-1)
522 # vertical effective stress at current depth
523 sigV = gamma * pyDepth
524 # procedure to define tult and z50
525 tzParam = get_tzParam(phi, diameter, sigV, eleSize)
526 tult = tzParam [0]
527 z50 = tzParam [1]
528 op.uniaxialMaterial('TzSimple1', i+100, 2, tult, z50, 0.0)
529
530
531# q-z spring material
532
533 # vertical effective stress at pile tip, no water table (depth is embedded pile length)
534sigVq = gamma * L2
535 # procedure to define qult and z50
536qzParam = get_qzParam (phi, diameter, sigVq, Gsoil)
537qult = qzParam [0]
538z50q = qzParam [1]
539
540#op.uniaxialMaterial('QzSimple1', 101, 2, qult, z50q) #, 0.0, 0.0
541op.uniaxialMaterial('TzSimple1', 101, 2, qult, z50q, 0.0)
542
543print("Finished creating all p-y, t-z, and z-z spring material objects...")
544
545
546#----------------------------------------------------------
547# create zero-length elements for springs
548#----------------------------------------------------------
549
550# element at the pile tip (has q-z spring)
551op.element('zeroLength', 1001, 1, 101, '-mat', 1, 101, '-dir', 1, 3)
552
553# remaining elements
554for i in range(2, nNodeEmbed+1):
555 op.element('zeroLength', 1000+i, i, 100+i, '-mat', i, 100+i, '-dir', 1, 3)
556
557print("Finished creating all zero-Length elements for springs...")
558
559#----------------------------------------------------------
560# create pile nodes
561#----------------------------------------------------------
562
563# pile nodes created with 3 dimensions, 6 degrees of freedom
564op.model('basic', '-ndm', 3, '-ndf', 6)
565
566# create pile nodes
567for i in range(1, nNodePile+1):
568 zCoord = eleSize * i
569 op.node(i+200, 0.0, 0.0, zCoord)
570
571print("Finished creating all pile nodes...")
572
573# create coordinate-transformation object
574op.geomTransf('Linear', 1, 0.0, -1.0, 0.0)
575
576
577# create fixity at pile head (location of loading)
578op.fix(200+nNodePile, 0, 1, 0, 1, 0, 1)
579
580
581# create fixities for remaining pile nodes
582for i in range(201, 200+nNodePile):
583 op.fix(i, 0, 1, 0, 1, 0, 1)
584
585print("Finished creating all pile node fixities...")
586
587#----------------------------------------------------------
588# define equal dof between pile and spring nodes
589#----------------------------------------------------------
590
591for i in range(1, nNodeEmbed+1):
592 op.equalDOF(200+i, 100+i, 1, 3)
593
594print("Finished creating all equal degrees of freedom...")
595
596#----------------------------------------------------------
597# pile section
598#----------------------------------------------------------
599##########################################################################################################################################################################
600
601#########################################################################################################################################################################
602
603#----------------------------------------------------------
604# create elastic pile section
605#----------------------------------------------------------
606
607secTag = 1
608E = 25000000.0
609A = 0.785
610Iz = 0.049
611Iy = 0.049
612G = 9615385.0
613J = 0.098
614
615matTag = 3000
616op.section('Elastic', 1, E, A, Iz, Iy, G, J)
617
618# elastic torsional material for combined 3D section
619op.uniaxialMaterial('Elastic', 3000, 1e10)
620
621# create combined 3D section
622secTag3D = 3
623op.section('Aggregator', secTag3D, 3000, 'T', '-section', 1)
624
625
626#########################################################################################################################################################################
627
628##########################################################################################################################################################################
629
630# elastic pile section
631#import elasticPileSection
632
633#----------------------------------------------------------
634# create pile elements
635#----------------------------------------------------------
636op.beamIntegration('Legendre', 1, secTag3D, 3) # we are using gauss-Legendre integration as it is the default integration scheme used in opensees tcl (check dispBeamColumn)
637
638for i in range(201, 201+nElePile):
639 op.element('dispBeamColumn', i, i, i+1, 1, 1)
640
641print("Finished creating all pile elements...")
642
643#----------------------------------------------------------
644# create recorders
645#----------------------------------------------------------
646
647# record information at specified increments
648timeStep = 0.5
649
650# record displacements at pile nodes
651op.recorder('Node', '-file', 'pileDisp.out','-time', '-dT', timeStep, '-nodeRange', 201, 200 + nNodePile, '-dof', 1,2,3, 'disp')
652
653# record reaction force in the p-y springs
654op.recorder('Node', '-file', 'reaction.out','-time', '-dT', timeStep, '-nodeRange', 1, nNodePile, '-dof', 1, 'reaction')
655
656# record element forces in pile elements
657op.recorder('Element', '-file', 'pileForce.out','-time', '-dT', timeStep, '-eleRange', 201, 200+nElePile, 'globalForce')
658
659print("Finished creating all recorders...")
660
661#----------------------------------------------------------
662# create the loading
663#----------------------------------------------------------
664
665op.setTime(10.0)
666
667# apply point load at the uppermost pile node in the x-direction
668values = [0.0, 0.0, 1.0, 1.0]
669time = [0.0, 10.0, 20.0, 10000.0]
670
671nodeTag = 200+nNodePile
672loadValues = [3500.0, 0.0, 0.0, 0.0, 0.0, 0.0]
673op.timeSeries('Path', 1, '-values', *values, '-time', *time, '-factor', 1.0)
674
675op.pattern('Plain', 10, 1)
676op.load(nodeTag, *loadValues)
677
678print("Finished creating loading object...")
679
680#----------------------------------------------------------
681# create the analysis
682#----------------------------------------------------------
683op.integrator('LoadControl', 0.05)
684op.numberer('RCM')
685op.system('SparseGeneral')
686op.constraints('Transformation')
687op.test('NormDispIncr', 1e-5, 20, 1)
688op.algorithm('Newton')
689op.analysis('Static')
690
691print("Starting Load Application...")
692op.analyze(201)
693
694print("Load Application finished...")
695#print("Loading Analysis execution time: [expr $endT-$startT] seconds.")
696
697#op.wipe
698
699op.reactions()
700Nodereactions = dict()
701Nodedisplacements = dict()
702for i in range(201,nodeTag+1):
703 Nodereactions[i] = op.nodeReaction(i)
704 Nodedisplacements[i] = op.nodeDisp(i)
705print('Node Reactions are: ', Nodereactions)
706print('Node Displacements are: ', Nodedisplacements)
707
708
709
710
711