Greetings everyone,
I am currently trying to do an acoustic modal analysis on a cavity with a mesh I created with gmsh. Here is the main code :
Cool imports
import dolfinx as dfx
from dolfinx.fem import Function, FunctionSpace, form
from dolfinx.fem.petsc import assemble_matrix
from slepc4py import SLEPc
from ufl import dx, grad, inner, TrialFunction, TestFunction
from dolfinx.io import XDMFFile, gmshio
from mpi4py import MPI
import numpy as np
Resolution parameters
N_ext = 20 #number of acoustic modes extracted
deg = 1 #degree of the polynomial resolution
Mesh creation
mesh, cell_tags, facet_tags = gmshio.read_from_msh(“Cavity_2D_round_gmsh.msh”, MPI.COMM_WORLD, 0, gdim = 2)
Air parameters
rho = 1.225 #air density
c = 341 #soundspeed
Test & trial functions
V = FunctionSpace(mesh, (“CG”, deg))
u = TrialFunction(V)
v = TestFunction(V)
k_form = form(inner(grad(u), grad(v))*dx)
m_form = form(inner(u,v)*dx)
K = assemble_matrix(k_form, ) #where we input the dirichlet bc if we had one in the bracket
M = assemble_matrix(m_form, )
K.assemble()
M.assemble()
Setup of the eigensolver
solver = SLEPc.EPS().create()
solver.setDimensions(N_ext) #number of acoustic modes extracted
solver.setProblemType(SLEPc.EPS.ProblemType.GHEP) #stetting up the type of problem to solve, here a General Hermitian Eigenvalue Problem
st = SLEPc.ST().create() #initializing spectral transformation module
st.setType(SLEPc.ST.Type.SINVERT) #setting up shift and invert algorithm
st.setShift(0.1) #setting up the shift
st.setFromOptions() #finalizing setting up from options
solver.setST(st) #assigning these settups to our solver
solver.setOperators(K,M) #adjusting to stiffness and mass matrixs
Solving the problem
solver.solve()
xr, xi = K.createVecs() #initializing real and complex part of our solution to vectors the same size as the stiffness matrix
tol, maxit = solver.getTolerances() #tolerance and maximum iteration
nconv = solver.getConverged() #get the number of converged modes
print(“Number of extracted modes : %i” %N_ext)
print(“Degree of the polynomial resolution method : %i” %deg)
print(“Number of iterations of the method : %i” %solver.getIterationNumber())
print(“Solution method : %s” % solver.getType())
print(“”)
print(“Stopping condition : tol = %4g, maxit = %d” %(tol, maxit))
eig_vect =
eig_freq =
if nconv > 0 :
for i in range(nconv) :
k = solver.getEigenpair(i, xr, xi) #put the i-th eigenpair and put it in the real and imaginary part
fn = np.sqrt(k.real)/(2*np.pi) * c #converting pulsation to frequency
eig_freq.append(fn)
print("%12f HZ" % fn)
vect = xr.getArray()
eig_vect.append(vect.copy())
Save the eigenmodes to visualize it in Paraview for exemple
if nconv > 0 :
for i in range(nconv) :
with XDMFFile(mesh.comm, “Mode_” + str(np.round(eig_freq[i])) + “_Hz.xdmf”, “w”) as xdmf:
u = Function(V)
u.vector[:] = eig_vect[i]
u.name = “u”
xdmf.write_mesh(mesh)
xdmf.write_function(u)
However, when I try to import it with gmshio I got this error :
Info : Reading 'Cavity_2D_round_gmsh.msh'...
Info : 29 entities
Info : 2515 nodes
Info : 4516 elements
Info : Done reading 'Cavity_2D_round_gmsh.msh'
Traceback (most recent call last):
File ~/anaconda3/envs/fenicsx-env/lib/python3.12/site-packages/spyder_kernels/py3compat.py:356 in compat_exec
exec(code, globals, locals)
File ~/Documents/Test1/Test modal analysis.py:46
mesh, cell_tags, facet_tags = gmshio.read_from_msh(“Cavity_2D_round_gmsh.msh”, MPI.COMM_WORLD, 0, gdim = 2)
File ~/anaconda3/envs/fenicsx-env/lib/python3.12/site-packages/dolfinx/io/gmshio.py:332 in read_from_msh
msh = model_to_mesh(gmsh.model, comm, rank, gdim=gdim, partitioner=partitioner)
File ~/anaconda3/envs/fenicsx-env/lib/python3.12/site-packages/dolfinx/io/gmshio.py:268 in model_to_mesh
local_entities, local_values = _cpp.io.distribute_entity_data(
MemoryError: std::bad_alloc
I tried only running the mesh import but all I obtained was a kernel restart.
Thanks for your help and have a nice day !