General element end fixity

[mpnewcombe]

Hi There,

I am attempting to make a general conversion program to transfer a 3D static model (from Space Gass) into an OS input file. The biggest problem I appear to be having is with the element end fixity conditions. I insert another node at the end of each element and constrain the dof between the end of the element and this node to simulate the appropriate end fixity.

I am trying to run a ground motion through a simple cantilever and the message I get is:
constraint at dof 0 already specified for constrained node in MP_Constraint at node x.

My Inputs are as follows:

master_input.tcl ########################################################################################################################

# SETUP
wipe; # Clear memory of all past model definisions
model BasicBuilder -ndm 3 -ndf 6 ; # Define model builder, ndm = number of dimensions, ndf = number of DOF
source DisplayModel3D_input.tcl ; # Display model
source DisplayPlane_input.tcl ; # Display a plane in a model
set dataDir Data; # Setup name of output data directory
file mkdir $dataDir;
set GMdir "EQS/"; # Ground-motion file directory

# UNITS
set mm 1
set N 1
set sec 1
set LunitTXT "mm"
set FunitTXT "N"
set TunitTXT "sec"
set m [expr 1000*$mm]
set kN [expr 1000*$N]
set mm2 [expr $mm*$mm]
set MPa [expr $N/$mm2]
set pi 3.141592654
set g [expr 9.81*$m/($sec*$sec)]

# NODES
source node_input.tcl ;
source node_end1_input.tcl ;
source node_end2_input.tcl ;

# NODE FIXITY
source node_fixity_input.tcl ;

# SECTION PROPERTIES
source section_input.tcl ;

# MATERIAL PROPERTIES
source material_input.tcl ;

# AGG SECTION PROPERTIES
source agg_section_input.tcl ;


# ELEMENTS
# Geometric transformation
# Beam and columns are treated seperately if P-delt analysis is performed
set IDBeamTransf 1 ; # All beams have transform ID of 1
set IDColTransf 2 ; # All columns have transform ID of 2
set ColTransType Linear ; # Can use either Linear, PDelta or Corotational.
geomTransf Linear $IDBeamTransf 0 0 1 ; # Specify geometric transformation for beams (along x -axis)
geomTransf $ColTransType $IDColTransf 0 0 1 ; # Specify geometric transformation for columns (along y axis).

source mod_element_input.tcl ;
# Element end constraints
source element_end1_constraint_input.tcl ;
source element_end2_constraint_input.tcl ;

# SEISMIC MASSES
source lumped_mass_input.tcl ;

# GRAVITY LOADS
# Subroutine for gravity loads
pattern Plain 1 Linear {
source loads_input.tcl ;
}

# Iterative solution options
constraints Plain ; # How boundary conditions are handled
numberer Plain ; # Renumbering of dof to optimize band-width
system BandGeneral ; # How to store and solve the system of equations in the analysis
test NormDispIncr 1.00E-08 6 ; # Determine if convergence has been achieved at end of interation step
algorithm KrylovNewton ; # Krylov-Newton's solution algorithum
integrator LoadControl 0.1 ; # Determine the next time step for an analysis - apply gravity in 10 steps
analysis Static ; # Type of analysis: Static or Transient.
analyze 10 ; # Perform gravity analysis (gravity code = 10)
loadConst -time 0 ; # Hold gravity constant and restart time

# RECORDERS
# Time-history data
set dataSubDir EQ1 ;
file mkdir data/$dataSubDir;

# Create regions
source regions_input.tcl ;

# Floor displacements
recorder Node -file $dataDir/$dataSubDir/Disp.out -time -node 1 2 -dof 1 disp ;

# Record drift histories
recorder Drift -file $dataDir/$dataSubDir/Drift.out -time -iNode 1 -jNode 2 -dof 1 -perpDirn 2

# Columns shears and Moments
recorder Element -file $dataDir/$dataSubDir/MemActions1.out -time -region 1 globalForce ;

# READ GM FILE
# Read PEER GM File format
proc ReadSMDFile {inFilename outFilename dt} {
# Reference for dt
upvar $dt DT ;
# Give error if file can't be read
if [catch {open $inFilename r} inFileID] {
puts stderr "Cannot open $inFilename for reading"
} else {
# Open output file for writing
set outFileID [open $outFilename w]
# Flag for finding dt - assumes second value in second column is time-step
set flag 0
# Look at each line in the file
foreach line [split [read $inFileID] \n] {
if {[llength $line]==0} { # If line is blank do nothing
continue
} elseif {$flag == 1} { # Echo ground motion values to output file
puts $outFileID $line
} else { # Search header lines for dt
foreach word [split $line] {
if {$flag == 1} {
set DT $word
break
}
if {[string match $word "DT="] == 1} { # Find dt and set the flag
set flag 1
}
}
}
}
# Close output and input file
close $outFileID ;
close $inFileID ;
}
}

