      SKP
*
*********************************************
****    SUBROUTINE  MATRIX  MULT         ****
*********************************************
*ROUTINE IS OF FORM B3(M,P)=B1(M,N)*B2(N,P) *
*B1,B2,B3 ARE ADDRESSES OF BASE ADDRESSES OF*
*THREE MATRICES                             *
*********************************************
*
MULT  NOP
      JSB LCHK1     TEST B1,B2 FOR UNASSIGNED TERMS
*                                   CHECK DIMENSIONS
      LDA B3+1      PARAMETERS OF B3
      AND MSK0      SAVE COLUMN COUNT
      STA T6
      LDA B2+1      PARAMETERS OF B2
      AND MSK0
      CPA T6        COLUMNS EQUAL
      RSS             IN NUMBER?
      JMP LERR      NO
*                                   COMBINE B3,B2 PARAMETERS
*                                   INTO (M,N) AND COMPARE
*                                   WITH THOSE OF B1
      LDA B3+1      PARAMETERS OF B3
      AND M256
      STA 1         STORE ROW IN MSP OF B
      LDA B2+1      PARAMETERS OF B2
      ALF,ALF       GET ROW COUNT
      AND MSK0        IN (A)
      ADA 1         COMBINE A AND B
      LDB B1+1      PARAMETERS OF B1
      JSB COMPR     COMPARE ROW AND COL
*                                   DIMENSIONS ARE COMPATIBLE
*                                   M,N ARE STORED IN T3,T4
*                                   SAVE B2 AS DESTROYED IN
      LDA B2        MULT
      STA T5
      LDA T3        SET
      CMA,INA         ROW
      STA T9            COUNTER
MULT4 LDA T6
      CMA,INA
      STA T10       SET COUNTER
      LDA T5
      STA B2        RESTORE BASE ADDRESS B2
MULT3 CLA
      STA T11       COUNTER FOR B2. INCR BY
*                                   2*P AND POINTS TO NEXT TERM
*                                   IN COL TO BE MULTIPLIED
      STA T12       COUNTER FOR B1. INCR BY 2
*                                   AND POINTS TO NEXT TERM
*                                   IN ROW TO BE MULTIPLIED
      CLB
      JSB .DST      CLEAR TO ZERO
      DEF B3,I
MULT2 LDB B1        COMPUTE PROD OF ONE TERM
      ADB T12       IN ROW BY ONE TERM IN COL
      STB T18
      LDB B2
      ADB T11
      JSB .DLD
      DEF 1,I
      JSB .FMPA,I
      DEF T18,I
      JSB .FADA,I   COMPUTES RUNNING SUM
      DEF B3,I
      JSB .DST
      DEF B3,I
      ISZ T12       SELECT NEXT TERM IN ROW
      ISZ T12
      LDA T6        SELECT NEXT TERM IN COL
      ALS
      ADA T11
      STA T11
*                                   TEST IF HAVE MULT ONE ROW
*                                   BY ONE COLUMN
      LDA T4
      ALS
      CPA T12
      JMP *+2
      JMP MULT2     MULT AND ADD IN NEXT TERM
*                                   SUMMATION OF PRODUCTS FOR
*                                   ONE TERM OF B3 IS DONE
*                                   MULT SAME ROW BY NEXT COL
      ISZ B3        INCR RECEIVING MAT
      ISZ B3
      ISZ B2        BASE ADDRESS OF NEXT COL
      ISZ B2
*                                   TEST IF HAVE MULT ONE ROW
*                                   BY ALL COLUMNS
      ISZ T10       SKIP IF INNERPRODUCT DONE
      JMP MULT3     COMPUTE SAME ROW*NEXT COL
*                                   SELECT NEXT ROW
      LDA T4
      ALS
      ADA B1
      STA B1        ADDRESS OF NEXT ROW
      ISZ T9
      JMP MULT4     MULT ROW BY ALL COLUMNS
      JMP MULT,I    EXIT TO MAIN PROGRAM
      SKP
*
*********************************************
****   SUBROUTINE  MATRIX  INVERT        ****
*********************************************
*OPERATION OF FORM  MAT B3 = INV B1         *
*B1,B3 ARE ADDRESSES OF BASE ADDRESS OF     *
*MATRIX TO BE INVERTED AND RECEIVING MATRIX *
*RESPECTIVELY. B2 IS REDUNDANT. METHOD USED *
*IS GAUSSIAN ELIMINATION WITH COLUMN        *
*PIVOTING                                   *
*********************************************
*
LINV  NOP           SUBROUTINE MATRIX INVERT
      JSB LCHK2     TEST B1 FOR UNASSIGNED TERMS
      LDA B1+1      DIMENSIONS OF MATRIX B1
      LDB B3+1      DIMENSIONS OF MATRIX B3
      JSB COMPR     CHECK DIMENSIONS
