home *** CD-ROM | disk | FTP | other *** search
/ Geek Gadgets 1 / ADE-1.bin / ade-dist / octave-1.1.1p1-src.tgz / tar.out / fsf / octave / libcruft / minpack / fdjac1.f < prev    next >
Text File  |  1996-09-28  |  5KB  |  152 lines
  1.       SUBROUTINE FDJAC1(FCN,N,X,FVEC,FJAC,LDFJAC,IFLAG,ML,MU,EPSFCN,
  2.      *                  WA1,WA2)
  3.       INTEGER N,LDFJAC,IFLAG,ML,MU
  4.       DOUBLE PRECISION EPSFCN
  5.       DOUBLE PRECISION X(N),FVEC(N),FJAC(LDFJAC,N),WA1(N),WA2(N)
  6.       EXTERNAL FCN
  7. C     **********
  8. C
  9. C     SUBROUTINE FDJAC1
  10. C
  11. C     THIS SUBROUTINE COMPUTES A FORWARD-DIFFERENCE APPROXIMATION
  12. C     TO THE N BY N JACOBIAN MATRIX ASSOCIATED WITH A SPECIFIED
  13. C     PROBLEM OF N FUNCTIONS IN N VARIABLES. IF THE JACOBIAN HAS
  14. C     A BANDED FORM, THEN FUNCTION EVALUATIONS ARE SAVED BY ONLY
  15. C     APPROXIMATING THE NONZERO TERMS.
  16. C
  17. C     THE SUBROUTINE STATEMENT IS
  18. C
  19. C       SUBROUTINE FDJAC1(FCN,N,X,FVEC,FJAC,LDFJAC,IFLAG,ML,MU,EPSFCN,
  20. C                         WA1,WA2)
  21. C
  22. C     WHERE
  23. C
  24. C       FCN IS THE NAME OF THE USER-SUPPLIED SUBROUTINE WHICH
  25. C         CALCULATES THE FUNCTIONS. FCN MUST BE DECLARED
  26. C         IN AN EXTERNAL STATEMENT IN THE USER CALLING
  27. C         PROGRAM, AND SHOULD BE WRITTEN AS FOLLOWS.
  28. C
  29. C         SUBROUTINE FCN(N,X,FVEC,IFLAG)
  30. C         INTEGER N,IFLAG
  31. C         DOUBLE PRECISION X(N),FVEC(N)
  32. C         ----------
  33. C         CALCULATE THE FUNCTIONS AT X AND
  34. C         RETURN THIS VECTOR IN FVEC.
  35. C         ----------
  36. C         RETURN
  37. C         END
  38. C
  39. C         THE VALUE OF IFLAG SHOULD NOT BE CHANGED BY FCN UNLESS
  40. C         THE USER WANTS TO TERMINATE EXECUTION OF FDJAC1.
  41. C         IN THIS CASE SET IFLAG TO A NEGATIVE INTEGER.
  42. C
  43. C       N IS A POSITIVE INTEGER INPUT VARIABLE SET TO THE NUMBER
  44. C         OF FUNCTIONS AND VARIABLES.
  45. C
  46. C       X IS AN INPUT ARRAY OF LENGTH N.
  47. C
  48. C       FVEC IS AN INPUT ARRAY OF LENGTH N WHICH MUST CONTAIN THE
  49. C         FUNCTIONS EVALUATED AT X.
  50. C
  51. C       FJAC IS AN OUTPUT N BY N ARRAY WHICH CONTAINS THE
  52. C         APPROXIMATION TO THE JACOBIAN MATRIX EVALUATED AT X.
  53. C
  54. C       LDFJAC IS A POSITIVE INTEGER INPUT VARIABLE NOT LESS THAN N
  55. C         WHICH SPECIFIES THE LEADING DIMENSION OF THE ARRAY FJAC.
  56. C
  57. C       IFLAG IS AN INTEGER VARIABLE WHICH CAN BE USED TO TERMINATE
  58. C         THE EXECUTION OF FDJAC1. SEE DESCRIPTION OF FCN.
  59. C
  60. C       ML IS A NONNEGATIVE INTEGER INPUT VARIABLE WHICH SPECIFIES
  61. C         THE NUMBER OF SUBDIAGONALS WITHIN THE BAND OF THE
  62. C         JACOBIAN MATRIX. IF THE JACOBIAN IS NOT BANDED, SET
  63. C         ML TO AT LEAST N - 1.
  64. C
  65. C       EPSFCN IS AN INPUT VARIABLE USED IN DETERMINING A SUITABLE
  66. C         STEP LENGTH FOR THE FORWARD-DIFFERENCE APPROXIMATION. THIS
  67. C         APPROXIMATION ASSUMES THAT THE RELATIVE ERRORS IN THE
  68. C         FUNCTIONS ARE OF THE ORDER OF EPSFCN. IF EPSFCN IS LESS
  69. C         THAN THE MACHINE PRECISION, IT IS ASSUMED THAT THE RELATIVE
  70. C         ERRORS IN THE FUNCTIONS ARE OF THE ORDER OF THE MACHINE
  71. C         PRECISION.
  72. C
  73. C       MU IS A NONNEGATIVE INTEGER INPUT VARIABLE WHICH SPECIFIES
  74. C         THE NUMBER OF SUPERDIAGONALS WITHIN THE BAND OF THE
  75. C         JACOBIAN MATRIX. IF THE JACOBIAN IS NOT BANDED, SET
  76. C         MU TO AT LEAST N - 1.
  77. C
  78. C       WA1 AND WA2 ARE WORK ARRAYS OF LENGTH N. IF ML + MU + 1 IS AT
  79. C         LEAST N, THEN THE JACOBIAN IS CONSIDERED DENSE, AND WA2 IS
  80. C         NOT REFERENCED.
  81. C
  82. C     SUBPROGRAMS CALLED
  83. C
  84. C       MINPACK-SUPPLIED ... DPMPAR
  85. C
  86. C       FORTRAN-SUPPLIED ... DABS,DMAX1,DSQRT
  87. C
  88. C     MINPACK. VERSION OF JUNE 1979.
  89. C     BURTON S. GARBOW, KENNETH E. HILLSTROM, JORGE J. MORE
  90. C
  91. C     **********
  92.       INTEGER I,J,K,MSUM
  93.       DOUBLE PRECISION EPS,EPSMCH,H,TEMP,ZERO
  94.       DOUBLE PRECISION DPMPAR
  95.       DATA ZERO /0.0D0/
  96. C
  97. C     EPSMCH IS THE MACHINE PRECISION.
  98. C
  99.       EPSMCH = DPMPAR(1)
  100. C
  101.       EPS = DSQRT(DMAX1(EPSFCN,EPSMCH))
  102.       MSUM = ML + MU + 1
  103.       IF (MSUM .LT. N) GO TO 40
  104. C
  105. C        COMPUTATION OF DENSE APPROXIMATE JACOBIAN.
  106. C
  107.          DO 20 J = 1, N
  108.             TEMP = X(J)
  109.             H = EPS*DABS(TEMP)
  110.             IF (H .EQ. ZERO) H = EPS
  111.             X(J) = TEMP + H
  112.             CALL FCN(N,X,WA1,IFLAG)
  113.             IF (IFLAG .LT. 0) GO TO 30
  114.             X(J) = TEMP
  115.             DO 10 I = 1, N
  116.                FJAC(I,J) = (WA1(I) - FVEC(I))/H
  117.    10          CONTINUE
  118.    20       CONTINUE
  119.    30    CONTINUE
  120.          GO TO 110
  121.    40 CONTINUE
  122. C
  123. C        COMPUTATION OF BANDED APPROXIMATE JACOBIAN.
  124. C
  125.          DO 90 K = 1, MSUM
  126.             DO 60 J = K, N, MSUM
  127.                WA2(J) = X(J)
  128.                H = EPS*DABS(WA2(J))
  129.                IF (H .EQ. ZERO) H = EPS
  130.                X(J) = WA2(J) + H
  131.    60          CONTINUE
  132.             CALL FCN(N,X,WA1,IFLAG)
  133.             IF (IFLAG .LT. 0) GO TO 100
  134.             DO 80 J = K, N, MSUM
  135.                X(J) = WA2(J)
  136.                H = EPS*DABS(WA2(J))
  137.                IF (H .EQ. ZERO) H = EPS
  138.                DO 70 I = 1, N
  139.                   FJAC(I,J) = ZERO
  140.                   IF (I .GE. J - MU .AND. I .LE. J + ML)
  141.      *               FJAC(I,J) = (WA1(I) - FVEC(I))/H
  142.    70             CONTINUE
  143.    80          CONTINUE
  144.    90       CONTINUE
  145.   100    CONTINUE
  146.   110 CONTINUE
  147.       RETURN
  148. C
  149. C     LAST CARD OF SUBROUTINE FDJAC1.
  150. C
  151.       END
  152.