shell element accuracy issue

[ArdalanNejat]

Dear fmk,
I'm modeling a simple 1 story wall using every possible element in OpenSees( ShellMITC4,quad,spsbrick,stdbrick,...) and using it correctly like verification example "PlanarShearWall" but results are so different from Sap2000(I couldn't reproduce answers you assign to sap results in "PlanarShearWall.tcl" file) and ANSYS,(mesh size in OpenSees, SAP and ANSYS is same and equal to 1 by 1 )would you please guide me a little it's so urgent, please....
here is the code:
wipe

set nxBay 1
set nyFloor 1
set floorHeight 120.
set bayWidth 120.
set numFloor 1
set numBay 1
set E 3000
set v 0.2
set t 12

model basic -ndm 3 -ndf 6

section ElasticMembranePlateSection 1 $E $v $t 0.0

set Plate shell

set eleArgs "1"

set nx [expr $numBay * $nxBay]; # number of elements along building length
set ny [expr $numFloor * $nyFloor]; # number of elements along building height

set nodeTop [expr ($nx+1)*$ny + 1]


set blockCmd "block2D $nx $ny 1 1 $Plate $eleArgs {
1 0. 0. 0.
2 [expr $bayWidth * $numBay] 0. 0.
3 [expr $bayWidth * $numBay] [expr $floorHeight * $numFloor] 0.
4 0. [expr $floorHeight * $numFloor] 0.
}"

eval $blockCmd

# add some loads
pattern Plain 1 Linear {
load $nodeTop 100 0 0 0 0 0
}

fixY 0.0 1 1 1 1 1 1


#floor constraints
for {set i 1} {$i <= $numFloor} {incr i 1} {
set mNode [expr ($nx+1)*$nyFloor*$i + 1]
for {set j 1} {$j <= $nx} {incr j 1} {
equalDOF $mNode [expr $mNode+$j] 1
}
}

integrator LoadControl 1.0
algorithm Linear
numberer RCM
constraints Plain
#system ProfileSPD
system UmfPack
analysis Static

analyze 1

set disp [nodeDisp $nodeTop 1]
puts $disp
wipe
#print -ele
#print -node


here is top node displacement: OpenSees 0.01109 SAP2000: 0.0145

[fmk]

the difference is simply that the element formulations are different and you need to understand a little finite element theory before you use any application (SAP, OpenSees, ....)

the quad and shell elements of OpenSees do not use incompatable modes and as such are not good for a coarse mesh, and a 1x1 is a very very coase mesh. SAP probably includes incompatible modes in their formulation. as you add more elements OpenSees and SAP will converge to something for this particular problem. In OpenSees there is the enhancedQuad element that will give more similar results (problem is that element does not take a thickness and so you have to play with the material properties til i make the change).

i modified the script to show what i mean, going from a 1x1 mesh to a 10x10 mesh and then a 100x100 mesh with different elements .. as you refine mesh for this particular problem you will converge:

nx: 1 ny: 1 quad: 0.01238095238095237791
nx: 1 ny: 1 enhancedQuad: 0.01499999999999999771
nx: 1 ny: 1 SSPquad: 0.02666666666666666491
nx: 1 ny: 1 shell: 0.01109333333333333005
nx: 10 ny: 10 quad: 0.01863382196876817859
nx: 10 ny: 10 enhancedQuad: 0.01870440890866922889
nx: 10 ny: 10 SSPquad: 0.01879851848985728338
nx: 10 ny: 10 shell: 0.01854763394700455550
nx: 100 ny: 100 quad: 0.01881530069271759864
nx: 100 ny: 100 enhancedQuad: 0.01881632230090994964
nx: 100 ny: 100 SSPquad: 0.01881770055460455296
nx: 100 ny: 100 shell: 0.01881021976067661847

i strongly encourage you to review basic text books that deal with the finite element method before you start using it!


set nxBay 10
set nyFloor 10
set floorHeight 120.
set bayWidth 120.
set numFloor 1
set numBay 1

foreach nxBay {1 10 100} nyFloor {1 10 100} {

set eleTypes {quad enhancedQuad SSPquad shell}

foreach eleType $eleTypes {

set E 3000
set v 0.2
set t 12

if {$eleType == "shell"} {
model basic -ndm 3 -ndf 6
} else {
model basic -ndm 2 -ndf 2
}


if {$eleType == "enhancedQuad" } {
set eleArgs "PlaneStress2D 1"
set E [expr $E*$t]
}

if {$eleType == "quad"} {
set eleArgs "$t PlaneStress2D 1"
}

if {$eleType == "SSPquad"} {
set eleArgs "1 PlaneStress $t"
}

if {$eleType == "shell"} {
set eleArgs "1"
}

section ElasticMembranePlateSection 1 $E $v $t 0.0
nDMaterial ElasticIsotropic 1 $E $v

set nx [expr $numBay * $nxBay]; # number of elements along building len\
gth
set ny [expr $numFloor * $nyFloor]; # number of elements along building\
height

set nodeTop [expr ($nx+1)*$ny + 1]

if {$eleType == "shell"} {
set blockCmd "block2D $nx $ny 1 1 $eleType \"$eleArgs\" {
1 0. 0. 0.
2 [expr $bayWidth * $numBay] 0. 0.
3 [expr $bayWidth * $numBay] [expr $floorHeight * $numFloor] 0.
4 0. [expr $floorHeight * $numFloor] 0.
}"

eval $blockCmd

# add some loads
pattern Plain 1 Linear {
load $nodeTop 100 0 0 0 0 0
}
eval $blockCmd

# add some loads
pattern Plain 1 Linear {
load $nodeTop 100 0 0 0 0 0
}

fixY 0.0 1 1 1 1 1 1
} else {
set blockCmd "block2D $nx $ny 1 1 $eleType \"$eleArgs\" {
1 0. 0.
2 [expr $bayWidth * $numBay] 0.
3 [expr $bayWidth * $numBay] [expr $floorHeight * $numFloo\
r]
4 0. [expr $floorHeight * $numFloor]
}"

eval $blockCmd

# add some loads
pattern Plain 1 Linear {
load $nodeTop 100 0
}

fixY 0.0 1 1
}

#floor constraints
for {set i 1} {$i <= $numFloor} {incr i 1} {
set mNode [expr ($nx+1)*$nyFloor*$i + 1]
for {set j 1} {$j <= $nx} {incr j 1} {
equalDOF $mNode [expr $mNode+$j] 1
}
}
integrator LoadControl 1.0
algorithm Linear
numberer RCM
constraints Plain
system ProfileSPD
#system UmfPack
analysis Static

analyze 1

set disp [nodeDisp $nodeTop 1]
puts "nx: $nxBay ny: $nyFloor $eleType: $disp"
wipe
}
}

[ArdalanNejat]

dear fmk thank you so much for quick and brilliant response of yours this meant a lot for me thanks for your precoius time, and off course i will review my text books about finite element.