Enter your search terms Submit search form
	Web pdp-10.trailing-edge.com

C	THE PROGRAM AMPBD.F4 IS A SAMPLE PROGRAM, WHEN USED IN
C	CONJUNCTION AMPBX(LOADED AND SAVED AS AMPBD.SAV) WILL
C	PRODUCE A RUNNING EXAMPLE OF THE USE OF AMPBX.
      SUBROUTINE AMPB2(IND,TEMP,X,DX,Y,F,N,ICOUNT,NITER,MTST)
      DIMENSION TEMP(500),Y(500),F(500),AK(2)
C  INTEGRATE TO NEXT POINT
401   M=N+1
      IF (ICOUNT-1)  5001, 406, 411
406   N1=NITER+1
      AK(1)=.92962963
      AK(2)=-.070370370
C  RESTORE PREVIOUS VALUES
411   X=TEMP(2*M+1)
      DO 421 I=2,M
      IP2M=I+2*M
      Y(I-1)=TEMP(IP2M)
      IP3M=IP2M+M
421   F(I-1)=TEMP(IP3M)
      IF (IND)  431, 5001, 501
431   IF (ICOUNT-4)  601, 701, 701
C  REDUCE DX
501   DX=DX/2.0
      ICOUNT=1
C  INTEGRATION BY RUNGE KUTTA
601   CALL RKPB2(TEMP,X,DX,Y,F,N)
      GO TO 1001
C  INTEGRATION BY ADAMS-MOULTON
C  INDEPENDENT VARIABLE
701   X=X+DX
      TEMP(M+1)=X
C  PREDICTOR VALUES
801   DO 821 I=2,M
      IPM=I+M
      IP2M=IPM+M
      IP3M=IP2M+M
      IP4M=IP3M+M
      IP5M=IP4M+M
      IP6M=IP5M+M
      Y(I-1)=TEMP(IP2M)+DX*(55.*TEMP(IP3M)-59.*TEMP(IP4M)+
     +     37.*TEMP(IP5M)-9.*TEMP(IP6M))/24.
      TEMP(IPM)=Y(I-1)
      IF (MTST)  811, 821,811
811	IF(ICOUNT-4) 5001,901,821
821     Y(I-1)=Y(I-1)+AK(1)/AK(2)*TEMP(1)
C  CORRECTOR VALUES
901   DO 991 J1=1,N1
 	CALL DERIV
      DO 921 I=2,M
	IPM=I+M
	IP2M=IPM+M
	IP3M=IP2M+M
	IP4M=IP3M+M
	IP5M=IP4M+M
      N=I-1
      Y(N)=TEMP(IP2M)+DX*(9.*F(I-1)+19.*TEMP(IP3M)-5.*TEMP(IP4M)
     +    +TEMP(IP5M))/24.
      TEMP(I)=AK(2)*(Y(N)-TEMP(IPM))
      IF (MTST)  911, 921,911
911   Y(N)=Y(N)+TEMP(I)
921	CONTINUE
991	CONTINUE
C  RESTORE NORMAL MODE
1001	IND=-1
5001	RETURN
	END
C     AMPB1************ADAMS-MOULTON INTEGRATION-1**********
      SUBROUTINE AMPB1(IND,TEMP,X,DX,Y,F,N,ICOUNT,NITER,MTST)
      DIMENSION TEMP(500),Y(500),F(500)
1     IF (IND)  101, 11, 101
11    ICOUNT=0
      GO TO 301
101   M=N+1
      DO 191 I=1,M
      DO 121 J1=1,ICOUNT
      I1=I+(ICOUNT+3-J1)*M
      I2=I1+M
121   TEMP(I2)=TEMP(I1)
191	CONTINUE
201	IF (IND) 301,301,302
302     IF (ICOUNT-6)  301, 211, 5001
211   DX=2.*DX
      DO 261 J1=1,3
      DO 241 I=1,M
      I1=I+(J1+3)*M
	IF (I-1) 5001,221,222
222     I2=I+(2*J1+3)*M
      TEMP(I1)=TEMP(I2)
      GO TO 241
221   TEMP(I1)=TEMP(I1)*2.
241	CONTINUE
261	CONTINUE
      ICOUNT=3
301   ICOUNT=  MIN0(ICOUNT+1,6)
      CALL RKPB1(TEMP,X,DX,Y,F,N)
1001  IND=-1
5001	RETURN
	END