# GROUND MOTION (uniform input for at all support nodes)
set GMdirection 1 ; # Ground-motion direction
set GMfile "EQ1" ; # Ground-motion filenames
set GMfact 0.93 ; # Ground-motion scaling factor

# GROUND-MOTION-ANALYSIS PARAMETERS
set DtAnalysis 0.005 ; # Time-step Dt for Lateral Analysis (seconds)
set TmaxAnalysis [expr 50*$sec]; # Maximum duration for ground-motion analysis
# Constraints handler - Determines how contraint equations are enforced in analysis
variable constraintsTypeDynamic Transformation;
constraints $constraintsTypeDynamic ;
# DOF Renumbering
variable numbererTypeDynamic RCM ; # Plain uses the numbering provided by user, RCM renumbers DOF to minimize band-width.
numberer $numbererTypeDynamic ;
# System Linear Equation Solvers
variable systemTypeDynamic BandGeneral;
system $systemTypeDynamic ;
# Convergance tests
variable TolDynamic 1.00E-08 ; # Convergance tolerance
variable maxNumIterDynamic 10 ; # Maximum number of interations (if exceeded "failure to converge" will appear)
variable pringFlagDynamic 0 ; # Option for printing information on convergance (1: to print).
variable testTypeDynamic EnergyIncr ; # Convergance test type
test $testTypeDynamic $TolDynamic $maxNumIterDynamic $pringFlagDynamic;
# For improved convergence procedure (not required here)
variable maxNumIterConvergeDynamic 2000 ;
variable pringFlagConvergeDynamic 0 ;
# Solution algorithum
variable algorithmTypeDynamic KrylovNewton; # Uses subspace accelerator to accelerate the convergence of the modified newton method.
algorithm $algorithmTypeDynamic ;
# Integrator Type
variable NewmarkGamma 0.5 ; # Newmark intergrator gamma parameter
variable NewmarkBeta 0.25 ; # Newmark intergrator beta parameter
variable integratorTypeDynamic Newmark ;
integrator $integratorTypeDynamic $NewmarkGamma $NewmarkBeta ;
# Analysis Type (Static solves KU = R without mass or damping matrices)
variable analysisTypeDynamic Transient ;
analysis $analysisTypeDynamic ;

# DISPLAY DEFORMED SHAPE
set ViewScale 10 ; # Amplify display of deformed shape
DisplayModel3D DeformedShape $ViewScale ; # Display deformed shape, the scaling factor needs to be adjusted for each model

# ELASTIC DAMPING
# Opensees only considers Rayleigh Damping: D=$alphaM*M + $betaKcurr*Kcurrent + $betaKcomm*KlastCommit + $beatKinit*$Kinitial
set xDamp 15 ; # Damping ratio
set MpropSwitch 1 ; # Switch on/off mass proportional damping
set KcurrSwitch 0 ; # Switch on/off tangent stiffness proportional damping
set KcommSwitch 1 ; # Switch on/off tangent stiffness proportional damping at last time-step
set KinitSwitch 0 ; # Switch on/off initial stiffness proportional damping
set nEigenI 1 ; # Eigenvalue for Mode 1
set nEigenJ 1 ; # Eigenvalue for Mode 3
set lambdaN [eigen [expr $nEigenJ]] ; # Define number of eigenvalue analyses required
set lambdaI [lindex $lambdaN [expr $nEigenI-1]] ; # Eigenvalue solution for mode I
set lambdaJ [lindex $lambdaN [expr $nEigenJ-1]] ; # Eigenvalue solution for mode J
set omegaI [expr pow($lambdaI,0.5)] ; # Angular frequency for Mode I
set omegaJ [expr pow($lambdaJ,0.5)] ; # Angular frequency for Mode J
set alphaM [expr $MpropSwitch*$xDamp*(2*$omegaI*$omegaJ)/($omegaI+$omegaJ)]; # M-prop damping; D = alphaM*M
set betaKcurr [expr $KcurrSwitch*2.*$xDamp/($omegaI+$omegaJ)]; # Current-K; +beatKcurr*KCurrent
set betaKcomm [expr $KcommSwitch*2*$xDamp/($omegaI+$omegaJ)]; # Last-committed K; +betaKcomm*KlastCommitt
set betaKinit [expr $KinitSwitch*2.*$xDamp/($omegaI+$omegaJ)]; # Initial-K; +beatKinit*Kini

