Google
 

Trailing-Edge - PDP-10 Archives - decuslib10-01 - 43,50110/bidyp2.eng
There are 2 other files named bidyp2.eng in the archive. Click here to see a list.
100'  NAME--BIDYP2
110'
120'  DESCRIPTION--GENERAL PURPOSE DYNAMIC PROGRAMMING PROGRAM.
130'  IT IS LIMITED TO ONE STATE VARIABLE.
140'
150'  SOURCE--A.O. CONVERSE, THAYER SCHOOL OF ENGINEERING
160'
170'  INSTRUCTIONS--THE SOLUTION IS IMBEDDED IN THE INLET STATE.
180'  THE USER MUST SUPPLY THE FOLLOWING FUCNTIONAL RELATIONSHIPS:
190'     1.770-809 TRANSFORMATION FUNCTION.  IN THE PRESENT EXAMPLE
200'       IT RELATES THE ALLOCATION VARIABLE TO THE INLET AND
210'       OUTLET STATES.  USUALLY IT WOULD RELATE EITHER THE OUTLET
220'       OR INLET STATES TO THE ALLOCATION VARIABLE AND THE OTHER
230'       STATE VARIABLE. IT IS CALLED FNA.
240'     2.810-849 RETURN FUNCTION.  IT RELATES THE RETURN ON THE
250'       CURRENT STAGE THE MOST CONVENIENT OF THE FOLLOWING THREE
260'       VARIABLES:INLET STATE, OUTLET STATE, AND ALLOCATION VAR.
270'     3.850-869SETS THE VALUES OF THE INLET STATE VARIABLE IN WHICH
280'       THE SOLUTION IS TO BE IMBEDDED.
290'     4.870-960, INDICATES THE MAXIMUM RETURN FROM THE SUBSEQUENT
300'       AS A FUNCTION OF THE OUTLET STATE OF THE CURRENT STAGE.
310'       AS IT NOW STANDS IT USES LINEAR INTERPOLATION.
320'     5.970-2140, EVALUATES MAXIMUM RETURN ON CURRENT AND 
330'       SUBSEQUENT STAGES AND OPTIMAL ALLOCATION TO CURRENT STAGE
340'        AS A FUNCTION OF: THE INLET STATE, THE LOWER LIMIT ON THE
350'        VARIABLE BEING TREATED AS THE INDEPENDENT VARIABLE IN
360'        SEARCH, THE UPPER LIMIT OF THE SAME, THE MINIMUM ALLOWED
370'        SEPARATION AND MAXIMUM SIZE OF THE FINAL INTERVAL. A
380'        FIBONACCI SEARCH IS USED.
390'     6.710-769, SPECIAL FUNCTION NEED IN THE EXAMPLE. IT RELATES
400'       SOLUTE CONCENTRATION IN THE RAFFINATE TO THE SOLVENT 
410'        SOLUTE CONCENTRATION IN THE EXTRACT.
420'
430'  NOTE:FOR A LISTING OF THE NOMENCLATURE CONTINUE LISTING THROUGH
440'  LINE 1460.
450'
460'
470'  *  *  *  *  *  *  *  MAIN PROGRAM  *  *  *  *  *  *  *  *  *  *  *
480'
490' ***** NOMENCLATURE *****
500'
510' F(I,J) = MATRIX RELATING OPT CONTROL(J=2), TOTAL MAX RETURN(J=3),
520'          OPT OUTLET STATE(J=4) TO THE VARIOUS VALUES(I=1 TO N)
530'          OF THE INLET STATE(J=1).
540' G(I,J) = F(I,J) FOR THE FOLLOWING STAGES
550' FNA = TRANSFORMATION FUNCTION
560' FNB = RETURN FUNCTION
570' FNC = MAX RETURN FROM SUBSEQUENT STAGES
580' FNE = MAXRETURN ON CURRENT AND SUBSEQUENT STAGES
590' FNH = VALUE OF INLET STATE
600' FNI = SOLUTE CONC IN EXTRACT AS A FUNC OF CONC IN RAFFINATE
610' A1 = SALE PRICE OF EXTRACT
620' B1 = COST OF SOLVENT
630' Q = FLOW RATE OF FEED
640' K1 = NUMBER OF STAGES
650' I1 = NUMBER OF VALUES OF THE STATE VARIABLE
660'
670 DIM F(21,4),G(21,4)
680'
690'
700 DEF FNI(X)
710' FNI(X) = SOLUTE CONC IN EXTRACT ASA FUNC OF SOLUTE CONC IF
720' RAFFINATE(X)
730 LET Z9 = -1.52E-3 + 2.78*X - 3.385*(X^2) - 107.7*(X^3)
740 LET FNI = Z9 + 404*(X^4)
750 FNEND
760'
770'FNA(X1,X2) = SOLVENT ALLOCATION AS A FUNCTION OF INLET SOLUTE
780'CONCENTRATION, X1, AND OUT SOLUTE CONC IN RAFFINATE,X2
790 DEF FNA(X1,X2) = (X1-X2)*Q/FNI(X2)
800'
810' FNB(X1,X2) = RETURN FROM CURRENT STAGE AS A FUNCTION OF INLET,X1,
820' AND OUTLET,X2, STATES
830 DEF FNB(X1,X2) = (A1*FNI(X2) - B1)*FNA(X1,X2)
840'
850 DEF FNH(I,K) = 0.2 - 0.01*(I-1)' I TH VALUE OF INLET STATE
860'     AT THE K TH STAGE 
870 DEF FNC(X)' MAX RETURN FROM SUBSEQUENT STAGES AS A FUNCTION OF THE
880'    OUTLET STATE, X.
890 IF K < K1 THEN 920' NOT AT LAST STAGE
900 LET FNC = 0' NO SUBSEQUENT STAGES ,HENCE THE VALUE IS 0
910 GO TO 960' GO TO FNEND
920 LET L = INT((0.2-X)/0.01 + 1)' IDENTIFIES ROW NUMBER IN G(I,J)
930'    OF THE NEXT SMALLEST CONCENTRATION
940 LET FNC = G(L,3)+((X-G(L,1))/(G(L+1,1)-G(L,1)))*(G(L+1,3)-G(L,3))
950' THE ABOVE IS LINEAR INTERPOLATION
960 FNEND
970 DEF FNE(X1,A,B,E,H)' MAX RETURN ON CURRENT AND SUBSEQUENT
980'     STAGES AS A FUNCTION OF INLET STATE,X1,LOWER LIMIT ON ALLOWED
990'     OUTLET STATES,A, UPPER LIMIT ON ALLOWED OUTLET STATES,B,
1000'     MIN SEPARATION OF EVALUATIONS,E, AND MAX SIZE OF FINAL
1010'     INTERVAL,H.
1020'  (REVISED 4/17/68 BY MALCOLM LEWIS, THAYER SCHOOL)
1030'
1040' THIS PROGRAM LOCATES AND EVALUATES THE MAXIMUM OF A UNIMODAL
1050'    FUNCTION OF ONE VARIABLE WITHIN A SPECIFIED INTERVAL
1060'    OF THAT VARIABLE BY FIBONACCI METHOD TAKEN FROM LIBRARY
1070'    PROGRAM ENTITLED "FIB***".
1080'
1090'
1100' GENERAL INFORMATION:
1110'     NOMENCLATURE:
1120'          A=LOWER LIMIT OF INITIAL INTERVAL
1130'          B=UPPER LIMIT OF INITIAL INTERVAL
1140'          C=WIDTH OF INITIAL INTERVAL
1150'          E=MINIMUM ALLOWED SEPARATION OF EVALUATIONS
1160'          F=CURRENT VALUE OF FUNCTION
1170'          F1=EVALUATION #1 FOR CURRENT INTERVAL
1180'          F2=EVALUATION #2 FOR CURRENT INTERVAL
1190'          G1,G2,G3=FIBONACCI NUMBERS IN CURRENT USE
1200'          H=MAXIMUM ALLOWABLE FINAL INTERVAL SIZE
1210'          H1=COMPUTED WIDTH OF FINAL INTERVAL
1220'          I2 = COUNTER IN THE "L2 COMPUTATION" SUBROUTINE
1230'          J=COMPUTED TOTAL NUMBER OF EXPERIMENTS REQUIRED
1240'          K$=INPUT SWITCH ON INTERMEDIATE PRINTOUT OPTION
1250'          L1,L2=INTERVAL FRACTION POSITION INDICES
1260'          N=COUNTER IN MAIN SEARCH ROUTINE
1270'          X=CURRENT VALUE OF INDEPENDENT VARIABLE
1280'          S1,S2 = LOCATIONS OF F1 AND F2
1290'          Z=INTERMEDIATE VARIABLE IN MAIN SEARCH ROUTINE
1300'     
1310'     COMPUTATION PROCEDURE:
1320'          F1,F2 ARE PLACED AT X1>X2.  F1,F2 ARE COMPARED AND
1330'            THE INTERVAL (A,B) IS DIMINISHED  BY REDEFINING
1340'            A OR B ACCORDING TO THE FOLLOWING SCHEME:
1350'
1360'   IF F1<F2 THEN             *  IF F1>F2 THEN
1370'       REDEFINE: B=S1       *       REDEFINE: A=S2
1380'                 S1=S2                        S2=S1
1390'                 F1=F2       *                F2=F1
1400'                 S2=B-L2*C   *                S1=A+L2*C
1410'
1420'          THE ENTIRE CYCLE IS REPEATED J TIMES, AT WHICH POINT
1430'            THE INTERVAL WIDTH IS LESS THAN H.
1440'
1450'
1460 REM  ****************MAIN PROGRAM********************
1470'****** INITIALIZATION ******
1480 LET C=B-A  ' SIZE OF INITIAL INTERVAL
1490 GOSUB 1990  ' GO COMPUTE L2
1500 LET L1=1
1510 LET X=A+L2*C  ' SET X AT DIST(L2*C) TO RIGHT OF POINT A.
1520 LET S1=X ' START S1 AND S2 AT THE SAME POSITION
1530 LET S2=X
1540 GO SUB 1910   '   EVALUATE F AT X=S1
1550 LET F1=F  ' SET FI = F(S1)
1560'
1570' ****** MAIN SEARCH ROUTINE *******
1580 FOR N=2 TO J  ' PERFORM (J-2) EXPERIMENTS
1590 IF A <S2 THEN 1660  ' FOR A > S2 DECREMENT X.(REMEMBER:S2>S1)
1600 LET X=A+L2*C ' FOR A>S2, INCREMENT X: X(I)=A+L(I)*(B-A)
1610 LET S2 = S1  ' MOVE S2 TO RIGHT TO S1
1620 LET S1 = X  '  MOVE S1 TO NEW HIGHER X-VALUE
1630 GO SUB 1910  ' TO EVALUATE F AT X=S1
1640 LET F1 = F   ' SET F1 = F(S1)
1650 GO TO 1710
1660 LET X=B-L2*C  ' FOR A<S2, DECREMENT X: X(I) =B-L(I)*(B-A)
1670 LET S1=S2  ' MOVE S1 TO LEFT TO S2
1680 LET S2 = X  ' MOVE S2 TO NEW LOWER VALUE OF X
1690 GO SUB 1910 ' TO EVALUATE F AT X=S2
1700 LET F2=F '  SET F2=F(S2)
1710'
1720 IF F1<F2 THEN 1760  ' FOR F1>F2, MOVE SEARCH REGION TO RIGHT
1730 LET A = S2 ' RAISE THE LOWER BOUND ON THE SEARCH REGION
1740 LET F2=F1  ' SAVE NEW MAX VALUE OF F(X)
1750 GO TO 1790
1760 LET F1=F2  ' FOR F1<F2, MOVE REGION TO LEFT; SAVE NEW F(MAX)
1770 LET B=S1  ' LOWER UPPER BOUND ON THE SEARCH REGION
1780' RECURSIVELY COMPUTE L2
1790 LET Z=L2  ' TEMPORARILY SAVE  OLD VALUE OF L2
1800 LET L2=L1-L2         ' SET L(N+1)=L(N-1)-L(N)
1810 LET L1=Z  ' LET L1= OLD L2-VALUE
1820 NEXT N
1830'
1840' ******* FINAL OUTPUT *******
1850 LET X=(A+B)/2  ' COMPUTE MIDPOINT OF FINAL INTERVAL
1860 GOSUB 1910  ' GO EVALUATE FINAL VALUE OF F(X)
1870 LET FNE = F' SETS FUNCTION TO MAX VALUE
1880 GO TO 2140' FNEND
1890 REM
1900 REM
1910 REM  *************DEFINE YOUR FUNCTION HERE**************
1920 LET F = FNB(X1,X) + FNC(X)' FUNCTION TO BE MAXIMIZED
1930 REM
1940REM
1950REM
1960REM
1970 RETURN
1980'
1990' ******* SUBR TO COMPUTE L2********
2000 LET G1=1  ' FIRST FIBONACCI NUMBER
2010 LET G2=1  ' SECOND FIBONACCI NUMBER
2020 LET I2 = 2' THE ITERATION STARTS AT I2 = 2
2030 LET G3=G2+G1  ' FIBONACCI RECURSION FORMULA
2040 LET H1=(C+G1*E)/G3  ' WIDTH OF FINAL INTERVAL, FOR I TRIALS
2050 IF H1<=H THEN 2110  ' IF H1 ALLOWABLE, COMPUTE L2
2060 LET I2 = I2 + 1' H1 TOO LARGE; TRY N = I2 + 1
2070 LET G1=G2  ' NOW G(I-2)=G(I-1)
2080 LET G2=G3  ' NOW G(I-1)=G(I)
2090 IF I2 > 50 THEN 2110' IF N IS REASONABLY LARGE(>50), CALC L2
2100 GO TO 2030'GO FURTHER OUT IN THE FIBONACCI SERIES
2110 LET L2=G2/G3+((-1)^I)*E/(G3*C)  ' COMPUTE L2
2120 LET J = I2' SETS TOTAL NUMBER OF EVALUATIONS = I2
2130 RETURN
2140 FNEND
2150'       ARE TO BE USED AT THE K TH STAGE.
2160' ***** MAIN PROGRAM *****
2170'
2180 PRINT"THIS PROGRAM MUST BE RUN IN EXPERIMENTAL BASIC"
2190 PRINT
2200 PRINT
2210 GO SUB 2420' TO READ DATA
2220 PRINT"INLET STATE","OPT ALLOCATION","MAX RETURN","OPT OUTLET"
2230 PRINT"-----------","--------------","------------","----------"
2240 FOR K = K1 TO 1 STEP -1' SETS STAGE NUMBER, PROCEEDS BACKWARD
2250 FOR I = 1 TO I1' SETS ROW IN F(I,J),I.E. ELEMENT OF INLET STATE
2260'       VECTOR
2270 LET F(I,1) = FNH(I,K)' SETS VALUE FOR INLET STATE VARIABLE
2280 LET X1 = F(I,1)
2290 LET F(I,3) = FNE(X1,1E-6,X1,.0001,.001)
2300 LET F(I,4) = X' OPT VALUE OF OUTLET STATE, EVALUATED IN FNE
2310'        CONSIDERED AT THIS POINT
2320 LET F(I,2) = FNA(F(I,1),F(I,4))' OPT SOLV ALLOCATION
2330 NEXT I
2340 MAT G= F
2350 PRINT "STAGE NUMBER: ";K
2360 PRINT
2370 MAT PRINT F
2380 PRINT
2390 PRINT
2400 NEXT K
2410 GO TO 2540' PROGRAM END
2420' *****DATA*****
2430 READ I1' NUMBER OF VALUES IN STATE VARIABLE
2440 DATA 20
2450 READ K1' NUMBER OF STAGES
2460 DATA 3
2470 READ A1' SALE P[ICE OF EXTRACT
2480 DATA 100
2490 READ B1' COST OF SOLVENT
2500 DATA 5
2510 READ Q' FLOW RATE OF FEED
2520 DATA 1
2530 RETURN
2540 END