Drucker Prager material doesn't consider Flow beyond yield

[tusharkantimandal]

Hi all,
I was trying to model concrete as nDMaterial drucker prager. Doing so, I found that:
(1) Beyond yield it doesn't consider the flow of material. That's why the stress-strain curve NEVER goes down, it is continuously increasing.
and The stress values get magnified when hardening parameters are assigned some non-zero values compared to zero hardening values.
(2) Upto yield stress-strain curve is considered as linearly elastic.

I tried for a single brick element, nodes at one face of the element are fixed and it is pushed at the other end (by uniaxial loading). Displacement controlled analysis is done using monitoring direction as "1"
I was trying assess the uniaxial behaviour of this material model by using poisson ratio =0.0. I actually tried both way: poisson's ratio 'zero'/'Non-zero'.
(1) when Poisson's ratio=0.0: Sigma_2, Sigma_3=0.0, upto yield and some Non-zero value of Sigma_2, Sigma_3 beyond yield and Sigma_1,Sigma_2, Sigma_3 increasing continuously even beyond 1000 MPa, (wheras fck=40MPa, fcr=4.427 MPa used for analysis);
(2) when Poission ratio= Non-zero same phenomenon with Sigma_2, Sigma_3=Non-zero before yield as well.

From these what can we conclude?
Drucker-Prager Material model in OpenSees cann't model the post-yield behaviour of the material, It can identify the stress-condition when the material will reach the yield situation in any general loading condition only??
Any material model how much complex it may be It should show the proper uniaxial behaviour when the loading and boundary condition is set like that. But D-P or any of the nonlinear nDmaterial model like J2Plasticty doesn't fulfill that criteria. So are these examples of very improper material model??

Hi friends, am I wrong what I just written??
waiting for your suggestions..............

#### here is the test problem i am trying with####
wipe

model BasicBuilder -ndm 3 -ndf 3

set ele_type 1;#1=std brick, 0=SSP brick

set mat_type 2;#0=elastic isotropic, 1= Vonmises, 2=D-P

set analysis_control 1;#0=load control analysis, 1=Displacement controlled analysis...

set fck 40

set ft [expr 0.7*pow($fck,0.5)];

set output out_t;

file mkdir $output;

set mon_dof 1;

set mon_node 2;

set dLoad -0.0001;

set P 1;

set maxU -5;#max limiting displacement

set dUmin -0.01; # Displacement increment

set dUmax -0.01;

set dU -0.01



#--create the nodes

set x 100;

set y 1;

set z 1;

node 1 0.0 0.0 0.0;

node 2 $x 0.0 0.0;

node 3 $x $y 0.0;

node 4 0.0 $y 0.0;

node 5 0.0 0.0 $z;

node 6 $x 0.0 $z;

node 7 $x $y $z;

node 8 0.0 $y $z;



#boundary conditions

fix 1 1 1 1

fix 4 1 1 1

fix 8 1 1 1

fix 5 1 1 1

#fix 2 0 1 1

#fix 3 0 1 1

#fix 7 0 1 1

#fix 6 0 1 1

#equaldof....

equalDOF 2 3 1

equalDOF 2 7 1

equalDOF 2 6 1

#nDmaterial....

set matTag 1;

set E [expr 2*$fck/0.002];

set v 0.0;

set G [expr $E/(2*(1+$v))];

set K [expr $E/(3*(1-2*$v))];

set density [expr 2500*pow(10,-9)];

if {$mat_type==0} {

nDMaterial ElasticIsotropic $matTag $E $v $density

}

if {$mat_type==1} {

set sig0 $fck;

set sigInf [expr 1.0*$fck];

set delta 0.0;

set H [expr 0.0];

nDMaterial J2Plasticity $matTag $K $G $sig0 $sigInf $delta $H

}

if {$mat_type==3} {

# Yield surface

set DPys "-DP"

# Potential surface

set DPps "-DP 0.1"



# Scalar evolution law: linear hardening coef=1.0

set ES1 "-Leq 1.10"



# Initial stress

set stressp "0.1 0 0 0 0.1 0 0 0 0.1"



# EPState

#------------E--------Eo--------v----rho-----------------alpha----k

set EPS " 30000.0 30000.0 0.18 2.4 -NOD 0 -Nos 2 0.315 0.0 -stressp $stressp"



# Create nDMaterial using Template Elastic-Plastic Model

nDMaterial Template3Dep $matTag -YS $DPys -PS $DPps -EPS $EPS -ELS1 $ES1

}