*                                   ROW AND COL VALUES T3,T4
*                                   MAKE COPY OF MATRIX B1
*                                       IN FREE CORE
      LDA B3        SAVE
      STA T13         B3
      LDA B1+1      COMPUTE SIZE
      JSB MDIM        OF MATRIX
      CMA,INA
      ARS           SAVE
      STA T2          ELEMENT
      ALS               COUNTER
      LDB LSTPT     SAVE
      INB             ADDRESS OF
      STB B2            FREE CORE
      STB B3
      CMB,INB       COMPUTE SIZE OF
      ADB HSTPT       FREE CORE AREA
      ADA 1         ENOUGH
      SSA             CORE LEFT?
      JMP E1        NO
      JSB REPLC     YES, COPY SOURCE MATRIX
      LDA T13       RESTORE
      STA B3          B3
      JSB LIDN      SET DESTINATION TO IDENTITY
      LDA T13       RESTORE ITS
      STA B3          BASE ADDRESS
      CLA           INITIALIZE
      STA T12         MAXIMUM
      STA T13           ELEMENT
      LDA B2        COPY B2 INTO B1 AS
      STA B1        B2 NEEDED LATER
LIN11 LDA B1,I      LOAD
      ISZ B1          NEXT
      LDB B1,I          ELEMENT
      ISZ B1
      SSA           GET ABSOLUTE VALUE
      JSB ARINV        IF NUMBER IS NEGATIVE
      STA T18       SAVE NUMBER
      STB T19
      JSB .FSBA,I   SUBTRACT EXISTING MAX.
      DEF T12           VALUE
      SSA           SKIP AND SWAP IF POSITIVE
      JMP LIN10
      LDA T18       SWAP
      LDB T19
      STA T12
      STB T13
LIN10 ISZ T2        ALL ELEMENTS EXHAUSTED?
      JMP LIN11     NO
      LDA T12       COMPUTE RELATIVE TOLERANCE
      LDB T13       TOL=ABSOLUTE TOL * MAX VALUE
      JSB .FMPA,I
      DEF T16       ABSOLUTE TOLERANCE
      STA MLBX1     RELATIVE
      STB MLBX1+1     TOLERANCE
      CLA           INITIALIZE PIVOT
      STA LPIV
      ISZ T4        REQUIRE CONSTANT (ROW+1)
LINV1 ISZ LPIV      SELECT NEXT PIVOT
      LDA LPIV      TEST IF HAVE PROCESSED
      CPA T4        LAST PIVOT
      JMP LINV,I    NORMAL EXIT TO MAIN PROG
*                                   SCAN PIVOTAL COLUMN FOR
*                                   LARGEST ELEMENT
      LDA LPIV      COMPUTE ADDRESS OF PIVOT
      LDB LPIV      COLUMN USING ROUTINE LWHR
      STA T2        ROW COUNTER
      JSB LWHR      ON RETURN, ADDRESS IN A
      STA T1
      CLA
      STA T12       T12,T13 IS STORE
      STA T13       FOR GREATEST VALUE
LINV2 JSB .DLD      LOAD FP NUMBER
      DEF T1,I
      SSA           OBTAIN ABSOLUTE VALUE
      JSB ARINV      IF NUMBER IS NEGATIVE
      STA T18       STORE VALUE OF FP NUMBER
      STB T19
      JSB .FSBA,I   SUBTR EXISTING LARGEST VALUE
      DEF T12
      SSA           SKIP AND SWAP IF POSITIVE
      JMP LINV7     T2 STILL CONTAINS MAX VALUE
      LDA T18       STORE NEW MAX VALUE
      LDB T19
      STA T12
      STB T13
      LDA T2        SET T5 TO POSITION IN
      STA T5        COLUMN OF MAX VALUE
LINV7 ISZ T2
      LDA T2        TEST FOR LAST TERM IN COL
      CPA T4
      JMP LINV8     SWAP ROWS
      LDA T3        COMPUTE
      ALS           NEXT ADDRESS
      ADA T1             IN PIVOT
      STA T1              COLUMN
      JMP LINV2     SELECT NEXT TERM
