code 34143; comment MCA 2305; procedure TFMSYMTRI1(A, N, D, B, BB, EM); value N;integer N; array A, B, BB, D, EM; begin integer I, J, R, R1, P, Q, TI, TJ; real S, W, X, A1, B0, BB0, D0, NORM, MACHTOL; NORM:= 0; TJ:= 0; for J:= 1 step 1 until N do begin W:= 0; for I:= 1 step 1 until J do W:= ABS(A[I + TJ]) +W; TJ:= TJ + J; TI:= TJ + J; for I:= J + 1 step 1 until N do begin W:= ABS(A[TI]) + W; TI:= TI + I end; if W > NORM then NORM:= W end; MACHTOL:= EM[0] * NORM; EM[1]:= NORM; Q:= (N + 1) * N // 2; R:= N; for R1:= N - 1 step -1 until 1 do begin P:= Q - R; D[R]:= A[Q]; X:= VECVEC(P + 1, Q - 2, 0, A, A); A1:= A[Q - 1]; if SQRT(X) <= MACHTOL then begin B0:= B[R1]:= A1; BB[R1]:= B0 * B0; A[Q]:= 1 end else begin BB0:= BB[R1]:= A1 * A1 + X; B0:= if A1 > 0 then -SQRT(BB0) else SQRT(BB0); A1:= A[Q - 1]:= A1 - B0; W:= A[Q]:= 1 / (A1 * B0); TJ:= 0; for J:= 1 step 1 until R1 do begin TI:= TJ + J; S:= VECVEC(TJ + 1, TI, P - TJ, A, A); TJ:= TI + J; B[J]:= (SEQVEC(J + 1, R1, TJ, P, A, A) + S) * W; TJ:= TI end; ELMVEC(1, R1, P, B, A, VECVEC(1,R1,P,B,A)* W *.5); TJ:= 0; for J:= 1 step 1 until R1 do begin TI:= TJ + J; ELMVEC(TJ + 1, TI, P - TJ, A, A, B[J]);ELMVEC(TJ + 1, TI, -TJ, A, B, A[J + P]); TJ:= TI end; B[R1]:= B0 end; Q:= P; R:= R1 end; D[1]:= A[1]; A[1]:= 1; B[N]:= BB[N]:= 0 end TFMSYMTRI1 eop