TITLE 'MAT-B00,08/22/73,DWG702985'
PAGE
*
*
* E X T E R N A L C O M M U N I C A T I O N
*
*
* DEFINITIONS
*
DEF DMATDIV DYADIC MATRIX DIVIDE OP ROUTINE
DEF MMATINV MONADIC MATRIX INVERT OP ROUTINE
DEF MAT@ START OF PROCEDURE
*
* REFERENCES
*
REF ALOCHNW ALLOCATE HEADER AND N WORDS
REF ALOCRS ALLOCATE RESULT DATA BLOCK
REF DREF DE-REF
REF DTYPEF DYADIC TYPE CHECK: FLOT RESULT
REF DXRETURN DYADIC EXECUTION DRIVER RETURN
REF ERLENGTH LENGTH ERROR
REF ERRANK RANK ERROR
REF ERSING SINGULAR MATRIX ERROR EXIT
REF EXECUTE EXECUTE XSEG
REF FLOT0 FLOATING POINT 0.0
REF FLOT1 FLOATING POINT 1.0
REF FSQRT F SQUARE ROOT EVAL
REF GENLOAD GEN LOAD BY RSTYPE
REF GIVEBACK SHRINK DATA BLOCK
REF GXSEGINI GEN XSEG INITIALIZATION
REF LFARG LEFT ARG PNTR
REF LOOPLOC LOOP LOCATION
REF MTYPEF MONADIC TYPE CHECK: FLOT RESULT
REF MXRETURN MONADIC EXECUTION DRIVER RETURN
REF OPBREAK OP BREAK HANDLER
REF RESULT RESULT PNTR
REF RETURN RETURN ADR CELL
REF RSADR RESULT ADDRESS
REF RSRANK RESULT RANK
REF RSSIZE RESULT SIZE
REF RTARG RIGHT ARG PNTR
REF SETADR SET UP ARG ADR CELL
REF XSEGBASE XSEG BASE
REF XSEGBRK XSEG BREAK FLAG
REF XSETUP SET ADR FOR INDEXED ACCESS
PAGE
*
*
* A S S E M B L Y P A R A M E T E R S
*
*
SYSTEM SIG5F
PROGSECT CSECT 1
*
* REGISTERS
*
N EQU 1 XSEG EXECUTION-TIME INDEX
N1 EQU 4 XSEG EXECUTION REG
X EQU 1 GENERAL INDEX REG
K EQU 2 SUBSCRIPT REG
T EQU 2 TYPE
XL EQU 3 XSEG LOC
A EQU 4 ARG ADR/INDEX
LX EQU 5 INDEX LINK REG
LX7 EQU 7 INDEX LINK REG
GIJ EQU 1 INDEX TO ACCESS G(I,J)
GIK EQU 2 INDEX TO ACCESS G(I,K)
GKI EQU 3 INDEX TO ACCESS G(K,I)
VI EQU 2 INDEX TO ACCESS V(I)
RI EQU 2 INDEX TO ACCESS R(I)
FJ EQU 3 INDEX TO ACCESS F(J)
FI EQU 3 INDEX TO ACCESS F(I)
X0 EQU 2 EVEN/ODD INDEX PAIR
X1 EQU 3 *
AIJ EQU 4 INDEX TO ACCESS A(I,K)
AF EQU 6 ACCUM FOR FLOT VALUES
AF1 EQU 7 *
BF EQU 8 2ND ACCUM FOR FLOT VALUES
BF1 EQU 9 *
BUF EQU 7 BUFFER FOR MOVING DATA/CODE GROUPS
R EQU 8 GENERAL WORK REG
S EQU 11 SIZE
NI EQU 12 I-LOOP COUNT
NJ EQU 13 J-LOOP COUNT
L3 EQU 12 LINK REG
L2 EQU 13 LINK REG
L1 EQU 14 LINK REG
PAGE
*
*
* P R O C S
*
*
* (THIS MODULE RESIDES IN AN OVERLAY, SO THE USUAL
* TEMP ALLOCATION PROCEDURE WILL NOT WORK).
*
TLOC SET 0
MATTEMPS EQU XSEGBASE+50
*
TEMP CNAME 1
DTEMP CNAME 2
PROC
DO1 NAME=2
TLOC SET TLOC+(TLOC&1)
DISP TLOC+50
LF EQU MATTEMPS+TLOC
TLOC SET TLOC+NAME
PEND
*
*
EVEN CNAME 0
ODD CNAME 1
PROC
LF EQU %
ERROR,1,(CF(2)+NAME)&1 'REGISTER HAS WRONG PARITY'
PEND
*
*
EQUAL CNAME
PROC
LF EQU %
ERROR,1,1-(CF(2)=CF(3)) 'REGISTERS MUST BE EQUAL'
PEND
*
*
NB CNAME X'680'
NBGE CNAME X'681'
NBLE CNAME X'682'
NBE CNAME X'683'
NBL CNAME X'691'
NBG CNAME X'692'
NBNE CNAME X'693'
PROC
ERROR,1,(AF>=0)+(NUM(AF)>1) 'AF MUST BE NEG CONST ADR'
LF GEN,12,20 NAME-1,AF
PEND
PAGE
*
*
* XSEG GEN PROCS
*
*
OPEN GEN
GEN CNAME
OPEN M,N,MN,I
PROC
LF EQU %
ERROR,1,1-(NUM(CF)=3) 'WRONG NUMBER OF CF ARGS'
M SET CF(2)
N SET CF(3)
MN SET M+N
ERROR,1,1-(NUM(AF)=(M>0)+(N>0)) 'WRONG NUMBER OF AF ARGS'
DO M>0
I DO N
LW,AF(1)+M+I-1 AF(2)+I-1
FIN
I DO MN*(MN<3)
STW,AF(1)+I-1 I-1,XL
ELSE
LCI MN
STM,AF(1) 0,XL
FIN
ELSE
I DO MN*(MN<3)
LW,BUF AF(1)+I-1
STW,BUF I-1,XL
ELSE
LCI N
LM,BUF AF(1)
STM,BUF 0,XL
FIN
FIN
AI,XL MN
PEND
CLOSE M,N,MN,I
PAGE
*
*
* M A T R I X I N V E R T / D I V I D E O P R O U T I N E S
*
*
USECT PROGSECT
MAT@ RES 0 START OF PROCEDURE
*
*
MMATINV EQU % MONADIC MATRIX INVERT
BAL,L1 MTYPEF ALLOW NUMERIC ARG; SET FLOT RESULT
LW,A RTARG GET ARG'S 1ST DIMEN (M)
LW,R 2,A CONTINUE LIKE MATRIX DIVIDE,
B MATOP1 PRETENDING THERE WAS AN M-BY-M
* LEFT ARG.
*
*
DMATDIV EQU % DYADIC MATRIX INVERT
BAL,L1 DTYPEF ALLOW NUMERIC ARGS; SET FLOT RESULT
LI,X 1
LB,R *LFARG,X GET RANK OF LEFT ARG
LW,A LFARG
CLM,R 1OR2 MUST BE 1 OR 2 (VECTOR OR MATRIX)
BCS,9 ERRANK
BCR,3 1Z1 IF RANK =1, SET P=1
LW,R 3,A RANK =2, P = 2'ND DIMEN
BEZ ERLENGTH
MATOP1 EQU %
1Z1 STW,R PVAL
LW,R 2,A M = 1'ST DIMEN (LFARG IS M-BY-P)
STW,R MVAL
LW,A RTARG
LI,X 1
LB,X *RTARG,X GET RIGHT ARG RANK
CI,X 2
BNE ERRANK MUST BE 2 (MATRIX)
STW,X RSRANK SET RESULT TO MATRIX
CW,R 2,A RIGHT ARG MUST BE M-BY-N
BNE ERLENGTH
LW,S 3,A WITH 0<N<=M.
BEZ ERLENGTH
CW,S MVAL
BG ERLENGTH
STW,S NVAL
AW,S PVAL COMPUTE ROW SIZE OF G
STW,S ROWSIZE (M BY N+P).
ODD,S
MW,S MVAL COMPUTE SIZE OF G
* NOTE: THIS SIZE IS >0 SINCE
* M*(N+P) > M*N >= N*N > 0.
STW,S RSSIZE WE'LL ALOCATE THE G DATA BLOCK
BAL,L1 ALOCRS IN RESULT POSITION. AT THE MOMENT,
* IT'S BIGGER THAN THE FINAL RESULT;
* WE'LL PARE IT DOWN LATER.
LW,R NVAL
STW,R 2,A SET RESULT DIMENS:
LW,R PVAL N-BY-P.
STW,R 3,A
BDR,R 1Z5 IF P=1,
LI,X 1 REDUCE RESULT RANK TO 1
STB,X *RESULT,X (2ND DIMEN BECOMES GAP WORD).
1Z5 LW,A LFARG
BEZ 1Z2 IF LFARG EXISTS (IE., MATRIX DIV),
LI,A 0 SET ITS ADR FOR INDEXED LOAD.
BAL,LX XSETUP
1Z2 LI,A 1 SET RTADR FOR INDEXED LOAD
BAL,LX XSETUP
LI,X -1 SET RSADR FOR SEQUENCIAL STORE
BAL,LX SETADR
BAL,LX GXSEGINI INIT XSEG CODE (RESULT ISN'T NULL)
GEN,0,3 CODE1 GEN XSEG TO BUILD G;
LI,A 1 EITHER G = B,A (CATENATED
BAL,L1 GENLOAD OR G = B,I COLUMN-WISE)
GEN,0,6 CODE2 WHERE A = LFARG,
LW,R RSADR B = RTARG,
AWM,R -6,XL I = IDENTITY MATRIX
STW,XL LOOPLOC
LW,A LFARG TEST DIVIDE/INVERT
BNEZ 1Z3
GEN,0,8 CODE3 INVERT: GEN CODE TO BUILD IDENTITY
AWM,XL -5,XL MATRIX, ROW BY ROW.
AWM,XL -7,XL
LW,R PVAL
STW,R TESTVAL
B 1Z4
1Z3 LI,A 0 DIVIDE: GEN CODE TO APPEND ROWS
BAL,L1 GENLOAD OF LFARG.
GEN,0,3 CODE5
1Z4 LW,R RSADR
AWM,R -3,XL
B EXECUTE EXECUTE XSEG
*
*
* XSEG CODE:
*
* LCW,N RSSIZE 0 INIT STORE INDEX / RESULT COUNT
CODE1 LI,K 0 1 INIT LEFT INDEX
STW,K KTEMP 2 INIT RIGHT INDEX
LW,N1 NVAL 3 INIT LEFT LOOP COUNT
* LOAD/CNV RTARG(K) 4 LEFT LOOP: MOVE A ROW OF 'B'
CODE2 STD,AF 0 RESULT(N) TO RESULT.
AI,N 1 BUMP STORE AND
AI,K 1 LEFT INDECES.
BDR,N1 XSEGBASE+4 COUNT LEFT LOOP
XW,K KTEMP SAVE LEFT INDEX, GET RIGHT INDEX
LW,N1 PVAL INIT RIGHT LOOP COUNT
*
* FOR INVERT:
*
*LOOPLOC EQU %
CODE3 CW,N1 TESTVAL RIGHT LOOP: COMPUTE A ROW OF 'I'
NBNE -4 (%+3)
LD,AF FLOT1 ONE ELEMENT IS 1.0,
NB -3 (%+2)
LD,AF FLOT0 ALL OTHERS ARE 0.0.
STD,AF 0 RESULT(N) MOVE TO RESULT (G)
BIR,N CODE4 BUMP RESULT INDEX
B XSEGDONE EXIT WHEN G FILLED
CODE4 BDR,N1 *LOOPLOC * COUNT RIGHT LOOP
MTW,-1 TESTVAL * END OF ROW, ADVANCE POSITION OF 1.0
XW,K KTEMP * SAVE RIGHT INDEX (UNUSED), GET LEFT
B XSEGBASE+3 * GO DO LEFT LOOP
*
* FOR DIVIDE:
*
*LOOPLOC LOAD/CNV LFARG(K) RIGHT LOOP: MOVE A ROW OF 'A'
CODE5 STD,AF 0 RESULT(N) TO RESULT.
BIR,N CODE6 BUMP STORE INDEX
B XSEGDONE EXIT WHEN G FILLED
CODE6 AI,K 1 * BUMP RIGHT INDEX
BDR,N1 *LOOPLOC * COUNT RIGHT LOOP
XW,K KTEMP * SAVE RIGHT INDEX, GET LEFT
B XSEGBASE+3 * GO DO LEFT LOOP
PAGE
*
*
XSEGDONE EQU % WE NOW HAVE OUR INITIAL G MATRIX
LI,A 0 DON'T ALLOW BREAKS DURING ALOC'S
STW,A XSEGBRK AND DREF'S. (FLAG WAS SET
* BY 'EXECUTE').
XW,A LFARG
BEZ 2Z1 IF LFARG EXITS,
BAL,LX7 DREF GET RID OF IT.
LI,A 0
2Z1 XW,A RTARG GET RID OF RTARG
BAL,LX7 DREF
LW,S MVAL ALOC DB FOR ROW-MAX VECTOR
SLS,S 1 R(1),...,R(M).
BAL,LX7 ALOCHNW
STW,A LFARG
LI,R X'0100' SET IT TO LEGIT TYPE
STH,R *LFARG
LW,S NVAL ALOC DB FOR SCALING FACTOR VECTOR,
SLS,S 1 F(1),...,F(N).
BAL,LX7 ALOCHNW
STW,A RTARG
LI,R X'0100' SET IT TO LEGIT TYPE
STH,R *RTARG
LI,R OPBREAK NOW THAT ALL MEMORY FIDDLING IS
STW,R XSEGBRK DONE, IT'S SAFE TO ALLOW BREAKS.
AI,A 2 GET PAST RTARG HEADER
SLS,A -1
STW,A F1ADR SET DBLWORD ADR OF F(1)
LW,A LFARG
AI,A 2 GET PAST LFARG HEADER
STW,A V1ADR SET WORD ADR OF V(1) (INTG.VECTOR)
SLS,A -1
STW,A R1ADR SET DBLWORD ADR OF R(1)
LW,A RESULT
AI,A 4 GET PAST RESULT HEADER AND DIMENS
SLS,A -1
STW,A G11ADR SET DBLWORD ADR OF G(1,1)
PAGE
*
*
* DO SCALING
*
* ROW SCALING: FOR I=1 TO M,
* SET R(I) = MAX OF ABS(G(I,J)) FOR J=1 TO N
*
LW,RI R1ADR
LW,GIJ G11ADR
LW,NI MVAL I LOOP COUNT = M
LD,AF FLOT0 INIT BIGGEST ELMT = 0.0
3Z1 STD,AF BIGGEST REMEMBER BIGGEST ELEMENT
LD,AF FLOT0 I LOOP: INIT R(I)=0.0
STD,AF 0,RI
LW,NJ NVAL J LOOP COUNT = N
3Z2 LAD,AF 0,GIJ J LOOP
CD,AF 0,RI SET R(I)=MAX(R(I),ABS(G(I,J)))
BLE 3Z3
STD,AF 0,RI
3Z3 AI,GIJ 1 MOVE GIJ ALONG ROW
BDR,NJ 3Z2 COUNT J LOOP
LD,AF BIGGEST IS THIS STILL THE BIGGEST ?
CD,AF 0,RI IF NOT,
BGE 3Z31 UPDATE IT.
LD,AF 0,RI
3Z31 AW,GIJ PVAL MOVE ACCROSS LAST P COLS, TO NXT ROW
AI,RI 1 MOVE RI TO NEXT ELMT
BDR,NI 3Z1 COUNT I LOOP
FML,AF EPSILON SET PIVOT TEST VALUE
STD,AF PIVOTMIN = BIGGEST ELEMENT * EPSILON.
*
* COLUMN SCALING: FOR J=1 TO N,
* SET F(J) = MAX OF ABS(G(I,J))/R(I)
* FOR I = 1 TO M.
*
LW,FJ F1ADR
LW,GIJ G11ADR
LW,NJ NVAL J LOOP COUNT = N
3Z4 STW,GIJ G1JADR J LOOP: REMEMBER ADR OF G(1,J)
LD,AF FLOT0
STD,AF 0,FJ INIT F(J)=0
LW,RI R1ADR INIT R TO R(1)
LW,NI MVAL I LOOP COUNT = M
3Z5 LAD,AF 0,GIJ I LOOP
FDL,AF 0,RI A = ABS(G(I,J))/R(I)
CD,AF 0,FJ SET F(J)= MAX(F(J),A)
BLE 3Z6
STD,AF 0,FJ
3Z6 AW,GIJ ROWSIZE BUMP G DOWN COLUMN
AI,RI 1 BUMP R
BDR,NI 3Z5 COUNT I LOOP
AI,FJ 1 BUMP TO NEXT F(J)
LW,GIJ G1JADR RECALL ADR OF G(1,J)
AI,GIJ 1 BUMP TO NEXT COL
BDR,NJ 3Z4 COUNT J LOOP
*
* APPLY SCALE FACTORS TO MATRIX COLUMNS:
* FOR I=1 TO M, AND J=1 TO N,
* SET G(I,J) = G(I,J)/F(J)
*
LW,GIJ G11ADR
LW,NI MVAL I LOOP COUNT = M
3Z7 LW,FJ F1ADR
LW,NJ NVAL J LOOP COUNT = N
3Z8 LD,AF 0,GIJ
FDL,AF 0,FJ SET G(I,J) = G(I,J)/F(J)
STD,AF 0,GIJ
AI,GIJ 1 BUMP G ALONG ROW
AI,FJ 1 BUMP TO NEXT F
BDR,NJ 3Z8 COUNT J LOOP
AW,GIJ PVAL BUMP G PAST LAST P COLS, TO NXT ROW
BDR,NI 3Z7 COUNT I LOOP
*
* INITIALIZE COLUMN PERMUTATION VECTOR:
* SET V(I) = I-1 FOR I = 1 TO N.
*
LW,X NVAL I LOOP COUNT = N
3Z9 STW,X *V1ADR,X V(I+1)=I
BDR,X 3Z9 COUNT I LOOP
STW,X *V1ADR V(1)=0
PAGE
*
*
* PERFORM HOUSEHOLDER TRANSFORMATIONS,
* LEADING TO UPPER TRIANGULAR SYSTEM.
*
* INITIALIZE K LOOP PARAMETERS
* K RUNS FROM 1 TO N.
*
LW,GIJ G11ADR SET INIT (K=1) VALUES FOR:
STW,GIJ GK1ADR ADR(G(K,1)),
STW,GIJ GKKADR ADR(G(K,K)),
LW,R MVAL
STW,R MK1VAL M-K+1,
LW,R NVAL
STW,R NK1VAL N-K+1,
LI,R 0
STW,R K1VAL K-1.
*
* SELECT PIVOT ELEMENT: FIND S, T SUCH THAT
* ABS(G(S,T)) = MAX(ABS(G(I,J)))
* FOR I = K TO M, J = K TO N.
*
KLOOP LW,GIJ GK1ADR START G AT G(K,1)
LD,AF FLOT0 INIT ABS(G(S,T)) = 0
STD,AF ABSGST
LW,NI MK1VAL I LOOP COUNT = M-K+1
4Z2 AW,GIJ K1VAL I LOOP: SET G TO G(I,K)
LW,NJ NK1VAL J LOOP COUNT = N-K+1
4Z3 LAD,AF 0,GIJ J LOOP:
CD,AF ABSGST IF ABS(G(I,J)) > ABSGST,
BLE 4Z4
STD,AF ABSGST SET ABSGST = ABS(G(I,J)),
STW,GIJ GSTADR AND REMEMBER S AND T.
4Z4 AI,GIJ 1 BUMP G ALONG ROW
BDR,NJ 4Z3 COUNT J LOOP
AW,GIJ PVAL BUMP G PAST LAST P COLS TO NEXT ROW
BDR,NI 4Z2 COUNT I LOOP
LD,AF ABSGST
CD,AF PIVOTMIN IF PIVOT ELMT DANGEROUSLY LOW,
BLE ERSING PROCLAIM THE MATRIX SINGULAR.
*
* PERMUTE COLUMNS K/T: FOR I = 1 TO M,
* EXCHANGE G(I,K) AND G(I,T),
* EXCHANGE V(K) AND V(T).
*
LW,X1 GSTADR GET ADR(G(S,T))
LI,X0 0
EVEN,X0 X0/X1 ARE AN EVEN/ODD PAIR
EQUAL,X0+1,X1
SW,X1 G11ADR SEPARATE S,T INDEX INTO DISTINCT
DW,X0 ROWSIZE ROW & COL INDECES:
STW,X0 T1VAL X0 = T-1,
* X1 = S-1.
CW,X0 K1VAL IF T=K, NO EXCHANGE NEEDED
BE 4Z6
LW,GIJ K1VAL GET K-1
LW,R *V1ADR,X0 EXCHANGE V(T) AND V(K)
XW,R *V1ADR,GIJ
STW,R *V1ADR,X0
AW,X0 G11ADR COMPUTE ADR(G(1,T)),
AW,GIJ G11ADR AND ADR(G(1,K)).
LW,NI MVAL I LOOP COUNT = M
4Z5 LD,AF 0,GIJ I LOOP
LD,BF 0,X0 EXCHANGE G(I,K) AND G(I,T)
STD,AF 0,X0
STD,BF 0,GIJ
AW,GIJ ROWSIZE BUMP G DOWN COLUMN
AW,X0 ROWSIZE *
BDR,NI 4Z5 COUNT I LOOP
4Z6 EQU %
*
* PERMUTE ROWS S/K: FOR J = 1 TO N+P,
* EXCHANGE G(S,J) AND G(K,J).
*
CW,X1 K1VAL IF S=K, NO EXCHANGE NECESSARY
BE 4Z8
LW,X1 GSTADR
SW,X1 T1VAL GET ADR(G(S,1)),
LW,GIJ GK1ADR AND ADR(G(K,1)).
LW,NJ ROWSIZE J LOOP COUNT = N+P
4Z7 LD,AF 0,GIJ J LOOP
LD,BF 0,X1 EXCHANGE G(S,J) AND G(K,J)
STD,AF 0,X1
STD,BF 0,GIJ
AI,GIJ 1 BUMP G ALONG ROW
AI,X1 1 *
BDR,NJ 4Z7 COUNT J LOOP
4Z8 EQU %
*
*
* APPLY K'TH TRANSFORMATION:
*
* S = SQRT(G(K,K)**2 + G(K+1,K)**2 +...+ G(M,K)**2)
* U = G(K,K) + S*SGN(G(K,K))
* B = 1/(S*S + S*ABS(G(K,K)))
* FOR J = K+1 TO N+P,
* Y = B*U*G(K,J)+B*(G(K+1,J)*G(K+1,K) +...
* ...+ G(M,J)*G(M,K))
* G(I,J):=G(I,J)-Y*G(I,K) FOR I=K+1 TO M
* G(K,J):=G(K,J)-U*Y
* G(K,K):=-S*SGN(G(K,K)) = G(K,K)-U
*
LD,AF FLOT0 INIT SUM (FOR S**2) TO 0.0
LW,NI MK1VAL I LOOP COUNT = M-K+1
LW,GIJ GKKADR START AT G(K,K)
5Z0 STD,AF S2VAL I LOOP
LD,AF 0,GIJ ADD IN G(I,K)**2
FML,AF 0,GIJ
FAL,AF S2VAL
AW,GIJ ROWSIZE BUMP G DOWN COLUMN
BDR,NI 5Z0 COUNT I LOOP
BAL,LX FSQRT COMPUTE S=SQRT(S**2)
STD,AF SVAL SAVE S
FAL,AF ABSGKK S+ABS(G(K,K))
FML,AF SVAL S**2 + S*ABS(G(K,K)) = 1/B
LD,BF FLOT1
FDL,BF AF B
STD,BF BVAL
LW,GIJ GKKADR GET ADR(G(K,K))
LD,AF 0,GIJ
BGEZ 5Z1 IF G(K,K)<0,
FSL,AF SVAL U = G(K,K)-S,
LD,BF SVAL G(K,K):=S.
B 5Z2 IF G(K,K)>0,
5Z1 FAL,AF SVAL U = G(K,K)+S,
LCD,BF SVAL G(K,K):=-S.
5Z2 STD,AF UVAL U = G(K,K)+S*SGN(G(K,K))
STD,BF 0,GIJ G(K,K) := -S*SGN(G(K,K))
LW,NJ NK1VAL
AW,NJ PVAL J LOOP COUNT = N+P-K
B 5Z8
5Z3 STW,GIJ GKJADR J LOOP: SET ADR(G(K,J)),
LW,GIK GKKADR AND ADR(G(K,K)).
LD,AF 0,GIJ START COMPUTATION OF Y:
FML,AF UVAL AF = U*G(K,J).
LW,NI MK1VAL I LOOP COUNT = M-K (POSSIBLY 0)
B 5Z5 START I LOOP WITH I=K+1
5Z4 STD,AF YVAL I LOOP:
LD,AF 0,GIK AF=AF+G(I,K)*G(I,J).
FML,AF 0,GIJ
FAL,AF YVAL
5Z5 AW,GIK ROWSIZE BUMP G DOWN COLUMN
AW,GIJ ROWSIZE
BDR,NI 5Z4 COUNT I LOOP
FML,AF BVAL Y = B*AF
STD,AF YVAL
LW,GIK GKKADR RESTORE INIT G ADR'S
LW,GIJ GKJADR
LCD,AF UVAL G(K,J):=G(K,J)-Y*U
LW,NI MK1VAL I LOOP COUNT = M-K (POSSIBLY 0)
B 5Z7 START I LOOP WITH I=K+1
5Z6 LCD,AF 0,GIK
5Z7 FML,AF YVAL G(I,J):=G(I,J)-Y*G(I,K)
FAL,AF 0,GIJ
STD,AF 0,GIJ
AW,GIK ROWSIZE BUMP G DOWN COLUMN
AW,GIJ ROWSIZE
BDR,NI 5Z6 COUNT I LOOP
LW,GIJ GKJADR UPDATE ADR(G(K,J)) BY
5Z8 AI,GIJ 1 BUMPING IT ALONG ROW K.
BDR,NJ 5Z3 COUNT J LOOP
*
*
* UPDATE K LOOP PARAMETERS
*
MTW,-1 NK1VAL UPDATE N-K+1 VALUE
* IF N-K+1=0 (K>N),
BEZ SOLVE WE'RE DONE.
MTW,-1 MK1VAL K<=N: UPDATE M-K+1,
MTW,1 K1VAL K-1,
LW,S ROWSIZE
AWM,S GK1ADR ADR(G(K,1)),
AI,S 1
AWM,S GKKADR AND ADR(G(K,K)).
B KLOOP CONTINUE K LOOP
PAGE
*
*
* SOLVE UPPER TRIANGULAR SYSTEM
*
*
* FOR K = N TO 1, AND J = N+1 TO N+P, SET
* G(K,J) := (G(K,J) - G(K,K+1)*G(K+1,J) - ...
* ... - G(K,N)*G(N,J))/G(K,K)
*
SOLVE MTW,1 NK1VAL RESTORE N-K+1 VALUE
LW,GIJ GKKADR START G(K,N) AT G(N,N)
6Z1 STW,GIJ GKNADR REMEMBER ADR(G(K,N))
AI,GIJ 1 START G(K,J) AT G(K,N+1)
LW,NJ PVAL J LOOP COUNT = P
6Z2 STW,GIJ GKJADR REMEMBER ADR(G(K,J))
LW,GKI GKKADR START G(K,I) AT G(K,K)
LD,AF 0,GIJ START SUM WITH G(K,J)
LW,NI NK1VAL I LOOP COUNT = N-K (POSSIBLY 0)
B 6Z4 START WITH I=K+1
6Z3 STD,AF SUM
LCD,AF 0,GKI SUM := SUM-G(K,I)*G(I,J)
FML,AF 0,GIJ
FAL,AF SUM
6Z4 AI,GKI 1 BUMP G(K,I) ALONG ROW
AW,GIJ ROWSIZE BUMP G(I,J) DOWN COLUMN
BDR,NI 6Z3 COUNT I LOOP
LW,GIJ GKKADR COMPUTE NEW G(K,J)
FDL,AF 0,GIJ = SUM/G(K,K).
LW,GIJ GKJADR
STD,AF 0,GIJ STORE NEW G(K,J)
AI,GIJ 1 UPDATE ADR(G(K,J))
BDR,NJ 6Z2 COUNT J LOOP
MTW,-1 K1VAL COUNT K LOOP (UPDATE K-1)
BLZ NORMALYZ IF K=0, WE'RE DONE
MTW,1 NK1VAL K>0: UPDATE N-K+1,
LCW,R ROWSIZE
AI,R -1
AWM,R GKKADR ADR(G(K,K)),
LW,GIJ GKNADR AND ADR(G(K,N)).
SW,GIJ ROWSIZE
B 6Z1
PAGE
*
*
* NORMALIZE THE SOLUTION TO CORRECT FOR
* COLUMN SCALING AND INTERCHANGES.
*
* FOR I = 1 TO N, AND J = 1 TO P,
* EXCHANGE G(I,N+J) WITH G(L,N+J),
* AND COPY V(I) TO V(L), WHERE
* V(L)=I.
*
NORMALYZ LW,AIJ G11ADR START 'A' AT A(1,1)
LW,NI NVAL I LOOP COUNT = N
LI,VI 0 INIT TO V(1)
7Z1 CW,VI *V1ADR,VI IF V(I)=I, NO EXCHANGE NEEDED
BNE 7Z10
AW,AIJ ROWSIZE (IN WHICH CASE, BUMP TO NEXT ROW)
B 7Z3
7Z10 LW,GIJ VI L:=I
7Z11 LW,GIJ *V1ADR,GIJ L:=V(L)
CW,VI *V1ADR,GIJ SEARCH V FOR I: IF I.NE.V(L),
BNE 7Z11 SET L:=V(L) AND TRY AGAIN.
LW,R *V1ADR,VI COPY V(I)
STW,R *V1ADR,GIJ TO V(L).
ODD,GIJ
MW,GIJ ROWSIZE COMPUTE ADR(G(L,N+1))
AW,GIJ NVAL
AW,GIJ G11ADR
AW,AIJ NVAL MOVE 'A' TO A(I,N+1)
LW,NJ PVAL J LOOP COUNT = P
7Z2 LD,AF 0,GIJ EXCHANGE G(L,N+J)
LD,BF 0,AIJ AND A(I,N+J).
STD,AF 0,AIJ
STD,BF 0,GIJ
AI,GIJ 1 BUMP G ALONG ROW
AI,AIJ 1 BUMP A ALONG ROW
BDR,NJ 7Z2 COUNT J LOOP
7Z3 AI,VI 1 BUMP TO NEXT V(I)
BDR,NI 7Z1 COUNT I LOOP
*
*
* FOR I = 1 TO N, AND J = 1 TO P,
* SET A(I,J) = G(I,N+J)/F(I).
* 'A' IS THE N-BY-P CORNER OF G IN
* WHICH THE FINAL RESULT APPEARS.
*
LW,AIJ G11ADR START A AT A(1,1)
LW,GIJ G11ADR START G AT G(1,1)
LW,FI F1ADR START F AT F(1)
LW,NI NVAL I LOOP COUNT = N
7Z4 AW,GIJ NVAL START G(I,J) AT G(I,N+1)
LW,NJ PVAL J LOOP COUNT = P
7Z5 LD,AF 0,GIJ GET G(I,N+J)
FDL,AF 0,FI DIVIDE BY F(I)
STD,AF 0,AIJ STORE AS A(I,J)
AI,GIJ 1 BUMP G AND A
AI,AIJ 1 ALONG ROW.
BDR,NJ 7Z5 COUNT J LOOP
AI,FI 1 BUMP TO NEXT F(I)
BDR,NI 7Z4 COUNT I LOOP
*
*
* NOW GIVE BACK THE M*(N+P)-N*P DOUBLEWORDS OF 'G' THAT
* ARE NO LONGER NEEDED.
*
LI,A 0 MUSTN'T ALLOW BREAKS DURING MEMORY
STW,A XSEGBRK MANAGEMENT OPERATIONS TO FOLLOW.
LW,S NVAL FINAL RESULT SIZE = N*P
MW,S PVAL
XW,S RSSIZE COMPUTE NR OF DBLWORDS TO GIVE BACK
SW,S RSSIZE = OLD G SIZE - RESULT SIZE.
SLS,S 1 CHANGE TO NR OF WORDS
LW,A RESULT
BAL,LX7 GIVEBACK GIVE THEM BACK
LI,A DXRETURN IF DYADIC ENTRY USED, SIMPLY EXIT
CW,A RETURN DYADICALLY; BOTH LFARG (V(I))
BE DXRETURN AND RTARG (F(J)) WILL BE DEREFFED.
LI,A 0 OTHERWISE, WE HAVE TO DEREF
XW,A LFARG LFARG OURSELVES.
BAL,LX7 DREF
B MXRETURN
PAGE
*
*
* DATA / TEMPS
*
*
BOUND 8
1OR2 DATA 1,2 RANGE 1 TO 2, FOR ARG RANK TEST
EPSILON DATA X'34100000',0 16**(-13) FOR SINGULARITY TEST
*
*
V1ADR TEMP WORD ADR OF V(1)
R1ADR TEMP DBLWORD ADR OF R(1)
F1ADR TEMP DBLWORD ADR OF F(1)
G11ADR TEMP DBLWORD ADR OF G(1,1)
GK1ADR TEMP DBLWORD ADR OF G(K,1)
GKKADR TEMP DBLWORD ADR OF G(K,K)
GKJADR TEMP DBLWORD ADR OF G(K,J)
GKNADR TEMP DBLWORD ADR OF G(K,N)
G1JADR TEMP DBLWORD ADR OF G(1,J)
GSTADR TEMP DBLWORD ADR OF G(S,T)
*
MVAL TEMP VALUE OF M
NVAL TEMP VALUE OF N
PVAL TEMP VALUE OF P
K1VAL TEMP VALUE OF K-1
T1VAL TEMP VALUE OF T-1
MK1VAL TEMP VALUE OF M-K+1
NK1VAL TEMP VALUE OF N-K+1
ROWSIZE TEMP VALUE OF N+P
TESTVAL TEMP INDEX OF NEXT DIAGONAL POS. (XSEG)
KTEMP TEMP LEFT/RIGHT ARG INDEX (XSEG)
*
ABSGST DTEMP VALUE OF ABS(G(S,T))
ABSGKK EQU ABSGST = ABS(G(K,K)) AFTER EXCHANGES
S2VAL DTEMP SUM(I=K,M) OF G(I,K)**2
SVAL EQU S2VAL SQRT OF ABOVE
BVAL DTEMP 1/(S**2 + S*ABS(G(K,K)))
UVAL DTEMP G(K,K) + S*SGN(G(K,K))
YVAL DTEMP B*U*G(K,J)+B*SUM(G(I,K)*G(I,J))
SUM DTEMP G(K,J)-SUM(G(K,I)*G(I,J))
PIVOTMIN DTEMP PIVOT TEST VALUE
BIGGEST DTEMP BIGGEST ELEMENT OF MATRIX
*
*
7Z END