*                                   SWAP ROWS LPIV AND T5
LINV8 LDA LPIV      COMPUTE ADDRESS
      CLB,INB       OF PIVOTAL ROW
      JSB LWHR
      STA T1        ADDRESS OF PIVOTAL ROW
      LDA T5
      CLB,INB
      JSB LWHR
      STA T2        ADDR OF ROW TO BE SWAPPED
      LDA LPIV
      CLB,INB       COMPUTE ADDRESS OF
      JSB LWHR2     PIVOTAL ROW IN I-MATRIX
      STA T9
      STA T10       KEEP COPY
      LDA T5
      CLB,INB       COMPUTE ADDR OF ROW TO
      JSB LWHR2     BE SWAPPED IN I-MATRIX
      STA T11
      LDA T3
      CMA,INA
      STA T12       COUNTER FOR TERMS IN A ROW
LINV3 JSB .DLD      SWAP ONE ELEMENT OF ROW
      DEF T1,I
      STA T18
      STB T19
      JSB .DLD
      DEF T2,I
      STA T1,I
      ISZ T1
      STB T1,I
      ISZ T1
      LDA T18
      LDB T19
      STA T2,I
      ISZ T2
      STB T2,I
      ISZ T2
      JSB .DLD      SWAP ONE ELEMENT IN A ROW
      DEF T9,I      OF I-MATRIX
      STA T18
      STB T19
      JSB .DLD
      DEF T11,I
      STA T9,I
      ISZ T9
      STB T9,I
      ISZ T9
      LDA T18
      LDB T19
      STA T11,I
      ISZ T11
      STB T11,I
      ISZ T11
      ISZ T12       SKIP IF DONE
      JMP LINV3     SWAP NEXT ELEMENT
*                                   HAVE LARGEST ELEMENT IN
*                                   PIVOTAL POSITION. FIND
*                                   VALUE AND TEST TO ZERO
*                                   FOR SINGULAR MATRIX
      LDA LPIV      COMPUTE
      LDB LPIV        ADDRESS OF
      JSB LWHR          PIVOT
      STA T1             ELEMENT
      JSB .DLD        PIVOT VALUE
      DEF T1,I
      SSA           OBTAIN ABSOLUTE VALUE
      JSB ARINV      IF NUMBER IS NEGATIVE
      JSB .FSBA,I   SUBTRACT TOLERANCE AND
      DEF MLBX1
      SSA           COMPARE TO ZERO
      JSB ERROR     PRINT'NEARLY SING MATRIX'
*                                   DIVIDE PIVOT ROW AND ROW
*                                   IN I-MAT BY PIVOT VALUE
LDUM1 LDA T1        ADDRESS OF PIOT ELEMENT
      STA T2
      LDA HONE      LOAD
      LDB .2          1.0
      JSB .FDVA,I
      DEF T1,I
      STA T18       INVERSE OF PIVOT
      STB T19
*                                   MULT ROW BY 1/PIVOT
*                                   STARTING AT PIVOT+1
      LDA LPIV
      STA T11       COUNTER FOR ROW
LINV6 ISZ T11       INCREMENT COUNTER
      LDA T11
      CPA T4        TEST FOR END OF ROW
      JMP LIN12
      ISZ T2        ADDRESS OF NEXT ELEMENT
      ISZ T2
      JSB .DLD
      DEF T2,I
      JSB .FMPA,I
      DEF T18
      JSB .DST
      DEF T2,I
      JMP LINV6
*                                   MULT ROW IN I-MATRIX BY
*                                   1/PIVOT. SKIP IF ELEMENT=0
LIN12 LDA T10       ADDRESS OF PIVOT ROW
      STA T5        IN I-MATRIX
      LDA T3
      CMA,INA  SET
      STA T11       ROW COUNTER
LIN13 JSB .DLD
      DEF T5,I
      SZA,RSS       SKIP MULTIPLICATION IF ZERO
      SZB
      JMP *+2       NOT ZERO
      JMP LIN14     ZERO
      JSB .FMPA,I
      DEF T18
      JSB .DST
      DEF T5,I
LIN14 ISZ T5        NEXT ELEMENT IN I-MATRIX
      ISZ T5
      ISZ T11  DONE?
      JMP LIN13     NO
*                                   PERFORM ROW MANIPULATIONS
*                                   AND SUBTRACTIONS TO REDUCE
*                                   PIVOT COLUMN TO ZERO
      CLA
      STA B1
LINV4 ISZ B1   SELECT NEXT ROW
      LDA B1
      CPA T4        TEST FOR LAST ROW
      JMP LINV1     SELECT NEXT PIVOT
      CPA LPIV      TEST TO SKIP PIVOTAL ROW
      JMP LINV4     SKIP PIVOTAL ROW
      LDA B1
      CLB,INB
      JSB LWHR2     ADDRESS OF ROW TO BE TRANSFORMED
      STA T11           IN I-MATRIX
