Selecting/Importing Physical Groups (from GMSH .msh files) in ANSYS/mapdl

[fzanst]

Hello,
I have imported a GMSH .msh file using the Tutorial that can be found in the Issues section (or here: https://mapdldocs.pyansys.com/examples/extended_examples/ex_01-gmsh_example/ex_01-gmsh_example.html) for a structural analysis. Everything is working fine so far, I was just wondering if there was a way to use Physical groups from GMSH in ANSYS, as it would be less complicated to select surfaces/nodes that way, or if the Physical Groups get "lost" during the file import/conversion from .msh to .cdb?

Thanks

[germa89]

Hi @fzanst

At the moment that feature (support physical groups) is not available.

Having said that, I would probably create a CDB file for each group, and then import them one by one and assign them component selections at the same time you import them.
I think it should work, but I haven't tried it.

[fzanst]

Hello @germa89, I have tried to solve my issue but not been able to.
It seems the automatic offsetting of nodes mentioned in the command reference does not work, maybe I am missing something? ("If a set of data exists prior to the CDREAD operation, that data set is offset upward to allow the new data to fit without overlap"). Even using NOOFFSET,CLEAR prior to importing the second CDB file does not offset the nodes the way it should. The nodes from the first import are overridden by the nodes of the second import (up to the max number of nodes of the second import).

I have created a minimal example where I import two CDB files of meshed square planes (created in GMSH) that showcase my problem. Importing them using the code below results in this:

Here is a minimal example that shows my problem:
CDB-Files.zip

from ansys.mapdl.core import launch_mapdl
import os
import sys

os.chdir(os.path.dirname(sys.argv[0]))
mapdl = launch_mapdl()
mapdl.clear()
mapdl.prep7()

mapdl.clocal("11", 0, "", "0.15")
mapdl.csys(0)
mapdl.cdread("db", "square1.cdb")
mapdl.transfer("11","","all")
mapdl.cm("square1", "elem")
mapdl.save()

mapdl.nooffset("clear")
mapdl.cdread("db", "square2.cdb")
mapdl.save()

mapdl.eplot()
mapdl.finish()

[germa89]

Hi @fzanst

I dont know why this conversation didn't pop up in my notifications until now.

As mentioned, physical groups are not supported.

About the numbering issue, @mikerife do you have any idea to do it?

[mikerife]

Hi @fzanst & @germa89 let me look into this cdb/numbering issue and I'll get back asap.
Mike

[akaszynski]

Here's a minimum working example:

Create a Mesh Using gmsh

import pyvista as pv
from ansys.mapdl.reader import save_as_archive
from ansys.mapdl.core import launch_mapdl
import gmsh
import math
import sys

gmsh.initialize()
gmsh.clear()

# Copied from t1.py...
lc = 1e-2
gmsh.model.geo.addPoint(0, 0, 0, lc, 1)
gmsh.model.geo.addPoint(0.1, 0, 0, lc, 2)
gmsh.model.geo.addPoint(0.1, 0.3, 0, lc, 3)
gmsh.model.geo.addPoint(0, 0.3, 0, lc, 4)
gmsh.model.geo.addLine(1, 2, 1)
gmsh.model.geo.addLine(3, 2, 2)
gmsh.model.geo.addLine(3, 4, 3)
gmsh.model.geo.addLine(4, 1, 4)
gmsh.model.geo.addCurveLoop([4, 1, -2, 3], 1)
gmsh.model.geo.addPlaneSurface([1], 1)
ps = gmsh.model.addPhysicalGroup(2, [1])  # this is just the first area

# create second area
gmsh.model.geo.addPoint(0.3, 0, 0, lc, 5)
gmsh.model.geo.addPoint(0.4, 0, 0, lc, 6)
gmsh.model.geo.addPoint(0.4, 0.3, 0, lc, 7)
gmsh.model.geo.addPoint(0.3, 0.3, 0, lc, 8)
gmsh.model.geo.addLine(5, 6, 5)
gmsh.model.geo.addLine(7, 6, 6)
gmsh.model.geo.addLine(7, 8, 7)
gmsh.model.geo.addLine(8, 5, 8)
gmsh.model.geo.addCurveLoop([8, 5, -6, 7], 2)
gmsh.model.geo.addPlaneSurface([2], 2)
ps = gmsh.model.addPhysicalGroup(2, [2])  # this is just the first area

# must sync geo before meshing
gmsh.model.geo.synchronize()
gmsh.model.mesh.generate(3)

# write the gmsh to disk
filename = "t1.msh"
gmsh.write(filename)

# read in the mesh using pyvista. plot it for fun just to show the physical maps
mesh = pv.read(filename)
mesh.plot(scalars="gmsh:physical", categories=True, show_edges=True, cmap="blues")

Convert physical groups to Ansys Components

Now, here's the trick. Node components need to be bool arrays to for them to be saved to the archive:

# the archive reader writer must have node components stored as bool arrays
for group in range(1, 3):
    mesh.cell_data[f'GROUP_{group}'] = mesh.cell_data['gmsh:physical'] == group

Write to an archive using pymapdl_reader

archive_filename = "archive.cdb"
save_as_archive("archive.cdb", mesh)

# load in mesh in MAPDL and plot the elements
mapdl = launch_mapdl()
mapdl.clear()
mapdl.cdread("db", archive_filename)
mapdl.eplot()

Show components are indeed loaded into MAPDL

>>> mapdl.mesh.n_elem
1448
>>> mapdl.cmsel('S', 'GROUP_1', 'ELEM')
722
>>> mapdl.cmsel('S', 'GROUP_2', 'ELEM')
722
>>> mapdl.allsel()
>>> mapdl.mesh.n_elem
1448

If you need ansys-mapdl-reader to support physical groups, ping me and I'll add them.

I'm starting to realize how awesome gmsh is and how we should support it better.

18 days late, sorry about that.

[fzanst]

Thanks for the help, the physical group selection works perfectly this way!

Now the only problem remaining is the nodes getting overridden when importing a second CDB file, as written in my last comment. I tried importing both files using the Ansys APDL GUI, but it also doesn't work. So far the only thing that prevents the second import from overriding the first import seems to be using the mapdl.numoff command for both the nodes and the elements of the first import, but this is a temporary and not a very clean solution as the needed offset value differs with every file.

[akaszynski]

Here's a minimum working example:

Create a Mesh Using gmsh

import pyvista as pv
from ansys.mapdl.reader import save_as_archive
from ansys.mapdl.core import launch_mapdl
import gmsh
import math
import sys

gmsh.initialize()
gmsh.clear()

# Copied from t1.py...
lc = 1e-2
gmsh.model.geo.addPoint(0, 0, 0, lc, 1)
gmsh.model.geo.addPoint(0.1, 0, 0, lc, 2)
gmsh.model.geo.addPoint(0.1, 0.3, 0, lc, 3)
gmsh.model.geo.addPoint(0, 0.3, 0, lc, 4)
gmsh.model.geo.addLine(1, 2, 1)
gmsh.model.geo.addLine(3, 2, 2)
gmsh.model.geo.addLine(3, 4, 3)
gmsh.model.geo.addLine(4, 1, 4)
gmsh.model.geo.addCurveLoop([4, 1, -2, 3], 1)
gmsh.model.geo.addPlaneSurface([1], 1)
ps = gmsh.model.addPhysicalGroup(2, [1])  # this is just the first area

# create second area
gmsh.model.geo.addPoint(0.3, 0, 0, lc, 5)
gmsh.model.geo.addPoint(0.4, 0, 0, lc, 6)
gmsh.model.geo.addPoint(0.4, 0.3, 0, lc, 7)
gmsh.model.geo.addPoint(0.3, 0.3, 0, lc, 8)
gmsh.model.geo.addLine(5, 6, 5)
gmsh.model.geo.addLine(7, 6, 6)
gmsh.model.geo.addLine(7, 8, 7)
gmsh.model.geo.addLine(8, 5, 8)
gmsh.model.geo.addCurveLoop([8, 5, -6, 7], 2)
gmsh.model.geo.addPlaneSurface([2], 2)
ps = gmsh.model.addPhysicalGroup(2, [2])  # this is just the first area

# must sync geo before meshing
gmsh.model.geo.synchronize()
gmsh.model.mesh.generate(3)

# write the gmsh to disk
filename = "t1.msh"
gmsh.write(filename)

# read in the mesh using pyvista. plot it for fun just to show the physical maps
mesh = pv.read(filename)
mesh.plot(scalars="gmsh:physical", categories=True, show_edges=True, cmap="blues")

Convert physical groups to Ansys Components

Now, here's the trick. Node components need to be bool arrays to for them to be saved to the archive:

# the archive reader writer must have node components stored as bool arrays
for group in range(1, 3):
    mesh.cell_data[f'GROUP_{group}'] = mesh.cell_data['gmsh:physical'] == group

Write to an archive using pymapdl_reader

archive_filename = "archive.cdb"
save_as_archive("archive.cdb", mesh)

# load in mesh in MAPDL and plot the elements
mapdl = launch_mapdl()
mapdl.clear()
mapdl.cdread("db", archive_filename)
mapdl.eplot()

Show components are indeed loaded into MAPDL

>>> mapdl.mesh.n_elem
1448
>>> mapdl.cmsel('S', 'GROUP_1', 'ELEM')
722
>>> mapdl.cmsel('S', 'GROUP_2', 'ELEM')
722
>>> mapdl.allsel()
>>> mapdl.mesh.n_elem
1448

If you need ansys-mapdl-reader to support physical groups, ping me and I'll add them.

I'm starting to realize how awesome gmsh is and how we should support it better.

18 days late, sorry about that.

[mikerife]

Hi @fzanst @germa89 @akaszynski

fzanst sorry the MAPDL help is not clear on this. Usually the cdb file will contain the necessary NUMOFF commands to offset the numbering of the existing FEM so as to fit the FEM contained in the cdb into the model with numbering starting at 1. Looks like gmsh does not write that command to the cdb file. So you would need to add it manually. The NOOFFSET command can bypass NUMOFF. The Clear option does not result in an offset being created. The cdb files have element types, nodes and elements. Whether to offset nodes and elements is usually straight-forward. But element types, materials, real constants, etc may not be. Anyway this will import the two cdb files properly:

from ansys.mapdl.core import launch_mapdl
import os
import sys
path = os.getcwd()
mapdl = launch_mapdl(run_location = path, additional_switches = '-smp')

mapdl.clear()
mapdl.prep7()
mapdl.clocal("11", 0, "", "0.15")
mapdl.csys(0)
mapdl.cdread("db", "square1.cdb")
mapdl.transfer("11","","all")
mapdl.cm("square1", "elem")

n_max = mapdl.get_value('node',0,'num','max')
e_max = mapdl.get_value('elem',0,'num','max')

mapdl.cmsel('u','square1')
mapdl.numoff('node',n_max)
mapdl.numoff('elem',e_max)
mapdl.cdread("db", "square2.cdb")
mapdl.cm('square2','elem')
mapdl.allsel()

mapdl.cmsel('s','square1')
mapdl.cmsel('a','square2')

mapdl.show('png','rev')
mapdl.cmplot('all','elem')
mapdl.show()

mapdl.exit()

The plot should show this:

Mike

[germa89]

Amazing, thank you very very much @mikerife

[fzanst]

Thanks for the quick reply,
I have also been using an almost identical fix using the numoff function, but the problem using this function is that the nodes / elements of the first import are offset by the amount of nodes / elements of the first import. This means that if the second import has a larger number nodes / elements, the offset doesn't work as it should.

Example: Import one has 10 nodes and is offset by 10, meaning its node numbers are 11-20. If the second import has 15 nodes, the nodes of the second import will be 1-15, overriding the first few nodes of the first import (11-15). You can try it by switching square1 and square2 in the cdread functions in the script you posted, the small square(square2) will be completely overridden by the large one(square1).

As a workaround I am currently reading the number of nodes and elements of the second import from GMSH and using these as offset values, but I was/am looking for a cleaner solution, if one exists.

Franz

[mikerife]

Hi @fzanst I'd not consider this a workaround, rather just working with the APDL commands as they are. We need to remember that this functionality goes back decades and one assumption was that the cdb files were being generated by MAPDL which auto writes the correct NUMOFF to the cdb file when using CDWRITE. Anyway parsing through a file to find a specific value is a PITA in MAPDL so here is a Python version. I think I could do better with respect to memory usage on accessing the file to find the two specific lines. But it works.

The 4th (starting counting at 1) field of the NBLOCK and EBLOCK cdb commands are the maximum node and element number in the block. We can search for these lines; get the 4th field value; then use that as the appropriate numoff. Note that this is still not fully generic as (generally speaking) cdb files can contain N and EN commands to define nodes and elements. A better way would be to:

1 Read in first cdb
2 Save db (geom.db)
3 Clear MAPDL db
4 Read in second cdb
5 *get max node/element etc
6 Clear MAPDL db and resume geom.db
7 Read in second cdb with numoff using values from step 5

Anyway here is the updated PyMAPDL script:

from ansys.mapdl.core import launch_mapdl
import os
import sys
path = os.getcwd()
mapdl = launch_mapdl(run_location = path, additional_switches = '-smp')

mapdl.clear()
mapdl.prep7()
mapdl.clocal("11", 0, "", "0.15")
mapdl.csys(0)
mapdl.cdread("db", "square1.cdb")
mapdl.transfer("11","","all")
mapdl.cm("square1", "elem")

def lines_that_start_with(string, fp):
    return [line for line in fp if line.startswith(string)]

with open("square2.cdb", "r") as fp:
    nblock_line = lines_that_start_with('NBLOCK',fp)
with open("square2.cdb", "r") as fp:    
    eblock_line = lines_that_start_with('EBLOCK',fp)
    
nblock = ''.join(nblock_line)
eblock = ''.join(eblock_line)

n_max = int(nblock.split(',')[3])
e_max = int(eblock.split(',')[3])

mapdl.cmsel('u','square1')
mapdl.numoff('node',n_max)
mapdl.numoff('elem',e_max)
mapdl.cdread("db", "square2.cdb")
mapdl.cm('square2','elem')
mapdl.allsel()

mapdl.cmsel('s','square1')
mapdl.cmsel('a','square2')

mapdl.show('png','rev')
mapdl.cmplot('all','elem')
mapdl.show()

mapdl.exit()

Mike

[fzanst]

Thanks a lot Mike! This is now working perfectly, sorry for bothering so much!

Franz