1*DECK DNLS1Q 2 SUBROUTINE DNLS1Q (LUN, KPRINT, IPASS) 3C***BEGIN PROLOGUE DNLS1Q 4C***PURPOSE Quick check for DNLS1E, DNLS1 and DCOV. 5C***LIBRARY SLATEC 6C***TYPE DOUBLE PRECISION (SNLS1Q-S, DNLS1Q-D) 7C***KEYWORDS QUICK CHECK 8C***AUTHOR (UNKNOWN) 9C***DESCRIPTION 10C 11C This subroutine performs a quick check on the subroutines DNLS1E 12C (and DNLS1) and DCOV. 13C 14C***ROUTINES CALLED DENORM, DFCN1, DFCN2, DFCN3, DFDJC3, PASS, D1MACH, 15C DCOV, DNLS1E 16C***REVISION HISTORY (YYMMDD) 17C ?????? DATE WRITTEN 18C 890911 Removed unnecessary intrinsics. (WRB) 19C 890911 REVISION DATE from Version 3.2 20C 891214 Prologue converted to Version 4.0 format. (BAB) 21C 930214 Declarations sections added, code revised to test error 22C returns for all values of KPRINT and code polished. (WRB) 23C***END PROLOGUE DNLS1Q 24C .. Scalar Arguments .. 25 INTEGER IPASS, KPRINT, LUN 26C .. Local Scalars .. 27 DOUBLE PRECISION FNORM, FNORMS, ONE, SIGMA, TEMP1, TEMP2, TEMP3, 28 + TOL, TOL2, ZERO 29 INTEGER I, IFLAG, INFO, INFOS, IOPT, KONTRL, LDFJAC, LWA, M, N, 30 + NERR, NPRINT 31 LOGICAL FATAL 32C .. Local Arrays .. 33 DOUBLE PRECISION FJAC(10,2), FJROW(2), FJTJ(3), FVEC(10), WA(40), 34 + X(2) 35 INTEGER IW(2) 36C .. External Functions .. 37 DOUBLE PRECISION D1MACH, DENORM 38 INTEGER NUMXER 39 EXTERNAL D1MACH, DENORM, NUMXER 40C .. External Subroutines .. 41 EXTERNAL DFCN1, DFCN2, DFCN3, DFDJC3, PASS, DCOV, DNLS1E, XGETF, 42 + XSETF 43C .. Intrinsic Functions .. 44 INTRINSIC ABS, SQRT 45C***FIRST EXECUTABLE STATEMENT DNLS1Q 46 IF (KPRINT .GE. 2) WRITE (LUN,9000) 47C 48 IPASS = 1 49 INFOS = 1 50 FNORMS = 1.1151779D+01 51 M = 10 52 N = 2 53 LWA = 40 54 LDFJAC = 10 55 NPRINT = -1 56 IFLAG = 1 57 ZERO = 0.0D0 58 ONE = 1.0D0 59 TOL = MAX(SQRT(40.0D0*D1MACH(4)),1.0D-12) 60 TOL2 = SQRT(TOL) 61C 62C OPTION=2, the full Jacobian is stored and the user provides the 63C Jacobian. 64C 65 IOPT = 2 66 X(1) = 3.0D-1 67 X(2) = 4.0D-1 68 CALL DNLS1E(DFCN2,IOPT,M,N,X,FVEC,TOL,NPRINT,INFO,IW,WA,LWA) 69 FNORM = DENORM(M,FVEC) 70 IF (INFO.EQ.INFOS .AND. ABS(FNORM-FNORMS)/FNORMS.LE.TOL2) THEN 71 FATAL = .FALSE. 72 IF (KPRINT .GE. 3) CALL PASS(LUN,1,1) 73 ELSE 74 IPASS = 0 75 FATAL = .TRUE. 76 IF (KPRINT .GE. 2) CALL PASS(LUN,1,0) 77 ENDIF 78 IF ((FATAL.AND.KPRINT.GE.2) .OR. KPRINT.GE.3) 79 + WRITE (LUN,9010) INFOS,FNORMS,INFO,FNORM 80C 81C Form JAC-transpose*JAC. 82C 83 SIGMA = FNORM*FNORM/(M-N) 84 IFLAG = 2 85 CALL DFCN2(IFLAG,M,N,X,FVEC,FJAC,LDFJAC) 86 DO 20 I = 1,3 87 FJTJ(I) = ZERO 88 20 CONTINUE 89 DO 30 I = 1,M 90 FJTJ(1) = FJTJ(1) + FJAC(I,1)**2 91 FJTJ(2) = FJTJ(2) + FJAC(I,1)*FJAC(I,2) 92 FJTJ(3) = FJTJ(3) + FJAC(I,2)**2 93 30 CONTINUE 94C 95C Calculate the covariance matrix. 96C 97 CALL DCOV(DFCN2,IOPT,M,N,X,FVEC,FJAC,LDFJAC,INFO,WA(1),WA(N+1), 98 + WA(2*N+1),WA(3*N+1)) 99C 100C Form JAC-transpose*JAC * covariance matrix (should = SIGMA*I). 101C 102 TEMP1 = (FJTJ(1)*FJAC(1,1)+FJTJ(2)*FJAC(1,2))/SIGMA 103 TEMP2 = (FJTJ(1)*FJAC(1,2)+FJTJ(2)*FJAC(2,2))/SIGMA 104 TEMP3 = (FJTJ(2)*FJAC(1,2)+FJTJ(3)*FJAC(2,2))/SIGMA 105 IF (INFO.EQ.INFOS .AND. ABS(TEMP1-ONE).LT.TOL2 .AND. 106 + ABS(TEMP2).LT.TOL2 .AND. ABS(TEMP3-ONE).LT.TOL2) THEN 107 FATAL = .FALSE. 108 IF (KPRINT .GE. 3) CALL PASS(LUN,2,1) 109 ELSE 110 IPASS = 0 111 FATAL = .TRUE. 112 IF (KPRINT .GE. 2) CALL PASS(LUN,2,0) 113 ENDIF 114 IF ((FATAL.AND.KPRINT.GE.2) .OR. KPRINT.GE.3) 115 + WRITE (LUN,9020) INFOS,INFO,TEMP1,TEMP2,TEMP3 116C 117C OPTION=1, the full Jacobian is stored and the code approximates 118C the Jacobian. 119C 120 IOPT = 1 121 X(1) = 3.0D-1 122 X(2) = 4.0D-1 123 CALL DNLS1E(DFCN1,IOPT,M,N,X,FVEC,TOL,NPRINT,INFO,IW,WA,LWA) 124 FNORM = DENORM(M,FVEC) 125 IF (INFO.EQ.INFOS .AND. ABS(FNORM-FNORMS)/FNORMS.LE.TOL2) THEN 126 FATAL = .FALSE. 127 IF (KPRINT .GE. 3) CALL PASS(LUN,3,1) 128 ELSE 129 IPASS = 0 130 FATAL = .TRUE. 131 IF (KPRINT .GE. 2) CALL PASS(LUN,3,0) 132 ENDIF 133 IF ((FATAL.AND.KPRINT.GE.2) .OR. KPRINT.GE.3) 134 + WRITE (LUN,9010) INFOS,FNORMS,INFO,FNORM 135C 136C Form JAC-transpose*JAC. 137C 138 SIGMA = FNORM*FNORM/(M-N) 139 IFLAG = 1 140 CALL DFDJC3(DFCN1,M,N,X,FVEC,FJAC,LDFJAC,IFLAG,ZERO,WA) 141 DO 60 I = 1,3 142 FJTJ(I) = ZERO 143 60 CONTINUE 144 DO 70 I = 1,M 145 FJTJ(1) = FJTJ(1) + FJAC(I,1)**2 146 FJTJ(2) = FJTJ(2) + FJAC(I,1)*FJAC(I,2) 147 FJTJ(3) = FJTJ(3) + FJAC(I,2)**2 148 70 CONTINUE 149C 150C Calculate the covariance matrix. 151C 152 CALL DCOV(DFCN1,IOPT,M,N,X,FVEC,FJAC,LDFJAC,INFO,WA(1),WA(N+1), 153 + WA(2*N+1),WA(3*N+1)) 154C 155C Form JAC-transpose*JAC * covariance matrix (should = SIGMA*I). 156C 157 TEMP1 = (FJTJ(1)*FJAC(1,1)+FJTJ(2)*FJAC(1,2))/SIGMA 158 TEMP2 = (FJTJ(1)*FJAC(1,2)+FJTJ(2)*FJAC(2,2))/SIGMA 159 TEMP3 = (FJTJ(2)*FJAC(1,2)+FJTJ(3)*FJAC(2,2))/SIGMA 160 IF (INFO.EQ.INFOS .AND. ABS(TEMP1-ONE).LT.TOL2 .AND. 161 + ABS(TEMP2).LT.TOL2 .AND. ABS(TEMP3-ONE).LT.TOL2) THEN 162 FATAL = .FALSE. 163 IF (KPRINT .GE. 3) CALL PASS(LUN,4,1) 164 ELSE 165 IPASS = 0 166 FATAL = .TRUE. 167 IF (KPRINT .GE. 2) CALL PASS(LUN,4,0) 168 ENDIF 169 IF ((FATAL.AND.KPRINT.GE.2) .OR. KPRINT.GE.3) 170 + WRITE (LUN,9020) INFOS,INFO,TEMP1,TEMP2,TEMP3 171C 172C OPTION=3, the full Jacobian is not stored. Only the product of 173C the Jacobian transpose and Jacobian is stored. The user provides 174C the Jacobian one row at a time. 175C 176 IOPT = 3 177 X(1) = 3.0D-1 178 X(2) = 4.0D-1 179 CALL DNLS1E(DFCN3,IOPT,M,N,X,FVEC,TOL,NPRINT,INFO,IW,WA,LWA) 180 FNORM = DENORM(M,FVEC) 181 IF (INFO.EQ.INFOS .AND. ABS(FNORM-FNORMS)/FNORMS.LE.TOL2) THEN 182 FATAL = .FALSE. 183 IF (KPRINT .GE. 3) CALL PASS(LUN,5,1) 184 ELSE 185 IPASS = 0 186 FATAL = .TRUE. 187 IF (KPRINT .GE. 2) CALL PASS(LUN,5,0) 188 ENDIF 189 IF ((FATAL.AND.KPRINT.GE.2) .OR. KPRINT.GE.3) 190 + WRITE (LUN,9010) INFOS,FNORMS,INFO,FNORM 191C 192C Form JAC-transpose*JAC. 193C 194 SIGMA = FNORM*FNORM/(M-N) 195 DO 100 I = 1,3 196 FJTJ(I) = ZERO 197 100 CONTINUE 198 IFLAG = 3 199 DO 110 I = 1,M 200 CALL DFCN3(IFLAG,M,N,X,FVEC,FJROW,I) 201 FJTJ(1) = FJTJ(1) + FJROW(1)**2 202 FJTJ(2) = FJTJ(2) + FJROW(1)*FJROW(2) 203 FJTJ(3) = FJTJ(3) + FJROW(2)**2 204 110 CONTINUE 205C 206C Calculate the covariance matrix. 207C 208 CALL DCOV(DFCN3,IOPT,M,N,X,FVEC,FJAC,LDFJAC,INFO,WA(1),WA(N+1), 209 + WA(2*N+1),WA(3*N+1)) 210C 211C Form JAC-transpose*JAC * covariance matrix (should = SIGMA*I). 212C 213 TEMP1 = (FJTJ(1)*FJAC(1,1)+FJTJ(2)*FJAC(1,2))/SIGMA 214 TEMP2 = (FJTJ(1)*FJAC(1,2)+FJTJ(2)*FJAC(2,2))/SIGMA 215 TEMP3 = (FJTJ(2)*FJAC(1,2)+FJTJ(3)*FJAC(2,2))/SIGMA 216 IF (INFO.EQ.INFOS .AND. ABS(TEMP1-ONE).LT.TOL2 .AND. 217 + ABS(TEMP2).LT.TOL2 .AND. ABS(TEMP3-ONE).LT.TOL2) THEN 218 FATAL = .FALSE. 219 IF (KPRINT .GE. 3) CALL PASS(LUN,6,1) 220 ELSE 221 IPASS = 0 222 FATAL = .TRUE. 223 IF (KPRINT .GE. 2) CALL PASS(LUN,6,0) 224 ENDIF 225 IF ((FATAL.AND.KPRINT.GE.2) .OR. KPRINT.GE.3) 226 + WRITE (LUN,9020) INFOS,INFO,TEMP1,TEMP2,TEMP3 227C 228C Test improper input parameters. 229C 230 CALL XGETF (KONTRL) 231 IF (KPRINT .LE. 2) THEN 232 CALL XSETF (0) 233 ELSE 234 CALL XSETF (1) 235 ENDIF 236 FATAL = .FALSE. 237 CALL XERCLR 238C 239 IF (KPRINT .GE. 3) WRITE (LUN, 9060) 240C 241 LWA = 35 242 IOPT = 2 243 X(1) = 3.0D-1 244 X(2) = 4.0D-1 245 CALL DNLS1E(DFCN2,IOPT,M,N,X,FVEC,TOL,NPRINT,INFO,IW,WA,LWA) 246 IF (INFO.NE.0 .OR. NUMXER(NERR).NE.2) FATAL = .TRUE. 247C 248 M = 0 249 CALL DCOV(DFCN2,IOPT,M,N,X,FVEC,FJAC,LDFJAC,INFO,WA(1),WA(N+1), 250 + WA(2*N+1),WA(3*N+1)) 251 IF (INFO.NE.0 .OR. NUMXER(NERR).NE.2) FATAL = .TRUE. 252C 253C Restore KONTRL and check to see if the tests of error detection 254C passed. 255C 256 CALL XSETF (KONTRL) 257 IF (FATAL) THEN 258 IPASS = 0 259 IF (KPRINT .GE. 2) THEN 260 WRITE (LUN, 9070) 261 ENDIF 262 ELSE 263 IF (KPRINT .GE. 3) THEN 264 WRITE (LUN, 9080) 265 ENDIF 266 ENDIF 267C 268C Print PASS/FAIL message. 269C 270 IF (IPASS.EQ.1 .AND. KPRINT.GE.2) WRITE (LUN,9100) 271 IF (IPASS.EQ.0 .AND. KPRINT.GE.1) WRITE (LUN,9110) 272C 273 130 RETURN 274C 275 9000 FORMAT ('1' / ' Test DNLS1E, DNLS1 and DCOV') 276 9010 FORMAT (' EXPECTED VALUE OF INFO AND RESIDUAL NORM', I5, D20.9 / 277 + ' RETURNED VALUE OF INFO AND RESIDUAL NORM', I5, D20.9 /) 278 9020 FORMAT (' EXPECTED AND RETURNED VALUE OF INFO', I5, 10X, I5 / 279 + ' RETURNED PRODUCT OF (J-TRANS*J)*COVARIANCE MATRIX/SIGMA' 280 + / ' (SHOULD = THE IDENTITY -- 1.0, 0.0, 1.0)' / 3D20.9 /) 281 9060 FORMAT (/ ' TRIGGER 2 ERROR MESSAGES',/) 282 9070 FORMAT (' AT LEAST ONE INCORRECT ARGUMENT TEST FAILED') 283 9080 FORMAT (' ALL INCORRECT ARGUMENT TESTS PASSED') 284 9100 FORMAT (/' *************DNLS1E PASSED ALL TESTS*****************') 285 9110 FORMAT (/' ************DNLS1E FAILED SOME TESTS*****************') 286 END 287