fdurastante/psblas3
1
0
Fork 0
mirror of https://github.com/sfilippone/psblas3.git synced 2026-10-20 04:42:08 +00:00
Code Issues Packages Projects Releases Wiki Activity

Modified interface to GEALL, GEINS et al.

Browse Source 
Changed ITRACE behaviour in solver routines.
Updated documentation.
This commit is contained in:
Salvatore Filippone
2006-04-24 12:02:04 +00:00
parent 41d7626b66
commit d003f55d95
48 changed files with 18464 additions and 17667 deletions
Download Patch File Download Diff File Expand all files Collapse all files
src/comm/psb_dovrl.f90
+18 -30
View File
@@ -40,10 +40,9 @@
! jx - integer(optional). The starting column of the global matrix.
! ik - integer(optional). The number of columns to gather.
! work - real(optional). A working area.
! choice - logical(optional). ???.
! update_type - integer(optional). ???.
! update - integer(optional). ???.
!
subroutine psb_dovrlm(x,desc_a,info,jx,ik,work,choice,update_type)
subroutine psb_dovrlm(x,desc_a,info,jx,ik,work,update)
use psb_descriptor_type
use psb_const_mod
use psi_mod
@@ -56,15 +55,14 @@ subroutine psb_dovrlm(x,desc_a,info,jx,ik,work,choice,update_type)
type(psb_desc_type), intent(in) :: desc_a
integer, intent(out) :: info
real(kind(1.d0)), optional, target :: work(:)
logical, intent(in), optional :: choice
integer, intent(in), optional :: update_type,jx,ik
integer, intent(in), optional :: update,jx,ik
! locals
integer :: int_err(5), icontxt, nprow, npcol, myrow, mycol,&
& err_act, m, n, iix, jjx, temp(2), ix, ijx, nrow, ncol, k, maxk, iupdate,&
& imode, err, liwork, i
real(kind(1.d0)),pointer :: iwork(:), xp(:,:)
logical :: ichoice
logical :: do_update
character(len=20) :: name, ch_err
name='psb_dovrlm'
@@ -111,17 +109,13 @@ subroutine psb_dovrlm(x,desc_a,info,jx,ik,work,choice,update_type)
k = maxk
end if
if (present(choice)) then
ichoice = choice
if (present(update)) then
iupdate = update
else
ichoice = .true.
endif
if (present(update_type)) then
iupdate = update_type
else
iupdate = psb_none_
iupdate = psb_avg_
endif
do_update = (iupdate /= psb_none_)
imode = IOR(psb_swap_send_,psb_swap_recv_)
! check vector correctness
@@ -166,7 +160,7 @@ subroutine psb_dovrlm(x,desc_a,info,jx,ik,work,choice,update_type)
end if
! exchange overlap elements
if(ichoice) then
if(do_update) then
xp => x(iix:size(x,1),jjx:jjx+k-1)
call psi_swapdata(imode,k,1.d0,xp,&
& desc_a,iwork,info,data=psb_comm_ovr_)
@@ -263,10 +257,9 @@ end subroutine psb_dovrlm
! desc_a - type(<psb_desc_type>). The communication descriptor.
! info - integer. Eventually returns an error code.
! work - real(optional). A working area.
! choice - logical(optional). ???.
! update_type - integer(optional). ???.
! update - integer(optional). ???.
!
subroutine psb_dovrlv(x,desc_a,info,work,choice,update_type)
subroutine psb_dovrlv(x,desc_a,info,work,update)
use psb_descriptor_type
use psi_mod
use psb_const_mod
@@ -279,15 +272,14 @@ subroutine psb_dovrlv(x,desc_a,info,work,choice,update_type)
type(psb_desc_type), intent(in) :: desc_a
integer, intent(out) :: info
real(kind(1.d0)), optional, target :: work(:)
logical, intent(in), optional :: choice
integer, intent(in), optional :: update_type
integer, intent(in), optional :: update
! locals
integer :: int_err(5), icontxt, nprow, npcol, myrow, mycol,&
& err_act, m, n, iix, jjx, temp(2), ix, ijx, nrow, ncol, k, maxk, iupdate,&
& imode, err, liwork, i
real(kind(1.d0)),pointer :: iwork(:)
logical :: ichoice
logical :: do_update
character(len=20) :: name, ch_err
name='psb_dovrlv'
@@ -320,17 +312,13 @@ subroutine psb_dovrlv(x,desc_a,info,work,choice,update_type)
k = 1
if (present(choice)) then
ichoice = choice
if (present(update)) then
iupdate = update
else
ichoice = .true.
endif
if (present(update_type)) then
iupdate = update_type
else
iupdate = psb_none_
iupdate = psb_avg_
endif
do_update = (iupdate /= psb_none_)
imode = IOR(psb_swap_send_,psb_swap_recv_)
! check vector correctness
@@ -375,7 +363,7 @@ subroutine psb_dovrlv(x,desc_a,info,work,choice,update_type)
end if
! exchange overlap elements
if(ichoice) then
if(do_update) then
call psi_swapdata(imode,1.d0,x(iix:size(x)),&
& desc_a,iwork,info,data=psb_comm_ovr_)
end if
src/comm/psb_zovrl.f90
+17 -29
View File
@@ -40,10 +40,9 @@
! jx - integer(optional). The starting column of the global matrix.
! ik - integer(optional). The number of columns to gather.
! work - real(optional). A working area.
! choice - logical(optional). ???.
! update_type - integer(optional). ???.
! update - integer(optional). ???.
!
subroutine psb_zovrlm(x,desc_a,info,jx,ik,work,choice,update_type)
subroutine psb_zovrlm(x,desc_a,info,jx,ik,work,update)
use psb_descriptor_type
use psb_const_mod
use psi_mod
@@ -56,15 +55,14 @@ subroutine psb_zovrlm(x,desc_a,info,jx,ik,work,choice,update_type)
type(psb_desc_type), intent(in) :: desc_a
integer, intent(out) :: info
complex(kind(1.d0)), optional, target :: work(:)
logical, intent(in), optional :: choice
integer, intent(in), optional :: update_type,jx,ik
integer, intent(in), optional :: update,jx,ik
! locals
integer :: int_err(5), icontxt, nprow, npcol, myrow, mycol,&
& err_act, m, n, iix, jjx, temp(2), ix, ijx, nrow, ncol, k, maxk, iupdate,&
& imode, err, liwork, i
complex(kind(1.d0)),pointer :: iwork(:), xp(:,:)
logical :: ichoice
logical :: do_update
character(len=20) :: name, ch_err
name='psb_zovrlm'
@@ -111,17 +109,13 @@ subroutine psb_zovrlm(x,desc_a,info,jx,ik,work,choice,update_type)
k = maxk
end if
if (present(choice)) then
ichoice = choice
if (present(update)) then
iupdate = update
else
ichoice = .true.
endif
if (present(update_type)) then
iupdate = update_type
else
iupdate = psb_none_
iupdate = psb_avg_
endif
do_update = (iupdate /= psb_none_)
imode = IOR(psb_swap_send_,psb_swap_recv_)
! check vector correctness
@@ -166,7 +160,7 @@ subroutine psb_zovrlm(x,desc_a,info,jx,ik,work,choice,update_type)
end if
! exchange overlap elements
if(ichoice) then
if(do_update) then
xp => x(iix:size(x,1),jjx:jjx+k-1)
call psi_swapdata(imode,k,zone,xp,&
& desc_a,iwork,info,data=psb_comm_ovr_)
@@ -263,10 +257,9 @@ end subroutine psb_zovrlm
! desc_a - type(<psb_desc_type>). The communication descriptor.
! info - integer. Eventually returns an error code.
! work - real(optional). A working area.
! choice - logical(optional). ???.
! update_type - integer(optional). ???.
! update - integer(optional). ???.
!
subroutine psb_zovrlv(x,desc_a,info,work,choice,update_type)
subroutine psb_zovrlv(x,desc_a,info,work,update)
use psb_descriptor_type
use psi_mod
use psb_const_mod
@@ -279,15 +272,14 @@ subroutine psb_zovrlv(x,desc_a,info,work,choice,update_type)
type(psb_desc_type), intent(in) :: desc_a
integer, intent(out) :: info
complex(kind(1.d0)), optional, target :: work(:)
logical, intent(in), optional :: choice
integer, intent(in), optional :: update_type
integer, intent(in), optional :: update
! locals
integer :: int_err(5), icontxt, nprow, npcol, myrow, mycol,&
& err_act, m, n, iix, jjx, temp(2), ix, ijx, nrow, ncol, k, maxk, iupdate,&
& imode, err, liwork, i
complex(kind(1.d0)),pointer :: iwork(:)
logical :: ichoice
logical :: do_update
character(len=20) :: name, ch_err
name='psb_zovrlv'
@@ -320,17 +312,13 @@ subroutine psb_zovrlv(x,desc_a,info,work,choice,update_type)
k = 1
if (present(choice)) then
ichoice = choice
else
ichoice = .true.
endif
if (present(update_type)) then
iupdate = update_type
if (present(update)) then
iupdate = update
else
iupdate = psb_none_
endif
do_update = (iupdate /= psb_none_)
imode = IOR(psb_swap_send_,psb_swap_recv_)
! check vector correctness
@@ -376,7 +364,7 @@ subroutine psb_zovrlv(x,desc_a,info,work,choice,update_type)
! exchange overlap elements
if(ichoice) then
if(do_update) then
call psi_swapdata(imode,zone,x(iix:size(x)),&
& desc_a,iwork,info,data=psb_comm_ovr_)
end if
src/methd/psb_dbicg.f90
+19 -19
View File
@@ -102,7 +102,7 @@ subroutine psb_dbicg(a,prec,b,x,eps,desc_a,info,&
& r(:), p(:), zt(:), pt(:), z(:), rt(:),qt(:)
integer, pointer :: iperm(:), ipnull(:), ipsave(:), int_err(:)
real(kind(1.d0)) ::rerr
integer ::litmax, liter, naux, m, mglob, it, itrac,&
integer ::litmax, liter, naux, m, mglob, it, itrace_,&
& nprows,npcols,me,mecol, n_row, n_col, istop_, err_act
character ::diagl, diagu
logical, parameter :: debug = .false.
@@ -159,7 +159,7 @@ subroutine psb_dbicg(a,prec,b,x,eps,desc_a,info,&
naux=4*n_col
allocate(aux(naux),stat=info)
call psb_geall(mglob,9,wwrk,desc_a,info)
call psb_geall(wwrk,desc_a,info,n=9)
call psb_geasb(wwrk,desc_a,info)
if(info.ne.0) then
info=4011
@@ -186,9 +186,9 @@ subroutine psb_dbicg(a,prec,b,x,eps,desc_a,info,&
endif
if (present(itrace)) then
itrac = itrace
itrace_ = itrace
else
itrac = -1
itrace_ = 0
end if
diagl = 'u'
@@ -241,15 +241,8 @@ subroutine psb_dbicg(a,prec,b,x,eps,desc_a,info,&
if (istop_ == 1) then
xni = psb_geamax(x,desc_a,info)
rerr = rni/(ani*xni+bni)
if (itrac /= -1) then
if (me.eq.0) write(itrac,'(a,i4,5(2x,es10.4))') 'bicg: ',itx,rerr,rni,bni,&
&xni,ani
endif
else if (istop_ == 2) then
rerr = rni/bn2
if (itrac /= -1) then
if (me.eq.0) write(itrac,'(a,i4,3(2x,es10.4))') 'bicg: ',itx,rerr,rni,bn2
endif
endif
if(info.ne.0) then
@@ -261,6 +254,11 @@ subroutine psb_dbicg(a,prec,b,x,eps,desc_a,info,&
if (rerr<=eps) then
exit restart
end if
If (itrace_ > 0) then
if ((mod(itx,itrace_)==0).and.(me == 0))&
& write(*,'(a,i4,3(2x,es10.4))') 'bicg: ',itx,rerr
end If
iteration: do
it = it + 1
@@ -315,23 +313,25 @@ subroutine psb_dbicg(a,prec,b,x,eps,desc_a,info,&
if (istop_ == 1) then
xni = psb_geamax(x,desc_a,info)
rerr = rni/(ani*xni+bni)
if (itrac /= -1) then
if (me.eq.0) write(itrac,'(a,i4,5(2x,es10.4))') 'bicg: ',itx,rerr,rni,bni,&
&xni,ani
endif
else if (istop_ == 2) then
rerr = rni/bn2
if (itrac /= -1) then
if (me.eq.0) write(itrac,'(a,i4,3(2x,es10.4))') 'bicg: ',itx,rerr,rni,bn2
endif
endif
if (rerr<=eps) then
exit restart
end if
if (itx.ge.litmax) exit restart
If (itrace_ > 0) then
if ((mod(itx,itrace_)==0).and.(me == 0))&
& write(*,'(a,i4,3(2x,es10.4))') 'bicg: ',itx,rerr
end If
end do iteration
end do restart
If (itrace_ > 0) then
if (me == 0) write(*,'(a,i4,3(2x,es10.4))') 'bicg: ',itx,rerr
end If
if (present(err)) err=rerr
if (present(iter)) iter = itx
if (rerr>eps) then
src/methd/psb_dcg.f90
+13 -12
View File
@@ -102,7 +102,7 @@ Subroutine psb_dcg(a,prec,b,x,eps,desc_a,info,&
real(kind(1.d0)) ::rerr
real(kind(1.d0)) ::alpha, beta, rho, rho_old, rni, xni, bni, ani,bn2,&
& sigma
integer :: litmax, liter, istop_, naux, m, mglob, it, itx, itrac,&
integer :: litmax, liter, istop_, naux, m, mglob, it, itx, itrace_,&
& nprows,npcols,me,mecol, n_col, isvch, ich, icontxt, n_row,err_act, int_err(5)
character ::diagl, diagu
logical, parameter :: exchange=.true., noexchange=.false.
@@ -148,7 +148,7 @@ Subroutine psb_dcg(a,prec,b,x,eps,desc_a,info,&
naux=4*n_col
allocate(aux(naux), stat=info)
call psb_geall(mglob,5,wwrk,desc_a,info)
call psb_geall(wwrk,desc_a,info,n=5)
call psb_geasb(wwrk,desc_a,info)
if (info.ne.0) then
info=4011
@@ -170,9 +170,9 @@ Subroutine psb_dcg(a,prec,b,x,eps,desc_a,info,&
endif
if (present(itrace)) then
itrac = itrace
itrace_ = itrace
else
itrac = -1
itrace_ = 0
end if
itx=0
@@ -243,23 +243,24 @@ Subroutine psb_dcg(a,prec,b,x,eps,desc_a,info,&
rni = psb_geamax(r,desc_a,info)
xni = psb_geamax(x,desc_a,info)
rerr = rni/(ani*xni+bni)
If (itrac /= -1) Then
If (me.Eq.0) Write(itrac,'(a,i4,5(2x,es10.4))') 'cg: ',itx,rerr,rni,bni,&
&xni,ani
Endif
Else If (istop_ == 2) Then
rni = psb_genrm2(r,desc_a,info)
rerr = rni/bn2
If (itrac /= -1) Then
If (me.Eq.0) Write(itrac,'(a,i4,3(2x,es10.4)))') 'cg: ',itx,rerr,rni,bn2
Endif
Endif
if (rerr<=eps) exit restart
if (itx>= litmax) exit restart
If (itrace_ > 0) then
if ((mod(itx,itrace_)==0).and.(me == 0))&
& write(*,'(a,i4,3(2x,es10.4))') 'cg: ',itx,rerr
end If
end do iteration
end do restart
If (itrace_ > 0) then
if (me == 0) write(*,'(a,i4,3(2x,es10.4))') 'cg: ',itx,rerr
end If
if (present(err)) err=rerr
if (present(iter)) iter = itx
src/methd/psb_dcgs.f90
+23 -23
View File
@@ -100,7 +100,7 @@ Subroutine psb_dcgs(a,prec,b,x,eps,desc_a,info,&
& r(:), p(:), v(:), s(:), t(:), z(:), f(:), rt(:),qt(:),uv(:)
Integer, Pointer :: iperm(:), ipnull(:), ipsave(:)
Real(Kind(1.d0)) ::rerr
Integer ::litmax, liter, naux, m, mglob, it, itrac,int_err(5),&
Integer ::litmax, liter, naux, m, mglob, it, itrace_,int_err(5),&
& nprows,npcols,me,mecol, n_row, n_col,istop_, err_act
Character ::diagl, diagu
Logical, Parameter :: exchange=.True., noexchange=.False.
@@ -153,7 +153,7 @@ Subroutine psb_dcgs(a,prec,b,x,eps,desc_a,info,&
naux=4*n_col
Allocate(aux(naux),stat=info)
Call psb_geall(mglob,11,wwrk,desc_a,info)
Call psb_geall(wwrk,desc_a,info,n=11)
Call psb_geasb(wwrk,desc_a,info)
if (info.ne.0) Then
info=4011
@@ -181,9 +181,9 @@ Subroutine psb_dcgs(a,prec,b,x,eps,desc_a,info,&
Endif
If (Present(itrace)) Then
itrac = itrace
itrace_ = itrace
Else
itrac = -1
itrace_ = 0
End If
! Ensure global coherence for convergence checks.
@@ -227,26 +227,23 @@ Subroutine psb_dcgs(a,prec,b,x,eps,desc_a,info,&
rni = psb_geamax(r,desc_a,info)
xni = psb_geamax(x,desc_a,info)
rerr = rni/(ani*xni+bni)
if (itrac /= -1) then
If (me == 0) Write(itrac,'(a,i4,5(2x,es10.4))') 'cgs: ',&
& itx,rerr,rni,bni,xni,ani
endif
else if (istop_ == 2) then
rni = psb_genrm2(r,desc_a,info)
rerr = rni/bn2
if (itrac /= -1) then
If (me == 0) Write(itrac,'(a,i4,3(2x,es10.4))') 'cgs: ',itx,rerr,rni,bn2
endif
endif
if(info/=0)then
info=4011
call psb_errpush(info,name)
goto 9999
end if
end if
If (rerr<=eps) Then
Exit restart
End If
If (itrace_ > 0) then
if ((mod(itx,itrace_)==0).and.(me == 0))&
& write(*,'(a,i4,3(2x,es10.4))') 'cgs: ',itx,rerr
end If
iteration: Do
it = it + 1
@@ -303,28 +300,31 @@ Subroutine psb_dcgs(a,prec,b,x,eps,desc_a,info,&
rni = psb_geamax(r,desc_a,info)
xni = psb_geamax(x,desc_a,info)
rerr = rni/(ani*xni+bni)
if (itrac /= -1) then
If (me == 0) Write(itrac,'(a,i4,5(2x,es10.4))') 'cgs: ',&
& itx,rerr,rni,bni,xni,ani
endif
else if (istop_ == 2) then
rni = psb_genrm2(r,desc_a,info)
rerr = rni/bn2
if (itrac /= -1) then
If (me == 0) Write(itrac,'(a,i4,3(2x,es10.4))') 'cgs: ',&
& itx,rerr,rni,bn2
endif
endif
If (rerr<=eps) Then
Exit restart
End If
If (itx.Ge.litmax) Exit restart
If (itrace_ > 0) then
if ((mod(itx,itrace_)==0).and.(me == 0))&
& write(*,'(a,i4,3(2x,es10.4))') 'cgs: ',itx,rerr
end If
End Do iteration
End Do restart
If (itrace_ > 0) then
if ((mod(itx,itrace_)==0).and.(me == 0))&
& write(*,'(a,i4,3(2x,es10.4))') 'cgs: ',itx,rerr
end If
If (Present(err)) err=rerr
If (Present(iter)) iter = itx
If (rerr>eps) Then
src/methd/psb_dcgstab.f90
+24 -25
View File
@@ -100,7 +100,7 @@ Subroutine psb_dcgstab(a,prec,b,x,eps,desc_a,info,&
& r(:), p(:), v(:), s(:), t(:), z(:), f(:)
Integer, Pointer :: iperm(:), ipnull(:), ipsave(:)
Real(Kind(1.d0)) ::rerr
Integer ::litmax, liter, naux, m, mglob, it,itrac,&
Integer ::litmax, liter, naux, m, mglob, it,itrace_,&
& nprows,npcols,myrow,mycol, n_row, n_col
Character ::diagl, diagu
Logical, Parameter :: debug = .false.
@@ -158,7 +158,7 @@ Subroutine psb_dcgstab(a,prec,b,x,eps,desc_a,info,&
naux=6*n_col
allocate(aux(naux),stat=info)
if (info == 0) call psb_geall(mglob,8,wwrk,desc_a,info)
if (info == 0) call psb_geall(wwrk,desc_a,info,n=8)
if (info == 0) call psb_geasb(wwrk,desc_a,info)
if (info /= 0) then
info=4011
@@ -182,9 +182,9 @@ Subroutine psb_dcgstab(a,prec,b,x,eps,desc_a,info,&
Endif
If (Present(itrace)) Then
itrac = itrace
itrace_ = itrace
Else
itrac = -1
itrace_ = 0
End If
diagl = 'U'
@@ -245,8 +245,8 @@ Subroutine psb_dcgstab(a,prec,b,x,eps,desc_a,info,&
rn0 = rni
End If
If (rn0 == 0.d0 ) Then
If (itrac /= -1) Then
If (myrow == 0) Write(itrac,*) 'BiCGSTAB: ',itx,rn0
If (itrace_ > 0 ) Then
If (myrow == 0) Write(*,*) 'BiCGSTAB: ',itx,rn0
Endif
Exit restart
End If
@@ -254,15 +254,8 @@ Subroutine psb_dcgstab(a,prec,b,x,eps,desc_a,info,&
If (istop_ == 1) Then
xni = psb_geamax(x,desc_a,info)
rerr = rni/(ani*xni+bni)
If (itrac /= -1) Then
If (myrow == 0) Write(itrac,'(a,i4,5(2x,es10.4))') 'bicgstab: ',itx,rerr,rni,bni,&
&xni,ani
Endif
Else If (istop_ == 2) Then
rerr = rni/bn2
If (itrac /= -1) Then
If (myrow == 0) Write(itrac,'(a,i4,3(2x,es10.4))') 'bicgstab: ',itx,rerr,rni,bn2
Endif
Endif
if (info /= 0) Then
info=4011
@@ -270,10 +263,14 @@ Subroutine psb_dcgstab(a,prec,b,x,eps,desc_a,info,&
goto 9999
End If
If (rerr<=eps) Then
Exit restart
End If
If (itrace_ > 0) then
if (((itx==0).or.(mod(itx,itrace_)==0)).and.(myrow == 0)) &
& write(*,'(a,i4,3(2x,es10.4))') 'bicgstab: ',itx,rerr
Endif
iteration: Do
it = it + 1
@@ -357,33 +354,35 @@ Subroutine psb_dcgstab(a,prec,b,x,eps,desc_a,info,&
rni = psb_geamax(r,desc_a,info)
xni = psb_geamax(x,desc_a,info)
rerr = rni/(ani*xni+bni)
If (itrac /= -1) Then
If (myrow == 0) Write(itrac,'(a,i4,5(2x,es10.4))') &
& 'bicgstab: ',itx,rerr,rni,bni,xni,ani
Endif
Else If (istop_ == 2) Then
rni = psb_genrm2(r,desc_a,info)
rerr = rni/bn2
If (itrac /= -1) Then
If (myrow == 0) Write(itrac,'(a,i4,3(2x,es10.4)))') &
& 'bicgstab: ',itx,rerr,rni,bn2
Endif
Endif
If (rerr<=eps) Then
Exit restart
End If
If (itx.Ge.litmax) Exit restart
If (itrace_ > 0) then
if ((mod(itx,itrace_)==0).and.(myrow == 0)) &
& write(*,'(a,i4,3(2x,es10.4)))') &
& 'bicgstab: ',itx,rerr
Endif
End Do iteration
End Do restart
If (itrace_ > 0) then
if (myrow == 0) write(*,'(a,i4,3(2x,es10.4))') 'bicgstab: ',itx,rerr
Endif
If (Present(err)) err=rerr
If (Present(iter)) iter = itx
If (rerr>eps) Then
Write(0,*) 'BI-CGSTAB FAILED TO CONVERGE TO ',EPS,&
& ' IN ',ITX,' ITERATIONS '
Write(0,*) 'BI-CGSTAB failed to converge to ',EPS,&
& ' in ',ITX,' iterations. '
End If
Deallocate(aux)
src/methd/psb_dcgstabl.f90
+19 -25
View File
@@ -109,7 +109,7 @@ Subroutine psb_dcgstabl(a,prec,b,x,eps,desc_a,info,&
&pv1(:), pv2(:), pm1(:,:), pm2(:,:)
Integer, Pointer :: iperm(:), ipnull(:), ipsave(:)
Real(Kind(1.d0)) ::rerr
Integer ::litmax, liter, naux, m, mglob, it, itrac,&
Integer ::litmax, liter, naux, m, mglob, it, itrace_,&
& nprows,npcols,me,mecol, n_row, n_col, nl, err_act
Character ::diagl, diagu
Logical, Parameter :: exchange=.True., noexchange=.False.
@@ -161,9 +161,9 @@ Subroutine psb_dcgstabl(a,prec,b,x,eps,desc_a,info,&
Endif
If (Present(itrace)) Then
itrac = itrace
itrace_ = itrace
Else
itrac = -1
itrace_ = 0
End If
If (Present(irst)) Then
@@ -183,9 +183,9 @@ Subroutine psb_dcgstabl(a,prec,b,x,eps,desc_a,info,&
call psb_errpush(info,name)
goto 9999
End If
Call psb_geall(mglob,10,wwrk,desc_a,info)
Call psb_geall(mglob,nl+1,uh,desc_a,info,js=0)
Call psb_geall(mglob,nl+1,rh,desc_a,info,js=0)
Call psb_geall(wwrk,desc_a,info,n=10)
Call psb_geall(uh,desc_a,info,n=nl+1)
Call psb_geall(rh,desc_a,info,n=nl+1)
Call psb_geasb(wwrk,desc_a,info)
Call psb_geasb(uh,desc_a,info)
Call psb_geasb(rh,desc_a,info)
@@ -255,17 +255,9 @@ Subroutine psb_dcgstabl(a,prec,b,x,eps,desc_a,info,&
rni = psb_geamax(r,desc_a,info)
xni = psb_geamax(x,desc_a,info)
rerr = rni/(ani*xni+bni)
if (itrac /= -1) then
If (me == 0) Write(itrac,'(a,i4,5(2x,es10.4))') 'bicgstab(l): ',&
& itx,rerr,rni,bni,xni,ani
endif
else if (istop_ == 2) then
rni = psb_genrm2(r,desc_a,info)
rerr = rni/bn2
if (itrac /= -1) then
If (me == 0) Write(itrac,'(a,i4,3(2x,es10.4))') 'bicgstab(l): ',&
& itx,rerr,rni,bn2
endif
endif
if (info.ne.0) Then
info=4011
@@ -276,6 +268,10 @@ Subroutine psb_dcgstabl(a,prec,b,x,eps,desc_a,info,&
If (rerr<=eps) Then
Exit restart
End If
If (itrace_ > 0) then
if ((mod(itx,itrace_)==0).and.(me == 0))&
& write(*,'(a,i4,3(2x,es10.4))') 'bicgstab(l): ',itx,rerr
end If
iteration: Do
it = it + nl
@@ -364,29 +360,27 @@ Subroutine psb_dcgstabl(a,prec,b,x,eps,desc_a,info,&
rni = psb_geamax(rh(:,0),desc_a,info)
xni = psb_geamax(x,desc_a,info)
rerr = rni/(ani*xni+bni)
if (itrac /= -1) then
If (me == 0) Write(itrac,'(a,i4,5(2x,es10.4))') 'bicgstab(l): ',&
& itx,rerr,rni,bni,xni,ani
endif
else if (istop_ == 2) then
rni = psb_genrm2(rh(:,0),desc_a,info)
rerr = rni/bn2
if (itrac /= -1) then
If (me == 0) Write(itrac,'(a,i4,3(2x,es10.4))') 'bicgstab(l): ',&
& itx,rerr,rni,bn2
endif
endif
If (rerr<=eps) Then
Exit restart
End If
If (itx.Ge.litmax) Exit restart
If (itrace_ > 0) then
if ((mod(itx,itrace_)==0).and.(me == 0))&
& write(*,'(a,i4,3(2x,es10.4))') 'bicgstab(l): ',itx,rerr
end If
End Do iteration
End Do restart
If (itrace_ > 0) then
if (me == 0) write(*,'(a,i4,3(2x,es10.4))') 'bicgstab(l): ',itx,rerr
end If
If (Present(err)) err=rerr
If (Present(iter)) iter = itx
If (rerr>eps) Then
src/methd/psb_dgmresr.f90
+16 -21
View File
@@ -111,7 +111,7 @@ Subroutine psb_dgmresr(a,prec,b,x,eps,desc_a,info,&
&pv1(:), pv2(:), pm1(:,:), rr(:,:)
Integer, Pointer :: iperm(:), ipnull(:), ipsave(:), ierrv(:)
Real(Kind(1.d0)) :: rerr, scal, gm
Integer ::litmax, liter, naux, m, mglob, it,k, itrac,&
Integer ::litmax, liter, naux, m, mglob, it,k, itrace_,&
& nprows,npcols,me,mecol, n_row, n_col, nl, int_err(5)
Character ::diagl, diagu
Logical, Parameter :: exchange=.True., noexchange=.False.
@@ -164,9 +164,9 @@ Subroutine psb_dgmresr(a,prec,b,x,eps,desc_a,info,&
Endif
If (Present(itrace)) Then
itrac = itrace
itrace_ = itrace
Else
itrac = -1
itrace_ = 0
End If
If (Present(irst)) Then
@@ -188,8 +188,8 @@ Subroutine psb_dgmresr(a,prec,b,x,eps,desc_a,info,&
goto 9999
End If
Call psb_geall(mglob,nl+1,v,desc_a,info)
Call psb_geall(mglob,w,desc_a,info)
Call psb_geall(v,desc_a,info,n=nl+1)
Call psb_geall(w,desc_a,info)
Call psb_geasb(v,desc_a,info)
Call psb_geasb(w,desc_a,info)
if (info.ne.0) Then
@@ -247,17 +247,9 @@ Subroutine psb_dgmresr(a,prec,b,x,eps,desc_a,info,&
rni = psb_geamax(v(:,1),desc_a,info)
xni = psb_geamax(x,desc_a,info)
rerr = rni/(ani*xni+bni)
if (itrac /= -1) then
If (me == 0) Write(itrac,'(a,i4,5(2x,es10.4))') 'gmresr(l): ',&
& itx,rerr,rni,bni,xni,ani
endif
else if (istop_ == 2) then
rni = psb_genrm2(v(:,1),desc_a,info)
rerr = rni/bn2
if (itrac /= -1) then
If (me == 0) Write(itrac,'(a,i4,3(2x,es10.4))') 'gmresr(l): ',&
& itx,rerr,rni,bn2
endif
endif
if (info.ne.0) Then
info=4011
@@ -268,6 +260,10 @@ Subroutine psb_dgmresr(a,prec,b,x,eps,desc_a,info,&
If (rerr<=eps) Then
Exit restart
End If
If (itrace_ > 0) then
if ((mod(itx,itrace_)==0).and.(me == 0))&
& write(*,'(a,i4,3(2x,es10.4))') 'gmres(l): ',itx,rerr
end If
If (itx.Ge.litmax) Exit restart
@@ -304,17 +300,9 @@ Subroutine psb_dgmresr(a,prec,b,x,eps,desc_a,info,&
rni = abs(rs(i+1))
xni = psb_geamax(x,desc_a,info)
rerr = rni/(ani*xni+bni)
if (itrac /= -1) then
If (me == 0) Write(itrac,'(a,i4,5(2x,es10.4))') 'gmresr(l): ',&
& itx,rerr,rni,bni,xni,ani
endif
else if (istop_ == 2) then
rni = abs(rs(i+1))
rerr = rni/bn2
if (itrac /= -1) then
If (me == 0) Write(itrac,'(a,i4,3(2x,es10.4))') 'gmresr(l): ',&
& itx,rerr,rni,bn2
endif
endif
if (rerr < eps ) then
@@ -325,6 +313,10 @@ Subroutine psb_dgmresr(a,prec,b,x,eps,desc_a,info,&
end do
exit restart
end if
If (itrace_ > 0) then
if ((mod(itx,itrace_)==0).and.(me == 0))&
& write(*,'(a,i4,3(2x,es10.4))') 'gmres(l): ',itx,rerr
end If
end Do inner
if (debug) write(0,*) 'Before DTRSM :',rs(1:nl)
@@ -335,6 +327,9 @@ Subroutine psb_dgmresr(a,prec,b,x,eps,desc_a,info,&
end do
End Do restart
If (itrace_ > 0) then
if (me == 0) write(*,'(a,i4,3(2x,es10.4))') 'gmres(l): ',itx,rerr
end If
If (Present(err)) err=rerr
If (Present(iter)) iter = itx
src/methd/psb_zcgs.f90
+17 -21
View File
@@ -100,7 +100,7 @@ Subroutine psb_zcgs(a,prec,b,x,eps,desc_a,info,&
& r(:), p(:), v(:), s(:), t(:), z(:), f(:), rt(:),qt(:),uv(:)
Integer, Pointer :: iperm(:), ipnull(:), ipsave(:)
Real(Kind(1.d0)) ::rerr
Integer ::litmax, liter, naux, m, mglob, it, itrac,int_err(5),&
Integer ::litmax, liter, naux, m, mglob, it, itrace_,int_err(5),&
& nprows,npcols,me,mecol, n_row, n_col,istop_, err_act
Character ::diagl, diagu
Logical, Parameter :: exchange=.True., noexchange=.False.
@@ -153,7 +153,7 @@ Subroutine psb_zcgs(a,prec,b,x,eps,desc_a,info,&
naux=4*n_col
Allocate(aux(naux),stat=info)
Call psb_geall(mglob,11,wwrk,desc_a,info)
Call psb_geall(wwrk,desc_a,info,n=11)
Call psb_geasb(wwrk,desc_a,info)
if (info.ne.0) Then
info=4011
@@ -181,9 +181,9 @@ Subroutine psb_zcgs(a,prec,b,x,eps,desc_a,info,&
Endif
If (Present(itrace)) Then
itrac = itrace
itrace_ = itrace
Else
itrac = -1
itrace_ = 0
End If
! Ensure global coherence for convergence checks.
@@ -227,16 +227,9 @@ Subroutine psb_zcgs(a,prec,b,x,eps,desc_a,info,&
rni = psb_geamax(r,desc_a,info)
xni = psb_geamax(x,desc_a,info)
rerr = rni/(ani*xni+bni)
if (itrac /= -1) then
If (me == 0) Write(itrac,'(a,i4,5(2x,es10.4))') 'cgs: ',&
& itx,rerr,rni,bni,xni,ani
endif
else if (istop_ == 2) then
rni = psb_genrm2(r,desc_a,info)
rerr = rni/bn2
if (itrac /= -1) then
If (me == 0) Write(itrac,'(a,i4,3(2x,es10.4))') 'cgs: ',itx,rerr,rni,bn2
endif
endif
if(info/=0)then
info=4011
@@ -247,6 +240,10 @@ Subroutine psb_zcgs(a,prec,b,x,eps,desc_a,info,&
If (rerr<=eps) Then
Exit restart
End If
If (itrace_ > 0) then
if ((mod(itx,itrace_)==0).and.(me == 0))&
& write(*,'(a,i4,3(2x,es10.4))') 'cgs: ',itx,rerr
end If
iteration: Do
it = it + 1
@@ -303,27 +300,26 @@ Subroutine psb_zcgs(a,prec,b,x,eps,desc_a,info,&
rni = psb_geamax(r,desc_a,info)
xni = psb_geamax(x,desc_a,info)
rerr = rni/(ani*xni+bni)
if (itrac /= -1) then
If (me == 0) Write(itrac,'(a,i4,5(2x,es10.4))') 'cgs: ',&
& itx,rerr,rni,bni,xni,ani
endif
else if (istop_ == 2) then
rni = psb_genrm2(r,desc_a,info)
rerr = rni/bn2
if (itrac /= -1) then
If (me == 0) Write(itrac,'(a,i4,3(2x,es10.4))') 'cgs: ',&
& itx,rerr,rni,bn2
endif
endif
If (rerr<=eps) Then
Exit restart
End If
If (itx.Ge.litmax) Exit restart
If (itrace_ > 0) then
if ((mod(itx,itrace_)==0).and.(me == 0))&
& write(*,'(a,i4,3(2x,es10.4))') 'cgs: ',itx,rerr
end If
End Do iteration
End Do restart
If (itrace_ > 0) then
if (me == 0) write(*,'(a,i4,3(2x,es10.4))') 'cgs: ',itx,rerr
end If
If (Present(err)) err=rerr
If (Present(iter)) iter = itx
src/methd/psb_zcgstab.f90
+19 -26
View File
@@ -100,7 +100,7 @@ Subroutine psb_zcgstab(a,prec,b,x,eps,desc_a,info,&
& r(:), p(:), v(:), s(:), t(:), z(:), f(:)
Integer, Pointer :: iperm(:), ipnull(:), ipsave(:)
Real(Kind(1.d0)) :: rerr
Integer :: litmax, liter, naux, m, mglob, it,itrac,&
Integer :: litmax, liter, naux, m, mglob, it,itrace_,&
& nprows,npcols,myrow,mycol, n_row, n_col
Character ::diagl, diagu
Logical, Parameter :: debug = .false.
@@ -158,7 +158,7 @@ Subroutine psb_zcgstab(a,prec,b,x,eps,desc_a,info,&
naux=6*n_col
allocate(aux(naux),stat=info)
call psb_geall(mglob,8,wwrk,desc_a,info)
call psb_geall(wwrk,desc_a,info,n=8)
call psb_geasb(wwrk,desc_a,info)
if (info /= 0) then
info=4011
@@ -182,9 +182,9 @@ Subroutine psb_zcgstab(a,prec,b,x,eps,desc_a,info,&
Endif
If (Present(itrace)) Then
itrac = itrace
itrace_ = itrace
Else
itrac = -1
itrace_ = 0
End If
diagl = 'U'
@@ -245,24 +245,14 @@ Subroutine psb_zcgstab(a,prec,b,x,eps,desc_a,info,&
rn0 = rni
End If
If (rn0 == 0.d0 ) Then
If (itrac /= -1) Then
If (myrow == 0) Write(itrac,*) 'BiCGSTAB: ',itx,rn0
Endif
Exit restart
End If
If (istop_ == 1) Then
xni = psb_geamax(x,desc_a,info)
rerr = rni/(ani*xni+bni)
If (itrac /= -1) Then
If (myrow == 0) Write(itrac,'(a,i4,5(2x,es10.4))') 'bicgstab: ',itx,rerr,rni,bni,&
&xni,ani
Endif
Else If (istop_ == 2) Then
rerr = rni/bn2
If (itrac /= -1) Then
If (myrow == 0) Write(itrac,'(a,i4,3(2x,es10.4))') 'bicgstab: ',itx,rerr,rni,bn2
Endif
Endif
if (info /= 0) Then
info=4011
@@ -270,11 +260,15 @@ Subroutine psb_zcgstab(a,prec,b,x,eps,desc_a,info,&
goto 9999
End If
If (rerr<=eps) Then
Exit restart
End If
If (itrace_ > 0) then
if ((mod(itx,itrace_)==0).and.(myrow == 0))&
& write(*,'(a,i4,3(2x,es10.4))') 'bicgstab: ',itx,rerr
end If
iteration: Do
it = it + 1
itx = itx + 1
@@ -356,18 +350,9 @@ Subroutine psb_zcgstab(a,prec,b,x,eps,desc_a,info,&
rni = psb_geamax(r,desc_a,info)
xni = psb_geamax(x,desc_a,info)
rerr = rni/(ani*xni+bni)
If (itrac /= -1) Then
If (myrow == 0) Write(itrac,'(a,i4,5(2x,es10.4))') &
& 'bicgstab: ',itx,rerr,rni,bni,xni,ani
Endif
Else If (istop_ == 2) Then
rni = psb_genrm2(r,desc_a,info)
rerr = rni/bn2
If (itrac /= -1) Then
If (myrow == 0) Write(itrac,'(a,i4,3(2x,es10.4)))') &
& 'bicgstab: ',itx,rerr,rni,bn2
Endif
Endif
If (rerr<=eps) Then
@@ -375,14 +360,22 @@ Subroutine psb_zcgstab(a,prec,b,x,eps,desc_a,info,&
End If
If (itx.Ge.litmax) Exit restart
If (itrace_ > 0) then
if ((mod(itx,itrace_)==0).and.(myrow == 0))&
& write(*,'(a,i4,3(2x,es10.4))') 'bicgstab: ',itx,rerr
end If
End Do iteration
End Do restart
If (itrace_ > 0) then
if (myrow == 0) write(*,'(a,i4,3(2x,es10.4))') 'bicgstab: ',itx,rerr
end If
If (Present(err)) err=rerr
If (Present(iter)) iter = itx
If (rerr>eps) Then
Write(0,*) 'BI-CGSTAB FAILED TO CONVERGE TO ',EPS,&
& ' IN ',ITX,' ITERATIONS '
Write(0,*) 'BI-cgstab failed to converge to ',eps,&
& ' in ',itx,' iterations. '
End If
Deallocate(aux)
src/prec/psb_dbaseprc_aply.f90
+2 -2
View File
@@ -234,11 +234,11 @@ subroutine psb_dbaseprc_aply(prec,x,beta,y,desc_data,trans,work,info)
case(psb_none_)
! Would work anyway, but since it's supposed to do nothing...
! call f90_psovrl(ty,prec%desc_data,update_type=prec%a_restrict)
! call f90_psovrl(ty,prec%desc_data,update=prec%a_restrict)
case(psb_sum_,psb_avg_)
call psb_ovrl(ty,prec%desc_data,info,&
& update_type=prec%iprcparm(prol_),work=aux)
& update=prec%iprcparm(prol_),work=aux)
if(info /=0) then
info=4010
ch_err='psb_ovrl'
src/prec/psb_dslu_bld.f90
+1 -1
View File
@@ -89,7 +89,7 @@ subroutine psb_dslu_bld(a,desc_a,p,info)
call blacs_barrier(icontxt,'All')
endif
call psb_dcsdp(a,atmp,info)
call psb_csdp(a,atmp,info)
if(info /= 0) then
info=4010
ch_err='psb_dcsdp'
src/prec/psb_zbaseprc_aply.f90
+2 -2
View File
@@ -234,11 +234,11 @@ subroutine psb_zbaseprc_aply(prec,x,beta,y,desc_data,trans,work,info)
case(psb_none_)
! Would work anyway, but since it's supposed to do nothing...
! call f90_psovrl(ty,prec%desc_data,update_type=prec%a_restrict)
! call f90_psovrl(ty,prec%desc_data,update=prec%a_restrict)
case(psb_sum_,psb_avg_)
call psb_ovrl(ty,prec%desc_data,info,&
& update_type=prec%iprcparm(prol_),work=aux)
& update=prec%iprcparm(prol_),work=aux)
if(info /=0) then
info=4010
ch_err='psb_ovrl'
src/prec/psb_zbldaggrmat.f90
+2 -2
View File
@@ -624,7 +624,7 @@ contains
if (debug) write(0,*) me,'Done NUMBMM 2'
if (p%iprcparm(smth_kind_) == smth_omg_) then
call psb_transp(am1,am2,fmt='COO')
call psb_transc(am1,am2,fmt='COO')
nzl = am2%infoa(psb_nnz_)
i=0
!
@@ -643,7 +643,7 @@ contains
am2%infoa(psb_nnz_) = i
call psb_ipcoo2csr(am2,info)
else
call psb_transp(am1,am2)
call psb_transc(am1,am2)
endif
if (debug) write(0,*) me,'starting sphalo/ rwxtd'
src/prec/psb_zslu_bld.f90
+4 -4
View File
@@ -51,7 +51,7 @@ subroutine psb_zslu_bld(a,desc_a,p,info)
character(len=5) :: fmt
character :: upd='F'
integer :: i,j,nza,nzb,nzt,icontxt, me,mycol,nprow,npcol,err_act
logical, parameter :: debug=.false.
logical, parameter :: debug=.true.
character(len=20) :: name, ch_err
interface psb_asmatbld
@@ -85,14 +85,14 @@ subroutine psb_zslu_bld(a,desc_a,p,info)
atmp%fida='COO'
if (Debug) then
write(0,*) me, 'SPLUBLD: Calling csdp'
write(0,*) me, 'ZSLUBLD: Calling csdp'
call blacs_barrier(icontxt,'All')
endif
call psb_zcsdp(a,atmp,info)
call psb_csdp(a,atmp,info)
if(info /= 0) then
info=4010
ch_err='psb_zcsdp'
ch_err='psb_csdp'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
src/psblas/psb_dspmm.f90
+12 -12
View File
@@ -188,26 +188,26 @@ subroutine psb_dspmm(alpha,a,x,beta,y,desc_a,info,&
liwork= 2*ncol
if (a%pr(1) /= 0) liwork = liwork + n * ik
if (a%pl(1) /= 0) liwork = liwork + m * ik
if (present(work)) then
if (present(work)) then
if (size(work) >= liwork) then
aliw =.false.
aliw =.false.
else
aliw=.true.
aliw=.true.
endif
else
aliw=.true.
aliw=.true.
end if
if (aliw) then
call psb_realloc(liwork,iwork,info)
if(info.ne.0) then
info=4010
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call psb_realloc(liwork,iwork,info)
if(info.ne.0) then
info=4010
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
else
iwork => work
iwork => work
endif
iwork(1)=dzero
src/psblas/psb_dspsm.f90
+26 -26
View File
@@ -196,26 +196,26 @@ subroutine psb_dspsm(alpha,a,x,beta,y,desc_a,info,&
liwork= 2*ncol
if (a%pr(1) /= 0) llwork = liwork + m * ik
if (a%pl(1) /= 0) llwork = llwork + m * ik
if (present(work)) then
if (present(work)) then
if (size(work) >= liwork) then
aliw =.false.
aliw =.false.
else
aliw=.true.
aliw=.true.
endif
else
aliw=.true.
aliw=.true.
end if
if (aliw) then
call psb_realloc(liwork,iwork,info)
if(info.ne.0) then
info=4010
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call psb_realloc(liwork,iwork,info)
if(info.ne.0) then
info=4010
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
else
iwork => work
iwork => work
endif
iwork(1)=0.d0
@@ -490,27 +490,27 @@ subroutine psb_dspsv(alpha,a,x,beta,y,desc_a,info,&
liwork= 2*ncol
if (a%pr(1) /= 0) llwork = liwork + m * ik
if (a%pl(1) /= 0) llwork = llwork + m * ik
if (present(work)) then
if (size(work) >= liwork) then
aliw =.false.
aliw =.false.
else
aliw=.true.
aliw=.true.
endif
else
aliw=.true.
aliw=.true.
end if
if (aliw) then
call psb_realloc(liwork,iwork,info)
if(info.ne.0) then
info=4010
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call psb_realloc(liwork,iwork,info)
if(info.ne.0) then
info=4010
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
else
iwork => work
iwork => work
endif
iwork(1)=0.d0
@@ -596,7 +596,7 @@ subroutine psb_dspsv(alpha,a,x,beta,y,desc_a,info,&
end select
end if
if (aliw) deallocate(iwork)
if (aliw) deallocate(iwork)
if(.not.present(d)) deallocate(id)
call psb_erractionrestore(err_act)
src/psblas/psb_zspmm.f90
+13 -13
View File
@@ -183,26 +183,26 @@ subroutine psb_zspmm(alpha,a,x,beta,y,desc_a,info,&
liwork= 2*ncol
if (a%pr(1) /= 0) liwork = liwork + n * ik
if (a%pl(1) /= 0) liwork = liwork + m * ik
if (present(work)) then
if (present(work)) then
if (size(work) >= liwork) then
aliw =.false.
aliw =.false.
else
aliw=.true.
aliw=.true.
endif
else
aliw=.true.
aliw=.true.
end if
if (aliw) then
call psb_realloc(liwork,iwork,info)
if(info.ne.0) then
info=4010
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call psb_realloc(liwork,iwork,info)
if(info.ne.0) then
info=4010
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
else
iwork => work
iwork => work
endif
iwork(1)=zzero
@@ -315,7 +315,7 @@ if (present(work)) then
y(iiy+nrow+1-1:iiy+ncol,1:ik)=zzero
! local Matrix-vector product
call psb_csmm(alpha,a,x(iix:lldx,jjx:jjx+ik-1),&
& beta,y(iiy:lldy,jjy:jjy+ik-1),info,trans=itrans)
src/psblas/psb_zspsm.f90
+86 -86
View File
@@ -199,26 +199,26 @@ subroutine psb_zspsm(alpha,a,x,beta,y,desc_a,info,&
liwork= 2*ncol
if (a%pr(1) /= 0) llwork = liwork + m * ik
if (a%pl(1) /= 0) llwork = llwork + m * ik
if (present(work)) then
if (present(work)) then
if (size(work) >= liwork) then
aliw =.false.
aliw =.false.
else
aliw=.true.
aliw=.true.
endif
else
aliw=.true.
aliw=.true.
end if
if (aliw) then
call psb_realloc(liwork,iwork,info)
if(info.ne.0) then
info=4010
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call psb_realloc(liwork,iwork,info)
if(info.ne.0) then
info=4010
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
else
iwork => work
iwork => work
endif
iwork(1)=0.d0
@@ -428,14 +428,14 @@ subroutine psb_zspsv(alpha,a,x,beta,y,desc_a,info,&
! check on blacs grid
call blacs_gridinfo(icontxt, nprow, npcol, myrow, mycol)
if (nprow == -1) then
info = 2010
call psb_errpush(info,name)
goto 9999
info = 2010
call psb_errpush(info,name)
goto 9999
else if (npcol /= 1) then
info = 2030
int_err(1) = npcol
call psb_errpush(info,name)
goto 9999
info = 2030
int_err(1) = npcol
call psb_errpush(info,name)
goto 9999
endif
! just this case right now
@@ -469,7 +469,7 @@ subroutine psb_zspsv(alpha,a,x,beta,y,desc_a,info,&
goto 9999
end if
else
itrans = 'N'
itrans = 'N'
endif
m = desc_a%matrix_data(psb_m_)
@@ -479,9 +479,9 @@ subroutine psb_zspsv(alpha,a,x,beta,y,desc_a,info,&
lldy = size(y)
if((lldx.lt.ncol).or.(lldy.lt.ncol)) then
info=3010
call psb_errpush(info,name)
goto 9999
info=3010
call psb_errpush(info,name)
goto 9999
end if
iwork => null()
@@ -489,38 +489,38 @@ subroutine psb_zspsv(alpha,a,x,beta,y,desc_a,info,&
liwork= 2*ncol
if (a%pr(1) /= 0) llwork = liwork + m * ik
if (a%pl(1) /= 0) llwork = llwork + m * ik
if (present(work)) then
if (size(work) >= liwork) then
aliw =.false.
aliw =.false.
else
aliw=.true.
aliw=.true.
endif
else
aliw=.true.
aliw=.true.
end if
if (aliw) then
call psb_realloc(liwork,iwork,info)
if(info.ne.0) then
info=4010
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call psb_realloc(liwork,iwork,info)
if(info.ne.0) then
info=4010
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
else
iwork => work
iwork => work
endif
iwork(1)=0.d0
if(present(d)) then
lld = size(d)
id => d
lld = size(d)
id => d
else
lld=1
allocate(id(1))
id=1.d0
lld=1
allocate(id(1))
id=1.d0
end if
! checking for matrix correctness
@@ -529,25 +529,25 @@ subroutine psb_zspsv(alpha,a,x,beta,y,desc_a,info,&
call psb_chkvect(m,ik,size(x),ix,jx,desc_a%matrix_data,info,iix,jjx)
call psb_chkvect(m,ik,size(y),iy,jy,desc_a%matrix_data,info,iiy,jjy)
if(info.ne.0) then
info=4010
ch_err='psb_chkvect/mat'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_chkvect/mat'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if(ja.ne.ix) then
! this case is not yet implemented
info = 3030
! this case is not yet implemented
info = 3030
end if
if((iix.ne.1).or.(iiy.ne.1)) then
! this case is not yet implemented
info = 3040
! this case is not yet implemented
info = 3040
end if
if(info.ne.0) then
call psb_errpush(info,name)
goto 9999
call psb_errpush(info,name)
goto 9999
end if
! Perform local triangular system solve
@@ -556,46 +556,46 @@ subroutine psb_zspsv(alpha,a,x,beta,y,desc_a,info,&
call psb_cssm(alpha,a,xp,beta,yp,info,unitd=lunitd,d=id,trans=itrans)
if(info.ne.0) then
info = 4010
ch_err='dcssm'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info = 4010
ch_err='dcssm'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
! update overlap elements
if(lchoice.gt.0) then
call psi_swapdata(ior(psb_swap_send_,psb_swap_recv_),&
& zone,yp,desc_a,iwork,info)
call psi_swapdata(ior(psb_swap_send_,psb_swap_recv_),&
& zone,yp,desc_a,iwork,info)
i=0
! switch on update type
select case (lchoice)
case(psb_square_root_)
do while(desc_a%ovrlap_elem(i).ne.-ione)
y(desc_a%ovrlap_elem(i+psb_ovrlp_elem_)) =&
& y(desc_a%ovrlap_elem(i+psb_ovrlp_elem_))/&
& sqrt(real(desc_a%ovrlap_elem(i+psb_n_dom_ovr_)))
i = i+2
end do
case(psb_avg_)
do while(desc_a%ovrlap_elem(i).ne.-ione)
y(desc_a%ovrlap_elem(i+psb_ovrlp_elem_)) =&
& y(desc_a%ovrlap_elem(i+psb_ovrlp_elem_))/&
& real(desc_a%ovrlap_elem(i+psb_n_dom_ovr_))
i = i+2
end do
case(psb_sum_)
! do nothing
case default
! wrong value for choice argument
info = 70
int_err=(/10,lchoice,0,0,0/)
call psb_errpush(info,name,i_err=int_err)
goto 9999
end select
i=0
! switch on update type
select case (lchoice)
case(psb_square_root_)
do while(desc_a%ovrlap_elem(i).ne.-ione)
y(desc_a%ovrlap_elem(i+psb_ovrlp_elem_)) =&
& y(desc_a%ovrlap_elem(i+psb_ovrlp_elem_))/&
& sqrt(real(desc_a%ovrlap_elem(i+psb_n_dom_ovr_)))
i = i+2
end do
case(psb_avg_)
do while(desc_a%ovrlap_elem(i).ne.-ione)
y(desc_a%ovrlap_elem(i+psb_ovrlp_elem_)) =&
& y(desc_a%ovrlap_elem(i+psb_ovrlp_elem_))/&
& real(desc_a%ovrlap_elem(i+psb_n_dom_ovr_))
i = i+2
end do
case(psb_sum_)
! do nothing
case default
! wrong value for choice argument
info = 70
int_err=(/10,lchoice,0,0,0/)
call psb_errpush(info,name,i_err=int_err)
goto 9999
end select
end if
if (aliw) deallocate(iwork)
if (aliw) deallocate(iwork)
if(.not.present(d)) deallocate(id)
call psb_erractionrestore(err_act)
@@ -605,8 +605,8 @@ subroutine psb_zspsv(alpha,a,x,beta,y,desc_a,info,&
call psb_erractionrestore(err_act)
if (err_act.eq.act_abort) then
call psb_error(icontxt)
return
call psb_error(icontxt)
return
end if
return
end subroutine psb_zspsv
src/serial/psb_dcsdp.f90
+9 -3
View File
@@ -118,8 +118,7 @@ subroutine psb_dcsdp(a, b,info,ifc,check,trans,unitd,upd,dupl)
if (check_=='R') then
allocate(work(max(size(a%aspk),size(b%aspk))+1000),stat=info)
else
allocate(work(max(size(a%ia1),size(b%ia1),&
& size(a%ia2),size(b%ia2))+max(a%m,b%m)+1000),stat=info)
allocate(work(max(size(a%ia1),size(a%ia2))+max(a%m,b%m)+1000),stat=info)
endif
if (info /= 0) then
@@ -196,13 +195,20 @@ subroutine psb_dcsdp(a, b,info,ifc,check,trans,unitd,upd,dupl)
call psb_errpush(info,name,a_err=ch_err)
goto 9999
endif
if ((size(b%aspk) < aspk_size) .or.&
if (.not.associated(b%aspk).or.&
&.not.associated(b%ia1).or.&
&.not.associated(b%ia2).or.&
&.not.associated(b%pl).or.&
&.not.associated(b%pr)) then
call psb_sp_reall(b,ia1_size,ia2_size,aspk_size,info)
else if ((size(b%aspk) < aspk_size) .or.&
&(size(b%ia1) < ia1_size) .or.&
&(size(b%ia2) < ia2_size) .or.&
&(size(b%pl) < b%m) .or.&
&(size(b%pr) < b%k )) then
call psb_sp_reall(b,ia1_size,ia2_size,aspk_size,info)
endif
if (info /= no_err) then
info=4010
ch_err='psb_sp_reall'
src/serial/psb_zcsdp.f90
+15 -6
View File
@@ -67,7 +67,7 @@ subroutine psb_zcsdp(a, b,info,ifc,check,trans,unitd,upd,dupl)
& ipc, i, count, err_act, ierrv(5), i1, i2, ia
character :: check_,trans_,unitd_, up
Integer, Parameter :: maxtry=8
logical, parameter :: debug=.false.
logical, parameter :: debug=.true.
character(len=20) :: name, ch_err
interface psb_cest
@@ -114,11 +114,11 @@ subroutine psb_zcsdp(a, b,info,ifc,check,trans,unitd,upd,dupl)
unitd_ = 'U'
endif
if (check_=='R') then
allocate(work(max(size(a%aspk),size(b%aspk))+1000),stat=info)
else
allocate(work(max(size(a%ia1),size(b%ia1),&
& size(a%ia2),size(b%ia2))+max(a%m,b%m)+1000),stat=info)
allocate(work(max(size(a%ia1),size(a%ia2))+max(a%m,b%m)+1000),stat=info)
endif
if (info /= 0) then
@@ -127,7 +127,7 @@ subroutine psb_zcsdp(a, b,info,ifc,check,trans,unitd,upd,dupl)
goto 9999
end if
if (ifc_<1) then
write(0,*) 'dcsdp90 Error: invalid ifc ',ifc_
write(0,*) 'csdp90 Error: invalid ifc ',ifc_
info = -4
call psb_errpush(info,name)
goto 9999
@@ -181,13 +181,15 @@ subroutine psb_zcsdp(a, b,info,ifc,check,trans,unitd,upd,dupl)
b%m=a%m
b%k=a%k
call psb_spinfo(psb_nztotreq_,a,size_req,info)
if (debug) write(0,*) 'DCSDP : size_req 1:',size_req
if (debug) write(0,*) 'DCSDP : size_req 1:',size_req,a%m,a%k
!
n_row=b%m
n_col=b%k
call psb_cest(b%fida, n_row,n_col,size_req,&
& ia1_size, ia2_size, aspk_size, upd_,info)
!!$ write(0,*) size(b%aspk),size(b%ia1),size(b%ia2),size(b%pl),size(b%pr),&
!!$ & ia1_size, ia2_size, aspk_size,b%fida,b%m,b%k
if (info /= no_err) then
info=4010
@@ -195,13 +197,20 @@ subroutine psb_zcsdp(a, b,info,ifc,check,trans,unitd,upd,dupl)
call psb_errpush(info,name,a_err=ch_err)
goto 9999
endif
if ((size(b%aspk) < aspk_size) .or.&
if (.not.associated(b%aspk).or.&
&.not.associated(b%ia1).or.&
&.not.associated(b%ia2).or.&
&.not.associated(b%pl).or.&
&.not.associated(b%pr)) then
call psb_sp_reall(b,ia1_size,ia2_size,aspk_size,info)
else if ((size(b%aspk) < aspk_size) .or.&
&(size(b%ia1) < ia1_size) .or.&
&(size(b%ia2) < ia2_size) .or.&
&(size(b%pl) < b%m) .or.&
&(size(b%pr) < b%k )) then
call psb_sp_reall(b,ia1_size,ia2_size,aspk_size,info)
endif
if (info /= no_err) then
info=4010
ch_err='psb_sp_reall'
src/serial/psb_ztransc.f90
+92
View File
@@ -0,0 +1,92 @@
!!$
!!$ Parallel Sparse BLAS v2.0
!!$ (C) Copyright 2006 Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
! File: psb_ztransc.f90
! Subroutine:
! Parameters:
subroutine psb_ztransc(a,b,c,fmt)
use psb_spmat_type
use psb_tools_mod
use psb_string_mod
use psb_serial_mod, only : psb_ipcoo2csr, psb_ipcsr2coo, psb_fixcoo
implicit none
type(psb_zspmat_type) :: a,b
integer, optional :: c
character(len=*), optional :: fmt
character(len=5) :: fmt_
integer ::c_, info, nz, i
integer, pointer :: itmp(:)
if (present(c)) then
c_=c
else
c_=1
endif
if (present(fmt)) then
fmt_ = toupper(fmt)
else
fmt_='CSR'
endif
if (associated(b%aspk)) call psb_sp_free(b,info)
call psb_sp_clone(a,b,info)
if (b%fida=='CSR') then
call psb_ipcsr2coo(b,info)
else if (b%fida=='COO') then
! do nothing
else
write(0,*) 'Unimplemented case in TRANSC '
endif
itmp => b%ia1
b%ia1 => b%ia2
b%ia2 => itmp
b%m = a%k
b%k = a%m
do i=1, b%infoa(psb_nnz_)
b%aspk(i) = conjg(b%aspk(i))
end do
!!$ write(0,*) 'Calling IPCOO2CSR from transc90 ',b%m,b%k
if (fmt_=='CSR') then
call psb_ipcoo2csr(b,info)
b%fida='CSR'
else if (fmt_=='COO') then
call psb_fixcoo(b,info)
b%fida='COO'
else
write(0,*) 'Unknown FMT in TRANSC : "',fmt_,'"'
endif
return
end subroutine psb_ztransc
src/serial/psb_ztransp.f90
+1 -1
View File
@@ -45,7 +45,7 @@ subroutine psb_ztransp(a,b,c,fmt)
character(len=5) :: fmt_
integer ::c_, info, nz
integer, pointer :: itmp(:)=>null()
integer, pointer :: itmp(:)
if (present(c)) then
c_=c
else
src/tools/psb_dallc.f90
+28 -91
View File
@@ -35,13 +35,11 @@
! Allocates dense matrix for PSBLAS routines
!
! Parameters:
! m - number of rows.
! n - number of columns.
! x - the matrix to be allocated.
! desc_a - the communication descriptor.
! info - eventually returns an error code
! js - (optional) the starting column
subroutine psb_dalloc(m, n, x, desc_a, info, js)
! info - possibly returns an error code
! n - optional number of columns.
subroutine psb_dalloc(x, desc_a, info, n)
!....allocate dense matrix for psblas routines.....
use psb_descriptor_type
use psb_const_mod
@@ -50,18 +48,17 @@ subroutine psb_dalloc(m, n, x, desc_a, info, js)
implicit none
!....parameters...
integer, intent(in) :: m,n
real(kind(1.d0)), pointer :: x(:,:)
type(psb_desc_type), intent(in) :: desc_a
integer :: info
integer, optional, intent(in) :: js
integer, optional, intent(in) :: n
!locals
integer :: nprow,npcol,myrow,mycol,err,n_col,n_row,i,j,jj,err_act
integer :: icontxt,dectype
integer :: icontxt,dectype,n_
integer :: int_err(5),temp(1),exch(3)
real(kind(1.d0)) :: real_err(5)
integer, allocatable:: prc_v(:)
character(len=20) :: name, ch_err
name='psb_dallc'
@@ -88,88 +85,56 @@ subroutine psb_dalloc(m, n, x, desc_a, info, js)
dectype=desc_a%matrix_data(psb_dec_type_)
!... check m and n parameters....
if (m.lt.0) then
info = 10
int_err(1) = 1
int_err(2) = m
call psb_errpush(info,name,int_err)
goto 9999
else if (n.lt.0) then
info = 10
int_err(1) = 2
int_err(2) = n
call psb_errpush(info,name,int_err)
else if (.not.psb_is_ok_dec(dectype)) then
if (.not.psb_is_ok_dec(dectype)) then
info = 3110
call psb_errpush(info,name)
goto 9999
else if (m.ne.desc_a%matrix_data(psb_n_)) then
info = 300
int_err(1) = 1
int_err(2) = m
int_err(3) = 4
int_err(4) = psb_n_
int_err(5) = desc_a%matrix_data(psb_n_)
call psb_errpush(info,name,int_err)
goto 9999
endif
if (present(js)) then
j=js
if (present(n)) then
n_ = n
else
j=1
n_ = 1
endif
!global check on m and n parameters
!global check on n parameters
if (myrow.eq.psb_root_) then
exch(1)=m
exch(2)=n
exch(3)=j
call igebs2d(icontxt,psb_all_,psb_topdef_, ithree,ione, exch, ithree)
exch(1)=n_
call igebs2d(icontxt,psb_all_,psb_topdef_, ione,ione, exch, ione)
else
call igebr2d(icontxt,psb_all_,psb_topdef_, ithree,ione, exch, ithree, psb_root_, 0)
if (exch(1).ne.m) then
call igebr2d(icontxt,psb_all_,psb_topdef_, ione,ione, exch, ione, psb_root_, 0)
if (exch(1).ne.n_) then
info=550
int_err(1)=1
call psb_errpush(info,name,int_err)
goto 9999
else if (exch(2).ne.n) then
info=550
int_err(1)=2
call psb_errpush(info,name,int_err)
goto 9999
else if (exch(3).ne.j) then
info=550
int_err(1)=3
call psb_errpush(info,name,int_err)
goto 9999
endif
endif
!....allocate x .....
if (psb_is_asb_dec(dectype).or.psb_is_upd_dec(dectype)) then
n_col = max(1,desc_a%matrix_data(psb_n_col_))
allocate(x(n_col,j:j+n-1),stat=info)
allocate(x(n_col,n_),stat=info)
if (info.ne.0) then
info=4010
ch_err='allocate'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
endif
do jj=j,j+n-1
do j=1,n_
do i=1,n_col
x(i,j) = 0.0d0
end do
end do
else if (psb_is_bld_dec(dectype)) then
n_row = max(1,desc_a%matrix_data(psb_n_row_))
allocate(x(n_row,j:j+n-1),stat=info)
allocate(x(n_row,n_),stat=info)
if (info.ne.0) then
info=4010
ch_err='allocate'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
endif
do jj=j,j+n-1
do j = 1, n_
do i=1,n_row
x(i,j) = 0.0d0
end do
@@ -224,11 +189,10 @@ end subroutine psb_dalloc
! Allocates dense matrix for PSBLAS routines
!
! Parameters:
! m - number of rows.
! x - the matrix to be allocated.
! desc_a - the communication descriptor.
! info - eventually returns an error code
subroutine psb_dallocv(m, x, desc_a,info)
! info - possibly returns an error code
subroutine psb_dallocv(x, desc_a,info,n)
!....allocate sparse matrix structure for psblas routines.....
use psb_descriptor_type
use psb_const_mod
@@ -238,14 +202,14 @@ subroutine psb_dallocv(m, x, desc_a,info)
implicit none
!....parameters...
integer, intent(in) :: m
real(kind(1.d0)), pointer :: x(:)
type(psb_desc_type), intent(in):: desc_a
integer :: info
real(kind(1.d0)), pointer :: x(:)
type(psb_desc_type), intent(in) :: desc_a
integer :: info
integer, optional, intent(in) :: n
!locals
integer :: nprow,npcol,myrow,mycol,err,n_col,n_row,dectype,i,err_act
integer :: icontxt
integer :: icontxt, n_
integer :: int_err(5),temp(1),exch
real(kind(1.d0)) :: real_err(5)
logical, parameter :: debug=.false.
@@ -275,40 +239,13 @@ subroutine psb_dallocv(m, x, desc_a,info)
if (debug) write(0,*) 'dall: dectype',dectype
if (debug) write(0,*) 'dall: is_ok? dectype',psb_is_ok_dec(dectype)
!... check m and n parameters....
if (m.lt.0) then
info = 10
int_err(1) = 1
int_err(2) = m
call psb_errpush(info,name,int_err)
goto 9999
else if (.not.psb_is_ok_dec(dectype)) then
if (.not.psb_is_ok_dec(dectype)) then
info = 3110
call psb_errpush(info,name)
goto 9999
else if (m.ne.desc_a%matrix_data(psb_n_)) then
info = 300
int_err(1) = 1
int_err(2) = m
int_err(3) = 4
int_err(4) = psb_n_
int_err(5) = desc_a%matrix_data(psb_n_)
call psb_errpush(info,name,int_err)
goto 9999
endif
!global check on m and n parameters
if (myrow.eq.psb_root_) then
exch = m
call igebs2d(icontxt,psb_all_,psb_topdef_, ione,ione, exch, ione)
else
call igebr2d(icontxt,psb_all_,psb_topdef_, ione,ione, exch, ione, psb_root_, 0)
if (exch .ne. m) then
info = 550
int_err(1) = 1
call psb_errpush(info,name,int_err)
goto 9999
endif
endif
! As this is a rank-1 array, optional parameter N is actually ignored.
!....allocate x .....
if (psb_is_asb_dec(dectype).or.psb_is_upd_dec(dectype)) then
src/tools/psb_dins.f90
+3 -3
View File
@@ -44,7 +44,7 @@
! info - integer. Eventually returns an error code
! iblck - integer(optional). First row of submatrix belonging to blck to be inserted.
! jblck - integer(optional). First col of submatrix belonging to blck to be inserted.
subroutine psb_dins(m, n, x, ix, jx, blck, desc_a, info,&
subroutine psb_dins(m, n, blck, x, ix, jx, desc_a, info,&
& iblck, jblck,dupl)
!....insert dense submatrix to dense matrix .....
use psb_descriptor_type
@@ -270,7 +270,7 @@ end subroutine psb_dins
! desc_a - type(<psb_desc_type>). The communication descriptor.
! info - integer. Eventually returns an error code
! iblck - integer(optional). First row of submatrix belonging to blck to be inserted.
subroutine psb_dinsvm(m, x, ix, jx, blck, desc_a,info,&
subroutine psb_dinsvm(m, blck, x, ix, jx, desc_a,info,&
& iblck,dupl)
!....insert dense submatrix to dense matrix .....
use psb_descriptor_type
@@ -485,7 +485,7 @@ end subroutine psb_dinsvm
! info - integer. Eventually returns an error code
! iblck - integer(optional). First row of submatrix belonging to blck to be inserted.
! insflag - integer(optional). ???
subroutine psb_dinsvv(m, x, ix, blck, desc_a, info,&
subroutine psb_dinsvv(m, blck, x, ix, desc_a, info,&
& iblck,insflag,dupl)
!....insert dense submatrix to dense matrix .....
use psb_descriptor_type
src/tools/psb_ialloc.f90
+88 -146
View File
@@ -34,13 +34,11 @@
! Allocates dense integer matrix for PSBLAS routines
!
! Parameters:
! m - integer. The number of rows.
! n - integer. The number of columns.
! x - integer,dimension(:,:). The matrix to be allocated.
! desc_a - type(<psb_desc_type>). The communication descriptor.
! info - integer. Eventually returns an error code
! js - integer(optional). The starting column
subroutine psb_ialloc(m, n, x, desc_a, info,js)
! x - the matrix to be allocated.
! desc_a - the communication descriptor.
! info - possibly returns an error code
! n - optional number of columns.
subroutine psb_ialloc(x, desc_a, info, n)
!....allocate dense matrix for psblas routines.....
use psb_descriptor_type
use psb_const_mod
@@ -48,19 +46,18 @@ subroutine psb_ialloc(m, n, x, desc_a, info,js)
implicit none
!....parameters...
integer, intent(in) :: m,n
integer, pointer :: x(:,:)
type(psb_desc_type), intent(inout) :: desc_a
integer, intent(out) :: info
integer, optional, intent(in) :: js
integer, pointer :: x(:,:)
type(psb_desc_type), intent(in) :: desc_a
integer, intent(out) :: info
integer, optional, intent(in) :: n
!locals
integer :: j,nprow,npcol,myrow,mypcol,&
& n_col,n_row, err_act
integer :: icontxt,dectype
integer :: nprow,npcol,myrow,mypcol,err,n_col,n_row,i,j,jj,err_act
integer :: icontxt,dectype,n_
integer :: int_err(5),temp(1),exch(3)
real(kind(1.d0)) :: real_err(5)
character(len=20) :: name, char_err
character(len=20) :: name, ch_err
if(psb_get_errstatus().ne.0) return
info=0
@@ -84,97 +81,70 @@ subroutine psb_ialloc(m, n, x, desc_a, info,js)
dectype=desc_a%matrix_data(psb_dec_type_)
!... check m and n parameters....
if (m.lt.0) then
info = 10
int_err(1) = 1
int_err(2) = m
call psb_errpush(info,name,int_err)
goto 9999
else if (n.lt.0) then
info = 10
int_err(1) = 2
int_err(2) = n
call psb_errpush(info,name,int_err)
goto 9999
else if (.not.psb_is_ok_dec(dectype)) then
info = 3110
call psb_errpush(info,name)
goto 9999
else if (m.ne.desc_a%matrix_data(psb_n_)) then
info = 300
int_err(1) = 1
int_err(2) = m
int_err(3) = 4
int_err(4) = psb_n_
int_err(5) = desc_a%matrix_data(psb_n_)
call psb_errpush(info,name,int_err)
goto 9999
if (.not.psb_is_ok_dec(dectype)) then
info = 3110
call psb_errpush(info,name)
goto 9999
endif
if (present(js)) then
j=js
if (present(n)) then
n_ = n
else
j=1
n_ = 1
endif
!global check on m and n parameters
!global check on n parameters
if (myrow.eq.psb_root_) then
exch(1)=m
exch(2)=n
exch(3)=j
call igebs2d(icontxt,psb_all_,psb_topdef_, ithree,ione, exch, ithree)
exch(1)=n_
call igebs2d(icontxt,psb_all_,psb_topdef_, ione,ione, exch, ione)
else
call igebr2d(icontxt,psb_all_,psb_topdef_, ithree,ione, exch, ithree, psb_root_, 0)
if (exch(1).ne.m) then
info=550
int_err(1)=1
call psb_errpush(info,name,int_err)
goto 9999
else if (exch(2).ne.n) then
info=550
int_err(1)=2
call psb_errpush(info,name,int_err)
goto 9999
else if (exch(3).ne.j) then
info=550
int_err(1)=3
call psb_errpush(info,name,int_err)
goto 9999
endif
call igebr2d(icontxt,psb_all_,psb_topdef_, ione,ione, exch, ione, psb_root_, 0)
if (exch(1).ne.n_) then
info=550
int_err(1)=1
call psb_errpush(info,name,int_err)
goto 9999
endif
endif
!....allocate x .....
if (psb_is_asb_dec(dectype).or.psb_is_upd_dec(dectype)) then
n_col = max(1,desc_a%matrix_data(psb_n_col_))
allocate(x(n_col,j:j+n-1),stat=info)
if (info.ne.0) then
info=2025
int_err(1)=n_col
call psb_errpush(info,name,int_err)
goto 9999
endif
n_col = max(1,desc_a%matrix_data(psb_n_col_))
allocate(x(n_col,n_),stat=info)
if (info.ne.0) then
info=4010
ch_err='allocate'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
endif
do j=1,n_
do i=1,n_col
x(i,j) = 0
end do
end do
else if (psb_is_bld_dec(dectype)) then
n_row = max(1,desc_a%matrix_data(psb_n_row_))
allocate(x(n_row,j:j+n-1),stat=info)
if (info.ne.0) then
info=2025
int_err(1)=n_row
call psb_errpush(info,name,int_err)
goto 9999
endif
n_row = max(1,desc_a%matrix_data(psb_n_row_))
allocate(x(n_row,n_),stat=info)
if (info.ne.0) then
info=4010
ch_err='allocate'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
endif
do j = 1, n_
do i=1,n_row
x(i,j) = 0
end do
end do
endif
x = 0
call psb_erractionrestore(err_act)
return
9999 continue
call psb_erractionrestore(err_act)
if (err_act.eq.act_ret) then
return
else
call psb_error(icontxt)
if (err_act.eq.act_abort) then
call psb_error(icontxt)
return
end if
return
@@ -220,26 +190,28 @@ end subroutine psb_ialloc
! x - integer,dimension(:). The matrix to be allocated.
! desc_a - type(<psb_desc_type>). The communication descriptor.
! info - integer. Eventually returns an error code
subroutine psb_iallocv(m, x, desc_a, info)
subroutine psb_iallocv(x, desc_a, info,n)
!....allocate sparse matrix structure for psblas routines.....
use psb_descriptor_type
use psb_const_mod
use psb_realloc_mod
use psb_error_mod
implicit none
!....parameters...
integer, intent(in) :: m
integer, pointer :: x(:)
type(psb_desc_type), intent(in) :: desc_a
integer, intent(out) :: info
integer, pointer :: x(:)
type(psb_desc_type), intent(in) :: desc_a
integer, intent(out) :: info
integer, optional, intent(in) :: n
!locals
integer :: nprow,npcol,myrow,mypcol,err,n_col,n_row,dectype,err_act
integer :: icontxt
integer :: int_err(5),temp(1),exch(2)
integer :: nprow,npcol,myrow,mycol,err,n_col,n_row,dectype,i,err_act
integer :: icontxt, n_
integer :: int_err(5),temp(1),exch
real(kind(1.d0)) :: real_err(5)
logical, parameter :: debug=.false.
character(len=20) :: name, char_err
character(len=20) :: name, ch_err
if(psb_get_errstatus().ne.0) return
info=0
@@ -247,59 +219,31 @@ subroutine psb_iallocv(m, x, desc_a, info)
call psb_erractionsave(err_act)
icontxt=desc_a%matrix_data(psb_ctxt_)
call blacs_gridinfo(icontxt, nprow, npcol, myrow, mypcol)
call blacs_gridinfo(icontxt, nprow, npcol, myrow, mycol)
! ....verify blacs grid correctness..
if (nprow.eq.-1) then
info = 2010
call psb_errpush(info,name)
goto 9999
info = 2010
call psb_errpush(info,name)
goto 9999
else if (npcol.ne.1) then
info = 2030
int_err(1) = npcol
call psb_errpush(info,name,int_err)
goto 9999
info = 2030
int_err(1) = npcol
call psb_errpush(info,name,int_err)
goto 9999
endif
dectype=desc_a%matrix_data(psb_dec_type_)
if (debug) write(0,*) 'dall: dectype',dectype
if (debug) write(0,*) 'dall: is_ok? dectype',psb_is_ok_dec(dectype)
!... check m and n parameters....
if (m.lt.0) then
info = 10
int_err(1) = 1
int_err(2) = m
call psb_errpush(info,name,int_err)
goto 9999
else if (.not.psb_is_ok_dec(dectype)) then
info = 3110
call psb_errpush(info,name)
goto 9999
else if (m.ne.desc_a%matrix_data(psb_n_)) then
info = 300
int_err(1) = 1
int_err(2) = m
int_err(3) = 4
int_err(4) = psb_n_
int_err(5) = desc_a%matrix_data(psb_n_)
call psb_errpush(info,name,int_err)
goto 9999
endif
!global check on m and n parameters
if (myrow.eq.psb_root_) then
exch(1) = m
call igebs2d(icontxt,psb_all_,psb_topdef_, ione,ione, exch, ione)
else
call igebr2d(icontxt,psb_all_,psb_topdef_, ione,ione, exch, ione, psb_root_, 0)
if (exch(1) .ne. m) then
info = 550
int_err(1) = 1
call psb_errpush(info,name,int_err)
goto 9999
endif
if (.not.psb_is_ok_dec(dectype)) then
info = 3110
call psb_errpush(info,name)
goto 9999
endif
! As this is a rank-1 array, optional parameter N is actually ignored.
!....allocate x .....
if (psb_is_asb_dec(dectype).or.psb_is_upd_dec(dectype)) then
@@ -329,11 +273,9 @@ subroutine psb_iallocv(m, x, desc_a, info)
9999 continue
call psb_erractionrestore(err_act)
if (err_act.eq.act_ret) then
return
else
call psb_error(icontxt)
if (err_act.eq.act_abort) then
call psb_error(icontxt)
return
end if
return
src/tools/psb_iins.f90
+3 -3
View File
@@ -44,7 +44,7 @@
! info - integer. Eventually returns an error code
! iblck - integer(optional). First row of submatrix belonging to blck to be inserted.
! jblck - integer(optional). First col of submatrix belonging to blck to be inserted.
subroutine psb_iins(m, n, x, ix, jx, blck, desc_a, info,&
subroutine psb_iins(m, n, blck, x, ix, jx, desc_a, info,&
& iblck, jblck,dupl)
!....insert dense submatrix to dense matrix .....
use psb_descriptor_type
@@ -266,7 +266,7 @@ end subroutine psb_iins
! desc_a - type(<psb_desc_type>). The communication descriptor.
! info - integer. Eventually returns an error code
! iblck - integer(optional). First row of submatrix belonging to blck to be inserted.
subroutine psb_iinsvm(m, x, ix, jx, blck, desc_a, info,&
subroutine psb_iinsvm(m, blck, x, ix, jx, desc_a, info,&
& iblck,dupl)
!....insert dense submatrix to dense matrix .....
use psb_descriptor_type
@@ -413,7 +413,7 @@ end subroutine psb_iinsvm
! desc_a - type(<psb_desc_type>). The communication descriptor.
! info - integer. Eventually returns an error code
! iblck - integer(optional). First row of submatrix belonging to blck to be inserted.
subroutine psb_iinsvv(m, x, ix, blck, desc_a, info,&
subroutine psb_iinsvv(m, blck, x, ix, desc_a, info,&
& iblck,dupl)
!....insert dense submatrix to dense matrix .....
use psb_descriptor_type
src/tools/psb_zallc.f90
+28 -92
View File
@@ -35,13 +35,11 @@
! Allocates dense matrix for PSBLAS routines
!
! Parameters:
! m - number of rows.
! n - number of columns.
! x - the matrix to be allocated.
! desc_a - the communication descriptor.
! info - eventually returns an error code
! js - (optional) the starting column
subroutine psb_zalloc(m, n, x, desc_a, info, js)
! info - possibly returns an error code
! n - optional number of columns.
subroutine psb_zalloc(x, desc_a, info, n)
!....allocate dense matrix for psblas routines.....
use psb_descriptor_type
use psb_const_mod
@@ -50,18 +48,16 @@ subroutine psb_zalloc(m, n, x, desc_a, info, js)
implicit none
!....parameters...
integer, intent(in) :: m,n
complex(kind(1.d0)), pointer :: x(:,:)
complex(kind(1.d0)), pointer :: x(:,:)
type(psb_desc_type), intent(in) :: desc_a
integer :: info
integer, optional, intent(in) :: js
integer, optional, intent(in) :: n
!locals
integer :: nprow,npcol,myrow,mycol,err,n_col,n_row,i,j,jj,err_act
integer :: icontxt,dectype
integer :: icontxt,dectype,n_
integer :: int_err(5),temp(1),exch(3)
real(kind(1.d0)) :: real_err(5)
integer, allocatable:: prc_v(:)
character(len=20) :: name, ch_err
name='psb_zallc'
@@ -88,88 +84,56 @@ subroutine psb_zalloc(m, n, x, desc_a, info, js)
dectype=desc_a%matrix_data(psb_dec_type_)
!... check m and n parameters....
if (m.lt.0) then
info = 10
int_err(1) = 1
int_err(2) = m
call psb_errpush(info,name,int_err)
goto 9999
else if (n.lt.0) then
info = 10
int_err(1) = 2
int_err(2) = n
call psb_errpush(info,name,int_err)
else if (.not.psb_is_ok_dec(dectype)) then
if (.not.psb_is_ok_dec(dectype)) then
info = 3110
call psb_errpush(info,name)
goto 9999
else if (m.ne.desc_a%matrix_data(psb_n_)) then
info = 300
int_err(1) = 1
int_err(2) = m
int_err(3) = 4
int_err(4) = psb_n_
int_err(5) = desc_a%matrix_data(psb_n_)
call psb_errpush(info,name,int_err)
goto 9999
endif
if (present(js)) then
j=js
if (present(n)) then
n_ = n
else
j=1
n_ = 1
endif
!global check on m and n parameters
!global check on n parameters
if (myrow.eq.psb_root_) then
exch(1)=m
exch(2)=n
exch(3)=j
call igebs2d(icontxt,psb_all_,psb_topdef_, ithree,ione, exch, ithree)
exch(1)=n_
call igebs2d(icontxt,psb_all_,psb_topdef_, ione,ione, exch, ione)
else
call igebr2d(icontxt,psb_all_,psb_topdef_, ithree,ione, exch, ithree, psb_root_, 0)
if (exch(1).ne.m) then
call igebr2d(icontxt,psb_all_,psb_topdef_, ione,ione, exch, ione, psb_root_, 0)
if (exch(1).ne.n_) then
info=550
int_err(1)=1
call psb_errpush(info,name,int_err)
goto 9999
else if (exch(2).ne.n) then
info=550
int_err(1)=2
call psb_errpush(info,name,int_err)
goto 9999
else if (exch(3).ne.j) then
info=550
int_err(1)=3
call psb_errpush(info,name,int_err)
goto 9999
endif
endif
!....allocate x .....
if (psb_is_asb_dec(dectype).or.psb_is_upd_dec(dectype)) then
n_col = max(1,desc_a%matrix_data(psb_n_col_))
allocate(x(n_col,j:j+n-1),stat=info)
allocate(x(n_col,n_),stat=info)
if (info.ne.0) then
info=4010
ch_err='allocate'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
endif
do jj=j,j+n-1
do j=1,n_
do i=1,n_col
x(i,j) = 0.0d0
end do
end do
else if (psb_is_bld_dec(dectype)) then
n_row = max(1,desc_a%matrix_data(psb_n_row_))
allocate(x(n_row,j:j+n-1),stat=info)
allocate(x(n_row,n_),stat=info)
if (info.ne.0) then
info=4010
ch_err='allocate'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
endif
do jj=j,j+n-1
do j = 1, n_
do i=1,n_row
x(i,j) = 0.0d0
end do
@@ -224,11 +188,10 @@ end subroutine psb_zalloc
! Allocates dense matrix for PSBLAS routines
!
! Parameters:
! m - number of rows.
! x - the matrix to be allocated.
! desc_a - the communication descriptor.
! info - eventually returns an error code
subroutine psb_zallocv(m, x, desc_a,info)
! info - possibly returns an error code
subroutine psb_zallocv(x, desc_a,info,n)
!....allocate sparse matrix structure for psblas routines.....
use psb_descriptor_type
use psb_const_mod
@@ -238,14 +201,14 @@ subroutine psb_zallocv(m, x, desc_a,info)
implicit none
!....parameters...
integer, intent(in) :: m
complex(kind(1.d0)), pointer :: x(:)
type(psb_desc_type), intent(in):: desc_a
integer :: info
complex(kind(1.d0)), pointer :: x(:)
type(psb_desc_type), intent(in) :: desc_a
integer :: info
integer, optional, intent(in) :: n
!locals
integer :: nprow,npcol,myrow,mycol,err,n_col,n_row,dectype,i,err_act
integer :: icontxt
integer :: icontxt, n_
integer :: int_err(5),temp(1),exch
real(kind(1.d0)) :: real_err(5)
logical, parameter :: debug=.false.
@@ -275,40 +238,13 @@ subroutine psb_zallocv(m, x, desc_a,info)
if (debug) write(0,*) 'dall: dectype',dectype
if (debug) write(0,*) 'dall: is_ok? dectype',psb_is_ok_dec(dectype)
!... check m and n parameters....
if (m.lt.0) then
info = 10
int_err(1) = 1
int_err(2) = m
call psb_errpush(info,name,int_err)
goto 9999
else if (.not.psb_is_ok_dec(dectype)) then
if (.not.psb_is_ok_dec(dectype)) then
info = 3110
call psb_errpush(info,name)
goto 9999
else if (m.ne.desc_a%matrix_data(psb_n_)) then
info = 300
int_err(1) = 1
int_err(2) = m
int_err(3) = 4
int_err(4) = psb_n_
int_err(5) = desc_a%matrix_data(psb_n_)
call psb_errpush(info,name,int_err)
goto 9999
endif
!global check on m and n parameters
if (myrow.eq.psb_root_) then
exch = m
call igebs2d(icontxt,psb_all_,psb_topdef_, ione,ione, exch, ione)
else
call igebr2d(icontxt,psb_all_,psb_topdef_, ione,ione, exch, ione, psb_root_, 0)
if (exch .ne. m) then
info = 550
int_err(1) = 1
call psb_errpush(info,name,int_err)
goto 9999
endif
endif
! As this is a rank-1 array, optional parameter N is actually ignored.
!....allocate x .....
if (psb_is_asb_dec(dectype).or.psb_is_upd_dec(dectype)) then
src/tools/psb_zins.f90
+3 -3
View File
@@ -44,7 +44,7 @@
! info - integer. Eventually returns an error code
! iblck - integer(optional). First row of submatrix belonging to blck to be inserted.
! jblck - integer(optional). First col of submatrix belonging to blck to be inserted.
subroutine psb_zins(m, n, x, ix, jx, blck, desc_a, info,&
subroutine psb_zins(m, n, blck, x, ix, jx, desc_a, info,&
& iblck, jblck,dupl)
!....insert dense submatrix to dense matrix .....
use psb_descriptor_type
@@ -268,7 +268,7 @@ end subroutine psb_zins
! desc_a - type(<psb_desc_type>). The communication descriptor.
! info - integer. Eventually returns an error code
! iblck - integer(optional). First row of submatrix belonging to blck to be inserted.
subroutine psb_zinsvm(m, x, ix, jx, blck, desc_a,info,&
subroutine psb_zinsvm(m, blck, x, ix, jx, desc_a,info,&
& iblck,dupl)
!....insert dense submatrix to dense matrix .....
use psb_descriptor_type
@@ -483,7 +483,7 @@ end subroutine psb_zinsvm
! info - integer. Eventually returns an error code
! iblck - integer(optional). First row of submatrix belonging to blck to be inserted.
! insflag - integer(optional). ???
subroutine psb_zinsvv(m, x, ix, blck, desc_a, info,&
subroutine psb_zinsvv(m, blck, x, ix, desc_a, info,&
& iblck,insflag,dupl)
!....insert dense submatrix to dense matrix .....
use psb_descriptor_type
src/tools/psb_zspasb.f90
+1
View File
@@ -162,6 +162,7 @@ subroutine psb_zspasb(a,desc_a, info, afmt, upd, dupl)
! Do the real conversion into the requested storage format
! result is put in A
write(0,*) 'Calling csdp from SPASB'
call psb_csdp(atemp,a,info,ifc=2,upd=upd_,dupl=dupl_)
IF (debug) WRITE (*, *) myrow,' ASB: From DCSDP',info,' ',A%FIDA