%include "itrimp_hostcodes"
%CONSTINTEGER TARGET=M88K
!%EXTERNALROUTINESPEC IOCP %ALIAS "s_iocp"(%INTEGER EP,N)
%EXTERNALINTEGERFN BITS(%INTEGER X)
%INTEGER I,J
      J=0
      J=J+X&1 %AND X=X>>1 %WHILE X#0
      %RESULT=J
%END
%EXTERNALINTEGERFN PARITY(%INTEGER I)
      %RESULT=1-(I&1)*2
%END
%EXTERNALROUTINE SOLVE LN EQ(%LONGREALARRAYNAME A,B, %INTEGER N,
    %LONGREALNAME DET)
%LONGREAL AMAX,CH
%INTEGER I,J,J MAX,S
      ->L3 %IF N>0
      PRINTSTRING("
      SOLVE LN EQ DATA FAULT: N=")
      WRITE(N,2); NEWLINE; %STOP
L3:   ->L1 %IF N>1
      DET=A(1,1)
      ->L2 %IF DET=0
      B(1)=B(1)/DET
      ->L2
L1:   DET=1
      %CYCLE I=1,1,N-1
         A MAX=0; J MAX=0
         %CYCLE J=I,1,N
            ->L4 %IF MOD(A(J,I))<=AMAX
            AMAX=MOD(A(J,I)); JMAX=J
L4:      %REPEAT
         ->L5 %IF J MAX=I
         DET=-DET
         ->L6 %IF J MAX#0
         DET=0; ->L2
L6:      %CYCLE J=I,1,N
            CH=A(I,J)
            A(I,J)=A(J MAX,J)
            A(J MAX,J)=CH
         %REPEAT
         CH=B(I)
         B(I)=B(J MAX)
         B(J MAX)=CH
L5:      CH=A(I,I)
         DET=DET*CH
         %CYCLE J=I+1,1,N
            A MAX=A(J,I)/CH
            %CYCLE S=I+1,1,N
               A(J,S)=A(J,S)-A(I,S)*A MAX
            %REPEAT
            B(J)=B(J)-B(I)*A MAX
         %REPEAT
      %REPEAT
      CH=A(N,N)
      DET=DET*CH
      ->L2 %IF DET=0
      B(N)=B(N)/CH
      %CYCLE I=N-1,-1,1
         CH=B(I)
         %CYCLE J=I+1,1,N
            CH=CH-A(I,J)*B(J)
         %REPEAT
         B(I)=CH/A(I,I)
      %REPEAT
L2:
%END
%EXTERNALROUTINE DIV MATRIX(%LONGREALARRAYNAME A,B, %INTEGER N,M,
    %LONGREALNAME DET)
                                        ! A=INV(B)A:BNXN, ANXM
%LONGREAL AMAX,CH
%INTEGER I,J,JMAX,S,K
      ->L3 %IF N>0
      PRINTSTRING("
      DIV MATRIX DATA FAULT N=")
      WRITE(N,2)
      NEWLINE; %STOP
L3:   ->L1 %IF N>1
      DET=B(1,1)
      ->L2 %IF DET=0
      %CYCLE I=1,1,M
         A(1,I)=A(1,I)/DET
      %REPEAT
      ->L2
L1:   DET=1
      %CYCLE I=1,1,N-1
         AMAX=0; JMAX=0
         %CYCLE J=I,1,N
            ->L4 %IF MOD(B(J,I))<=AMAX
            AMAX=MOD(B(J,I)); JMAX=J
L4:      %REPEAT
         ->L5 %IF J MAX=I
         DET=-DET
         ->L6 %IF JMAX#0
         DET=0; ->L2
L6:      %CYCLE J=I,1,N
            CH=B(I,J)
            B(I,J)=B(JMAX,J)
            B(JMAX,J)=CH
         %REPEAT
         %CYCLE K=1,1,M
            CH=A(I,K)
            A(I,K)=A(JMAX,K)
            A(JMAX,K)=CH
         %REPEAT
L5:      CH=B(I,I)
         DET=DET*CH
         %CYCLE J=I+1,1,N
            AMAX=B(J,I)/CH
            %CYCLE S=I+1,1,N
               B(J,S)=B(J,S)-B(I,S)*AMAX
            %REPEAT
            %CYCLE K=1,1,M
               A(J,K)=A(J,K)-A(I,K)*AMAX
            %REPEAT
         %REPEAT
      %REPEAT
      CH=B(N,N)
      DET=DET*CH
      ->L2 %IF DET=0
      %CYCLE K=1,1,M
         A(N,K)=A(N,K)/CH
      %REPEAT
      %CYCLE I=N-1,-1,1
         AMAX=B(I,I)
         %CYCLE K=1,1,M
            CH=A(I,K)
            %CYCLE J=I+1,1,N
               CH=CH-B(I,J)*A(J,K)
            %REPEAT
            A(I,K)=CH/AMAX
         %REPEAT
      %REPEAT
L2:
%END
%EXTERNALROUTINE UNIT(%LONGREALARRAYNAME A, %INTEGER N)
%INTEGER I,J
      ->L10 %IF N>0
      PRINTSTRING("
MATRIX BOUND ZERO OR NEGATIVE")
      NEWLINE
      %MONITOR
      %STOP
L10:
      %CYCLE I=1,1,N
         %CYCLE J=1,1,N
            A(I,J)=0
         %REPEAT
         A(I,I)=1
      %REPEAT
%END
%EXTERNALROUTINE INVERT(%LONGREALARRAYNAME A,B, %INTEGER N, %LONGREALNAME DET)
                                        ! A=INV B USING DIV MATRIX
      ->L3 %IF N>0
      PRINTSTRING("
      INVERT DATA FAULT N=")
      WRITE(N,2); NEWLINE; %STOP
L3:   UNIT(A,N)
      DIV MATRIX(A,B,N,N,DET)
%END
%EXTERNALLONGREALFN DET(%LONGREALARRAYNAME A, %INTEGER N)
%LONGREALARRAY B(1:N); %LONGREAL DETVAL
%INTEGER I
      ->L10 %IF N>0
      PRINTSTRING("
MATRIX BOUND ZERO OR NEGATIVE")
      NEWLINE
      %MONITOR
      %STOP
L10:
      %CYCLE I=1,1,N
         B(I)=0
      %REPEAT
      SOLVE LN EQ(A,B,N,DETVAL)
      %RESULT=DETVAL
%END
%EXTERNALROUTINE NULL(%LONGREALARRAYNAME A, %INTEGER N,M)
%INTEGER I,J
      ->L10 %IF N>0 %AND M>0
      PRINTSTRING("
MATRIX BOUND ZERO OR NEGATIVE")
      NEWLINE
      %MONITOR
      %STOP
L10:
      %CYCLE I=1,1,N
         %CYCLE J=1,1,M
            A(I,J)=0
         %REPEAT
      %REPEAT
%END
%EXTERNALROUTINE ADD MATRIX(%LONGREALARRAYNAME A,B,C, %INTEGER N,M)
%INTEGER I,J
      ->L10 %IF N>0 %AND M>0
      PRINTSTRING("
MATRIX BOUND ZERO OR NEGATIVE")
      NEWLINE
      %MONITOR
      %STOP
L10:
      %CYCLE I=1,1,N
         %CYCLE J=1,1,M
            A(I,J)=B(I,J)+C(I,J)
         %REPEAT
      %REPEAT
%END
%EXTERNALROUTINE SUB MATRIX(%LONGREALARRAYNAME A,B,C, %INTEGER N,M)
%INTEGER I,J
      ->L10 %IF N>0 %AND M>0
      PRINTSTRING("
MATRIX BOUND ZERO OR NEGATIVE")
      NEWLINE
      %MONITOR
      %STOP
L10:
      %CYCLE I=1,1,N
         %CYCLE J=1,1,M
            A(I,J)=B(I,J)-C(I,J)
         %REPEAT
      %REPEAT
%END
%EXTERNALROUTINE COPY MATRIX(%LONGREALARRAYNAME A,B, %INTEGER N,M)
%INTEGER I,J
      ->L10 %IF N>0 %AND M>0
      PRINTSTRING("
MATRIX BOUND ZERO OR NEGATIVE")
      NEWLINE
      %MONITOR
      %STOP
L10:
      %CYCLE I=1,1,N
         %CYCLE J=1,1,M
            A(I,J)=B(I,J)
         %REPEAT
      %REPEAT
%END
%EXTERNALROUTINE MULT MATRIX(%LONGREALARRAYNAME A,B,C, %INTEGER N,P,M)
                                        ! A=B*C, A IS N X M
%INTEGER I,J,K
%LONGREAL R
      ->L10 %IF N>0 %AND M>0 %AND P>0
      PRINTSTRING("
MATRIX BOUND ZERO OR NEGATIVE")
      NEWLINE
      %MONITOR
      %STOP
L10:
      %CYCLE I=1,1,N
         %CYCLE J=1,1,M
            R=0
            %CYCLE K=1,1,P
               R=R+B(I,K)*C(K,J)
            %REPEAT
            A(I,J)=R
         %REPEAT
      %REPEAT
%END
%EXTERNALROUTINE MULT TR MATRIX(%LONGREALARRAYNAME A,B,C, %INTEGER N,P,M)
%LONGREAL R
                                        ! A=B*C, A IS N X M
%INTEGER I,J,K
      ->L10 %IF N>0 %AND M>0 %AND P>0
      PRINTSTRING("
MATRIX BOUND ZERO OR NEGATIVE")
      NEWLINE
      %MONITOR
      %STOP
L10:
      %CYCLE I=1,1,N
         %CYCLE J=1,1,M
            R=0
            %CYCLE K=1,1,P
               R=R+B(I,K)*C(J,K)
            %REPEAT
            A(I,J)=R
         %REPEAT
      %REPEAT
%END
%EXTERNALROUTINE TRANS MATRIX(%LONGREALARRAYNAME A,B, %INTEGER N,M)
                                        ! AN X M, B M X N
%INTEGER I,J
      ->L10 %IF N>0 %AND M>0
      PRINTSTRING("
MATRIX BOUND ZERO OR NEGATIVE")
      NEWLINE
      %MONITOR
      %STOP
L10:
      %CYCLE I=1,1,N
         %CYCLE J=1,1,M
            A(I,J)=B(J,I)
         %REPEAT
      %REPEAT
%END
%EXTERNALREALFN RANDOM(%INTEGERNAME I, %INTEGER N)
%REAL X
%INTEGER J
      ->L2 %IF N<0; ->L1 %IF N>1
      %IF target=emas %THENSTART
         J=I
         *LSS_65539
         *IMYD_J
         *AND_X'000000007FFFFFFF'
         *STUH_ %B
         *ST_J
         I=J
      %FINISHELSE I=(65539*I)&X'7FFFFFFF'
      %RESULT=0.0000000004656613*I
L1:   X=0
      %CYCLE N=1,1,N
         X=X+RANDOM(I,1)
      %REPEAT
      %RESULT=X
L2:   PRINTSTRING("NEGATIVE ARGUMENT IN RANDOM")
      %MONITOR
      %STOP
%END
%EXTERNALLONGREALFN LRANDOM( %INTEGERNAME I, %INTEGER N)
%LONGREAL X
%constinteger magic=65539
%INTEGER J
      ->L2 %IF N < 0; ->L1 %IF N > 1
      %if target= emas %Start
      J = I
      *LSS_magic
      *IMYD_J
      *AND_X'000000007FFFFFFF'
      *STUH_ %B
       *ST_J
      I = J
      %finish %else %if target=ibm %or target=ibmxa %or target=amdahl %Start
      j=i
        *l_1,magic
        *m_0,j
         *sll_1,1
      *srl_1,1
        *st_1,j
      i=j
     %finish %else I=(magic*I)&X'7FFFFFFF'
      %RESULT = 0.0000000004656613 * I
L1:   X = 0
      %CYCLE N=1,1,N; X = X + LRANDOM(I,1); %REPEAT; %RESULT = X
L2:   PRINTSTRING("NEGATIVE ARGUMENT IN LRANDOM"); %MONITOR
      %STOP
%END;                                    !OF LRANDOM
%ENDOFFILE