*                                   COMPUTE MULTIPLIER WHICH
*                                   IS THAT ELEMENT IN ROW TO
*                                   BE TRANSFORMED WHICH LIES
*                                   IN THE PIVOTAL COLUMN
      LDA B1
      LDB LPIV
      JSB LWHR
      STA T9        SAVE ADDRESS
      JSB .DLD
      DEF 0,I
      STA T7        VALUE OF MULTIPLIER
      STB T8
*                                   DO ELIMINATION OF ROWS IN
*                                   ORIGINAL MATRIX. START AT
*                                   COLUMN LPIV+1
      LDA LPIV
      STA T13       COUNTER
      LDA T1
      STA T2
LINV5 ISZ T13
      LDA T13
      CPA T4        TEST FOR LAST TERM IN ROW
      JMP LIN15
      ISZ T9        T9 IS ADDRESS OF
      ISZ T9        ELEMENT TO BE CHANGED
      ISZ T2        T2 IS ADDR OF CORRESPONDING
      ISZ T2        ELEMENT IN PIVOTAL ROW
      LDA T7
      LDB T8
      JSB .FMPA,I
      DEF T2,I
      STA T18       MULTIPLIER*VALUE IN
      STB T19           PIVOT ROW
      JSB .DLD
      DEF T9,I
      JSB .FSBA,I
      DEF T18
      JSB .DST      TRANSFORMED ELEMENT
      DEF T9,I
      JMP LINV5     SELECT NEXT TERM
*                                   DO ELIMINATION OF ROWS IN
*                                   IDENTITY MATRIX. START AT
*                                   BEGINNING OF ROW AND LEAVE
*                                   ELEMENT UNCHANGED WHEN ZERO
*                                   IN PIVOTAL ROW.
LIN15 LDA T10       ADDRESS OF
      STA T5           PIVOTAL ROW
      LDA T3
      CMA,INA       SET
      STA T13         COUNTER
LIN18 LDA T5,I
      ISZ T5
      LDB T5,I
      ISZ T5
      SZA,RSS       SKIP IF ZERO
      SZB
      JMP *+2       NOT ZERO
      JMP LIN17     ZERO
      JSB .FMPA,I   MULTIPLY BY
      DEF T7             MULTIPLIER
      STA T18
      STB T19
      JSB .DLD
      DEF T11,I
      JSB .FSBA,I
      DEF T18
      JSB .DST
      DEF T11,I
LIN17 ISZ T11
      ISZ T11
      ISZ T13
      JMP LIN18     SELECT NEXT TERM
      JMP LINV4     ELIMINATE NEXT ROW
*
*
*********************************************
*****    SUBROUTINE LWHR                *****
*********************************************
*SUBROUTINE COMPUTES ADDRESS OF AN ELEMENT  *
*IN MATRIX GIVEN BY B2. ROW AND COL VALUES  *
*ARE SUPPLIED IN A,B. ADDRESS IS LEFT IN A  *
*ENTRY LWHR2 COMPUTES ADDR IN MAT B3        *
*********************************************
*
LWHR  NOP
      STB T7        SAVE COLUMN #
      ADA M1
      JSB MPY
      DEF T3        (A-1)*T3
      ADA T7
      ADA M1        +(B-1)
      ALS
      ADA B2        DDR=B2+2((A-1)*T3+(B-1))
      JMP LWHR,I
LWHR2 NOP
      STB T7
      ADA M1
      JSB MPY
      DEF T3
      ADA T7
      ADA M1
      ALS
      ADA B3
      JMP LWHR2,I
*
*
*********************************************
*              CONSTANTS                    *
*********************************************
*
T1    BSS 1         TEMPORARY CONSTANTS
T2    BSS 1
T3    BSS 1
T4    BSS 1
T5    BSS 1
T6    BSS 1
T7    BSS 1
T8    BSS 1
T9    BSS 1
T10   BSS 1
T11   BSS 1
T12   BSS 1
T13   BSS 1
T16   DEC +1E-6     ABSOLUTE TOLERANCE
T18   BSS 1
T19   BSS 1
LPIV  BSS 1
LPLUS JSB .FADA,I   GENERATES CODE
      DEF B2,I
LMIN  JSB .FSBA,I   GENERATES CODE
LTIME JSB .FMPA,I   GENERATES CODE
INCB2 ISZ B2        GENERATES CODE
FINIS EQU *
      END
