|
|
|
The software which is named "fossils" here is based on a much more general finite element code built with the help of the [Gmsh SDK](https://gmsh.info/). This code is able to solve any kind of linear statics problems. In practise, the code can be used through a python module named `cxxfem`.
|
|
|
|
|
|
|
|
Two examples of mechanical problems are given in the `models/others` folder of the installed version of fossils. They correspond to:
|
|
|
|
* [beam2d.py](https://gitlab.uliege.be/rboman/fossils/-/blob/master/cxxfem/tests/beam2d.py): bending of a beam (2D plane stress),
|
|
|
|
* [beam3d.py](https://gitlab.uliege.be/rboman/fossils/-/blob/master/cxxfem/tests/beam3d.py): combined bending and torsion of a beam (3D).
|
|
|
|
|
|
|
|
Let us have a look at [beam2d.py](https://gitlab.uliege.be/rboman/fossils/-/blob/master/cxxfem/tests/beam2d.py).
|
|
|
|
|
|
|
|
The file starts with a function named `build_square` which creates the problem geometry using a series Gmsh calls. The syntax is described in the Gmsh manual. It is also possible to load a `.geo` file from Gmsh.
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
The result of this function is a rectangle with several geometrical entities named with explicit names ("top" is the top curve, "interior" is the interior surface of the beam, etc.). These names are also called "physical groups" in Gmsh. They are used to assign properties and boundary conditions to a list of geometrical entities.
|
|
|
|
|
|
|
|
Let us have a look at the `__main__` part of the script.
|
|
|
|
|
|
|
|
The first thing to do is to create a `Problem` object.
|
|
|
|
```
|
|
|
|
pbl = fem.Problem()
|
|
|
|
```
|
|
|
|
This is the main object which contains all the data related to the problem to be solved.
|
|
|
|
|
|
|
|
Then, the `build_square` function is called with particular parameters (beam size and mesh size)
|
|
|
|
```
|
|
|
|
build_square(0.0, 0.0, 10.0, 1.0, 80, 8, True, True)
|
|
|
|
```
|
|
|
|
The next step is to define a material. This is done by adding a `Medium` object to the `Problem`. The target physical group ("interior") is also given, as well as the values of the Young's modulus and the Poisson coefficient:
|
|
|
|
```
|
|
|
|
m = fem.Medium(pbl, "interior", 100.0, 0.3)
|
|
|
|
```
|
|
|
|
Next, we can define boundary conditions:
|
|
|
|
```
|
|
|
|
d1 = fem.Dirichlet(pbl, "left", "X", 0.0)
|
|
|
|
d2 = fem.Dirichlet(pbl, "left", "Y", 0.0)
|
|
|
|
n3 = fem.Neumann(pbl, "top", "Y", -1e-3)
|
|
|
|
d4 = fem.Dirichlet(pbl, "interior", "Z", 0.0)
|
|
|
|
```
|
|
|
|
Their type can be either `Dirichlet` (the value of the unknown is prescribed) or `Neumann` (the derivative of the unknown is prescribed).
|
|
|
|
|
|
|
|
The syntax is rather simple. For example, the first line means that we want to prescribe the value of the X component of the displacement to be 0 and the "left" side of the beam. Similarly, the Y component of the displacement is fixed to 0 too, so that the left side is fully clamped. The Neumann condition means that we want to apply a uniformly distributed load of -1e-3 along Y on the top side of the beam. Finally, the Z component of the dispacement field is fixed to 0 since it is a 2D problem.
|
|
|
|
|
|
|
|
At this stage, the `Problem` is fully defined. It should be solved by a `Solver` object.
|
|
|
|
|
|
|
|
The solver is defined by the following line. It takes the target problem as argument:
|
|
|
|
```
|
|
|
|
solver = fem.Solver(pbl)
|
|
|
|
```
|
|
|
|
|
|
|
|
The solver is then executed:
|
|
|
|
```
|
|
|
|
solver.solve()
|
|
|
|
```
|
|
|
|
This command starts the calculations. When they are finished, the results can be displayed and analysed. This part is managed by a `Post` object, which takes the target solver as argument:
|
|
|
|
```
|
|
|
|
post = fem.Post(solver)
|
|
|
|
```
|
|
|
|
The post-processing views (stress/strain fields calculation from the displacements previously computed) are built with the following command:
|
|
|
|
```
|
|
|
|
post.build_views()
|
|
|
|
```
|
|
|
|
The views that will be displayed later in Gmsh can be selected with:
|
|
|
|
```
|
|
|
|
post.show_views( [ "stress_tensor", "force_vector" ] )
|
|
|
|
```
|
|
|
|
It is also possible (but not mandatory) to write the results to disk in the current workspace folder:
|
|
|
|
```
|
|
|
|
post.write()
|
|
|
|
```
|
|
|
|
Since the FE code solves small-strain problems, the deformation of the mesh should be magnified in order to be visible. This is done with the following command:
|
|
|
|
```
|
|
|
|
post.deform(5)
|
|
|
|
```
|
|
|
|
The `Post` object is also able to extract results on given physical groups. For example, we can extract the values of the `force_vector` on the nodes lying of the `left` side of the beam:
|
|
|
|
```
|
|
|
|
forces = post.probe("force_vector", "left")
|
|
|
|
```
|
|
|
|
The values are returned in a python list such as [`v1x`, `v1y`, `v1z`, `v2x`, `v2y`, `v2z`, ..., `vnx`, `vny`, `vnZ`] where [`v1x`, `v1y`, `v1z`] are the values of the 3 components of the force vector at the first node of the physical group; [`v1x`, `v1y`, `v1z`] are the values related to the second node, and so on...
|
|
|
|
|
|
|
|
It is thus easy to calculate the resultant vector whith the following python commands:
|
|
|
|
```
|
|
|
|
print(f'Rx = {sum(forces[::3])} N')
|
|
|
|
print(f'Ry = {sum(forces[1::3])} N')
|
|
|
|
```
|
|
|
|
Similar commands can be used to extract the vertical displacement of the beam:
|
|
|
|
```
|
|
|
|
disp = post.probe("Y", "pt2")
|
|
|
|
print(f'tip displacement = {disp[0]} mm')
|
|
|
|
```
|
|
|
|
The last command of the script starts Gmsh and displays the results:
|
|
|
|
```
|
|
|
|
post.view()
|
|
|
|
```
|
|
|
|
|
|
|
|
 |
