A1. What does the software do in a finite element calculation? The example of beam structures

The purpose of this introduction is to show, in a simple example, what lies behind a static calculation carried out with a finite element beam software.

In the wiki, by misuse of language, the term "bar" is often used to refer to "beam" type elements and not exclusively to "bar" type elements in the strict sense (i.e. working only under normal stress).

A few reminders:

Computer calculation is based on a representation of the structure by a set of beams whose intersections are ~~nodes.~~nodes. Its purpose is to determine the displacements of the nodes of the structure, i.e. the displacements of the ends of the beams, under the applied loads.

Figure 1 - Sign convention - forces and displacements

The forces at the ends of the beam depend linearly on the displacements at these points. It is shown that there is a matrix relationship between forces and displacements for a beam connecting two nodes i and j, such that, for a given coordinate system, [K_ij][q_ij]=[F_ij] (Figures 1 and 2), allowing the components of the vectors [F_ij] to be expressed as a function of those of the vectors [q_ij]. [K_ij] is called the stiffness matrix of the beam, [F_ij] the force vector at the nodes, and [q_ij] the vector of nodal displacements.

Figure 2 - Stiffness matrix of a bi-clamped beam

The matrix [Kij*], associated with the global coordinate system, is deduced from the matrix [Kij] expressed in the local coordinate system of the beam by applying a transformation [Λ]T [Kij] [Λ] (figure 3), making the displacements compatible with those of the global coordinate system. [Kij*] = [Λ]T [Kij] [Λ] is then the stiffness matrix of the beam expressed in the global coordinate system.

Figure 3 – Rotation matrix

To calculate the displacements of all the nodes, the software:

"assembles" the matrices [Kij*] of each of the beams to form the overall stiffness matrix [K] of the structure (figure 4),
inverses [K] → [K]-1,
multiplies the matrix [K]-1 by the matrix of the external forces [F], previously brought back to the nodes (if they are applied to a beam other than at the ends), in the global coordinate system, to calculate the displacement vector [q].

Figure 4 - Principle of assembly

Finally, to calculate the internal forces of the structure, the software multiplies the stiffness matrix [Kij] of each beam by the vector of displacements [qij], in the local reference frame, at nodes i and j. Thus, the forces are calculated at the nodes and added, if necessary, to the reactions of the bi-clamped beam (i.e. the forces introduced in the global model, but here projected into the local reference frame of the beam) - see the second example below. The forces and displacements along the beams can easily be deduced from those at the nodes using the Strength of Materials formulae.

Application to a simple example:

Figure 5 - Calculated structure

This example is processed using a spreadsheet, the source file is provided in Part 3 of this guide - it can be downloaded by clicking on the link below. Source file in Excel format of the example above.

This spreadsheet covers two examples: the example developed below for a 2-beam structure (to remain easily readable in paper format) and a second example with 4 beams. These examples allow you to visualize and understand the sequence of tasks.

The required beam characteristics for a plane problem are the section S, inertia I, Young's modulus E, length L, and orientation α about the horizontal for instance.

E = 36 000 MPa for all beams.

The numerical values of the stiffness matrix for each beam are easily determined from the literal values in figure 2.

Beam 1: Stiffness matrix

Beam 2: Stiffness matrix

The transformation is applied to these two matrices to make them compatible with the global coordinate system.

Rotation matrix for α = 1.571 Rad:

And its transpose:

Rotation matrix for α = 0.197 Rad:

And its transpose:

Beam 1 - Result of the operation [K*] = [Λ]T [K] [Λ]:

Beam 2 - Result of the operation [K*] = [Λ]T [K] [Λ]:

All these unit matrices are then assembled, i.e. all the matrix values relating to the same nodal displacement are cumulated.

For our example, it means summing the values of the two matrices corresponding to the common nodes (in the general case, the beams do not always connect two nodes whose numbering follow each other as shown below):

The blocked displacement components can be removed from the matrix since we already know that their value is equal to zero:

This matrix is easily inverted using the spreadsheet:

Load 1: External forces applied on a node

Figure 6 - Load definition

The nodal displacements are calculated (in the global coordinate system):

The forces in the column (beam 1) can be determined after transforming the displacements to the local coordinate system of the beam:

We have a compression of 12.68 kN, a shear effort of 1000 kN, and a bending moment that changes sign along the beam (we check that 1000 kN * 8 m - 3654.89 = 4345.11 kN.m in the base of the column - ok).

Figure 7 - Diagrams N, V, M - example 1

Loading 2: External forces applied to a beam (i.e. other than at the nodes).

Figure 8 – Loading definition

In this case, it is necessary to calculate the support reactions of the bi-clamped beam under these forces beforehand, as they are the ones that have to be incorporated into the matrix of external loads (the software only knows the nodes!). A bi-clamped beam reaction table is enough. Of course, the sign conventions must be observed.

Application to the horizontal distributed load p = -1000 kN/m on the column (figure 8). In this case, it is known that the clamp moments are -pL²/12 and the shear forces at the restraints are +/-pL/2; therefore, with L=8 m and taking into account the sign convention:

Multiplying [K*]^-1 by these efforts provides the values of the displacements of all the degrees of freedom. It is then enough to multiply the stiffness matrix of beam 1 by these displacements to recover the forces at the nodes.

The efforts in the column (beam 1) can be calculated, after converting the displacements from the global coordinate system to the local coordinate system of the beam...

... and add them together with the forces of perfect clamping, with the appropriate sign convention:

We have a compression effort of 1298 kN, a shear effort that increases from 0 to 8000 kN (pL=1000 kN/m * 8 m = 8000 kN), and a bending moment at the base of the column which is much higher than the perfect clamping moment.

The same methodology applies to beam 2 to obtain the efforts at the nodes.

Figure 9 - Diagrams N, V, M - example 2

Cases where the beams present important differences in stiffness:

For example, the vertical beam is stiffened very strongly. We decide, keeping a length of 8 m, to increase its cross-section to 1m * 1012m (b*ht). The matrix [K*] is shown below. The difference between the larger and smaller values can be easily observed ... which can lead in some cases to numerical instabilities (the resolution leads to dividing the matrix terms by one another).

Figure 10 - Matrix with contrasting values