Source code for pymech.neksuite.map
import struct
import numpy as np
from ..log import logger
[docs]def readma2(fname):
"""A function for reading binary map files (``*.ma2``) for nek5000.
The map file comtains, for each element in the mesh, the id of the MPI rank that owns it
followed by the ids of the vertices of the element in the global address space.
The partitioning is determined by generating the undirected graph formed by the mesh, then
repeatedly computing and deviding the graph using the Fiedler vector of the graph Laplacian.
Parameters
----------
fname : str
file name
Returns
-------
cell : 2d int ndarray
list of the vertices for each element of the mesh (in global address space)
procmap : 1d int ndarray
processor map (0-based) indicating ownership for each element of the mesh
"""
try:
infile = open(fname, "rb")
except OSError as e:
logger.critical(f"I/O error ({e.errno}): {e.strerror}")
return -1
#
# read header
header = infile.read(132).split()
nel = int(header[1])
# number of active elements (nrank - noutflow)
# nactive = int(header[2])
# total number of element vertices in the mesh
# npts = (2**ldim)*nel
npts = int(header[5])
# number of unique element vertices in the mesh
# nrank = int(header[6])
# number of points on outflow boundaries ('o ')
# noutflow = int(header[7])
# NOTE: these values can be computed from the others but are included in the
# header
# number of levels in the binary partition tree
# depth = log2(nel)
# depth = int(header[3])
# maximum number of elements in partition tree of depth d
# d2 = 2**d
# d2 = int(header[4])
# always double precision
wdsz = 4
inttype = "i"
# detect endianness
etagb = infile.read(4)
etagL = struct.unpack("<f", etagb)[0]
etagL = int(etagL * 1e5) / 1e5
etagB = struct.unpack(">f", etagb)[0]
etagB = int(etagB * 1e5) / 1e5
if etagL == 6.54321:
logger.debug("Reading little-endian file\n")
emode = "<"
# endian = "little"
elif etagB == 6.54321:
logger.debug("Reading big-endian file\n")
emode = ">"
# endian = "big"
else:
logger.error("Could not interpret endianness")
return -3
# read the entire contents of the file
# for each element, there are nvert vertices and a processor id
nvert = int(npts / nel) # 2**ldim
buf = infile.read((nvert + 1) * wdsz * nel)
# processor map (0-based)
procmap = np.empty((nel,))
# list of vertices for each element (in global address space)
cell = np.empty((nel, nvert))
for iel in range(nel):
fi = np.frombuffer(
buf,
dtype=emode + inttype,
count=nvert + 1,
offset=(nvert + 1) * wdsz * iel,
)
procmap[iel] = fi[0]
cell[iel, :] = fi[1:]
# close file
infile.close()
#
# output
return cell, procmap