rayleigh $alphaM $betaKcurr $betaKinit $betaKcomm ; # Rayleigh damping code

# PERFORMING DYNAMIC GROUND-MOTION ANALYSIS
# Uniform earthquake excitation
set IDloadTag 400 ; # Tag for uniform support excitation
set inFile $GMdir$GMfile.txt ; # Input name
set outFile $GMdir$GMfile.g3 ; # Set variable holding new filename
ReadSMDFile $inFile $outFile dt ; # Call procedure to convert the GM file
set GMfatt [expr $g*$GMfact] ; # Data in input file is in units of 'g' and is scaled by GMfact
set AccelSeries "Series -dt $dt -filePath $outFile -factor $GMfatt" ; # Time series information
pattern UniformExcitation $IDloadTag $GMdirection -accel $AccelSeries ; # Create uniform excitation
# Perform analysis
set Nsteps [expr int($TmaxAnalysis/$DtAnalysis)]; # Number of steps in analysis
set ok [analyze $Nsteps $DtAnalysis]; # Actually perform analysis; returns ok = 0 if analysis was successful.
# If analysis was not successful change some parameters to achieve convergance
# Performance is slower inside this loop
# Time bounded analysis
if {$ok != 0} { ;
set ok 0 ; # Reset analysis
set controlTime [getTime] ; # Find time at which convergance is not achieved
while {$controlTime < $TmaxAnalysis && $ok == 0} { # While time is less than maximum time and analysis is ok...
set controlTime [getTime]
set ok [analyze 1 $DtAnalysis] # Set initial values
if {$ok != 0} { # If not converged...
puts "Trying Newton with Initial Tangent..."
test NormDispIncr $Tol 1000 0
algorithm Newton -initial
set ok [analyze 1 $DtAnalysis]
test $testTypeDynamic $TolDynamic $maxNumIterDynamic 0
algorithm $algorithmTypeDynamic
}
if {$ok != 0} { # If still not converged...
puts "Trying Broyden..."
algorithm Broyden 8
set ok [analyze 1 $DtAnalysis]
algorithm $algorithmTypeDynamic
}
if {$ok != 0} { # If still not converged...
puts "Trying Newton with line search..."
algorithm NewtonLineSearch 0.8
set ok [analyze 1 $DtAnalysis]
algorithm $algorithmTypeDynamic
}
}
}; # End analysis if ok does not equal zero (convergence is not achieved).


node_input.tcl #################################################################################################################

# NODES
# Tag X(mm) Y(mm) Z(mm)
node 1 0 0 0
node 2 0 1000 0

node_end1_input.tcl ##############################################################################################################

# NODES AT ITH MEMBER END

# Tag X(mm) Y(mm) Z(mm)
node 3 0 0 0

node_end2_input.tcl ###############################################################################################################

# NODES AT jTH MEMBER END

# Tag X(mm) Y(mm) Z(mm)
node 4 0 1000 0

node_fixity_input.tcl #################################################################################################################

# NODE FIXITY

# Tag DX DY DZ RX RY RZ
fix 1 1 1 1 1 1 1

section_input.tcl #############################################################################################################

# SECTION PROPERTIES

# Sec Tag E A Iz Iy G J
section Elastic 1 10000 18225 27679000 27679000 666.6667 46778000

material_input.tcl #############################################################################################################

# MATERIAL PROPERTIES
# SHEAR COMPONENT FOR ELASTIC SECTIONS
# Mat Tag K
uniaxialMaterial Elastic 1 9720000 ;

agg_section_input.tcl ##########################################################################################################

# AGGREGATED SECTIONS
# COMBINING FLEXURAL/AXIAL AND SHEAR COMPONENT FOR ELASTIC SECTIONS
# Sec Tag Mat Tag1 String1 Mat Tag2 String2 Sec Tag
section Aggregator 2 1 Vy 1 Vz -section 1 ;

mod_element_input.tcl ###########################################################################################################

# ELEMENTS
# Elements are modelled between end nodes to ensure appropriate end fixity is achieved.
# Tag iNode jNode Int Pnts/Mat Tag sec tag transTag
element forceBeamColumn 1 3 4 5 1 1

