! $Id: main.f90 189 2011-02-08 19:20:24Z pamestoy $
! Solution of a sparse linear system with Preconditioned Conjugate Gradient.
! The preconditioner is an incomplete ("fill level") Cholesky factorization.

program ic_test
  use utils
  use icfacto
  use pcg
  implicit none

  character*72 :: mat_file
  integer :: want_amd, level, trace
  integer :: n, anz, lnz
  integer :: i, ierror, allocok, itmax, prec
  real(kind(1.d0)) :: eps
  logical :: doprec

  type(csc_mat) :: A, &   ! the input matrix A
       &           L      ! the incomplete factor of A

  real(kind(1.d0)), allocatable :: b(:), &  ! RHS
       &                           x(:)     ! solution

  ! Read parameters file
  call read_parameters(mat_file, want_amd, itmax, eps, prec, level, trace)

  ! Read matrix with general matrix-market format
  call csc_read_matrix(mat_file, A, Ierror, Trace)
  if(ierror.gt.0) then
     write(*,*)'Problem when reading matrix', Ierror
     stop
  endif

  ! if wand_amd is true permute matrix with approximate 
  ! minimum degree ordering
  call csc_ordering(want_amd, A, trace)

  ! compute incomplete factorization (forms the preconditioner)
  call csc_incomplete_facto(A, L, level, trace)

  ! build the right-hand side
  allocate(b(A%n), stat=allocok)
  allocate(x(A%n), stat=allocok)
  if (allocok.ne.0) then
    write(*,*)'Problem with allocation of rhs.', allocok
    stop
  end if
  b = 1.d0
  x = 0.d0

  ! use preconditioned conjugate gradient to solves the system
  ! and print statistics
  doprec = prec .eq. 1
  call conjugate_gradient(A, b, x, doprec, L, itmax, eps, trace)

  ! clean everything before exit
  deallocate(b, x, A%jptr, A%irn, A%val, L%jptr, L%irn, L%val)

  stop

end program ic_test
