1 SUBROUTINE CDNTL (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 CDNTL 5C***REFER TO CDRIV3 6C Subroutine CDNTL is called to set parameters on the first call 7C to CDSTP, 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, CDCST 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 CDCST,CDSCL,CGEFA,CGESL,CGBFA,CGBSL,SCNRM2 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 CDNTL 29 COMPLEX A(MATDIM,*), SAVE1(*), SAVE2(*), Y(*), YH(N,*), YWT(*) 30 REAL EL(13,12), EPS, H, HMAX, HNEW, HOLD, OLDL0, RC, RH, RMAX, 31 8 RMINIT, SCNRM2, SMAX, SMIN, SUM, SUM0, T, TQ(3,12), TREND 32 INTEGER IPVT(*) 33 LOGICAL CONVRG, IER 34 PARAMETER(RMINIT = 10000.E0) 35C***FIRST EXECUTABLE STATEMENT CDNTL 36 IER = .FALSE. 37 IF (JTASK .GE. 0) THEN 38 IF (JTASK .EQ. 0) CALL CDCST (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 CGEFA (A, MATDIM, N, IPVT, INFO) 56 IF (INFO .NE. 0) THEN 57 IER = .TRUE. 58 RETURN 59 END IF 60 CALL CGESL (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 CGBFA (A, MATDIM, N, ML, MU, IPVT, INFO) 64 IF (INFO .NE. 0) THEN 65 IER = .TRUE. 66 RETURN 67 END IF 68 CALL CGBSL (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. CMPLX(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) = CMPLX(0.E0) 82 END IF 83 END IF 84 DO 170 I = 1,NDE 85 170 SAVE1(I) = SAVE2(I)/YWT(I) 86 SUM = SCNRM2(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 CDCST (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 CDSCL (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