|
|
1.1 root 1: C PROGRAM MAIN
2: C
3: C ******** THIS DRIVER WAS SET UP FOR THE PRIME 400 AT
4: C ******** THE CENTER FOR SCIENTIFIC COMPUTATION AND INTERACTIVE GRAPHICS
5: C ******** DREXEL UNIVERSITY
6: C ******** PHILADELPHIA, PA 19104
7: C
8: C THIS PROGRAM AND ASSOCIATED SUBROUTINES SHOULD BE COMPILED
9: C USING THE LONG INTEGER OPTION.
10: C
11: INTEGER*2 I
12: LOGICAL L,CLOS$A
13: CALL MATLAB(0)
14: I4=4
15: DO 100 I=1,12
16: L = CLOS$A(I)
17: 100 CONTINUE
18: CALL EXIT
19: END
20: SUBROUTINE FILES(LUNIT,NAME,IOSTAT)
21: INTEGER LUNIT,NAME(32)
22: INTEGER*2 I1,I2,I3,I4,I5,I8,ITYPE,ICODE,IUNIT,I16
23: INTEGER*4 NAM8(2)
24: DATA I1,I3,I4,I5,I8,I16/1,3,4,5,8,16/
25: C
26: C SYSTEM DEPENDENT ROUTINE TO ALLOCATE FILES
27: C LUNIT = LOGICAL UNIT NUMBER
28: C NAME = FILE NAME, 1 CHARACTER PER WORD
29: C
30: C MODIFY SUBROUTINE MATLAB SO THAT RTE = 1 AND WTE = 1
31: C ATTACH UNIT 9 TO HELP FILE
32: C
33: L = -LUNIT
34: IF (LUNIT .LT. 0) REWIND L
35: IF (LUNIT .LT. 0) RETURN
36: C
37: ENCODE(8,10,NAM8) (NAME(I),I=1,8)
38: 10 FORMAT(8A1)
39: IUNIT = LUNIT-4
40: CALL SRCH$$(I4,NAM8,0,IUNIT,ITYPE,ICODE)
41: CALL SRCH$$(I3,NAM8,I8,IUNIT,ITYPE,ICODE)
42: RETURN
43: END
44: SUBROUTINE SAVLOD(LUNIT,ID,M,N,IMG,JOB,XREAL,XIMAG)
45: INTEGER LUNIT,ID(4),M,N,IMG,JOB
46: DOUBLE PRECISION XREAL(1),XIMAG(1)
47: C
48: C IMPLEMENT SAVE AND LOAD
49: C LUNIT = LOGICAL UNIT NUMBER
50: C ID = NAME, FORMAT 4A1
51: C M, N = DIMENSIONS
52: C IMG = NONZERO IF XIMAG IS NONZERO
53: C JOB = 0 FOR SAVE
54: C = SPACE AVAILABLE FOR LOAD
55: C XREAL, XIMAG = REAL AND OPTIONAL IMAGINARY PARTS
56: C
57: C SYSTEM DEPENDENT FORMATS
58: 101 FORMAT(4A1,3I4)
59: 102 FORMAT(4D25.18)
60: C
61: IF (JOB .GT. 0) GO TO 20
62: C
63: C SAVE
64: 10 WRITE(LUNIT,101) ID,M,N,IMG
65: DO 15 J = 1, N
66: K = (J-1)*M+1
67: L = J*M
68: WRITE(LUNIT,102) (XREAL(I),I=K,L)
69: IF (IMG .NE. 0) WRITE(LUNIT,102) (XIMAG(I),I=K,L)
70: 15 CONTINUE
71: RETURN
72: C
73: C LOAD
74: 20 READ(LUNIT,101,END=30) ID,M,N,IMG
75: IF (M*N .GT. JOB) GO TO 30
76: DO 25 J = 1, N
77: K = (J-1)*M+1
78: L = J*M
79: READ(LUNIT,102,END=30) (XREAL(I),I=K,L)
80: IF (IMG .NE. 0) READ(LUNIT,102,END=30) (XIMAG(I),I=K,L)
81: 25 CONTINUE
82: RETURN
83: C
84: C END OF FILE
85: 30 M = 0
86: N = 0
87: RETURN
88: END
89: SUBROUTINE FORMZ(LUNIT,X,Y)
90: DOUBLE PRECISION X,Y,XHEX(2),YHEX(2)
91: INTEGER*2 I4
92: DATA I4/4/
93: C
94: C SYSTEM DEPENDENT ROUTINE TO PRINT WITH Z FORMAT
95: C
96: CALL HEXCHR(X,I4,XHEX)
97: IF(Y.EQ.0.D0)GO TO 20
98: CALL HEXCHR(Y,I4,YHEX)
99: WRITE(LUNIT,10)XHEX,YHEX
100: RETURN
101: 20 WRITE(LUNIT,10)XHEX
102: RETURN
103: 10 FORMAT(2A8,4X,2A8)
104: END
105: SUBROUTINE HEXCHR(VALUE,VALWDS,STRING)
106: C* VALUE IS THE VALUE TO BE CONVERTED TO HEXADECIMAL. ANY NUMBER OF WORDS.
107: C* VALWDS IS THE NUMBER OF WORDS IN THE VALUE TO BE TRANSLATED:
108: C*
109: C* INTEGER*2 -- 1
110: C* INTEGER*4 -- 2
111: C* REAL*4 -- 2
112: C* REAL*8 -- 4
113: C*
114: C**===> ETC.
115: C* STRING IS AN ARRAY INTO WHICH THE CHARACTERS WILL BE PACKED.
116: C* (IT DOESN'T MATTER WHAT KIND OF ARRAY IT IS).
117: C*
118: C*
119: C* WARNING: IF THE STRING PROVIDED IS NOT LONG ENOUGH (IT HAS TO BE TWICE
120: C* AS LONG IN TERMS OF BITS OR WORDS AS THE VALUE ARRAY), THERE
121: C* IS NO WAY TO DETECT IF IT RUNS OFF AND DESTROYS OTHER MEMORY.
122: C* (THIS IS PRIME STANDARD, BY THE WAY.)
123: C* ALSO, ANYTHING BEYOND THE 2*VALWDS ELEMENT OF THE STRING ARRAY
124: C* IS MEANINGLESS.
125: C*
126: C*
127:
128:
129:
130: INTEGER*2 VALUE(1),STRING(1),VALWDS,ONE,TWO,THREE,FOUR,I
131: DO 10 I=1,VALWDS
132: STRING(2*I-1) = 0
133: STRING(2*I) = 0
134: ONE = RS(LT(VALUE(I),4),12)
135: TWO = RS(LT(RT(VALUE(I),12),8),8)
136: THREE = RS(LT(RT(VALUE(I),8),12),4)
137: FOUR = RT(VALUE(I),4)
138: IF (ONE .LE. :011) STRING(2*I-1) =
139: + OR(LS(OR(:260,ONE),8),STRING(2*I-1))
140: IF (ONE .GT. :011) STRING(2*I-1) =
141: + OR(LS(ONE + :267,8),STRING(2*I-1))
142: IF (TWO .LE. :011) STRING(2*I-1) =
143: + OR(TWO,:260,STRING(2*I-1))
144: IF (TWO .GT. :011) STRING(2*I-1) =
145: + OR(TWO+ :267,STRING(2*I-1))
146: IF (THREE .LE. :011) STRING(2*I) =
147: + OR(LS(OR(:260,THREE),8),STRING(2*I))
148: IF (THREE .GT. :011) STRING(2*I) =
149: + OR(LS(THREE + :267,8),STRING(2*I))
150: IF (FOUR .LE. :011) STRING(2*I) =
151: + OR(FOUR,:260,STRING(2*I))
152: IF (FOUR .GT. :011) STRING(2*I) =
153: + OR(FOUR+ :267,STRING(2*I))
154: 10 CONTINUE
155: RETURN
156: END
157: DOUBLE PRECISION FUNCTION FLOP(X)
158: DOUBLE PRECISION X
159: C SYSTEM DEPENDENT FUNCTION
160: C COUNT AND POSSIBLY CHOP EACH FLOATING POINT OPERATION
161: C FLP(1) IS FLOP COUNTER
162: C FLP(2) IS NUMBER OF PLACES TO BE CHOPPED
163: C
164: INTEGER SYM,SYN(4),BUF(256),CHAR,FLP(2),FIN,FUN,LHS,RHS,RAN(2)
165: COMMON /COM/ SYM,SYN,BUF,CHAR,FLP,FIN,FUN,LHS,RHS,RAN
166: C
167: DOUBLE PRECISION MASK(11),XX,MM
168: INTEGER*4 MAS(2,11),LX(2),LM(2)
169: C LOGICAL LX(4),LM(4)
170: EQUIVALENCE (LX(1),XX),(LM(1),MM),(MASK,MAS)
171: C
172: C *** THESE ARE THE MASKS TO CHOP HEX DIGITS ON THE PRIME 400 ***
173: C
174: DATA MAS /-:00000000001,-:00003600001,
175: $ -:00000000001,-:00077600001,
176: $ -:00000000001,-:01777600001,
177: $ -:00000000001,-:37777600001,
178: $ -:00000000020,:00000177777,
179: $ -:00000000400,:00000177777,
180: $ -:00000010000,:00000177777,
181: $ -:00000200000,:00000177777,
182: $ -:00004000000,:00000177777,
183: $ -:00100000000,:00000177777,
184: $ -:02000000000,:00000177777/
185: FLP(1) = FLP(1) + 1
186: K = FLP(2)
187: FLOP = X
188: IF (K .LE. 0) RETURN
189: FLOP = 0.0D0
190: IF (K .GE. 12) RETURN
191: XX = X
192: MM = MASK(K)
193: LX(1) = AND(LX(1),LM(1))
194: LX(2) = AND(LX(2),LM(2))
195: FLOP = XX
196: RETURN
197: END
198: SUBROUTINE XCHAR(BUF,K)
199: INTEGER BUF(1),K
200: C
201: C SYSTEM DEPENDENT ROUTINE TO HANDLE SPECIAL CHARACTERS
202: C
203: WRITE(6,10) BUF(1)
204: 10 FORMAT(1X,A1,' is not a MATLAB character.')
205: RETURN
206: END
207: SUBROUTINE USER(A,M,N,S,T)
208: DOUBLE PRECISION A(M,N),S,T
209: C
210: INTEGER A3(9)
211: DATA A3 /-149,537,-27,-50,180,-9,-154,546,-25/
212: IF (A(1,1) .NE. 3.0D0) RETURN
213: DO 10 I = 1, 9
214: A(I,1) = A3(I)
215: 10 CONTINUE
216: M = 3
217: N = 3
218: RETURN
219: END
220: SUBROUTINE PROMPT(PAUSE)
221: INTEGER PAUSE
222: C
223: C ISSUE MATLAB PROMPT WITH OPTIONAL PAUSE
224: C
225: INTEGER DDT,ERR,FMT,LCT(4),LIN(1024),LPT(6),RIO,WIO,RTE,WTE,HIO
226: COMMON /IOP/ DDT,ERR,FMT,LCT,LIN,LPT,RIO,WIO,RTE,WTE,HIO
227: WRITE(WTE,10)
228: IF (WIO .NE. 0) WRITE(WIO,10)
229: 10 FORMAT(/1X,'<>')
230: IF (PAUSE .EQ. 1) READ(RTE,20) DUMMY
231: 20 FORMAT(A1)
232: RETURN
233: END
234: SUBROUTINE PLOT(LUNIT,X,Y,N,P,K,BUF)
235: DOUBLE PRECISION X(N),Y(N),P(1)
236: INTEGER BUF(79)
237: C
238: C PLOT X VS. Y ON LUNIT
239: C IF K IS NONZERO, THEN P(1),...,P(K) ARE EXTRA PARAMETERS
240: C BUF IS WORK SPACE
241: C
242: DOUBLE PRECISION XMIN,YMIN,XMAX,YMAX,DY,DX,Y1,Y0
243: INTEGER AST,BLANK,H,W
244: DATA AST/1H*/,BLANK/1H /,H/20/,W/79/
245: C
246: C H = HEIGHT, W = WIDTH
247: C
248: XMIN = X(1)
249: XMAX = X(1)
250: YMIN = Y(1)
251: YMAX = Y(1)
252: DO 10 I = 1, N
253: XMIN = DMIN1(XMIN,X(I))
254: XMAX = DMAX1(XMAX,X(I))
255: YMIN = DMIN1(YMIN,Y(I))
256: YMAX = DMAX1(YMAX,Y(I))
257: 10 CONTINUE
258: DX = XMAX - XMIN
259: IF (DX .EQ. 0.0D0) DX = 1.0D0
260: DY = YMAX - YMIN
261: WRITE(LUNIT,35)
262: DO 40 L = 1, H
263: DO 20 J = 1, W
264: BUF(J) = BLANK
265: 20 CONTINUE
266: Y1 = YMIN + (H-L+1)*DY/H
267: Y0 = YMIN + (H-L)*DY/H
268: JMAX = 1
269: DO 30 I = 1, N
270: IF (Y(I) .GT. Y1) GO TO 30
271: IF (L.NE.H .AND. Y(I).LE.Y0) GO TO 30
272: J = 1 + (W-1)*(X(I) - XMIN)/DX
273: BUF(J) = AST
274: JMAX = MAX0(JMAX,J)
275: 30 CONTINUE
276: WRITE(LUNIT,35) (BUF(J),J=1,JMAX)
277: 35 FORMAT(1X,79A1)
278: 40 CONTINUE
279: RETURN
280: END
281: SUBROUTINE EDIT(BUF,N)
282: INTEGER BUF(N)
283: C
284: C CALLED AFTER INPUT OF A SINGLE BACKSLASH
285: C BUF CONTAINS PREVIOUS INPUT LINE, ONE CHAR PER WORD
286: C ENTER LOCAL EDITOR IF AVAILABLE
287: C OTHERWISE JUST
288: RETURN
289: END
290: DOUBLE PRECISION FUNCTION DSINH(X)
291: IMPLICIT DOUBLE PRECISION (A-Z)
292: DATA A0/-.139005324307231D6/,
293: $ A1/-.180697428204813D5/,
294: $ A2/-.404150536907816D3/,
295: $ B0/-.139005324307231D6/,
296: $ B1/+.509781123072247D4/,
297: $ B2/-.954080394016126D2/
298: IF (DABS(X) .LT. 0.5D0) GO TO 10
299: E = DEXP(X)
300: DSINH = 0.5D0*(E-1.0D0/E)
301: RETURN
302: 10 Y = X*X
303: P = ((A2*Y)+A1)*Y+A0
304: Q = ((Y+B2)*Y+B1)*Y+B0
305: DSINH = X*P/Q
306: RETURN
307: END
308: DOUBLE PRECISION FUNCTION DCOSH(X)
309: IMPLICIT DOUBLE PRECISION (A-Z)
310: E = DEXP(X)
311: DCOSH = 0.5D0*(E+1.0D0/E)
312: RETURN
313: END
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.