/*  連立方程式についてのプログラムを書いてテストにする. 
*/
def gaussE(A,Y) {
  N = size(A)[0];
  for (K=1; K<=N-1; K++) {
    for (I=K+1; I<=N; I++) {
      P = A[I-1][K-1]/A[K-1][K-1];
      for (J=K; J<=N; J++) {
        A[I-1][J-1] = A[I-1][J-1]-P*A[K-1][J-1];
      }
      Y[I-1] = Y[I-1]-P*Y[K-1];
    } 
    print(A); print(Y); print("----------------");
  }
}

def gtest1() {
  A=newmat(3,3,[[1,-1,2],[-1,2,-3],[3,1,1]]);
  Y=newvect(3,[5,-6,8]);
  print(A); print(Y); print("----------------");
  gaussE(A,Y);
  return [A,Y];
}

def laplacian(N) {
   A = newmat(N,N);
   for (I=0; I<N; I++) {
	 if (I > 0) A[I][I-1] = 1;
     	A[I][I] = -2;
          if  (l < N-1) A[I][I+1] = 1;
   }
   return A;
}
[gtest1(), laplacian(3)];