Reading CalculiX mesh with Python

This article explains how to make CalculiX nodes and elements accessible with Python.

Reading the mesh

When one uses the command send all abq in CalculiX GraphiX, this produces a file all.msh. This contains the definition of all the nodes and elements.

For starters, the file is read into Python and each line has extra whitespace on the left and right removed.

In [3]: with open("all.msh") as mf:
...:     lines = [ln.strip() for ln in mf]
...:

In [4]: len(lines)
Out[4]: 59901

This produces a list of lines. Generally, there is one set of nodes, and one set of elements. To process the lines, it needs to be clear where the nodes and elements begin.

In [6]: [n for n, ln in enumerate(lines) if '*ELEMENT' in ln]
Out[6]: [44340]

In [7]: [n for n, ln in enumerate(lines) if '*NODE' in ln]
Out[7]: [0]

In [8]: lines[0]
Out[8]: '*NODE, NSET=Nall'

In [9]: lines[44340]
Out[9]: '*ELEMENT, TYPE=C3D20R, ELSET=Eall'

This mesh contains one set of nodes and one set of C3D20R elements.

A node looks like this.

In [10]: lines[1]
Out[10]: '1,0.000000000000e+00,2.100000000000e-02,0.000000000000e+00'

Node 1 is located on coordinates x=0, y=0.021 and z=0. The nodes are converted to a dict of 3-tuples of floats as follows.

In [17]: nodedata = (ln.split(",") for ln in lines[1:44340])
Out[17]: <generator object <genexpr> at 0x87f095270>

In [18]: nodes = {int(n): (float(x), float(y), float(z)) for n, x, y, z in nodedata}

In [19]: len(nodes)
Out[19]: 44339

Note the use of a generator expression to simplify the processing without creating a long intermediate list.

The data for each C3D20(R) element is spread out over two lines.

In [21]: lines[44341:44343]
Out[21]:
['1,     1,     2,     3,     4,     5,     6,     7,     8,     9,    10,',
'11,    12,    17,    18,    19,    20,    13,    14,    15,    16']

Note how the first line ends with a comma and the second line doesn’t. For processing, it’s most convenient to combine the lines. The first number on the first line is the element number. The remaining numbers are the numbers of the nodes that make up the element.

In [25]: combined = ((a+b).split(",") for a, b in zip(lines[44341::2], lines[44342::2]))
Out[25]: <generator object <genexpr> at 0x89a19e970>

In [26]: elements = {int(items[0]): tuple(int(j) for j in items[1:]) for items in combined};

In [27]: len(elements)
Out[27]: 7780

In [28]: elements[1]
Out[28]: (1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 17, 18, 19, 20, 13, 14, 15, 16)

In [30]: (len(lines) - 44341)/2
Out[30]: 7780.0

Now both the nodes and elements are available in Python.

Reading a surface set

A surface set contains a list of numbered element faces. It starts with reading the file.

In [11]: with open("loadsur.sur") as sf:
    ...:     lines = [ln.strip() for ln in sf]
    ...:

In [12]: [n for n, ln in enumerate(lines) if '*SURFACE' in ln]
Out[12]: [1]

There is only one surface set in this case

Then the relevant lines are split and converted in a dictionary that maps element numbers to surface numbers.

In [16]: surfgen = (ln.split(",") for ln in lines[2:])
Out[16]: <generator object <genexpr> at 0x8a1cfdcf0>

In [17]: surface = {int(elnum): int(surfnum[2]) for elnum, surfnum in surfgen};

For comments, please send me an e-mail.

← From python script to executable with cython Profiling with pyinstrument →

Roland's homepage

Reading CalculiX mesh with Python

Reading the mesh

Reading a surface set

Related articles