formulating an eigenvalue problem

Question

Hi every one!

Please consider the following weak formulation:
a(u,v)=omega*L(u,v)
where the omegas are the eigenvalues of the problem. My question is:
How can I convert the mentioned weak formulation to an eigenvalue problem in FEniCS and solve it? I read the simple eigenvalue problem demo provided by FeniCS but I don't know how to apply the method to my problem.

Thanks!

Answer 1

There's a brilliant demonstration of the Maxwell eigenvalue problem here

http://fenics.readthedocs.io/projects/dolfin/en/latest/demos/maxwell-eigenvalues/python/demo_maxwell-eigenvalues.py.html

Comments

Great! Thank you very much Nate.
Another question: How should I figure out whether to use assemble() or assemble_system()?
also what is dummy in the following part:

dummy = v[0]*dx
A = PETScMatrix()
assemble_system(a, dummy, bcs, A_tensor=A)

Thanks

---

assemble_system() will assemble A and b simultaneously with boundary conditions preserving symmetry. This means that you preserve the important property that A is Hermitian when solving the generalised Hermitian eigenvalue problem. The dummy argument is required just as a place filler for right hand side assembly. It's not used in the actual problem calculation.

Using assemble() and calling bc.apply() for each bc in bcs will alter the rows of A to satisfy the Dirichlet BCs, however, not the columns. This means that A will not be Hermitian.

---

Thanks Nate. I tried to follow the instruction in the demo. However, as my case contains complex numbers, I don't know how to use this part for my case:

dummy = v[0]*dx
A = PETScMatrix()
assemble_system(a, dummy, bcs, A_tensor=A)
B = PETScMatrix()
assemble_system(b, dummy, bcs, A_tensor=B)

Also I should note that in my problem there are only periodic boundary conditions (so there is no bcs argument). my code up to here is as follows:

from __future__ import print_function
from fenics import *
import numpy as np
import matplotlib.pyplot as plt

if not has_linear_algebra_backend("PETSc"):
    print("DOLFIN has not been configured with PETSc. Exiting.")
    exit()
if not has_slepc():
    print("DOLFIN has not been configured with SLEPc. Exiting.")
    exit()

###Parameters

k=as_vector([0.5,0.5,0.5]);
a=5.e-2; #cube size


####Material properties
E, rho, nu = 300E+8, 8000, 0.3
lam=E*nu/((1+nu)*(1-2*nu))
mu=E/(2*(1+nu))

# Sub domain for Periodic boundary condition
class PeriodicBoundary(SubDomain):

    # Left boundary is "target domain" G
    def inside(self, x, on_boundary):
        # return True if on left or bottom boundary AND NOT on one of the two slave edges
        return bool(( near(x[0], -a/2.) or near(x[1], -a/2.) or near(x[2], -a/2.)) and
            (not ((near(x[0], a/2.) and near(x[2], -a/2.)) or
                  (near(x[0], -a/2.) and near(x[2], a/2.)) or
                  (near(x[1], a/2.) and near(x[2], -a/2.))or
                  (near(x[1], -a/2.) and near(x[2], a/2.)))) and on_boundary)


# Map right boundary (H) to left boundary (G)
    def map(self, x, y):

        if near(x[0], a/2.) and near(x[2], a/2.):
            y[0] = x[0] - a/2
            y[1] = x[1]
            y[2] = x[2] - a/2
        elif near(x[1], a/2.) and near(x[2], a/2.):
            y[0] = x[0]
            y[1] = x[1] - a/2
            y[2] = x[2] - a/2
        elif near(x[0], a/2):
            y[0] = x[0] - a/2
            y[1] = x[1]
            y[2] = x[2]
        elif near(x[1], a/2):
            y[0] = x[0]
            y[1] = x[1] - a/2
            y[2] = x[2]
        elif near(x[2], a/2):
            y[0] = x[0]
            y[1] = x[1]
            y[2] = x[2] - a/2
        else:
            y[0] = -1000
            y[1] = -1000
            y[2] = -1000

###Mesh information
mesh = BoxMesh(Point(-a/2, -a/2, -a/2),Point(a/2, a/2, a/2),36,36,36)

V = VectorElement("Lagrange", mesh.ufl_cell(), 1)
W = MixedElement([V, V])
Vcomplex = FunctionSpace(mesh, W, constrained_domain=PeriodicBoundary())

###Variational Problem
uR, uI = TrialFunctions(Vcomplex)
wR, wI = TestFunctions(Vcomplex)

