EM eigenmodes (vectorial)

Hello all,

I am writing here again as the community has proven to be so helpful for me.
I am trying to find the eigenvalues and -vectors of the waveequation for vectorial quantities,
but it seems like there are some details in how to properly define the problem that I do not understand.
Find code excerpts below:

    N=10
    mesh = UnitCubeMesh(N,N,N)
    V = VectorFunctionSpace(mesh, family='CG', degree=1, dim=3)
    bc = DirichletBC(V,[0, 0, 0],boundary)

    materials = MeshFunction("size_t", mesh, mesh.topology().dim(), 0)
    eps = Permittivity(materials, 2, 1, degree=0)
    subdomain_0 = Omega_0()
    subdomain_1 = Omega_1()
    subdomain_0.mark(materials, 0)
    subdomain_1.mark(materials, 1)

    space = helmholtz_wellposed(mesh, V, bc, 1E7, eps)

The method helmholtz_wellposed is slightly adjusted from one that is on the website:

def helmholtz_wellposed(mesh, V, bc, approxEig, eps):
    u, v = TrialFunction(V), TestFunction(V)
    a = inner(grad(u), grad(v))*dx
    L = inner(as_vector((0,0,0)),v)   (1)
    m = eps*u*v*dx                          (2)
    A, _ = assemble_system(a, L, bc)
    B = assemble(m)

I mainly have questions about the lines marked as (1) and (2). First of all, for some reason (1) throws an error (“ufl.log.UFLException: This integral is missing an integration domain.”) if I multiply with dx. Why?

Removing the dx then throws an error at line (2): “ufl.log.UFLException: Invalid ranks 1 and 1 in product.”. Same question: Where is this error coming from and how do I solve it?

Because the representation is reduced to a zero integral, the compiler does not get information about the mesh from v, and it has to be supplied through the measure, i.e. dx(domain=mesh).

As V is a VectorFunctionSpace you need to use eps * inner(u, v) * dx to reduce the m to a matrix.

Thank you very much. Now the problem is “just” a runtime problem:

RuntimeError: 

*** -------------------------------------------------------------------------
*** DOLFIN encountered an error. If you are not able to resolve this issue
*** using the information listed below, you can ask for help at
***
***     fenics-support@googlegroups.com
***
*** Remember to include the error message listed below and, if possible,
*** include a *minimal* running example to reproduce the error.
***
*** -------------------------------------------------------------------------
*** Error:   Unable to assemble system.
*** Reason:  expected a linear form for L.
*** Where:   This error was encountered inside SystemAssembler.cpp.
*** Process: 0
*** 
*** DOLFIN version: 2019.1.0
*** Git changeset:  74d7efe1e84d65e9433fd96c50f1d278fa3e3f3f
*** -------------------------------------------------------------------------

Do you guys have any ideas how to go further?

The problem is that your right-hand side has been reduced to 0*dx, and thus is no longer linear, and cannot be assembled with assemble_system.
To by-pass this, you can call:

    L = inner(Constant((0, 0, 0)), v)*dx(domain=mesh)

and your code will work.

Please note that for further enquires, please supply a minimal code that can reproduce the issue.

1 Like

Changing the line

L = inner(as_vector((0,0,0)),v)*dx(domain=mesh)

to

L = inner(as_vector((1,1,1)),v)*dx(domain=mesh)

resulted in correct assembly, but this is not the way to solve for eigenmodes I guess :frowning:
Has to have something to do with it compiling to 0