search

connortannahill / bacoli_py

Python package for the error controlled numerical solution of 1D time-dependent PDEs
Other
3 stars 2 forks source link

Abort trap 6; segmentation faults exit python interpreter rather than raise exception #9

Open mrclary opened 4 years ago

mrclary commented 4 years ago

So, I have a minimal example (attached) that results in an Abort trap 6.

  1. Is it possible to raise an exception rather than kill the interpreter?
  2. Maybe you can help me understand why this fails in the first place?

I have a simple ode du/dt = nu * u(t, x) + source(x) with u(0, t) = u(1, t) = 0 and u(x, t) = 0

If I use an "analytic" expression for source(x), no problems. If I use an interpolating function for source(x), there-be-dragons (?!). Linear interpolation does not work and 2-order spline interpolation only works with fewer than 161 spatial points.

You can see in the definition of source(x) various commented options with comments describing which ones fail and which ones work. test_bacoli_py.py.zip

connortannahill commented 4 years ago

Hi Ryan,

I replicated your issue where a segmentation fault occurs when NX is to large. I don't know what's causing this exactly and will have to dig deeper when I have some time. One thing that seems to have fixed this issue is setting the t_int='r' when creating the Solver object.

I did get your code to run by changing a few things, I've attached it below, let me know if this works for you.

test_bacoli_py.py.zip

There are a couple of things to note about what your attempting to do, and reasons why the solver could fail:

u_t = f(x,t,u,u_x,u_xx),

where the dependence on u_xx is not optional. That is, the right hand side of the partial differential equation must include a u_xx dependence. I've added a very small u_xx term to the PDE in the example you provided, you can remove it if you wish, but again, this may result in occasional non-convergence.

mrclary commented 4 years ago

Thank you so much for your prompt response.

I replicated your issue where a segmentation fault occurs when NX is to large. I don't know what's causing this exactly and will have to dig deeper when I have some time. One thing that seems to have fixed this issue is setting the t_int='r' when creating the Solver object.

Good to know. Would changing nint_max or maxord also have an affect here? I haven't tried those yet.

I did get your code to run by changing a few things, I've attached it below, let me know if this works for you.

test_bacoli_py.py.zip

I'll take a look and play with it some more. Thanks.

There are a couple of things to note about what your attempting to do, and reasons why the solver could fail:

* The solver is not designed to solve an ODE. The differential equation must be of the form

u_t = f(x,t,u,u_x,u_xx),

where the dependence on u_xx is not optional. That is, the right hand side of the partial differential equation must include a u_xx dependence. I've added a very small u_xx term to the PDE in the example you provided, you can remove it if you wish, but again, this may result in occasional non-convergence.

Ahh, I see. This is good to know. This will be satisfied by my real problem. As I just started exploring bacoli_py, I was intentionally starting real simple just to become comfortable with the package and gradually add in more elements until I reach my final problem.

* The reason you may see non-convergence when using the linear/quadratic splines is likely because the spline function has discontinuous first/second derivatives. This causes massive issues for the underlying solver. To remedy this, I used Hermite cubic splines (setting kind='cubic') in your source terms and now both of the suggested source terms work for any NX I've tried.

Again, thank you, and good to know.

It looks like I had a convolution of three problems: 1. too small nint_max, 2. discontinuous 1st/2nd derivative spline, and 3. u_xx = 0