Bug in elemental loads

[TsarpalisD]

There is a bug in the way the elemental loads are taken into account for the construction of element's stiffness matrix. When calling getTangentStiff() the element's resisting forces in the natural system are calculated by q = kv + q0. Adding q0 to q leads to the negation of the elemental load, as it has been implictly taken into account in kv. As an example:

wipe;
model BasicBuilder -ndm 2 -ndf 3
node 1 0.0 0.0
node 2 0.0 5.0
fix 1 1 1 1
set E [expr 2.1e+8]
set A [expr 1.9e-3]
set Iz [expr 2.866e-6]
#geomTransf Corotational 1
#geomTransf Linear 1
geomTransf PDelta 1
element elasticBeamColumn 1 1 2 $A $E $Iz 1
pattern Plain 1 "Linear" {
eleLoad -ele 1 -type -beamUniform 0.0 -5.0
load 2 5.0 0.0 0.0
}
recorder Node -file displ.txt -time -node 2 -dof 1 disp
initialize
system BandGeneral
constraints Plain
numberer RCM
test NormUnbalance 1.0e-3 15 0
algorithm Newton
set nsteps 10
integrator LoadControl [expr 1./$nsteps]
analysis Static
analyze $nsteps

Using the script above the displacement of the top node is equal to 0.346149 while the correct value is 0.4373 as calculated by SAP2000 or OpenSees with equivalent nodal forces. This is due to the absence of p-delta effects as they are negated by q=kv+q0. Thus, keep in mind that elemental loads do not contribute to geometric nonlinearity of the structure.