From 05a5d8fa37015864f5f605cb4e85b666c80835b3 Mon Sep 17 00:00:00 2001 From: Salvatore Filippone Date: Thu, 12 Apr 2007 11:14:29 +0000 Subject: [PATCH] Fixed application of Block-Jacobi preconditioner, folding diagonal scale back into the serial sparse BLAS where it belongs. --- base/modules/psb_psblas_mod.f90 | 8 ++++---- base/psblas/psb_dspsm.f90 | 24 ++++++++++++------------ base/psblas/psb_zspsm.f90 | 24 ++++++++++++------------ prec/psb_dbjac_aply.f90 | 13 ++++--------- prec/psb_zbjac_aply.f90 | 12 ++++-------- 5 files changed, 36 insertions(+), 45 deletions(-) diff --git a/base/modules/psb_psblas_mod.f90 b/base/modules/psb_psblas_mod.f90 index cd8b5f7c..7d80ae90 100644 --- a/base/modules/psb_psblas_mod.f90 +++ b/base/modules/psb_psblas_mod.f90 @@ -417,7 +417,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(inout),target :: work(:), diag(:) + real(kind(1.d0)), optional, intent(in),target :: work(:), diag(:) integer, intent(out) :: info end subroutine psb_dspsm subroutine psb_dspsv(alpha, t, x, beta, y,& @@ -432,7 +432,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(inout),target :: work(:), diag(:) + real(kind(1.d0)), optional, intent(in),target :: work(:), diag(:) integer, intent(out) :: info end subroutine psb_dspsv subroutine psb_zspsm(alpha, t, x, beta, y,& @@ -448,7 +448,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(inout),target :: work(:), diag(:) + complex(kind(1.d0)), optional, intent(in),target :: work(:), diag(:) integer, intent(out) :: info end subroutine psb_zspsm subroutine psb_zspsv(alpha, t, x, beta, y,& @@ -463,7 +463,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(inout),target :: work(:), diag(:) + complex(kind(1.d0)), optional, intent(in),target :: work(:), diag(:) integer, intent(out) :: info end subroutine psb_zspsv end interface diff --git a/base/psblas/psb_dspsm.f90 b/base/psblas/psb_dspsm.f90 index aae113c0..9c21a974 100644 --- a/base/psblas/psb_dspsm.f90 +++ b/base/psblas/psb_dspsm.f90 @@ -72,7 +72,7 @@ ! work - real,dimension(:)(optional). Working area. ! 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_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_desc_type), intent(in) :: desc_a 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(:) character, intent(in), optional :: trans, unitd integer, intent(in), optional :: choice @@ -216,9 +216,9 @@ subroutine psb_dspsm(alpha,a,x,beta,y,desc_a,info,& iwork(1)=0.d0 - if(present(d)) then - lld = size(d) - id => d + if(present(diag)) then + lld = size(diag) + id => diag else lld=1 allocate(id(1)) @@ -301,7 +301,7 @@ subroutine psb_dspsm(alpha,a,x,beta,y,desc_a,info,& end if if(aliw) deallocate(iwork) - if(.not.present(d)) deallocate(id) + if(.not.present(diag)) deallocate(id) call psb_erractionrestore(err_act) return @@ -379,7 +379,7 @@ end subroutine psb_dspsm ! work - real,dimension(:)(optional). Working area. ! 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_serial_mod 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_desc_type), intent(in) :: desc_a 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(:) character, intent(in), optional :: trans, unitd integer, intent(in), optional :: choice @@ -510,9 +510,9 @@ subroutine psb_dspsv(alpha,a,x,beta,y,desc_a,info,& iwork(1)=0.d0 - if (present(d)) then - lld = size(d) - id => d + if (present(diag)) then + lld = size(diag) + id => diag else lld=1 allocate(id(1)) @@ -594,7 +594,7 @@ subroutine psb_dspsv(alpha,a,x,beta,y,desc_a,info,& end if if (aliw) deallocate(iwork) - if(.not.present(d)) deallocate(id) + if(.not.present(diag)) deallocate(id) call psb_erractionrestore(err_act) return diff --git a/base/psblas/psb_zspsm.f90 b/base/psblas/psb_zspsm.f90 index 6b0bd811..02319c46 100644 --- a/base/psblas/psb_zspsm.f90 +++ b/base/psblas/psb_zspsm.f90 @@ -72,7 +72,7 @@ ! work - real,dimension(:)(optional). Working area. ! 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_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_desc_type), intent(in) :: desc_a 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(:) character, intent(in), optional :: trans, unitd integer, intent(in), optional :: choice @@ -219,9 +219,9 @@ subroutine psb_zspsm(alpha,a,x,beta,y,desc_a,info,& iwork(1)=0.d0 - if(present(d)) then - lld = size(d) - id => d + if(present(diag)) then + lld = size(diag) + id => diag else lld=1 allocate(id(1)) @@ -304,7 +304,7 @@ subroutine psb_zspsm(alpha,a,x,beta,y,desc_a,info,& end if if(aliw) deallocate(iwork) - if(.not.present(d)) deallocate(id) + if(.not.present(diag)) deallocate(id) call psb_erractionrestore(err_act) return @@ -382,7 +382,7 @@ end subroutine psb_zspsm ! work - real,dimension(:)(optional). Working area. ! 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_serial_mod 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_desc_type), intent(in) :: desc_a 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(:) character, intent(in), optional :: trans, unitd integer, intent(in), optional :: choice @@ -508,9 +508,9 @@ subroutine psb_zspsv(alpha,a,x,beta,y,desc_a,info,& iwork(1)=0.d0 - if(present(d)) then - lld = size(d) - id => d + if(present(diag)) then + lld = size(diag) + id => diag else lld=1 allocate(id(1)) @@ -592,7 +592,7 @@ subroutine psb_zspsv(alpha,a,x,beta,y,desc_a,info,& end if if (aliw) deallocate(iwork) - if(.not.present(d)) deallocate(id) + if(.not.present(diag)) deallocate(id) call psb_erractionrestore(err_act) return diff --git a/prec/psb_dbjac_aply.f90 b/prec/psb_dbjac_aply.f90 index 845c7572..9ae5d292 100644 --- a/prec/psb_dbjac_aply.f90 +++ b/prec/psb_dbjac_aply.f90 @@ -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) call psb_info(ictxt, me, np) - diagl='U' - diagu='U' select case(trans) case('N','n') @@ -109,20 +107,18 @@ subroutine psb_dbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info) case('N','n') 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 - 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,& - & trans='N',unit=diagu,choice=psb_none_, work=aux) + & trans='N',unit='U',choice=psb_none_, work=aux) if(info /=0) goto 9999 case('T','t','C','c') 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 - 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,& - & trans=trans,unit=diagl,choice=psb_none_,work=aux) + & trans=trans,unit='U',choice=psb_none_,work=aux) if(info /=0) goto 9999 end select @@ -157,4 +153,3 @@ subroutine psb_dbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info) return end subroutine psb_dbjac_aply - diff --git a/prec/psb_zbjac_aply.f90 b/prec/psb_zbjac_aply.f90 index 25131ea1..767db45a 100644 --- a/prec/psb_zbjac_aply.f90 +++ b/prec/psb_zbjac_aply.f90 @@ -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) call psb_info(ictxt, me, np) - diagl='U' - diagu='U' select case(trans) case('N','n') @@ -109,20 +107,18 @@ subroutine psb_zbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info) case('N','n') 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 - 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,& - & trans='N',unit=diagu,choice=psb_none_, work=aux) + & trans='N',unit='U',choice=psb_none_, work=aux) if(info /=0) goto 9999 case('T','t','C','c') 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 - 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,& - & trans=trans,unit=diagl,choice=psb_none_,work=aux) + & trans=trans,unit='U',choice=psb_none_,work=aux) if(info /=0) goto 9999 end select