!!$ !!$ 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_dcsdp.f90 ! ! Subroutine: psb_dcsdp ! This subroutine performs the assembly of ! the local part of a sparse distributed matrix ! ! Arguments: ! a - type(). The input matrix to be assembled. ! b - type(). The assembled output matrix. ! info - integer. Eventually returns an error code. ! ifc - integer(optional). ??? ! check - character(optional). ??? ! trans - character(optional). ??? ! unitd - character(optional). ??? ! subroutine psb_zspcnv2(a, b,info,afmt,upd,dupl) use psb_const_mod use psb_error_mod use psb_spmat_type use psb_string_mod use psb_serial_mod, psb_protect_name => psb_zspcnv2 implicit none !....Parameters... Type(psb_zspmat_type), intent(in) :: A Type(psb_zspmat_type), intent(out) :: B Integer, intent(out) :: info Integer, intent(in), optional :: upd,dupl character(len=*), optional, intent(in) :: afmt !...Locals... complex(kind(1.d0)) :: d(1) complex(kind(1.d0)), allocatable :: work(:) type(psb_zspmat_type) :: temp_a Integer :: nzr, ntry, ifc_, ia1_size,& & ia2_size, aspk_size,size_req,n_row,n_col,upd_,dupl_ integer :: ip1, ip2, nnz, iflag, ichk, nnzt,& & ipc, i, count, err_act, i1, i2, ia character :: check_,trans_,unitd_ character(len=5) :: afmt_ Integer, Parameter :: maxtry=8 logical, parameter :: debug=.false. character(len=20) :: name, ch_err name='psb_spcnv' info = 0 call psb_erractionsave(err_act) ntry=0 ifc_ = 2 check_ = 'N' trans_ = 'N' unitd_ = 'U' allocate(work(max(size(a%ia1),size(a%ia2))+a%m+1000),stat=info) if (info /= 0) then info=2040 call psb_errpush(info,name) goto 9999 end if if (present(upd)) then call psb_sp_setifld(upd,psb_upd_,b,info) end if if (present(dupl)) then call psb_sp_setifld(dupl,psb_dupl_,b,info) end if if (present(afmt)) then afmt_ = afmt else afmt_ = psb_fidef_ end if upd_ = psb_sp_getifld(psb_upd_,b,info) select case(upd_) case(psb_upd_srch_,psb_upd_perm_) ! Legal value, do nothing case default ! Fix bad value upd_ = psb_upd_dflt_ call psb_sp_setifld(upd_,psb_upd_,b,info) end select dupl_ = psb_sp_getifld(psb_dupl_,b,info) select case(dupl_) case(psb_dupl_ovwrt_,psb_dupl_add_,psb_dupl_err_) ! Legal value, do nothing case default ! Fix bad value dupl_ = psb_dupl_def_ call psb_sp_setifld(dupl_,psb_dupl_,b,info) end select ! ...matrix conversion... b%m=a%m b%k=a%k b%fida=afmt_ size_req = psb_sp_get_nnzeros(a) if (debug) write(0,*) 'DCSDP : size_req 1:',size_req ! n_row=b%m n_col=b%k call psb_cest(afmt_,n_row,n_col,size_req,& & ia1_size, ia2_size, aspk_size, upd_,info) b%fida=afmt_ if (info /= psb_no_err_) then info=4010 ch_err='psb_cest' call psb_errpush(info,name,a_err=ch_err) goto 9999 endif if (.not.allocated(b%aspk).or.& &.not.allocated(b%ia1).or.& &.not.allocated(b%ia2).or.& &.not.allocated(b%pl).or.& &.not.allocated(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 /= psb_no_err_) then info=4010 ch_err='psb_sp_reall' call psb_errpush(info,name,a_err=ch_err) goto 9999 endif b%pl(:) = 0 b%pr(:) = 0 b%descra = a%descra select case (tolower(a%fida)) case ('csr') select case (tolower(b%fida)) case ('csr') call zcrcr(trans_, a%m, a%k, unitd_, d, a%descra, a%aspk,& & a%ia1, a%ia2, a%infoa, b%pl, b%descra, b%aspk, b%ia1,& & b%ia2, b%infoa, b%pr, size(b%aspk), size(b%ia1),& & size(b%ia2), work, size(work), info) if (info/=0) then info=4010 ch_err='dcrcr' call psb_errpush(info,name,a_err=ch_err) goto 9999 end if case ('jad') !...converting to JAD !...output matrix may not be big enough do call zcrjd(trans_, a%m, a%k, unitd_, d, a%descra, a%aspk,& & a%ia1, a%ia2, a%infoa, b%pl, b%descra, b%aspk, b%ia1,& & b%ia2, b%infoa, b%pr, size(b%aspk), size(b%ia1),& & size(b%ia2), work, size(work), nzr, info) if (info /= 0) then call psb_errpush(4010,name,a_err='dcrjd') goto 9999 endif ntry = ntry + 1 if (debug) then write(0,*) 'On out from dcrjad ',nzr,info end if if (nzr == 0) exit if (ntry > maxtry ) then write(0,*) 'Tried reallocating for DCRJAD for ',maxtry,': giving up now.' info=2040 call psb_errpush(info,name) goto 9999 endif call psb_cest(afmt_,n_row,n_col,nzr,& & ia1_size, ia2_size, aspk_size, upd_,info) call psb_sp_reall(b,ia1_size,ia2_size,aspk_size,info) if (info /= 0) then info=2040 call psb_errpush(info,name) goto 9999 endif end do if (info/=0) then call psb_errpush(info,name) goto 9999 end if case ('coo') if (debug) write(0,*) 'Calling CRCO ',a%descra call zcrco(trans_, a%m, a%k, unitd_, d, a%descra, a%aspk,& & a%ia1, a%ia2, a%infoa, b%pl, b%descra, b%aspk, b%ia1,& & b%ia2, b%infoa, b%pr, size(b%aspk), size(b%ia1),& & size(b%ia2), work, size(work), info) if (info/=0) then call psb_errpush(4010,name,a_err='dcrco') goto 9999 end if case default info=4010 call psb_errpush(info,name) goto 9999 end select case ('coo','coi') select case (tolower(b%fida)) case ('csr') call zcocr(trans_, a%m, a%k, unitd_, d, a%descra, a%aspk,& & a%ia2, a%ia1, a%infoa, b%pl, b%descra, b%aspk, b%ia1,& & b%ia2, b%infoa, b%pr, size(b%aspk), size(b%ia1),& & size(b%ia2), work, 2*size(work), info) if (info/=0) then call psb_errpush(4010,name,a_err='dcocr') goto 9999 end if case ('jad') call psb_sp_all(temp_a, size(b%ia1),size(b%ia2),size(b%aspk),info) if (info /= 0) then info=2040 call psb_errpush(info,name) goto 9999 endif temp_a%m = a%m temp_a%k = a%k temp_a%infoa=b%infoa !...Dirty trick: converting to CSR and then to JAD call zcocr(trans_, a%m, a%k, unitd_, d, a%descra, a%aspk,& & a%ia2, a%ia1, a%infoa, temp_a%pl, temp_a%descra, & & temp_a%aspk, temp_a%ia1, temp_a%ia2, temp_a%infoa, temp_a%pr, & & size(temp_a%aspk), size(temp_a%ia1),& & size(temp_a%ia2), work, 2*size(work), info) if (info/=0) then call psb_errpush(4010,name,a_err='dcocr') goto 9999 end if do call zcrjd(trans_, temp_a%m, temp_a%k, unitd_, d, temp_a%descra, & & temp_a%aspk, temp_a%ia1, temp_a%ia2, temp_a%infoa, & & b%pl, b%descra, b%aspk, b%ia1, b%ia2, b%infoa, b%pr, & & size(b%aspk), size(b%ia1),& & size(b%ia2), work, size(work), nzr, info) if (info/=0) then call psb_errpush(4010,name,a_err='dcrjd') goto 9999 end if ntry = ntry + 1 if (debug) then write(0,*) 'On out from dcrjad ',nzr,info end if if (nzr == 0) exit if (ntry > maxtry ) then write(0,*) 'Tried reallocating for DCRJAD for ',maxtry,& & ': giving up now.' info=2040 call psb_errpush(info,name) goto 9999 endif call psb_sp_reall(b,nzr,info,ifc=ifc_) if (info /= 0) then info=2040 call psb_errpush(info,name) goto 9999 endif end do case ('coo') call zcoco(trans_, a%m, a%k, unitd_, d, a%descra, a%aspk,& & a%ia1, a%ia2, a%infoa, b%pl, b%descra, b%aspk, b%ia1,& & b%ia2, b%infoa, b%pr, size(b%aspk), size(b%ia1),& & size(b%ia2), work, 2*size(work), info) if (info/=0) then call psb_errpush(4010,name,a_err='dcoco') goto 9999 end if case default info=4010 call psb_errpush(info,name) goto 9999 end select case default info=4010 call psb_errpush(info,name) goto 9999 end select if (psb_sp_getifld(psb_upd_,b,info) /= psb_upd_perm_) & & call psb_sp_trim(b,info) call psb_sp_setifld(psb_spmat_asb_,psb_state_,b,info) call psb_erractionrestore(err_act) return 9999 continue call psb_erractionrestore(err_act) if (err_act.eq.psb_act_abort_) then call psb_error() return end if return end subroutine psb_zspcnv2 subroutine psb_zspcnv1(a, info, afmt, upd, dupl) use psb_spmat_type use psb_regen_mod use psb_serial_mod, psb_protect_name => psb_zspcnv1 use psb_const_mod use psi_mod use psb_error_mod use psb_string_mod !use psb_penv_mod implicit none !...Parameters.... type(psb_zspmat_type), intent (inout) :: a integer, intent(out) :: info integer,optional, intent(in) :: dupl, upd character(len=*), optional, intent(in) :: afmt integer :: int_err(5) type(psb_zspmat_type) :: atemp integer :: np,me,n_col,iout, err_act integer :: spstate integer :: upd_, dupl_ integer :: ictxt,n_row logical, parameter :: debug=.false., debugwrt=.false. character(len=20) :: name, ch_err info = 0 int_err(1)=0 name = 'psb_spcnv' call psb_erractionsave(err_act) if (present(upd)) then call psb_sp_setifld(upd,psb_upd_,a,info) end if if (present(dupl)) then call psb_sp_setifld(dupl,psb_dupl_,a,info) end if upd_ = psb_sp_getifld(psb_upd_,a,info) select case(upd_) case(psb_upd_srch_,psb_upd_perm_) ! Legal value, do nothing case default ! Fix bad value upd_ = psb_upd_dflt_ call psb_sp_setifld(upd_,psb_upd_,a,info) end select dupl_ = psb_sp_getifld(psb_dupl_,a,info) select case(dupl_) case(psb_dupl_ovwrt_,psb_dupl_add_,psb_dupl_err_) ! Legal value, do nothing case default ! Fix bad value dupl_ = psb_dupl_def_ call psb_sp_setifld(dupl_,psb_dupl_,a,info) end select spstate = psb_sp_getifld(psb_state_,a,info) if (info /= 0) then goto 9999 endif if (debug) write(0,*) 'Sparse matrix state:',spstate,psb_spmat_bld_,psb_spmat_upd_ if (spstate /= psb_spmat_upd_) then ! Should we first figure out if we can do it in place? if (debug) write(0,*) 'Update:',upd_,psb_upd_srch_,psb_upd_perm_ if (upd_ == psb_upd_srch_) then if (present(afmt)) then select case (tolower(a%fida)) case('coo') select case(tolower(afmt)) case('coo') call psb_fixcoo(a,info) goto 9998 case('csr') call psb_ipcoo2csr(a,info) goto 9998 case('csc') call psb_ipcoo2csc(a,info) goto 9998 end select case('csr') select case(tolower(afmt)) case('coo') call psb_ipcsr2coo(a,info) goto 9998 end select end select end if end if call psb_sp_clone(a,atemp,info) if(info /= psb_no_err_) then info=4010 ch_err='psb_sp_clone' call psb_errpush(info,name,a_err=ch_err) goto 9999 ! convert to user requested format after the temp copy end if if (debugwrt) then iout = 30+me open(iout) call psb_csprt(iout,atemp,head='Input mat') close(iout) endif ! Do the real conversion into the requested storage format ! result is put in A call psb_spcnv(atemp,a,info,afmt=afmt,upd=upd,dupl=dupl) IF (debug) WRITE (*, *) me,' ASB: From SPCNV',info,' ',A%FIDA if (info /= psb_no_err_) then info=4010 ch_err='psb_csdp' call psb_errpush(info,name,a_err=ch_err) goto 9999 endif if (debugwrt) then iout = 60+me open(iout) call psb_csprt(iout,a,head='Output mat') close(iout) endif call psb_sp_free(atemp,info) else if (spstate == psb_spmat_upd_) then ! ! Second case: we come from an update loop. ! select case(tolower(a%fida)) case('csr') call csr_regen(a,info) case ('coo','coi') call coo_regen(a,info) case('jad') call jad_regen(a,info) case default info=136 ch_err=a%fida call psb_errpush(info,name,a_err=ch_err) goto 9999 end select ! check on error retuned by dcsdp if (info /= psb_no_err_) then info = 4010 ch_err='xx_regen' call psb_errpush(info,name,a_err=ch_err) goto 9999 end if else info = 600 call psb_errpush(info,name) goto 9999 if (debug) write(0,*) 'Sparse matrix state:',spstate,psb_spmat_bld_,psb_spmat_upd_ endif 9998 continue call psb_sp_setifld(psb_spmat_asb_,psb_state_,a,info) call psb_erractionrestore(err_act) return 9999 continue call psb_erractionrestore(err_act) if (err_act.eq.psb_act_abort_) then call psb_error() return end if return end subroutine psb_zspcnv1