if {$mat_type==2} {

#--define material parameters for the model

set mf_E 1.0;

set mf_fy 1.0;

set mf_rho 1.0;

#---bulk modulus

set k $K;

#---shear modulus

set G $G;

#---yield stress

set sigY [expr $fck*$ft/($fck+$ft)*pow(2,0.5)*$mf_fy];

#---failure surface and associativity

set rho [expr ($fck-$ft)/($fck+$ft)*0.8165*$mf_rho];

set rhoBar [expr $rho*0.0];

#---isotropic hardening

set Kinf 0.0

set Ko 0.0

set delta1 0.0

#---kinematic hardening

set H 0.0

set theta 1.0

#---tension softening

set delta2 0.0

#---mass density

set mDen $density;#in kg/m^3

set matTag 1;

#--material models

# type tag k G sigY rho rhoBar Kinf Ko delta1 delta2 H theta density

nDMaterial DruckerPrager $matTag $K $G $sigY $rho $rhoBar $Kinf $Ko $delta1 $delta2 $H $theta $mDen

}

#--create the element

# type tag nodes matID bforce1 bforce2 bforce3

if {$ele_type==1} {

element stdBrick 1 1 2 3 4 5 6 7 8 $matTag 0.0 0.0 0.0

}

if {$ele_type==0} {

element SSPbrick 1 1 2 3 4 5 6 7 8 $matTag 0.0 0.0 0.0

}



puts "model Built..."



#-------------------------------------------------------

# create the recorders

#-------------------------------------------------------

recorder Node -node 2 -time -file $output/node2.out -dof 1 disp;

if {$ele_type==0} {

recorder Element -ele 1 -time -file $output/stress1.out stress

recorder Element -ele 1 -time -file $output/strain1.out strain

}

if {$ele_type==1} {

recorder Element -ele 1 -time -file $output/stress1.out material $matTag stress

recorder Element -ele 1 -time -file $output/strain1.out material $matTag strain

}

puts "recorders set..."



#-------------------------------------------------------

# create the loading

#-------------------------------------------------------



#Ref load...

set P1 [expr $P/4.0];

set patternTag 1;

pattern Plain $patternTag Linear {

# Create nodal loads at node 2

# nd FX FY FZ

load 2 $P1 0.0 0.0

load 3 $P1 0.0 0.0

load 7 $P1 0.0 0.0

load 6 $P1 0.0 0.0

}

#-------------------------------------------------------

# create the analysis

#-------------------------------------------------------

system BandGeneral



# Create the constraint handler, the transformation method

constraints Transformation



# Create the DOF numberer, the reverse Cuthill-McKee algorithm

numberer RCM

if {$analysis_control==1} {

# Change the integration scheme to be displacement control

# node dof init Jd min max

integrator DisplacementControl $mon_node $mon_dof $dU 1 $dUmin $dUmax

}

if {$analysis_control==0} {

integrator LoadControl $dLoad;

}

set Tol 1e-06;

set maxiter 100;#max. no of iterations

set printFlagDynamic 0;

test NormDispIncr $Tol $maxiter $printFlagDynamic;

algorithm ModifiedNewton -initial;

analysis Static;

if {$analysis_control==0} {

set ok 0;

set currentDisp 0;

while {$ok == 0 && [expr {abs($currentDisp)}] < [expr {abs($maxU)}]} {

puts "displacement= $currentDisp";

set ok [analyze 1];

if {$ok!=0} {

set dLoad [expr $dLoad/10.0];

integrator LoadControl $dLoad;

set ok [analyze 1];

source algos.tcl;

}

set currentDisp [nodeDisp $mon_node $mon_dof];

}

}

if {$analysis_control==1} {

set ok 0;

set currentDisp 0.0;

puts "max. displacement= $maxU"

# perform the analysis

while {$ok == 0 && [expr {abs($currentDisp)}] < [expr {abs($maxU)}]} {

set ok [analyze 1]

if {$ok!=0} {

set dU [expr $dU/10.0];

set dUmin [expr $dUmin/10.0];

integrator DisplacementControl $mon_node $mon_dof $dU 1 $dUmin $dUmax;

set ok [analyze 1];

source algos.tcl;

}

set currentDisp [nodeDisp $mon_node $mon_dof];

puts "displacement= $currentDisp"

}

if {$ok == 0} {

puts "Pushover analysis completed SUCCESSFULLY";

} else {

puts "Pushover analysis FAILED";

}

}

puts "$sigY";

puts "$rho";

[gerber]

I also try to use Drucker Prager material to simulate concrete behavior.
I made a uniaxial loading test in a conrete cube and I noticed an elastic totally plastic behavior.
The strength was right ( 20 MPa) but the strain was very small in the end of the elastic area (0,07%).
I am not sure if we can continue our effort by using D-P !

[tusharkantimandal]

gerber wrote:
> I also try to use Drucker Prager material to simulate concrete behavior.
> I made a uniaxial loading test in a conrete cube and I noticed an elastic
> totally plastic behavior.
> The strength was right ( 20 MPa) but the strain was very small in the end
> of the elastic area (0,07%).
> I am not sure if we can continue our effort by using D-P !
Can I get your code?..I need to see where I am doing mistake..it's not flowing at all in my case..I am getting continuously increasing stress value..
yea the strain value is not correct in my case also...we can modify the K,G values to get modified strain value since it's taking linearly elastic in-stead-of parabolic shape...so changing the E-value (actually it should be half of the E-value considered for parabolic shape for same fck-value)..

