# Achieving high relative accuracy (vs. absolute accuracy) using spectral methods

## Problem

I have a PDE that I'm trying to solve with spectral methods. The solution $y$ is always positive, and decays as $y \propto e^{-ax}$ for large $x$. The domain is $[0, \infty)$. (There are actually two independent variables but let's not worry about that right now.)

Unfortunately, although whatever numerical solution I obtain (using the Galerkin method) seems to have good absolute accuracy, it has poor relative accuracy towards larger $x$ due to the large dynamic range of the solution (coming from the exponential decay). I'm wondering how best to improve the relative accuracy of my solution. Any suggestions would be appreciated!

## Example test problem

This is not the actual PDE I'm trying to solve (which is actually 2d), but it illustrates the difficulty I'm having:

$$(1 + e^{-x})y' = -y + 2xe^{-x},\quad y(0) = 0.$$

The analytic solution is $x^2/(1+e^x)$ (plotted as dashed red below), which is always positive. However, solving this with a Laguerre function basis set (recombined so that each basis function is $0$ at $x = 0$) of order 20 yields the blue curve below, which diverges from the analytic solution (and even becomes negative) when it becomes small.

## What I've tried/considered

• I considered making the transformation $y = e^{-A(x)}$ and then solving for $A$. This would probably work, but unfortunately makes the equations nonlinear, which I was hoping not to deal with.
• I tried using Laguerre functions as a basis set, which have exponential decay built in. This didn't really work, though, as you can see above -- solutions still have low relative accuracy towards larger $x$.
• I'm planning to substitute $y = e^{-ax} z$ and then solve for $z$. I don't know $a$ a priori, but can estimate it. I think this will work all right, but I'm wondering if there's a more general method that I'm unaware of.
• Are you certain you don't have any errors in your code? I just solved this example problem with an Adaptive Runge-Kutta method and it solves it fine, so I don't know why the spectral approach would have difficulties unless there's implementation error. – spektr May 31 '16 at 4:48
• Thanks for looking into it! I used a fixed spectral order of 20, as I mentioned, so imperfect accuracy is expected. If I increase the order, the solution certainly becomes better. My issue is wanting to go for roughly uniform relative accuracy, rather than what seems to be uniform absolute accuracy. I see your point that switching to finite differences would probably help, but was looking specifically for advice on spectral methods. – Josh Burkart May 31 '16 at 6:45
• Your problem has a known exact solution. You can get arbitrarily good accuracy by using that. What are you really interested in solving? – David Ketcheson May 31 '16 at 11:39
• @JoshBurkart The FAQ says to please only post real problems you face, and for good reason. You're likely to get perfectly good answers (like "use the exact solution") that don't apply to your real problem. – David Ketcheson Jun 1 '16 at 15:43
• @DavidKetcheson I did post about a real problem I am facing. I then reduced it down to a demonstrative test problem that others could reproduce, which I labeled as "example test problem". I think this is a generally useful practice when requesting assistance, since my actual problem has lots of extra complexity that distracts from the issue I'm wondering about. Since there was confusion, I'll try to edit my post to make it clearer. – Josh Burkart Jun 1 '16 at 16:50

One thing I did do that you might not have was define the basis to exist on a domain that might have better numerical stability (especially as x increases). I made it so my basis was defined on the domain $\zeta \in [-1,1]$. Then I used the mapping between $\zeta$ and $x$ to modify the integrals and do the necessary computation needed to find the coefficients for the basis.