Simple symbolical and numerical FE solver optimized for Politecnico di Milano's course of Space Structures. This model can:
- Compute stiffness matrix from non-constant beam stiffness
- Include concentrated/distributed loads, prescribed displacements and springs
- Compute reaction forces and strain energy
- Recover internal actions and local displacements along the beam
- Plot the model and its deformed version, along with relevant info (i.e, displacements, node number)
Figure 1. Example of result produced by the model (see example 1).
An example is reported hereafter, see the full version here.
%% Build structure
str = Structure([ ...
% Nodes
Node("uv"), ...
Node("v"), ...
Node(""), ...
Node("uv") ...
], [ ...
% Beams
Beam([1 2], l1, EA, 0, 330), ...
Beam([2 3], l1, EA, 0, 330), ...
Beam([3 4], l2, EA, 0, 180), ...
Beam([2 4], l1, EA, 0, 210) ...
] ...
);
- This software is not meant to replace the student's personal calculations;
- The code may be subject to bugs or yield wrong results;
- The author of the code assumes no responsibility in this sense;
- The code can be used at the risk and benefit of the individual;
It is highly recommended to often download the newly released versions as they may fix bugs or add new features.
You can download the latest release here.
Alternatively, you can import the project directly to MATLAB following these steps:
- From MATLAB, go to New > Project > From Git.
- In the Repository path field, paste the following link:
https://github.com/whitehole07/space-structure-fem.git
. - Finally, click on Retrieve.
You are all set.
This solver is based on three basic elements:
- Nodes
- Beams
- Springs
After defining them, these elements are combined together to form a compound element, the Structure.
To define a node the following syntax holds:
node = Node(constraint)
where:
- constraint, string cointaining the letters associated to the constrained degrees of freedom (order doesn't matter):
- u, if the horizontal displacement is prevented;
- v, if the vertical displacement is prevented;
- t, if the rotation is prevented.
For example:
node = Node("u") % only the horizontal movement is prevented
node = Node("uv") % the horizontal and vertical movements are prevented
...
node = Node("uvt") % all the movements are prevented
node = Node("") % none of the movements is prevented
For a multi-node structure, the nodes can be collected in an array, this will be later referred to as nodes array:
nodes = [Node("uvt"), Node(""), Node("uv"), ...]
To define a beam the following syntax holds:
beam = Beam(nodes, l, EA, EJ, alpha)
where:
- nodes, array cointaining the position in the nodes array of the beam's nodes, in ascending order (e.g. [1, 3]);
- l, beam length, can be symbolical;
- EA, beam axial stiffness function, can be symbolical;
- EJ, beam bending stiffness function, can be symbolical;
- alpha, counter-clockwise angle between the horizon and the beam, starting from the first node in nodes (see Figure 21).
Figure 2. Adopted convention.
Similarly, for a multi-beam structure, the beams can be collected in an array:
beams = [Beam([1 2], l1, EA, 0, 330), Beam([3 4], l2, EA, 0, 180), ...]
To define a spring the following syntax holds:
spring = Spring(node_num, dir, stiff)
where:
- node_num, associated node's index in the nodes array;
- dir, spring orientation:
- u, horizontal;
- v, vertical;
- t, rotational.
- stiff, spring stiffness, can be symbolical.
Springs can be collected in an array too, the code is omitted.
The structure can now be assembled, the following syntax holds:
str = Structure(nodes, beams, springs)
For example:
nodes = [Node(""), Node(""), Node("")]
beams = [Beam([1 2], l, EA1 + (EA2 - EA1)*(x/l), 0, 0), Beam([2 3], l, EA3, 0, 0)]
springs = [Spring(1, "u", k1), Spring(3, "u", k2)]
str = Structure(nodes, beams, springs)
% or
str = Structure([ ...
% Nodes
Node(""), ...
Node(""), ...
Node("") ...
], [ ...
% Beams
Beam([1 2], l, EA1 + (EA2 - EA1)*(x/l), 0, 0), ...
Beam([2 3], l, EA3, 0, 0), ...
], [ ...
% Springs
Spring(1, "u", k1), ...
Spring(3, "u", k2) ...
])
The methods available to apply loads to the structure are reported here.
str.add_concentrated_load(node_num, dir, load)
where:
- node_num, associated node's index in the nodes array;
- dir, load orientation:
- u, horizontal force;
- v, vertical force;
- t, torque.
- load, load value, can be symbolical.
str.add_distributed_load(beam_num, dir, load)
where:
- beam_num, associated beam's index in the beams array;
- dir, distributed load direction:
- axial, beam axial direction;
- transversal, beam transversal direction.
- load, distributed load function, can be symbolical.
To apply a prescribed displacement the following method is available:
str.add_prescribed_displacement(node_num, dir, displ)
where:
- node_num, associated node's position in the nodes array;
- dir, displacement direction:
- u, horizontal;
- v, vertical;
- t, rotation.
- displ, prescribed displacement, can be symbolical.
To solve the problem the following syntax can be used:
str.solve(var, val)
where:
- var, array of symbolic variables used during the definition of the problem;
- val, array of respective values (same order).
Once the problem is solved, all the involved quantities are stored inside the Structure object. The latter can be retrieved from MATLAB's Command Window, all the accessible properties are visible:
>> str
str =
Structure with properties:
nodes: [...] % Array of nodes
beams: [...] % Array of beams
springs: [...] % Array of springs
kf: [...] % Symbolic full stiffness matrix
uf: [...] % Symbolic full nodal displacements
ff: [...] % Symbolic full load array
Rf: [...] % Symbolic full reaction forces
kf_num: [...] % Numerical full stiffness matrix
uf_num: [...] % Numerical full nodal displacements
ff_num: [...] % Numerical full load array
Rf_num: [...] % Numerical full reaction forces
k: [...] % Symbolic reduced stiffness matrix
u: [...] % Symboluc solution for nodal displacements
f: [...] % Symbolic reduced load array
R: [...] % Symbolic reduced reaction forces
k_num: [...] % Numerical reduced stiffness matrix
u_num: [...] % Numerical solution for nodal displacements
f_num: [...] % Numerical reduced load array
R_num: [...] % Numerical reduced reaction forces
strain_energy: [...] % Numerical strain energy
var: [...] % Array of symbolic variables
val: [...] % Array of values
To get the numerical solution for the nodal displacements:
str.u_num
and similarly for the other properties.
Instead, single beam's stiffness matrix can be recovered from each beam object:
str.beams(beam_num).k_cont_red % Beam contribution to global problem
Once the problem is solved, different post processing methods to recover derived quantities are available.
To recover the displacements at any point of a beam, expressed in the beam's local reference frame, the following method can be used:
[u, v, t] = str.get_displ_local(beam_num, xi)
Instead, to recover the displacements expressed in the global reference frame:
[u, v, t] = str.get_displ_global(beam_num, xi)
where:
- [u, v, t], horizontal and vertical displacements and rotation;
- beam_num, associated beam's index in the beams array;
- xi, normalized position along the beam (xi = x/l), starting from the first node.
To recover the internal actions at any point of a beam the following method can be used:
[N, M] = str.get_internal_actions(beam_num, xi)
where:
- [N, M], axial force and bending moment;
- beam_num, associated beam's index in the beams array;
- xi, normalized position along the beam (xi = x/l), starting from the first node.
Full examples are collected in the example directory.
Footnotes
-
Credits to prof. Riccardo Vescovini, course of Space Structures at Politecnico di Milano. ↩