Hi,
The eigenfunctions for a particle in a 1D infinite well are the usual sin functions of the form
U_1 = \sin(\pi x)
U_2 = \sin(2 \pi x)
U_3 = \sin(3 \pi x)
etc.
However, when I solve this problem numerically FEniCS will sometimes produce eigenfunctions with the wrong sign
The code below, for example, produces the eigenfunctions
U_1 = -\sin(\pi x)
U_2 = \sin(2 \pi x)
even though it reports positive eigenvalues. Is there a way to insist on eigenfunctions that produce positive eigenvalues?
from dolfin import *
mesh = IntervalMesh(200,-5,5)
V = FunctionSpace(mesh, 'Lagrange', 1)
def u0_boundary(x, on_boundary):
return on_boundary
bc = DirichletBC(V,Constant(0.0) , u0_boundary)
u = TrialFunction(V)
v = TestFunction(V)
a = (inner(grad(u), grad(v)) \
+ Constant('0.0')*u*v)*dx
m = u*v*dx
A = PETScMatrix()
M = PETScMatrix()
_ = PETScVector()
L = Constant(0.)*v*dx
assemble_system(a, L, bc, A_tensor=A, b_tensor=_)
assemble_system(m, L, A_tensor=M, b_tensor=_)
eigensolver = SLEPcEigenSolver(A,M)
eigensolver.parameters['spectrum'] = 'smallest magnitude'
eigensolver.parameters['tolerance'] = 1.e-15
neig = 2
eigensolver.solve(neig)
u = Function(V)
for i in range(neig):
r, c, rx, cx = eigensolver.get_eigenpair(i)
print 'eigenvalue:', r
u.vector()[:] = rx
plot(u,interactive=True)
Is there a way to insist only one form is calculated?