From help-octave-request at bevo dot che dot wisc dot edu  Fri Jan 12 09:55:00 2001
Subject: correction to leasqr algo. : in dfdp.m
From: courtois at ensia dot inra dot fr
To: help-octave at bevo dot che dot wisc dot edu
Date: Fri, 12 Jan 2001 16:52:38 +0100 (CET)

Hi

I found that using leasqr.m there is a problem when matrix x has more
columns than rows. You'll find an easy solution below. I replace length
by size(x,1) which return the number of rows and not the max (row,
cols).

BTW I found that those routines are both in matcompat and optim which is
potentially a problem.

Anyway it is great ! I use it often and love it.

Thanx to Paul and other contributors.

--------------------dfdp.m-----------------
function prt=dfdp(x,f,p,dp,func)
% numerical partial derivatives (Jacobian) df/dp for use with leasqr
% --------INPUT VARIABLES---------
% x=vec or matrix of indep var(used as arg to func) x=[x0 x1 ....]
% f=func(x,p) vector initialsed by user before each call to dfdp
% p= vec of current parameter values
% dp= fractional increment of p for numerical derivatives
%      dp(j)>0 central differences calculated
%      dp(j)<0 one sided differences calculated
%      dp(j)=0 sets corresponding partials to zero; i.e. holds p(j) fixed
% func=string naming the function (.m) file
%      e.g. to calc Jacobian for function expsum prt=dfdp(x,f,p,dp,'expsum')
%----------OUTPUT VARIABLES-------
% prt= Jacobian Matrix prt(i,j)=df(i)/dp(j)
%================================
n=length(p);      %dimensions
m=size(x,1); % Used to be : m=length(x);
ps=p; prt=zeros(m,n);del=zeros(n,1);       % initialise Jacobian to Zero
for j=1:n
      del(j)=dp(j) .*p(j);    %cal delx=fract(dp)*param value(p)
           if p(j)==0
           del(j)=dp(j);     %if param=0 delx=fraction
           end
      p(j)=ps(j) + del(j);
      if del(j)~=0, f1=feval(func,x,p);
           if dp(j) < 0, prt(:,j)=(f1-f)./del(j);
           else
           p(j)=ps(j)- del(j);
           prt(:,j)=(f1-feval(func,x,p))./(2 .*del(j));
           end
      end
      p(j)=ps(j);     %restore p(j)
end
return

-- 

Francis

 =========================================================================
 Dr. Francis COURTOIS               I  ENSIA - INRA
                                    I  1, avenue des Olympiades
 Food Process Control Lab.          I  91744 MASSY Cedex
 Food Engineering Department        I  FRANCE
 -----------------------------------I  Tel : (+33) (0) 1 69 93 51 29
 mailto:courtois at ensia dot inra dot fr      I  Fax : (+33) (0) 1 69 93 51 85
 http://ensia.inra.fr/~courtois     I  Cell: (+33) (0) 6 84 13 73 92
==========================================================================



-------------------------------------------------------------
Octave is freely available under the terms of the GNU GPL.

Octave's home on the web:  http://www.octave.org
How to fund new projects:  http://www.octave.org/funding.html
Subscription information:  http://www.octave.org/archive.html
-------------------------------------------------------------