Topic: Nonlinear beams and condensing dofs.
I've tried to make beam2d and beam3d possibly support nonlinear behavior. The way condensing dofs works makes this very hard.
I'd like to replace dof condensing in beam2d and beam3d (the only place where it is used) with slave nodes. As far as I understand it, the goal is to obtain a moment free joint, which is easily achieved with a slave node that locks the displacements, placed on top of an existing node. I think this is much easier to understand what is going on, but mainly, it will enable nonlinear support.
---------------
Also, Beam2D and Beam3D use a hardcoded stiffness matrix, and doens't use the B-matrix ever.
I've tried to replace the computations for the internal forces by integrating in the form "B^T * D * B"
The manual claims it to be a Timoshenko beam, and the header file claims "cubic lateral displacement, quadratic rotations"
For example Beam2D. With B defined as
B.at(1, 1) = -1. / l;
B.at(1, 4) = 1. / l;
B.at(2, 2) = ( 6. - 12. * ksi ) / ( l * l * ( 1. + 2. * kappa ) );
B.at(2, 3) = ( -2. * ( 2. + kappa ) + 6. * ksi ) / ( l * ( 1. + 2. * kappa ) );
B.at(2, 5) = ( -6. + 12. * ksi ) / ( l * l * ( 1. + 2. * kappa ) );
B.at(2, 6) = ( -2. * ( 1. - kappa ) + 6. * ksi ) / ( l * ( 1. + 2. * kappa ) );
B.at(3, 2) = ( -2. * kappa ) / ( l * ( 1. + 2. * kappa ) );
B.at(3, 3) = kappa / ( l * ( 1. + 2. * kappa ) );
B.at(3, 5) = 2. * kappa / ( l * ( 1. + 2. * kappa ) );
B.at(3, 6) = kappa / ( l * ( 1. + 2. * kappa ) );
then integral B^T * D * B dV will not produce the same stiffness matrix as the hardcoded one.
Perhaps some trick to avoid shearlocking? Reverse-engineering this is very hard.
If kappa is close to zero, then this difference is barely noticeable.
For Beam3D i get some bigger trouble though. Again, if kappa is close to zero, then I mostly get the same matrix, except for a few components (the same value in a few places, changes sign):
Tangent = FloatMatrix with dimensions : 12 12
3.938e+06 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00 -3.938e+06 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00
0.000e+00 4.823e+05 0.000e+00 0.000e+00 0.000e+00 2.412e+05 0.000e+00 -4.823e+05 0.000e+00 0.000e+00 0.000e+00 2.412e+05
0.000e+00 0.000e+00 7.973e+05 0.000e+00 -3.987e+05 0.000e+00 0.000e+00 0.000e+00 -7.973e+05 0.000e+00 -3.987e+05 0.000e+00
0.000e+00 0.000e+00 0.000e+00 3.920e+04 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00 -3.920e+04 0.000e+00 0.000e+00
0.000e+00 0.000e+00 -3.987e+05 0.000e+00 2.658e+05 0.000e+00 0.000e+00 0.000e+00 3.987e+05 0.000e+00 1.329e+05 0.000e+00
0.000e+00 2.412e+05 0.000e+00 0.000e+00 0.000e+00 1.608e+05 0.000e+00 -2.412e+05 0.000e+00 0.000e+00 0.000e+00 8.039e+04
-3.938e+06 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00 3.938e+06 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00
0.000e+00 -4.823e+05 0.000e+00 0.000e+00 0.000e+00 -2.412e+05 0.000e+00 4.823e+05 0.000e+00 0.000e+00 0.000e+00 -2.412e+05
0.000e+00 0.000e+00 -7.973e+05 0.000e+00 3.987e+05 0.000e+00 0.000e+00 0.000e+00 7.973e+05 0.000e+00 3.987e+05 0.000e+00
0.000e+00 0.000e+00 0.000e+00 -3.920e+04 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00 3.920e+04 0.000e+00 0.000e+00
0.000e+00 0.000e+00 -3.987e+05 0.000e+00 1.329e+05 0.000e+00 0.000e+00 0.000e+00 3.987e+05 0.000e+00 2.658e+05 0.000e+00
0.000e+00 2.412e+05 0.000e+00 0.000e+00 0.000e+00 8.039e+04 0.000e+00 -2.412e+05 0.000e+00 0.000e+00 0.000e+00 1.608e+05
B^T * D * B = FloatMatrix with dimensions : 12 12
3.937e+06 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00 -3.937e+06 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00
0.000e+00 4.823e+05 0.000e+00 0.000e+00 0.000e+00 -2.412e+05 0.000e+00 -4.823e+05 0.000e+00 0.000e+00 0.000e+00 -2.412e+05 <----- Sign changed!
0.000e+00 0.000e+00 7.973e+05 0.000e+00 -3.987e+05 0.000e+00 0.000e+00 0.000e+00 -7.973e+05 0.000e+00 -3.987e+05 0.000e+00
0.000e+00 0.000e+00 0.000e+00 3.920e+04 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00 -3.920e+04 0.000e+00 0.000e+00
0.000e+00 0.000e+00 -3.987e+05 0.000e+00 2.658e+05 0.000e+00 0.000e+00 0.000e+00 3.987e+05 0.000e+00 1.329e+05 0.000e+00
0.000e+00 -2.412e+05 0.000e+00 0.000e+00 0.000e+00 1.608e+05 0.000e+00 2.412e+05 0.000e+00 0.000e+00 0.000e+00 8.039e+04
-3.937e+06 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00 3.937e+06 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00
0.000e+00 -4.823e+05 0.000e+00 0.000e+00 0.000e+00 2.412e+05 0.000e+00 4.823e+05 0.000e+00 0.000e+00 0.000e+00 2.412e+05
0.000e+00 0.000e+00 -7.973e+05 0.000e+00 3.987e+05 0.000e+00 0.000e+00 0.000e+00 7.973e+05 0.000e+00 3.987e+05 0.000e+00
0.000e+00 0.000e+00 0.000e+00 -3.920e+04 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00 3.920e+04 0.000e+00 0.000e+00
0.000e+00 0.000e+00 -3.987e+05 0.000e+00 1.329e+05 0.000e+00 0.000e+00 0.000e+00 3.987e+05 0.000e+00 2.658e+05 0.000e+00
0.000e+00 -2.412e+05 0.000e+00 0.000e+00 0.000e+00 8.039e+04 0.000e+00 2.412e+05 0.000e+00 0.000e+00 0.000e+00 1.608e+05
Which must mean a serious bug in either the B matrix or in the stiffness matrix.