[m,n,x,fvec,fjac,ldfjac,xp,fvecp,mode,err]=chkder(m,n,x,fvec,fjac,ldfjac,xp,fvecp,mode,err);

function [m,n,x,fvec,fjac,ldfjac,xp,fvecp,mode,err]=chkder(m,n,x,fvec,fjac,ldfjac,xp,fvecp,mode,err);
%***BEGIN PROLOGUE  CHKDER
%***PURPOSE  Check the gradients of M nonlinear functions in N
%            variables, evaluated at a point X, for consistency
%            with the functions themselves.
%***LIBRARY   SLATEC
%***CATEGORY  F3, G4C
%***TYPE      SINGLE PRECISION (CHKDER-S, DCKDER-D)
%***KEYWORDS  GRADIENTS, JACOBIAN, MINPACK, NONLINEAR
%***AUTHOR  Hiebert, K. L. (SNLA)
%***DESCRIPTION
%
%   This subroutine is a companion routine to SNLS1,SNLS1E,SNSQ,and
%   SNSQE which may be used to check the calculation of the Jacobian.
%
%     subroutine CHKDER
%
%     This subroutine checks the gradients of M nonlinear functions
%     in N variables, evaluated at a point X, for consistency with
%     the functions themselves. The user must call CKDER twice,
%     first with MODE = 1 and then with MODE = 2.
%
%     MODE = 1. On input, X must contain the point of evaluation.
%               On output, XP is set to a neighboring point.
%
%     MODE = 2. On input, FVEC must contain the functions and the
%                         rows of FJAC must contain the gradients
%                         of the respective functions each evaluated
%                         at X, and FVECP must contain the functions
%                         evaluated at XP.
%               On output, ERR contains measures of correctness of
%                          the respective gradients.
%
%     The subroutine does not perform reliably if cancellation or
%     rounding errors cause a severe loss of significance in the
%     evaluation of a function. Therefore, none of the components
%     of X should be unusually small (in particular, zero) or any
%     other value which may cause loss of significance.
%
%     The SUBROUTINE statement is
%
%       subroutine CHKDER(M,N,X,FVEC,FJAC,LDFJAC,XP,FVECP,MODE,ERR)
%
%     where
%
%       M is a positive integer input variable set to the number
%         of functions.
%
%       N is a positive integer input variable set to the number
%         of variables.
%
%       X is an input array of length N.
%
%       FVEC is an array of length M. On input when MODE = 2,
%         FVEC must contain the functions evaluated at X.
%
%       FJAC is an M by N array. On input when MODE = 2,
%         the rows of FJAC must contain the gradients of
%         the respective functions evaluated at X.
%
%       LDFJAC is a positive integer input parameter not less than M
%         which specifies the leading dimension of the array FJAC.
%
%       XP is an array of length N. On output when MODE = 1,
%         XP is set to a neighboring point of X.
%
%       FVECP is an array of length M. On input when MODE = 2,
%         FVECP must contain the functions evaluated at XP.
%
%       MODE is an integer input variable set to 1 on the first call
%         and 2 on the second. Other values of MODE are equivalent
%         to MODE = 1.
%
%       ERR is an array of length M. On output when MODE = 2,
%         ERR contains measures of correctness of the respective
%         gradients. If there is no severe loss of significance,
%         then if ERR(I) is 1.0 the I-th gradient is correct,
%         while if ERR(I) is 0.0 the I-th gradient is incorrect.
%         For values of ERR between 0.0 and 1.0, the categorization
%         is less certain. In general, a value of ERR(I) greater
%         than 0.5 indicates that the I-th gradient is probably
%         correct, while a value of ERR(I) less than 0.5 indicates
%         that the I-th gradient is probably incorrect.
%
%***REFERENCES  M. J. D. Powell, A hybrid method for nonlinear equa-
%                 tions. In Numerical Methods for Nonlinear Algebraic
%                 Equations, P. Rabinowitz, Editor.  Gordon and Breach,
%                 1988.
%***ROUTINES CALLED  R1MACH
%***REVISION HISTORY  (YYMMDD)
%   800301  DATE WRITTEN
%   890531  Changed all specific intrinsics to generic.  (WRB)
%   890831  Modified array declarations.  (WRB)
%   890831  REVISION DATE from Version 3.2
%   891214  Prologue converted to Version 4.0 format.  (BAB)
%   900326  Removed duplicate information from DESCRIPTION section.
%           (WRB)
%   920501  Reformatted the REFERENCES section.  (WRB)
%***end PROLOGUE  CHKDER
persistent eps epsf epslog epsmch factor firstCall i j one temp zero ; if isempty(firstCall),firstCall=1;end; 

