# Weak Formulation of Helmholtz equation with a complex coefficient and complex source term

I'm trying to solve

$$\nabla^2u + ku = f\text{ in } \Omega,$$ $$u=0 \text{ on }\partial\Omega,$$ where $\Omega$ is a unit square and $k,f$ have real and imaginary components as $k=k_r+k_ij$ and $f(x,y)=f_r(x,y) + j\cdot f_i(x,y)$, respectively.

My approach to solve this problem is to assume the solution has both a real and imaginary component as well: $u=u_r + ju_i$. Substituting this int the PDE above, we get:

$$\nabla^2u_r + j\nabla^2u_i+ (k_ru_r-k_iu_i)+j(k_iu_r+k_ru_i) = f_r+jf_i$$

Equating reals with reals and imaginaries with imaginaries, we obtain a system of two coupled PDE's:

$$\nabla^2u_r + (k_ru_r-k_iu_i) = f_r$$ $$\nabla^2u_i + (k_iu_r+k_ru_i) = f_i$$

Multipling each equation by a separate test function $v_r$, $v_i$, respectively, and integrating to get the weak forms, I get

$$-(\nabla u_r,\nabla v_r) + k_r(u_r,v_r)-k_i(u_i,v_r) = \int_\Omega f_rv_r$$ $$-(\nabla u_i,\nabla v_i) + k_i(u_r,v_i)+k_r(u_i,v_i) = \int_\Omega f_iv_i$$

Adding the two equations I get a single equation: $$-(\nabla u_r,\nabla v_r) + k_r(u_r,v_r)-k_i(u_i,v_r) -(\nabla u_i,\nabla v_i) + k_i(u_r,v_i)+k_r(u_i,v_i) = \int_\Omega f_rv_r+\int_\Omega f_iv_i$$

I then choose a finite dimensional space of piecewise polynomial functions V such that $u_r,v_r,u_i,v_i\in V$. That is, I choose all trial and test functions to come from the same space.

I've tried using Method of Manufactured Solutions to test my implementation in two cases. In case 1, I assume that the solution is real and of the form $u=sin(\pi x)sin(\pi y)$, which yields a source term $f = (-2\pi^2+2j)sin(\pi x)sin(\pi y)$. In case 2, I assume the solution has both real and imaginary parts and is of the form $u=sin(\pi x)sin(\pi y) + j sin(\pi x)sin(\pi y)$, which yields a source term $f= [(-2\pi^2-2) + j(-2\pi^2+2) ]sin(\pi x)sin(\pi y)$.

For case 1, the error decreases as I refine both the mesh size and the polynomial order, as expected. However, in case 2, I do not get convergence at all. Clearly, it could be a programming error, but I want to make sure that my approach is sound. In particular:

1. Is it ok to separate the complex equation into a system of real PDE's?

2. Am I wrong for multiplying real test functions into the separated equations? Should I have multiplied the complex equation by a complex test function of the form $v_r+j_vi$?

3. Is it ok to add the two weak forms together to come up with a single variational expression? I've seen this done when using mixed finite element spaces for the stokes equation, but I'm not entirely sure that this is ok for this problem.

Any help would be greatly appreciated.

Let $A$ and $B$ be real matrices, and $b$ and $c$ be real vectors. Then if we define $z$ as the complex solution to a complex system of equations, $$(A+jB)z = (b+jc)$$ then $z$ is equivalently given as $z=x+jy$, where $$\begin{bmatrix}A & -B \\ B & A\end{bmatrix}\begin{bmatrix}x \\ y\end{bmatrix} = \begin{bmatrix}b \\ c\end{bmatrix}$$ These simple facts extend to linear operator equations. That is, suppose that $A+jB$ were the complex Helmholtz differential operator that you are trying to invert. Then you can follow these identical steps to obtain two sets of real, coupled operator equations.
3. This is where you've made a mistake. To give a finite-dimensional analogy, you have effectively collapsed the 2x2 system shown above (with $2n$ real equations and $2n$ real unknowns) into $(A+B)x+(A-B)y = b+c$, which is an underdetermined system with $n$ real equations and $2n$ real unknowns. You then made this problem well-defined by setting $y=0$. So you're really solving $(A+B)x = (b+c)$, and this is right only when $B=0$ and $c=0$.