MODULE R3;

    (* 
    
         Linear algebra in R^3
         
         Author: Jorge Stolfi
         
         Last edited: M. Carrard - 18/04/1993 - Transcricao para modula-3
         
     *)

IMPORT Math;

(* Return z[i] := x[i] + y[i] *)
PROCEDURE add( x, y: T ): T  =
VAR z: T;
BEGIN
  FOR i := 0 TO 2 DO  z[i] := x[i] + y[i]; END;
  RETURN z;
END add;


(* Return z[i] := x[i] + y[i] *)
PROCEDURE sub( x, y: T ): T  =
VAR z: T;
BEGIN
  FOR i := 0 TO 2 DO  z[i] := x[i] - y[i]; END;
  RETURN z;
END sub;


(* Return z[i] := x[i] * a *)
PROCEDURE scale( x: T; a: REAL ): T = 
VAR z: T;
BEGIN
  FOR i := 0 TO 2 DO z[i] := x[i] * a;  END;
  RETURN z;
END scale;


(* Scales x to unit L_infinity norm *)  
PROCEDURE reduce_inf( x: T ): T =
VAR m: REAL;
    z: T;
BEGIN
  m := 0.0;
  FOR i := 0 TO 2 DO 
     IF ABS( x[i] ) > m THEN m := ABS( x[i] ); END;
  END;
  FOR i := 0 TO 2 DO z[i] := x[i] / m;  END;
  RETURN z;
END reduce_inf;  


(* Scale x to unit Euclidean norm *)
PROCEDURE reduce( x: T ): T =
VAR z: T;
    m: REAL;
BEGIN
  m := 0.0;
  FOR i:= 0 TO 2 DO m := m + x[i] * x[i]; END;
  m := FLOAT( Math.sqrt( FLOAT( m, LONGREAL) ) );
  FOR i:= 0 TO 2 DO z[i] := x[i] / m; END;
  RETURN z;
END reduce; 


(* Dot product: sum of x[i] * y[i] *)
PROCEDURE dot( x, y: T ): REAL =
VAR m: REAL;
BEGIN
  m := x[0]*y[0] + x[1]*y[1] + x[2]*y[2];
  RETURN m;
END dot;


(* Cross-product: p x q *)
PROCEDURE cross( p, q: T ): T =
VAR d01, d02, d12: REAL;
    z: T;
BEGIN
  d01 := p[0]*q[1] - p[1]*q[0];
  d02 := p[0]*q[2] - p[2]*q[0];
  d12 := p[1]*q[2] - p[2]*q[1];

  z[0] := + d12;
  z[1] := - d02;
  z[2] := + d01;  
  RETURN z;
END cross;


(* Return the component of x that is orthogonal tu u *)
PROCEDURE orthize( x, u: T ): T =
VAR z: T;
    xi, ui, sxu, suu, c: REAL;
BEGIN
   suu := 0.0;
   sxu := 0.0;
   FOR i := 0 TO 2 DO 
      xi := x[i]; ui := u[i];
      suu := suu + ui*ui;
      sxu := sxu + xi*ui; 
   END;
   c := sxu/suu;
   FOR i := 0 TO 2 DO 
      xi := x[i]; ui := u[i];
      z[i] := xi - c*ui;
   END;
   RETURN z;
END orthize;

BEGIN 
END R3.