Meshes and floating bodies¶
Importing a mesh with Meshmagick¶
To load an existing mesh file, use the following syntax:
import capytaine as cpt
mesh = cpt.load_mesh('path/to/mesh.dat', file_format='nemoh')
body = cpt.FloatingBody(mesh=mesh)
The above example uses Nemoh’s mesh format.
Thanks to Meshmagick, numerous other mesh format can be imported.
The file format can be given with the file_format
optional argument.
If no format is given, the code will try to infer it from the file extension:
mesh = cpt.load_mesh('path/to/mesh.msh') # gmsh file
The formats currently supported in reading are listed in the following table (adapted from the documentation of Meshmagick).
File extension |
Software |
Keywords |
Extra features |
---|---|---|---|
.mar |
NEMOH [1] |
nemoh, mar |
Symmetries |
.nem |
NEMOH [1] |
nemoh_mesh, nem |
|
.gdf |
WAMIT [2] |
wamit, gdf |
|
.inp |
DIODORE [3] |
diodore-inp, inp |
|
.DAT |
DIODORE [3] |
diodore-dat |
|
.pnl |
HAMS |
pnl, hams |
Symmetries |
.hst |
HYDROSTAR [4] |
hydrostar, hst |
Symmetries |
.nat |
natural, nat |
||
.msh |
GMSH 2 [5] |
gmsh, msh |
|
.rad |
RADIOSS |
rad, radioss |
|
.stl |
stl |
||
.vtu |
PARAVIEW [6] |
vtu |
|
.vtp |
PARAVIEW [6] |
vtp |
|
.vtk |
PARAVIEW [6] |
paraview-legacy, vtk |
|
.tec |
TECPLOT [7] |
tecplot, tec |
|
.med |
SALOME [8] |
med, salome |
Not all metadata is taken into account when reading the mesh file. For instance, the body symmetry is taken into account only for the .mar and .hst file formats. Feel free to open an issue on Github to suggest improvements.
Importing a mesh with Meshio¶
Mesh can also be imported using the meshio library. Unlike the Meshmagick mesh readers mentioned above, this library is not packaged with Capytaine and need to be installed independently:
pip install meshio
A meshio mesh object can be used directly to initialize a FloatingBody
:
import meshio
mesh = meshio.read("myfile.stl")
body = cpt.FloatingBody(mesh=mesh, dofs=...)
Alternatively, the meshio mesh object can converted to Capytaine’s mesh
format with the load_from_meshio()
function:
from capytaine.io.meshio import load_from_meshio
cpt_mesh = cpt.load_from_meshio(mesh, name="My floating body")
This features allows to use pygmsh to generate the mesh, since this library returns mesh in the same format as meshio. Below is an example of a mesh generation with pygmsh (which also needs to be installed independently):
import pygmsh
offset = 1e-2
T1 = 0.16
T2 = 0.37
r1 = 0.88
r2 = 0.35
with pygmsh.occ.Geometry() as geom:
cyl = geom.add_cylinder([0, 0, 0], [0, 0, -T1], r1)
cone = geom.add_cone([0, 0, -T1], [0, 0, -T2], r1, r2)
geom.translate(cyl, [0, 0, offset])
geom.translate(cone, [0, 0, offset])
geom.boolean_union([cyl, cone])
gmsh_mesh = geom.generate_mesh(dim=2)
body = cpt.FloatingBody(gmsh_mesh)
Display and animation¶
Use the show
method to display the mesh in 3D using VTK (if installed):
mesh.show()
Once a FloatingBody
with dofs has been defineds, the animate
method can be used to visualize a given motion of the body:
anim = body.animate(motion={"Heave": 0.1, "Surge": 0.1j}, loop_duration=1.0)
anim.run()
The above example will present an interactive animation of the linear combination of heave and surge.
Jupyter notebooks can also include a (non-interactive) video of the animation:
anim.embed_in_notebook(camera_position=(-1.0, -1.0, 1.0), resolution=(400, 300))
Geometric transformations¶
Several functions are available to transform existing bodies and meshes.
Most transformation methods exist in two versions:
one, named as a infinitive verb (translate, rotate, …), is an in-place transformation;
the other, named as a past participle (translated, rotated, …), is the same transformation but returning a new object.
In most cases, performance is not significant and the latter method should be preferred since it makes code slightly easier to debug.
Below is a list of most of the available methods. All of them can be applied to both meshes or to floating bodies, in which case the degrees of freedom will also be transformed:
# TRANSLATIONS
mesh.translated_x(10.0)
mesh.translated_y(10.0)
mesh.translated_z(10.0)
mesh.translated([10.0, 5.0, 2.0])
# Translation such that point_a would become equal to point_b
mesh.translated_point_to_point(point_a=[5, 6, 7], point_b=[4, 3, 2])
# ROTATIONS
mesh.rotated_x(3.14/5) # Rotation of pi/5 around the Ox axis
mesh.rotated_y(3.14/5) # Rotation of pi/5 around the Oy axis
mesh.rotated_z(3.14/5) # Rotation of pi/5 around the Oz axis
# Rotation of pi/5 around an arbitrary axis.
from capytaine import Axis
my_axis = Axis(vector=[1, 1, 1], point=[3, 4, 5])
mesh.rotated(axis=my_axis, angle=3.14/5)
# Rotation around a point such that vec1 would become equal to vec2
mesh.rotated_around_center_to_align_vector(
center=(0, 0, 0),
vec1=(1, 4, 7),
vec2=(9, 2, 1)
)
# REFLECTIONS
from capytaine import Plane
mesh.mirrored(Plane(normal=[1, 2, 1], point=[0, 4, 5]))
All the above method can also be applied to Plane
and Axis
objects.
Joining¶
Meshes and bodies can be merged together with the +
operator:
both_bodies = body_1 + body_2
The +
operation is associative, that is (body_1 + body_2) + body_3
is equivalent to body_1 + (body_2 + body_3)
.
It is also commutative, up to some internal details which are usually not relevant.
However for more than two bodies, it is recommended to use instead the
join_bodies
method:
all_bodies = body_1.join_bodies(body_2, body_3, body_4)
When two floating bodies with dofs are merged, the resulting body inherits from
the dofs of the individual bodies with the new name body_name__dof_name
.
For instance:
body_1.add_translation_dof(name="Heave")
body_2.add_translation_dof(name="Heave")
both_bodies = body_1 + body_2
assert 'body_1__Heave' in both_bodies.dofs
assert 'body_2__Heave' in both_bodies.dofs
Clipping¶
Meshes and bodies can be clipped with the clip
and clipped
methods.
As for the geometric transformations, the former is in-place whereas the second
returns a new object.
These methods take a Plane
object as argument. The plane is defined by a point belonging to it and a normal
vector:
xOy_Plane = Plane(point=(0, 0, 0), normal=(0, 0, 1))
clipped_body = body.clipped(xOy_Plane)
Beware that the orientation of the normal vector of the Plane
will
determine which part of the mesh will be returned:
higher_part = body.clipped(Plane(point=(0, 0, 0), normal=(0, 0, -1)))
lower_part = body.clipped(Plane(point=(0, 0, 0), normal=(0, 0, 1)))
# body = lower_part + higher_part
The method immersed_part
will clip the body with respect to two
horizontal planes at \(z=0\) and \(z=-h\):
clipped_body = body.immersed_part(water_depth=10)
Center of mass and rotation dofs¶
The center of gravity of the body can be defined by assigning a vector of 3
elements to the center_of_mass
attribute:
body.center_of_mass = np.array([0.0, -1.0, -1.0])
The center of mass is used in some hydrostatics computation.
It is not required for hydrodynamical coefficients, except for the definition of the rotation degrees of freedom.
When defining a rotation dof, the code looks for attributes called
rotation_center
, center_of_mass
or * geometric_center
(in that order),
and use them to define the rotation axis.
If none of them are define, the rotation is defined around the origin of the domain \((0, 0, 0)\).
Defining an integration quadrature¶
During the resolution of the BEM problem, the Green function has to be integrated on each panel of the mesh. Parts of the Green function (such as the \(1/r\) Rankine terms) are integrated using an exact analytical expression for the integral. Other parts of the Green function rely on numerical integration. By default, this numerical integration is done by taking the value at the center of the panel and multiplying by its area. For a more accurate intagration, an higher order quadrature can be defined.
To define a quadrature scheme for a mesh, run the following command:
body.mesh.compute_quadrature(method="Gauss-Legendre 2")
The quadrature data can then be accessed at:
body.mesh.quadrature_points
and will be used automatically when needed.
Warning
Transformations of the mesh (merging, clipping, …) may reset the quadrature. Compute it only on your final mesh.
Warning
Quadratures schemes have been designed with quadrilateral panels. They work on triangular panels, but might not be as optimal then.
Alternatively, the compute_quadrature()
method also accepts methods from the Quadpy package:
import quadpy
body.mesh.compute_quadrature(method=quadpy.c2.get_good_scheme(8))