~ubuntu-branches/ubuntu/saucy/python-scipy/saucy

« back to all changes in this revision

Viewing changes to Lib/optimize/minpack/lmstr1.f

  • Committer: Bazaar Package Importer
  • Author(s): Ondrej Certik
  • Date: 2008-06-16 22:58:01 UTC
  • mfrom: (2.1.24 intrepid)
  • Revision ID: james.westby@ubuntu.com-20080616225801-irdhrpcwiocfbcmt
Tags: 0.6.0-12
* The description updated to match the current SciPy (Closes: #489149).
* Standards-Version bumped to 3.8.0 (no action needed)
* Build-Depends: netcdf-dev changed to libnetcdf-dev

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
 
      subroutine lmstr1(fcn,m,n,x,fvec,fjac,ldfjac,tol,info,ipvt,wa,
2
 
     *                  lwa)
3
 
      integer m,n,ldfjac,info,lwa
4
 
      integer ipvt(n)
5
 
      double precision tol
6
 
      double precision x(n),fvec(m),fjac(ldfjac,n),wa(lwa)
7
 
      external fcn
8
 
c     **********
9
 
c
10
 
c     subroutine lmstr1
11
 
c
12
 
c     the purpose of lmstr1 is to minimize the sum of the squares of
13
 
c     m nonlinear functions in n variables by a modification of
14
 
c     the levenberg-marquardt algorithm which uses minimal storage.
15
 
c     this is done by using the more general least-squares solver
16
 
c     lmstr. the user must provide a subroutine which calculates
17
 
c     the functions and the rows of the jacobian.
18
 
c
19
 
c     the subroutine statement is
20
 
c
21
 
c       subroutine lmstr1(fcn,m,n,x,fvec,fjac,ldfjac,tol,info,
22
 
c                         ipvt,wa,lwa)
23
 
c
24
 
c     where
25
 
c
26
 
c       fcn is the name of the user-supplied subroutine which
27
 
c         calculates the functions and the rows of the jacobian.
28
 
c         fcn must be declared in an external statement in the
29
 
c         user calling program, and should be written as follows.
30
 
c
31
 
c         subroutine fcn(m,n,x,fvec,fjrow,iflag)
32
 
c         integer m,n,iflag
33
 
c         double precision x(n),fvec(m),fjrow(n)
34
 
c         ----------
35
 
c         if iflag = 1 calculate the functions at x and
36
 
c         return this vector in fvec.
37
 
c         if iflag = i calculate the (i-1)-st row of the
38
 
c         jacobian at x and return this vector in fjrow.
39
 
c         ----------
40
 
c         return
41
 
c         end
42
 
c
43
 
c         the value of iflag should not be changed by fcn unless
44
 
c         the user wants to terminate execution of lmstr1.
45
 
c         in this case set iflag to a negative integer.
46
 
c
47
 
c       m is a positive integer input variable set to the number
48
 
c         of functions.
49
 
c
50
 
c       n is a positive integer input variable set to the number
51
 
c         of variables. n must not exceed m.
52
 
c
53
 
c       x is an array of length n. on input x must contain
54
 
c         an initial estimate of the solution vector. on output x
55
 
c         contains the final estimate of the solution vector.
56
 
c
57
 
c       fvec is an output array of length m which contains
58
 
c         the functions evaluated at the output x.
59
 
c
60
 
c       fjac is an output n by n array. the upper triangle of fjac
61
 
c         contains an upper triangular matrix r such that
62
 
c
63
 
c                t     t           t
64
 
c               p *(jac *jac)*p = r *r,
65
 
c
66
 
c         where p is a permutation matrix and jac is the final
67
 
c         calculated jacobian. column j of p is column ipvt(j)
68
 
c         (see below) of the identity matrix. the lower triangular
69
 
c         part of fjac contains information generated during
70
 
c         the computation of r.
71
 
c
72
 
c       ldfjac is a positive integer input variable not less than n
73
 
c         which specifies the leading dimension of the array fjac.
74
 
c
75
 
c       tol is a nonnegative input variable. termination occurs
76
 
c         when the algorithm estimates either that the relative
77
 
c         error in the sum of squares is at most tol or that
78
 
c         the relative error between x and the solution is at
79
 
c         most tol.
80
 
c
81
 
c       info is an integer output variable. if the user has
82
 
c         terminated execution, info is set to the (negative)
83
 
c         value of iflag. see description of fcn. otherwise,
84
 
c         info is set as follows.
85
 
c
86
 
c         info = 0  improper input parameters.
87
 
c
88
 
c         info = 1  algorithm estimates that the relative error
89
 
c                   in the sum of squares is at most tol.
90
 
c
91
 
c         info = 2  algorithm estimates that the relative error
92
 
c                   between x and the solution is at most tol.
93
 
c
94
 
c         info = 3  conditions for info = 1 and info = 2 both hold.
95
 
c
96
 
c         info = 4  fvec is orthogonal to the columns of the
97
 
c                   jacobian to machine precision.
98
 
c
99
 
c         info = 5  number of calls to fcn with iflag = 1 has
100
 
c                   reached 100*(n+1).
101
 
c
102
 
c         info = 6  tol is too small. no further reduction in
103
 
c                   the sum of squares is possible.
104
 
c
105
 
c         info = 7  tol is too small. no further improvement in
106
 
c                   the approximate solution x is possible.
107
 
c
108
 
c       ipvt is an integer output array of length n. ipvt
109
 
c         defines a permutation matrix p such that jac*p = q*r,
110
 
c         where jac is the final calculated jacobian, q is
111
 
c         orthogonal (not stored), and r is upper triangular.
112
 
c         column j of p is column ipvt(j) of the identity matrix.
113
 
c
114
 
c       wa is a work array of length lwa.
115
 
c
116
 
c       lwa is a positive integer input variable not less than 5*n+m.
117
 
c
118
 
c     subprograms called
119
 
c
120
 
c       user-supplied ...... fcn
121
 
c
122
 
c       minpack-supplied ... lmstr
123
 
c
124
 
c     argonne national laboratory. minpack project. march 1980.
125
 
c     burton s. garbow, dudley v. goetschel, kenneth e. hillstrom,
126
 
c     jorge j. more
127
 
c
128
 
c     **********
129
 
      integer maxfev,mode,nfev,njev,nprint
130
 
      double precision factor,ftol,gtol,xtol,zero
131
 
      data factor,zero /1.0d2,0.0d0/
132
 
      info = 0
133
 
c
134
 
c     check the input parameters for errors.
135
 
c
136
 
      if (n .le. 0 .or. m .lt. n .or. ldfjac .lt. n .or. tol .lt. zero
137
 
     *    .or. lwa .lt. 5*n + m) go to 10
138
 
c
139
 
c     call lmstr.
140
 
c
141
 
      maxfev = 100*(n + 1)
142
 
      ftol = tol
143
 
      xtol = tol
144
 
      gtol = zero
145
 
      mode = 1
146
 
      nprint = 0
147
 
      call lmstr(fcn,m,n,x,fvec,fjac,ldfjac,ftol,xtol,gtol,maxfev,
148
 
     *           wa(1),mode,factor,nprint,info,nfev,njev,ipvt,wa(n+1),
149
 
     *           wa(2*n+1),wa(3*n+1),wa(4*n+1),wa(5*n+1))
150
 
      if (info .eq. 8) info = 4
151
 
   10 continue
152
 
      return
153
 
c
154
 
c     last card of subroutine lmstr1.
155
 
c
156
 
      end