C*+*+*+*+*+* C This program was produced by the ATOMFT translator version 3.11 C Copyright(c) 1979-93, by Y. F. Chang C*+*+*+*+*+*+*+*+*+*+*+*+*+*+*+*+*+*+*+*+*+*+*+*+*+*+*+*+*+*+ C This is for the inputs below. C COPTION DOUBLE,DUMP=1 C DIFF ( Y , X , 1 ) = SQRT(SQRT(0.1 * X + 0.2)) C-------- C no instructions in Second input block C-------- DIMENSION TMPS( 38, 1),TMPV( 39) DIMENSION IPASS(20),RPASS(20) EQUIVALENCE (IPASS(1),NUMEQS),(IPASS(2),LENSER), A (IPASS(3),LENVAR),(IPASS(4),MPRINT),(IPASS(5),LIST), A (IPASS(6),MSTIFF),(IPASS(7),LRUN),(IPASS(8),KTRDCV), A (IPASS(9),KNTSTP),(IPASS(10),KTSTIF),(IPASS(11),KXPNUM), A (IPASS(12),KDIGS),(IPASS(13),KENDFG),(IPASS(14),NTERMS), A (IPASS(15),KOVER),(RPASS(1),RADIUS),(RPASS(2),H), A (RPASS(3),HNEW),(RPASS(4),ERRLIM),(RPASS(5),ADJSTF), A (RPASS(6),XPRINT),(RPASS(7),DLTXPT),(TMPS(1,1),TMPV(1)) A,(RPASS(8),START),(RPASS(9),END),(RPASS(10),ORDER) DOUBLE PRECISION RADIUS,H,HNEW,ERRLIM,ADJSTF, A XPRINT,DLTXPT,RPASS,AL,START,END,ORDER DOUBLE PRECISION TMPS,TMPV CHARACTER*6 NAMES EQUIVALENCE (TMPS(1,1),Y(1)) DOUBLE PRECISION CNDITN, Y(38), X(2), TMPAAC(30), TMPAAB(30), A TMPAAA(2) DOUBLE PRECISION DSQRT, SHIFT DIMENSION NAMES(2) DATA NAMES(1)/' X'/ DATA NAMES(2)/' Y'/ Y(33) = 2.11 C-------- C Initialize variables to default values. C-------- NTERMS = 1 NSTEPS = 40 MPRINT = 4 LIST = 0 DETUNE = 1.0 NUMEQS = 1 LENVAR = 38 ERRLIM = 1.D-12 LENSER = 30 KTRDCV = 2 ADJSTF = 1.E-2 H = 1.4131D0 START = 0.D0 END = 0.D0 MSTIFF = 0 KENDFG = 3 LRUN = 1 KXPNUM = 3070038 KDIGS = 15 CALL DHEAD(5,TMPV,NAMES,IPASS,RPASS) DLTXPT = 0.D0 C-------- C start of Third input block C-------- MPRINT = 2 NSTEPS = 1000 START = 0.1d0 END = 0.5d0 DLTXPT = 0.01d0 Y(1) = (0.8d0 * (START + 2.0d0) C * DSQRT(DSQRT(0.1d0 * START + 0.2d0))) WRITE(LIST,120) START,END,DLTXPT,Y(1) 120 FORMAT(8F16.10) C-------- C end of Third input block C-------- NUMEQS = 1 LENVAR = 38 KOVER = 1 CALL DHEAD(1,TMPV,NAMES,IPASS,RPASS) C More initializations C-------- IF(MSTIFF.GE.20) LENSER = MIN0(15,LENSER) IF(MSTIFF.EQ.21) LENSER = MIN0(10,LENSER) DLTXPT = DSIGN(DLTXPT,(END-START)) H = DSIGN(H,(END-START)) XPRINT = START + DLTXPT LRUN = 1 KTSTIF = 0 IF(LENSER.GT.(LENVAR- 8)) LENSER = LENVAR - 8 C-------- C Loop for integration steps. Inside the loop, print the desired output C-------- 17 DO 27 KINTS=1,NSTEPS KNTSTP = KINTS 19 CONTINUE X(1) = START X(2) = H C-------- C Preliminary series calculations C-------- TMPAAA(1) = (1.D-1*X(1)) + 2.D-1 TMPAAB(1) = DSQRT(TMPAAA(1)) TMPAAC(1) = DSQRT(TMPAAB(1)) Y(2) = TMPAAC(1)*(H) TMPAAA(2) = (1.D-1*X(2)) TMPAAB(2) = TMPAAA(2)/TMPAAB(1) *5.D-1 TMPAAC(2) = TMPAAB(2)/TMPAAC(1) *5.D-1 Y(3) = TMPAAC(2)*(H/2.D0) C-------- C Loop for series calculations C-------- DO 23 K= 4,LENSER KA = K - 1 KB = K - 2 KC = K - 3 TMPAAB(KA) = 0.D0 SHIFT = 0.D0 KZ = 1 + KA DO 1000 N=2, KA L = KZ - N 1000 SHIFT = SHIFT + TMPAAB(N)*TMPAAB(L) TMPAAB(KA) = SHIFT TMPAAB(KA) = - TMPAAB(KA)/TMPAAB(1) *5.D-1 TMPAAC(KA) = 0.D0 SHIFT = 0.D0 DO 1001 N=2, KA L = KZ - N 1001 SHIFT = SHIFT + TMPAAC(N)*TMPAAC(L) TMPAAC(KA) = SHIFT TMPAAC(KA) = (TMPAAB(KA) - TMPAAC(KA))/TMPAAC(1) *5.D-1 Y(LENVAR) = (K) Y(K) = TMPAAC(KA)*(H/(KA)) C-------- C Test and adjust H to avoid over/under flow. C-------- CALL DHEAD(2,TMPV,NAMES,IPASS,RPASS) IF(LRUN.EQ.0) GO TO 19 23 CONTINUE C-------- C Calculate radius of convergence and take optimum step. C-------- CALL DRDCV(TMPV,NAMES,IPASS,RPASS) 24 CALL DRSET(TMPV,NAMES,IPASS,RPASS) C-------- C no instructions in Fourth input block C-------- 25 GO TO (26,28,24), KENDFG 26 AL = RADIUS*DETUNE H = DSIGN(AL,H) START = START + HNEW 27 CONTINUE CALL DHEAD(3,TMPV,NAMES,IPASS,RPASS) 28 CONTINUE 29 STOP END