Folded in a bunch of fixes coming from testing on IBM SP. A few

wrinkles are still out there.
psblas3-type-indexed
Salvatore Filippone 19 years ago
parent 2c21b017b3
commit 0b428f4c7d

@ -1,6 +1,9 @@
Changelog. A lot less detailed than usual, at least for past
history.
2006/04/21: A bunch of fixes related to various scratch spmat(s) initialization
problems that were revealed while testing on SP5.
2006/04/18: Changed interface to spasb and csdp: better handling of
regeneration. To be tested further for sophisticated uses.

@ -19,8 +19,8 @@ Make.inc file; we have tested with AIX XLF, Intel ifc/Linux, Lahey
F95/Linux, Nag f95/Linux, GNU Fortran/Linux. If you succeed in compiling with
other compiler/operating systems please let us know.
IBM SP2.
The library has been tested on an IBM SP2 with XLC and XLF
IBM SP.
The library has been tested on an IBM SP2, SP4 and SP5, with XLC and XLF
compilers, and a version of the BLACS based on MPI.
The rather baroque setting
F90=xlf90 -qsuffix=f=f90
@ -117,8 +117,8 @@ Credits for version 2.0:
Salvatore Filippone
Alfredo Buttari
The MPcube preconditioners contained in directory src/prec were
developed with the contribution of:
The multilevel parallel preconditioners contained in directory
src/prec were developed with the contribution of:
Pasqua D'Ambra
Daniela Di Serafino

