out aa7;
off nat;
%on rounded;

procedure delta(ii,jj);
begin 
if ii eq jj then return 1 else return 0; 
end;

procedure addone(x); rrr:=x+1;

% random numbers with uniform distribution generated in Mathematica 
%no=Table[Random[Real],{i,1,100 }]

noise:=
{0.7611189572266488, 0.809932949360867, 
  0.32754272053822, 0.3711809230728502, 
  0.2198155852281288, 0.2810910296348176, 
  0.5994259353572776, 0.987082750593675, 
  0.6503306993502929, 0.6922800675720912, 
  0.7549523349478609, 0.445705988068098}$

n:=2;
m:=2;
%-----------------
AA:=TP(MAT( (A(1), A(2) )));
%-----------------
matrix KM(n,n);
for tt:=1:n do for pp:=1:m do km(tt,pp):=k(tt,pp);

%-----------------
RR:=TP(MAT( (R(1), R(2) )));
RR - KM*AA;

NN:=TP(MAT( (NO(1), NO(2) )));

RMN:=TP(MAT( (d(1), d(2) )));

%-----------------
matrix SM(n,n);
for tt:=1:n do for pp:=1:m do sm(tt,pp):=s(tt,pp);


SM*RMN - RMN;


S(1,2):=-(D(1)*S(1,1) - D(1))/D(2);
S(2,1):=-(D(2)*S(2,2) - D(2) )/D(1);

SM;

SM*RMN - RMN;


SM*(KM*AA + NN) - (KM*AA + NN); 
%---------------------------------------
k(1,1):=4;k(1,2):=3;
k(2,1):=3;k(2,2):=4;
%a(1):=1; a(2):=2;
%no(1):=8;no(2):=5;

d(1):=18; d(2):=16;
km;
aa;
nn;


KM*AA;
KM*AA + NN;

%--------
%s(1,1):=2;s(2,2):=3;
%NO(1):=(43*A(2) - 113*A(1) + 16*ARBCOMPLEX(1))/23;
%NO(2):=ARBCOMPLEX(1);
%SM;

%--------
%s(1,1):=10;s(2,2):=3;
%NO(1):=(43*A(2) - 113*A(1) + 16*ARBCOMPLEX(1))/23;
%NO(2):=ARBCOMPLEX(1);
SM;


%-----------------------------
factor s;
SM;
SM*RMN - RMN;
eqq:=SM*(KM*AA + NN) - (KM*AA + NN); 
xx:=solve({eqq(1,1), eqq(2,1)},{NO(1), NO(2)});
eqq(1,1);
eqq(2,1);
%---------------------------




%S(1,2):=-(D(1)*S(1,1) - D(1))/D(2);
%S(2,1):=-(D(2)*S(2,2) - D(2) )/D(1);

%eqq:=SM*(R+Qq) - R;
qq:=TP(MAT( (q(1), q(2) )));

matrix SXM(n,n);
for tt:=1:n do for pp:=1:m do sxm(tt,pp):=sx(tt,pp);


sx(1,2):=( - D(1)*sx(1,1) + D(1) - Q(1)*sx(1,1))/(D(2) + Q(2));
sx(2,1):=( - D(2)*sx(2,2) + D(2) - Q(2)*sx(2,2))/(D(1) + Q(1));

sxm;

eqqx:=sxM*(KM*AA + NN+Qq) - (KM*AA + NN);

xx1:=solve({eqqx(1,1), eqqx(2,1)},{NO(1), NO(2)});



%a(1):=50; a(2):=100;
%xx;

%no(2):=300;
%xx;

%ab:=300;
%NO(2):=AB;
%NO(1):=(16*AB - 1350)/23;

d(1);
d(2);

shut aa7;
bye;
end;