waiting for your reply....

[saumyashah]

hello everyone
can anyone please explain me the tension softening parameter in drucker prager material what should be it's value in case of masonry structure having tensile stress of 4.6N/mm2
i have got literature in which it is mention in terms of fracture energy can anybody help me in relating this with opensees parameters

[saumyashah]

dear all

please someone can explain me softening parameter in drucker prager material
i have tried different values of delta 2 but there is no change in stress - strain curve i am posting my script please give me suggestions

After maximum tensile stress the stress strain curve remains constant and there is no tension softening
wipe

model BasicBuilder -ndm 3 -ndf 3

#--create the 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

#--uniaxial test boundary conditions
fix 1 0 1 1
fix 2 0 0 1
fix 3 1 0 1
fix 4 1 1 1
fix 5 0 1 0
fix 6 0 0 0
fix 7 1 0 0
fix 8 1 1 0

#--define material parameters for the model
set E 30760;

set v 0.2;

set G [expr $E/(2*(1+$v))];

set K [expr $E/(3*(1-2*$v))];


#---yield stress
set sigY 8.36
#---failure surface and associativity
set rho 0.67
set rhoBar 0.67
#---isotropic hardening
set Kinf 0.0
set Ko 0.0
set delta1 0.0
#---kinematic hardening
set H 0.0
set theta 1.0
#---tension softening
set delta2 3000
#---mass density
set mDen 2.5

#--material models
# type tag k G sigY rho rhoBar Kinf Ko delta1 delta2 H theta density
nDMaterial DruckerPrager 2 $K $G $sigY $rho $rhoBar $Kinf $Ko $delta1 $delta2 $H $theta $mDen

#--create the element
# type tag nodes matID bforce1 bforce2 bforce3
element stdBrick 1 1 2 3 4 5 6 7 8 2 0.0 0.0 0.0

puts "model Built..."

#-------------------------------------------------------
# create the recorders
#-------------------------------------------------------
set step 0.1

# record nodal displacements
recorder Node -file displacementst.out -time -node 5 -dof 1 2 3 disp

# record the element stress, strain, and state at one of the Gauss points
recorder Element -ele 1 -time -file stress1tec.out material 2 stress
recorder Element -ele 1 -time -file strain1tec.out material 2 strain
recorder Element -ele 1 -time -file state1tec.out material 2 state


puts "recorders set..."

#-------------------------------------------------------
# create the loading
#-------------------------------------------------------

#--pressure magnitude
set p 10.0
set pNode [expr $p/4]

#--loading object for hydrostatic pressure
pattern Plain 1 "Linear" {
load 1 $pNode 0.0 0.0
load 2 $pNode 0.0 0.0

load 5 $pNode 0.0 0.0
load 6 $pNode 0.0 0.0

}

#--loading object deviator stress
#pattern Plain 2 {Series -time {0 10 100} -values {0 1 5} -factor 1} {
# load 5 0.0 0.0 $pNode
# load 6 0.0 0.0 $pNode
# load 7 0.0 0.0 $pNode
# load 8 0.0 0.0 $pNode
#}

#-------------------------------------------------------
# create the analysis
#-------------------------------------------------------
numberer RCM
system SparseGeneral
constraints Transformation
test NormDispIncr 1e-5 10 1
algorithm Newton
analysis Static



set dU 0.00001; # Displacement increment

# Change the integration scheme to be displacement control
# node dof init Jd min max
integrator DisplacementControl 5 1 $dU 1 $dU $dU

recorder Node -node 5 -time -file node5s.out -dof 1 disp;

set maxU 0.0050000; # Max displacement
set numSteps [expr int($maxU/$dU)]

# Perform the analysis
set ok [analyze $numSteps]

if {$ok != 0} {

set currentDisp [nodeDisp 5 1]
print currentDisp
set ok 0
while {$ok == 0 && $currentDisp < $maxU} {

set ok [analyze 1]

# if the analysis fails try initial tangent iteration
if {$ok != 0} {
puts "regular newton failed .. lets try an initail stiffness for this step"
test NormDispIncr 1.0e-12 2000 1
algorithm ModifiedNewton -initial
set ok [analyze 1]
if {$ok == 0} {puts "that worked .. back to regular newton"}
test NormDispIncr 1.0e-8 100
algorithm Newton
}

set currentDisp [nodeDisp 5 1]
}
}

puts "";
if {$ok == 0} {
puts "Pushover analysis completed SUCCESSFULLY";
} else {
puts "Pushover analysis FAILED";
}