Methodology Suggestion for Wave-propagation Problem using Finite Elements

I want to simulate the propagation of a sinusoidal plane wave in a rectangular domain using Finite Elements Method. First, the wave should propagate through a fluid medium, then it will encounter a wall of solid material with some thickness, and then continue its propagation in the fluid. I am aware I have to consider some layer of absorbing material to avoid the spurious effects of wave reflection. What I found in the literature to be most suitable for this type of problem is the PML.

I can illustrate my objective with the following image

Ideally, I’d use an open-source Finite Element software. I have been trying to make this work with FEniCS, and I managed to get some good results in a time-domain simulation when there’s no PML. However, I couldn’t find a way to implement the PML. Many examples use a frequency-domain approach, which I’m not sure if can be used in this type of problem.

After some research I found a Thesis in which the author developed a PML mixed-formulation of time-domain equations suitable for nonlinear problems. The author also gave the weak form of these equations. But as far as I know, FEniCS demands the equations to be in the bilinear-linear form $a(u,v) = L(v)$, or $F(u,v)=0$ – no mixed-formulation. Also, the coding seems to be a bit complicated for a non-expert.

I’ve been trying to rewrite his equations into this different Weak Form ever since, but I haven’t a formal background in FE analysis, so up to now I’ve been unsuccessful. I can share my progress if you want me to.

The point is, I’m persuaded that there must be a simpler way to get this done. So I’d like to get from you some advice of methodology I should be following, or maybe point me to some meaningful reference that will help me achieve my goal.

Thanks in advance!

Edit

The general wave equation is $partial ^2 u / partial t ^2 = c^2 nabla ^2 u + f$.

The weak form is

$a(u,v) = int_{Omega} (vu^{n+1} + c^2 Delta t^2 nabla v nabla u^{n+1})d Omega$

$L_{n+1}(v) = int_{Omega}(2u^{n} – u^{n-1} + Delta t^2 f^{n+1})vd Omega$

where the exponent $n$ denotes the time step.

I have successfully implemented this into FEniCS by assembling a linear system with $a$, $L$ and the boundary conditions.

The next step is to implement the PML. I’ve seen lots of work using frequency-domain formulations, but as I mentioned, I’m not sure how to work with them. So I found this Thesis, where the author developed a set of mixed-formulation equations that account for the PML and are in time-domain. Consider the following for the fluid domain:

(eq. 5.9)

and its weak form

(eq. 5.10a)

(eq. 5.10c)

The problem is I was unable to translate those into FEniCS syntax. There are two equations and FEniCS can only solve one of the following:

• Linear systems $Ax=b$
• Linear variational problems $a=L$
• Nonlinear variational problems $F=0$

That means I must, somehow, rewrite the full set of governing equations into FEniCS syntax. I couldn’t find a way to this. I tried adding eqs. 5.10a and c together and making them equal 0, in order to use the nonlinear variational solver.

It didn’t work though. Not only the math became over complicated, but also the coding got messy. I felt I had to do a PhD in FEA just to know how code that (I’m a grad student). Instead I’m looking for something easier to understand and to code.