fdurastante/amg4psblas
1
0
Fork 0
mirror of https://github.com/sfilippone/amg4psblas.git synced 2026-10-06 14:44:55 +00:00
Code Issues Packages Projects Releases Wiki Activity

Reorganized USE hierarchy.

Browse Source 
Reorganized naming scheme preliminary to applying MLD_ prefix.
This commit is contained in:
Salvatore Filippone
2007-02-21 14:32:23 +00:00
parent 8040237a77
commit 3173937349
41 changed files with 1069 additions and 1547 deletions
Download Patch File Download Diff File Expand all files Collapse all files
Makefile
+2 -2
View File
@@ -9,13 +9,13 @@ INCDIRS=-I. -I$(LIBDIR)
MODOBJS= psb_prec_type.o psb_prec_mod.o
MPFOBJS=psb_dbldaggrmat.o psb_zbldaggrmat.o
F90OBJS=psb_dasmatbld.o psb_dslu_bld.o psb_dumf_bld.o psb_dilu_fct.o\
psb_dmlprc_bld.o psb_dsp_renum.o psb_dilu_bld.o \
psb_dmlprc_bld.o psb_dsp_renum.o psb_dbjac_bld.o psb_dilu_bld.o \
psb_dprecbld.o psb_dprecfree.o psb_dprecset.o \
psb_dbaseprc_bld.o psb_ddiagsc_bld.o psb_dgenaggrmap.o \
psb_dprc_aply.o psb_dmlprc_aply.o \
psb_dbaseprc_aply.o psb_dbjac_aply.o\
psb_zasmatbld.o psb_zslu_bld.o psb_zumf_bld.o psb_zilu_fct.o\
psb_zmlprc_bld.o psb_zsp_renum.o psb_zilu_bld.o \
psb_zmlprc_bld.o psb_zsp_renum.o psb_zbjac_bld.o psb_zilu_bld.o \
psb_zprecbld.o psb_zprecfree.o psb_zprecset.o \
psb_zbaseprc_bld.o psb_zdiagsc_bld.o psb_zgenaggrmap.o \
psb_zprc_aply.o psb_zmlprc_aply.o \
psb_dasmatbld.f90
+2 -1
View File
@@ -54,7 +54,7 @@
Subroutine psb_dasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_dasmatbld
Implicit None
! .. Array Arguments ..
@@ -77,6 +77,7 @@ Subroutine psb_dasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
& tot_recv, ircode, n_row,nhalo, nrow_a,err_act
Logical,Parameter :: debug=.false., debugprt=.false.
character(len=20) :: name, ch_err
name='psb_dasmatbld'
if(psb_get_errstatus().ne.0) return
info=0
psb_dbaseprc_aply.f90
+2 -15
View File
@@ -41,7 +41,8 @@ subroutine psb_dbaseprc_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
!
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_dbaseprc_aply
implicit none
type(psb_desc_type),intent(in) :: desc_data
@@ -61,20 +62,6 @@ subroutine psb_dbaseprc_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
logical,parameter :: debug=.false., debugprt=.false.
character(len=20) :: name, ch_err
interface psb_bjac_aply
subroutine psb_dbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
use psb_base_mod
use psb_prec_type
type(psb_desc_type), intent(in) :: desc_data
type(psb_dbaseprc_type), intent(in) :: prec
real(kind(0.d0)),intent(inout) :: x(:), y(:)
real(kind(0.d0)),intent(in) :: alpha,beta
character(len=1) :: trans
real(kind(0.d0)),target :: work(:)
integer, intent(out) :: info
end subroutine psb_dbjac_aply
end interface
name='psb_dbaseprc_aply'
info = 0
call psb_erractionsave(err_act)
psb_dbaseprc_bld.f90
+9 -108
View File
@@ -37,7 +37,8 @@
subroutine psb_dbaseprc_bld(a,desc_a,p,info,upd)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_dbaseprc_bld
Implicit None
type(psb_dspmat_type), target :: a
@@ -46,56 +47,6 @@ subroutine psb_dbaseprc_bld(a,desc_a,p,info,upd)
integer, intent(out) :: info
character, intent(in), optional :: upd
interface psb_diagsc_bld
subroutine psb_ddiagsc_bld(a,desc_data,p,upd,info)
use psb_base_mod
use psb_prec_type
integer, intent(out) :: info
type(psb_dspmat_type), intent(in), target :: a
type(psb_desc_type),intent(in) :: desc_data
type(psb_dbaseprc_type), intent(inout) :: p
character, intent(in) :: upd
end subroutine psb_ddiagsc_bld
end interface
interface psb_ilu_bld
subroutine psb_dilu_bld(a,desc_data,p,upd,info)
use psb_base_mod
use psb_prec_type
integer, intent(out) :: info
type(psb_dspmat_type), intent(in), target :: a
type(psb_desc_type),intent(in) :: desc_data
type(psb_dbaseprc_type), intent(inout) :: p
character, intent(in) :: upd
end subroutine psb_dilu_bld
end interface
interface psb_slu_bld
subroutine psb_dslu_bld(a,desc_a,p,info)
use psb_base_mod
use psb_prec_type
implicit none
type(psb_dspmat_type), intent(inout) :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_dbaseprc_type), intent(inout) :: p
integer, intent(out) :: info
end subroutine psb_dslu_bld
end interface
interface psb_umf_bld
subroutine psb_dumf_bld(a,desc_a,p,info)
use psb_base_mod
use psb_prec_type
implicit none
type(psb_dspmat_type), intent(inout) :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_dbaseprc_type), intent(inout) :: p
integer, intent(out) :: info
end subroutine psb_dumf_bld
end interface
! Local scalars
Integer :: err, nnzero, n_row, n_col,I,j,k,ictxt,&
& me,mycol,np,npcol,mglob,lw, mtype, nrg, nzg, err_act
@@ -112,7 +63,7 @@ subroutine psb_dbaseprc_bld(a,desc_a,p,info,upd)
info=0
err=0
call psb_erractionsave(err_act)
name = 'psb_baseprc_bld'
name = 'psb_dbaseprc_bld'
if (debug) write(0,*) 'Entering baseprc_bld'
info = 0
@@ -139,17 +90,10 @@ subroutine psb_dbaseprc_bld(a,desc_a,p,info,upd)
!
! Should add check to ensure all procs have the same...
!
! ALso should define symbolic names for the preconditioners.
!
call psb_check_def(p%iprcparm(p_type_),'base_prec',&
& diag_,is_legal_base_prec)
!!$ allocate(p%desc_data,stat=info)
!!$ if (info /= 0) then
!!$ call psb_errpush(4010,name,a_err='Allocate')
!!$ goto 9999
!!$ end if
call psb_nullify_desc(p%desc_data)
@@ -175,7 +119,7 @@ subroutine psb_dbaseprc_bld(a,desc_a,p,info,upd)
goto 9999
end if
case (bjac_,asm_)
case(bjac_,asm_)
call psb_check_def(p%iprcparm(n_ovr_),'overlap',&
& 0,is_legal_n_ovr)
@@ -188,59 +132,16 @@ subroutine psb_dbaseprc_bld(a,desc_a,p,info,upd)
call psb_check_def(p%iprcparm(f_type_),'fact',&
& f_ilu_n_,is_legal_ml_fact)
if (debug) write(0,*)me, ': Calling PSB_ILU_BLD'
if (debug) write(0,*)me, ': Calling PSB_BJAC_BLD'
if (debug) call psb_barrier(ictxt)
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 psb_barrier(ictxt)
if(info /= 0) then
info=4010
ch_err='psb_ilu_bld'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
case(f_slu_)
if(debug) write(0,*)me,': calling slu_bld'
call psb_slu_bld(a,desc_a,p,info)
if(info /= 0) then
info=4010
ch_err='slu_bld'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
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'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
case(f_none_)
write(0,*) 'Fact=None in BASEPRC_BLD Bjac/ASM??'
call psb_bjac_bld(a,desc_a,p,iupd,info)
if(info /= 0) then
info=4010
ch_err='Inconsistent prec f_none_'
call psb_errpush(info,name,a_err=ch_err)
call psb_errpush(info,name,a_err='psb_bjac_bld')
goto 9999
end if
case default
write(0,*) 'Unknown factor type in baseprc_bld bjac/asm: ',&
&p%iprcparm(f_type_)
info=4010
ch_err='Unknown f_type_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end select
case default
info=4010
ch_err='Unknown p_type_'
psb_dbjac_aply.f90
+3 -2
View File
@@ -43,7 +43,8 @@ subroutine psb_dbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
!
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_dbjac_aply
implicit none
type(psb_desc_type), intent(in) :: desc_data
@@ -63,7 +64,7 @@ subroutine psb_dbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
logical,parameter :: debug=.false., debugprt=.false.
character(len=20) :: name, ch_err
name='psb_bjac_aply'
name='psb_dbjac_aply'
info = 0
call psb_erractionsave(err_act)
psb_dbjac_bld.f90
+368
View File
@@ -0,0 +1,368 @@
!!$
!!$
!!$ MD2P4
!!$ Multilevel Domain Decomposition Parallel Preconditioner Package for PSBLAS
!!$ for
!!$ Parallel Sparse BLAS v2.0
!!$ (C) Copyright 2006 Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$ Daniela di Serafino Second University of Naples
!!$ Pasqua D'Ambra ICAR-CNR
!!$
!!$ 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 MD2P4 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 MD2P4 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.
!!$
!!$
!*****************************************************************************
!* *
!* This is where the action takes place. *
!* ASMATBLD does the setup: building the prec descriptor plus retrieving *
!* matrix rows if needed *
!* *
!* *
!* *
!* *
!* some open code does the renumbering *
!* *
!* *
!* *
!* *
!*****************************************************************************
subroutine psb_dbjac_bld(a,desc_a,p,upd,info)
use psb_base_mod
use psb_prec_mod, mld_protect_name => psb_dbjac_bld
implicit none
!
! .. Scalar Arguments ..
integer, intent(out) :: info
! .. array Arguments ..
type(psb_dspmat_type), intent(in), target :: a
type(psb_dbaseprc_type), intent(inout) :: p
type(psb_desc_type), intent(in) :: desc_a
character, intent(in) :: upd
! .. Local Scalars ..
integer :: i, j, jj, k, kk, m
integer :: int_err(5)
character :: trans, unitd
type(psb_dspmat_type) :: blck, atmp
real(kind(1.d0)) :: t1,t2,t3,t4,t5,t6, t7, t8
logical, parameter :: debugprt=.false., debug=.false., aggr_dump=.false.
integer nztota, nztotb, nztmp, nzl, nnr, ir, err_act,&
& n_row, nrow_a,n_col, nhalo, ind, iind, i1,i2,ia
integer :: ictxt,np,me
character(len=20) :: name, ch_err
character(len=5), parameter :: coofmt='COO'
if(psb_get_errstatus().ne.0) return
info=0
name='psb_dbjac_bld'
call psb_erractionsave(err_act)
ictxt=psb_cd_get_context(desc_a)
call psb_info(ictxt, me, np)
m = a%m
if (m < 0) then
info = 10
int_err(1) = 1
int_err(2) = m
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
trans = 'N'
unitd = 'U'
if (p%iprcparm(n_ovr_) < 0) then
info = 11
int_err(1) = 1
int_err(2) = p%iprcparm(n_ovr_)
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
call psb_nullify_sp(blck)
call psb_nullify_sp(atmp)
t1= psb_wtime()
if(debug) write(0,*)me,': calling psb_asmatbld',p%iprcparm(p_type_),p%iprcparm(n_ovr_)
if (debug) call psb_barrier(ictxt)
call psb_asmatbld(p%iprcparm(p_type_),p%iprcparm(n_ovr_),a,&
& blck,desc_a,upd,p%desc_data,info,outfmt=coofmt)
if(info/=0) then
call psb_errpush(4010,name,a_err='psb_asmatbld')
goto 9999
end if
t2= psb_wtime()
if (debug) write(0,*)me,': out of psb_asmatbld'
if (debug) call psb_barrier(ictxt)
select case(p%iprcparm(iren_))
case (1:)
!
! Here we allocate a full copy to hold local A and received BLK
! Done inside sp_renum.
!
call psb_sp_renum(a,desc_a,blck,p,atmp,info)
if (info/=0) then
call psb_errpush(4010,name,a_err='psb_sp_renum')
goto 9999
end if
t3 = psb_wtime()
if (debugprt) then
call psb_barrier(ictxt)
open(40+me)
call psb_csprt(40+me,atmp,head='% Local matrix')
close(40+me)
endif
if (debug) write(0,*) me,' Factoring rows ',&
&atmp%m,a%m,blck%m,atmp%ia2(atmp%m+1)-1
select case(p%iprcparm(f_type_))
case(f_ilu_n_,f_ilu_e_)
call psb_ipcoo2csr(atmp,info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='psb_ipcoo2csr')
goto 9999
end if
call psb_ilu_bld(atmp,p%desc_data,p,upd,info)
if (info/=0) then
call psb_errpush(4010,name,a_err='psb_ilu_bld')
goto 9999
end if
if (debugprt) then
open(80+me)
call psb_csprt(80+me,p%av(l_pr_),head='% Local L factor')
write(80+me,*) '% Diagonal: ',p%av(l_pr_)%m
do i=1,p%av(l_pr_)%m
write(80+me,*) i,i,p%d(i)
enddo
call psb_csprt(80+me,p%av(u_pr_),head='% Local U factor')
close(80+me)
endif
case(f_slu_)
call psb_ipcoo2csr(atmp,info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='psb_ipcoo2csr')
goto 9999
end if
call psb_slu_bld(atmp,p%desc_data,p,info)
if(info /= 0) then
call psb_errpush(4010,name,a_err='slu_bld')
goto 9999
end if
case(f_umf_)
call psb_ipcoo2csc(atmp,info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='psb_ipcoo2csc')
goto 9999
end if
call psb_umf_bld(atmp,p%desc_data,p,info)
if(debug) write(0,*)me,': Done umf_bld ',info
if (info /= 0) then
call psb_errpush(4010,name,a_err='umf_bld')
goto 9999
end if
case(f_none_)
info=4010
call psb_errpush(info,name,a_err='Inconsistent prec f_none_')
goto 9999
case default
info=4010
call psb_errpush(info,name,a_err='Unknown f_type_')
goto 9999
end select
call psb_sp_free(atmp,info)
if(info/=0) then
call psb_errpush(4010,name,a_err='psb_sp_free')
goto 9999
end if
case(0) ! No renumbering
select case(p%iprcparm(f_type_))
case(f_ilu_n_,f_ilu_e_)
call psb_ipcoo2csr(blck,info,rwshr=.true.)
if(info/=0) then
call psb_errpush(4010,name,a_err='psb_ipcoo2csr')
goto 9999
end if
call psb_ilu_bld(a,desc_a,p,upd,info,blck=blck)
if(info/=0) then
call psb_errpush(4010,name,a_err='psb_ilu_bld')
goto 9999
end if
if (debugprt) then
open(80+me)
call psb_csprt(80+me,p%av(l_pr_),head='% Local L factor')
write(80+me,*) '% Diagonal: ',p%av(l_pr_)%m
do i=1,p%av(l_pr_)%m
write(80+me,*) i,i,p%d(i)
enddo
call psb_csprt(80+me,p%av(u_pr_),head='% Local U factor')
close(80+me)
endif
case(f_slu_)
atmp%fida='COO'
call psb_csdp(a,atmp,info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='psb_csdp')
goto 9999
end if
call psb_rwextd(atmp%m+blck%m,atmp,info,blck,rowscale=.true.)
if (info == 0) call psb_ipcoo2csr(atmp,info)
if (info == 0) call psb_slu_bld(atmp,p%desc_data,p,info)
if(info /= 0) then
call psb_errpush(4010,name,a_err='slu_bld')
goto 9999
end if
call psb_sp_free(atmp,info)
if(info/=0) then
call psb_errpush(4010,name,a_err='psb_sp_free')
goto 9999
end if
case(f_umf_)
atmp%fida='COO'
call psb_csdp(a,atmp,info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='psb_csdp')
goto 9999
end if
call psb_rwextd(atmp%m+blck%m,atmp,info,blck,rowscale=.true.)
if (info == 0) call psb_ipcoo2csc(atmp,info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='psb_ipcoo2csc')
goto 9999
end if
call psb_umf_bld(atmp,p%desc_data,p,info)
if(debug) write(0,*)me,': Done umf_bld ',info
if (info /= 0) then
call psb_errpush(4010,name,a_err='umf_bld')
goto 9999
end if
call psb_sp_free(atmp,info)
if(info/=0) then
call psb_errpush(4010,name,a_err='psb_sp_free')
goto 9999
end if
case(f_none_)
info=4010
call psb_errpush(info,name,a_err='Inconsistent prec f_none_')
goto 9999
case default
info=4010
call psb_errpush(info,name,a_err='Unknown f_type_')
goto 9999
end select
case default
info=4010
call psb_errpush(info,name,a_err='Invalid renum_')
goto 9999
end select
t6 = psb_wtime()
call psb_sp_free(blck,info)
if(info/=0) then
call psb_errpush(4010,name,a_err='psb_sp_free')
goto 9999
end if
if (debug) write(0,*) me,'End of ilu_bld'
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_dbjac_bld
psb_dbldaggrmat.f90
+3 -1
View File
@@ -36,7 +36,8 @@
!!$
subroutine psb_dbldaggrmat(a,desc_a,ac,desc_ac,p,info)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_dbldaggrmat
implicit none
type(psb_dspmat_type), intent(in), target :: a
@@ -49,6 +50,7 @@ subroutine psb_dbldaggrmat(a,desc_a,ac,desc_ac,p,info)
logical, parameter :: aggr_dump=.false.
integer ::ictxt,np,me, err_act
character(len=20) :: name, ch_err
name='psb_dbldaggrmat'
if(psb_get_errstatus().ne.0) return
info=0
psb_ddiagsc_bld.f90
+2 -2
View File
@@ -36,7 +36,7 @@
!!$
subroutine psb_ddiagsc_bld(a,desc_a,p,upd,info)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_ddiagsc_bld
Implicit None
type(psb_dspmat_type), target :: a
@@ -61,7 +61,7 @@ subroutine psb_ddiagsc_bld(a,desc_a,p,upd,info)
info=0
err=0
call psb_erractionsave(err_act)
name = 'psb_diagsc_bld'
name = 'psb_ddiagsc_bld'
if (debug) write(0,*) 'Entering diagsc_bld'
info = 0
psb_dgenaggrmap.f90
+3 -2
View File
@@ -36,7 +36,8 @@
!!$
subroutine psb_dgenaggrmap(aggr_type,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_dgenaggrmap
implicit none
integer, intent(in) :: aggr_type
type(psb_dspmat_type), intent(in) :: a
@@ -56,7 +57,7 @@ subroutine psb_dgenaggrmap(aggr_type,a,desc_a,nlaggr,ilaggr,info)
if(psb_get_errstatus().ne.0) return
info=0
name = 'psb_bldaggrmat'
name = 'psb_dgenaggrmap'
call psb_erractionsave(err_act)
!
! Note. At the time being we are ignoring aggr_type
psb_dilu_bld.f90
+45 -203
View File
@@ -49,141 +49,80 @@
!* *
!* *
!*****************************************************************************
subroutine psb_dilu_bld(a,desc_a,p,upd,info)
subroutine psb_dilu_bld(a,desc_a,p,upd,info,blck)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_dilu_bld
implicit none
!
! .. Scalar Arguments ..
integer, intent(out) :: info
! .. array Arguments ..
type(psb_dspmat_type), intent(in), target :: a
type(psb_dbaseprc_type), intent(inout) :: p
type(psb_desc_type), intent(in) :: desc_a
character, intent(in) :: upd
integer, intent(out) :: info
! .. array Arguments ..
type(psb_dspmat_type), intent(in), target :: a
type(psb_dbaseprc_type), intent(inout) :: p
type(psb_desc_type), intent(in) :: desc_a
character, intent(in) :: upd
type(psb_dspmat_type), intent(in), optional :: blck
! .. Local Scalars ..
integer :: i, j, jj, k, kk, m
integer :: int_err(5)
character :: trans, unitd
type(psb_dspmat_type) :: blck, atmp
real(kind(1.d0)) :: t1,t2,t3,t4,t5,t6, t7, t8
logical, parameter :: debugprt=.false., debug=.false., aggr_dump=.false.
integer nztota, nztotb, nztmp, nzl, nnr, ir, err_act,&
& n_row, nrow_a,n_col, nhalo, ind, iind, i1,i2,ia
& n_row, nrow_a, ind, iind, i1,i2,ia
integer :: ictxt,np,me
character(len=20) :: name, ch_err
interface psb_ilu_fct
subroutine psb_dilu_fct(a,l,u,d,info,blck)
use psb_base_mod
integer, intent(out) :: info
type(psb_dspmat_type),intent(in) :: a
type(psb_dspmat_type),intent(inout) :: l,u
type(psb_dspmat_type),intent(in), optional, target :: blck
real(kind(1.d0)), intent(inout) :: d(:)
end subroutine psb_dilu_fct
end interface
interface psb_asmatbld
Subroutine psb_dasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
use psb_base_mod
use psb_prec_type
integer, intent(in) :: ptype,novr
Type(psb_dspmat_type), Intent(in) :: a
Type(psb_dspmat_type), Intent(inout) :: blk
Type(psb_desc_type), Intent(inout) :: desc_p
Type(psb_desc_type), Intent(in) :: desc_data
Character, Intent(in) :: upd
integer, intent(out) :: info
character(len=5), optional :: outfmt
end Subroutine psb_dasmatbld
end interface
interface psb_sp_renum
subroutine psb_dsp_renum(a,desc_a,blck,p,atmp,info)
use psb_base_mod
use psb_prec_type
implicit none
type(psb_dspmat_type), intent(in) :: a,blck
type(psb_dspmat_type), intent(inout) :: atmp
type(psb_dbaseprc_type), intent(inout) :: p
type(psb_desc_type), intent(in) :: desc_a
integer, intent(out) :: info
end subroutine psb_dsp_renum
end interface
if(psb_get_errstatus().ne.0) return
info=0
name='psb_ilu_bld'
name='psb_dilu_bld'
call psb_erractionsave(err_act)
ictxt=psb_cd_get_context(desc_a)
call psb_info(ictxt, me, np)
m = a%m
if (m < 0) then
info = 10
int_err(1) = 1
int_err(2) = m
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
trans = 'N'
unitd = 'U'
if (p%iprcparm(n_ovr_) < 0) then
info = 11
int_err(1) = 1
int_err(2) = p%iprcparm(n_ovr_)
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
call psb_nullify_sp(blck)
call psb_nullify_sp(atmp)
t1= psb_wtime()
if(debug) write(0,*)me,': calling psb_asmatbld',p%iprcparm(p_type_),p%iprcparm(n_ovr_)
if (debug) call psb_barrier(ictxt)
call psb_asmatbld(p%iprcparm(p_type_),p%iprcparm(n_ovr_),a,&
& blck,desc_a,upd,p%desc_data,info)
if(info/=0) then
info=4010
ch_err='psb_asmatbld'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
t2= psb_wtime()
if (debug) write(0,*)me,': out of psb_asmatbld'
if (debug) call psb_barrier(ictxt)
if (allocated(p%av)) then
if (size(p%av) < bp_ilu_avsz) then
call psb_errpush(4010,name,a_err='Insufficient av size')
goto 9999
do i=1,size(p%av)
call psb_sp_free(p%av(i),info)
if (info /= 0) then
! Actually, we don't care here about this.
! Just let it go.
! return
end if
enddo
deallocate(p%av,stat=info)
endif
else
call psb_errpush(4010,name,a_err='AV not associated')
goto 9999
end if
if (.not.allocated(p%av)) then
allocate(p%av(max_avsz),stat=info)
if (info /= 0) then
call psb_errpush(4000,name)
goto 9999
end if
endif
!!$ call psb_csprt(50+me,a,head='% (A)')
nrow_a = psb_cd_get_local_rows(desc_a)
nztota = psb_sp_get_nnzeros(a)
nztotb = psb_sp_get_nnzeros(blck)
if (debug) write(0,*)me,': out get_nnzeros',nztota
if (present(blck)) then
nztota = nztota + psb_sp_get_nnzeros(blck)
end if
if (debug) write(0,*)me,': out get_nnzeros',nztota,a%m,a%k
if (debug) call psb_barrier(ictxt)
n_col = psb_cd_get_local_cols(desc_a)
nhalo = n_col-nrow_a
n_row = p%desc_data%matrix_data(psb_n_row_)
p%av(l_pr_)%m = n_row
p%av(l_pr_)%k = n_row
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+nztotb,info)
if (info == 0) call psb_sp_all(n_row,n_row,p%av(u_pr_),nztota+nztotb,info)
call psb_sp_all(n_row,n_row,p%av(l_pr_),nztota,info)
if (info == 0) call psb_sp_all(n_row,n_row,p%av(u_pr_),nztota,info)
if(info/=0) then
info=4010
ch_err='psb_sp_all'
@@ -206,100 +145,17 @@ subroutine psb_dilu_bld(a,desc_a,p,upd,info)
endif
if (debug) then
write(0,*) me,'Done psb_asmatbld'
call psb_barrier(ictxt)
endif
if (p%iprcparm(iren_) > 0) then
!
! Here we allocate a full copy to hold local A and received BLK
!
nztota = psb_sp_get_nnzeros(a)
nztotb = psb_sp_get_nnzeros(blck)
call psb_sp_all(atmp,nztota+nztotb,info)
if(info/=0) then
info=4011
call psb_errpush(info,name)
goto 9999
end if
! write(0,*) 'ILU_BLD ',nztota,nztotb,a%m
call psb_sp_renum(a,desc_a,blck,p,atmp,info)
if(info/=0) then
info=4010
ch_err='psb_sp_renum'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
t3 = psb_wtime()
if (debugprt) then
call psb_barrier(ictxt)
open(40+me)
call psb_csprt(40+me,atmp,head='% Local matrix')
close(40+me)
endif
if (debug) write(0,*) me,' Factoring rows ',&
&atmp%m,a%m,blck%m,atmp%ia2(atmp%m+1)-1
!
! Ok, factor the matrix.
!
t5 = psb_wtime()
blck%m=0
call psb_ilu_fct(atmp,p%av(l_pr_),p%av(u_pr_),p%d,info,blck=blck)
if(info/=0) then
info=4010
ch_err='psb_ilu_fct'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
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
end if
else if (p%iprcparm(iren_) == 0) then
t3 = psb_wtime()
! This is where we have mo renumbering, thus no need
! for ATMP
if (debugprt) then
open(40+me)
call psb_barrier(ictxt)
call psb_csprt(40+me,a,iv=p%desc_data%loc_to_glob,&
& head='% Local matrix')
if (p%iprcparm(p_type_)==asm_) then
call psb_csprt(40+me,blck,iv=p%desc_data%loc_to_glob,&
& irs=a%m,head='% Received rows')
endif
close(40+me)
endif
t5= psb_wtime()
if (debug) write(0,*) me,' Going for ilu_fct'
if (debug) call psb_barrier(ictxt)
call psb_ilu_fct(a,p%av(l_pr_),p%av(u_pr_),p%d,info,blck=blck)
if(info/=0) then
info=4010
ch_err='psb_ilu_fct'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if (debug) write(0,*) me,' Done dilu_fct'
endif
!
! Ok, factor the matrix.
!
t5 = psb_wtime()
call psb_ilu_fct(a,p%av(l_pr_),p%av(u_pr_),p%d,info,blck=blck)
if(info/=0) then
info=4010
ch_err='psb_ilu_fct'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if (debugprt) then
@@ -327,28 +183,14 @@ subroutine psb_dilu_bld(a,desc_a,p,upd,info)
! write(0,'(i3,1x,a,3(1x,g18.9))') me,'renum/factor time',t3-t2,t6-t5
! if (me==0) write(0,'(a,3(1x,g18.9))') 'renum/factor time',t3-t2,t6-t5
call psb_sp_free(blck,info)
if(info/=0) then
info=4010
ch_err='psb_sp_free'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if (psb_sp_getifld(psb_upd_,p%av(u_pr_),info) /= psb_upd_perm_) then
call psb_sp_trimsize(p%av(u_pr_),i1,i2,ia,info)
if (info == 0) call psb_sp_reall(p%av(u_pr_),i1,i2,ia,info)
if (info /=0) then
write(0,*) 'Error from trimsize 1',info
endif
endif
if (psb_sp_getifld(psb_upd_,p%av(l_pr_),info) /= psb_upd_perm_) then
call psb_sp_trimsize(p%av(l_pr_),i1,i2,ia,info)
if (info == 0) call psb_sp_reall(p%av(l_pr_),i1,i2,ia,info)
if (info /=0) then
write(0,*) 'Error from trimsize 2',info
endif
endif
psb_dilu_fct.f90
+2 -1
View File
@@ -57,7 +57,8 @@ subroutine psb_dilu_fct(a,l,u,d,info,blck)
type(psb_dspmat_type), pointer :: blck_
character(len=20) :: name, ch_err
logical, parameter :: debug=.false.
name='psb_dcsrlu'
name='psb_dilu_fct'
info = 0
call psb_erractionsave(err_act)
! .. Executable Statements ..
psb_dmlprc_aply.f90
+2 -16
View File
@@ -84,7 +84,7 @@ subroutine psb_dmlprc_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info)
!
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_dmlprc_aply
implicit none
type(psb_desc_type),intent(in) :: desc_data
@@ -112,21 +112,7 @@ subroutine psb_dmlprc_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info)
end type psb_mlprec_wrk_type
type(psb_mlprec_wrk_type), allocatable :: mlprec_wrk(:)
interface psb_baseprc_aply
subroutine psb_dbaseprc_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
use psb_base_mod
use psb_prec_type
type(psb_desc_type),intent(in) :: desc_data
type(psb_dbaseprc_type), intent(in) :: prec
real(kind(0.d0)),intent(inout) :: x(:), y(:)
real(kind(0.d0)),intent(in) :: alpha,beta
character(len=1) :: trans
real(kind(0.d0)),target :: work(:)
integer, intent(out) :: info
end subroutine psb_dbaseprc_aply
end interface
name='psb_mlprc_aply'
name='psb_dmlprc_aply'
info = 0
call psb_erractionsave(err_act)
psb_dmlprc_bld.f90
+2 -41
View File
@@ -37,7 +37,7 @@
subroutine psb_dmlprc_bld(a,desc_a,p,info)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_dmlprc_bld
implicit none
type(psb_dspmat_type), intent(in), target :: a
@@ -51,48 +51,9 @@ subroutine psb_dmlprc_bld(a,desc_a,p,info)
character(len=20) :: name, ch_err
logical, parameter :: debug=.false.
type(psb_dspmat_type) :: ac
interface psb_baseprc_bld
subroutine psb_dbaseprc_bld(a,desc_a,p,info,upd)
use psb_base_mod
use psb_prec_type
type(psb_dspmat_type), target :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_dbaseprc_type),intent(inout) :: p
integer, intent(out) :: info
character, intent(in), optional :: upd
end subroutine psb_dbaseprc_bld
end interface
interface psb_genaggrmap
subroutine psb_dgenaggrmap(aggr_type,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod
use psb_prec_type
implicit none
integer, intent(in) :: aggr_type
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer, allocatable :: ilaggr(:),nlaggr(:)
integer, intent(out) :: info
end subroutine psb_dgenaggrmap
end interface
interface psb_bldaggrmat
subroutine psb_dbldaggrmat(a,desc_a,ac,desc_ac,p,info)
use psb_base_mod
use psb_prec_type
type(psb_dspmat_type), intent(in), target :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_dspmat_type), intent(out),target :: ac
type(psb_desc_type), intent(inout) :: desc_ac
type(psb_dbaseprc_type), intent(inout), target :: p
integer, intent(out) :: info
end subroutine psb_dbldaggrmat
end interface
integer :: ictxt, np, me
name='psb_mlprec_bld'
name='psb_dmlprec_bld'
if (psb_get_errstatus().ne.0) return
info = 0
ictxt = psb_cd_get_context(desc_a)
psb_dprc_aply.f90
+7 -47
View File
@@ -37,7 +37,8 @@
subroutine psb_dprc_aply(prec,x,y,desc_data,info,trans, work)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_dprc_aply
implicit none
type(psb_desc_type),intent(in) :: desc_data
@@ -54,35 +55,7 @@ subroutine psb_dprc_aply(prec,x,y,desc_data,info,trans, work)
logical,parameter :: debug=.false., debugprt=.false.
character(len=20) :: name
interface psb_baseprc_aply
subroutine psb_dbaseprc_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
use psb_base_mod
use psb_prec_type
type(psb_desc_type),intent(in) :: desc_data
type(psb_dbaseprc_type), intent(in) :: prec
real(kind(0.d0)),intent(inout) :: x(:), y(:)
real(kind(0.d0)),intent(in) :: alpha,beta
character(len=1) :: trans
real(kind(0.d0)),target :: work(:)
integer, intent(out) :: info
end subroutine psb_dbaseprc_aply
end interface
interface psb_mlprc_aply
subroutine psb_dmlprc_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info)
use psb_base_mod
use psb_prec_type
type(psb_desc_type),intent(in) :: desc_data
type(psb_dbaseprc_type), intent(in) :: baseprecv(:)
real(kind(0.d0)),intent(in) :: alpha,beta
real(kind(0.d0)),intent(inout) :: x(:), y(:)
character :: trans
real(kind(0.d0)),target :: work(:)
integer, intent(out) :: info
end subroutine psb_dmlprc_aply
end interface
name='psb_dprc_aply'
name='psb_dprec_aply'
info = 0
call psb_erractionsave(err_act)
@@ -181,7 +154,8 @@ end subroutine psb_dprc_aply
subroutine psb_dprc_aply1(prec,x,desc_data,info,trans)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_dprc_aply1
implicit none
type(psb_desc_type),intent(in) :: desc_data
@@ -191,27 +165,13 @@ subroutine psb_dprc_aply1(prec,x,desc_data,info,trans)
character(len=1), optional :: trans
logical,parameter :: debug=.false., debugprt=.false.
interface
subroutine psb_dprc_aply(prec,x,y,desc_data,info,trans, work)
use psb_base_mod
use psb_prec_type
implicit none
type(psb_desc_type),intent(in) :: desc_data
type(psb_dprec_type), intent(in) :: prec
real(kind(0.d0)),intent(inout) :: x(:), y(:)
integer, intent(out) :: info
character(len=1), optional :: trans
real(kind(0.d0)), optional, target :: work(:)
end subroutine psb_dprc_aply
end interface
! Local variables
character :: trans_
integer :: ictxt,np,me,i, err_act
real(kind(1.d0)), pointer :: WW(:), w1(:)
character(len=20) :: name
name='psb_dprec1'
name='psb_dprec_aply1'
info = 0
call psb_erractionsave(err_act)
psb_dprecbld.f90
+3 -3
View File
@@ -1,3 +1,4 @@
!!$
!!$
!!$ MD2P4
@@ -37,8 +38,7 @@
subroutine psb_dprecbld(a,desc_a,p,info,upd)
use psb_base_mod
use psb_prec_type
use psb_prec_mod
use psb_prec_mod, protect => psb_dprecbld
Implicit None
type(psb_dspmat_type), target :: a
@@ -60,7 +60,7 @@ subroutine psb_dprecbld(a,desc_a,p,info,upd)
info=0
err=0
call psb_erractionsave(err_act)
name = 'psb_precbld'
name = 'psb_dprecbld'
if (debug) write(0,*) 'Entering precbld',P%prec,desc_a%matrix_data(:)
info = 0
psb_dprecfree.f90
+3 -2
View File
@@ -36,7 +36,8 @@
!!$
subroutine psb_dprecfree(p,info)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_dprecfree
implicit none
type(psb_dprec_type), intent(inout) :: p
integer, intent(out) :: info
@@ -47,7 +48,7 @@ subroutine psb_dprecfree(p,info)
if(psb_get_errstatus().ne.0) return
info=0
name = 'psdprecfree'
name = 'psb_dprecfree'
call psb_erractionsave(err_act)
me=-1
psb_dprecset.f90
+2 -1
View File
@@ -37,7 +37,8 @@
subroutine psb_dprecset(p,ptype,info,iv,rs,rv,ilev,nlev)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_dprecset
implicit none
type(psb_dprec_type), intent(inout) :: p
character(len=*), intent(in) :: ptype
psb_dslu_bld.f90
+13 -118
View File
@@ -36,142 +36,46 @@
!!$
subroutine psb_dslu_bld(a,desc_a,p,info)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_dslu_bld
implicit none
type(psb_dspmat_type), intent(inout) :: a
type(psb_dspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_dbaseprc_type), intent(inout) :: p
integer, intent(out) :: info
type(psb_dspmat_type) :: blck, atmp
character(len=5) :: fmt
character :: upd='F'
integer :: i,j,nza,nzb,nzt,ictxt,me,np,err_act
logical, parameter :: debug=.false.
character(len=20) :: name, ch_err
interface psb_asmatbld
Subroutine psb_dasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
use psb_base_mod
use psb_prec_type
integer, intent(in) :: ptype,novr
Type(psb_dspmat_type), Intent(in) :: a
Type(psb_dspmat_type), Intent(inout) :: blk
Type(psb_desc_type), Intent(inout) :: desc_p
Type(psb_desc_type), Intent(in) :: desc_data
Character, Intent(in) :: upd
integer, intent(out) :: info
character(len=5), optional :: outfmt
end Subroutine psb_dasmatbld
end interface
if(psb_get_errstatus().ne.0) return
info=0
name='psb_slu_bld'
name='psb_dslu_bld'
call psb_erractionsave(err_act)
ictxt = desc_a%matrix_data(psb_ctxt_)
ictxt = psb_cd_get_context(desc_a)
call psb_info(ictxt, me, np)
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 psb_barrier(ictxt)
endif
call psb_csdp(a,atmp,info)
if(info /= 0) then
info=4010
ch_err='psb_csdp'
call psb_errpush(info,name,a_err=ch_err)
if (toupper(a%fida) /= 'CSR') then
write(0,*) 'Unimplemented input to SLU_BLD'
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 psb_barrier(ictxt)
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
end if
nzb = blck%infoa(psb_nnz_)
if (Debug) then
write(0,*) me, 'SPLUBLD: Done asmatbld',info,nzb,blck%fida
call psb_barrier(ictxt)
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 psb_barrier(ictxt)
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
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
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
end if
nzt = psb_sp_get_nnzeros(atmp)
nzt = psb_sp_get_nnzeros(a)
if (Debug) then
write(0,*) me,'Calling psb_slu_factor ',nzt,atmp%m,&
& atmp%k,p%desc_data%matrix_data(psb_n_row_)
write(0,*) me,'Calling psb_slu_factor ',nzt,a%m,&
& a%k,p%desc_data%matrix_data(psb_n_row_)
call psb_barrier(ictxt)
endif
call psb_dslu_factor(atmp%m,nzt,&
& atmp%aspk,atmp%ia2,atmp%ia1,p%iprcparm(slu_ptr_),info)
if(info /= 0) then
call psb_dslu_factor(a%m,nzt,&
& a%aspk,a%ia2,a%ia1,p%iprcparm(slu_ptr_),info)
if (info /= 0) then
ch_err='psb_slu_fact'
call psb_errpush(4110,name,a_err=ch_err,i_err=(/info,0,0,0,0/))
goto 9999
@@ -182,15 +86,6 @@ subroutine psb_dslu_bld(a,desc_a,p,info)
call psb_barrier(ictxt)
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
end if
call psb_erractionrestore(err_act)
return
psb_dsp_renum.f90
+36 -193
View File
@@ -36,7 +36,8 @@
!!$
subroutine psb_dsp_renum(a,desc_a,blck,p,atmp,info)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_dsp_renum
implicit none
! .. array Arguments ..
@@ -57,24 +58,32 @@ subroutine psb_dsp_renum(a,desc_a,blck,p,atmp,info)
if (psb_get_errstatus().ne.0) return
info=0
name='apply_renum'
name='psb_dsp_renum'
call psb_erractionsave(err_act)
ictxt=psb_cd_get_context(desc_a)
call psb_info(ictxt, me, np)
!!!!!!!!!!!!!!!! CHANGE FOR NON-CSR A
!
! CHANGE: Start with a COO atmp. Then change if/when necessary.
! Exit with a COO atmp.
!
! Renumbering type:
! 1. Global column indices
! (2. GPS band reduction disabled for the time being)
nztota=psb_sp_get_nnzeros(a)
nztotb=psb_sp_get_nnzeros(blck)
call psb_sp_reall(atmp,nztota+nztotb,info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='psb_sp_reall')
goto 9999
end if
atmp%fida='COO'
call psb_csdp(a,atmp,info)
call psb_rwextd(a%m+blck%m,atmp,info,blck,rowscale=.false.)
if (p%iprcparm(iren_)==renum_glb_) then
atmp%m = a%m + blck%m
atmp%k = a%k
atmp%fida='CSR'
atmp%descra = 'GUN'
! This is the renumbering coherent with global indices..
mglob = psb_cd_get_global_rows(desc_a)
@@ -104,112 +113,10 @@ subroutine psb_dsp_renum(a,desc_a,blck,p,atmp,info)
enddo
t3 = psb_wtime()
! Build ATMP with new numbering.
nztmp=size(atmp%aspk)
allocate(itmp(max(8,atmp%m+2,nztmp+2)),rtmp(atmp%m),stat=info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='Allocate')
goto 9999
end if
j = 1
atmp%ia2(1) = 1
do i=1, atmp%m
ir = p%perm(i)
if (ir <= a%m ) then
nzl = a%ia2(ir+1) - a%ia2(ir)
if (nzl > size(rtmp)) then
call psb_realloc(nzl,rtmp,info)
if(info/=0) then
info=4010
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
endif
jj = a%ia2(ir)
k=0
do kk=1, nzl
if (a%ia1(jj+kk-1)<=atmp%m) then
k = k + 1
rtmp(k) = a%aspk(jj+kk-1)
atmp%ia1(j+k-1) = p%invperm(a%ia1(jj+kk-1))
endif
enddo
call psb_msort(atmp%ia1(j:j+k-1),ix=itmp2)
do kk=1,k
atmp%aspk(j+kk-1) = rtmp(itmp2(kk))
enddo
else if (ir <= atmp%m ) then
ir = ir - a%m
nzl = blck%ia2(ir+1) - blck%ia2(ir)
if (nzl > size(rtmp)) then
call psb_realloc(nzl,rtmp,info)
if(info/=0) then
info=4010
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
endif
jj = blck%ia2(ir)
k=0
do kk=1, nzl
if (blck%ia1(jj+kk-1)<=atmp%m) then
k = k + 1
rtmp(k) = blck%aspk(jj+kk-1)
atmp%ia1(j+k-1) = p%invperm(blck%ia1(jj+kk-1))
endif
enddo
call psb_msort(atmp%ia1(j:j+k-1),ix=itmp2)
do kk=1,k
atmp%aspk(j+kk-1) = rtmp(itmp2(kk))
enddo
else
write(0,*) 'Row index error 1 :',i,ir
endif
j = j + k
atmp%ia2(i+1) = j
enddo
t4 = psb_wtime()
deallocate(itmp,itmp2,rtmp)
else if (p%iprcparm(iren_)==renum_gps_) then
atmp%m = a%m + blck%m
atmp%k = a%k
atmp%fida='CSR'
atmp%descra = 'GUN'
do i=1, a%m
atmp%ia2(i) = a%ia2(i)
do j= a%ia2(i), a%ia2(i+1)-1
atmp%ia1(j) = a%ia1(j)
enddo
enddo
atmp%ia2(a%m+1) = a%ia2(a%m+1)
nztota = atmp%ia2(a%m+1) -1
if (blck%m>0) then
do i=1, blck%m
atmp%ia2(a%m+i) = nztota+blck%ia2(i)
do j= blck%ia2(i), blck%ia2(i+1)-1
atmp%ia1(nztota+j) = blck%ia1(j)
enddo
enddo
atmp%ia2(atmp%m+1) = nztota+blck%ia2(blck%m+1)
endif
nztmp = atmp%ia2(atmp%m+1) - 1
call psb_ipcoo2csr(atmp,info)
nztmp = psb_sp_get_nnzeros(atmp)
! This is a renumbering with Gibbs-Poole-Stockmeyer
! band reduction. Switched off for now. To be fixed,
! gps_reduction should get p%perm.
@@ -263,90 +170,26 @@ subroutine psb_dsp_renum(a,desc_a,blck,p,atmp,info)
p%invperm(p%perm(k)) = k
enddo
t3 = psb_wtime()
! Build ATMP with new numbering.
allocate(itmp2(max(8,atmp%m+2,nztmp+2)),rtmp(atmp%m),stat=info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='Allocate')
goto 9999
end if
j = 1
atmp%ia2(1) = 1
do i=1, atmp%m
ir = p%perm(i)
if (ir <= a%m ) then
nzl = a%ia2(ir+1) - a%ia2(ir)
if (nzl > size(rtmp)) then
call psb_realloc(nzl,rtmp,info)
if(info/=0) then
info=4010
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
endif
jj = a%ia2(ir)
k=0
do kk=1, nzl
if (a%ia1(jj+kk-1)<=atmp%m) then
k = k + 1
rtmp(k) = a%aspk(jj+kk-1)
atmp%ia1(j+k-1) = p%invperm(a%ia1(jj+kk-1))
endif
enddo
call psb_msort(atmp%ia1(j:j+k-1),ix=itmp2)
do kk=1,k
atmp%aspk(j+kk-1) = rtmp(itmp2(kk))
enddo
else if (ir <= atmp%m ) then
ir = ir - a%m
nzl = blck%ia2(ir+1) - blck%ia2(ir)
if (nzl > size(rtmp)) then
call psb_realloc(nzl,rtmp,info)
if(info/=0) then
info=4010
ch_err='psb_realloc'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
endif
jj = blck%ia2(ir)
k=0
do kk=1, nzl
if (blck%ia1(jj+kk-1)<=atmp%m) then
k = k + 1
rtmp(k) = blck%aspk(jj+kk-1)
atmp%ia1(j+k-1) = p%invperm(blck%ia1(jj+kk-1))
endif
enddo
call psb_msort(atmp%ia1(j:j+k-1),ix=itmp2)
do kk=1,k
atmp%aspk(j+kk-1) = rtmp(itmp2(kk))
enddo
else
write(0,*) 'Row index error 1 :',i,ir
endif
j = j + k
atmp%ia2(i+1) = j
enddo
t4 = psb_wtime()
deallocate(itmp,itmp2,rtmp)
call psb_ipcsr2coo(atmp,info)
end if
! Rebuild ATMP with new numbering.
nztmp=psb_sp_get_nnzeros(atmp)
do i=1,nztmp
atmp%ia1(i) = p%perm(a%ia1(i))
atmp%ia2(i) = p%invperm(a%ia2(i))
end do
call psb_fixcoo(atmp,info)
if (info /= 0) then
call psb_errpush(4010,name,a_err='psb_fixcoo')
goto 9999
end if
t4 = psb_wtime()
call psb_erractionrestore(err_act)
return
psb_dumf_bld.f90
+13 -117
View File
@@ -36,7 +36,8 @@
!!$
subroutine psb_dumf_bld(a,desc_a,p,info)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_dumf_bld
implicit none
type(psb_dspmat_type), intent(inout) :: a
@@ -45,140 +46,43 @@ subroutine psb_dumf_bld(a,desc_a,p,info)
integer, intent(out) :: info
type(psb_dspmat_type) :: blck, atmp
character(len=5) :: fmt
character :: upd='F'
integer :: i,j,nza,nzb,nzt,ictxt,me,np,err_act
integer :: i_err(5)
logical, parameter :: debug=.false.
character(len=20) :: name, ch_err
interface psb_asmatbld
Subroutine psb_dasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
use psb_base_mod
use psb_prec_type
integer, intent(in) :: ptype,novr
Type(psb_dspmat_type), Intent(in) :: a
Type(psb_dspmat_type), Intent(inout) :: blk
Type(psb_desc_type), Intent(inout) :: desc_p
Type(psb_desc_type), Intent(in) :: desc_data
Character, Intent(in) :: upd
integer, intent(out) :: info
character(len=5), optional :: outfmt
end Subroutine psb_dasmatbld
end interface
info=0
name='psb_umf_bld'
name='psb_dumf_bld'
call psb_erractionsave(err_act)
ictxt = desc_A%matrix_data(psb_ctxt_)
call psb_info(ictxt, me, np)
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 psb_barrier(ictxt)
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)
if (toupper(a%fida) /= 'CSC') then
write(0,*) 'Unimplemented input to UMF_BLD'
goto 9999
end if
nza = psb_sp_get_nnzeros(atmp)
nzb = psb_sp_get_nnzeros(a)
if (Debug) then
write(0,*) me, 'UMFBLD: Done csdp',info,nza,atmp%m,atmp%k,nzb
call psb_barrier(ictxt)
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
end if
nzb = psb_sp_get_nnzeros(blck)
if (Debug) then
write(0,*) me, 'UMFBLD: Done asmatbld',info,nzb,blck%fida
call psb_barrier(ictxt)
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 psb_barrier(ictxt)
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
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
enddo
atmp%infoa(psb_nnz_) = i
nzt = psb_sp_get_nnzeros(a)
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
end if
nzt = psb_sp_get_nnzeros(atmp)
if (Debug) then
write(0,*) me,'Calling psb_umf_factor ',nzt,atmp%m,&
& atmp%k,p%desc_data%matrix_data(psb_n_row_)
write(0,*) me,'Calling psb_umf_factor ',nzt,a%m,&
& a%k,p%desc_data%matrix_data(psb_n_row_)
open(80+me)
call psb_csprt(80+me,atmp)
call psb_csprt(80+me,a)
close(80+me)
call psb_barrier(ictxt)
endif
call psb_dumf_factor(atmp%m,nzt,&
& atmp%aspk,atmp%ia1,atmp%ia2,&
call psb_dumf_factor(a%m,nzt,&
& a%aspk,a%ia1,a%ia2,&
& p%iprcparm(umf_symptr_),p%iprcparm(umf_numptr_),info)
if (info /= 0) then
i_err(1) = info
info=4110
ch_err='psb_umf_fact'
call psb_errpush(info,name,a_err=ch_err,i_err=i_err)
call psb_errpush(info,name,a_err='psb_umf_fact',i_err=i_err)
goto 9999
end if
@@ -186,14 +90,6 @@ subroutine psb_dumf_bld(a,desc_a,p,info)
write(0,*) me, 'UMFBLD: Done umf_Factor',info,p%iprcparm(umf_numptr_)
call psb_barrier(ictxt)
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
end if
call psb_erractionrestore(err_act)
return
psb_prec_mod.f90
+28 -6
View File
@@ -249,7 +249,6 @@ module psb_prec_mod
real(kind(0.d0)),target :: work(:)
integer, intent(out) :: info
end subroutine psb_dbjac_aply
subroutine psb_zbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
use psb_base_mod
use psb_prec_type
@@ -285,8 +284,8 @@ module psb_prec_mod
end subroutine psb_zdiagsc_bld
end interface
interface psb_ilu_bld
subroutine psb_dilu_bld(a,desc_data,p,upd,info)
interface psb_bjac_bld
subroutine psb_dbjac_bld(a,desc_data,p,upd,info)
use psb_base_mod
use psb_prec_type
integer, intent(out) :: info
@@ -294,8 +293,8 @@ module psb_prec_mod
type(psb_desc_type),intent(in) :: desc_data
type(psb_dbaseprc_type), intent(inout) :: p
character, intent(in) :: upd
end subroutine psb_dilu_bld
subroutine psb_zilu_bld(a,desc_data,p,upd,info)
end subroutine psb_dbjac_bld
subroutine psb_zbjac_bld(a,desc_data,p,upd,info)
use psb_base_mod
use psb_prec_type
integer, intent(out) :: info
@@ -303,6 +302,29 @@ module psb_prec_mod
type(psb_desc_type),intent(in) :: desc_data
type(psb_zbaseprc_type), intent(inout) :: p
character, intent(in) :: upd
end subroutine psb_zbjac_bld
end interface
interface psb_ilu_bld
subroutine psb_dilu_bld(a,desc_data,p,upd,info,blck)
use psb_base_mod
use psb_prec_type
integer, intent(out) :: info
type(psb_dspmat_type), intent(in), target :: a
type(psb_desc_type),intent(in) :: desc_data
type(psb_dbaseprc_type), intent(inout) :: p
character, intent(in) :: upd
type(psb_dspmat_type), intent(in), optional :: blck
end subroutine psb_dilu_bld
subroutine psb_zilu_bld(a,desc_data,p,upd,info,blck)
use psb_base_mod
use psb_prec_type
integer, intent(out) :: info
type(psb_zspmat_type), intent(in), target :: a
type(psb_desc_type),intent(in) :: desc_data
type(psb_zbaseprc_type), intent(inout) :: p
character, intent(in) :: upd
type(psb_zspmat_type), intent(in), optional :: blck
end subroutine psb_zilu_bld
end interface
@@ -364,7 +386,7 @@ module psb_prec_mod
end subroutine psb_zilu_fct
end interface
interface psb_as_matbld
interface psb_asmatbld
Subroutine psb_dasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
use psb_base_mod
use psb_prec_type
psb_prec_type.f90
+24 -1
View File
@@ -216,16 +216,19 @@ module psb_prec_type
contains
subroutine psb_out_prec_descr(p)
use psb_base_mod
type(psb_dprec_type), intent(in) :: p
call psb_file_prec_descr(6,p)
end subroutine psb_out_prec_descr
subroutine psb_zout_prec_descr(p)
use psb_base_mod
type(psb_zprec_type), intent(in) :: p
call psb_zfile_prec_descr(6,p)
end subroutine psb_zout_prec_descr
subroutine psb_file_prec_descr(iout,p)
use psb_base_mod
integer, intent(in) :: iout
type(psb_dprec_type), intent(in) :: p
integer :: ilev
@@ -305,6 +308,7 @@ contains
end subroutine psb_file_prec_descr
function psb_prec_short_descr(p)
use psb_base_mod
type(psb_dprec_type), intent(in) :: p
character(len=20) :: psb_prec_short_descr
psb_prec_short_descr = ' '
@@ -373,6 +377,7 @@ contains
subroutine psb_zfile_prec_descr(iout,p)
use psb_base_mod
integer, intent(in) :: iout
type(psb_zprec_type), intent(in) :: p
@@ -447,6 +452,7 @@ contains
end subroutine psb_zfile_prec_descr
function psb_zprec_short_descr(p)
use psb_base_mod
type(psb_zprec_type), intent(in) :: p
character(len=20) :: psb_zprec_short_descr
psb_zprec_short_descr = ' '
@@ -517,6 +523,7 @@ contains
function is_legal_base_prec(ip)
use psb_base_mod
integer, intent(in) :: ip
logical :: is_legal_base_prec
@@ -524,6 +531,7 @@ contains
return
end function is_legal_base_prec
function is_legal_n_ovr(ip)
use psb_base_mod
integer, intent(in) :: ip
logical :: is_legal_n_ovr
@@ -531,6 +539,7 @@ contains
return
end function is_legal_n_ovr
function is_legal_renum(ip)
use psb_base_mod
integer, intent(in) :: ip
logical :: is_legal_renum
! For the time being we are disabling renumbering options.
@@ -538,6 +547,7 @@ contains
return
end function is_legal_renum
function is_legal_jac_sweeps(ip)
use psb_base_mod
integer, intent(in) :: ip
logical :: is_legal_jac_sweeps
@@ -545,19 +555,21 @@ contains
return
end function is_legal_jac_sweeps
function is_legal_prolong(ip)
use psb_base_mod
integer, intent(in) :: ip
logical :: is_legal_prolong
is_legal_prolong = ((ip>=psb_none_).and.(ip<=psb_square_root_))
return
end function is_legal_prolong
function is_legal_restrict(ip)
use psb_base_mod
integer, intent(in) :: ip
logical :: is_legal_restrict
is_legal_restrict = ((ip==psb_nohalo_).or.(ip==psb_halo_))
return
end function is_legal_restrict
function is_legal_ml_type(ip)
use psb_base_mod
integer, intent(in) :: ip
logical :: is_legal_ml_type
@@ -565,6 +577,7 @@ contains
return
end function is_legal_ml_type
function is_legal_ml_aggr_kind(ip)
use psb_base_mod
integer, intent(in) :: ip
logical :: is_legal_ml_aggr_kind
@@ -572,6 +585,7 @@ contains
return
end function is_legal_ml_aggr_kind
function is_legal_ml_smooth_pos(ip)
use psb_base_mod
integer, intent(in) :: ip
logical :: is_legal_ml_smooth_pos
@@ -579,6 +593,7 @@ contains
return
end function is_legal_ml_smooth_pos
function is_legal_ml_smth_kind(ip)
use psb_base_mod
integer, intent(in) :: ip
logical :: is_legal_ml_smth_kind
@@ -586,6 +601,7 @@ contains
return
end function is_legal_ml_smth_kind
function is_legal_ml_coarse_mat(ip)
use psb_base_mod
integer, intent(in) :: ip
logical :: is_legal_ml_coarse_mat
@@ -593,6 +609,7 @@ contains
return
end function is_legal_ml_coarse_mat
function is_legal_ml_fact(ip)
use psb_base_mod
integer, intent(in) :: ip
logical :: is_legal_ml_fact
@@ -600,6 +617,7 @@ contains
return
end function is_legal_ml_fact
function is_legal_ml_lev(ip)
use psb_base_mod
integer, intent(in) :: ip
logical :: is_legal_ml_lev
@@ -607,6 +625,7 @@ contains
return
end function is_legal_ml_lev
function is_legal_omega(ip)
use psb_base_mod
real(kind(1.d0)), intent(in) :: ip
logical :: is_legal_omega
@@ -614,6 +633,7 @@ contains
return
end function is_legal_omega
function is_legal_ml_eps(ip)
use psb_base_mod
real(kind(1.d0)), intent(in) :: ip
logical :: is_legal_ml_eps
@@ -623,6 +643,7 @@ contains
subroutine psb_icheck_def(ip,name,id,is_legal)
use psb_base_mod
integer, intent(inout) :: ip
integer, intent(in) :: id
character(len=*), intent(in) :: name
@@ -640,6 +661,7 @@ contains
end subroutine psb_icheck_def
subroutine psb_dcheck_def(ip,name,id,is_legal)
use psb_base_mod
real(kind(1.d0)), intent(inout) :: ip
real(kind(1.d0)), intent(in) :: id
character(len=*), intent(in) :: name
@@ -823,6 +845,7 @@ contains
function pr_to_str(iprec)
use psb_base_mod
integer, intent(in) :: iprec
character(len=10) :: pr_to_str
psb_zasmatbld.f90
+3 -1
View File
@@ -54,7 +54,8 @@
Subroutine psb_zasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_zasmatbld
Implicit None
! .. Array Arguments ..
@@ -77,6 +78,7 @@ Subroutine psb_zasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
& tot_recv, ircode, n_row,nhalo, nrow_a,err_act
Logical,Parameter :: debug=.false., debugprt=.false.
character(len=20) :: name, ch_err
name='psb_zasmatbld'
if(psb_get_errstatus().ne.0) return
info=0
psb_zbaseprc_aply.f90
+2 -15
View File
@@ -40,7 +40,8 @@ subroutine psb_zbaseprc_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
! where K is a a basic preconditioner stored in prec
!
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_zbaseprc_aply
implicit none
type(psb_desc_type),intent(in) :: desc_data
@@ -59,20 +60,6 @@ subroutine psb_zbaseprc_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
real(kind(1.d0)) :: t1, t2, t3, t4, t5, t6, t7
logical,parameter :: debug=.false., debugprt=.false.
character(len=20) :: name, ch_err
interface psb_bjac_aply
subroutine psb_zbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
use psb_base_mod
use psb_prec_type
type(psb_desc_type), intent(in) :: desc_data
type(psb_zbaseprc_type), intent(in) :: prec
complex(kind(0.d0)),intent(inout) :: x(:), y(:)
complex(kind(0.d0)),intent(in) :: alpha,beta
character(len=1) :: trans
complex(kind(0.d0)),target :: work(:)
integer, intent(out) :: info
end subroutine psb_zbjac_aply
end interface
name='psb_zbaseprc_aply'
info = 0
psb_zbaseprc_bld.f90
+10 -103
View File
@@ -37,7 +37,8 @@
subroutine psb_zbaseprc_bld(a,desc_a,p,info,upd)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_zbaseprc_bld
Implicit None
type(psb_zspmat_type), target :: a
@@ -46,52 +47,6 @@ subroutine psb_zbaseprc_bld(a,desc_a,p,info,upd)
integer, intent(out) :: info
character, intent(in), optional :: upd
interface psb_diagsc_bld
subroutine psb_zdiagsc_bld(a,desc_data,p,upd,info)
use psb_base_mod
use psb_prec_type
integer, intent(out) :: info
type(psb_zspmat_type), intent(in), target :: a
type(psb_desc_type),intent(in) :: desc_data
type(psb_zbaseprc_type), intent(inout) :: p
character, intent(in) :: upd
end subroutine psb_zdiagsc_bld
end interface
interface psb_ilu_bld
subroutine psb_zilu_bld(a,desc_data,p,upd,info)
use psb_base_mod
use psb_prec_type
integer, intent(out) :: info
type(psb_zspmat_type), intent(in), target :: a
type(psb_desc_type),intent(in) :: desc_data
type(psb_zbaseprc_type), intent(inout) :: p
character, intent(in) :: upd
end subroutine psb_zilu_bld
end interface
interface psb_slu_bld
subroutine psb_zslu_bld(a,desc_a,p,info)
use psb_base_mod
use psb_prec_type
type(psb_zspmat_type), intent(inout) :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_zbaseprc_type), intent(inout) :: p
integer, intent(out) :: info
end subroutine psb_zslu_bld
end interface
interface psb_umf_bld
subroutine psb_zumf_bld(a,desc_a,p,info)
use psb_base_mod
use psb_prec_type
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_zbaseprc_type), intent(inout) :: p
integer, intent(out) :: info
end subroutine psb_zumf_bld
end interface
! Local scalars
Integer :: err, nnzero, n_row, n_col,I,j,k,ictxt,&
& me,mycol,np,npcol,mglob,lw, mtype, nrg, nzg, err_act
@@ -108,7 +63,7 @@ subroutine psb_zbaseprc_bld(a,desc_a,p,info,upd)
info=0
err=0
call psb_erractionsave(err_act)
name = 'psb_baseprc_bld'
name = 'psb_zbaseprc_bld'
if (debug) write(0,*) 'Entering baseprc_bld'
info = 0
@@ -140,13 +95,8 @@ subroutine psb_zbaseprc_bld(a,desc_a,p,info,upd)
call psb_check_def(p%iprcparm(p_type_),'base_prec',&
& diag_,is_legal_base_prec)
!!$ allocate(p%desc_data,stat=info)
!!$ if (info /= 0) then
!!$ call psb_errpush(4010,name,a_err='Allocate')
!!$ goto 9999
!!$ end if
!!$
!!$ call psb_nullify_desc(p%desc_data)
call psb_nullify_desc(p%desc_data)
select case(p%iprcparm(p_type_))
case (noprec_)
@@ -183,59 +133,16 @@ subroutine psb_zbaseprc_bld(a,desc_a,p,info,upd)
call psb_check_def(p%iprcparm(f_type_),'fact',&
& f_ilu_n_,is_legal_ml_fact)
if (debug) write(0,*)me, ': Calling PSB_ILU_BLD'
if (debug) write(0,*)me, ': Calling PSB_BJAC_BLD'
if (debug) call psb_barrier(ictxt)
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 psb_barrier(ictxt)
if(info /= 0) then
info=4010
ch_err='psb_ilu_bld'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
case(f_slu_)
if(debug) write(0,*)me,': calling slu_bld'
call psb_slu_bld(a,desc_a,p,info)
if(info /= 0) then
info=4010
ch_err='slu_bld'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
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'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
case(f_none_)
write(0,*) 'Fact=None in BASEPRC_BLD Bja/ASM??'
call psb_bjac_bld(a,desc_a,p,iupd,info)
if(info /= 0) then
info=4010
ch_err='Inconsistent prec f_none_'
call psb_errpush(info,name,a_err=ch_err)
call psb_errpush(info,name,a_err='psb_bjac_bld')
goto 9999
end if
case default
write(0,*) 'Unknown factor type in baseprc_bld bja/asm: ',&
&p%iprcparm(f_type_)
info=4010
ch_err='Unknown f_type_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end select
case default
info=4010
ch_err='Unknown p_type_'
psb_zbjac_aply.f90
+3 -2
View File
@@ -43,7 +43,8 @@ subroutine psb_zbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
!
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_zbjac_aply
implicit none
type(psb_desc_type), intent(in) :: desc_data
@@ -63,7 +64,7 @@ subroutine psb_zbjac_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
logical,parameter :: debug=.false., debugprt=.false.
character(len=20) :: name, ch_err
name='psb_bjac_aply'
name='psb_zbjac_aply'
info = 0
call psb_erractionsave(err_act)
psb_zbjac_bld.f90
+374
View File
@@ -0,0 +1,374 @@
!!$
!!$
!!$ MD2P4
!!$ Multilevel Domain Decomposition Parallel Preconditioner Package for PSBLAS
!!$ for
!!$ Parallel Sparse BLAS v2.0
!!$ (C) Copyright 2006 Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$ Daniela di Serafino Second University of Naples
!!$ Pasqua D'Ambra ICAR-CNR
!!$
!!$ 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 MD2P4 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 MD2P4 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.
!!$
!!$
!*****************************************************************************
!* *
!* This is where the action takes place. *
!* ASMATBLD does the setup: building the prec descriptor plus retrieving *
!* matrix rows if needed *
!* *
!* *
!* *
!* *
!* some open code does the renumbering *
!* *
!* *
!* *
!* *
!*****************************************************************************
subroutine psb_zbjac_bld(a,desc_a,p,upd,info)
use psb_base_mod
use psb_prec_mod, mld_protect_name => psb_zbjac_bld
implicit none
!
! .. Scalar Arguments ..
integer, intent(out) :: info
! .. array Arguments ..
type(psb_zspmat_type), intent(in), target :: a
type(psb_zbaseprc_type), intent(inout) :: p
type(psb_desc_type), intent(in) :: desc_a
character, intent(in) :: upd
! .. Local Scalars ..
integer :: i, j, jj, k, kk, m
integer :: int_err(5)
character :: trans, unitd
type(psb_zspmat_type) :: blck, atmp
real(kind(1.d0)) :: t1,t2,t3,t4,t5,t6, t7, t8
logical, parameter :: debugprt=.false., debug=.false., aggr_dump=.false.
integer nztota, nztotb, nztmp, nzl, nnr, ir, err_act,&
& n_row, nrow_a,n_col, nhalo, ind, iind, i1,i2,ia
integer :: ictxt,np,me
character(len=20) :: name, ch_err
character(len=5), parameter :: coofmt='COO'
if(psb_get_errstatus().ne.0) return
info=0
name='psb_bjac_bld'
call psb_erractionsave(err_act)
ictxt=psb_cd_get_context(desc_a)
call psb_info(ictxt, me, np)
m = a%m
if (m < 0) then
info = 10
int_err(1) = 1
int_err(2) = m
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
trans = 'N'
unitd = 'U'
if (p%iprcparm(n_ovr_) < 0) then
info = 11
int_err(1) = 1
int_err(2) = p%iprcparm(n_ovr_)
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
call psb_nullify_sp(blck)
call psb_nullify_sp(atmp)
t1= psb_wtime()
if(debug) write(0,*)me,': calling psb_asmatbld',p%iprcparm(p_type_),p%iprcparm(n_ovr_)
if (debug) call psb_barrier(ictxt)
call psb_asmatbld(p%iprcparm(p_type_),p%iprcparm(n_ovr_),a,&
& blck,desc_a,upd,p%desc_data,info,outfmt=coofmt)
if(info/=0) then
info=4010
ch_err='psb_asmatbld'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
t2= psb_wtime()
if (debug) write(0,*)me,': out of psb_asmatbld'
if (debug) call psb_barrier(ictxt)
if (p%iprcparm(iren_) > 0) then
!
! Here we allocate a full copy to hold local A and received BLK
! Done inside sp_renum.
!
call psb_sp_renum(a,desc_a,blck,p,atmp,info)
if(info/=0) then
info=4010
ch_err='psb_sp_renum'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
t3 = psb_wtime()
if (debugprt) then
call psb_barrier(ictxt)
open(40+me)
call psb_csprt(40+me,atmp,head='% Local matrix')
close(40+me)
endif
if (debug) write(0,*) me,' Factoring rows ',&
&atmp%m,a%m,blck%m,atmp%ia2(atmp%m+1)-1
select case(p%iprcparm(f_type_))
case(f_ilu_n_,f_ilu_e_)
call psb_ipcoo2csr(atmp,info)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_ipcoo2csr')
goto 9999
end if
call psb_ilu_bld(atmp,p%desc_data,p,upd,info)
if(info/=0) then
info=4010
call psb_errpush(info,name,a_err='psb_ilu_bld')
goto 9999
end if
if (debugprt) then
open(80+me)
call psb_csprt(80+me,p%av(l_pr_),head='% Local L factor')
write(80+me,*) '% Diagonal: ',p%av(l_pr_)%m
do i=1,p%av(l_pr_)%m
write(80+me,*) i,i,p%d(i)
enddo
call psb_csprt(80+me,p%av(u_pr_),head='% Local U factor')
close(80+me)
endif
case(f_slu_)
call psb_ipcoo2csr(atmp,info)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_ipcoo2csr')
goto 9999
end if
call psb_slu_bld(atmp,p%desc_data,p,info)
if(info /= 0) then
info=4010
call psb_errpush(info,name,a_err='slu_bld')
goto 9999
end if
case(f_umf_)
call psb_ipcoo2csc(atmp,info)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_ipcoo2csc')
goto 9999
end if
call psb_umf_bld(atmp,p%desc_data,p,info)
if(debug) write(0,*)me,': Done umf_bld ',info
if (info /= 0) then
info = 4010
call psb_errpush(info,name,a_err='umf_bld')
goto 9999
end if
case(f_none_)
info=4010
call psb_errpush(info,name,a_err='Inconsistent prec f_none_')
goto 9999
case default
info=4010
call psb_errpush(info,name,a_err='Unknown f_type_')
goto 9999
end select
call psb_sp_free(atmp,info)
if(info/=0) then
info=4010
call psb_errpush(info,name,a_err='psb_sp_free')
goto 9999
end if
else if (p%iprcparm(iren_) == 0) then
select case(p%iprcparm(f_type_))
case(f_ilu_n_,f_ilu_e_)
call psb_ipcoo2csr(blck,info,rwshr=.true.)
if (info==0) call psb_ilu_bld(a,desc_a,p,upd,info,blck=blck)
if(info/=0) then
info=4010
call psb_errpush(info,name,a_err='psb_ilu_bld')
goto 9999
end if
if (debugprt) then
open(80+me)
call psb_csprt(80+me,p%av(l_pr_),head='% Local L factor')
write(80+me,*) '% Diagonal: ',p%av(l_pr_)%m
do i=1,p%av(l_pr_)%m
write(80+me,*) i,i,p%d(i)
enddo
call psb_csprt(80+me,p%av(u_pr_),head='% Local U factor')
close(80+me)
endif
case(f_slu_)
atmp%fida='COO'
call psb_csdp(a,atmp,info)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_csdp')
goto 9999
end if
call psb_rwextd(atmp%m+blck%m,atmp,info,blck,rowscale=.true.)
call psb_ipcoo2csr(atmp,info)
call psb_slu_bld(atmp,p%desc_data,p,info)
if(info /= 0) then
info=4010
call psb_errpush(info,name,a_err='slu_bld')
goto 9999
end if
call psb_sp_free(atmp,info)
if(info/=0) then
info=4010
call psb_errpush(info,name,a_err='psb_sp_free')
goto 9999
end if
case(f_umf_)
atmp%fida='COO'
atmp%fida='COO'
call psb_csdp(a,atmp,info)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_csdp')
goto 9999
end if
call psb_rwextd(atmp%m+blck%m,atmp,info,blck,rowscale=.true.)
if (info == 0) call psb_ipcoo2csc(atmp,info)
if (info /= 0) then
info=4010
call psb_errpush(info,name,a_err='psb_ipcoo2csc')
goto 9999
end if
call psb_umf_bld(atmp,p%desc_data,p,info)
if(debug) write(0,*)me,': Done umf_bld ',info
if (info /= 0) then
info = 4010
call psb_errpush(info,name,a_err='umf_bld')
goto 9999
end if
call psb_sp_free(atmp,info)
if(info/=0) then
info=4010
call psb_errpush(info,name,a_err='psb_sp_free')
goto 9999
end if
case(f_none_)
info=4010
call psb_errpush(info,name,a_err='Inconsistent prec f_none_')
goto 9999
case default
info=4010
call psb_errpush(info,name,a_err='Unknown f_type_')
goto 9999
end select
endif
t6 = psb_wtime()
call psb_sp_free(blck,info)
if(info/=0) then
info=4010
ch_err='psb_sp_free'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if (debug) write(0,*) me,'End of ilu_bld'
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_zbjac_bld
psb_zbldaggrmat.f90
+3 -1
View File
@@ -36,7 +36,8 @@
!!$
subroutine psb_zbldaggrmat(a,desc_a,ac,desc_ac,p,info)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_zbldaggrmat
implicit none
type(psb_zspmat_type), intent(in), target :: a
@@ -49,6 +50,7 @@ subroutine psb_zbldaggrmat(a,desc_a,ac,desc_ac,p,info)
logical, parameter :: aggr_dump=.false.
integer ::ictxt,np,me, err_act
character(len=20) :: name, ch_err
name='psb_zbldaggrmat'
if(psb_get_errstatus().ne.0) return
info=0
psb_zdiagsc_bld.f90
+3 -2
View File
@@ -37,7 +37,8 @@
subroutine psb_zdiagsc_bld(a,desc_a,p,upd,info)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_zdiagsc_bld
Implicit None
type(psb_zspmat_type), target :: a
@@ -62,7 +63,7 @@ subroutine psb_zdiagsc_bld(a,desc_a,p,upd,info)
info=0
err=0
call psb_erractionsave(err_act)
name = 'psb_diagsc_bld'
name = 'psb_zdiagsc_bld'
if (debug) write(0,*) 'Entering diagsc_bld'
info = 0
psb_zgenaggrmap.f90
+3 -2
View File
@@ -36,7 +36,8 @@
!!$
subroutine psb_zgenaggrmap(aggr_type,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_zgenaggrmap
implicit none
integer, intent(in) :: aggr_type
type(psb_zspmat_type), intent(in) :: a
@@ -56,7 +57,7 @@ subroutine psb_zgenaggrmap(aggr_type,a,desc_a,nlaggr,ilaggr,info)
if(psb_get_errstatus().ne.0) return
info=0
name = 'psb_bldaggrmat'
name = 'psb_zgenaggrmap'
call psb_erractionsave(err_act)
!
! Note. At the time being we are ignoring aggr_type
psb_zilu_bld.f90
+42 -193
View File
@@ -49,9 +49,10 @@
!* *
!* *
!*****************************************************************************
subroutine psb_zilu_bld(a,desc_a,p,upd,info)
subroutine psb_zilu_bld(a,desc_a,p,upd,info,blck)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_zilu_bld
implicit none
!
! .. Scalar Arguments ..
@@ -61,128 +62,66 @@ subroutine psb_zilu_bld(a,desc_a,p,upd,info)
type(psb_zbaseprc_type), intent(inout) :: p
type(psb_desc_type), intent(in) :: desc_a
character, intent(in) :: upd
type(psb_zspmat_type), intent(in), optional :: blck
! .. Local Scalars ..
integer :: i, j, jj, k, kk, m
integer :: int_err(5)
character :: trans, unitd
type(psb_zspmat_type) :: blck, atmp
real(kind(1.d0)) :: t1,t2,t3,t4,t5,t6, t7, t8
logical, parameter :: debugprt=.false., debug=.false., aggr_dump=.false.
integer nztota, nztotb, nztmp, nzl, nnr, ir, err_act,&
& n_row, nrow_a,n_col, nhalo, ind, iind, i1,i2,ia
& n_row, nrow_a, ind, iind, i1,i2,ia
integer :: ictxt,np,me
character(len=20) :: name, ch_err
interface psb_ilu_fct
subroutine psb_zilu_fct(a,l,u,d,info,blck)
use psb_base_mod
use psb_prec_type
integer, intent(out) :: info
type(psb_zspmat_type),intent(in) :: a
type(psb_zspmat_type),intent(inout) :: l,u
type(psb_zspmat_type),intent(in), optional, target :: blck
complex(kind(1.d0)), intent(inout) :: d(:)
end subroutine psb_zilu_fct
end interface
interface psb_asmatbld
Subroutine psb_zasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
use psb_base_mod
use psb_prec_type
integer, intent(in) :: ptype,novr
Type(psb_zspmat_type), Intent(in) :: a
Type(psb_zspmat_type), Intent(inout) :: blk
Type(psb_desc_type), Intent(inout) :: desc_p
Type(psb_desc_type), Intent(in) :: desc_data
Character, Intent(in) :: upd
integer, intent(out) :: info
character(len=5), optional :: outfmt
end Subroutine psb_zasmatbld
end interface
interface psb_sp_renum
subroutine psb_zsp_renum(a,desc_a,blck,p,atmp,info)
use psb_base_mod
use psb_prec_type
type(psb_zspmat_type), intent(in) :: a,blck
type(psb_zspmat_type), intent(inout) :: atmp
type(psb_zbaseprc_type), intent(inout) :: p
type(psb_desc_type), intent(in) :: desc_a
integer, intent(out) :: info
end subroutine psb_zsp_renum
end interface
if(psb_get_errstatus().ne.0) return
info=0
name='psb_ilu_bld'
name='psb_zilu_bld'
call psb_erractionsave(err_act)
ictxt=psb_cd_get_context(desc_a)
call psb_info(ictxt, me, np)
m = a%m
if (m < 0) then
info = 10
int_err(1) = 1
int_err(2) = m
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
trans = 'N'
unitd = 'U'
if (p%iprcparm(n_ovr_) < 0) then
info = 11
int_err(1) = 1
int_err(2) = p%iprcparm(n_ovr_)
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
call psb_nullify_sp(blck)
call psb_nullify_sp(atmp)
t1= psb_wtime()
if(debug) write(0,*)me,': calling psb_asmatbld',p%iprcparm(p_type_),p%iprcparm(n_ovr_)
if (debug) call psb_barrier(ictxt)
call psb_asmatbld(p%iprcparm(p_type_),p%iprcparm(n_ovr_),a,&
& blck,desc_a,upd,p%desc_data,info)
if(info/=0) then
info=4010
ch_err='psb_asmatbld'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
t2= psb_wtime()
if (debug) write(0,*)me,': out of psb_asmatbld'
if (debug) call psb_barrier(ictxt)
if (allocated(p%av)) then
if (size(p%av) < bp_ilu_avsz) then
call psb_errpush(4010,name,a_err='Insufficient av size')
goto 9999
do i=1,size(p%av)
call psb_sp_free(p%av(i),info)
if (info /= 0) then
! Actually, we don't care here about this.
! Just let it go.
! return
end if
enddo
deallocate(p%av,stat=info)
endif
else
call psb_errpush(4010,name,a_err='AV not associated')
goto 9999
end if
if (.not.allocated(p%av)) then
allocate(p%av(max_avsz),stat=info)
if (info /= 0) then
call psb_errpush(4000,name)
goto 9999
end if
endif
!!$ call psb_csprt(50+me,a,head='% (A)')
nrow_a = psb_cd_get_local_rows(desc_a)
nztota = psb_sp_get_nnzeros(a)
nztotb = psb_sp_get_nnzeros(blck)
if (debug) write(0,*)me,': out get_nnzeros',nztota
if (present(blck)) then
nztota = nztota + psb_sp_get_nnzeros(blck)
end if
if (debug) write(0,*)me,': out get_nnzeros',nztota,a%m,a%k
if (debug) call psb_barrier(ictxt)
n_col = psb_cd_get_local_cols(desc_a)
nhalo = n_col-nrow_a
n_row = p%desc_data%matrix_data(psb_n_row_)
p%av(l_pr_)%m = n_row
p%av(l_pr_)%k = n_row
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+nztotb,info)
if (info == 0) call psb_sp_all(n_row,n_row,p%av(u_pr_),nztota+nztotb,info)
call psb_sp_all(n_row,n_row,p%av(l_pr_),nztota,info)
if (info == 0) call psb_sp_all(n_row,n_row,p%av(u_pr_),nztota,info)
if(info/=0) then
info=4010
ch_err='psb_sp_all'
@@ -205,101 +144,17 @@ subroutine psb_zilu_bld(a,desc_a,p,upd,info)
endif
if (debug) then
write(0,*) me,'Done psb_asmatbld'
call psb_barrier(ictxt)
endif
if (p%iprcparm(iren_) > 0) then
!
! Here we allocate a full copy to hold local A and received BLK
!
nztota = psb_sp_get_nnzeros(a)
nztotb = psb_sp_get_nnzeros(blck)
call psb_sp_all(atmp,nztota+nztotb,info)
if(info/=0) then
info=4011
call psb_errpush(info,name)
goto 9999
end if
! write(0,*) 'ILU_BLD ',nztota,nztotb,a%m
call psb_sp_renum(a,desc_a,blck,p,atmp,info)
if(info/=0) then
info=4010
ch_err='psb_sp_renum'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
t3 = psb_wtime()
if (debugprt) then
call psb_barrier(ictxt)
open(40+me)
call psb_csprt(40+me,atmp,head='% Local matrix')
close(40+me)
endif
if (debug) write(0,*) me,' Factoring rows ',&
&atmp%m,a%m,blck%m,atmp%ia2(atmp%m+1)-1
!
! Ok, factor the matrix.
!
t5 = psb_wtime()
blck%m=0
call psb_ilu_fct(atmp,p%av(l_pr_),p%av(u_pr_),p%d,info,blck=blck)
if(info/=0) then
info=4010
ch_err='psb_ilu_fct'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
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
end if
else if (p%iprcparm(iren_) == 0) then
t3 = psb_wtime()
! This is where we have mo renumbering, thus no need
! for ATMP
if (debugprt) then
open(40+me)
call psb_barrier(ictxt)
call psb_csprt(40+me,a,iv=p%desc_data%loc_to_glob,&
& head='% Local matrix')
if (p%iprcparm(p_type_)==asm_) then
call psb_csprt(40+me,blck,iv=p%desc_data%loc_to_glob,&
& irs=a%m,head='% Received rows')
endif
close(40+me)
endif
t5= psb_wtime()
if (debug) write(0,*) me,' Going for ilu_fct'
if (debug) call psb_barrier(ictxt)
call psb_ilu_fct(a,p%av(l_pr_),p%av(u_pr_),p%d,info,blck=blck)
if(info/=0) then
info=4010
ch_err='psb_ilu_fct'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if (debug) write(0,*) me,' Done dilu_fct'
endif
!
! Ok, factor the matrix.
!
t5 = psb_wtime()
call psb_ilu_fct(a,p%av(l_pr_),p%av(u_pr_),p%d,info,blck=blck)
if(info/=0) then
info=4010
ch_err='psb_ilu_fct'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if (debugprt) then
@@ -318,6 +173,8 @@ subroutine psb_zilu_bld(a,desc_a,p,upd,info)
close(80+me)
endif
!!$ call psb_csprt(60+me,a,head='% (A)')
! ierr = MPE_Log_event( ifcte, 0, "st SIMPLE" )
t6 = psb_wtime()
@@ -325,14 +182,6 @@ subroutine psb_zilu_bld(a,desc_a,p,upd,info)
! write(0,'(i3,1x,a,3(1x,g18.9))') me,'renum/factor time',t3-t2,t6-t5
! if (me==0) write(0,'(a,3(1x,g18.9))') 'renum/factor time',t3-t2,t6-t5
call psb_sp_free(blck,info)
if(info/=0) then
info=4010
ch_err='psb_sp_free'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if (psb_sp_getifld(psb_upd_,p%av(u_pr_),info) /= psb_upd_perm_) then
call psb_sp_trimsize(p%av(u_pr_),i1,i2,ia,info)
if (info == 0) call psb_sp_reall(p%av(u_pr_),i1,i2,ia,info)
psb_zilu_fct.f90
+2 -2
View File
@@ -42,7 +42,6 @@ subroutine psb_zilu_fct(a,l,u,d,info,blck)
!
!
use psb_base_mod
use psb_prec_type
implicit none
! .. Scalar Arguments ..
integer, intent(out) :: info
@@ -55,7 +54,8 @@ subroutine psb_zilu_fct(a,l,u,d,info,blck)
integer :: i, j, jj, k, kk, l1, l2, ll, low1, low2,m,ma,err_act
type(psb_zspmat_type), pointer :: blck_
character(len=20) :: name, ch_err
name='psb_zcsrlu'
name='psb_zilu_fct'
info = 0
call psb_erractionsave(err_act)
! .. Executable Statements ..
psb_zmlprc_aply.f90
+3 -16
View File
@@ -84,7 +84,8 @@ subroutine psb_zmlprc_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info)
!
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_zmlprc_aply
implicit none
type(psb_desc_type),intent(in) :: desc_data
@@ -111,21 +112,7 @@ subroutine psb_zmlprc_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info)
end type psb_mlprec_wrk_type
type(psb_mlprec_wrk_type), allocatable :: mlprec_wrk(:)
interface psb_baseprc_aply
subroutine psb_zbaseprc_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
use psb_base_mod
use psb_prec_type
type(psb_desc_type),intent(in) :: desc_data
type(psb_zbaseprc_type), intent(in) :: prec
complex(kind(1.d0)),intent(inout) :: x(:), y(:)
complex(kind(1.d0)),intent(in) :: alpha,beta
character(len=1) :: trans
complex(kind(1.d0)),target :: work(:)
integer, intent(out) :: info
end subroutine psb_zbaseprc_aply
end interface
name='psb_mlprc_aply'
name='psb_zmlprc_aply'
info = 0
call psb_erractionsave(err_act)
psb_zmlprc_bld.f90
+3 -41
View File
@@ -37,7 +37,8 @@
subroutine psb_zmlprc_bld(a,desc_a,p,info)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_zmlprc_bld
implicit none
type(psb_zspmat_type), intent(in), target :: a
@@ -51,48 +52,9 @@ subroutine psb_zmlprc_bld(a,desc_a,p,info)
character(len=20) :: name, ch_err
logical, parameter :: debug=.false.
type(psb_zspmat_type) :: ac
interface psb_baseprc_bld
subroutine psb_zbaseprc_bld(a,desc_a,p,info,upd)
use psb_base_mod
use psb_prec_type
type(psb_zspmat_type), target :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_zbaseprc_type),intent(inout) :: p
integer, intent(out) :: info
character, intent(in), optional :: upd
end subroutine psb_zbaseprc_bld
end interface
interface psb_genaggrmap
subroutine psb_zgenaggrmap(aggr_type,a,desc_a,nlaggr,ilaggr,info)
use psb_base_mod
use psb_prec_type
implicit none
integer, intent(in) :: aggr_type
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
integer, allocatable :: ilaggr(:),nlaggr(:)
integer, intent(out) :: info
end subroutine psb_zgenaggrmap
end interface
interface psb_bldaggrmat
subroutine psb_zbldaggrmat(a,desc_a,ac,desc_ac,p,info)
use psb_base_mod
use psb_prec_type
type(psb_zspmat_type), intent(in), target :: a
type(psb_zbaseprc_type), intent(inout),target :: p
type(psb_zspmat_type), intent(out),target :: ac
type(psb_desc_type), intent(in) :: desc_a
type(psb_desc_type), intent(inout) :: desc_ac
integer, intent(out) :: info
end subroutine psb_zbldaggrmat
end interface
integer :: ictxt, np, me
name='psb_mlprec_bld'
name='psb_zmlprec_bld'
if (psb_get_errstatus().ne.0) return
info = 0
ictxt = psb_cd_get_context(desc_a)
psb_zprc_aply.f90
+5 -44
View File
@@ -37,7 +37,8 @@
subroutine psb_zprc_aply(prec,x,y,desc_data,info,trans, work)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_zprc_aply
implicit none
type(psb_desc_type),intent(in) :: desc_data
@@ -53,34 +54,6 @@ subroutine psb_zprc_aply(prec,x,y,desc_data,info,trans, work)
integer :: ictxt,np,me,err_act
logical,parameter :: debug=.false., debugprt=.false.
character(len=20) :: name
interface psb_baseprc_aply
subroutine psb_zbaseprc_aply(alpha,prec,x,beta,y,desc_data,trans,work,info)
use psb_base_mod
use psb_prec_type
type(psb_desc_type),intent(in) :: desc_data
type(psb_zbaseprc_type), intent(in) :: prec
complex(kind(0.d0)),intent(inout) :: x(:), y(:)
complex(kind(0.d0)),intent(in) :: alpha,beta
character(len=1) :: trans
complex(kind(0.d0)),target :: work(:)
integer, intent(out) :: info
end subroutine psb_zbaseprc_aply
end interface
interface psb_mlprc_aply
subroutine psb_zmlprc_aply(alpha,baseprecv,x,beta,y,desc_data,trans,work,info)
use psb_base_mod
use psb_prec_type
type(psb_desc_type),intent(in) :: desc_data
type(psb_zbaseprc_type), intent(in) :: baseprecv(:)
complex(kind(0.d0)),intent(in) :: alpha,beta
complex(kind(0.d0)),intent(inout) :: x(:), y(:)
character :: trans
complex(kind(0.d0)),target :: work(:)
integer, intent(out) :: info
end subroutine psb_zmlprc_aply
end interface
name='psb_zprc_aply'
info = 0
@@ -180,7 +153,8 @@ end subroutine psb_zprc_aply
!!$
subroutine psb_zprc_aply1(prec,x,desc_data,info,trans)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_zprc_aply1
implicit none
type(psb_desc_type),intent(in) :: desc_data
@@ -190,25 +164,12 @@ subroutine psb_zprc_aply1(prec,x,desc_data,info,trans)
character(len=1), optional :: trans
logical,parameter :: debug=.false., debugprt=.false.
interface
subroutine psb_zprc_aply(prec,x,y,desc_data,info,trans, work)
use psb_base_mod
use psb_prec_type
type(psb_desc_type),intent(in) :: desc_data
type(psb_zprec_type), intent(in) :: prec
complex(kind(0.d0)),intent(inout) :: x(:), y(:)
integer, intent(out) :: info
character(len=1), optional :: trans
complex(kind(0.d0)), optional, target :: work(:)
end subroutine psb_zprc_aply
end interface
! Local variables
character :: trans_
integer :: ictxt,np,me,i, isz, err_act, int_err(5)
complex(kind(1.d0)), pointer :: WW(:), w1(:)
character(len=20) :: name, ch_err
name='psb_zprec1'
name='psb_zprec_aply1'
info = 0
call psb_erractionsave(err_act)
psb_zprecbld.f90
+2 -3
View File
@@ -37,8 +37,7 @@
subroutine psb_zprecbld(a,desc_a,p,info,upd)
use psb_base_mod
use psb_prec_type
use psb_prec_mod
use psb_prec_mod, mld_protect_name => psb_zprecbld
Implicit None
type(psb_zspmat_type), target :: a
@@ -61,7 +60,7 @@ subroutine psb_zprecbld(a,desc_a,p,info,upd)
info=0
err=0
call psb_erractionsave(err_act)
name = 'psb_precbld'
name = 'psb_zprecbld'
if (debug) write(0,*) 'Entering precbld',P%prec,desc_a%matrix_data(:)
info = 0
psb_zprecfree.f90
+3 -2
View File
@@ -36,7 +36,8 @@
!!$
subroutine psb_zprecfree(p,info)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_zprecfree
implicit none
!....parameters...
@@ -49,7 +50,7 @@ subroutine psb_zprecfree(p,info)
if(psb_get_errstatus().ne.0) return
info=0
name = 'pszprecfree'
name = 'psb_zprecfree'
call psb_erractionsave(err_act)
me=-1
psb_zprecset.f90
+2 -1
View File
@@ -37,7 +37,8 @@
subroutine psb_zprecset(p,ptype,info,iv,rs,rv,ilev,nlev)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_zprecset
implicit none
type(psb_zprec_type), intent(inout) :: p
psb_zslu_bld.f90
+15 -118
View File
@@ -36,140 +36,46 @@
!!$
subroutine psb_zslu_bld(a,desc_a,p,info)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_zslu_bld
implicit none
type(psb_zspmat_type), intent(inout) :: a
type(psb_zspmat_type), intent(in) :: a
type(psb_desc_type), intent(in) :: desc_a
type(psb_zbaseprc_type), intent(inout) :: p
integer, intent(out) :: info
type(psb_zspmat_type) :: blck, atmp
character(len=5) :: fmt
character :: upd='F'
integer :: i,j,nza,nzb,nzt,ictxt, me,np,err_act
integer :: i,j,nza,nzb,nzt,ictxt,me,np,err_act
logical, parameter :: debug=.false.
character(len=20) :: name, ch_err
interface psb_asmatbld
Subroutine psb_zasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
use psb_base_mod
use psb_prec_type
integer, intent(in) :: ptype,novr
Type(psb_zspmat_type), Intent(in) :: a
Type(psb_zspmat_type), Intent(inout) :: blk
Type(psb_desc_type), Intent(inout) :: desc_p
Type(psb_desc_type), Intent(in) :: desc_data
Character, Intent(in) :: upd
integer, intent(out) :: info
character(len=5), optional :: outfmt
end Subroutine psb_zasmatbld
end interface
if(psb_get_errstatus().ne.0) return
info=0
name='psb_slu_bld'
name='psb_zslu_bld'
call psb_erractionsave(err_act)
ictxt = desc_A%matrix_data(psb_ctxt_)
ictxt = psb_cd_get_context(desc_a)
call psb_info(ictxt, me, np)
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 psb_barrier(ictxt)
endif
call psb_csdp(a,atmp,info)
if(info /= 0) then
info=4010
ch_err='psb_csdp'
call psb_errpush(info,name,a_err=ch_err)
if (toupper(a%fida) /= 'CSR') then
write(0,*) 'Unimplemented input to SLU_BLD'
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 psb_barrier(ictxt)
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
end if
nzb = blck%infoa(psb_nnz_)
if (Debug) then
write(0,*) me, 'SPLUBLD: Done asmatbld',info,nzb,blck%fida
call psb_barrier(ictxt)
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 psb_barrier(ictxt)
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
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
enddo
atmp%infoa(psb_nnz_) = i
nzt = psb_sp_get_nnzeros(a)
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
end if
nzt = psb_sp_get_nnzeros(atmp)
if (Debug) then
write(0,*) me,'Calling psb_slu_factor ',nzt,atmp%m,&
& atmp%k,p%desc_data%matrix_data(psb_n_row_)
write(0,*) me,'Calling psb_slu_factor ',nzt,a%m,&
& a%k,p%desc_data%matrix_data(psb_n_row_)
call psb_barrier(ictxt)
endif
call psb_zslu_factor(atmp%m,nzt,&
& atmp%aspk,atmp%ia1,atmp%ia2,p%iprcparm(slu_ptr_),info)
if(info /= 0) then
call psb_zslu_factor(a%m,nzt,&
& a%aspk,a%ia2,a%ia1,p%iprcparm(slu_ptr_),info)
if (info /= 0) then
ch_err='psb_slu_fact'
call psb_errpush(4110,name,a_err=ch_err,i_err=(/info,0,0,0,0/))
goto 9999
@@ -180,15 +86,6 @@ subroutine psb_zslu_bld(a,desc_a,p,info)
call psb_barrier(ictxt)
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
end if
call psb_erractionrestore(err_act)
return
psb_zsp_renum.f90
+3 -2
View File
@@ -36,7 +36,8 @@
!!$
subroutine psb_zsp_renum(a,desc_a,blck,p,atmp,info)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_zsp_renum
implicit none
! .. array Arguments ..
@@ -57,7 +58,7 @@ subroutine psb_zsp_renum(a,desc_a,blck,p,atmp,info)
if (psb_get_errstatus().ne.0) return
info=0
name='apply_renum'
name='psb_zsp_renum'
call psb_erractionsave(err_act)
ictxt=psb_cd_get_context(desc_a)
psb_zumf_bld.f90
+14 -117
View File
@@ -36,7 +36,8 @@
!!$
subroutine psb_zumf_bld(a,desc_a,p,info)
use psb_base_mod
use psb_prec_type
use psb_prec_mod, mld_protect_name => psb_zumf_bld
implicit none
type(psb_zspmat_type), intent(inout) :: a
@@ -45,139 +46,43 @@ subroutine psb_zumf_bld(a,desc_a,p,info)
integer, intent(out) :: info
type(psb_zspmat_type) :: blck, atmp
character(len=5) :: fmt
character :: upd='F'
integer :: i,j,nza,nzb,nzt,ictxt, me,np,err_act
integer :: i,j,nza,nzb,nzt,ictxt,me,np,err_act
integer :: i_err(5)
logical, parameter :: debug=.false.
character(len=20) :: name, ch_err
interface psb_asmatbld
Subroutine psb_zasmatbld(ptype,novr,a,blk,desc_data,upd,desc_p,info,outfmt)
use psb_base_mod
use psb_prec_type
integer, intent(in) :: ptype,novr
Type(psb_zspmat_type), Intent(in) :: a
Type(psb_zspmat_type), Intent(inout) :: blk
Type(psb_desc_type), Intent(inout) :: desc_p
Type(psb_desc_type), Intent(in) :: desc_data
Character, Intent(in) :: upd
integer, intent(out) :: info
character(len=5), optional :: outfmt
end Subroutine psb_zasmatbld
end interface
info=0
name='psb_umf_bld'
name='psb_zumf_bld'
call psb_erractionsave(err_act)
ictxt = desc_A%matrix_data(psb_ctxt_)
call psb_info(ictxt, me, np)
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 psb_barrier(ictxt)
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)
if (toupper(a%fida) /= 'CSC') then
write(0,*) 'Unimplemented input to UMF_BLD'
goto 9999
end if
nza = psb_sp_get_nnzeros(atmp)
nzb = psb_sp_get_nnzeros(a)
if (Debug) then
write(0,*) me, 'UMFBLD: Done csdp',info,nza,atmp%m,atmp%k,nzb
call psb_barrier(ictxt)
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
end if
nzb = psb_sp_get_nnzeros(blck)
if (Debug) then
write(0,*) me, 'UMFBLD: Done asmatbld',info,nzb,blck%fida
call psb_barrier(ictxt)
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 psb_barrier(ictxt)
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
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
enddo
atmp%infoa(psb_nnz_) = i
nzt = psb_sp_get_nnzeros(a)
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
end if
nzt = psb_sp_get_nnzeros(atmp)
if (Debug) then
write(0,*) me,'Calling psb_umf_factor ',nzt,atmp%m,&
& atmp%k,p%desc_data%matrix_data(psb_n_row_)
write(0,*) me,'Calling psb_umf_factor ',nzt,a%m,&
& a%k,p%desc_data%matrix_data(psb_n_row_)
open(80+me)
call psb_csprt(80+me,atmp)
call psb_csprt(80+me,a)
close(80+me)
call psb_barrier(ictxt)
endif
call psb_zumf_factor(atmp%m,nzt,&
& atmp%aspk,atmp%ia1,atmp%ia2,&
call psb_zumf_factor(a%m,nzt,&
& a%aspk,a%ia1,a%ia2,&
& p%iprcparm(umf_symptr_),p%iprcparm(umf_numptr_),info)
if (info /= 0) then
i_err(1) = info
info=4110
ch_err='psb_umf_fact'
call psb_errpush(info,name,a_err=ch_err,i_err=i_err)
call psb_errpush(info,name,a_err='psb_umf_fact',i_err=i_err)
goto 9999
end if
@@ -185,14 +90,6 @@ subroutine psb_zumf_bld(a,desc_a,p,info)
write(0,*) me, 'UMFBLD: Done umf_Factor',info,p%iprcparm(umf_numptr_)
call psb_barrier(ictxt)
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
end if
call psb_erractionrestore(err_act)
return