fdurastante/psblas3
1
0
Fork 0
mirror of https://github.com/sfilippone/psblas3.git synced 2026-10-16 02:49:03 +00:00
Code Issues Packages Projects Releases Wiki Activity

psblas3:

Browse Source 
base/Makefile
 base/newserial/psbn_base_mat_mod.f03
 base/newserial/psbn_d_base_mat_mod.f03
 base/newserial/psbn_d_coo_impl.f03
 base/newserial/psbn_d_csr_impl.f03
 base/newserial/psbn_d_csr_mat_mod.f03
 base/newserial/psbn_mat_impl.f03
 base/newserial/psbn_mat_mod.f03
 test/serial/Makefile
 test/serial/d_matgen.f03
 test/serial/psbn_d_cxx_impl.f03
 test/serial/psbn_d_cxx_mat_mod.f03

Improved error messages. 
Various cosmetic fixes.
This commit is contained in:
Salvatore Filippone
2009-09-10 16:43:24 +00:00
parent 7ee565006c
commit 607a0aa949
12 changed files with 1874 additions and 1648 deletions
Download Patch File Download Diff File Expand all files Collapse all files
base/Makefile
+1
View File
@@ -22,6 +22,7 @@ clean:
(cd tools; make clean)
(cd serial; make clean)
(cd psblas; make clean)
(cd newserial; make clean)
veryclean: clean
/bin/rm -f $(HERE)/$(LIBNAME) $(LIBMOD) *$(.mod)
base/newserial/psbn_base_mat_mod.f03
+34 -18
View File
@@ -30,6 +30,35 @@ module psbn_base_mat_mod
integer, private :: state, duplicate
logical, private :: triangle, unitd, upper, sorted
contains
! ====================================
!
! Getters
!
!
! ====================================
procedure, pass(a) :: get_nrows
procedure, pass(a) :: get_ncols
procedure, pass(a) :: get_nzeros
procedure, pass(a) :: get_size
procedure, pass(a) :: get_state
procedure, pass(a) :: get_dupl
procedure, pass(a) :: get_fmt
procedure, pass(a) :: is_null
procedure, pass(a) :: is_bld
procedure, pass(a) :: is_upd
procedure, pass(a) :: is_asb
procedure, pass(a) :: is_sorted
procedure, pass(a) :: is_upper
procedure, pass(a) :: is_lower
procedure, pass(a) :: is_triangle
procedure, pass(a) :: is_unit
! ====================================
!
! Setters
!
! ====================================
procedure, pass(a) :: set_nrows
procedure, pass(a) :: set_ncols
procedure, pass(a) :: set_dupl
@@ -45,31 +74,18 @@ module psbn_base_mat_mod
procedure, pass(a) :: set_unit
procedure, pass(a) :: get_nrows
procedure, pass(a) :: get_ncols
procedure, pass(a) :: get_nzeros
procedure, pass(a) :: get_size
procedure, pass(a) :: get_state
procedure, pass(a) :: get_dupl
procedure, pass(a) :: get_fmt
procedure, pass(a) :: is_null
procedure, pass(a) :: is_bld
procedure, pass(a) :: is_upd
procedure, pass(a) :: is_asb
procedure, pass(a) :: is_sorted
procedure, pass(a) :: is_upper
procedure, pass(a) :: is_lower
procedure, pass(a) :: is_triangle
procedure, pass(a) :: is_unit
! ====================================
!
! Data management
!
! ====================================
procedure, pass(a) :: get_neigh
procedure, pass(a) :: allocate_mnnz
procedure, pass(a) :: reallocate_nz
procedure, pass(a) :: free
generic, public :: allocate => allocate_mnnz
generic, public :: reallocate => reallocate_nz
procedure, pass(a) :: print => sparse_print
end type psbn_base_sparse_mat
base/newserial/psbn_d_base_mat_mod.f03
+817 -678
View File
File diff suppressed because it is too large Load Diff
base/newserial/psbn_d_coo_impl.f03
+210 -420
View File
@@ -2,7 +2,7 @@
subroutine d_coo_cssm_impl(alpha,a,x,beta,y,info,trans)
use psb_error_mod
use psbn_d_base_mat_mod, psb_protect_name => d_coo_cssm_impl
implicit none
implicit none
class(psbn_d_coo_sparse_mat), intent(in) :: a
real(psb_dpk_), intent(in) :: alpha, beta, x(:,:)
real(psb_dpk_), intent(inout) :: y(:,:)
@@ -852,6 +852,54 @@ subroutine d_coo_csmm_impl(alpha,a,x,beta,y,info,trans)
end subroutine d_coo_csmm_impl
function d_coo_csnmi_impl(a) result(res)
use psb_error_mod
use psbn_d_base_mat_mod, psb_protect_name => d_coo_csnmi_impl
implicit none
class(psbn_d_coo_sparse_mat), intent(in) :: a
real(psb_dpk_) :: res
integer :: i,j,k,m,n, nnz, ir, jc, nc
real(psb_dpk_) :: acc
real(psb_dpk_), allocatable :: tmp(:,:)
logical :: tra
Integer :: err_act
character(len=20) :: name='d_base_csnmi'
logical, parameter :: debug=.false.
res = dzero
nnz = a%get_nzeros()
i = 1
j = i
do while (i<=nnz)
do while ((a%ia(j) == a%ia(i)).and. (j <= nnz))
j = j+1
enddo
acc = dzero
do k=i, j-1
acc = acc + abs(a%val(k))
end do
res = max(res,acc)
end do
end function d_coo_csnmi_impl
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!
!
!
! Data management
!
!
!
!
!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
subroutine d_coo_csput_impl(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl)
use psb_error_mod
@@ -950,80 +998,6 @@ subroutine d_coo_csput_impl(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl)
contains
!!$
!!$ subroutine psb_inner_upd(nz,ia,ja,val,nza,aspk,maxsz,&
!!$ & imin,imax,jmin,jmax,info,gtl)
!!$ implicit none
!!$
!!$ integer, intent(in) :: nz, imin,imax,jmin,jmax,maxsz
!!$ integer, intent(in) :: ia(:),ja(:)
!!$ integer, intent(inout) :: nza
!!$ real(psb_dpk_), intent(in) :: val(:)
!!$ real(psb_dpk_), intent(inout) :: aspk(:)
!!$ integer, intent(out) :: info
!!$ integer, intent(in), optional :: gtl(:)
!!$ integer :: i,ir,ic, ng,nzl
!!$ character(len=20) :: name, ch_err
!!$
!!$
!!$ name='psb_inner_upd'
!!$ nzl = 0
!!$ if (present(gtl)) then
!!$ ng = size(gtl)
!!$ if ((nza > nzl)) then
!!$ do i=1, nz
!!$ nza = nza + 1
!!$ if (nza>maxsz) then
!!$ call psb_errpush(50,name,i_err=(/7,maxsz,5,0,nza /))
!!$ info = -71
!!$ return
!!$ endif
!!$ aspk(nza) = val(i)
!!$ end do
!!$ else
!!$ do i=1, nz
!!$ ir = ia(i)
!!$ ic = ja(i)
!!$ if ((ir >=1).and.(ir<=ng).and.(ic>=1).and.(ic<=ng)) then
!!$ ir = gtl(ir)
!!$ ic = gtl(ic)
!!$ if ((ir >=imin).and.(ir<=imax).and.(ic>=jmin).and.(ic<=jmax)) then
!!$ nza = nza + 1
!!$ if (nza>maxsz) then
!!$ info = -72
!!$ return
!!$ endif
!!$ aspk(nza) = val(i)
!!$ end if
!!$ end if
!!$ end do
!!$ end if
!!$ else
!!$ if ((nza >= nzl)) then
!!$ do i=1, nz
!!$ nza = nza + 1
!!$ if (nza>maxsz) then
!!$ info = -73
!!$ return
!!$ endif
!!$ aspk(nza) = val(i)
!!$ end do
!!$ else
!!$ do i=1, nz
!!$ ir = ia(i)
!!$ ic = ja(i)
!!$ if ((ir >=imin).and.(ir<=imax).and.(ic>=jmin).and.(ic<=jmax)) then
!!$ nza = nza + 1
!!$ if (nza>maxsz) then
!!$ info = -74
!!$ return
!!$ endif
!!$ aspk(nza) = val(i)
!!$ end if
!!$ end do
!!$ end if
!!$ end if
!!$ end subroutine psb_inner_upd
subroutine psb_inner_ins(nz,ia,ja,val,nza,ia1,ia2,aspk,maxsz,&
& imin,imax,jmin,jmax,info,gtl)
@@ -1082,151 +1056,56 @@ contains
end subroutine psb_inner_ins
subroutine d_coo_srch_upd(nz,ia,ja,val,a,&
& imin,imax,jmin,jmax,info,gtl)
use psb_const_mod
use psb_realloc_mod
use psb_string_mod
use psb_serial_mod
implicit none
class(psbn_d_coo_sparse_mat), intent(inout) :: a
integer, intent(in) :: nz, imin,imax,jmin,jmax
integer, intent(in) :: ia(:),ja(:)
real(psb_dpk_), intent(in) :: val(:)
integer, intent(out) :: info
integer, intent(in), optional :: gtl(:)
integer :: i,ir,ic, ilr, ilc, ip, &
& i1,i2,nc,nnz,dupl,ng
integer :: debug_level, debug_unit
character(len=20) :: name='d_coo_srch_upd'
info = 0
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
subroutine d_coo_srch_upd(nz,ia,ja,val,a,&
& imin,imax,jmin,jmax,info,gtl)
dupl = a%get_dupl()
if (.not.a%is_sorted()) then
info = -4
return
end if
ilr = -1
ilc = -1
nnz = a%get_nzeros()
if (present(gtl)) then
ng = size(gtl)
select case(dupl)
case(psbn_dupl_ovwrt_,psbn_dupl_err_)
! Overwrite.
! Cannot test for error, should have been caught earlier.
do i=1, nz
ir = ia(i)
ic = ja(i)
if ((ir >=1).and.(ir<=ng).and.(ic>=1).and.(ic<=ng)) then
ir = gtl(ir)
if ((ir > 0).and.(ir <= a%m)) then
ic = gtl(ic)
if (ir /= ilr) then
i1 = psb_ibsrch(ir,nnz,a%ia)
i2 = i1
do
if (i2+1 > nnz) exit
if (a%ia(i2+1) /= a%ia(i2)) exit
i2 = i2 + 1
end do
do
if (i1-1 < 1) exit
if (a%ia(i1-1) /= a%ia(i1)) exit
i1 = i1 - 1
end do
ilr = ir
else
i1 = 1
i2 = 1
end if
nc = i2-i1+1
ip = psb_issrch(ic,nc,a%ja(i1:i2))
if (ip>0) then
a%val(i1+ip-1) = val(i)
else
info = i
return
end if
else
if (debug_level >= psb_debug_serial_) &
& write(debug_unit,*) trim(name),&
& ': Discarding row that does not belong to us.'
endif
end if
end do
case(psbn_dupl_add_)
! Add
do i=1, nz
ir = ia(i)
ic = ja(i)
if ((ir >=1).and.(ir<=ng).and.(ic>=1).and.(ic<=ng)) then
ir = gtl(ir)
ic = gtl(ic)
if ((ir > 0).and.(ir <= a%m)) then
if (ir /= ilr) then
i1 = psb_ibsrch(ir,nnz,a%ia)
i2 = i1
do
if (i2+1 > nnz) exit
if (a%ia(i2+1) /= a%ia(i2)) exit
i2 = i2 + 1
end do
do
if (i1-1 < 1) exit
if (a%ia(i1-1) /= a%ia(i1)) exit
i1 = i1 - 1
end do
ilr = ir
else
i1 = 1
i2 = 1
end if
nc = i2-i1+1
ip = psb_issrch(ic,nc,a%ja(i1:i2))
if (ip>0) then
a%val(i1+ip-1) = a%val(i1+ip-1) + val(i)
else
info = i
return
end if
else
if (debug_level >= psb_debug_serial_) &
& write(debug_unit,*) trim(name),&
& ': Discarding row that does not belong to us.'
end if
end if
end do
case default
info = -3
if (debug_level >= psb_debug_serial_) &
& write(debug_unit,*) trim(name),&
& ': Duplicate handling: ',dupl
end select
else
select case(dupl)
case(psbn_dupl_ovwrt_,psbn_dupl_err_)
! Overwrite.
! Cannot test for error, should have been caught earlier.
do i=1, nz
ir = ia(i)
ic = ja(i)
use psb_const_mod
use psb_realloc_mod
use psb_string_mod
use psb_serial_mod
implicit none
class(psbn_d_coo_sparse_mat), intent(inout) :: a
integer, intent(in) :: nz, imin,imax,jmin,jmax
integer, intent(in) :: ia(:),ja(:)
real(psb_dpk_), intent(in) :: val(:)
integer, intent(out) :: info
integer, intent(in), optional :: gtl(:)
integer :: i,ir,ic, ilr, ilc, ip, &
& i1,i2,nc,nnz,dupl,ng
integer :: debug_level, debug_unit
character(len=20) :: name='d_coo_srch_upd'
info = 0
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
dupl = a%get_dupl()
if (.not.a%is_sorted()) then
info = -4
return
end if
ilr = -1
ilc = -1
nnz = a%get_nzeros()
if (present(gtl)) then
ng = size(gtl)
select case(dupl)
case(psbn_dupl_ovwrt_,psbn_dupl_err_)
! Overwrite.
! Cannot test for error, should have been caught earlier.
do i=1, nz
ir = ia(i)
ic = ja(i)
if ((ir >=1).and.(ir<=ng).and.(ic>=1).and.(ic<=ng)) then
ir = gtl(ir)
if ((ir > 0).and.(ir <= a%m)) then
ic = gtl(ic)
if (ir /= ilr) then
i1 = psb_ibsrch(ir,nnz,a%ia)
i2 = i1
@@ -1253,14 +1132,21 @@ contains
info = i
return
end if
end if
end do
case(psbn_dupl_add_)
! Add
do i=1, nz
ir = ia(i)
ic = ja(i)
else
if (debug_level >= psb_debug_serial_) &
& write(debug_unit,*) trim(name),&
& ': Discarding row that does not belong to us.'
endif
end if
end do
case(psbn_dupl_add_)
! Add
do i=1, nz
ir = ia(i)
ic = ja(i)
if ((ir >=1).and.(ir<=ng).and.(ic>=1).and.(ic<=ng)) then
ir = gtl(ir)
ic = gtl(ic)
if ((ir > 0).and.(ir <= a%m)) then
if (ir /= ilr) then
@@ -1289,19 +1175,107 @@ contains
info = i
return
end if
else
if (debug_level >= psb_debug_serial_) &
& write(debug_unit,*) trim(name),&
& ': Discarding row that does not belong to us.'
end if
end do
end if
end do
case default
info = -3
if (debug_level >= psb_debug_serial_) &
& write(debug_unit,*) trim(name),&
& ': Duplicate handling: ',dupl
end select
case default
info = -3
if (debug_level >= psb_debug_serial_) &
& write(debug_unit,*) trim(name),&
& ': Duplicate handling: ',dupl
end select
end if
else
end subroutine d_coo_srch_upd
select case(dupl)
case(psbn_dupl_ovwrt_,psbn_dupl_err_)
! Overwrite.
! Cannot test for error, should have been caught earlier.
do i=1, nz
ir = ia(i)
ic = ja(i)
if ((ir > 0).and.(ir <= a%m)) then
if (ir /= ilr) then
i1 = psb_ibsrch(ir,nnz,a%ia)
i2 = i1
do
if (i2+1 > nnz) exit
if (a%ia(i2+1) /= a%ia(i2)) exit
i2 = i2 + 1
end do
do
if (i1-1 < 1) exit
if (a%ia(i1-1) /= a%ia(i1)) exit
i1 = i1 - 1
end do
ilr = ir
else
i1 = 1
i2 = 1
end if
nc = i2-i1+1
ip = psb_issrch(ic,nc,a%ja(i1:i2))
if (ip>0) then
a%val(i1+ip-1) = val(i)
else
info = i
return
end if
end if
end do
case(psbn_dupl_add_)
! Add
do i=1, nz
ir = ia(i)
ic = ja(i)
if ((ir > 0).and.(ir <= a%m)) then
if (ir /= ilr) then
i1 = psb_ibsrch(ir,nnz,a%ia)
i2 = i1
do
if (i2+1 > nnz) exit
if (a%ia(i2+1) /= a%ia(i2)) exit
i2 = i2 + 1
end do
do
if (i1-1 < 1) exit
if (a%ia(i1-1) /= a%ia(i1)) exit
i1 = i1 - 1
end do
ilr = ir
else
i1 = 1
i2 = 1
end if
nc = i2-i1+1
ip = psb_issrch(ic,nc,a%ja(i1:i2))
if (ip>0) then
a%val(i1+ip-1) = a%val(i1+ip-1) + val(i)
else
info = i
return
end if
end if
end do
case default
info = -3
if (debug_level >= psb_debug_serial_) &
& write(debug_unit,*) trim(name),&
& ': Duplicate handling: ',dupl
end select
end if
end subroutine d_coo_srch_upd
end subroutine d_coo_csput_impl
@@ -1518,7 +1492,7 @@ subroutine d_fix_coo_impl(a,info,idir)
if (nza < 2) return
dupl_ = a%get_dupl()
!!$
call d_fix_coo_inner(nza,dupl_,a%ia,a%ja,a%val,i,info,idir_)
call a%set_sorted()
@@ -1526,188 +1500,6 @@ subroutine d_fix_coo_impl(a,info,idir)
call a%set_asb()
!!$ allocate(iaux(nza+2),stat=info)
!!$ if (info /= 0) return
!!$
!!$
!!$ select case(idir_)
!!$
!!$ case(0) ! Row major order
!!$
!!$ call msort_up(nza,a%ia(1),iaux(1),iret)
!!$ if (iret == 0) &
!!$ & call psb_ip_reord(nza,a%val,a%ia,a%ja,iaux)
!!$ i = 1
!!$ j = i
!!$ do while (i <= nza)
!!$ do while ((a%ia(j) == a%ia(i)))
!!$ j = j+1
!!$ if (j > nza) exit
!!$ enddo
!!$ nzl = j - i
!!$ call msort_up(nzl,a%ja(i),iaux(1),iret)
!!$ if (iret == 0) &
!!$ & call psb_ip_reord(nzl,a%val(i:i+nzl-1),&
!!$ & a%ia(i:i+nzl-1),a%ja(i:i+nzl-1),iaux)
!!$ i = j
!!$ enddo
!!$
!!$ i = 1
!!$ irw = a%ia(i)
!!$ icl = a%ja(i)
!!$ j = 1
!!$
!!$ select case(dupl_)
!!$ case(psbn_dupl_ovwrt_)
!!$
!!$ do
!!$ j = j + 1
!!$ if (j > nza) exit
!!$ if ((a%ia(j) == irw).and.(a%ja(j) == icl)) then
!!$ a%val(i) = a%val(j)
!!$ else
!!$ i = i+1
!!$ a%val(i) = a%val(j)
!!$ a%ia(i) = a%ia(j)
!!$ a%ja(i) = a%ja(j)
!!$ irw = a%ia(i)
!!$ icl = a%ja(i)
!!$ endif
!!$ enddo
!!$
!!$ case(psbn_dupl_add_)
!!$
!!$ do
!!$ j = j + 1
!!$ if (j > nza) exit
!!$ if ((a%ia(j) == irw).and.(a%ja(j) == icl)) then
!!$ a%val(i) = a%val(i) + a%val(j)
!!$ else
!!$ i = i+1
!!$ a%val(i) = a%val(j)
!!$ a%ia(i) = a%ia(j)
!!$ a%ja(i) = a%ja(j)
!!$ irw = a%ia(i)
!!$ icl = a%ja(i)
!!$ endif
!!$ enddo
!!$
!!$ case(psbn_dupl_err_)
!!$ do
!!$ j = j + 1
!!$ if (j > nza) exit
!!$ if ((a%ia(j) == irw).and.(a%ja(j) == icl)) then
!!$ call psb_errpush(130,name)
!!$ goto 9999
!!$ else
!!$ i = i+1
!!$ a%val(i) = a%val(j)
!!$ a%ia(i) = a%ia(j)
!!$ a%ja(i) = a%ja(j)
!!$ irw = a%ia(i)
!!$ icl = a%ja(i)
!!$ endif
!!$ enddo
!!$ case default
!!$ write(0,*) 'Error in fix_coo: unsafe dupl',dupl_
!!$
!!$ end select
!!$
!!$
!!$ if(debug_level >= psb_debug_serial_)&
!!$ & write(debug_unit,*) trim(name),': end second loop'
!!$
!!$ case(1) ! Col major order
!!$
!!$ call msort_up(nza,a%ja(1),iaux(1),iret)
!!$ if (iret == 0) &
!!$ & call psb_ip_reord(nza,a%val,a%ia,a%ja,iaux)
!!$ i = 1
!!$ j = i
!!$ do while (i <= nza)
!!$ do while ((a%ja(j) == a%ja(i)))
!!$ j = j+1
!!$ if (j > nza) exit
!!$ enddo
!!$ nzl = j - i
!!$ call msort_up(nzl,a%ia(i),iaux(1),iret)
!!$ if (iret == 0) &
!!$ & call psb_ip_reord(nzl,a%val(i:i+nzl-1),&
!!$ & a%ia(i:i+nzl-1),a%ja(i:i+nzl-1),iaux)
!!$ i = j
!!$ enddo
!!$
!!$ i = 1
!!$ irw = a%ia(i)
!!$ icl = a%ja(i)
!!$ j = 1
!!$
!!$
!!$ select case(dupl_)
!!$ case(psbn_dupl_ovwrt_)
!!$ do
!!$ j = j + 1
!!$ if (j > nza) exit
!!$ if ((a%ia(j) == irw).and.(a%ja(j) == icl)) then
!!$ a%val(i) = a%val(j)
!!$ else
!!$ i = i+1
!!$ a%val(i) = a%val(j)
!!$ a%ia(i) = a%ia(j)
!!$ a%ja(i) = a%ja(j)
!!$ irw = a%ia(i)
!!$ icl = a%ja(i)
!!$ endif
!!$ enddo
!!$
!!$ case(psbn_dupl_add_)
!!$ do
!!$ j = j + 1
!!$ if (j > nza) exit
!!$ if ((a%ia(j) == irw).and.(a%ja(j) == icl)) then
!!$ a%val(i) = a%val(i) + a%val(j)
!!$ else
!!$ i = i+1
!!$ a%val(i) = a%val(j)
!!$ a%ia(i) = a%ia(j)
!!$ a%ja(i) = a%ja(j)
!!$ irw = a%ia(i)
!!$ icl = a%ja(i)
!!$ endif
!!$ enddo
!!$
!!$ case(psbn_dupl_err_)
!!$ do
!!$ j = j + 1
!!$ if (j > nza) exit
!!$ if ((a%ia(j) == irw).and.(a%ja(j) == icl)) then
!!$ call psb_errpush(130,name)
!!$ goto 9999
!!$ else
!!$ i = i+1
!!$ a%val(i) = a%val(j)
!!$ a%ia(i) = a%ia(j)
!!$ a%ja(i) = a%ja(j)
!!$ irw = a%ia(i)
!!$ icl = a%ja(i)
!!$ endif
!!$ enddo
!!$ case default
!!$ write(0,*) 'Error in fix_coo: unsafe dupl',dupl_
!!$ end select
!!$ if (debug_level >= psb_debug_serial_)&
!!$ & write(debug_unit,*) trim(name),': end second loop'
!!$ case default
!!$ write(debug_unit,*) trim(name),': unknown direction ',idir_
!!$ end select
!!$
!!$
!!$ call a%set_sorted()
!!$ call a%set_nzeros(i)
!!$ call a%set_asb()
!!$
!!$ deallocate(iaux)
!!$
call psb_erractionrestore(err_act)
return
@@ -1719,8 +1511,6 @@ subroutine d_fix_coo_impl(a,info,idir)
end if
return
end subroutine d_fix_coo_impl
base/newserial/psbn_d_csr_impl.f03
+52
View File
@@ -1,4 +1,17 @@
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!
!
!
! Computational routines
!
!
!
!
!
!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
subroutine d_csr_csmv_impl(alpha,a,x,beta,y,info,trans)
use psb_error_mod
use psbn_d_csr_mat_mod, psb_protect_name => d_csr_csmv_impl
@@ -940,6 +953,45 @@ contains
end subroutine d_csr_cssm_impl
function d_csr_csnmi_impl(a) result(res)
use psb_error_mod
use psbn_d_csr_mat_mod, psb_protect_name => d_csr_csnmi_impl
implicit none
class(psbn_d_csr_sparse_mat), intent(in) :: a
real(psb_dpk_) :: res
integer :: i,j,k,m,n, nr, ir, jc, nc
real(psb_dpk_) :: acc
logical :: tra
Integer :: err_act
character(len=20) :: name='d_csnmi'
logical, parameter :: debug=.false.
res = dzero
do i = 1, a%get_nrows()
acc = dzero
do j=a%irp(i),a%irp(i+1)-1
acc = acc + abs(a%val(j))
end do
res = max(res,acc)
end do
end function d_csr_csnmi_impl
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!
!
!
! Data management
!
!
!
!
!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
subroutine d_csr_csput_impl(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl)
base/newserial/psbn_d_csr_mat_mod.f03
+284 -217
View File
@@ -8,29 +8,30 @@ module psbn_d_csr_mat_mod
real(psb_dpk_), allocatable :: val(:)
contains
procedure, pass(a) :: get_nzeros => d_csr_get_nzeros
procedure, pass(a) :: d_base_csmm => d_csr_csmm
procedure, pass(a) :: d_base_csmv => d_csr_csmv
procedure, pass(a) :: d_base_cssm => d_csr_cssm
procedure, pass(a) :: d_base_cssv => d_csr_cssv
procedure, pass(a) :: reallocate_nz => d_csr_reallocate_nz
procedure, pass(a) :: csput => d_csr_csput
procedure, pass(a) :: allocate_mnnz => d_csr_allocate_mnnz
procedure, pass(a) :: cp_to_coo => d_cp_csr_to_coo
procedure, pass(a) :: cp_from_coo => d_cp_csr_from_coo
procedure, pass(a) :: cp_to_fmt => d_cp_csr_to_fmt
procedure, pass(a) :: cp_from_fmt => d_cp_csr_from_fmt
procedure, pass(a) :: mv_to_coo => d_mv_csr_to_coo
procedure, pass(a) :: mv_from_coo => d_mv_csr_from_coo
procedure, pass(a) :: mv_to_fmt => d_mv_csr_to_fmt
procedure, pass(a) :: mv_from_fmt => d_mv_csr_from_fmt
procedure, pass(a) :: free => d_csr_free
procedure, pass(a) :: print => d_csr_print
procedure, pass(a) :: get_fmt => d_csr_get_fmt
procedure, pass(a) :: get_nzeros => d_csr_get_nzeros
procedure, pass(a) :: get_fmt => d_csr_get_fmt
procedure, pass(a) :: d_base_csmm => d_csr_csmm
procedure, pass(a) :: d_base_csmv => d_csr_csmv
procedure, pass(a) :: d_base_cssm => d_csr_cssm
procedure, pass(a) :: d_base_cssv => d_csr_cssv
procedure, pass(a) :: csnmi => d_csr_csnmi
procedure, pass(a) :: reallocate_nz => d_csr_reallocate_nz
procedure, pass(a) :: csput => d_csr_csput
procedure, pass(a) :: allocate_mnnz => d_csr_allocate_mnnz
procedure, pass(a) :: cp_to_coo => d_cp_csr_to_coo
procedure, pass(a) :: cp_from_coo => d_cp_csr_from_coo
procedure, pass(a) :: cp_to_fmt => d_cp_csr_to_fmt
procedure, pass(a) :: cp_from_fmt => d_cp_csr_from_fmt
procedure, pass(a) :: mv_to_coo => d_mv_csr_to_coo
procedure, pass(a) :: mv_from_coo => d_mv_csr_from_coo
procedure, pass(a) :: mv_to_fmt => d_mv_csr_to_fmt
procedure, pass(a) :: mv_from_fmt => d_mv_csr_from_fmt
procedure, pass(a) :: free => d_csr_free
procedure, pass(a) :: print => d_csr_print
end type psbn_d_csr_sparse_mat
private :: d_csr_get_nzeros, d_csr_csmm, d_csr_csmv, d_csr_cssm, d_csr_cssv, &
& d_csr_csput, d_csr_reallocate_nz, d_csr_allocate_mnnz, &
& d_csr_free, d_csr_print, d_csr_get_fmt, &
& d_csr_free, d_csr_print, d_csr_get_fmt, d_csr_csnmi, &
& d_cp_csr_to_coo, d_cp_csr_from_coo, &
& d_mv_csr_to_coo, d_mv_csr_from_coo, &
& d_cp_csr_to_fmt, d_cp_csr_from_fmt, &
@@ -181,10 +182,31 @@ module psbn_d_csr_mat_mod
end subroutine d_csr_csmm_impl
end interface
interface d_csr_csnmi_impl
function d_csr_csnmi_impl(a) result(res)
use psb_const_mod
import psbn_d_csr_sparse_mat
class(psbn_d_csr_sparse_mat), intent(in) :: a
real(psb_dpk_) :: res
end function d_csr_csnmi_impl
end interface
contains
!=====================================
!
!
!
! Getters
!
!
!
!
!
!=====================================
function d_csr_get_fmt(a) result(res)
implicit none
class(psbn_d_csr_sparse_mat), intent(in) :: a
@@ -192,6 +214,26 @@ contains
res = 'CSR'
end function d_csr_get_fmt
function d_csr_get_nzeros(a) result(res)
implicit none
class(psbn_d_csr_sparse_mat), intent(in) :: a
integer :: res
res = a%irp(a%m+1)-1
end function d_csr_get_nzeros
!=====================================
!
!
!
! Data management
!
!
!
!
!
!=====================================
subroutine d_csr_reallocate_nz(nz,a)
use psb_error_mod
@@ -227,14 +269,6 @@ contains
end subroutine d_csr_reallocate_nz
function d_csr_get_nzeros(a) result(res)
implicit none
class(psbn_d_csr_sparse_mat), intent(in) :: a
integer :: res
res = a%irp(a%m+1)-1
end function d_csr_get_nzeros
subroutine d_csr_csput(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl)
use psb_const_mod
use psb_error_mod
@@ -299,195 +333,6 @@ contains
end subroutine d_csr_csput
subroutine d_csr_csmv(alpha,a,x,beta,y,info,trans)
use psb_error_mod
implicit none
class(psbn_d_csr_sparse_mat), intent(in) :: a
real(psb_dpk_), intent(in) :: alpha, beta, x(:)
real(psb_dpk_), intent(inout) :: y(:)
integer, intent(out) :: info
character, optional, intent(in) :: trans
character :: trans_
integer :: i,j,k,m,n, nnz, ir, jc
real(psb_dpk_) :: acc
logical :: tra
Integer :: err_act
character(len=20) :: name='d_csr_csmv'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
if (.not.a%is_asb()) then
info = 1121
call psb_errpush(info,name)
goto 9999
endif
call d_csr_csmm_impl(alpha,a,x,beta,y,info,trans)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_csr_csmv
subroutine d_csr_csmm(alpha,a,x,beta,y,info,trans)
use psb_error_mod
implicit none
class(psbn_d_csr_sparse_mat), intent(in) :: a
real(psb_dpk_), intent(in) :: alpha, beta, x(:,:)
real(psb_dpk_), intent(inout) :: y(:,:)
integer, intent(out) :: info
character, optional, intent(in) :: trans
character :: trans_
integer :: i,j,k,m,n, nnz, ir, jc, nc
real(psb_dpk_), allocatable :: acc(:)
logical :: tra
Integer :: err_act
character(len=20) :: name='d_csr_csmm'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
call d_csr_csmm_impl(alpha,a,x,beta,y,info,trans)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_csr_csmm
subroutine d_csr_cssv(alpha,a,x,beta,y,info,trans)
use psb_error_mod
implicit none
class(psbn_d_csr_sparse_mat), intent(in) :: a
real(psb_dpk_), intent(in) :: alpha, beta, x(:)
real(psb_dpk_), intent(inout) :: y(:)
integer, intent(out) :: info
character, optional, intent(in) :: trans
character :: trans_
integer :: i,j,k,m,n, nnz, ir, jc
real(psb_dpk_) :: acc
real(psb_dpk_), allocatable :: tmp(:)
logical :: tra
Integer :: err_act
character(len=20) :: name='d_csr_cssv'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
if (.not.a%is_asb()) then
info = 1121
call psb_errpush(info,name)
goto 9999
endif
if (.not. (a%is_triangle())) then
write(0,*) 'Called SM on a non-triangular mat!'
info = 1121
call psb_errpush(info,name)
goto 9999
end if
call d_csr_cssm_impl(alpha,a,x,beta,y,info,trans)
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_csr_cssv
subroutine d_csr_cssm(alpha,a,x,beta,y,info,trans)
use psb_error_mod
implicit none
class(psbn_d_csr_sparse_mat), intent(in) :: a
real(psb_dpk_), intent(in) :: alpha, beta, x(:,:)
real(psb_dpk_), intent(inout) :: y(:,:)
integer, intent(out) :: info
character, optional, intent(in) :: trans
character :: trans_
integer :: i,j,k,m,n, nnz, ir, jc, nc
real(psb_dpk_) :: acc
real(psb_dpk_), allocatable :: tmp(:,:)
logical :: tra
Integer :: err_act
character(len=20) :: name='d_csr_csmm'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
if (.not.a%is_asb()) then
info = 1121
call psb_errpush(info,name)
goto 9999
endif
if (.not. (a%is_triangle())) then
write(0,*) 'Called SM on a non-triangular mat!'
info = 1121
call psb_errpush(info,name)
goto 9999
end if
call d_csr_cssm_impl(alpha,a,x,beta,y,info,trans)
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_csr_cssm
subroutine d_csr_free(a)
implicit none
@@ -907,4 +752,226 @@ contains
end subroutine d_csr_print
!=====================================
!
!
!
! Computational routines
!
!
!
!
!
!
!=====================================
subroutine d_csr_csmv(alpha,a,x,beta,y,info,trans)
use psb_error_mod
implicit none
class(psbn_d_csr_sparse_mat), intent(in) :: a
real(psb_dpk_), intent(in) :: alpha, beta, x(:)
real(psb_dpk_), intent(inout) :: y(:)
integer, intent(out) :: info
character, optional, intent(in) :: trans
character :: trans_
integer :: i,j,k,m,n, nnz, ir, jc
real(psb_dpk_) :: acc
logical :: tra
Integer :: err_act
character(len=20) :: name='d_csr_csmv'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
if (.not.a%is_asb()) then
info = 1121
call psb_errpush(info,name)
goto 9999
endif
call d_csr_csmm_impl(alpha,a,x,beta,y,info,trans)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_csr_csmv
subroutine d_csr_csmm(alpha,a,x,beta,y,info,trans)
use psb_error_mod
implicit none
class(psbn_d_csr_sparse_mat), intent(in) :: a
real(psb_dpk_), intent(in) :: alpha, beta, x(:,:)
real(psb_dpk_), intent(inout) :: y(:,:)
integer, intent(out) :: info
character, optional, intent(in) :: trans
character :: trans_
integer :: i,j,k,m,n, nnz, ir, jc, nc
real(psb_dpk_), allocatable :: acc(:)
logical :: tra
Integer :: err_act
character(len=20) :: name='d_csr_csmm'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
call d_csr_csmm_impl(alpha,a,x,beta,y,info,trans)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_csr_csmm
subroutine d_csr_cssv(alpha,a,x,beta,y,info,trans)
use psb_error_mod
implicit none
class(psbn_d_csr_sparse_mat), intent(in) :: a
real(psb_dpk_), intent(in) :: alpha, beta, x(:)
real(psb_dpk_), intent(inout) :: y(:)
integer, intent(out) :: info
character, optional, intent(in) :: trans
character :: trans_
integer :: i,j,k,m,n, nnz, ir, jc
real(psb_dpk_) :: acc
real(psb_dpk_), allocatable :: tmp(:)
logical :: tra
Integer :: err_act
character(len=20) :: name='d_csr_cssv'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
if (.not.a%is_asb()) then
info = 1121
call psb_errpush(info,name)
goto 9999
endif
if (.not. (a%is_triangle())) then
write(0,*) 'Called SM on a non-triangular mat!'
info = 1121
call psb_errpush(info,name)
goto 9999
end if
call d_csr_cssm_impl(alpha,a,x,beta,y,info,trans)
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_csr_cssv
subroutine d_csr_cssm(alpha,a,x,beta,y,info,trans)
use psb_error_mod
implicit none
class(psbn_d_csr_sparse_mat), intent(in) :: a
real(psb_dpk_), intent(in) :: alpha, beta, x(:,:)
real(psb_dpk_), intent(inout) :: y(:,:)
integer, intent(out) :: info
character, optional, intent(in) :: trans
character :: trans_
integer :: i,j,k,m,n, nnz, ir, jc, nc
real(psb_dpk_) :: acc
real(psb_dpk_), allocatable :: tmp(:,:)
logical :: tra
Integer :: err_act
character(len=20) :: name='d_csr_csmm'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
if (.not.a%is_asb()) then
info = 1121
call psb_errpush(info,name)
goto 9999
endif
if (.not. (a%is_triangle())) then
write(0,*) 'Called SM on a non-triangular mat!'
info = 1121
call psb_errpush(info,name)
goto 9999
end if
call d_csr_cssm_impl(alpha,a,x,beta,y,info,trans)
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_csr_cssm
function d_csr_csnmi(a) result(res)
use psb_error_mod
use psb_const_mod
implicit none
class(psbn_d_csr_sparse_mat), intent(in) :: a
real(psb_dpk_) :: res
Integer :: err_act
character(len=20) :: name='csnmi'
logical, parameter :: debug=.false.
res = d_csr_csnmi_impl(a)
return
end function d_csr_csnmi
end module psbn_d_csr_mat_mod
base/newserial/psbn_mat_impl.f03
+8 -3
View File
@@ -125,7 +125,7 @@ subroutine psbn_d_spcnv(a,b,info,type,mold,upd,dupl)
call psb_errpush(info,name)
goto 9999
end if
call altmp%cp_from_fmt(a%a, info)
if (info /= 0) then
@@ -208,15 +208,20 @@ subroutine psbn_d_spcnv_ip(a,info,type,mold,dupl)
else
allocate(psbn_d_csr_sparse_mat :: altmp, stat=info)
end if
if (info /= 0) then
info = 4000
call psb_errpush(info,name)
goto 9999
end if
call altmp%mv_from_fmt(a%a, info)
if (allocated(altmp)) then
call altmp%mv_from_fmt(a%a, info)
else
write(0,*) 'Unallocated altmp??'
info = -1
end if
if (info /= 0) then
info = 4010
call psb_errpush(info,name,a_err="mv_from")
base/newserial/psbn_mat_mod.f03
+55 -25
View File
@@ -54,6 +54,7 @@ module psbn_d_mat_mod
procedure, pass(a) :: print => sparse_print
! Computational routines
procedure, pass(a) :: csnmi
procedure, pass(a) :: d_csmv
procedure, pass(a) :: d_csmm
generic, public :: csmm => d_csmm, d_csmv
@@ -71,12 +72,12 @@ module psbn_d_mat_mod
& reallocate_nz, free, d_csmv, d_csmm, d_cssv, d_cssm, sparse_print, &
& set_nrows, set_ncols, set_dupl, set_state, set_null, set_bld, &
& set_upd, set_asb, set_sorted, set_upper, set_lower, set_triangle, &
& set_unit
& set_unit, csnmi
contains
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!=====================================
!
!
!
@@ -86,7 +87,7 @@ contains
!
!
!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!=====================================
function sparse_get_fmt(a) result(res)
@@ -342,7 +343,7 @@ contains
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!=====================================
!
!
!
@@ -353,7 +354,7 @@ contains
!
!
!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!=====================================
subroutine set_nrows(m,a)
@@ -762,7 +763,7 @@ contains
end subroutine set_upper
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!=====================================
!
!
!
@@ -772,7 +773,7 @@ contains
!
!
!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!=====================================
@@ -793,7 +794,7 @@ contains
logical, parameter :: debug=.false.
info = 0
call psb_erractionsave(err_act)
call psb_get_erraction(err_act)
if (.not.allocated(a%a)) then
info = 1121
call psb_errpush(info,name)
@@ -802,12 +803,9 @@ contains
call a%a%print(iout,iv,eirs,eics,head,ivr,ivc)
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
@@ -873,8 +871,8 @@ contains
character(len=20) :: name='csall'
logical, parameter :: debug=.false.
call psb_get_erraction(err_act)
call psb_erractionsave(err_act)
info = 0
allocate(psbn_d_coo_sparse_mat :: a%a, stat=info)
if (info /= 0) then
@@ -885,11 +883,9 @@ contains
call a%a%allocate(nr,nc,nz)
call a%set_bld()
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
@@ -916,9 +912,6 @@ contains
call a%a%reallocate(nz)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
@@ -940,7 +933,6 @@ contains
character(len=20) :: name='free'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
if (.not.allocated(a%a)) then
info = 1121
call psb_errpush(info,name)
@@ -949,7 +941,6 @@ contains
call a%a%free()
call psb_erractionrestore(err_act)
return
9999 continue
@@ -1113,6 +1104,9 @@ contains
goto 9999
end if
if (debug) write(0,*) 'Converting from ',&
& a%get_fmt(),' to ',altmp%get_fmt()
call altmp%cp_from_fmt(a%a, info)
if (info /= 0) then
@@ -1202,6 +1196,9 @@ contains
goto 9999
end if
if (debug) write(0,*) 'Converting in-place from ',&
& a%get_fmt(),' to ',altmp%get_fmt()
call altmp%mv_from_fmt(a%a, info)
if (info /= 0) then
@@ -1229,9 +1226,7 @@ contains
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!=====================================
!
!
!
@@ -1242,7 +1237,7 @@ contains
!
!
!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!=====================================
subroutine d_csmm(alpha,a,x,beta,y,info,trans)
@@ -1266,7 +1261,7 @@ contains
endif
call a%a%csmm(alpha,x,beta,y,info,trans)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
@@ -1302,7 +1297,7 @@ contains
endif
call a%a%csmm(alpha,x,beta,y,info,trans)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
@@ -1338,6 +1333,7 @@ contains
endif
call a%a%cssm(alpha,x,beta,y,info,trans)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
@@ -1375,6 +1371,7 @@ contains
call a%a%cssm(alpha,x,beta,y,info,trans)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
@@ -1390,5 +1387,38 @@ contains
end subroutine d_cssv
function csnmi(a) result(res)
use psb_error_mod
use psb_const_mod
implicit none
class(psbn_d_sparse_mat), intent(in) :: a
real(psb_dpk_) :: res
Integer :: err_act, info
character(len=20) :: name='csnmi'
logical, parameter :: debug=.false.
if (.not.allocated(a%a)) then
info = 1121
call psb_errpush(info,name)
goto 9999
endif
res = a%a%csnmi()
return
9999 continue
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end function csnmi
end module psbn_d_mat_mod
test/serial/Makefile
+2 -1
View File
@@ -38,7 +38,8 @@ spde: spde.o
clean:
/bin/rm -f d_coo_matgen.o d_matgen.o tpg.o ppde.o spde.o $(EXEDIR)/ppde
/bin/rm -f d_coo_matgen.o d_matgen.o tpg.o ppde.o spde.o \
psbn_d_cxx_mat_mod.o psbn_d_cxx_impl.o $(EXEDIR)/ppde
verycleanlib:
(cd ../..; make veryclean)
lib:
test/serial/d_matgen.f03
+7 -3
View File
@@ -161,7 +161,7 @@ contains
type(psbn_d_cxx_sparse_mat) :: acxx
! deltah dimension of each grid cell
! deltat discretization time
real(psb_dpk_) :: deltah
real(psb_dpk_) :: deltah, anorm
real(psb_dpk_),parameter :: rhs=0.d0,one=1.d0,zero=0.d0
real(psb_dpk_) :: t0, t1, t2, t3, tasb, talc, ttot, tgen, tcpy, tmov
real(psb_dpk_) :: a1, a2, a3, a4, b1, b2, b3
@@ -172,7 +172,7 @@ contains
info = 0
name = 'create_matrix'
call psb_erractionsave(err_act)
!!$ call psb_erractionsave(err_act)
call psb_info(ictxt, iam, np)
@@ -380,6 +380,8 @@ contains
end if
tasb = psb_wtime()-t1
call a_n%print(20)
anorm = a_n%csnmi()
write(0,*) 'Nrm infinity ',anorm
!!$
t1 = psb_wtime()
call a_n%cscnv(info,mold=acxx)
@@ -392,6 +394,8 @@ contains
end if
tmov = psb_wtime()-t1
call a_n%print(21)
anorm = a_n%csnmi()
write(0,*) 'Nrm infinity ',anorm
!!$
if(iam == psb_root_) then
@@ -404,7 +408,7 @@ contains
!!$ write(*,'("-total time : ",es12.5)') ttot
end if
call psb_erractionrestore(err_act)
!!$ call psb_erractionrestore(err_act)
return
9999 continue
test/serial/psbn_d_cxx_impl.f03
+87 -33
View File
@@ -1,4 +1,17 @@
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!
!
!
! Computational routines
!
!
!
!
!
!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
subroutine d_cxx_csmv_impl(alpha,a,x,beta,y,info,trans)
use psb_error_mod
use psbn_d_cxx_mat_mod, psb_protect_name => d_cxx_csmv_impl
@@ -940,6 +953,45 @@ contains
end subroutine d_cxx_cssm_impl
function d_cxx_csnmi_impl(a) result(res)
use psb_error_mod
use psbn_d_cxx_mat_mod, psb_protect_name => d_cxx_csnmi_impl
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
real(psb_dpk_) :: res
integer :: i,j,k,m,n, nr, ir, jc, nc
real(psb_dpk_) :: acc
logical :: tra
Integer :: err_act
character(len=20) :: name='d_csnmi'
logical, parameter :: debug=.false.
res = dzero
do i = 1, a%get_nrows()
acc = dzero
do j=a%irp(i),a%irp(i+1)-1
acc = acc + abs(a%val(j))
end do
res = max(res,acc)
end do
end function d_cxx_csnmi_impl
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!
!
!
! Data management
!
!
!
!
!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
subroutine d_cxx_csput_impl(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl)
@@ -1449,7 +1501,6 @@ subroutine d_mv_cxx_to_fmt_impl(a,b,info)
end subroutine d_mv_cxx_to_fmt_impl
subroutine d_cp_cxx_to_fmt_impl(a,b,info)
use psb_const_mod
use psb_realloc_mod
@@ -1471,6 +1522,7 @@ subroutine d_cp_cxx_to_fmt_impl(a,b,info)
info = 0
select type (b)
class is (psbn_d_coo_sparse_mat)
call a%cp_to_coo(b,info)
@@ -1478,9 +1530,43 @@ subroutine d_cp_cxx_to_fmt_impl(a,b,info)
call tmp%cp_from_fmt(a,info)
if (info == 0) call b%mv_from_coo(tmp,info)
end select
end subroutine d_cp_cxx_to_fmt_impl
subroutine d_mv_cxx_from_fmt_impl(a,b,info)
use psb_const_mod
use psb_realloc_mod
use psbn_d_base_mat_mod
use psbn_d_cxx_mat_mod, psb_protect_name => d_mv_cxx_from_fmt_impl
implicit none
class(psbn_d_cxx_sparse_mat), intent(inout) :: a
class(psbn_d_base_sparse_mat), intent(inout) :: b
integer, intent(out) :: info
!locals
type(psbn_d_coo_sparse_mat) :: tmp
logical :: rwshr_
Integer :: nza, nr, i,j,irw, idl,err_act, nc
Integer, Parameter :: maxtry=8
integer :: debug_level, debug_unit
character(len=20) :: name
info = 0
select type (b)
class is (psbn_d_coo_sparse_mat)
call a%mv_from_coo(b,info)
class default
call tmp%mv_from_fmt(b,info)
if (info == 0) call a%mv_from_coo(tmp,info)
end select
end subroutine d_mv_cxx_from_fmt_impl
subroutine d_cp_cxx_from_fmt_impl(a,b,info)
use psb_const_mod
use psb_realloc_mod
@@ -1510,35 +1596,3 @@ subroutine d_cp_cxx_from_fmt_impl(a,b,info)
if (info == 0) call a%mv_from_coo(tmp,info)
end select
end subroutine d_cp_cxx_from_fmt_impl
subroutine d_mv_cxx_from_fmt_impl(a,b,info)
use psb_const_mod
use psb_realloc_mod
use psbn_d_base_mat_mod
use psbn_d_cxx_mat_mod, psb_protect_name => d_mv_cxx_from_fmt_impl
implicit none
class(psbn_d_cxx_sparse_mat), intent(inout) :: a
class(psbn_d_base_sparse_mat), intent(inout) :: b
integer, intent(out) :: info
!locals
type(psbn_d_coo_sparse_mat) :: tmp
logical :: rwshr_
Integer :: nza, nr, i,j,irw, idl,err_act, nc
Integer, Parameter :: maxtry=8
integer :: debug_level, debug_unit
character(len=20) :: name
info = 0
select type (b)
class is (psbn_d_coo_sparse_mat)
call a%mv_from_coo(b,info)
class default
call tmp%mv_from_fmt(b,info)
if (info == 0) call a%mv_from_coo(tmp,info)
end select
end subroutine d_mv_cxx_from_fmt_impl
test/serial/psbn_d_cxx_mat_mod.f03
+317 -250
View File
@@ -8,29 +8,30 @@ module psbn_d_cxx_mat_mod
real(psb_dpk_), allocatable :: val(:)
contains
procedure, pass(a) :: get_nzeros => d_cxx_get_nzeros
procedure, pass(a) :: d_base_csmm => d_cxx_csmm
procedure, pass(a) :: d_base_csmv => d_cxx_csmv
procedure, pass(a) :: d_base_cssm => d_cxx_cssm
procedure, pass(a) :: d_base_cssv => d_cxx_cssv
procedure, pass(a) :: reallocate_nz => d_cxx_reallocate_nz
procedure, pass(a) :: csput => d_cxx_csput
procedure, pass(a) :: allocate_mnnz => d_cxx_allocate_mnnz
procedure, pass(a) :: cp_to_coo => d_cp_cxx_to_coo
procedure, pass(a) :: cp_from_coo => d_cp_cxx_from_coo
procedure, pass(a) :: cp_to_fmt => d_cp_cxx_to_fmt
procedure, pass(a) :: cp_from_fmt => d_cp_cxx_from_fmt
procedure, pass(a) :: mv_to_coo => d_mv_cxx_to_coo
procedure, pass(a) :: mv_from_coo => d_mv_cxx_from_coo
procedure, pass(a) :: mv_to_fmt => d_mv_cxx_to_fmt
procedure, pass(a) :: mv_from_fmt => d_mv_cxx_from_fmt
procedure, pass(a) :: free => d_cxx_free
procedure, pass(a) :: print => d_cxx_print
procedure, pass(a) :: get_fmt => d_cxx_get_fmt
procedure, pass(a) :: get_nzeros => d_cxx_get_nzeros
procedure, pass(a) :: get_fmt => d_cxx_get_fmt
procedure, pass(a) :: d_base_csmm => d_cxx_csmm
procedure, pass(a) :: d_base_csmv => d_cxx_csmv
procedure, pass(a) :: d_base_cssm => d_cxx_cssm
procedure, pass(a) :: d_base_cssv => d_cxx_cssv
procedure, pass(a) :: csnmi => d_cxx_csnmi
procedure, pass(a) :: reallocate_nz => d_cxx_reallocate_nz
procedure, pass(a) :: csput => d_cxx_csput
procedure, pass(a) :: allocate_mnnz => d_cxx_allocate_mnnz
procedure, pass(a) :: cp_to_coo => d_cp_cxx_to_coo
procedure, pass(a) :: cp_from_coo => d_cp_cxx_from_coo
procedure, pass(a) :: cp_to_fmt => d_cp_cxx_to_fmt
procedure, pass(a) :: cp_from_fmt => d_cp_cxx_from_fmt
procedure, pass(a) :: mv_to_coo => d_mv_cxx_to_coo
procedure, pass(a) :: mv_from_coo => d_mv_cxx_from_coo
procedure, pass(a) :: mv_to_fmt => d_mv_cxx_to_fmt
procedure, pass(a) :: mv_from_fmt => d_mv_cxx_from_fmt
procedure, pass(a) :: free => d_cxx_free
procedure, pass(a) :: print => d_cxx_print
end type psbn_d_cxx_sparse_mat
private :: d_cxx_get_nzeros, d_cxx_csmm, d_cxx_csmv, d_cxx_cssm, d_cxx_cssv, &
& d_cxx_csput, d_cxx_reallocate_nz, d_cxx_allocate_mnnz, &
& d_cxx_free, d_cxx_print, d_cxx_get_fmt, &
& d_cxx_free, d_cxx_print, d_cxx_get_fmt, d_cxx_csnmi, &
& d_cp_cxx_to_coo, d_cp_cxx_from_coo, &
& d_mv_cxx_to_coo, d_mv_cxx_from_coo, &
& d_cp_cxx_to_fmt, d_cp_cxx_from_fmt, &
@@ -181,10 +182,31 @@ module psbn_d_cxx_mat_mod
end subroutine d_cxx_csmm_impl
end interface
interface d_cxx_csnmi_impl
function d_cxx_csnmi_impl(a) result(res)
use psb_const_mod
import psbn_d_cxx_sparse_mat
class(psbn_d_cxx_sparse_mat), intent(in) :: a
real(psb_dpk_) :: res
end function d_cxx_csnmi_impl
end interface
contains
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!
!
!
! Getters
!
!
!
!
!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
function d_cxx_get_fmt(a) result(res)
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
@@ -192,6 +214,26 @@ contains
res = 'CXX'
end function d_cxx_get_fmt
function d_cxx_get_nzeros(a) result(res)
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
integer :: res
res = a%irp(a%m+1)-1
end function d_cxx_get_nzeros
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!
!
!
! Data management
!
!
!
!
!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
subroutine d_cxx_reallocate_nz(nz,a)
use psb_error_mod
@@ -227,14 +269,6 @@ contains
end subroutine d_cxx_reallocate_nz
function d_cxx_get_nzeros(a) result(res)
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
integer :: res
res = a%irp(a%m+1)-1
end function d_cxx_get_nzeros
subroutine d_cxx_csput(nz,val,ia,ja,a,imin,imax,jmin,jmax,info,gtl)
use psb_const_mod
use psb_error_mod
@@ -299,195 +333,6 @@ contains
end subroutine d_cxx_csput
subroutine d_cxx_csmv(alpha,a,x,beta,y,info,trans)
use psb_error_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
real(psb_dpk_), intent(in) :: alpha, beta, x(:)
real(psb_dpk_), intent(inout) :: y(:)
integer, intent(out) :: info
character, optional, intent(in) :: trans
character :: trans_
integer :: i,j,k,m,n, nnz, ir, jc
real(psb_dpk_) :: acc
logical :: tra
Integer :: err_act
character(len=20) :: name='d_cxx_csmv'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
if (.not.a%is_asb()) then
info = 1121
call psb_errpush(info,name)
goto 9999
endif
call d_cxx_csmm_impl(alpha,a,x,beta,y,info,trans)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_cxx_csmv
subroutine d_cxx_csmm(alpha,a,x,beta,y,info,trans)
use psb_error_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
real(psb_dpk_), intent(in) :: alpha, beta, x(:,:)
real(psb_dpk_), intent(inout) :: y(:,:)
integer, intent(out) :: info
character, optional, intent(in) :: trans
character :: trans_
integer :: i,j,k,m,n, nnz, ir, jc, nc
real(psb_dpk_), allocatable :: acc(:)
logical :: tra
Integer :: err_act
character(len=20) :: name='d_cxx_csmm'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
call d_cxx_csmm_impl(alpha,a,x,beta,y,info,trans)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_cxx_csmm
subroutine d_cxx_cssv(alpha,a,x,beta,y,info,trans)
use psb_error_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
real(psb_dpk_), intent(in) :: alpha, beta, x(:)
real(psb_dpk_), intent(inout) :: y(:)
integer, intent(out) :: info
character, optional, intent(in) :: trans
character :: trans_
integer :: i,j,k,m,n, nnz, ir, jc
real(psb_dpk_) :: acc
real(psb_dpk_), allocatable :: tmp(:)
logical :: tra
Integer :: err_act
character(len=20) :: name='d_cxx_cssv'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
if (.not.a%is_asb()) then
info = 1121
call psb_errpush(info,name)
goto 9999
endif
if (.not. (a%is_triangle())) then
write(0,*) 'Called SM on a non-triangular mat!'
info = 1121
call psb_errpush(info,name)
goto 9999
end if
call d_cxx_cssm_impl(alpha,a,x,beta,y,info,trans)
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_cxx_cssv
subroutine d_cxx_cssm(alpha,a,x,beta,y,info,trans)
use psb_error_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
real(psb_dpk_), intent(in) :: alpha, beta, x(:,:)
real(psb_dpk_), intent(inout) :: y(:,:)
integer, intent(out) :: info
character, optional, intent(in) :: trans
character :: trans_
integer :: i,j,k,m,n, nnz, ir, jc, nc
real(psb_dpk_) :: acc
real(psb_dpk_), allocatable :: tmp(:,:)
logical :: tra
Integer :: err_act
character(len=20) :: name='d_cxx_csmm'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
if (.not.a%is_asb()) then
info = 1121
call psb_errpush(info,name)
goto 9999
endif
if (.not. (a%is_triangle())) then
write(0,*) 'Called SM on a non-triangular mat!'
info = 1121
call psb_errpush(info,name)
goto 9999
end if
call d_cxx_cssm_impl(alpha,a,x,beta,y,info,trans)
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_cxx_cssm
subroutine d_cxx_free(a)
implicit none
@@ -601,39 +446,6 @@ contains
return
end subroutine d_cp_cxx_to_fmt
subroutine d_mv_cxx_to_coo(a,b,info)
use psb_error_mod
use psb_realloc_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(inout) :: a
class(psbn_d_coo_sparse_mat), intent(out) :: b
integer, intent(out) :: info
Integer :: err_act
character(len=20) :: name='to_coo'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
info = 0
call d_mv_cxx_to_coo_impl(a,b,info)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
call psb_errpush(info,name)
if (err_act /= psb_act_ret_) then
call psb_error()
end if
return
end subroutine d_mv_cxx_to_coo
subroutine d_cp_cxx_from_fmt(a,b,info)
use psb_error_mod
@@ -667,6 +479,39 @@ contains
end subroutine d_cp_cxx_from_fmt
subroutine d_mv_cxx_to_coo(a,b,info)
use psb_error_mod
use psb_realloc_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(inout) :: a
class(psbn_d_coo_sparse_mat), intent(out) :: b
integer, intent(out) :: info
Integer :: err_act
character(len=20) :: name='to_coo'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
info = 0
call d_mv_cxx_to_coo_impl(a,b,info)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
call psb_errpush(info,name)
if (err_act /= psb_act_ret_) then
call psb_error()
end if
return
end subroutine d_mv_cxx_to_coo
subroutine d_mv_cxx_from_coo(a,b,info)
use psb_error_mod
use psb_realloc_mod
@@ -907,4 +752,226 @@ contains
end subroutine d_cxx_print
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!
!
!
! Computational routines
!
!
!
!
!
!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
subroutine d_cxx_csmv(alpha,a,x,beta,y,info,trans)
use psb_error_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
real(psb_dpk_), intent(in) :: alpha, beta, x(:)
real(psb_dpk_), intent(inout) :: y(:)
integer, intent(out) :: info
character, optional, intent(in) :: trans
character :: trans_
integer :: i,j,k,m,n, nnz, ir, jc
real(psb_dpk_) :: acc
logical :: tra
Integer :: err_act
character(len=20) :: name='d_cxx_csmv'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
if (.not.a%is_asb()) then
info = 1121
call psb_errpush(info,name)
goto 9999
endif
call d_cxx_csmm_impl(alpha,a,x,beta,y,info,trans)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_cxx_csmv
subroutine d_cxx_csmm(alpha,a,x,beta,y,info,trans)
use psb_error_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
real(psb_dpk_), intent(in) :: alpha, beta, x(:,:)
real(psb_dpk_), intent(inout) :: y(:,:)
integer, intent(out) :: info
character, optional, intent(in) :: trans
character :: trans_
integer :: i,j,k,m,n, nnz, ir, jc, nc
real(psb_dpk_), allocatable :: acc(:)
logical :: tra
Integer :: err_act
character(len=20) :: name='d_cxx_csmm'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
call d_cxx_csmm_impl(alpha,a,x,beta,y,info,trans)
if (info /= 0) goto 9999
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_cxx_csmm
subroutine d_cxx_cssv(alpha,a,x,beta,y,info,trans)
use psb_error_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
real(psb_dpk_), intent(in) :: alpha, beta, x(:)
real(psb_dpk_), intent(inout) :: y(:)
integer, intent(out) :: info
character, optional, intent(in) :: trans
character :: trans_
integer :: i,j,k,m,n, nnz, ir, jc
real(psb_dpk_) :: acc
real(psb_dpk_), allocatable :: tmp(:)
logical :: tra
Integer :: err_act
character(len=20) :: name='d_cxx_cssv'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
if (.not.a%is_asb()) then
info = 1121
call psb_errpush(info,name)
goto 9999
endif
if (.not. (a%is_triangle())) then
write(0,*) 'Called SM on a non-triangular mat!'
info = 1121
call psb_errpush(info,name)
goto 9999
end if
call d_cxx_cssm_impl(alpha,a,x,beta,y,info,trans)
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_cxx_cssv
subroutine d_cxx_cssm(alpha,a,x,beta,y,info,trans)
use psb_error_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
real(psb_dpk_), intent(in) :: alpha, beta, x(:,:)
real(psb_dpk_), intent(inout) :: y(:,:)
integer, intent(out) :: info
character, optional, intent(in) :: trans
character :: trans_
integer :: i,j,k,m,n, nnz, ir, jc, nc
real(psb_dpk_) :: acc
real(psb_dpk_), allocatable :: tmp(:,:)
logical :: tra
Integer :: err_act
character(len=20) :: name='d_cxx_csmm'
logical, parameter :: debug=.false.
call psb_erractionsave(err_act)
if (.not.a%is_asb()) then
info = 1121
call psb_errpush(info,name)
goto 9999
endif
if (.not. (a%is_triangle())) then
write(0,*) 'Called SM on a non-triangular mat!'
info = 1121
call psb_errpush(info,name)
goto 9999
end if
call d_cxx_cssm_impl(alpha,a,x,beta,y,info,trans)
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act == psb_act_abort_) then
call psb_error()
return
end if
return
end subroutine d_cxx_cssm
function d_cxx_csnmi(a) result(res)
use psb_error_mod
use psb_const_mod
implicit none
class(psbn_d_cxx_sparse_mat), intent(in) :: a
real(psb_dpk_) :: res
Integer :: err_act
character(len=20) :: name='csnmi'
logical, parameter :: debug=.false.
res = d_cxx_csnmi_impl(a)
return
end function d_cxx_csnmi
end module psbn_d_cxx_mat_mod