the following is a library file that I use which defines a number of procs for static analysis. it is used for a pushover or cyclic load reversals. please see it and use it as needed. please let me know if you find anything wrong or questionable with it. (The last proc is extra)
Code: Select all
############################################################################
############################################################################
# LibStatic.tcl: static gravity, lateral pushover or cyclic procedures
# Silvia Mazzoni, 2006
############################################################################
############################################################################
############################################################################
# procApplyGravity.tcl -- apply gravity load, set it constant and reset time to zero.
# Silvia Mazzoni, 2006
############################################################################
proc procApplyGravity { } {
# apply gravity load, set it constant and reset time to zero.
global constraintsTypeGravity;
global alphaSP; # Penalty/Lagrange constraints -- factor used when adding the single-point constraint into the system of equations
global alphaMP; # Penalty/Lagrange constraints -- factor used when adding the multi-point constraint into the system of equations
global numbererType;
global systemTypeGravity;
global testTypeGravity; # convergence-test type
global TolGravity; # Convergence Test: tolerance
global maxNumIterGravity; # Convergence Test: maximum number of iterations that will be performed before "failure to converge" is returned
global printFlagGravity; # Convergence Test: flag used to print information on convergence (optional)
global NstepGravity; # number of steps to apply gravity
global algorithmTypeGravity;
global algorithmCount;
if {$constraintsTypeGravity == "Plain" | $constraintsTypeGravity == "Transformation"} {
constraints $constraintsTypeGravity ; # Plain & Transformation constraints
} else {
constraints $constraintsTypeGravity $alphaSP $alphaMP ; # Penalty & Lagrange congs
}
if {$systemTypeGravity == "SparseGeneralPivot"} {
system SparseGeneral -piv; # optional pivoting for SparseGeneral system
} else {
system $systemTypeGravity ;
}
numberer $numbererType; # renumber dof's to minimize band-width (optimization), if you want to
test $testTypeGravity $TolGravity $maxNumIterGravity ; # tolerance, max no. of iterations, and print code , 1: every iteration
if {$algorithmTypeGravity == "BFGS" | $algorithmTypeGravity == "Broyden" } {
algorithm $algorithmTypeGravity $algorithmCount ; # use Newton's solution algorithm: updates tangent stiffness at every iteration
} elseif {$algorithmTypeGravity == "NewtonLineSearch"} {
algorithm $algorithmTypeGravity $NewtonLineSearchRatio; # use Newton's solution algorithm: updates tangent stiffness at every iteration
} else {
algorithm $algorithmTypeGravity ; # use Newton's solution algorithm: updates tangent stiffness at every iteration
}
set DGravity [expr 1./$NstepGravity]; # first load increment;
integrator LoadControl $DGravity
analysis Static
analyze $NstepGravity
loadConst -time 0.0
puts doneGravity
}; #
###########################################################################
######################################################################################
## procFlateral $iLcol $iFloorWeight
######################################################################################
# calculate distribution of lateral load based on mass/weight distributions along building height
# Fj = WjHj/sum(WiHi) * Weight
# by Silvia Mazzoni, 2006
#
proc procFlateral {iLcol iFloorWeight } {
set Y0 0.
set iYlevel ""
foreach Lcol $iLcol {
set Y0 [expr $Y0+$Lcol];
lappend iYlevel $Y0
}
set denom 0.
set Wtot 0.0
foreach FloorWeight $iFloorWeight Ylevel $iYlevel {
set denom [expr $denom + $FloorWeight*$Ylevel]
set Wtot [expr $Wtot+$FloorWeight]; # calculate total weight
}
set iFj ""
foreach FloorWeight $iFloorWeight Ylevel $iYlevel {
set Fj [expr ($Ylevel*$FloorWeight)/$denom*$Wtot]
lappend iFj $Fj
}
return $iFj
}; #
######################################################################################
######################################################################################
## procGenPeaks $Dmax $CycleType $Fact
######################################################################################
# generate incremental disps for Dmax
# by Silvia Mazzoni, 2006
#
proc procGenPeaks {Dmax {DincrStatic 0.01} {CycleType "Full"} {Fact 1} } {; # generate incremental disps for Dmax
file mkdir data
set outFileID [open data/tmpDsteps.tcl w]
set Disp 0.
puts $outFileID "set iDstep { "
puts $outFileID $Disp
puts $outFileID $Disp
set Dmax [expr $Dmax*$Fact]; # scale value
if {$Dmax<0} {; # avoid the divide by zero
set dx [expr -$DincrStatic]
} else {
set dx $DincrStatic;
}
set NstepsPeak [expr int(abs($Dmax)/$DincrStatic)]
for {set i 1} {$i <= $NstepsPeak} {incr i 1} {; # zero to one
set Disp [expr $Disp + $dx]
puts $outFileID $Disp
}
if {$CycleType !="Push"} {
for {set i 1} {$i <= $NstepsPeak} {incr i 1} {; # one to zero
set Disp [expr $Disp - $dx]
puts $outFileID $Disp
}
if {$CycleType !="HalfCycle"} {
for {set i 1} {$i <= $NstepsPeak} {incr i 1} {; # zero to minus one
set Disp [expr $Disp - $dx]
puts $outFileID $Disp
}
for {set i 1} {$i <= $NstepsPeak} {incr i 1} {; # minus one to zero
set Disp [expr $Disp + $dx]
puts $outFileID $Disp
}
}
}
puts $outFileID " }"
close $outFileID
source data/tmpDsteps.tcl
return $iDstep
}; #
######################################################################################
######################################################################################
## procConvergeStatic $Tol
######################################################################################
# modify analysis parameters and perform analysis to find convergence for displacement-driven analysis
# by Silvia Mazzoni, 2006
#
proc procConvergeStatic { Tol } {
# if analysis fails, we try some other stuff
# performance is slower inside this loop
global maxNumIterConvergeStatic; # maximum number of iterations that will be performed before "failure to converge" is returned
global maxNumIterStatic; # maximum number of iterations that will be performed before "failure to converge" is returned
global NewtonLineSearchRatio;
global printFlagStatic;
global printFlagConvergeStatic;
global algorithmTypeStatic
puts "Trying Newton with Initial Tangent .."
test NormDispIncr $Tol $maxNumIterConvergeStatic $printFlagConvergeStatic
algorithm Newton -initial
set ok [analyze 1]
test NormDispIncr $Tol $maxNumIterStatic $printFlagStatic
algorithm $algorithmTypeStatic
if {$ok != 0} {
puts "Trying NewtonWithLineSearch .."
algorithm NewtonLineSearch $NewtonLineSearchRatio
set ok [analyze 1]
algorithm $algorithmTypeStatic
}
return $ok
}; #
######################################################################################
######################################################################################
## procAnalysisStatic $iDmax $IDctrlNode $IDctrlDOF $CycleType $Ncycles $Fact
######################################################################################
# perform displacement-controlled static analysis
# by Silvia Mazzoni, 2006
#
proc procAnalysisStatic { iDmax IDctrlNode IDctrlDOF {DincrStatic 0.01} {CycleType "Full"} {Ncycles 1} {Fact 1} } {
global LunitTXT
global in
global constraintsTypeStatic;
global alphaSP; # Penalty/Lagrange constraints -- factor used when adding the single-point constraint into the system of equations
global alphaMP; # Penalty/Lagrange constraints -- factor used when adding the multi-point constraint into the system of equations
global numbererType;
global systemTypeStatic ;
global testTypeStatic; # convergence-test type
global TolStatic; # Convergence Test: tolerance
global maxNumIterStatic; # Convergence Test: maximum number of iterations that will be performed before "failure to converge" is returned
global printFlagStatic; # Convergence Test: flag used to print information on convergence (optional)
global NstepStatic; # number of steps to apply Static
global algorithmTypeStatic;
global algorithmCount;
global SaveDatabase; # save to database trigger
set fmt1 "%s analysis: CtrlNode %.3i, dof %.1i, Disp=%.4f %s"
foreach Dmax $iDmax {
set iDstep [procGenPeaks $Dmax $DincrStatic $CycleType $Fact]; # this proc is defined above
for {set i 1} {$i <= $Ncycles} {incr i 1} {
set zeroD 0
set D0 0.0
foreach Dstep $iDstep {
set D1 $Dstep
set Dincr [expr $D1 - $D0]
if {$constraintsTypeStatic == "Plain" | $constraintsTypeStatic == "Transformation"} {
constraints $constraintsTypeStatic ; # Plain & Transformation constraints
} else {
constraints $constraintsTypeStatic $alphaSP $alphaMP ; # Penalty & Lagrange congs
}
if {$systemTypeStatic == "SparseGeneralPivot"} {
system SparseGeneral -piv; # optional pivoting for SparseGeneral system
} else {
system $systemTypeStatic ;
}
numberer $numbererType; # renumber dof's to minimize band-width (optimization), if you want to
test $testTypeStatic $TolStatic $maxNumIterStatic ; # tolerance, max no. of iterations, and print code , 1: every iteration
if {$algorithmTypeStatic == "BFGS" | $algorithmTypeStatic == "Broyden" } {
algorithm $algorithmTypeStatic; # use Newton's solution algorithm: updates tangent stiffness at every iteration
} elseif { $algorithmTypeStatic == "NewtonLineSearch"} {
algorithm $algorithmTypeStatic $NewtonLineSearchRatio; # use Newton's solution algorithm: updates tangent stiffness at every iteration
} else {
algorithm $algorithmTypeStatic $algorithmCount; # use Newton's solution algorithm: updates tangent stiffness at every iteration
}
integrator DisplacementControl $IDctrlNode $IDctrlDOF $Dincr
analysis Static
set ok [analyze 1]
if {$ok != 0} {
set Nk 10;
for {set k 1} {$k <= $Nk} {incr k 1} {
integrator DisplacementControl $IDctrlNode $IDctrlDOF [expr $Dincr/$Nk]
set ok [analyze 1]
if {$ok != 0} {
set Nm 2;
for {set m 1} {$m <= $Nm} {incr m 1} {
integrator DisplacementControl $IDctrlNode $IDctrlDOF [expr $Dincr/$Nk/$Nm]
set ok [analyze 1]
if {$ok != 0} {
set ok [procConvergeStatic $TolStatic]; # this proc is defined above
if {$ok != 0} {
set ok [procConvergeStatic [expr $TolStatic*1000]]; # increase tolerance a little before giving up!
if {$ok != 0} {
puts [format $fmt1 "PROBLEM" $IDctrlNode $IDctrlDOF [nodeDisp $IDctrlNode $IDctrlDOF] $LunitTXT]
return -1
}
}
}; # end if statement
}; # end for m
}; # end if ok statement
}; # for k
}; # end if ok statement
set D0 $D1
}; # end Dstep
}; # end i
if {$SaveDatabase=="on" } {
save 1; # save to database the state at the end of last peak cycle
}
}; # end of iDmaxCycl
# announce if we made if over the convergence problems:
if {$ok != 0 } {
puts [format $fmt1 "PROBLEM" $IDctrlNode $IDctrlDOF [nodeDisp $IDctrlNode $IDctrlDOF] $LunitTXT]
} else {
puts [format $fmt1 "DONE" $IDctrlNode $IDctrlDOF [nodeDisp $IDctrlNode $IDctrlDOF] $LunitTXT]
}
}; #
######################################################################################
######################################################################################
## procLoadCtrlStaticAnalysis
######################################################################################
# perform force-controlled static analysis -- modify this script as needed
# by Silvia Mazzoni, 2006
# For some reason, you should have 0. 0. as the first two points, to make sure your output starts at zero.
# And also repeat the last load (7.0 7.0), as otherwise it only goes to one step before it.
proc procLoadCtrlStaticAnalysis {} {
set Tol 1e-6; # convergence tolerance
set dt .1
set Fseries "Series -dt $dt -values {0. 0. +5.0 7.0 7.0}"; # load series information
set LoadIDstatic 900
pattern Plain $LoadIDstatic $Fseries {
load 31 1.5 0 0
load 41 3.2 0 0
}
constraints Plain; # how it handles boundary conditions
numberer Plain; # renumber dof's to minimize band-width (optimization), if you want to
system BandGeneral
test NormDispIncr $Tol 6 ; # tolerance, max no. of iterations, and print code , 1: every iteration
algorithm Newton;# use Newton's solution algorithm: updates tangent stiffness at every iteration
set Nstep 100; # apply gravity in Nstep steps
set Dt [expr $dt*([llength $Fseries]-2)/$Nstep]; # first load increment;
integrator LoadControl $Dt
analysis Static
analyze $Nstep
}; #
######################################################################################
Code: Select all
# some general information:
variable problemSize Large; # option, Large or Small (less then 10 nodes)
# ANALYZE
variable NstepGravity 10; # number of integration steps for applying gravity.
# variable DincrStatic [expr 0.01*$in]; # displacement increment for static analysis
variable SaveDatabase "off"; # save to database trigger "on" saves at the end of each cycle
# CONSTRAINTS variables: determines how the constraint equations are enforced in the analysis
variable constraintsTypeGravity Plain; # options: Plain, Penalty, Lagrange, Transformation
variable constraintsTypeStatic Plain; # options: Plain, Penalty, Lagrange, Transformation
variable constraintsTypeDynamic Transformation; # options: Plain, Penalty, Lagrange, Transformation
variable alphaSP 1e6 ; # Penalty/Lagrange constraints -- factor used when adding the single-point constraint into the system of equations
variable alphaMP 1e6 ;# Penalty/Lagrange constraints -- factor used when adding the multi-point constraint into the system of equations
# NUMBERER: determines the mapping between equation numbers and degrees-of-freedom
variable numbererType Plain; # options: Plain, RCM
# SYSTEM: how to store and solve the system of equations in the analysis
variable systemTypeGravity BandGeneral; # options: BandGeneral, BandSPD, ProfileSPD, SparseGeneral, SparseGeneralPivot, UmfPack, SparseSPD
variable systemTypeStatic BandGeneral; # options: BandGeneral, BandSPD, ProfileSPD, SparseGeneral, SparseGeneralPivot, UmfPack, SparseSPD
variable systemTypeDynamic SparseGeneralPivot; # options: BandGeneral, BandSPD, ProfileSPD, SparseGeneral, SparseGeneralPivot, UmfPack, SparseSPD
# TEST: # convergence test to determine if convergence has been achieved at the end of an iteration step
variable testTypeGravity NormDispIncr; # options: NormUnbalance, NormDispIncr, EnergyIncr, RelativeNormUnbalance, RelativeNormDispIncr, RelativeEnergyIncr
# variable testTypeStatic RelativeNormDispIncr; # options: NormUnbalance, NormDispIncr, EnergyIncr, RelativeNormUnbalance, RelativeNormDispIncr, RelativeEnergyIncr
variable testTypeStatic NormDispIncr; # options: NormUnbalance, NormDispIncr, EnergyIncr, RelativeNormUnbalance, RelativeNormDispIncr, RelativeEnergyIncr
variable testTypeDynamic RelativeNormDispIncr; # options: NormUnbalance, NormDispIncr, EnergyIncr, RelativeNormUnbalance, RelativeNormDispIncr, RelativeEnergyIncr
variable TolGravity 1.e-6; # Convergence Test: tolerance
variable TolStatic 1.e-8; # Convergence Test: tolerance
# variable TolStatic 1.e-6; # Convergence Test: tolerance
variable TolDynamic 1.e-8; # Convergence Test: tolerance
variable maxNumIterGravity 6; # Convergence Test: maximum number of iterations that will be performed before "failure to converge" is returned
variable maxNumIterStatic 6; # Convergence Test: maximum number of iterations that will be performed before "failure to converge" is returned
variable maxNumIterDynamic 6; # Convergence Test: maximum number of iterations that will be performed before "failure to converge" is returned
# procConvergeStatic -- assist in achieving convergence
if {$problemSize=="Large"} {
variable maxNumIterConvergeStatic 5000; # Convergence Test: maximum number of iterations that will be performed before "failure to converge" is returned
variable maxNumIterConvergeDynamic 5000; # Convergence Test: maximum number of iterations that will be performed before "failure to converge" is returned
} else {
variable maxNumIterConvergeStatic 200; # Convergence Test: maximum number of iterations that will be performed before "failure to converge" is returned
variable maxNumIterConvergeDynamic 200; # Convergence Test: maximum number of iterations that will be performed before "failure to converge" is returned
}
variable printFlagGravity 0; # Convergence Test: flag used to print information on convergence (optional) # 0: no print output (default);
variable printFlagStatic 0; # Convergence Test: flag used to print information on convergence (optional) # 1: print information on each step;
variable printFlagDynamic 0; # Convergence Test: flag used to print information on convergence (optional) # 2: print information when convergence has been achieved;
variable printFlagConvergeStatic 0; # Convergence Test: flag used to print information on convergence (optional) 1: on, 0: off # 4: print norm, dU and dR vectors;
variable printFlagConvergeDynamic 0; # Convergence Test: flag used to print information on convergence (optional) 1: on, 0: off # 5: at convergence failure, carry on, print error message, but do not stop analysis
# ALGORITHM:
variable algorithmTypeGravity Newton; # options: Linear, Newton, NewtonLineSearch, ModifiedNewton, KrylovNewton, BFGS, Broyden,
variable algorithmTypeStatic Newton; # options: Linear, Newton, NewtonLineSearch, ModifiedNewton, KrylovNewton, BFGS, Broyden,
variable algorithmTypeDynamic Newton; # options: Linear, Newton, NewtonLineSearch, ModifiedNewton, KrylovNewton, BFGS, Broyden,
variable NewtonLineSearchRatio 0.8; # Algorithm: NewtonLineSearch, limiting ratio between the residuals before and after the incremental update (between 0.5 and 0.8)
variable algorithmCount 5; # Algorithm: BFGS/Broyden, number of iterations within a time step until a new tangent is formed
# INTEGRATOR: defined locally for static analysis.
variable integratorTypeDynamic Newmark; # options: Newmark, HHT
variable NewmarkGamma 0.5; # Newmark-integrator gamma parameter (also HHT)
variable NewmarkBeta 0.25; # Newmark-integrator gamma parameter
# ANALYSIS
variable analysisTypeDynamic Transient; # options: Transient, VariableTransient
# RAYLEIGH damping parameters, Where to put M/K-prop damping, switches
# D=$alphaM*M + $betaKcurr*Kcurrent + $betaKcomm*KlastCommit + $beatKinit*$Kinitial
variable MpropSwitch 1.0;
variable KcurrSwitch 0.0;
variable KcommSwitch 0.0;
variable KinitSwitch 1.0;
# DISPLAY variables
variable xPixels 600; # height of graphical window in pixels
variable yPixels 400; # height of graphical window in pixels
variable xLoc1 10; # horizontal location of graphical window (0=upper left-most corner)
variable yLoc1 10; # vertical location of graphical window (0=upper left-most corner)
variable xLoc2 $xLoc1; # horizontal location of graphical window (0=upper left-most corner)
variable yLoc2 [expr $yLoc1+$yPixels]; # vertical location of graphical window (0=upper left-most corner)
variable xLoc3 [expr $xLoc1+$xPixels]; # horizontal location of graphical window (0=upper left-most corner)
variable yLoc3 $yLoc1; # vertical location of graphical window (0=upper left-most corner)
# LOAD-TAG variables
variable loadIDgravity 100; # gravity load ID tag
variable loadIDstatic 200; # static load ID tag
variable IDloadTagGMA 310; # ground-motion load ID tag
variable IDloadTagGMB 320; # ground-motion load ID tag
variable IDloadTagGMC 330; # ground-motion load ID tag
variable IDgmSeries 350; # ground-motion series ID tag
Code: Select all
# view model -- node numbers
set xPixels 600; # height of graphical window in pixels
set yPixels 400; # height of graphical window in pixels
set xLoc1 10; # horizontal location of graphical window (0=upper left-most corner)
set yLoc1 10; # vertical location of graphical window (0=upper left-most corner)
set xLoc2 $xLoc1; # horizontal location of graphical window (0=upper left-most corner)
set yLoc2 [expr $yLoc1+$yPixels]; # vertical location of graphical window (0=upper left-most corner)
set xLoc3 [expr $xLoc1+$xPixels]; # horizontal location of graphical window (0=upper left-most corner)
set yLoc3 $yLoc1; # vertical location of graphical window (0=upper left-most corner)
set dAmp 5; # relative amplification factor for deformations
set viewEigen 3; # eigenmode to be viewed
procDisplayShape2D DeformedShape $dAmp $xLoc1 $yLoc1 $xPixels $yPixels
procDisplayShape2D NodeNumbers $dAmp $xLoc2 $yLoc2 $xPixels $yPixels
procDisplayShape2D ModeShape $dAmp $xLoc3 $yLoc3 $xPixels $yPixels $viewEigen
# GRAVITY LOADS # define gravity load applied to beams and columns -- eleLoad applies loads in local coordinate axis
pattern Plain 101 Linear {
eleLoad -ele 221 222 223 -type -beamUniform -$QdlBeam; ; # beams level 2 (in -ydirection)
eleLoad -ele 231 232 233 -type -beamUniform -$QdlBeam;
eleLoad -ele 241 242 243 -type -beamUniform -$QdlBeam
eleLoad -ele 111 112 113 114 -type -beamUniform 0 -$QdlCol; # columns level 1-2 (in -xdirection)
eleLoad -ele 121 122 123 124 -type -beamUniform 0 -$QdlCol;
eleLoad -ele 131 132 133 134 -type -beamUniform 0 -$QdlCol;
}
procApplyGravity ; # apply gravity load, set it constant and reset time to zero.
recorder Node -file RNode.out -time -node 11 12 13 14 -dof 1 reaction
set InputType EQ
set InputType Static
if {$InputType == "Static" } {
# Define Lateral Loads for Static Pushover Analysis, distributed along building height
# Fj = WjHj/sum(WiHi) * Weight
set LoadIDstatic 200; # assign load tag to static pushover load (optional, default 200)
set iFj [procFlateral $iLcol $iFloorWeight]
set Y0 0.
set iYlevel ""
foreach Lcol $iLcol {
set Y0 [expr $Y0+$Lcol];
lappend iYlevel $Y0
}
set iIDpushNode ""
set level 1
foreach Ylevel $iYlevel {
set level [expr $level+1]
set pier 1
set IDpushNode [expr $level*10+$pier]
lappend iIDpushNode $IDpushNode
}
# create load pattern for static pushover loads
pattern Plain $LoadIDstatic Linear {
foreach IDpushNode $iIDpushNode Fj $iFj {
load $IDpushNode $Fj 0 0
}
}
set AnalysisTypeTXT Push
set Lcol [expr $Lcol1+$Lcol2+$Lcol3];
if {$AnalysisTypeTXT=="Push" } {
set iDmax 0.11
set fact $Lcol
set CycleType Push
set Ncycles 1
} elseif {$AnalysisTypeTXT=="Cycl" } {
set iDmax "[expr 0.001] [expr 0.005] [expr 0.01] [expr 0.025] [expr 0.05] [expr 0.075] [expr 0.09] [expr 0.11]"
set fact $Lcol
set CycleType Full
set Ncycles 1
} else {
puts "No Analysis Type Specified"
return
}
set IDctrlNode 41
set IDctrlDOF 1
procAnalysisStatic $iDmax $IDctrlNode $IDctrlDOF $CycleType $Ncycles $fact
} else {
puts "No Input Type Specified"
}
okey, here is my libDisplay.tcl file, too -- the display command only works on the latest release of OpenSees 1.7.1:
Code: Select all
######################################################################################
## procDisplayPlane $ShapeType $dAmp $viewPlane $nEigen $quadrant
######################################################################################
proc procDisplayPlane {ShapeType dAmp viewPlane {nEigen 0} {quadrant 0}} {
######################################################################################
## setup display parameters for specified viewPlane
######################################################################################
#### -- by Silvia Mazzoni, 2006 (mazzoni@berkeley_NO_SPAM_.edu)
####
#### ShapeType : type of shape to display. # options: ModeShape , NodeNumbers , DeformedShape
#### dAmp : relative amplification factor for deformations
#### viewPlane : set local xy axes in global coordinates (XY,YX,XZ,ZX,YZ,ZY)
#### nEigen : if nEigen not=0, show mode shape for nEigen eigenvalue
#### quadrant: quadrant where to show this figure (0=full figure)
####
######################################################################################
set Xmin [lindex [nodeBounds] 0]; # view bounds in global coords - proc will add padding on the sides
set Ymin [lindex [nodeBounds] 1];
set Zmin [lindex [nodeBounds] 2];
set Xmax [lindex [nodeBounds] 3];
set Ymax [lindex [nodeBounds] 4];
set Zmax [lindex [nodeBounds] 5];
set Xo 0; # center of local viewing system
set Yo 0;
set Zo 0;
set uLocal [string index $viewPlane 0]; # viewPlane local-x axis in global coordinates
set vLocal [string index $viewPlane 1]; # viewPlane local-y axis in global coordinates
if {$viewPlane =="YZ" } {
set uMin $Ymin
set uMax $Ymax
set vMin $Zmin
set vMax $Zmax
set wMin $Xmin
set wMax $Xmax
} elseif {$viewPlane =="ZY" } {
set uMin $Zmin
set uMax $Zmax
set vMin $Ymin
set vMax $Ymax
set wMin $Xmin
set wMax $Xmax
} elseif {$viewPlane =="XY" } {
set uMin $Xmin
set uMax $Xmax
set vMin $Ymin
set vMax $Ymax
set wMin $Zmin
set wMax $Zmax
} elseif {$viewPlane =="YX" } {
set uMin $Ymin
set uMax $Ymax
set vMin $Xmin
set vMax $Xmax
set wMin $Zmin
set wMax $Zmax
} elseif {$viewPlane =="XZ" } {
set uMin $Xmin
set uMax $Xmax
set vMin $Zmin
set vMax $Zmax
set wMin $Ymin
set wMax $Ymax
} elseif {$viewPlane =="ZX" } {
set uMin $Zmin
set uMax $Zmax
set vMin $Xmin
set vMax $Xmax
set wMin $Ymin
set wMax $Ymax
} elseif {$viewPlane =="3D" } {
set uMin $Zmin+$Xmin
set uMax $Zmax+$Xmax
set vMin $Ymin
set vMax $Ymax
set wMin -10000
set wMax 10000
vup 0 1 0; # dirn defining up direction of view plane
} else {
return -1
}
set epsilon 1e-3; # make windows width or height not zero when the Max and Min values of a coordinate are the same
set uWide [expr $uMax - $uMin+$epsilon];
set vWide [expr $vMax - $vMin+$epsilon];
set uSide [expr 0.25*$uWide];
set vSide [expr 0.25*$vWide];
set uMin [expr $uMin - $uSide];
set uMax [expr $uMax + $uSide];
set vMin [expr $vMin - $vSide];
set vMax [expr $vMax + $vSide];
set uWide [expr $uMax - $uMin+$epsilon];
set vWide [expr $vMax - $vMin+$epsilon];
set uMid [expr ($uMin+$uMax)/2];
set vMid [expr ($vMin+$vMax)/2];
# keep the following general, as change the X and Y and Z for each view plane
# next three commmands define viewing system, all values in global coords
vrp $Xo $Yo $Zo; # point on the view plane in global coord, center of local viewing system
if {$vLocal == "X"} {
vup 1 0 0; # dirn defining up direction of view plane
} elseif {$vLocal == "Y"} {
vup 0 1 0; # dirn defining up direction of view plane
} elseif {$vLocal == "Z"} {
vup 0 0 1; # dirn defining up direction of view plane
}
if {$viewPlane =="YZ" } {
vpn 1 0 0; # direction of outward normal to view plane
prp 10000. $uMid $vMid ; # eye location in local coord sys defined by viewing system
plane 10000 -10000; # distance to front and back clipping planes from eye
} elseif {$viewPlane =="ZY" } {
vpn -1 0 0; # direction of outward normal to view plane
prp -10000. $vMid $uMid ; # eye location in local coord sys defined by viewing system
plane 10000 -10000; # distance to front and back clipping planes from eye
} elseif {$viewPlane =="XY" } {
vpn 0 0 1; # direction of outward normal to view plane
prp $uMid $vMid 10000; # eye location in local coord sys defined by viewing system
plane 10000 -10000; # distance to front and back clipping planes from eye
} elseif {$viewPlane =="YX" } {
vpn 0 0 -1; # direction of outward normal to view plane
prp $uMid $vMid -10000; # eye location in local coord sys defined by viewing system
plane 10000 -10000; # distance to front and back clipping planes from eye
} elseif {$viewPlane =="XZ" } {
vpn 0 -1 0; # direction of outward normal to view plane
prp $uMid -10000 $vMid ; # eye location in local coord sys defined by viewing system
plane 10000 -10000; # distance to front and back clipping planes from eye
} elseif {$viewPlane =="ZX" } {
vpn 0 1 0; # direction of outward normal to view plane
prp $uMid 10000 $vMid ; # eye location in local coord sys defined by viewing system
plane 10000 -10000; # distance to front and back clipping planes from eye
} elseif {$viewPlane =="3D" } {
vpn 1 0.25 1; # direction of outward normal to view plane
prp -100 $vMid 10000; # eye location in local coord sys defined by viewing system
plane 10000 -10000; # distance to front and back clipping planes from eye
} else {
return -1
}
# next three commands define view, all values in local coord system
if {$viewPlane =="3D" } {
viewWindow [expr $uMin-$uWide/4] [expr $uMax/2] [expr $vMin-0.25*$vWide] [expr $vMax]
} else {
viewWindow $uMin $uMax $vMin $vMax
}
projection 1; # projection mode, 0:prespective, 1: parallel
fill 1; # fill mode; needed only for solid elements
if {$quadrant == 0} {
port -1 1 -1 1 # area of window that will be drawn into (uMin,uMax,vMin,vMax);
} elseif {$quadrant == 1} {
port 0 1 0 1 # area of window that will be drawn into (uMin,uMax,vMin,vMax);
} elseif {$quadrant == 2} {
port -1 0 0 1 # area of window that will be drawn into (uMin,uMax,vMin,vMax);
} elseif {$quadrant == 3} {
port -1 0 -1 0 # area of window that will be drawn into (uMin,uMax,vMin,vMax);
} elseif {$quadrant == 4} {
port 0 1 -1 0 # area of window that will be drawn into (uMin,uMax,vMin,vMax);
}
if {$ShapeType == "ModeShape" } {
display -$nEigen 0 [expr 5.*$dAmp]; # display mode shape for mode $nEigen
} elseif {$ShapeType == "NodeNumbers" } {
display 1 -1 0 ; # display node numbers
} elseif {$ShapeType == "DeformedShape" } {
display 1 5 $dAmp; # display deformed shape
}
}; #
######################################################################################
######################################################################################
## procDisplayShape3D $ShapeType $dAmp $xLoc $yLoc $xPixels $yPixels $nEigen
######################################################################################
proc procDisplayShape3D { ShapeType {dAmp 5} {xLoc 10} {yLoc 10} {xPixels 750} {yPixels 600} {nEigen 1} } {
######################################################################################
## display Node Numbers, Deformed or Mode Shape in all 3 planes
######################################################################################
#### -- by Silvia Mazzoni, 2006 (mazzoni@berkeley_NO_SPAM_.edu)
####
#### ShapeType : type of shape to display. # options: ModeShape , NodeNumbers , DeformedShape
#### dAmp : relative amplification factor for deformations
#### xLoc,yLoc : horizontal & vertical location in pixels of graphical window (0,0=upper left-most corner)
#### xPixels,yPixels : width & height of graphical window in pixels
#### nEigen : if nEigen not=0, show mode shape for nEigen eigenvalue
####
#######################################################################################
global TunitTXT ; # load global unit variable
set Xmin [lindex [nodeBounds] 0]; # view bounds in global coords - proc will add padding on the sides
set Ymin [lindex [nodeBounds] 1];
set Zmin [lindex [nodeBounds] 2];
set Xmax [lindex [nodeBounds] 3];
set Ymax [lindex [nodeBounds] 4];
set Zmax [lindex [nodeBounds] 5];
if {$ShapeType == "ModeShape" } {
set lambdaN [eigen $nEigen]; # perform eigenvalue analysis for ModeShape
set lambda [lindex $lambdaN [expr $nEigen-1]];
set omega [expr pow($lambda,0.5)]
set PI [expr 2*asin(1.0)]; # define constant
set Tperiod [expr 2*$PI/$omega]; # period
set fmt1 "Mode Shape, Mode=%.1i Period=%.3f %s "
set windowTitle [format $fmt1 $nEigen $Tperiod $TunitTXT ]
} elseif {$ShapeType == "NodeNumbers" } {
set windowTitle "Node Numbers"
} elseif {$ShapeType == "DeformedShape" } {
set windowTitle0 "Deformed Shape "
}
if {$ShapeType == "DeformedShape" } {
set xPixels [expr int($xPixels/2)]
set yPixels [expr int($yPixels/2)]
set xLoc1 [expr $xLoc+$xPixels]
set yLoc1 [expr $yLoc+$yPixels]
set planeTXT "-Plane"
set viewPlane XY
set windowTitle $windowTitle0$viewPlane$planeTXT
recorder display $windowTitle $xLoc1 $yLoc $xPixels $yPixels -wipe ; # display recorder
procDisplayPlane $ShapeType $dAmp $viewPlane
set viewPlane ZY
set windowTitle $windowTitle0$viewPlane$planeTXT
recorder display $windowTitle $xLoc $yLoc $xPixels $yPixels -wipe ; # display recorder
procDisplayPlane $ShapeType $dAmp $viewPlane
set viewPlane ZX
set windowTitle $windowTitle0$viewPlane$planeTXT
recorder display $windowTitle $xLoc $yLoc1 $xPixels $yPixels -wipe ; # display recorder
procDisplayPlane $ShapeType $dAmp $viewPlane
set viewPlane 3D
set windowTitle $windowTitle0$viewPlane
recorder display $windowTitle $xLoc1 $yLoc1 $xPixels $yPixels -wipe ; # display recorder
procDisplayPlane $ShapeType $dAmp $viewPlane
} else {
recorder display $windowTitle $xLoc $yLoc $xPixels $yPixels -file $windowTitle -nowipe; # display recorder
set viewPlane XY
procDisplayPlane $ShapeType $dAmp $viewPlane $nEigen 1
set viewPlane ZY
procDisplayPlane $ShapeType $dAmp $viewPlane $nEigen 2
set viewPlane ZX
procDisplayPlane $ShapeType $dAmp $viewPlane $nEigen 3
set viewPlane 3D
procDisplayPlane $ShapeType $dAmp $viewPlane $nEigen 4
}
}; #
######################################################################################
######################################################################################
## procDisplayShape2D $ShapeType $dAmp $xLoc $yLoc $xPixels $yPixels $nEigen
######################################################################################
proc procDisplayShape2D { ShapeType {dAmp 5} {xLoc 10} {yLoc 10} {xPixels 750} {yPixels 600} {nEigen 0} } {
######################################################################################
## display Node Numbers, Deformed or Mode Shape in 2D problem
######################################################################################
#### -- by Silvia Mazzoni, 2006 (mazzoni@berkeley_NO_SPAM_.edu)
####
#### ShapeType : type of shape to display. # options: ModeShape , NodeNumbers , DeformedShape
#### dAmp : relative amplification factor for deformations
#### xLoc,yLoc : horizontal & vertical location in pixels of graphical window (0,0=upper left-most corner)
#### xPixels,yPixels : width & height of graphical window in pixels
#### nEigen : if nEigen not=0, show mode shape for nEigen eigenvalue
####
#######################################################################################
global TunitTXT
set Xmin [lindex [nodeBounds] 0]; # view bounds in global coords - proc will add padding on the sides
set Ymin [lindex [nodeBounds] 1];
set Zmin [lindex [nodeBounds] 2];
set Xmax [lindex [nodeBounds] 3];
set Ymax [lindex [nodeBounds] 4];
set Zmax [lindex [nodeBounds] 5];
if {$ShapeType == "ModeShape" } {
set lambdaN [eigen $nEigen]; # perform eigenvalue analysis for ModeShape
set lambda [lindex $lambdaN [expr $nEigen-1]];
set omega [expr pow($lambda,0.5)]
set PI [expr 2*asin(1.0)]; # define constant
set Tperiod [expr 2*$PI/$omega]; # period (sec.)
set fmt1 "Mode Shape, Mode=%.1i Period=%.3f %s "
set windowTitle [format $fmt1 $nEigen $Tperiod $TunitTXT]
} elseif {$ShapeType == "NodeNumbers" } {
set windowTitle "Node Numbers"
} elseif {$ShapeType == "DeformedShape" } {
set windowTitle "Deformed Shape"
}
set viewPlane XY
recorder display $windowTitle $xLoc $yLoc $xPixels $yPixels -wipe ; # display recorder
procDisplayPlane $ShapeType $dAmp $viewPlane $nEigen 0
}; #
######################################################################################
######################################################################################
## procDisplayDefaults $dAmp
######################################################################################
proc procDisplayAll { {dAmp 5} } {
######################################################################################
## display Node Numbers, Deformed AND Mode Shape using default values.
######################################################################################
# view model -- node numbers
global xPixels;
global yPixels;
global xLoc1;
global yLoc1;
global xLoc2;
global yLoc2;
global xLoc3;
global yLoc3;
# set dAmp 5; # relative amplification factor for deformations
set viewEigen 1; # eigenmode to be viewed
procDisplayShape2D DeformedShape $dAmp $xLoc1 $yLoc1 $xPixels $yPixels
procDisplayShape2D NodeNumbers $dAmp $xLoc2 $yLoc2 $xPixels $yPixels
procDisplayShape2D ModeShape $dAmp $xLoc3 $yLoc3 $xPixels $yPixels $viewEigen
}; #
######################################################################################
######################################################################################
## procDisplayDefaults $dAmp
######################################################################################
proc procDisplayDefoShape { {dAmp 5} } {
######################################################################################
## display Node Numbers, Deformed AND Mode Shape using default values.
######################################################################################
# view model -- node numbers
global xPixels;
global yPixels;
global xLoc1;
global yLoc1;
global xLoc2;
global yLoc2;
global xLoc3;
global yLoc3;
# set dAmp 5; # relative amplification factor for deformations
set viewEigen 1; # eigenmode to be viewed
procDisplayShape2D DeformedShape $dAmp $xLoc1 $yLoc1 800 600
}; #
######################################################################################as you can see, it is a bit "complicated", but pretty clean and self-explanatory.