femwind/fortran/lin_alg_test.f90

   1 program lin_alg_test
   2 use module_lin_alg
   3 implicit none
   4 integer, parameter :: n=3
   5 real, dimension(n,n) :: A, Ainv, Id, res
   6 integer i,j
   7
   8 A  = reshape([1., 2., 3., 4., 5., 6., 7., 8., 10.],[3, 3])
   9 Id = reshape([1., 0., 0., 0., 1., 0., 0., 0., 1.],[3, 3])
  10
  11 !do i=1,n
  12 !    do j=1,n
  13 !        call random_number(A(i,j))
  14 !        if(i .eq. j)then
  15 !            id(i,j)=1.
  16 !        else
  17 !            id(i,j)=0.
  18 !        endif
  19 !    enddo
  20 !enddo
  21
  22 call print_matrix('A',A)
  23 call inv3( A, Ainv)
  24 call print_matrix('Ainv',Ainv)
  25 res = matmul(A, Ainv)
  26 call print_matrix('A*Ainv',res)
  27 res = res - Id
  28 call print_matrix('A*Ainv - I',res)
  29
  30 print *,'lin_alg_test: error=',maxval(abs(res))
  31 end program lin_alg_test