.. |fenics| replace:: FEniCS
.. _convert-cerfacs-mesh:
How to convert a CERFACS mesh?
==============================
We can create simple meshes inside |fenics|, such as rectangles (``dolfin.cpp.generation.RectangleMesh``) and 3D boxes (``dolfin.cpp.generation.BoxMesh``). Nevertheless, non-academic applications normally require more complex geometries. For such cases, you are advised to use your own favorite CAD tool and mesh generator and rely on mesh conversion tools (such as `pyhip `_, `yamio `_/`meshio `_ and/or `dolfin-convert `_) to convert it to the desired format (``dolfin`` ``.xml``).
To show the simplest way of doing it, let's start by creating a mesh using ``pyhip``:
.. code-block:: python
from pyhip.commands.operations import generate_mesh_2d_3d
from pyhip.commands.writers import write_hdf5
lower_corner = (0., 0.)
upper_corner = (1., 1.)
resolution = (10, 10)
convert2tri = False
extru_range = (0., 1.)
extru_axis = 'z'
extru_res = 10
generate_mesh_2d_3d(lower_corner, upper_corner, resolution,
convert2tri=convert2tri, extru_range=extru_range,
extru_axis=extru_axis, extru_res=extru_res)
write_hdf5(mesh_filename)
hip_exit();
.. image:: images/cube.png
:align: center
|
``yamio`` (a wrapper of ``meshio`` with added functionality) allows the conversion of the mesh to a format readable by ``dolfin``:
.. code-block:: python
import yamio
mesh = yamio.read(mesh_filename)
mesh.write(dolfin_filename, file_format='dolfin-yamio')
Currently, this command creates two files, one containing information about mesh points and connectivities, and another containing information about the boundaries (in the future we will probably merge both files).
The mesh can now be read in |fenics|:
.. code-block:: python
from dolfin.cpp.mesh import Mesh
from dolfin.cpp.io import XDMFFile
mesh = Mesh()
with XDMFFile(dolfin_filename) as file:
file.read(mesh)
Boundary information can be loaded in the following way:
.. code-block:: python
import h5py
from dolfin.cpp.mesh import MeshFunctionSizet
from dolfin.mesh.meshvaluecollection import MeshValueCollection
# define markers
mvc = MeshValueCollection("size_t", mesh)
with XDMFFile(dolfin_bnd_filename) as file:
file.read(mvc, "bnd_patches")
boundary_markers = MeshFunctionSizet(mesh, mvc)
# load patch names
with h5py.File(h5_dolfin_bnd_filename, 'r') as h5_file:
patch_labels = [name.decode('utf-8').strip() for name in h5_file['PatchLabels'][()]]
.. note::
`fenics_utils `_ contains a function that abstracts these last two blocks of code:
.. code-block:: python
from fenics_utils.mesh.external import load_dolfin_mesh
mesh, boundary_markers, patch_labels = load_dolfin_mesh(dolfin_filename)
.. note::
Check :ref:`this ` out to understand how boundary patch information can be used to apply boundary conditions.