From a42b4b95aca5e429344a6ac0086ed2370f31445c Mon Sep 17 00:00:00 2001 From: Salvatore Filippone Date: Wed, 10 May 2006 15:47:38 +0000 Subject: [PATCH] *** empty log message *** --- src/modules/psb_realloc_mod.f90 | 91 ++++++++++++++++++++++++++++++--- src/tools/psb_dasb.f90 | 79 ++++++++++++---------------- src/tools/psb_iasb.f90 | 56 +++++++++----------- src/tools/psb_zasb.f90 | 45 +++++++--------- 4 files changed, 159 insertions(+), 112 deletions(-) diff --git a/src/modules/psb_realloc_mod.f90 b/src/modules/psb_realloc_mod.f90 index 5cdae29f..2802201d 100644 --- a/src/modules/psb_realloc_mod.f90 +++ b/src/modules/psb_realloc_mod.f90 @@ -38,6 +38,7 @@ module psb_realloc_mod module procedure psb_dreallocate2i1d module procedure psb_dreallocate1d module procedure psb_dreallocated2 + module procedure psb_dreallocatei2 module procedure psb_dreallocate2i1z module procedure psb_dreallocate1z module procedure psb_dreallocatez2 @@ -285,7 +286,7 @@ Contains ! ...Local Variables Real(kind(1.d0)),Pointer :: tmp(:,:) - Integer :: dim,err_act,err,i, m + Integer :: dim,err_act,err,i, m, dim2 character(len=20) :: name name='psb_dreallocated2' @@ -294,6 +295,7 @@ Contains if (associated(rrax)) then dim=size(rrax,1) + dim2=size(rrax,2) If (dim /= len1) Then Allocate(tmp(len1,len2),stat=info) @@ -306,10 +308,10 @@ Contains !!$ write(0,*) 'DA: copying ',min(len,dim) if (.true.) then do i=1,m - tmp(i,:) = rrax(i,:) + tmp(i,1:min(len2,dim2)) = rrax(i,1:min(len2,dim2)) end do else - tmp(1:m,:) = rrax(1:m,:) + tmp(1:m,1:min(len2,dim2)) = rrax(1:m,1:min(len2,dim2)) end if !!$ write(0,*) 'DA: copying done ',m Deallocate(rrax,stat=info) @@ -331,6 +333,7 @@ Contains endif if (present(pad)) then rrax(dim+1:len1,:) = pad + rrax(:,dim2+1:len2) = pad endif call psb_erractionrestore(err_act) return @@ -357,7 +360,7 @@ Contains ! ...Local Variables complex(kind(1.d0)),Pointer :: tmp(:,:) - Integer :: dim,err_act,err,i, m + Integer :: dim,err_act,err,i, m, dim2 character(len=20) :: name name='psb_dreallocatez2' @@ -366,6 +369,7 @@ Contains if (associated(rrax)) then dim=size(rrax,1) + dim2=size(rrax,2) If (dim /= len1) Then Allocate(tmp(len1,len2),stat=info) @@ -378,10 +382,10 @@ Contains !!$ write(0,*) 'DA: copying ',min(len,dim) if (.true.) then do i=1,m - tmp(i,:) = rrax(i,:) + tmp(i,1:min(len2,dim2)) = rrax(i,1:min(len2,dim2)) end do else - tmp(1:m,:) = rrax(1:m,:) + tmp(1:m,1:min(len2,dim2)) = rrax(1:m,1:min(len2,dim2)) end if !!$ write(0,*) 'DA: copying done ',m Deallocate(rrax,stat=info) @@ -403,6 +407,7 @@ Contains endif if (present(pad)) then rrax(dim+1:len1,:) = pad + rrax(:,dim2+1:len2) = pad endif call psb_erractionrestore(err_act) return @@ -419,6 +424,80 @@ Contains End Subroutine psb_dreallocatez2 + Subroutine psb_dreallocatei2(len1,len2,rrax,info,pad) + use psb_error_mod + ! ...Subroutine Arguments + Integer,Intent(in) :: len1,len2 + integer,pointer :: rrax(:,:) + integer :: info + integer, optional, intent(in) :: pad + + ! ...Local Variables + integer,Pointer :: tmp(:,:) + Integer :: dim,err_act,err,i, m, dim2 + character(len=20) :: name + + name='psb_dreallocatei2' + call psb_erractionsave(err_act) + info = 0 + + if (associated(rrax)) then + dim=size(rrax,1) + dim2=size(rrax,2) + + If (dim /= len1) Then + Allocate(tmp(len1,len2),stat=info) + if (info /= 0) then + err=4000 + call psb_errpush(err,name) + goto 9999 + end if + m = min(dim,len1) +!!$ write(0,*) 'DA: copying ',min(len,dim) + if (.true.) then + do i=1,m + tmp(i,1:min(len2,dim2)) = rrax(i,1:min(len2,dim2)) + end do + else + tmp(1:m,1:min(len2,dim2)) = rrax(1:m,1:min(len2,dim2)) + end if +!!$ write(0,*) 'DA: copying done ',m + Deallocate(rrax,stat=info) + if (info /= 0) then + err=4000 + call psb_errpush(err,name) + goto 9999 + end if + rrax=>tmp + End If + else + dim = 0 + Allocate(rrax(len1,len2),stat=info) + if (info /= 0) then + err=4000 + call psb_errpush(err,name) + goto 9999 + end if + endif + if (present(pad)) then + rrax(dim+1:len1,:) = pad + rrax(:,dim2+1:len2) = pad + endif + 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() + end if + return + + End Subroutine psb_dreallocatei2 + Subroutine psb_dreallocate2i(len,rrax,y,info,pad) use psb_error_mod diff --git a/src/tools/psb_dasb.f90 b/src/tools/psb_dasb.f90 index 7b359b23..63e27a63 100644 --- a/src/tools/psb_dasb.f90 +++ b/src/tools/psb_dasb.f90 @@ -41,8 +41,9 @@ subroutine psb_dasb(x, desc_a, info) !....assembly dense matrix x ..... use psb_descriptor_type use psb_const_mod - use psb_psblas_mod + use psb_comm_mod use psb_error_mod + use psb_realloc_mod implicit none type(psb_desc_type), intent(in) :: desc_a @@ -51,7 +52,6 @@ subroutine psb_dasb(x, desc_a, info) ! local variables integer :: err, icontxt,nprow,npcol,me,mypcol,temp,lwork,nrow,ncol, err_act - real(kind(1.d0)),pointer :: dtemp(:,:) integer :: int_err(5), i1sz, i2sz, dectype, i,j double precision :: real_err(5) real(kind(1.d0)),parameter :: one=1 @@ -62,18 +62,18 @@ subroutine psb_dasb(x, desc_a, info) info=0 name='psb_dasb' call psb_erractionsave(err_act) - + icontxt=desc_a%matrix_data(psb_ctxt_) dectype=desc_a%matrix_data(psb_dec_type_) - + call blacs_gridinfo(icontxt, nprow, npcol, me, mypcol) - + if ((.not.associated(desc_a%matrix_data))) then - info=3110 - call psb_errpush(info,name) - goto 9999 + info=3110 + call psb_errpush(info,name) + goto 9999 endif - + if (debug) write(*,*) 'asb start: ',nprow,npcol,me,& &desc_a%matrix_data(psb_dec_type_) ! ....verify blacs grid correctness.. @@ -93,7 +93,7 @@ subroutine psb_dasb(x, desc_a, info) call psb_errpush(info,name,i_err=int_err) goto 9999 endif - + ! check size icontxt=desc_a%matrix_data(psb_ctxt_) nrow=desc_a%matrix_data(psb_n_row_) @@ -102,22 +102,13 @@ subroutine psb_dasb(x, desc_a, info) i2sz = size(x,dim=2) if (debug) write(*,*) 'asb: ',i1sz,i2sz,nrow,ncol if (i1sz.lt.ncol) then - allocate(dtemp(ncol,i2sz),stat=info) - if (info.ne.0) then - info=2025 - int_err(1)=ncol - call psb_errpush(info,name,i_err=int_err) - goto 9999 - endif - - do j=1,size(x,2) - do i=1,nrow - dtemp(i,j) = x(i,j) - end do - end do - - deallocate(x) - x => dtemp + call psb_realloc(ncol,i2sz,x,info) + if (info.ne.0) then + info=2025 + int_err(1)=ncol + call psb_errpush(info,name,i_err=int_err) + goto 9999 + endif endif ! ..update halo elements.. @@ -186,8 +177,9 @@ subroutine psb_dasbv(x, desc_a, info) !....assembly dense matrix x ..... use psb_descriptor_type use psb_const_mod - use psb_psblas_mod + use psb_comm_mod use psb_error_mod + use psb_realloc_mod implicit none type(psb_desc_type), intent(in) :: desc_a @@ -197,7 +189,6 @@ subroutine psb_dasbv(x, desc_a, info) ! local variables integer :: err, icontxt,nprow,npcol,me,mypcol,temp,lwork integer :: int_err(5), i1sz,nrow,ncol, dectype, i, err_act - real(kind(1.d0)),pointer :: dtemp(:) double precision :: real_err(5) real(kind(1.d0)),parameter :: one=1 logical, parameter :: debug=.false. @@ -206,7 +197,7 @@ subroutine psb_dasbv(x, desc_a, info) info = 0 int_err(1) = 0 name = 'psb_dasbv' - + icontxt=desc_a%matrix_data(psb_ctxt_) dectype=desc_a%matrix_data(psb_dec_type_) @@ -214,18 +205,18 @@ subroutine psb_dasbv(x, desc_a, info) ! ....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,i_err=int_err) - goto 9999 + info = 2030 + int_err(1) = npcol + call psb_errpush(info,name,i_err=int_err) + goto 9999 else if (.not.psb_is_asb_dec(dectype)) then - info = 3110 - call psb_errpush(info,name) - goto 9999 + info = 3110 + call psb_errpush(info,name) + goto 9999 endif nrow=desc_a%matrix_data(psb_n_row_) @@ -234,20 +225,14 @@ subroutine psb_dasbv(x, desc_a, info) i1sz = size(x) if (debug) write(*,*) 'dasb: sizes ',i1sz,ncol if (i1sz.lt.ncol) then - allocate(dtemp(ncol),stat=info) + call psb_realloc(ncol,x,info) if (info.ne.0) then info=2025 int_err(1)=ncol call psb_errpush(info,name,i_err=int_err) goto 9999 endif - - do i=1,nrow - dtemp(i) = x(i) - end do - deallocate(x) - x => dtemp - endif + endif ! ..update halo elements.. call psb_halo(x,desc_a,info) diff --git a/src/tools/psb_iasb.f90 b/src/tools/psb_iasb.f90 index b0d026d7..c0593b21 100644 --- a/src/tools/psb_iasb.f90 +++ b/src/tools/psb_iasb.f90 @@ -41,8 +41,9 @@ subroutine psb_iasb(x, desc_a, info) !....assembly dense matrix x ..... use psb_descriptor_type use psb_const_mod - use psb_psblas_mod + use psb_comm_mod use psb_error_mod + use psb_realloc_mod implicit none type(psb_desc_type), intent(in) :: desc_a @@ -51,7 +52,6 @@ subroutine psb_iasb(x, desc_a, info) ! local variables integer :: icontxt,nprow,npcol,me,mypcol,temp,lwork,nrow,ncol,err_act - integer, pointer :: itemp(:,:) integer :: int_err(5), i1sz, i2sz, dectype, i real(kind(1.d0)) :: real_err(5) logical, parameter :: debug=.false. @@ -63,24 +63,24 @@ subroutine psb_iasb(x, desc_a, info) call psb_erractionsave(err_act) if ((.not.associated(desc_a%matrix_data))) then - info=3110 - call psb_errpush(info,name) - return + info=3110 + call psb_errpush(info,name) + return endif - + icontxt=desc_a%matrix_data(psb_ctxt_) dectype=desc_a%matrix_data(psb_dec_type_) - + call blacs_gridinfo(icontxt, nprow, npcol, me, mypcol) 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 ! check size @@ -91,17 +91,13 @@ subroutine psb_iasb(x, desc_a, info) i2sz = size(x,dim=2) if (debug) write(*,*) 'asb: ',i1sz,i2sz,nrow,ncol if (i1sz.lt.ncol) then - allocate(itemp(ncol,i2sz),stat=info) - if (info.ne.0) then - info=2025 - int_err(1)=ncol - call psb_errpush(info,name,int_err) - goto 9999 - endif - itemp(nrow+1:,:) = 0 - itemp(1:nrow,:) = x(1:nrow,:) - deallocate(x) - x => itemp + call psb_realloc(ncol,i2sz,x,info) + if (info.ne.0) then + info=2025 + int_err(1)=ncol + call psb_errpush(info,name,int_err) + goto 9999 + endif endif ! ..update halo elements.. @@ -163,8 +159,9 @@ subroutine psb_iasbv(x, desc_a, info) !....assembly dense matrix x ..... use psb_descriptor_type use psb_const_mod - use psb_psblas_mod + use psb_comm_mod use psb_error_mod + use psb_realloc_mod implicit none type(psb_desc_type), intent(in) :: desc_a @@ -174,7 +171,6 @@ subroutine psb_iasbv(x, desc_a, info) ! local variables integer :: icontxt,nprow,npcol,me,mypcol,temp,lwork, err_act integer :: int_err(5), i1sz,nrow,ncol, dectype, i - integer, pointer :: itemp(:) real(kind(1.d0)) :: real_err(5) logical, parameter :: debug=.false. character(len=20) :: name, ch_err @@ -206,17 +202,13 @@ subroutine psb_iasbv(x, desc_a, info) i1sz = size(x) if (debug) write(*,*) 'dasb: sizes ',i1sz,ncol if (i1sz.lt.ncol) then - allocate(itemp(ncol),stat=info) + call psb_realloc(ncol,x,info) if (info.ne.0) then info=2025 int_err(1)=ncol call psb_errpush(info,name,int_err) goto 9999 endif - itemp(nrow+1:) = 0 - itemp(1:nrow) = x(1:nrow) - deallocate(x) - x => itemp endif ! ..update halo elements.. diff --git a/src/tools/psb_zasb.f90 b/src/tools/psb_zasb.f90 index 6b3028cd..5c3a8618 100644 --- a/src/tools/psb_zasb.f90 +++ b/src/tools/psb_zasb.f90 @@ -41,8 +41,9 @@ subroutine psb_zasb(x, desc_a, info) !....assembly dense matrix x ..... use psb_descriptor_type use psb_const_mod - use psb_psblas_mod + use psb_comm_mod use psb_error_mod + use psb_realloc_mod implicit none type(psb_desc_type), intent(in) :: desc_a @@ -51,7 +52,6 @@ subroutine psb_zasb(x, desc_a, info) ! local variables integer :: err, icontxt,nprow,npcol,me,mypcol,temp,lwork,nrow,ncol, err_act - complex(kind(1.d0)),pointer :: ztemp(:,:) integer :: int_err(5), i1sz, i2sz, dectype, i,j double precision :: real_err(5) logical, parameter :: debug=.false. @@ -61,18 +61,18 @@ subroutine psb_zasb(x, desc_a, info) info=0 name='psb_zasb' call psb_erractionsave(err_act) - + icontxt=desc_a%matrix_data(psb_ctxt_) dectype=desc_a%matrix_data(psb_dec_type_) - + call blacs_gridinfo(icontxt, nprow, npcol, me, mypcol) - + if ((.not.associated(desc_a%matrix_data))) then - info=3110 - call psb_errpush(info,name) - goto 9999 + info=3110 + call psb_errpush(info,name) + goto 9999 endif - + if (debug) write(*,*) 'asb start: ',nprow,npcol,me,& &desc_a%matrix_data(psb_dec_type_) ! ....verify blacs grid correctness.. @@ -92,7 +92,7 @@ subroutine psb_zasb(x, desc_a, info) call psb_errpush(info,name,i_err=int_err) goto 9999 endif - + ! check size icontxt=desc_a%matrix_data(psb_ctxt_) nrow=desc_a%matrix_data(psb_n_row_) @@ -101,22 +101,13 @@ subroutine psb_zasb(x, desc_a, info) i2sz = size(x,dim=2) if (debug) write(*,*) 'asb: ',i1sz,i2sz,nrow,ncol if (i1sz.lt.ncol) then - allocate(ztemp(ncol,i2sz),stat=info) - if (info.ne.0) then - info=2025 - int_err(1)=ncol - call psb_errpush(info,name,i_err=int_err) - goto 9999 - endif - - do j=1,size(x,2) - do i=1,nrow - ztemp(i,j) = x(i,j) - end do - end do - - deallocate(x) - x => ztemp + call psb_realloc(ncol,i2sz,x,info) + if (info.ne.0) then + info=2025 + int_err(1)=ncol + call psb_errpush(info,name,i_err=int_err) + goto 9999 + endif endif ! ..update halo elements.. @@ -185,7 +176,7 @@ subroutine psb_zasbv(x, desc_a, info) !....assembly dense matrix x ..... use psb_descriptor_type use psb_const_mod - use psb_psblas_mod + use psb_comm_mod use psb_error_mod use psb_realloc_mod implicit none