@ -348,7 +348,7 @@ module psb_psblas_mod
real(kind(1.d0)), intent(in) :: alpha, beta
type(psb_desc_type), intent(in) :: desc_a
character, optional, intent(in) :: trans
real(kind(1.d0)), optional, intent(inout) :: work(:)
real(kind(1.d0)), optional, intent(inout),target :: work(:)
integer, optional, intent(in) :: k, jx, jy,doswap
integer, intent(out) :: info
end subroutine psb_dspmm
@ -362,7 +362,7 @@ module psb_psblas_mod
real(kind(1.d0)), intent(in) :: alpha, beta
type(psb_desc_type), intent(in) :: desc_a
character, optional, intent(in) :: trans
real(kind(1.d0)), optional, intent(inout) :: work(:)
real(kind(1.d0)), optional, intent(inout),target :: work(:)
integer, optional, intent(in) :: doswap
integer, intent(out) :: info
end subroutine psb_dspmv
@ -376,7 +376,7 @@ module psb_psblas_mod
complex(kind(1.d0)), intent(in) :: alpha, beta
type(psb_desc_type), intent(in) :: desc_a
character, optional, intent(in) :: trans
complex(kind(1.d0)), optional, intent(inout) :: work(:)
complex(kind(1.d0)), optional, intent(inout),target :: work(:)
integer, optional, intent(in) :: k, jx, jy,doswap
integer, intent(out) :: info
end subroutine psb_zspmm
@ -390,7 +390,7 @@ module psb_psblas_mod
complex(kind(1.d0)), intent(in) :: alpha, beta
type(psb_desc_type), intent(in) :: desc_a
character, optional, intent(in) :: trans
complex(kind(1.d0)), optional, intent(inout) :: work(:)
complex(kind(1.d0)), optional, intent(inout),target :: work(:)
integer, optional, intent(in) :: doswap
integer, intent(out) :: info
end subroutine psb_zspmv
@ -410,7 +410,7 @@ module psb_psblas_mod
character, optional, intent(in) :: trans, unit
integer, optional, intent(in) :: n, jx, jy
integer, optional, intent(in) :: choice
real(kind(1.d0)), optional, intent(in) :: work(:), diag(:)
real(kind(1.d0)), optional, intent(inout),target :: work(:), diag(:)
integer, intent(out) :: info
end subroutine psb_dspsm
subroutine psb_dspsv(alpha, t, x, beta, y,&
@ -425,7 +425,7 @@ module psb_psblas_mod
type(psb_desc_type), intent(in) :: desc_a
character, optional, intent(in) :: trans, unit
integer, optional, intent(in) :: choice
real(kind(1.d0)), optional, intent(in) :: work(:), diag(:)
real(kind(1.d0)), optional, intent(inout),target :: work(:), diag(:)
integer, intent(out) :: info
end subroutine psb_dspsv
subroutine psb_zspsm(alpha, t, x, beta, y,&
@ -441,7 +441,7 @@ module psb_psblas_mod
character, optional, intent(in) :: trans, unit
integer, optional, intent(in) :: n, jx, jy
integer, optional, intent(in) :: choice
complex(kind(1.d0)), optional, intent(in) :: work(:), diag(:)
complex(kind(1.d0)), optional, intent(inout),target :: work(:), diag(:)
integer, intent(out) :: info
end subroutine psb_zspsm
subroutine psb_zspsv(alpha, t, x, beta, y,&
@ -456,7 +456,7 @@ module psb_psblas_mod
type(psb_desc_type), intent(in) :: desc_a
character, optional, intent(in) :: trans, unit
integer, optional, intent(in) :: choice
complex(kind(1.d0)), optional, intent(in) :: work(:), diag(:)
complex(kind(1.d0)), optional, intent(inout),target :: work(:), diag(:)
integer, intent(out) :: info
end subroutine psb_zspsv
end interface

@ -61,13 +61,15 @@ Contains
Integer,Pointer :: tmp(:)
Integer :: dim, err_act, err,i
character(len=20) :: name
logical, parameter :: debug=.false.
name='psb_dreallocate1i'
call psb_erractionsave(err_act)
if(psb_get_errstatus().ne.0) return
info=0
if (associated(rrax)) then
if (debug) write(0,*) 'reallocate I',len
if (associated(rrax)) then a
dim=size(rrax)
If (dim /= len) Then
Allocate(tmp(len),stat=info)
@ -137,9 +139,12 @@ Contains
Real(kind(1.d0)),Pointer :: tmp(:)
Integer :: dim,err_act,err,i, m
character(len=20) :: name
logical, parameter :: debug=.false.
name='psb_dreallocate1d'
call psb_erractionsave(err_act)
info = 0
if (debug) write(0,*) 'reallocate D',len
if (associated(rrax)) then
dim=size(rrax)
@ -210,9 +215,12 @@ Contains
complex(kind(1.d0)),Pointer :: tmp(:)
Integer :: dim,err_act,err,i, m
character(len=20) :: name
logical, parameter :: debug=.false.
name='psb_dreallocate1z'
call psb_erractionsave(err_act)
info = 0
if (debug) write(0,*) 'reallocate Z',len
if (associated(rrax)) then
dim=size(rrax)
@ -286,6 +294,7 @@ Contains
name='psb_dreallocated2'
call psb_erractionsave(err_act)
info = 0
if (associated(rrax)) then
dim=size(rrax,1)
@ -357,6 +366,7 @@ Contains
name='psb_dreallocatez2'
call psb_erractionsave(err_act)
info = 0
if (associated(rrax)) then
dim=size(rrax,1)

@ -98,7 +98,7 @@ module psb_spmat_type
interface psb_sp_reall
module procedure psb_dspreallocate, psb_dspreall3, &
& psb_zspreallocate, psb_zspreall3
& psb_zspreall3, psb_zspreallocate
end interface
interface psb_sp_all
@ -757,6 +757,32 @@ contains
End Subroutine psb_zspall3
subroutine psb_zspreall3(a, ni1,ni2,nz,info)
implicit none
!....Parameters...
Type(psb_zspmat_type), intent(inout) :: A
Integer, intent(in) :: ni1,ni2,nz
Integer, intent(inout) :: info
!locals
logical, parameter :: debug=.false.
info = 0
call psb_realloc(nz,a%aspk,info)
if (info /= 0) return
call psb_realloc(ni2,a%ia2,info)
if (info /= 0) return
call psb_realloc(ni1,a%ia1,info)
if (info /= 0) return
call psb_realloc(max(1,a%m),a%pl,info)
if (info /= 0) return
call psb_realloc(max(1,a%k),a%pr,info)
if (info /= 0) return
Return
End Subroutine psb_zspreall3
subroutine psb_zspreallocate(a, nnz,info,ifc)
implicit none
@ -807,33 +833,6 @@ contains
End Subroutine psb_zspreallocate
subroutine psb_zspreall3(a, ni1,ni2,nd,info)
implicit none
!....Parameters...
Type(psb_zspmat_type), intent(inout) :: A
Integer, intent(in) :: ni1,ni2,nd
Integer, intent(inout) :: info
!locals
logical, parameter :: debug=.false.
info = 0
call psb_realloc(nd,a%aspk,info)
if (info /= 0) return
call psb_realloc(ni2,a%ia2,info)
if (info /= 0) return
call psb_realloc(ni1,a%ia1,info)
if (info /= 0) return
call psb_realloc(max(1,a%m),a%pl,info)
if (info /= 0) return
call psb_realloc(max(1,a%k),a%pr,info)
if (info /= 0) return
Return
End Subroutine psb_zspreall3
subroutine psb_zspclone(a, b,info)
implicit none
!....Parameters...

@ -102,6 +102,7 @@ Subroutine psb_dasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
! Block Jacobi. Copy the descriptor, just in case we want to
! do the renumbering.
!
If(debug) Write(0,*)' asmatbld calling allocate '
call psb_sp_all(0,0,blk,1,info)
if(info /= 0) then
info=4010
@ -111,9 +112,10 @@ Subroutine psb_dasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
end if
blk%fida = 'COO'
blk%infoa(psb_nnz_) = 0
If(debug) Write(0,*)' asmatbld done spallocate'
If (upd == 'F') Then
call psb_cdcpy(desc_data,desc_p,info)
If(debug) Write(0,*)' asmatbld done cdcpy'
if(info /= 0) then
info=4010
ch_err='psb_cdcpy'
@ -143,6 +145,7 @@ Subroutine psb_dasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
!
! This is really just Block Jacobi.....
!
If(debug) Write(0,*)' asmatbld calling allocate novr=0'
call psb_sp_all(0,0,blk,1,info)
if(info /= 0) then
info=4010
@ -155,6 +158,7 @@ Subroutine psb_dasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
if (debug) write(0,*) 'Calling desccpy'
if (upd == 'F') then
call psb_cdcpy(desc_data,desc_p,info)
If(debug) Write(0,*)' asmatbld done cdcpy'
if(info /= 0) then
info=4010
ch_err='psb_cdcpy'

@ -223,6 +223,7 @@ subroutine psb_dbaseprc_bld(a,desc_a,p,info,upd)
case(f_umf_)
if(debug) write(0,*)me,': calling umf_bld'
call psb_umf_bld(a,desc_a,p,info)
if(debug) write(0,*)me,': Done umf_bld ',info
if(info /= 0) then
info=4010
ch_err='umf_bld'

@ -64,7 +64,7 @@ subroutine psb_dbldaggrmat(a,desc_a,ac,p,desc_p,info)
select case (p%iprcparm(smth_kind_))
case (no_smth_)
call raw_aggregate(info)
call raw_aggregate(info)
if(info /= 0) then
call psb_errpush(4010,name,a_err='raw_aggregate')
@ -120,6 +120,7 @@ contains
call psb_erractionsave(err_act)
bg => ac
call psb_nullify_sp(b)
icontxt = desc_a%matrix_data(psb_ctxt_)
call blacs_gridinfo(icontxt,nprows,npcols,myprow,mypcol)
@ -347,7 +348,7 @@ contains
type(psb_dspmat_type) :: am3,am4
logical :: ml_global_nmb
logical, parameter :: test_dump=.false.
logical, parameter :: test_dump=.false.,debug=.false.
integer, parameter :: ncmax=16
real(kind(1.d0)) :: omega, anorm, tmp, dg
character(len=20) :: name, ch_err
@ -362,6 +363,9 @@ contains
call blacs_gridinfo(icontxt,nprows,npcols,myprow,mypcol)
bg => ac
call psb_nullify_sp(b)
call psb_nullify_sp(am3)
call psb_nullify_sp(am4)
am2 => p%av(sm_pr_t_)
am1 => p%av(sm_pr_)
@ -441,7 +445,15 @@ contains
! 1. Allocate Ptilde in sparse matrix form
call psb_sp_all(am4,ncol,info)
am4%fida='COO'
am4%m=ncol
if (ml_global_nmb) then
am4%k=ntaggr
call psb_sp_all(ncol,ntaggr,am4,ncol,info)
else
am4%k=naggr
call psb_sp_all(ncol,naggr,am4,ncol,info)
endif
if(info /= 0) then
call psb_errpush(4010,name,a_err='spall')
goto 9999
@ -462,14 +474,6 @@ contains
end do
am4%infoa(psb_nnz_) = nrow
endif
am4%fida='COO'
am4%m=ncol
if (ml_global_nmb) then
am4%k=ntaggr
else
am4%k=naggr
endif
if (test_dump) call &
@ -560,6 +564,7 @@ contains
if (test_dump) call psb_csprt(40+me,am3,head='% (I-wDA)',ivr=desc_a%loc_to_glob,&
& ivc=desc_a%loc_to_glob)
if (debug) write(0,*) me,'Done gather, going for SYMBMM 1'
!
! Symbmm90 does the allocation for its result.
!
@ -570,6 +575,7 @@ contains
call psb_symbmm(am3,am4,am1)
call psb_numbmm(am3,am4,am1)
if (debug) write(0,*) me,'Done NUMBMM 1'
call psb_sp_free(am4,info)
if(info /= 0) then
@ -615,6 +621,7 @@ contains
call psb_symbmm(a,am1,am3)
call psb_numbmm(a,am1,am3)
if (debug) write(0,*) me,'Done NUMBMM 2'
if (p%iprcparm(smth_kind_) == smth_omg_) then
call psb_transp(am1,am2,fmt='COO')
@ -638,6 +645,7 @@ contains
else
call psb_transp(am1,am2)
endif
if (debug) write(0,*) me,'starting sphalo/ rwxtd'
if (p%iprcparm(smth_kind_) == smth_omg_) then
! am2 = ((i-wDA)Ptilde)^T
@ -667,8 +675,11 @@ contains
end if
endif
if (debug) write(0,*) me,'starting symbmm 3'
call psb_symbmm(am2,am3,b)
if (debug) write(0,*) me,'starting numbmm 3'
call psb_numbmm(am2,am3,b)
if (debug) write(0,*) me,'Done NUMBMM 3'
!!$ if (aggr_dump) call csprt(50+me,am1,head='% Operator PTrans.')
call psb_sp_free(am3,info)
@ -731,6 +742,7 @@ contains
goto 9999
end if
if (debug) write(0,*) me,'Created aux descr. distr.'
call psb_cdasb(desc_p,info)
if(info /= 0) then
call psb_errpush(4010,name,a_err='psb_cdasb')
@ -738,6 +750,7 @@ contains
end if
if (debug) write(0,*) me,'Asmbld aux descr. distr.'
call psb_glob_to_loc(bg%ia1(1:nzl),desc_p,info,iact='I')
if(info /= 0) then

@ -157,6 +157,8 @@ subroutine psb_dilu_bld(a,desc_a,p,upd,info)
icontxt=desc_a%matrix_data(psb_ctxt_)
call psb_nullify_sp(blck)
call psb_nullify_sp(atmp)
t1= mpi_wtime()
if(debug) write(0,*)me,': calling psb_asmatbld',p%iprcparm(p_type_),p%iprcparm(n_ovr_)

@ -151,7 +151,7 @@ contains
if(psb_get_errstatus().ne.0) return
info=0
call psb_erractionsave(err_act)
call psb_nullify_sp(trw)
trw%m=0
trw%k=0
if(debug) write(0,*)'LUINT Allocating TRW'
@ -300,7 +300,7 @@ contains
!
info = 2
int_err(1) = i
write(ch_err,'(g20.10)'),dia
write(ch_err,'(g20.10)') dia
call psb_errpush(info,name,i_err=int_err,a_err=ch_err)
goto 9999
else
@ -439,7 +439,7 @@ contains
! Pivot too small: unstable factorization
!
int_err(1) = i
write(ch_err,'(g20.10)'),dia
write(ch_err,'(g20.10)') dia
info = 2
call psb_errpush(info,name,i_err=int_err,a_err=ch_err)
goto 9999

@ -150,6 +150,7 @@ subroutine psb_dmlprc_bld(a,desc_a,p,info)
goto 9999
end if
if (debug) write(0,*) 'Out from genaggrmap',p%nlaggr
nullify(desc_p)
allocate(desc_p)
call psb_nullify_desc(desc_p)
@ -165,6 +166,7 @@ subroutine psb_dmlprc_bld(a,desc_a,p,info)
call psb_baseprc_bld(ac,desc_p,p,info)
if (debug) write(0,*) 'Out from basaeprcbld',info
!
! We have used a separate ac because:

@ -68,7 +68,7 @@ subroutine psb_dslu_bld(a,desc_a,p,info)
integer, intent(out) :: info
character(len=5), optional :: outfmt
end Subroutine psb_dasmatbld
end interface
end interface
if(psb_get_errstatus().ne.0) return
info=0
@ -82,120 +82,121 @@ subroutine psb_dslu_bld(a,desc_a,p,info)
fmt = 'COO'
call psb_nullify_sp(blck)
call psb_nullify_sp(atmp)
atmp%fida='COO'
if (Debug) then
write(0,*) me, 'SPLUBLD: Calling csdp'
call blacs_barrier(icontxt,'All')
write(0,*) me, 'SPLUBLD: Calling csdp'
call blacs_barrier(icontxt,'All')
endif
call psb_dcsdp(a,atmp,info)
if(info /= 0) then
info=4010
ch_err='psb_dcsdp'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_dcsdp'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
nza = atmp%infoa(psb_nnz_)
if (Debug) then
write(0,*) me, 'SPLUBLD: Done csdp',info,nza,atmp%m,atmp%k
call blacs_barrier(icontxt,'All')
write(0,*) me, 'SPLUBLD: Done csdp',info,nza,atmp%m,atmp%k
call blacs_barrier(icontxt,'All')
endif
call psb_asmatbld(p%iprcparm(p_type_),p%iprcparm(n_ovr_),a,&
& blck,desc_a,upd,p%desc_data,info,outfmt=fmt)
if(info /= 0) then
info=4010
ch_err='psb_asmatbld'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_asmatbld'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
nzb = blck%infoa(psb_nnz_)
if (Debug) then
write(0,*) me, 'SPLUBLD: Done asmatbld',info,nzb,blck%fida
call blacs_barrier(icontxt,'All')
write(0,*) me, 'SPLUBLD: Done asmatbld',info,nzb,blck%fida
call blacs_barrier(icontxt,'All')
endif
if (nzb > 0 ) then
if (size(atmp%aspk)<nza+nzb) then
call psb_sp_reall(atmp,nza+nzb,info)
if(info /= 0) then
info=4010
ch_err='psb_sp_reall'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
endif
if (Debug) then
write(0,*) me, 'SPLUBLD: Done realloc',info,nza+nzb,atmp%fida
call blacs_barrier(icontxt,'All')
endif
do j=1,nzb
atmp%aspk(nza+j) = blck%aspk(j)
atmp%ia1(nza+j) = blck%ia1(j)
atmp%ia2(nza+j) = blck%ia2(j)
end do
atmp%infoa(psb_nnz_) = nza+nzb
atmp%m = atmp%m + blck%m
atmp%k = max(a%k,blck%k)
if (size(atmp%aspk)<nza+nzb) then
call psb_sp_reall(atmp,nza+nzb,info)
if(info /= 0) then
info=4010
ch_err='psb_sp_reall'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
endif
if (Debug) then
write(0,*) me, 'SPLUBLD: Done realloc',info,nza+nzb,atmp%fida
call blacs_barrier(icontxt,'All')
endif
do j=1,nzb
atmp%aspk(nza+j) = blck%aspk(j)
atmp%ia1(nza+j) = blck%ia1(j)
atmp%ia2(nza+j) = blck%ia2(j)
end do
atmp%infoa(psb_nnz_) = nza+nzb
atmp%m = atmp%m + blck%m
atmp%k = max(a%k,blck%k)
else
atmp%infoa(psb_nnz_) = nza
atmp%m = a%m
atmp%k = a%k
atmp%infoa(psb_nnz_) = nza
atmp%m = a%m
atmp%k = a%k
endif
i=0
do j=1, atmp%infoa(psb_nnz_)
if (atmp%ia2(j) <= atmp%m) then
i = i + 1
atmp%aspk(i) = atmp%aspk(j)
atmp%ia1(i) = atmp%ia1(j)
atmp%ia2(i) = atmp%ia2(j)
endif
if (atmp%ia2(j) <= atmp%m) then
i = i + 1
atmp%aspk(i) = atmp%aspk(j)
atmp%ia1(i) = atmp%ia1(j)
atmp%ia2(i) = atmp%ia2(j)
endif
enddo
atmp%infoa(psb_nnz_) = i
call psb_ipcoo2csr(atmp,info)
if(info /= 0) then
info=4010
ch_err='psb_ipcoo2csr'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_ipcoo2csr'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call psb_spinfo(psb_nztotreq_,atmp,nzt,info)
if(info /= 0) then
info=4010
ch_err='psb_spinfo'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_spinfo'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if (Debug) then
write(0,*) me,'Calling psb_slu_factor ',nzt,atmp%m,&
& atmp%k,p%desc_data%matrix_data(psb_n_row_)
call blacs_barrier(icontxt,'All')
write(0,*) me,'Calling psb_slu_factor ',nzt,atmp%m,&
& atmp%k,p%desc_data%matrix_data(psb_n_row_)
call blacs_barrier(icontxt,'All')
endif
call psb_dslu_factor(atmp%m,nzt,&
& atmp%aspk,atmp%ia2,atmp%ia1,p%iprcparm(slu_ptr_),info)
if(info /= 0) then
info=4010
ch_err='psb_slu_fact'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_slu_fact'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if (Debug) then
write(0,*) me, 'SPLUBLD: Done slu_Factor',info,p%iprcparm(slu_ptr_)
call blacs_barrier(icontxt,'All')
write(0,*) me, 'SPLUBLD: Done slu_Factor',info,p%iprcparm(slu_ptr_)
call blacs_barrier(icontxt,'All')
endif
call psb_sp_free(blck,info)
call psb_sp_free(atmp,info)
if(info /= 0) then
info=4010
ch_err='psb_sp_free'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_sp_free'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call psb_erractionrestore(err_act)
@ -204,8 +205,8 @@ subroutine psb_dslu_bld(a,desc_a,p,info)
9999 continue
call psb_erractionrestore(err_act)
if (err_act.eq.act_abort) then
call psb_error()
return
call psb_error()
return
end if
return

@ -82,6 +82,7 @@ subroutine psb_dumf_bld(a,desc_a,p,info)
fmt = 'COO'
call psb_nullify_sp(blck)
call psb_nullify_sp(atmp)
atmp%fida='COO'
if (Debug) then
write(0,*) me, 'UMFBLD: Calling csdp'

@ -102,6 +102,7 @@ Subroutine psb_zasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
! Block Jacobi. Copy the descriptor, just in case we want to
! do the renumbering.
!
If(debug) Write(0,*)' asmatbld calling allocate '
call psb_sp_all(0,0,blk,1,info)
if(info /= 0) then
info=4010
@ -111,9 +112,10 @@ Subroutine psb_zasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
end if
blk%fida = 'COO'
blk%infoa(psb_nnz_) = 0
If(debug) Write(0,*)' asmatbld done spallocate'
If (upd == 'F') Then
call psb_cdcpy(desc_data,desc_p,info)
If(debug) Write(0,*)' asmatbld done cdcpy'
if(info /= 0) then
info=4010
ch_err='psb_cdcpy'
@ -143,6 +145,7 @@ Subroutine psb_zasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
!
! This is really just Block Jacobi.....
!
If(debug) Write(0,*)' asmatbld calling allocate novr=0'
call psb_sp_all(0,0,blk,1,info)
if(info /= 0) then
info=4010
@ -155,6 +158,7 @@ Subroutine psb_zasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
if (debug) write(0,*) 'Calling desccpy'
if (upd == 'F') then
call psb_cdcpy(desc_data,desc_p,info)
If(debug) Write(0,*)' asmatbld done cdcpy'
if(info /= 0) then
info=4010
ch_err='psb_cdcpy'

@ -194,13 +194,14 @@ subroutine psb_zbaseprc_bld(a,desc_a,p,info,upd)
& f_ilu_n_,is_legal_ml_fact)
if (debug) write(0,*)me, ': Calling PSB_ILU_BLD'
if (debug) call blacs_barrier(icontxt,'All')
select case(p%iprcparm(f_type_))
case(f_ilu_n_,f_ilu_e_)
call psb_ilu_bld(a,desc_a,p,iupd,info)
if(debug) write(0,*)me,': out of psb_ilu_bld'
if (debug) call blacs_barrier(icontxt,'All')
if(info /= 0) then
info=4010
ch_err='psb_ilu_bld'
@ -222,6 +223,7 @@ subroutine psb_zbaseprc_bld(a,desc_a,p,info,upd)
case(f_umf_)
if(debug) write(0,*)me,': calling umf_bld'
call psb_umf_bld(a,desc_a,p,info)
if(debug) write(0,*)me,': Done umf_bld ',info
if(info /= 0) then
info=4010
ch_err='umf_bld'

@ -64,7 +64,7 @@ subroutine psb_zbldaggrmat(a,desc_a,ac,p,desc_p,info)
select case (p%iprcparm(smth_kind_))
case (no_smth_)
call raw_aggregate(info)
call raw_aggregate(info)
if(info /= 0) then
call psb_errpush(4010,name,a_err='raw_aggregate')
@ -120,6 +120,7 @@ contains
call psb_erractionsave(err_act)
bg => ac
call psb_nullify_sp(b)
icontxt = desc_a%matrix_data(psb_ctxt_)
call blacs_gridinfo(icontxt,nprows,npcols,myprow,mypcol)
@ -362,6 +363,9 @@ contains
call blacs_gridinfo(icontxt,nprows,npcols,myprow,mypcol)
bg => ac
call psb_nullify_sp(b)
call psb_nullify_sp(am3)
call psb_nullify_sp(am4)
am2 => p%av(sm_pr_t_)
am1 => p%av(sm_pr_)
@ -441,7 +445,15 @@ contains
! 1. Allocate Ptilde in sparse matrix form
call psb_sp_all(am4,ncol,info)
am4%fida='COO'
am4%m=ncol
if (ml_global_nmb) then
am4%k=ntaggr
call psb_sp_all(ncol,ntaggr,am4,ncol,info)
else
am4%k=naggr
call psb_sp_all(ncol,naggr,am4,ncol,info)
endif
if(info /= 0) then
call psb_errpush(4010,name,a_err='spall')
goto 9999
@ -462,14 +474,6 @@ contains
end do
am4%infoa(psb_nnz_) = nrow
endif
am4%fida='COO'
am4%m=ncol
if (ml_global_nmb) then
am4%k=ntaggr
else
am4%k=naggr
endif
if (test_dump) call &
@ -560,6 +564,7 @@ contains
if (test_dump) call psb_csprt(40+me,am3,head='% (I-wDA)',ivr=desc_a%loc_to_glob,&
& ivc=desc_a%loc_to_glob)
if (debug) write(0,*) me,'Done gather, going for SYMBMM 1'
!
! Symbmm90 does the allocation for its result.
!
@ -570,6 +575,7 @@ contains
call psb_symbmm(am3,am4,am1)
call psb_numbmm(am3,am4,am1)
if (debug) write(0,*) me,'Done NUMBMM 1'
call psb_sp_free(am4,info)
if(info /= 0) then
@ -615,6 +621,7 @@ contains
call psb_symbmm(a,am1,am3)
call psb_numbmm(a,am1,am3)
if (debug) write(0,*) me,'Done NUMBMM 2'
if (p%iprcparm(smth_kind_) == smth_omg_) then
call psb_transp(am1,am2,fmt='COO')
@ -638,6 +645,7 @@ contains
else
call psb_transp(am1,am2)
endif
if (debug) write(0,*) me,'starting sphalo/ rwxtd'
if (p%iprcparm(smth_kind_) == smth_omg_) then
! am2 = ((i-wDA)Ptilde)^T
@ -667,8 +675,11 @@ contains
end if
endif
if (debug) write(0,*) me,'starting symbmm 3'
call psb_symbmm(am2,am3,b)
if (debug) write(0,*) me,'starting numbmm 3'
call psb_numbmm(am2,am3,b)
if (debug) write(0,*) me,'Done NUMBMM 3'
!!$ if (aggr_dump) call csprt(50+me,am1,head='% Operator PTrans.')
call psb_sp_free(am3,info)
@ -731,6 +742,7 @@ contains
goto 9999
end if
if (debug) write(0,*) me,'Created aux descr. distr.'
call psb_cdasb(desc_p,info)
if(info /= 0) then
call psb_errpush(4010,name,a_err='psb_cdasb')
@ -738,6 +750,7 @@ contains
end if
if (debug) write(0,*) me,'Asmbld aux descr. distr.'
call psb_glob_to_loc(bg%ia1(1:nzl),desc_p,info,iact='I')
if(info /= 0) then

@ -156,9 +156,12 @@ subroutine psb_zilu_bld(a,desc_a,p,upd,info)
icontxt=desc_a%matrix_data(psb_ctxt_)
call psb_nullify_sp(blck)
call psb_nullify_sp(atmp)
t1= mpi_wtime()
if(debug) write(0,*)me,': calling psb_asmatbld',p%iprcparm(p_type_),p%iprcparm(n_ovr_)
if (debug) call blacs_barrier(icontxt,'All')
call psb_asmatbld(p%iprcparm(p_type_),p%iprcparm(n_ovr_),a,&
& blck,desc_a,upd,p%desc_data,info)
if(info/=0) then
@ -168,7 +171,8 @@ subroutine psb_zilu_bld(a,desc_a,p,upd,info)
goto 9999
end if
t2= mpi_wtime()
if(debug) write(0,*)me,': out of psb_asmatbld'
if (debug) write(0,*)me,': out of psb_asmatbld'
if (debug) call blacs_barrier(icontxt,'All')
if (associated(p%av)) then
if (size(p%av) < bp_ilu_avsz) then
@ -188,6 +192,9 @@ subroutine psb_zilu_bld(a,desc_a,p,upd,info)
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if (debug) write(0,*)me,': out spinfo',nztota
if (debug) call blacs_barrier(icontxt,'All')
n_col = desc_a%matrix_data(psb_n_col_)
nhalo = n_col-nrow_a
n_row = p%desc_data%matrix_data(psb_n_row_)
@ -197,7 +204,7 @@ subroutine psb_zilu_bld(a,desc_a,p,upd,info)
p%av(u_pr_)%m = n_row
p%av(u_pr_)%k = n_row
call psb_sp_all(n_row,n_row,p%av(l_pr_),nztota+lovr,info)
call psb_sp_all(n_row,n_row,p%av(u_pr_),nztota+lovr,info)
if (info == 0) call psb_sp_all(n_row,n_row,p%av(u_pr_),nztota+lovr,info)
if(info/=0) then
info=4010
ch_err='psb_sp_all'
@ -303,7 +310,8 @@ subroutine psb_zilu_bld(a,desc_a,p,upd,info)
endif
t5= mpi_wtime()
if (debug) write(0,*) me,' Going for dilu_fct'
if (debug) write(0,*) me,' Going for ilu_fct'
if (debug) call blacs_barrier(icontxt,'All')
call psb_ilu_fct(a,p%av(l_pr_),p%av(u_pr_),p%d,info,blck=blck)
if(info/=0) then
info=4010

@ -149,7 +149,7 @@ contains
if(psb_get_errstatus().ne.0) return
info=0
call psb_erractionsave(err_act)
call psb_nullify_sp(trw)
trw%m=0
trw%k=0
if(debug) write(0,*)'LUINT Allocating TRW'
@ -298,7 +298,7 @@ contains
!
info = 2
int_err(1) = i
write(ch_err,'(g20.10)'),abs(dia)
write(ch_err,'(g20.10)') abs(dia)
call psb_errpush(info,name,i_err=int_err,a_err=ch_err)
goto 9999
else
@ -437,7 +437,7 @@ contains
! Pivot too small: unstable factorization
!
int_err(1) = i
write(ch_err,'(g20.10)'),abs(dia)
write(ch_err,'(g20.10)') abs(dia)
info = 2
call psb_errpush(info,name,i_err=int_err,a_err=ch_err)
goto 9999

@ -150,6 +150,7 @@ subroutine psb_zmlprc_bld(a,desc_a,p,info)
goto 9999
end if
if (debug) write(0,*) 'Out from genaggrmap',p%nlaggr
nullify(desc_p)
allocate(desc_p)
call psb_nullify_desc(desc_p)
@ -165,6 +166,7 @@ subroutine psb_zmlprc_bld(a,desc_a,p,info)
call psb_baseprc_bld(ac,desc_p,p,info)
if (debug) write(0,*) 'Out from basaeprcbld',info
!
! We have used a separate ac because:

@ -68,7 +68,7 @@ subroutine psb_zslu_bld(a,desc_a,p,info)
integer, intent(out) :: info
character(len=5), optional :: outfmt
end Subroutine psb_zasmatbld
end interface
end interface
if(psb_get_errstatus().ne.0) return
info=0
@ -82,120 +82,121 @@ subroutine psb_zslu_bld(a,desc_a,p,info)
fmt = 'COO'
call psb_nullify_sp(blck)
call psb_nullify_sp(atmp)
atmp%fida='COO'
if (Debug) then
write(0,*) me, 'SPLUBLD: Calling csdp'
call blacs_barrier(icontxt,'All')
write(0,*) me, 'SPLUBLD: Calling csdp'
call blacs_barrier(icontxt,'All')
endif
call psb_zcsdp(a,atmp,info)
if(info /= 0) then
info=4010
ch_err='psb_zcsdp'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_zcsdp'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
nza = atmp%infoa(psb_nnz_)
if (Debug) then
write(0,*) me, 'SPLUBLD: Done csdp',info,nza,atmp%m,atmp%k
call blacs_barrier(icontxt,'All')
write(0,*) me, 'SPLUBLD: Done csdp',info,nza,atmp%m,atmp%k
call blacs_barrier(icontxt,'All')
endif
call psb_asmatbld(p%iprcparm(p_type_),p%iprcparm(n_ovr_),a,&
& blck,desc_a,upd,p%desc_data,info,outfmt=fmt)
if(info /= 0) then
info=4010
ch_err='psb_asmatbld'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_asmatbld'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
nzb = blck%infoa(psb_nnz_)
if (Debug) then
write(0,*) me, 'SPLUBLD: Done asmatbld',info,nzb,blck%fida
call blacs_barrier(icontxt,'All')
write(0,*) me, 'SPLUBLD: Done asmatbld',info,nzb,blck%fida
call blacs_barrier(icontxt,'All')
endif
if (nzb > 0 ) then
if (size(atmp%aspk)<nza+nzb) then
call psb_sp_reall(atmp,nza+nzb,info)
if(info /= 0) then
info=4010
ch_err='psb_sp_reall'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
endif
if (Debug) then
write(0,*) me, 'SPLUBLD: Done realloc',info,nza+nzb,atmp%fida
call blacs_barrier(icontxt,'All')
endif
do j=1,nzb
atmp%aspk(nza+j) = blck%aspk(j)
atmp%ia1(nza+j) = blck%ia1(j)
atmp%ia2(nza+j) = blck%ia2(j)
end do
atmp%infoa(psb_nnz_) = nza+nzb
atmp%m = atmp%m + blck%m
atmp%k = max(a%k,blck%k)
if (size(atmp%aspk)<nza+nzb) then
call psb_sp_reall(atmp,nza+nzb,info)
if(info /= 0) then
info=4010
ch_err='psb_sp_reall'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
endif
if (Debug) then
write(0,*) me, 'SPLUBLD: Done realloc',info,nza+nzb,atmp%fida
call blacs_barrier(icontxt,'All')
endif
do j=1,nzb
atmp%aspk(nza+j) = blck%aspk(j)
atmp%ia1(nza+j) = blck%ia1(j)
atmp%ia2(nza+j) = blck%ia2(j)
end do
atmp%infoa(psb_nnz_) = nza+nzb
atmp%m = atmp%m + blck%m
atmp%k = max(a%k,blck%k)
else
atmp%infoa(psb_nnz_) = nza
atmp%m = a%m
atmp%k = a%k
atmp%infoa(psb_nnz_) = nza
atmp%m = a%m
atmp%k = a%k
endif
i=0
do j=1, atmp%infoa(psb_nnz_)
if (atmp%ia2(j) <= atmp%m) then
i = i + 1
atmp%aspk(i) = atmp%aspk(j)
atmp%ia1(i) = atmp%ia1(j)
atmp%ia2(i) = atmp%ia2(j)
endif
if (atmp%ia2(j) <= atmp%m) then
i = i + 1
atmp%aspk(i) = atmp%aspk(j)
atmp%ia1(i) = atmp%ia1(j)
atmp%ia2(i) = atmp%ia2(j)
endif
enddo
atmp%infoa(psb_nnz_) = i
call psb_ipcoo2csr(atmp,info)
if(info /= 0) then
info=4010
ch_err='psb_ipcoo2csr'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_ipcoo2csr'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call psb_spinfo(psb_nztotreq_,atmp,nzt,info)
if(info /= 0) then
info=4010
ch_err='psb_spinfo'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_spinfo'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if (Debug) then
write(0,*) me,'Calling psb_slu_factor ',nzt,atmp%m,&
& atmp%k,p%desc_data%matrix_data(psb_n_row_)
call blacs_barrier(icontxt,'All')
write(0,*) me,'Calling psb_slu_factor ',nzt,atmp%m,&
& atmp%k,p%desc_data%matrix_data(psb_n_row_)
call blacs_barrier(icontxt,'All')
endif
call psb_zslu_factor(atmp%m,nzt,&
& atmp%aspk,atmp%ia2,atmp%ia1,p%iprcparm(slu_ptr_),info)
if(info /= 0) then
info=4010
ch_err='psb_slu_fact'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_slu_fact'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if (Debug) then
write(0,*) me, 'SPLUBLD: Done slu_Factor',info,p%iprcparm(slu_ptr_)
call blacs_barrier(icontxt,'All')
write(0,*) me, 'SPLUBLD: Done slu_Factor',info,p%iprcparm(slu_ptr_)
call blacs_barrier(icontxt,'All')
endif
call psb_sp_free(blck,info)
call psb_sp_free(atmp,info)
if(info /= 0) then
info=4010
ch_err='psb_sp_free'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_sp_free'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call psb_erractionrestore(err_act)
@ -204,8 +205,8 @@ subroutine psb_zslu_bld(a,desc_a,p,info)
9999 continue
call psb_erractionrestore(err_act)
if (err_act.eq.act_abort) then
call psb_error()
return
call psb_error()
return
end if
return

@ -69,7 +69,7 @@ subroutine psb_zumf_bld(a,desc_a,p,info)
integer, intent(out) :: info
character(len=5), optional :: outfmt
end Subroutine psb_zasmatbld
end interface
end interface
info=0
name='psb_umf_bld'
@ -82,97 +82,98 @@ subroutine psb_zumf_bld(a,desc_a,p,info)
fmt = 'COO'
call psb_nullify_sp(blck)
call psb_nullify_sp(atmp)
atmp%fida='COO'
if (Debug) then
write(0,*) me, 'UMFBLD: Calling csdp'
call blacs_barrier(icontxt,'All')
write(0,*) me, 'UMFBLD: Calling csdp'
call blacs_barrier(icontxt,'All')
endif
call psb_zcsdp(a,atmp,info)
if(info /= 0) then
info=4010
ch_err='psb_zcsdp'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_zcsdp'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call psb_spinfo(psb_nztotreq_,atmp,nza,info)
if (Debug) then
write(0,*) me, 'UMFBLD: Done csdp',info,nza,atmp%m,atmp%k
call blacs_barrier(icontxt,'All')
write(0,*) me, 'UMFBLD: Done csdp',info,nza,atmp%m,atmp%k
call blacs_barrier(icontxt,'All')
endif
call psb_asmatbld(p%iprcparm(p_type_),p%iprcparm(n_ovr_),a,&
& blck,desc_a,upd,p%desc_data,info,outfmt=fmt)
if(info /= 0) then
info=4010
ch_err='psb_asmatbld'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_asmatbld'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call psb_spinfo(psb_nztotreq_,blck,nzb,info)
if (Debug) then
write(0,*) me, 'UMFBLD: Done asmatbld',info,nzb,blck%fida
call blacs_barrier(icontxt,'All')
write(0,*) me, 'UMFBLD: Done asmatbld',info,nzb,blck%fida
call blacs_barrier(icontxt,'All')
endif
if (nzb > 0 ) then
if (size(atmp%aspk)<nza+nzb) then
call psb_sp_reall(atmp,nza+nzb,info)
if(info /= 0) then
info=4010
ch_err='psb_sp_reall'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
endif
if (Debug) then
write(0,*) me, 'UMFBLD: Done realloc',info,nza+nzb,atmp%fida
call blacs_barrier(icontxt,'All')
endif
do j=1,nzb
atmp%aspk(nza+j) = blck%aspk(j)
atmp%ia1(nza+j) = blck%ia1(j)
atmp%ia2(nza+j) = blck%ia2(j)
end do
atmp%infoa(psb_nnz_) = nza+nzb
atmp%m = atmp%m + blck%m
atmp%k = max(a%k,blck%k)
if (size(atmp%aspk)<nza+nzb) then
call psb_sp_reall(atmp,nza+nzb,info)
if(info /= 0) then
info=4010
ch_err='psb_sp_reall'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
endif
if (Debug) then
write(0,*) me, 'UMFBLD: Done realloc',info,nza+nzb,atmp%fida
call blacs_barrier(icontxt,'All')
endif
do j=1,nzb
atmp%aspk(nza+j) = blck%aspk(j)
atmp%ia1(nza+j) = blck%ia1(j)
atmp%ia2(nza+j) = blck%ia2(j)
end do
atmp%infoa(psb_nnz_) = nza+nzb
atmp%m = atmp%m + blck%m
atmp%k = max(a%k,blck%k)
else
atmp%infoa(psb_nnz_) = nza
atmp%m = a%m
atmp%k = a%k
atmp%infoa(psb_nnz_) = nza
atmp%m = a%m
atmp%k = a%k
endif
i=0
do j=1, atmp%infoa(psb_nnz_)
if (atmp%ia2(j) <= atmp%m) then
i = i + 1
atmp%aspk(i) = atmp%aspk(j)
atmp%ia1(i) = atmp%ia1(j)
atmp%ia2(i) = atmp%ia2(j)
endif
if (atmp%ia2(j) <= atmp%m) then
i = i + 1
atmp%aspk(i) = atmp%aspk(j)
atmp%ia1(i) = atmp%ia1(j)
atmp%ia2(i) = atmp%ia2(j)
endif
enddo
atmp%infoa(psb_nnz_) = i
call psb_ipcoo2csc(atmp,info)
if(info /= 0) then
info=4010
ch_err='psb_ipcoo2csc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_ipcoo2csc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call psb_spinfo(psb_nztotreq_,atmp,nzt,info)
if(info /= 0) then
info=4010
ch_err='psb_spinfo'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_spinfo'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if (Debug) then
write(0,*) me,'Calling psb_umf_factor ',nzt,atmp%m,&
& atmp%k,p%desc_data%matrix_data(psb_n_row_)
call blacs_barrier(icontxt,'All')
write(0,*) me,'Calling psb_umf_factor ',nzt,atmp%m,&
& atmp%k,p%desc_data%matrix_data(psb_n_row_)
call blacs_barrier(icontxt,'All')
endif
call psb_zumf_factor(atmp%m,nzt,&
@ -180,23 +181,23 @@ subroutine psb_zumf_bld(a,desc_a,p,info)
& p%iprcparm(umf_symptr_),p%iprcparm(umf_numptr_),info)
if(info /= 0) then
info=4010
ch_err='psb_umf_fact'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_umf_fact'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if (Debug) then
write(0,*) me, 'UMFBLD: Done umf_Factor',info,p%iprcparm(umf_numptr_)
call blacs_barrier(icontxt,'All')
write(0,*) me, 'UMFBLD: Done umf_Factor',info,p%iprcparm(umf_numptr_)
call blacs_barrier(icontxt,'All')
endif
call psb_sp_free(blck,info)
call psb_sp_free(atmp,info)
if(info /= 0) then
info=4010
ch_err='psb_sp_free'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
info=4010
ch_err='psb_sp_free'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
call psb_erractionrestore(err_act)
@ -205,8 +206,8 @@ subroutine psb_zumf_bld(a,desc_a,p,info)
9999 continue
call psb_erractionrestore(err_act)
if (err_act.eq.act_abort) then
call psb_error()
return
call psb_error()
return
end if
return

@ -96,7 +96,7 @@ subroutine psb_dspmm(alpha,a,x,beta,y,desc_a,info,&
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer, intent(out) :: info
real(kind(1.d0)), optional, pointer :: work(:)
real(kind(1.d0)), optional, target :: work(:)
character, intent(in), optional :: trans
integer, intent(in), optional :: k, jx, jy,doswap
@ -109,6 +109,7 @@ subroutine psb_dspmm(alpha,a,x,beta,y,desc_a,info,&
real(kind(1.d0)),pointer :: tmpx(:), xp(:,:), yp(:,:), iwork(:)
character :: itrans
character(len=20) :: name, ch_err
logical :: aliw
name='psb_dspmm'
if(psb_get_errstatus().ne.0) return
@ -187,26 +188,28 @@ 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(size(work).lt.liwork) then
call psb_realloc(liwork,work,info)
if (present(work)) then
if (size(work) >= liwork) then
aliw =.false.
else
aliw=.true.
endif
else
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
end if
iwork => work
else
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
end if
iwork => work
endif
iwork(1)=dzero
! checking for matrix correctness
@ -342,7 +345,7 @@ subroutine psb_dspmm(alpha,a,x,beta,y,desc_a,info,&
end if
if(.not.present(work)) deallocate(iwork)
if(aliw) deallocate(iwork)
nullify(iwork)
call psb_erractionrestore(err_act)
@ -433,7 +436,7 @@ subroutine psb_dspmv(alpha,a,x,beta,y,desc_a,info,&
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer, intent(out) :: info
real(kind(1.d0)), optional, pointer :: work(:)
real(kind(1.d0)), optional, target :: work(:)
character, intent(in), optional :: trans
integer, intent(in), optional :: doswap
@ -446,6 +449,7 @@ subroutine psb_dspmv(alpha,a,x,beta,y,desc_a,info,&
real(kind(1.d0)),pointer :: tmpx(:), iwork(:), xp(:), yp(:)
character :: itrans
character(len=20) :: name, ch_err
logical :: aliw
name='psb_dspmv'
if(psb_get_errstatus().ne.0) return
@ -505,16 +509,24 @@ subroutine psb_dspmv(alpha,a,x,beta,y,desc_a,info,&
lldx = size(x)
lldy = size(y)
iwork => null()
! check for presence/size of a work area
liwork= 2*ncol
if (a%pr(1) /= 0) liwork = liwork + n * ik
if (a%pl(1) /= 0) liwork = liwork + m * ik
! write(0,*)'---->>>',work(1)
if (present(work)) then
if(size(work).ge.liwork) then
iwork => work
liwork=size(work)
if (present(work)) then
if (size(work) >= liwork) then
aliw =.false.
else
aliw=.true.
endif
else
aliw=.true.
end if
aliw=.true.
if (aliw) then
call psb_realloc(liwork,iwork,info)
if(info.ne.0) then
info=4010
@ -522,16 +534,10 @@ subroutine psb_dspmv(alpha,a,x,beta,y,desc_a,info,&
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
end if
else
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
end if
iwork => work
endif
! checking for matrix correctness
call psb_chkmat(m,n,ia,ja,desc_a%matrix_data,info,iia,jja)
@ -644,9 +650,7 @@ subroutine psb_dspmv(alpha,a,x,beta,y,desc_a,info,&
end if
if(.not.present(work)) then
deallocate(iwork)
end if
if(aliw) deallocate(iwork)
nullify(iwork)
call psb_erractionrestore(err_act)

@ -91,7 +91,7 @@ subroutine psb_dspsm(alpha,a,x,beta,y,desc_a,info,&
type(psb_desc_type), intent(in) :: desc_a
integer, intent(out) :: info
real(kind(1.d0)), intent(in), optional, target :: d(:)
real(kind(1.d0)), optional, pointer :: work(:)
real(kind(1.d0)), optional, target :: work(:)
character, intent(in), optional :: trans, unitd
integer, intent(in), optional :: choice
integer, intent(in), optional :: k, jx, jy
@ -107,6 +107,7 @@ subroutine psb_dspsm(alpha,a,x,beta,y,desc_a,info,&
real(kind(1.d0)),pointer :: iwork(:), xp(:,:), yp(:,:), id(:)
character :: itrans
character(len=20) :: name, ch_err
logical :: aliw
name='psb_dspsm'
if(psb_get_errstatus().ne.0) return
@ -195,29 +196,32 @@ subroutine psb_dspsm(alpha,a,x,beta,y,desc_a,info,&
end if
! check for presence/size of a work area
iwork => null()
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).lt.liwork) then
call psb_realloc(liwork,work,info)
if (present(work)) then
if (size(work) >= liwork) then
aliw =.false.
else
aliw=.true.
endif
else
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
end if
iwork => work
else
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
end if
iwork => work
endif
iwork(1)=0.d0
if(present(d)) then
@ -302,7 +306,7 @@ subroutine psb_dspsm(alpha,a,x,beta,y,desc_a,info,&
end select
end if
if(.not.present(work)) deallocate(iwork)
if(aliw) deallocate(iwork)
if(.not.present(d)) deallocate(id)
call psb_erractionrestore(err_act)
@ -398,7 +402,7 @@ subroutine psb_dspsv(alpha,a,x,beta,y,desc_a,info,&
type(psb_desc_type), intent(in) :: desc_a
integer, intent(out) :: info
real(kind(1.d0)), intent(in), optional, target :: d(:)
real(kind(1.d0)), optional, pointer :: work(:)
real(kind(1.d0)), optional, target :: work(:)
character, intent(in), optional :: trans, unitd
integer, intent(in), optional :: choice
@ -413,6 +417,7 @@ subroutine psb_dspsv(alpha,a,x,beta,y,desc_a,info,&
real(kind(1.d0)),pointer :: iwork(:), xp(:), yp(:), id(:)
character :: itrans
character(len=20) :: name, ch_err
logical :: aliw
name='psb_dspsv'
if(psb_get_errstatus().ne.0) return
@ -484,30 +489,35 @@ subroutine psb_dspsv(alpha,a,x,beta,y,desc_a,info,&
goto 9999
end if
iwork => null()
! check for presence/size of a work area
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).lt.liwork) then
call psb_realloc(liwork,work,info)
if (size(work) >= liwork) then
aliw =.false.
else
aliw=.true.
endif
else
aliw=.true.
end if
aliw=.true.
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
end if
iwork => work
else
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
end if
iwork => work
endif
iwork(1)=0.d0
if (present(d)) then
@ -591,7 +601,7 @@ subroutine psb_dspsv(alpha,a,x,beta,y,desc_a,info,&
end select
end if
if(.not.present(work)) deallocate(iwork)
if (aliw) deallocate(iwork)
if(.not.present(d)) deallocate(id)
call psb_erractionrestore(err_act)

@ -96,6 +96,7 @@ Subroutine psb_dcdovrbld(n_ovr,desc_p,desc_a,a,&
icontxt = desc_a%matrix_data(psb_ctxt_)
!!$ call blacs_barrier(icontxt,'All')
Call blacs_gridinfo(icontxt,np,npcol,myrow,mycol)
call psb_nullify_sp(blk)
Allocate(brvindx(np+1),rvsz(np),sdsz(np),bsdindx(np+1),stat=info)
tl = 0.0

@ -50,6 +50,8 @@ zhb2mm: $(ZH2MOBJS)
zmm2hb: $(ZM2HOBJS)
$(MPF90) -o zmm2hb $(ZM2HOBJS) $(PSBLAS_LIB) $(BLACS)
srctst: srctst.o
$(MPF90) -o srctst srctst.o $(PSBLAS_LIB) $(BLACS)
.f90.o:
$(MPF90) $(F90COPT) $(INCDIRS) -c $<

@ -248,6 +248,15 @@ program df_sample
call psb_precset(pre,'asm',iv=(/novr,halo_,none_/))
case(rash_)
call psb_precset(pre,'asm',iv=(/novr,nohalo_,none_/))
case(7)
call psb_precset(pre,'asm',iv=(/ml,halo_,none_/))
call psb_precset(pre,'ml',&
& iv=(/mult_ml_prec_,loc_aggr_,smth_omg_,mat_distr_,post_smooth_,1,f_ilu_n_,4/))
case(8)
call psb_precset(pre,'asm',iv=(/ml,halo_,none_/))
call psb_precset(pre,'ml',&
& iv=(/mult_ml_prec_,loc_aggr_,smth_omg_,mat_distr_,post_smooth_,1,f_umf_,4/))
case default
call psb_precset(pre,'ilu')
end select

Loading…
Cancel
Save