Stupid question: how to do y' = - lambda y

[ethrelfall]

Dear community,

I want to solve $y'=-\lambda y$, subject to $y(0)=1$. Obviously the ans. is $y = \exp{(-\lambda x)}$.

I've written a script attempting to solve this; the first-order element version looks reasonable but higher orders always seem to contain artifacts i.e. spurious wiggles.

What have I done wrong? (Did search for matching.)

Script:

from firedrake import *

mesh = IntervalMesh(20,0,1)

V = FunctionSpace(mesh, "CG", 3)  # why don't higher orders work seem to properly?
u = TrialFunction(V)
v = TestFunction(V)
lam = Constant(1.0)
zee = Constant(0.0)

a = (grad(u)[0]*v + lam*u*v)*dx
L = zee*v*dx

g = Function(V)

bc = DirichletBC(V, 1.0, 1)

solve(a==L, g, bcs=bc)

File("first_order_ode.pvd").write(g)

[colinjcotter]

This is a hyperbolic system so in general a standard Galerkin method will lead to wiggles. However, I think that the problem here is treatment of the outflow boundary. You need to integrate by parts and drop the surface term at the outflow.

[ethrelfall]

Thanks for the reply Colin ... unfortunately I think the problem was the way I was visualizing the data (plot over line in Paraview) which seems to give the same artifacts even for functions I interpolated (i.e. not solved from ODE). Dumping the data into a csv gives something that looks OK so I think the script actually works fine ...

Sorry, guess it really was a stupid question!