Fixed application of Block-Jacobi preconditioner, folding diagonal

scale back into the serial sparse BLAS where it belongs.
psblas3-type-indexed
Salvatore Filippone 18 years ago
parent fd22d34830
commit 05a5d8fa37

@ -417,7 +417,7 @@ module psb_psblas_mod
character, optional, intent(in) :: trans, unit character, optional, intent(in) :: trans, unit
integer, optional, intent(in) :: n, jx, jy integer, optional, intent(in) :: n, jx, jy
integer, optional, intent(in) :: choice integer, optional, intent(in) :: choice
real(kind(1.d0)), optional, intent(inout),target :: work(:), diag(:) real(kind(1.d0)), optional, intent(in),target :: work(:), diag(:)
integer, intent(out) :: info integer, intent(out) :: info
end subroutine psb_dspsm end subroutine psb_dspsm
subroutine psb_dspsv(alpha, t, x, beta, y,& subroutine psb_dspsv(alpha, t, x, beta, y,&
@ -432,7 +432,7 @@ module psb_psblas_mod
type(psb_desc_type), intent(in) :: desc_a type(psb_desc_type), intent(in) :: desc_a
character, optional, intent(in) :: trans, unit character, optional, intent(in) :: trans, unit
integer, optional, intent(in) :: choice integer, optional, intent(in) :: choice
real(kind(1.d0)), optional, intent(inout),target :: work(:), diag(:) real(kind(1.d0)), optional, intent(in),target :: work(:), diag(:)
integer, intent(out) :: info integer, intent(out) :: info
end subroutine psb_dspsv end subroutine psb_dspsv
subroutine psb_zspsm(alpha, t, x, beta, y,& subroutine psb_zspsm(alpha, t, x, beta, y,&
@ -448,7 +448,7 @@ module psb_psblas_mod
character, optional, intent(in) :: trans, unit character, optional, intent(in) :: trans, unit
integer, optional, intent(in) :: n, jx, jy integer, optional, intent(in) :: n, jx, jy
integer, optional, intent(in) :: choice integer, optional, intent(in) :: choice
complex(kind(1.d0)), optional, intent(inout),target :: work(:), diag(:) complex(kind(1.d0)), optional, intent(in),target :: work(:), diag(:)
integer, intent(out) :: info integer, intent(out) :: info
end subroutine psb_zspsm end subroutine psb_zspsm
subroutine psb_zspsv(alpha, t, x, beta, y,& subroutine psb_zspsv(alpha, t, x, beta, y,&
@ -463,7 +463,7 @@ module psb_psblas_mod
type(psb_desc_type), intent(in) :: desc_a type(psb_desc_type), intent(in) :: desc_a
character, optional, intent(in) :: trans, unit character, optional, intent(in) :: trans, unit
integer, optional, intent(in) :: choice integer, optional, intent(in) :: choice
complex(kind(1.d0)), optional, intent(inout),target :: work(:), diag(:) complex(kind(1.d0)), optional, intent(in),target :: work(:), diag(:)
integer, intent(out) :: info integer, intent(out) :: info
end subroutine psb_zspsv end subroutine psb_zspsv
end interface end interface

@ -72,7 +72,7 @@
! work - real,dimension(:)(optional). Working area. ! work - real,dimension(:)(optional). Working area.
! !
subroutine psb_dspsm(alpha,a,x,beta,y,desc_a,info,& subroutine psb_dspsm(alpha,a,x,beta,y,desc_a,info,&
& trans, unitd, choice, d, k, jx, jy, work) & trans, unitd, choice, diag, k, jx, jy, work)
use psb_spmat_type use psb_spmat_type
use psb_serial_mod use psb_serial_mod
@ -91,7 +91,7 @@ subroutine psb_dspsm(alpha,a,x,beta,y,desc_a,info,&
type (psb_dspmat_type), intent(in) :: a type (psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a type(psb_desc_type), intent(in) :: desc_a
integer, intent(out) :: info integer, intent(out) :: info
real(kind(1.d0)), intent(in), optional, target :: d(:) real(kind(1.d0)), intent(in), optional, target :: diag(:)
real(kind(1.d0)), optional, target :: work(:) real(kind(1.d0)), optional, target :: work(:)
character, intent(in), optional :: trans, unitd character, intent(in), optional :: trans, unitd
integer, intent(in), optional :: choice integer, intent(in), optional :: choice
@ -216,9 +216,9 @@ subroutine psb_dspsm(alpha,a,x,beta,y,desc_a,info,&
iwork(1)=0.d0 iwork(1)=0.d0
if(present(d)) then if(present(diag)) then
lld = size(d) lld = size(diag)
id => d id => diag
else else
lld=1 lld=1
allocate(id(1)) allocate(id(1))
@ -301,7 +301,7 @@ subroutine psb_dspsm(alpha,a,x,beta,y,desc_a,info,&
end if end if
if(aliw) deallocate(iwork) if(aliw) deallocate(iwork)
if(.not.present(d)) deallocate(id) if(.not.present(diag)) deallocate(id)
call psb_erractionrestore(err_act) call psb_erractionrestore(err_act)
return return
@ -379,7 +379,7 @@ end subroutine psb_dspsm
! work - real,dimension(:)(optional). Working area. ! work - real,dimension(:)(optional). Working area.
! !
subroutine psb_dspsv(alpha,a,x,beta,y,desc_a,info,& subroutine psb_dspsv(alpha,a,x,beta,y,desc_a,info,&
& trans, unitd, choice, d, work) & trans, unitd, choice, diag, work)
use psb_spmat_type use psb_spmat_type
use psb_serial_mod use psb_serial_mod
use psb_descriptor_type use psb_descriptor_type
@ -397,7 +397,7 @@ subroutine psb_dspsv(alpha,a,x,beta,y,desc_a,info,&
type(psb_dspmat_type), intent(in) :: a type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a type(psb_desc_type), intent(in) :: desc_a
integer, intent(out) :: info integer, intent(out) :: info
real(kind(1.d0)), intent(in), optional, target :: d(:) real(kind(1.d0)), intent(in), optional, target :: diag(:)
real(kind(1.d0)), optional, target :: work(:) real(kind(1.d0)), optional, target :: work(:)
character, intent(in), optional :: trans, unitd character, intent(in), optional :: trans, unitd
integer, intent(in), optional :: choice integer, intent(in), optional :: choice
@ -510,9 +510,9 @@ subroutine psb_dspsv(alpha,a,x,beta,y,desc_a,info,&
iwork(1)=0.d0 iwork(1)=0.d0
if (present(d)) then if (present(diag)) then
lld = size(d) lld = size(diag)
id => d id => diag
else else
lld=1 lld=1
allocate(id(1)) allocate(id(1))
@ -594,7 +594,7 @@ subroutine psb_dspsv(alpha,a,x,beta,y,desc_a,info,&
end if end if
if (aliw) deallocate(iwork) if (aliw) deallocate(iwork)
if(.not.present(d)) deallocate(id) if(.not.present(diag)) deallocate(id)
call psb_erractionrestore(err_act) call psb_erractionrestore(err_act)
return return

@ -72,7 +72,7 @@
! work - real,dimension(:)(optional). Working area. ! work - real,dimension(:)(optional). Working area.
! !
subroutine psb_zspsm(alpha,a,x,beta,y,desc_a,info,& subroutine psb_zspsm(alpha,a,x,beta,y,desc_a,info,&
& trans, unitd, choice, d, k, jx, jy, work) & trans, unitd, choice, diag, k, jx, jy, work)
use psb_spmat_type use psb_spmat_type
use psb_serial_mod use psb_serial_mod
@ -91,7 +91,7 @@ subroutine psb_zspsm(alpha,a,x,beta,y,desc_a,info,&
type (psb_zspmat_type), intent(in) :: a type (psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a type(psb_desc_type), intent(in) :: desc_a
integer, intent(out) :: info integer, intent(out) :: info
complex(kind(1.d0)), intent(in), optional, target :: d(:) complex(kind(1.d0)), intent(in), optional, target :: diag(:)
complex(kind(1.d0)), optional, target :: work(:) complex(kind(1.d0)), optional, target :: work(:)
character, intent(in), optional :: trans, unitd character, intent(in), optional :: trans, unitd
integer, intent(in), optional :: choice integer, intent(in), optional :: choice
@ -219,9 +219,9 @@ subroutine psb_zspsm(alpha,a,x,beta,y,desc_a,info,&
iwork(1)=0.d0 iwork(1)=0.d0
if(present(d)) then if(present(diag)) then
lld = size(d) lld = size(diag)
id => d id => diag
else else
lld=1 lld=1
allocate(id(1)) allocate(id(1))
@ -304,7 +304,7 @@ subroutine psb_zspsm(alpha,a,x,beta,y,desc_a,info,&
end if end if
if(aliw) deallocate(iwork) if(aliw) deallocate(iwork)
if(.not.present(d)) deallocate(id) if(.not.present(diag)) deallocate(id)
call psb_erractionrestore(err_act) call psb_erractionrestore(err_act)
return return
@ -382,7 +382,7 @@ end subroutine psb_zspsm
! work - real,dimension(:)(optional). Working area. ! work - real,dimension(:)(optional). Working area.
! !
subroutine psb_zspsv(alpha,a,x,beta,y,desc_a,info,& subroutine psb_zspsv(alpha,a,x,beta,y,desc_a,info,&
& trans, unitd, choice, d, work) & trans, unitd, choice, diag, work)
use psb_spmat_type use psb_spmat_type
use psb_serial_mod use psb_serial_mod
use psb_descriptor_type use psb_descriptor_type
@ -400,7 +400,7 @@ subroutine psb_zspsv(alpha,a,x,beta,y,desc_a,info,&
type(psb_zspmat_type), intent(in) :: a type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a type(psb_desc_type), intent(in) :: desc_a
integer, intent(out) :: info integer, intent(out) :: info
complex(kind(1.d0)), intent(in), optional, target :: d(:) complex(kind(1.d0)), intent(in), optional, target :: diag(:)
complex(kind(1.d0)), optional, target :: work(:) complex(kind(1.d0)), optional, target :: work(:)
character, intent(in), optional :: trans, unitd character, intent(in), optional :: trans, unitd
integer, intent(in), optional :: choice integer, intent(in), optional :: choice
@ -508,9 +508,9 @@ subroutine psb_zspsv(alpha,a,x,beta,y,desc_a,info,&
iwork(1)=0.d0 iwork(1)=0.d0
if(present(d)) then if(present(diag)) then
lld = size(d) lld = size(diag)
id => d id => diag
else else
lld=1 lld=1
allocate(id(1)) allocate(id(1))
@ -592,7 +592,7 @@ subroutine psb_zspsv(alpha,a,x,beta,y,desc_a,info,&
end if end if
if (aliw) deallocate(iwork) if (aliw) deallocate(iwork)
if(.not.present(d)) deallocate(id) if(.not.present(diag)) deallocate(id)
call psb_erractionrestore(err_act) call psb_erractionrestore(err_act)
return return

@ -64,8 +64,6 @@ subroutine psb_dbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
ictxt=psb_cd_get_context(desc_data) ictxt=psb_cd_get_context(desc_data)
call psb_info(ictxt, me, np) call psb_info(ictxt, me, np)
diagl='U'
diagu='U'
select case(trans) select case(trans)
case('N','n') case('N','n')
@ -109,20 +107,18 @@ subroutine psb_dbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
case('N','n') case('N','n')
call psb_spsm(done,prec%av(l_pr_),x,dzero,ww,desc_data,info,& call psb_spsm(done,prec%av(l_pr_),x,dzero,ww,desc_data,info,&
& trans='N',unit=diagl,choice=psb_none_,work=aux) & trans='N',unit='L',diag=prec%d,choice=psb_none_,work=aux)
if(info /=0) goto 9999 if(info /=0) goto 9999
ww(1:n_row) = ww(1:n_row)*prec%d(1:n_row)
call psb_spsm(alpha,prec%av(u_pr_),ww,beta,y,desc_data,info,& call psb_spsm(alpha,prec%av(u_pr_),ww,beta,y,desc_data,info,&
& trans='N',unit=diagu,choice=psb_none_, work=aux) & trans='N',unit='U',choice=psb_none_, work=aux)
if(info /=0) goto 9999 if(info /=0) goto 9999
case('T','t','C','c') case('T','t','C','c')
call psb_spsm(done,prec%av(u_pr_),x,dzero,ww,desc_data,info,& call psb_spsm(done,prec%av(u_pr_),x,dzero,ww,desc_data,info,&
& trans=trans,unit=diagu,choice=psb_none_, work=aux) & trans=trans,unit='L',diag=prec%d,choice=psb_none_, work=aux)
if(info /=0) goto 9999 if(info /=0) goto 9999
ww(1:n_row) = ww(1:n_row)*prec%d(1:n_row)
call psb_spsm(alpha,prec%av(l_pr_),ww,beta,y,desc_data,info,& call psb_spsm(alpha,prec%av(l_pr_),ww,beta,y,desc_data,info,&
& trans=trans,unit=diagl,choice=psb_none_,work=aux) & trans=trans,unit='U',choice=psb_none_,work=aux)
if(info /=0) goto 9999 if(info /=0) goto 9999
end select end select
@ -157,4 +153,3 @@ subroutine psb_dbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
return return
end subroutine psb_dbjac_aply end subroutine psb_dbjac_aply

@ -64,8 +64,6 @@ subroutine psb_zbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
ictxt=psb_cd_get_context(desc_data) ictxt=psb_cd_get_context(desc_data)
call psb_info(ictxt, me, np) call psb_info(ictxt, me, np)
diagl='U'
diagu='U'
select case(trans) select case(trans)
case('N','n') case('N','n')
@ -109,20 +107,18 @@ subroutine psb_zbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
case('N','n') case('N','n')
call psb_spsm(zone,prec%av(l_pr_),x,zzero,ww,desc_data,info,& call psb_spsm(zone,prec%av(l_pr_),x,zzero,ww,desc_data,info,&
& trans='N',unit=diagl,choice=psb_none_,work=aux) & trans='N',unit='L',diag=prec%d,choice=psb_none_,work=aux)
if(info /=0) goto 9999 if(info /=0) goto 9999
ww(1:n_row) = ww(1:n_row)*prec%d(1:n_row)
call psb_spsm(alpha,prec%av(u_pr_),ww,beta,y,desc_data,info,& call psb_spsm(alpha,prec%av(u_pr_),ww,beta,y,desc_data,info,&
& trans='N',unit=diagu,choice=psb_none_, work=aux) & trans='N',unit='U',choice=psb_none_, work=aux)
if(info /=0) goto 9999 if(info /=0) goto 9999
case('T','t','C','c') case('T','t','C','c')
call psb_spsm(zone,prec%av(u_pr_),x,zzero,ww,desc_data,info,& call psb_spsm(zone,prec%av(u_pr_),x,zzero,ww,desc_data,info,&
& trans=trans,unit=diagu,choice=psb_none_, work=aux) & trans=trans,unit='L',diag=prec%d,choice=psb_none_, work=aux)
if(info /=0) goto 9999 if(info /=0) goto 9999
ww(1:n_row) = ww(1:n_row)*prec%d(1:n_row)
call psb_spsm(alpha,prec%av(l_pr_),ww,beta,y,desc_data,info,& call psb_spsm(alpha,prec%av(l_pr_),ww,beta,y,desc_data,info,&
& trans=trans,unit=diagl,choice=psb_none_,work=aux) & trans=trans,unit='U',choice=psb_none_,work=aux)
if(info /=0) goto 9999 if(info /=0) goto 9999
end select end select

Loading…
Cancel
Save