C+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C
C*************************
C**       LSQMQL        **
C*************************
C
C
      SUBROUTINE LSQMQL(M,N,X,ACC,H,CALFUN,FMIN,FVECC,IFLAG,W)
C
      IMPLICIT NONE
C
C Purpose: Minimise the sum of squares of 'm' given functions of 'n'
C          parameters
C
C          The function minimised is:
C
C                      M
C             F(x) = Sigma [f(j)(X(1),X(2),...X(N))]**2         (M>=N)
C                     j=1
C
C          The user must supply an initial approximation to the required
C          solution in X and a subroutine CALFUN to calculate the values
C          of the functions f(j), j=1,2,...M for any values of the parameters
C          in X. (see below for parameters of CALFUN routine)
C
C          This is an interface routine to the SLATEC public domain
C          routines using a modification of the Levenberg-Marquardt algorithm
C          That code is a combination of the MINPACK codes (Argonne) LMDER, 
C          LMDIF, and LMSTR. The routine used is a modified version of
C          the DNLS1 routine from the SLATEC library. The modification
C          allows the same version of CALFUN to be used for this refinement
C          as for the other Laue library refinement routines LSQMIN and
C          SSQMIN and removes options which cannot be used with such a
C          function.
C
C          This interface is designed to be called with very similar
C          parameters as those of the Laue libraries SSQMIN or LSQMIN
C          which perform similar functions. The user supplied subroutine
C          for calculating the required function values is the same in
C          each case.
C          
C
C Parameters:
C
C INTEGER               M (I)    The number of functions f(j) occuring in the 
C                                sum of squares.
C INTEGER               N (I)    The number of variables X(1) ... X(N)
C DOUBLE PRECISION      X (I/O)  An array of at least length N containing,
C                                on input, the user's estimate of the
C                                required parameters X(1) ... X(N)
C                                On output the refined parameter values
C                                are returned.
C DOUBLE PRECISION    ACC (I)    Tolerance in the Function value such that 
C                                failure to decrease by more than
C                                this amount in an iteration signals
C                                completness. (Or, more precisely, Used as 
C                                the value for FTOL and XTOL in call to the 
C                                SLATEC library routine DNLS1)
C DOUBLE PRECISION      H (I)    The step size to be used in calculating the
C                                partial derivatives (EPSFCN parameter in
C                                call to the SLATEC library routine DNLS1).
C EXTERNAL         CALFUN (I)    Subroutine for calculating the function
C                                values for a given set of parameters (see
C                                below)
C DOUBLE PRECISION   FMIN (O)    The returned minimum value sum of squares
C DOUBLE PRECISION  FVECC (O)    An array, dimensioned to at least M, returning
C                                the M values of f(j), j=1...M calculated using
C                                the refined parameter values.
C INTEGER           IFLAG (O)    Output flag =-1, error in parameters
 
C                                = 1  algorithm estimates that the relative 
C                                     error in the sum of squares is at most 
C                                     ACC.
C
C                                = 2  algorithm estimates that the relative 
C                                     error between X and the solution is at 
C                                     most ACC.
C
C                                = 3  conditions for IFLAG = 1 and IFLAG = 2 
C                                     both hold.
C
C                                = 4  FVECC is orthogonal to the columns of 
C                                     the Jacobian to machine precision.
C
C                                = 5  number of calls to CALFUN has reached 
C                                     200*(N+1).
C
C                                = 6  ACC is too small. No further reduction 
C                                     in the sum of squares is possible.
C
C                                = 7  ACC is too small.  No further 
C                                     improvement in the approximate solution 
C                                     X is possible.
C
C DOUBLE PRECISION      W (*)    An array for work space dimensioned to at
C                                least (M+6)*N + M.
C
C     The user written function CALFUN has the call
C
C                      CALFUN (M,N,FVECC,X)
C
C          Where the parameters have the same meaning as for LSQMQL; this
C          routine must not alter the values of M, N, or X.
C  
      INTEGER M, N, IFLAG
      DOUBLE PRECISION FVECC(M), X(N)
      DOUBLE PRECISION ACC, H, FMIN

      DOUBLE PRECISION W(N*N+4*N)
C
C====== FUNCTIONS:
C
      EXTERNAL CALFUN
C
C====== LOCALS:
C
      INTEGER INDEX
      INTEGER INDEX2
      INTEGER MODE
      INTEGER MAXFEV
      INTEGER NPRINT
      INTEGER IOPT
      INTEGER INFO
      INTEGER NFEV
      INTEGER NJEV
      INTEGER I
      DOUBLE PRECISION FTOL
      DOUBLE PRECISION XTOL
      DOUBLE PRECISION GTOL
      DOUBLE PRECISION EPSFCN
      DOUBLE PRECISION FACTOR
      INTEGER IXI, IPT, IPTT, IXIT, IXTEMP
      INTEGER IBIG, JFAIL
      DOUBLE PRECISION FP, FPTT, DEL, T
C
C-------------------------------------------------------------------------------
C
      IFLAG = 0
      MAXFEV = 2*100*(N + 1)
      IOPT = 1
      FTOL = ACC
      XTOL = ACC
      GTOL = 0.0D0
      FACTOR = 1.0D2
      MODE = 1
      NPRINT = 0
      INDEX = 5*N+M
      INDEX2 = INDEX + M*N
      CALL DNLS1_QQQ (CALFUN, IOPT, M, N, X, FVECC, W(INDEX+1), M,
     +           FTOL, XTOL, GTOL, MAXFEV, H, W(1), MODE,
     +           FACTOR, NPRINT, IFLAG, NFEV, NJEV,
     +           W(INDEX2+1), W(N+1), W(2*N+1), W(3*N+1),
     +           W(4*N+1), W(5*N+1))
      IF (IFLAG.EQ.0) IFLAG = -1
      FMIN = 0.0
      DO 10 I = 1, M
         FMIN = FMIN + FVECC(I)*FVECC(I)
 10   CONTINUE
      RETURN
      END
