Ordinary Differential Equations

1
Ordinary Differential Equations
2
Outline

• Announcements
• Homework II Solutions on web
• Homework III due Wed. by 5, by e-mail
• Homework II
• Differential Equations
• Numerical solution of ODEs
• Matlabs ODE solvers
• Optimization problems
• Matlabs optimization routines

3
Homework II

• Emphasis was on writing functions
• P.2 Only way to get R into ll2xy is through
  function argument
• function x,yll2xy(lat,lon,ref,R)
• P.8 Most confusion stems from definition of the
  running mean. A few others were less clear on
  programming

4
Running Mean

• Running mean takes the mean of successive data
  blocks
• Ex n3, x is length 6

y
x












y is m-n1
5
Running Mean

• for j1m-n1 loop over length of y
• end

y(j)0 for k1n y(j)y(j)x(jk-1) end y(j)y(j)/n y(j)mean(x(jjn-1))
6
RunningMean

• Another option
• yzeros(m-n1,1)
• for k1n
• yyx(kkm-n)
• end
• yy/n
• Which is faster, looping over m or n?

7
Differential Equations

• Ordinary differential equations (ODEs) arise in
  almost every field
• ODEs describe a function y in terms of its
  derivatives
• The goal is to solve for y

8
Example Logistic Growth

• N(t) is the function we want (number of animals)

9
Numerical Solution to ODEs

• In general, only simple (linear) ODEs can be
  solved analytically
• Most interesting ODEs are nonlinear, must solve
  numerically
• The idea is to approximate the derivatives by
  subtraction

10
Euler Method
11
Euler Method

• Simplest ODE scheme, but not very good
• 1st order, explicit, multi-step solver
• General multi-step solvers

(weighted mean of f evaluated at lots of ts)
12
Runge-Kutta Methods

• Multi-step solvers--each N is computed from N at
  several times
• can store previous Ns, so only one evaluation of
  f/iteration
• Runge-Kutta Methods multiple evaluations of
  f/iteration

13
Matlabs ODE solvers

• Matlab has several ODE solvers
• ode23 and ode45 are standard RK solvers
• ode15s and ode23s are specialized for stiff
  problems
• several others, check help ode23 or book

14
Matlabs ODE solvers

• All solvers use the same syntax
• t,Node23(_at_odefile, t, N0, options, params )
• odefile is the name of a function that implements
  f
• function fodefile(t, N, params), f is a column
  vector
• t is either start time, end time or a vector of
  times where you want N
• N0 initial conditions
• options control how solver works (defaults or
  okay)
• params parameters to be passed to odefile

15
Parameters

• To accept parameters, your ode function must be
  polymorphic
• fodefile(t,x)
• fodefile(t,x,params)
• function fodefile(t,x,params)
• if(narginlt3)
• paramsdefaults
• end
• f

16
nargin

• Matlab functions can be polymorphic--the same
  function can be called with different numbers and
  types of arguments
• Example plot(y), plot(x,y), plot(x,y,rp)
• In a function, nargin and nargout are the number
  of inputs provided and outputs requested by the
  caller.
• Even more flexibility using vargin and vargout

17
Example Lorenz equations

• Simplified model of convection cells
• In this case, N is a vector x,y,z and f must
  return a vector x,y,z

18
Optimization

• For a function f(x), we might want to know
• at what value of x is f smallest
• largest?
• equal to some value?
• All of these are optimization problems
• Matlab has several functions to perform
  optimization

19
Optimization

• Simplest functions are
• xfzero(_at_fun, x0)Finds f(x)0 starting from x0
• xfminbnd(_at_fun, xL, xR) Finds minimum in
  interval xL,xR
• xfminsearch(_at_fun,x0) Finds minimum starting at
  x0
• All of these functions require you to write a
  function implementing f
• function ffun(x)computes f(x)

20
Optimization Toolbox

• More optimization techniques, but used in the
  same way
• Some of these functions can make use of gradient
  information
• in higher dimensions, gradient tells you which
  way to go, can speed up procedure

f
df/dx 0
21
Optimization Toolbox

• To return gradient information, your function
  must be polymorphic
• xfun(x) return value
• x,dxfun(x) return value and gradient
• function x,dxfun(x)
• if(nargoutgt1)
• dx
• end
• x