1 SUBROUTINE SDNTL (EPS,F,FA,HMAX,HOLD,IMPL,JTASK,MATDIM,MAXORD, 2 8 MINT,MITER,ML,MU,N,NDE,SAVE1,T,Y,YWT,H,MNTOLD,MTROLD,NFE,RC,YH, 3 8 A,CONVRG,EL,IER,IPVT,NQ,NWAIT,RH,RMAX,SAVE2,TQ,TREND,ISWFLG) 4C***BEGIN PROLOGUE SDNTL 5C***REFER TO SDRIV3 6C Subroutine SDNTL is called to set parameters on the first call 7C to SDSTP, on an internal restart, or when the user has altered 8C MINT, MITER, and/or H. 9C On the first call, the order is set to 1 and the initial derivatives 10C are calculated. RMAX is the maximum ratio by which H can be 11C increased in one step. It is initially RMINIT to compensate 12C for the small initial H, but then is normally equal to RMNORM. 13C If a failure occurs (in corrector convergence or error test), RMAX 14C is set at RMFAIL for the next increase. 15C If the caller has changed MINT, or if JTASK = 0, SDCST is called 16C to set the coefficients of the method. If the caller has changed H, 17C YH must be rescaled. If H or MINT has been changed, NWAIT is 18C reset to NQ + 2 to prevent further increases in H for that many 19C steps. Also, RC is reset. RC is the ratio of new to old values of 20C the coefficient L(0)*H. If the caller has changed MITER, RC is 21C set to 0 to force the partials to be updated, if partials are used. 22C***ROUTINES CALLED SDCST,SDSCL,SGEFA,SGESL,SGBFA,SGBSL,SNRM2 23C***DATE WRITTEN 790601 (YYMMDD) 24C***REVISION DATE 850320 (YYMMDD) 25C***CATEGORY NO. I1A2,I1A1B 26C***AUTHOR KAHANER, D. K., NATIONAL BUREAU OF STANDARDS, 27C SUTHERLAND, C. D., LOS ALAMOS NATIONAL LABORATORY 28C***END PROLOGUE SDNTL 29 REAL A(MATDIM,*), EL(13,12), EPS, H, HMAX, HNEW, HOLD, 30 8 OLDL0, RC, RH, RMAX, RMINIT, SAVE1(*), SAVE2(*), SMAX, SMIN, 31 8 SNRM2, SUM, SUM0, T, TQ(3,12), TREND, Y(*), YH(N,*), YWT(*) 32 INTEGER IPVT(*) 33 LOGICAL CONVRG, IER 34 PARAMETER(RMINIT = 10000.E0) 35C***FIRST EXECUTABLE STATEMENT SDNTL 36 IER = .FALSE. 37 IF (JTASK .GE. 0) THEN 38 IF (JTASK .EQ. 0) CALL SDCST (MAXORD, MINT, ISWFLG, EL, TQ) 39 RC = 0.E0 40 CONVRG = .FALSE. 41 TREND = 1.E0 42 RMAX = RMINIT 43 NQ = 1 44 NWAIT = 3 45 CALL F (N, T, Y, SAVE2) 46 NFE = NFE + 1 47 IF (IMPL .NE. 0) THEN 48 IF (MITER .EQ. 3) THEN 49 IFLAG = 0 50 CALL USERS (Y, YH, YWT, SAVE1, SAVE2, T, H, EL, IMPL, N, 51 8 NDE, IFLAG) 52 ELSE IF (IMPL .EQ. 1) THEN 53 IF (MITER .EQ. 1 .OR. MITER .EQ. 2) THEN 54 CALL FA (N, T, Y, A, MATDIM, ML, MU, NDE) 55 CALL SGEFA (A, MATDIM, N, IPVT, INFO) 56 IF (INFO .NE. 0) THEN 57 IER = .TRUE. 58 RETURN 59 END IF 60 CALL SGESL (A, MATDIM, N, IPVT, SAVE2, 0) 61 ELSE IF (MITER .EQ. 4 .OR. MITER .EQ. 5) THEN 62 CALL FA (N, T, Y, A(ML+1,1), MATDIM, ML, MU, NDE) 63 CALL SGBFA (A, MATDIM, N, ML, MU, IPVT, INFO) 64 IF (INFO .NE. 0) THEN 65 IER = .TRUE. 66 RETURN 67 END IF 68 CALL SGBSL (A, MATDIM, N, ML, MU, IPVT, SAVE2, 0) 69 END IF 70 ELSE IF (IMPL .EQ. 2) THEN 71 CALL FA (N, T, Y, A, MATDIM, ML, MU, NDE) 72 DO 150 I = 1,NDE 73 IF(A(I,1) .EQ. 0.E0) THEN 74 IER = .TRUE. 75 RETURN 76 ELSE 77 SAVE2(I) = SAVE2(I)/A(I,1) 78 END IF 79 150 CONTINUE 80 DO 155 I = NDE+1,N 81 155 A(I,1) = 0.E0 82 END IF 83 END IF 84 DO 170 I = 1,NDE 85 170 SAVE1(I) = SAVE2(I)/YWT(I) 86 SUM = SNRM2(NDE, SAVE1, 1) 87 SUM0 = 1.E0/MAX(1.E0, ABS(T)) 88 SMAX = MAX(SUM0, SUM) 89 SMIN = MIN(SUM0, SUM) 90 SUM = SMAX*SQRT(1.E0 + (SMIN/SMAX)**2)/SQRT(REAL(NDE)) 91 H = SIGN(MIN(2.E0*EPS/SUM, ABS(H)), H) 92 DO 180 I = 1,N 93 180 YH(I,2) = H*SAVE2(I) 94 ELSE 95 IF (MITER .NE. MTROLD) THEN 96 MTROLD = MITER 97 RC = 0.E0 98 CONVRG = .FALSE. 99 END IF 100 IF (MINT .NE. MNTOLD) THEN 101 MNTOLD = MINT 102 OLDL0 = EL(1,NQ) 103 CALL SDCST (MAXORD, MINT, ISWFLG, EL, TQ) 104 RC = RC*EL(1,NQ)/OLDL0 105 NWAIT = NQ + 2 106 END IF 107 IF (H .NE. HOLD) THEN 108 NWAIT = NQ + 2 109 HNEW = H 110 RH = H/HOLD 111 H = HOLD 112 CALL SDSCL (HMAX, N, NQ, RMAX, H, RC, RH, YH) 113 H = SIGN(MIN(ABS(H), ABS(HNEW)), H) 114 END IF 115 END IF 116 END 117