Complex Eigenvalue problem with field variable

Linear Algebra
1

Hello, I’m implementing a eigenvalue problem with complex numbers. So far and with help of:
Complex inner product,
Coupled linear system as a workaround to complex numbers,
Complex eigenvalue problem
I’ve implemented the sesquilinear formulation. Here the code:

def inner_c(u_r,u_i,v_r,v_i):          #inner product for sesquilinear form
    r_r = inner(u_r,v_r) + inner(u_i,v_i)
    r_i = inner(u_i,v_r) - inner(u_r,v_i)
    return r_r + r_i

mesh = UnitSquareMesh(10,10)
# mixed function space:
V = FunctionSpace(mesh, MixedElement([FiniteElement("CG",mesh.ufl_cell(),1),
                                      FiniteElement("CG",mesh.ufl_cell(),1)]))
u = TrialFunction(V)
v = TestFunction(V)

u_r, u_j = split(u)
v_r, v_j  = split(v)

a = inner_c(grad(u_r),grad(u_j),grad(v_r), grad(v_j))*dx
A = assemble(a)

Now, the stiffness matrix has to be multiplied by a fiel variable c(x)=e^{j\cdot w} which may also be complex. Here my definition of c:

tol=1E-14
w_0=0
w_1=1

w = Expression('x[0]<=0.5+tol ? w_0:w_1', degree=0, 
               tol=tol, w_0=w_0, w_1=w_1)
c_real = cos(w)
c_img = sin(w)

If c were only real, a would look like:

a = (c*inner_c(grad(u_r),grad(u_j),grad(v_r), grad(v_j)))*dx

but since c(x) is complex, it appears to me that the multiplication may not be as “easy”. My assumption is that In general the multiplication would look like this:
(c_r+j\cdot c_j) (a_r+j\cdot a_j)=c_r(a_r+j\cdot a_j)+c_j(j\cdot a_r-a_j)
\Rightarrow Re=c_r\cdot a_r-c_j\cdot a_j , \qquad Im=c_r\cdot a_j+c_j\cdot a_r

Since the Matrix A has mixed real and imaginary columns/rows the multiplication would have to be respectively for each column/row. Here is some code I’ve tried, where for each real column/row it computs the real part of the multiplication and for the imaginary columns/row the imaginary part:

A_arr = np.array(A.array())
 
        
        for row in range(dim_A):            # real part A
            if row%2 == 0:                      # for even rows
                for clm in range (int(dim_A/2)):
                    A_arr[row,2*clm] = c_real* A_arr[row,2*clm]-c_img*A_arr[row,2*clm+1]
            else :                                   # for odd rows
                for clm in range (int(dim_A/2)):
                    A_arr[row,2*clm+1] = c_real*A_arr[row,2*clm+1]-c_img*A_arr[row,2*clm]
                    
                    
        for row in range(dim_A):            #imaginary part A
            if row%2 == 0:                      # for even rows 
                for clm in range (int(dim_A/2)):
                    A_arr[row,2*clm+1] =c_real* A_arr[row,2*clm+1]+c_img*A_arr[row,2*clm]
            else :                                   # for odd rows
                for clm in range (int(dim_A/2)):
                    A_arr[row,2*clm] = c_real*A_arr[row,2*clm]+c_img*A_arr[row,2*clm+1]

The implementation is just an idea so far, but i would like to know if the idea of the multiplication is correct.
It would be also helpful if anyone had some literature recommendation for the sesquilinar formulation and its FEM implementation.