element_end1_constraint_input.tcl########################################################################################################

# ELEMENT ITH END FIXITY

# Master Tag Slave Tag CONSTRAINT 1 CONSTRAINT 2 CONSTRAINT 3 CONSTRAINT 4 CONSTRAINT 5 CONSTRAINT 6
equalDOF 1 3 1 1 1 1 1 1

element_end2_constraint_input.tcl ########################################################################################################

# ELEMENT JTH END FIXITY

# Master Tag Slave Tag CONSTRAINT 1 CONSTRAINT 2 CONSTRAINT 3 CONSTRAINT 4 CONSTRAINT 5 CONSTRAINT 6
equalDOF 2 4 1 1 1 1 1 1

lumped_mass_input.tcl ##################################################################################################################

# SEISMIC MASSES
# Lumped masses excluding member self-weight
# Node Mx My Mz MRx MRy MRz
mass 1 4.55625 4.55625 4.55625 0 0 0
mass 2 104.55625 4.55625 4.55625 0 0 0

loads_input.tcl #########################################################################################################################

# GRAVITY LOADS
# Under EQ loading combination G+?lQ
# Node FX FY FZ Mx My Mz
load 1 0 -144.69681 0 0 0 0
load 2 0 -144.69681 0 0 0 0

regions_input.tcl #######################################################################################################################

# REGIONS
# Used for computation of interstorey shear and moment.
# Tag Type Type Numbers...
region 1 -ele 1 ;

##################################################################################################################################


Any assisstance would be greatly appreciated.

Best regards,

Mike Newcombe

[brag006]

Your model is set-up wrong. What's the purpose of node 3 and 4? There are no elements between node 1 and 3 and none between 2 and 4. It's enough to have 2 nodes 1 and 2. I did that and the gravity analysis ran fine.

[fmk]

if using Transformation 3 needs to be fixed and no equalDOF. Penalty doesn't mind multiple constraints.

now as for your model, brag006 has a point 3 and 4 are rather pointless as the fully constrain every node to move the same as 1 and 2!

[mpnewcombe]

Thank you both for your response. The reason I have the other nodes is to simulate end fixity for more complex models... In this model I have made it a simple as possible. I realize they are pointless, but it should still run right?

Please explain your comment "if using Transformation 3 needs to be fixed and no equalDOF. Penalty doesn't mind multiple constraints." Sorry, I am new to Opensees and don't understand.

Cheers,

Mike Newcombe

[brag006]

When I took those extra nodes your code seems to run fine.

He is saying if you choose the constraint called Transformation then your node 3 needs to be fixed with no equalDOF. But if you use the Penalty constraint then not a problem

[mpnewcombe]

Ok... Definately using Transformation contraints.

I've discovered that programming the member end fixity for conversion program to transfer a 3D static model (from Space Gass) into an OS input file is not going be easy anyway. Space Gass uses local member axes for end fixity while equalDOF is a global constraint condition.

I can forget about the extra nodes and model end fixity with global constraints which are compatible with equalDOF.

So I have made progress and can model 2D frame (in X and Y plane) with "model BasicBuilder -ndm 3 -ndf 6". However when I try to run a simple 3D frame OS crashes and I don't get an error message. Any idea what may cause this?

I thought it may be the way I have used section aggregator to capture shear deformation: e.g.:

section_input.tcl########################################################################
# SECTION PROPERTIES

# Sec Tag E A Iz Iy G J
section Elastic 1 200000 9320 114000000 38800000 80000 586000
section Elastic 2 200000 5720 121000000 8100000 80000 161000

agg_section_input.tcl########################################################################
# AGGREGATED SECTIONS
# COMBINING FLEXURAL/AXIAL AND SHEAR COMPONENT FOR ELASTIC SECTIONS
# Sec Tag Mat Tag1 String1 Mat Tag2 String2 Sec Tag
section Aggregator 3 1 Vy 1 Vz -section 1 ;
section Aggregator 4 2 Vy 2 Vz -section 2 ;
########################################################################

But when I use just section 1 or 2, it still crashes.

Thoughts? I can sent files though if required.

Many thanks,

Mike Newcombe

[brag006]

Best approach would be to output something to screen after every few lines that way you know were the programme is getting before it fails. If it fails quickly before you see were just add the command 'after $numOfmiliSecs'. That will pause the program for the set number of miliseconds giving you enough time to identify at which line the code is failing.

[mpnewcombe]

Teach a man to fish

Thank you, this will be very useful.