How to use FEniCS to handle a 3D Helmholtz equations

Question

I recently want to solve the three-dimensional Helmholtz equations with ABCs via the edge element method. But I am a newcomer with the program of FEniCS language, so I want to obtain some helps for achieving the program. Meanwhile, How to store the coefficient matrices A and right-hand side vector b of the discretized linear systems Au = b, i.e. coming from the FEM discretization of 3D Helmholtz equations by FEniCS? I already know that the discretized coefficient matrices A can be written as (nearly)

A = (K - k^2M + ikB); k is the wave-number.

K is the stiffness matrix, M is the mass matrix, and B is matrix from the discretized boundary conditions.

How to obtain and store the data of K, M, B, and b? and then I use them in MATLAB;

Thanks.

Answer 1

Hi, I don't know how B matrix is defined but the rest is shown below

from dolfin import *
from scipy.sparse import csr_matrix
from scipy.io import savemat

# Convert PETScMatrix to csr_matrix
to_csr = lambda matrix: csr_matrix(tuple(as_backend_type(matrix).mat().getValuesCSR()[::-1]),
                                   shape=(A.size(0), A.size(1)))

parameters['linear_algebra_backend'] = 'PETSc'

mesh = UnitSquareMesh(3, 3)
V = FunctionSpace(mesh, 'CG', 1)
u = TrialFunction(V)
v = TestFunction(V)
bc = DirichletBC(V, Constant(0.), 'on_boundary')

a = inner(grad(u), grad(v))*dx
m = inner(u, v)*dx
L = inner(Constant(1), v)*dx

A, b = assemble_system(a, L, bc)
M, _ = assemble_system(m, L, bc)

savemat('A.mat', dict(A=to_csr(A)))
savemat('M.mat', dict(M=to_csr(M))) 
savemat('b.mat', dict(b=b.array()))

Comments

@Mirok, many thanks for your answer. I am sorry that I fail to give the details of matrix B, now I give the description of matrix B via the following picture:

FEM for Helmholtz equations

or we refer to pp. 10-11, pp. 16-17

I want to solve the similar case. i.e., linear systems A u = f.

Then how to obtain the matrix B?

At last but not least, when I write the "from scipy.sparse import csr_matrix", but it will give the error in Command Shell. so what happen?

ImportError: No module named scipy.sparse

Computer: Windows 7, 64bit,
FEniCS 1.0.0.
Python 3.4