Solving eigenvalue problem on a nanometer scale mesh.

Question

Hi,

Below is a toy script for solving an eigenvalue problem on a 2D circle. It fails due to the small radius of the mesh. I.e. if radius and mesh are made larger (e.g. 5.0 and 0.7), the script runs correctly.

Is there any way to get FEniCs to accept such small scales explicitly? Unit conversion would be easy enough in the toy model but the real thing is much more complicated.

Thanks

from dolfin import *

radius = 5e-9
grain = 0.7e-9
#define mesh and function space 
mesh = CircleMesh(Point(0,0),radius,grain)
V = FunctionSpace(mesh, 'Lagrange', 1)

#build essential boundary conditions
def u0_boundary(x, on_boundary):
        return on_boundary
bc = DirichletBC(V,Constant(0.0) , u0_boundary)

#define functions
u = TrialFunction(V)
v = TestFunction(V)

Pot = Expression('1.0')

#define problem
A = PETScMatrix() 
M = PETScMatrix()
_ = PETScVector()
L = Constant(1.)*v*dx
m = u*v*dx
assemble_system(m, L, A_tensor=M, b_tensor=_)

a = (inner(grad(u), grad(v)) \
             + Pot*u*v)*dx
assemble_system(a, L, bc, A_tensor=A, b_tensor=_)

#create eigensolver
eigensolver = SLEPcEigenSolver(A,M)
eigensolver.parameters['spectrum'] = 'smallest real'
eigensolver.parameters['tolerance'] = 1.e-15

#solve for eigenvalues
eigensolver.solve(5)

uh = Function(V)
for j in range(0,5):
    #extract next eigenpair
    r, c, rx, cx = eigensolver.get_eigenpair(j)
    print 'eigenvalue:', r

    #assign eigenvector to function
    uh.vector()[:] = rx

Answer 1

You should try to non-dimensionalize the problem analytically first so that your radius only varies on a much more computationally friendly radius interval (e.g. radius ranging from 0 to 1). This will enable you to solve (essentially) the same problem on a mesh that this not vulnerable to small number arithmetic. There is no unique way of non-dimensionalizing, but I highly recommend reading a simple paper on the matter, such as Yovanovich & Muzychka's paper "Solutions of Poisson Equations within Singly and Doubly Connected Prismatic Domains", which illustrates a simple procedure to non-dimensionalize a domain which you can easily adapt to your problem.