### Weak Form Fuctions
def F1(u,w):
    f1=-(inner((mu*nabla_grad(u)+mu*nabla_grad(u).T+lam*(div(u)*Identity(3))),nabla_grad(w)))
    return f1
def F2(u,w):
    f2= -dot(mu*(div(outer(u,k))+div(outer(k,u))+lam*(div(dot(k,u)*Identity(3)))),w)
    return f2
def B1(u,w):
    b=-rho*inner(u,w)
    return b
### Weak Form
F_Real=F1(uR,wR)*dx-F1(uI,wI)*dx-F2(uR,wI)*dx-F2(uI,wR)*dx
F_Imag=F2(uR,wR)*dx-F2(uI,wI)*dx+F1(uR,wI)*dx+F1(uI,wR)*dx
F_Total=F_Real+F_Imag
B_Real=B1(uR,wR)*dx-B1(uI,wI)*dx
B_Imag=B1(uR,wI)*dx+B1(uI,wR)*dx
B_Total=B_Real+B_Imag


### Assemble
dummy = wR[0]*dx
A = PETScMatrix()
assemble_system(F_Total, dummy, A_tensor=A)
B = PETScMatrix()
assemble_system(B_Total, dummy, A_tensor=B)

So I would be thankful if you could instruct me how to modify the last part of the code (###Assemble) part. When I run the current code it gives the following error:

    Calling FFC just-in-time (JIT) compiler, this may take some time.
Moving new file over differing existing file:
src: /home/amokhtar/Fenics/jitfailure-ffc_form_018a10df52f6df3cfd97841b43a64d71a0bc00fb/error.log.4de1fdcb562d4a01bd945ce2473761d6
dst: /home/amokhtar/Fenics/jitfailure-ffc_form_018a10df52f6df3cfd97841b43a64d71a0bc00fb/error.log
backup: /home/amokhtar/Fenics/jitfailure-ffc_form_018a10df52f6df3cfd97841b43a64d71a0bc00fb/error.log.old
Compilation failed! Sources, command, and errors have been written to: /home/amokhtar/Fenics/jitfailure-ffc_form_018a10df52f6df3cfd97841b43a64d71a0bc00fb
Traceback (most recent call last):
  File "Alessandro.py", line 108, in <module>
    assemble(F_Total, tensor=A)
  File "/share/apps/miniconda2/envs/fenics/lib/python2.7/site-packages/dolfin/fem/assembling.py", line 193, in assemble
    dolfin_form = _create_dolfin_form(form, form_compiler_parameters)
  File "/share/apps/miniconda2/envs/fenics/lib/python2.7/site-packages/dolfin/fem/assembling.py", line 67, in _create_dolfin_form
    function_spaces=function_spaces)
  File "/share/apps/miniconda2/envs/fenics/lib/python2.7/site-packages/dolfin/fem/form.py", line 89, in __init__
    mpi_comm=mesh.mpi_comm())
  File "/share/apps/miniconda2/envs/fenics/lib/python2.7/site-packages/dolfin/compilemodules/jit.py", line 70, in mpi_jit
    return local_jit(*args, **kwargs)
  File "/share/apps/miniconda2/envs/fenics/lib/python2.7/site-packages/dolfin/compilemodules/jit.py", line 147, in jit
    "ffc.jit failed with message:\n%s" % (tb_text,))
  File "/share/apps/miniconda2/envs/fenics/lib/python2.7/site-packages/dolfin/cpp/common.py", line 2342, in dolfin_error
    return _common.dolfin_error(location, task, reason)
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 perform just-in-time compilation of form.
*** Reason:  ffc.jit failed with message:
Traceback (most recent call last):
  File "/share/apps/miniconda2/envs/fenics/lib/python2.7/site-packages/dolfin/compilemodules/jit.py", line 142, in jit
    result = ffc.jit(ufl_object, parameters=p)
  File "/share/apps/miniconda2/envs/fenics/lib/python2.7/site-packages/ffc/jitcompiler.py", line 210, in jit
    raise FFCJitError("A directory with files to reproduce the jit build failure has been created.")
FFCJitError: A directory with files to reproduce the jit build failure has been created.
.
*** Where:   This error was encountered inside jit.py.
*** Process: 0
***
*** DOLFIN version: 2017.1.0
*** Git changeset:  d90efb69cac2db97c88c66a2189d43b3f91d3c80
*** -------------------------------------------------------------------------

Aborted