VELACS model-1 Modelling

A forum dedicated to users with questions regarding soil materials and elements.

forum currently locked

Moderator: Moderators

omidpaul
Posts: 2
Joined: Wed Jul 11, 2012 10:47 pm
Location: Sharif University of Technology

VELACS model-1 Modelling

Hello all,
I am trying to model VELACS 1, it is one of my first attempts to use OpenSees so I have a little bit of problems.
I have modeled the geometry of the problem, and in gravity loading (soil self weight) I have faced a problem and get errors while running the code.
I have distributed the soil weight forces to the nodes.
I hope somebody can help me with it,

here's the code:

# Units: m, sec, ton, kN

wipe all
set startT [clock seconds]

################################### Soil Properties ##########################################

set Gs 2.68 ;# Soil Grain mass density (Gs)
set nou 0.05 ;# Poisson's Ratio
set rhof 1.0 ;# Fluid mass density
set Phi 31.4 ;# Friction angle at peak shear strength, in degree
set Dr 40.0 ;# Relative Density
set e0 [expr 0.894-(0.894-0.516)*\$Dr/100]
set poro0 [expr \$e0/(1+\$e0)]
set rhos [expr ((1.0-\$poro0)*\$Gs+\$poro0)*\$rhof] ;# Saturated Soil mass density

# Permeability (m/s)
set kx 2.1e-5
set ky 2.1e-5
set kz 2.1e-5

set gravX 0.0000
set gravY 0.0000
set gravZ -9.8100

################################### Define Number of Analysis steps ###########################

set NumIncr1 50 ; # Number of Steps in Static Analysis

################################### Define Soil Nodes ##########################################

model basic -ndm 2 -ndf 3

set xincr 2.54
set yincr 2.00

for {set i 1} {\$i<=[expr int(22.86/\$xincr)+1]} {incr i 1} {
for {set j 1} {\$j<=[expr int(10.00/\$yincr)+1]} {incr j 1} {
set xdim [expr (\$i-1)*\$xincr ]
set ydim [expr (\$j-1)*\$yincr ]
set NodeNum [expr \$j+ (int(10.00/\$yincr)+1)*(\$i-1)]
node \$NodeNum \$xdim \$ydim
puts "\$NodeNum \$xdim \$ydim"
}
}

########################## Define Constitutive Model for soil (Dafalias-Manzari) ###############

set Gr 7.5e4
set Br 2.0e5
set PeakShearStrain 0.1
set Pr 80
set pressDependCoe 0.5
set PTAng 27
set contrac 0.07
set dilat1 0.4
set dilat2 2
set liquefac1 10
set liquefac2 0.01
set liquefac3 1

nDMaterial PressureDependMultiYield 1 2 \$rhos \$Gr \$Br \$Phi \$PeakShearStrain \$Pr \$pressDependCoe \$PTAng \$contrac \$dilat1 \$dilat2 \$liquefac1 \$liquefac2 \$liquefac3

################################### Define Soil Elements #######################################

set m 1
for {set i 1} {\$i <= [expr int(22.86/\$xincr)]} {incr i 1} {
for {set j 1} {\$j <= [expr int(10.00/\$yincr)]} {incr j 1} {
set n1 [expr 6*\$i-(6-\$j)]
set n2 [expr 6*\$i-(6-\$j)+1]
set n3 [expr 6*\$i+\$j]
set n4 [expr 6*\$i+\$j+1]
set EleNum [expr (\$m)]
element quadUP \$EleNum \$n1 \$n2 \$n3 \$n4 1.0 1 [expr 2.2e6/\$poro0] \$rhof \$kx \$ky 0 -9.81 0
set m [expr (\$m+1)]
puts "\$EleNum \$n1 \$n2 \$n3 \$n4 \$j"
}
}

### update material stage to zero for elastic behavior

updateMaterialStage -material 1 -stage 0

################################### Boundary Conditions #######################################

### Lateral Boundaries

fixX 0.0 1 0 0 -tol 1e-10
fixY 0.0 1 1 0 -tol 1e-10
fixX 22.86 1 0 0 -tol 1e-10

for {set m 1} {\$m<=5} {incr m 1} {
equalDOF \$m [expr \$m+1] 2
equalDOF [expr \$m+54] [expr \$m+55] 2
}

### Surface Pore Pressure

for {set m 1} {\$m<=10} {incr m 1} {
fix [expr 6*\$m] 0 0 1
}

########################### Applying Gravity Loads to Nodes ###########################

set m 1
for {set i 1} {\$i <= [expr int(22.86/\$xincr)]} {incr i 1} {
for {set j 1} {\$j <= [expr int(10.00/\$yincr)]} {incr j 1} {
set n1 [expr 6*\$i-(6-\$j)]
set n2 [expr 6*\$i-(6-\$j)+1]
set n3 [expr 6*\$i+\$j]
set n4 [expr 6*\$i+\$j+1]
pattern Plain \$m Linear {
load \$n1 0 [expr \$xincr*\$yincr*9.81*(\$rhos-\$rhof)/4] 0
load \$n2 0 [expr \$xincr*\$yincr*9.81*(\$rhos-\$rhof)/4] 0
load \$n3 0 [expr \$xincr*\$yincr*9.81*(\$rhos-\$rhof)/4] 0
load \$n4 0 [expr \$xincr*\$yincr*9.81*(\$rhos-\$rhof)/4] 0
set m [expr (\$m+1)]
}
}
}
################################# Gravity Analysis #################################

system UmfPack
constraints Transformation
test NormDispIncr 0.01 50 1
algorithm Newton
numberer RCM
analysis Static

analyze 10

updateMaterialStage -material 1 -stage 1

analyze 40