Skip to content
This repository has been archived by the owner on Sep 9, 2024. It is now read-only.

Commit

Permalink
Fix indentation (again).
Browse files Browse the repository at this point in the history
  • Loading branch information
acroy committed Mar 26, 2015
1 parent 1857442 commit 0ae6410
Showing 1 changed file with 21 additions and 22 deletions.
43 changes: 21 additions & 22 deletions src/ODE.jl
Original file line number Diff line number Diff line change
Expand Up @@ -30,25 +30,25 @@ function hinit(F, x0, t0, p, reltol, abstol)
f0 = F(t0, x0)
d1 = norm(f0, Inf)/tau

if d0 < 1e-5 || d1 < 1e-5
h0 = 1e-6
else
h0 = 0.01*d0/d1
end
if d0 < 1e-5 || d1 < 1e-5
h0 = 1e-6
else
h0 = 0.01*d0/d1
end

# perform Euler step
x1 = x0 + h0*f0
f1 = F(t0 + h0, x1)

# estimate second derivative
# estimate second derivative
d2 = norm(f1 - f0, Inf)/(tau*h0)

if max(d1, d2) <= 1e-15
h1 = max(1e-6, 1e-3*h0)
else
pow = -(2. + log10(max(d1, d2)))/(p + 1.)
h1 = 10.^pow
end
if max(d1, d2) <= 1e-15
h1 = max(1e-6, 1e-3*h0)
else
pow = -(2. + log10(max(d1, d2)))/(p + 1.)
h1 = 10.^pow
end
return min(100*h0, h1), f0
end

Expand Down Expand Up @@ -231,7 +231,7 @@ function oderkf(F, x0, tspan, p, a, bs, bp; reltol = 1.0e-5, abstol = 1.0e-8,
pow = 1/p # use the higher order to estimate the next step size

c = sum(a, 2) # consistency condition
k = Array(typeof(x0), length(c)) # storage for intermediate steps
k = Array(typeof(x0), length(c)) # storage for intermediate steps

# Initialization
t = tspan[1]
Expand All @@ -240,10 +240,10 @@ function oderkf(F, x0, tspan, p, a, bs, bp; reltol = 1.0e-5, abstol = 1.0e-8,

h = initstep
if h == 0.
# initial guess at a step size
h, k[1] = hinit(F, x0, t, p, reltol, abstol)
else
k[1] = F(t,x0) # first stage
# initial guess at a step size
h, k[1] = hinit(F, x0, t, p, reltol, abstol)
else
k[1] = F(t,x0) # first stage
end
h = tdir*min(h, maxstep)

Expand Down Expand Up @@ -511,10 +511,10 @@ function ode23s(F, y0, tspan; reltol = 1.0e-5, abstol = 1.0e-8,

h = initstep
if h == 0.
# initial guess at a step size
h, F0 = hinit(F, y0, t, 3, reltol, abstol)
else
F0 = F(t,y)
# initial guess at a step size
h, F0 = hinit(F, y0, t, 3, reltol, abstol)
else
F0 = F(t,y)
end
h = tdir*min(h, maxstep)

Expand All @@ -524,7 +524,6 @@ function ode23s(F, y0, tspan; reltol = 1.0e-5, abstol = 1.0e-8,
yout = Array(typeof(y0), 1)
yout[1] = copy(y) # first output solution


J = jac(t,y) # get Jacobian of F wrt y

while abs(t) < abs(tfinal) && minstep < abs(h)
Expand Down

0 comments on commit 0ae6410

Please sign in to comment.