[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
User supplied jacobian makes ode23s very slow
From: |
Alexander Bessman |
Subject: |
User supplied jacobian makes ode23s very slow |
Date: |
Fri, 06 Nov 2015 13:01:39 +0100 |
Hello,
I'm using ode23s from odepkg to solve a simple 1D PDE, ∂u/∂t = ∂²u/dx².
I have discretized it via the forward-time central-space scheme,
yielding a tridiagonal system of ODEs, which ode23s has no problem
solving. However, when I supply the system's Jacobian via odeset, the
solution takes much longer to complete. The solution is still correct,
but I would have expected it to be faster, not slower. I've included my
code below. The result of running it is thus:
|Timesteps |CPU-time |Max stepsize |
|N |with Jac |w/o Jac |with Jac |w/o Jac |with Jac |w/o Jac |
|10 |432 |79 |1.3801 s |1.1953 s |1.9914e-2 |2.0000e-1|
|20 |1633 |98 |5.4225 s |2.6037 s |5.2711e-3 |1.9447e-1|
|40 |6442 |119 |26.7256 s |5.9159 s |8.4338e-4 |1.7272e-1|
I'd be grateful of someone could help me understand this phenomenon.
//Alexander Bessman
-----------------------------------------------------------------------
pkg load odepkg
cpu = zeros(2,3);
steps = zeros(2,3);
htmax = zeros(2,3);
iter = 0;
for N = [10, 20, 40]
% Jacobian
e = ones (N, 1);
jac = spdiags ([e, -2*e, e], [-1:1], N, N);
jac(end-1, end) = 2;
conf = odeset ('Jacobian', jac);
u0 = zeros(N,1);
iter = iter + 1;
tic
[t, y] = ode23s (@heat_eq, [0, 2], u0);
cpu(1, iter) = toc;
steps(1, iter) = length(t);
htmax(1, iter) = max(diff(t));
tic
[t, y] = ode23s (@heat_eq, [0, 2], u0, conf);
cpu(2, iter) = toc;
steps(2, iter) = length(t);
htmax(2, iter) = max(diff(t));
end
function dudt = heat_eq (t, u, kappa=1)
N = length(u);
e = ones (N, 1);
A = spdiags ([e, -2*e, e], [-1:1], N, N);
% Boundary conditions
A(end-1, end) = 2;
b = zeros(N, 1);
if t <= 1
b(1) = 1;
end
dudt = kappa * N^2 * (A*u + b);
- User supplied jacobian makes ode23s very slow,
Alexander Bessman <=