Manzari Dafalias Material: Difference between revisions
No edit summary |
No edit summary |
||
Line 70: | Line 70: | ||
* Valid Element Recorder queries are | * Valid Element Recorder queries are | ||
** '''stress''', '''strain''' | ** '''stress''', '''strain''' | ||
** '''alpha''' (or '''backstressratio''') for <math> | ** '''alpha''' (or '''backstressratio''') for <math>0</math> | ||
** '''fabric''' for <math>\mathbf{z}</math> | ** '''fabric''' for <math>\mathbf{z}</math> | ||
** '''alpha_in''' (or '''alphain''') for <math>\mathbf{\alpha_{in}}</math> | ** '''alpha_in''' (or '''alphain''') for <math>\mathbf{\alpha_{in}}</math> |
Revision as of 01:31, 17 December 2013
- Command_Manual
- Tcl Commands
- Modeling_Commands
- model
- uniaxialMaterial
- ndMaterial
- frictionModel
- section
- geometricTransf
- element
- node
- sp commands
- mp commands
- timeSeries
- pattern
- mass
- block commands
- region
- rayleigh
- Analysis Commands
- Output Commands
- Misc Commands
- DataBase Commands
This command is used to construct a multi-dimensional Manzari-Dafalias(2004) material.
nDmaterial ManzariDafalias $matTag $G0 $nu $e_init $Mc $c $lambda_c $e0 $ksi $P_atm $m $h0 $ch $nb $A0 $nd $z_max $cz $Den <$intScheme $TanType $JacoType $TolF $TolR> |
$matTag | integer tag identifying material |
$G0 | bulk modulus constant |
$nu | poisson ratio |
$e_init | initial void ratio |
$Mc | critical state stress ratio |
$c | ratio of critical state stress ratio in extension and compression |
$lambda_c | critical state line constant |
$e0 | critical void ratio at p = 0 |
$ksi | critical state line constant |
$P_atm | atmospheric pressure |
$m | yield surface constant (radius of yield surface in stress ratio space) |
$h0 | constant parameter |
$ch | constant parameter |
$nb | bounding surface parameter, $nb ≥ 0 |
$A0 | dilatancy parameter |
$nd | dilatancy surface parameter $nd ≥ 0 |
$z_max | fabric-dilatancy tensor parameter |
$cz | fabric-dilatancy tensor parameter |
$Den | mass density of the material |
<$intScheme> | plastic integration scheme, 0: forward explicit, 1: backward implicit, 2: improved backward implicit (default = 2) |
<$TanType> | material modulus matrix, 0: elastic stiffness, 1: continuum elastoplastic stiffness, 2: consistent elastoplastic stiffness (default = 2) |
<$JacoType> | jacobian matrix used for newton iterations, 0: finite difference jacobian, 1: analytical jacobian (default = 1) |
<$TolF> | drift from yield surface tolerance (default = 1.0e-7) |
<$TolR> | newton iterations residuals tolerance (default = 1.0e-7) |
The material formulations for the Manzari-Dafalias object are "ThreeDimensional" and "PlaneStrain"
Code Developed by: Alborz Ghofrani, Pedro Arduino, U Washington
Notes
- Valid Element Recorder queries are
- stress, strain
- alpha (or backstressratio) for <math>0</math>
- fabric for <math>\mathbf{z}</math>
- alpha_in (or alphain) for <math>\mathbf{\alpha_{in}}</math>
e.g. recorder Element -eleRange 1 $numElem -time -file stress.out stress
- Elastic or Elastoplastic response could be enforced by
Elastic: updateMaterialStage -material $matTag -stage 0 Elastoplastic: updateMaterialStage -material $matTag -stage 1
- Integration scheme ($intScheme):
0: Forward Euler integration scheme with no error control on the yield surface drift 1: Backward Euler integration scheme (This is the same as Closest Point Projection Method), Newton-Raphson algorithm is used to solve the implicit equation. May not converge due to sensitivity of Newton-Raphson's method to the initial guess. 2: The same as 'option 1' but in case of failure in Newton-Raphson, it uses different strategies to insure convergence.
- Jacobian ($jacoType):
0: Uses a finite difference jacobian <math> \mathbf{J}(x) = \frac{\mathbf{R}(x+h) - \mathbf{R}(x)}{h} </math> 1: Uses analytical jacobian <math> \mathbf{J}(x) = \frac{\partial{\mathbf{R}}}{\partial{x}} </math>
Theory
- <math> p = \frac{1}{3} \mathrm{tr}(\mathbf{\sigma}) </math>
- <math> \mathbf{s} = \mathrm{dev} (\mathbf{\sigma}) = \mathbf{\sigma} - \frac{1}{3} p \mathbf{1}
</math>
Elasticity
Elastic moduli are considered to be functions of p and current void ratio:
- <math> G = G_0 p_{atm}\frac{\left(2.97-e\right)^2}{1+e}\left(\frac{p}{p_{atm}}\right)^{1/2}
</math>
- <math> K = \frac{2(1+\nu)}{3(1-2\nu)} G </math>
The elastic stress-strain relationship is:
- <math> d\mathbf{e}^\mathrm{e} = \frac{d\mathbf{s}}{2G}
</math>
- <math> d\varepsilon^\mathrm{e}_v = \frac{dp}{K}
</math>
Critical State Line
A power relationship is assumed for the critical state line:
- <math> e_c = e_0 - \lambda_c\left(\frac{p_c}{p_{atm}}\right)^\xi </math>
where <math>e_0</math> is the void ratio at <math> p_c = 0 </math> and <math> \lambda_c </math> and <math> \xi </math> constants.
Yield Surface
Yield surface is a stress-ratio dependent surface in this model and is defined as
- <math> \left\| \mathbf{s} - p \mathbf{\alpha} \right\| - \sqrt\frac{2}{3}pm = 0 </math>
with <math> \mathbf{\alpha} </math> being the deviatoric back stress-ratio.
Plastic Strain Increment
The increment of the plastic strain tensor is given by
- <math> d\mathbf{\varepsilon}^p = \langle L \rangle \mathbf{R}
</math>
where
- <math> \mathbf{R} = \mathbf{R'} + \frac{1}{3} D \mathbf{1} </math>
therefore
- <math> d\mathbf{e}^p = \langle L \rangle \mathbf{R'}</math> and <math> d\varepsilon^p_v = \langle L \rangle D </math>
The hardening modulus in this model is defined as
- <math> K_p = \frac{2}{3} p h (\mathbf{\alpha}^b_{\theta} - \mathbf{\alpha}): \mathbf{n}
</math> where <math>\mathbf{n}</math> is the deviatoric part of the gradient to yield surface.
- <math> \mathbf{\alpha}^b_{\theta} = \sqrt{\frac{2}{3}} \left[g(\theta,c) M_c exp(-n^b\Psi) - m\right] \mathbf{n}
</math>, <math>\Psi</math> being the state parameter.
the hardening parameter <math> h </math> is defined as
- <math> h = \frac{b_0}{(\mathbf{\alpha}-\mathbf{\alpha_{in}}):\mathbf{n}}
</math>, <math>\mathbf{\alpha_{in}}</math> is the value of <math>\mathbf{\alpha}</math> at initiation of loading cycle.
- <math>b_0 = G_0 h_0 (1-c_h e) \left(\frac{p}{p_{atm}}\right)^{-1/2} </math>
Also the dilation parameters are defined as
- <math> D = A_d (\mathbf{\alpha}^d_{\theta}-\mathbf{\alpha}) : \mathbf{n}
</math>
- <math> \mathbf{\alpha}^d_{\theta} = \sqrt{\frac{2}{3}} \left[g(\theta,c) M_c exp(n^d\Psi) - m\right] \mathbf{n}
</math>
- <math> A_d = A_0 (1+\langle \mathbf{z : n}\rangle) </math>, where <math> \mathbf{z} </math> is the fabric tensor.
The evolution of fabric and the back stress-ratio tensors are defined as
- <math> d\mathbf{z} = - c_z \langle -d\varepsilon^p_v \rangle (z_{max}\mathbf{n}+\mathbf{z}) </math>
- <math> d\mathbf{\alpha} = \langle L \rangle (2/3) h (\mathbf{\alpha}^b_{\theta} - \mathbf{\alpha})
</math>
Example
This example, provides an undrained confined triaxial compression test using one 8-node SSPBrichUP element and ManzariDafalias material model.
# HHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHH #
# 3D Undrained Conventional Triaxial Compression Test Using One Element #
# University of Washington, Department of Civil and Environmental Eng #
# Geotechnical Eng Group, A. Ghofrani, P. Arduino - Dec 2013 #
# Basic units are m, Ton(metric), s #
# HHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHH #
wipe
# ------------------------ #
# Test Specific parameters #
# ------------------------ #
# Confinement Stress
set pConf -300.0
# Deviatoric strain
set devDisp -0.3
# Permeablity
set perm 1.0e-10
# Initial void ratio
set vR 0.8
# Rayleigh damping parameter
set damp 0.1
set omega1 0.0157
set omega2 64.123
set a1 [expr 2.0*$damp/($omega1+$omega2)]
set a0 [expr $a1*$omega1*$omega2]
# HHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHH
# HHHHHHHHHHHHHHHHHHHHHHHHHHHCreate ModelHHHHHHHHHHHHHHHHHHHHHHHHHHHHHH
# HHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHHH
# Create a 3D model with 4 Degrees of Freedom
model BasicBuilder -ndm 3 -ndf 4
# Create nodes
node 1 1.0 0.0 0.0
node 2 1.0 1.0 0.0
node 3 0.0 1.0 0.0
node 4 0.0 0.0 0.0
node 5 1.0 0.0 1.0
node 6 1.0 1.0 1.0
node 7 0.0 1.0 1.0
node 8 0.0 0.0 1.0
# Create Fixities
fix 1 0 1 1 1
fix 2 0 0 1 1
fix 3 1 0 1 1
fix 4 1 1 1 1
fix 5 0 1 0 1
fix 6 0 0 0 1
fix 7 1 0 0 1
fix 8 1 1 0 1
# Create material
# ManzariDafalias tag G0 nu e_init Mc c lambda_c e0 ksi P_atm m h0 ch nb A0 nd z_max cz Den intScheme TanType JacoType TolF TolR
nDMaterial ManzariDafalias 1 125 0.05 $vR 1.25 0.712 0.019 0.934 0.7 100 0.01 7.05 0.968 1.1 0.704 3.5 4 600 1.42 2 2 1 1.0e-10 1.0e-10
# Create element
# SSPbrickUP tag i j k l m n p q matTag fBulk fDen k1 k2 k3 void alpha <b1 b2 b3>
element SSPbrickUP 1 1 2 3 4 5 6 7 8 1 2.2e6 1.0 $perm $perm $perm $vR 1.5e-9
# Create recorders
recorder Node -file disp.out -time -nodeRange 1 8 -dof 1 2 3 disp
recorder Node -file press.out -time -nodeRange 1 8 -dof 4 vel
recorder Element -file stress.out -time stress
recorder Element -file strain.out -time strain
recorder Element -file alpha.out -time alpha
recorder Element -file fabric.out -time fabric
# Create analysis
constraints Penalty 1.0e18 1.0e18
test NormDispIncr 1.0e-5 20 1
algorithm Newton
numberer RCM
system BandGeneral
integrator Newmark 0.5 0.25
rayleigh $a0 0. $a1 0.0
analysis Transient
# Apply confinement pressure
set pNode [expr $pConf / 4.0]
pattern Plain 1 {Series -time {0 10000 1e10} -values {0 1 1} -factor 1} {
load 1 $pNode 0.0 0.0 0.0
load 2 $pNode $pNode 0.0 0.0
load 3 0.0 $pNode 0.0 0.0
load 4 0.0 0.0 0.0 0.0
load 5 $pNode 0.0 $pNode 0.0
load 6 $pNode $pNode $pNode 0.0
load 7 0.0 $pNode $pNode 0.0
load 8 0.0 0.0 $pNode 0.0
}
analyze 100 100
# Let the model rest and waves damp out
analyze 50 100
# Close drainage valves
for {set x 1} {$x<9} {incr x} {
remove sp $x 4
}
analyze 50 100
# Read vertical displacement of top plane
set vertDisp [nodeDisp 5 3]
# Apply deviatoric strain
set lValues [list 1 [expr 1+$devDisp/$vertDisp] [expr 1+$devDisp/$vertDisp]]
set ts "{Series -time {20000 1020000 10020000} -values {$lValues} -factor 1}"
# loading object deviator stress
eval "pattern Plain 2 $ts {
sp 5 3 $vertDisp
sp 6 3 $vertDisp
sp 7 3 $vertDisp
sp 8 3 $vertDisp
}"
# Set number and length of (pseudo)time steps
set dT 100
set numStep 10000
# Analyze and use substepping if needed
set remStep $numStep
set success 0
proc subStepAnalyze {dT subStep} {
if {$subStep > 10} {
return -10
}
for {set i 1} {$i < 3} {incr i} {
puts "Try dT = $dT"
set success [analyze 1 $dT]
if {$success != 0} {
set success [subStepAnalyze [expr $dT/2.0] [expr $subStep+1]]
if {$success == -10} {
puts "Did not converge."
return success
}
} else {
if {$i==1} {
puts "Substep $subStep : Left side converged with dT = $dT"
} else {
puts "Substep $subStep : Right side converged with dT = $dT"
}
}
}
return success
}
puts "Start analysis"
set startT [clock seconds]
while {$success != -10} {
set subStep 0
set success [analyze $remStep $dT]
if {$success == 0} {
puts "Analysis Finished"
break
} else {
set curTime [getTime]
puts "Analysis failed at $curTime . Try substepping."
set success [subStepAnalyze [expr $dT/2.0] [incr subStep]]
set curStep [expr int(($curTime-20000)/$dT + 1)]
set remStep [expr int($numStep-$curStep)]
puts "Current step: $curStep , Remaining steps: $remStep"
}
}
set endT [clock seconds]
puts "loading analysis execution time: [expr $endT-$startT] seconds."
wipe
References
Dafalias YF, Manzari MT. "Simple plasticity sand model accounting for fabric change effects". Journal of Engineering Mechanics 2004