x_shape=size(x);x=reshape(x,1,[]);
fvec_shape=size(fvec);fvec=reshape(fvec,1,[]);
fjac_shape=size(fjac);fjac=reshape([fjac(:).',zeros(1,ceil(numel(fjac)./prod([ldfjac])).*prod([ldfjac])-numel(fjac))],ldfjac,[]);
xp_shape=size(xp);xp=reshape(xp,1,[]);
fvecp_shape=size(fvecp);fvecp=reshape(fvecp,1,[]);
err_shape=size(err);err=reshape(err,1,[]);
if isempty(i), i=0; end;
if isempty(j), j=0; end;
if isempty(eps), eps=0; end;
if isempty(epsf), epsf=0; end;
if isempty(epslog), epslog=0; end;
if isempty(epsmch), epsmch=0; end;
if isempty(factor), factor=0; end;
if isempty(one), one=0; end;
if isempty(temp), temp=0; end;
if isempty(zero), zero=0; end;
%
if firstCall,   factor =[1.0e2];  end;
if firstCall,  one =[1.0e0];  end;
if firstCall,  zero=[0.0e0];  end;
firstCall=0;
%***FIRST EXECUTABLE STATEMENT  CHKDER
[epsmch ]=r1mach(4);
%
eps = sqrt(epsmch);
%
if( mode==2 )
%
%        MODE = 2.
%
epsf = factor.*epsmch;
epslog = log10(eps);
for i = 1 : m;
err(i) = zero;
end; i = fix(m+1);
for j = 1 : n;
temp = abs(x(j));
if( temp==zero )
temp = one;
end;
for i = 1 : m;
err(i) = err(i) + temp.*fjac(i,j);
end; i = fix(m+1);
end; j = fix(n+1);
for i = 1 : m;
temp = one;
if( fvec(i)~=zero && fvecp(i)~=zero &&abs(fvecp(i)-fvec(i))>=epsf.*abs(fvec(i)) )
temp = eps.*abs((fvecp(i)-fvec(i))./eps-err(i))./(abs(fvec(i))+abs(fvecp(i)));
end;
err(i) = one;
if( temp>epsmch && temp<eps )
err(i) =(log10(temp)-epslog)./epslog;
end;
if( temp>=eps )
err(i) = zero;
end;
end; i = fix(m+1);
else;
%
%        MODE = 1.
%
for j = 1 : n;
temp = eps.*abs(x(j));
if( temp==zero )
temp = eps;
end;
xp(j) = x(j) + temp;
end; j = fix(n+1);
end;
%
%
%     LAST CARD OF SUBROUTINE CHKDER.
%
x_shape=zeros(x_shape);x_shape(:)=x(1:numel(x_shape));x=x_shape;
fvec_shape=zeros(fvec_shape);fvec_shape(:)=fvec(1:numel(fvec_shape));fvec=fvec_shape;
fjac_shape=zeros(fjac_shape);fjac_shape(:)=fjac(1:numel(fjac_shape));fjac=fjac_shape;
xp_shape=zeros(xp_shape);xp_shape(:)=xp(1:numel(xp_shape));xp=xp_shape;
fvecp_shape=zeros(fvecp_shape);fvecp_shape(:)=fvecp(1:numel(fvecp_shape));fvecp=fvecp_shape;
err_shape=zeros(err_shape);err_shape(:)=err(1:numel(err_shape));err=err_shape;
end
%DECK CHKPR4