fdurastante
/
amg4psblas
mirror of https://github.com/sfilippone/amg4psblas.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

psblas3:

Browse Source 
mlprec/mld_caggrmat_nosmth_asb.F90
 mlprec/mld_caggrmat_smth_asb.F90
 mlprec/mld_daggrmat_nosmth_asb.F90
 mlprec/mld_daggrmat_smth_asb.F90
 mlprec/mld_saggrmat_nosmth_asb.F90
 mlprec/mld_saggrmat_smth_asb.F90
 mlprec/mld_zaggrmat_nosmth_asb.F90
 mlprec/mld_zaggrmat_smth_asb.F90
 tests/fileread/cf_sample.f90
 tests/fileread/df_sample.f90
 tests/fileread/sf_sample.f90
 tests/fileread/zf_sample.f90
 tests/pdegen/runs/ppde.inp


Fixup descriptor for replicated index space construction.
stopcriterion
Salvatore Filippone 14 years ago
parent f89f15b162
commit 44c29297ad
13 changed files with 263 additions and 208 deletions
Show Stats Download Patch File Download Diff File

1
mlprec/mld_caggrmat_nosmth_asb.F90
Unescape Escape View File

@ -183,6 +183,7 @@ subroutine mld_caggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,p,info)
if (p%iprcparm(mld_coarse_mat_) == mld_repl_mat_) then if (p%iprcparm(mld_coarse_mat_) == mld_repl_mat_) then
call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.) call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.)
if (info == psb_success_) call psb_cdasb(p%desc_ac,info)
if(info /= psb_success_) then if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_cdall') call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_cdall')
goto 9999 goto 9999

10
mlprec/mld_caggrmat_smth_asb.F90
Unescape Escape View File

@ -572,9 +572,11 @@ subroutine mld_caggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,p,info)
! !
! !
call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.) call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.)
call psb_gather(p%ac,b,p%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.) if (info == psb_success_) call psb_cdasb(p%desc_ac,info)
if(info /= psb_success_) goto 9999 if (info == psb_success_) &
if(info /= psb_success_) goto 9999 & call psb_gather(p%ac,b,p%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.)
if (info /= psb_success_) goto 9999
deallocate(nzbr,idisp,stat=info) deallocate(nzbr,idisp,stat=info)
if (info /= psb_success_) then if (info /= psb_success_) then
@ -608,6 +610,7 @@ subroutine mld_caggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,p,info)
! !
! !
call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.) call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.)
if (info == psb_success_) call psb_cdasb(p%desc_ac,info)
if(info /= psb_success_) then if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_cdall') call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_cdall')
goto 9999 goto 9999
@ -615,7 +618,6 @@ subroutine mld_caggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,p,info)
call psb_gather(p%ac,b,p%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.) call psb_gather(p%ac,b,p%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.)
if(info /= psb_success_) goto 9999 if(info /= psb_success_) goto 9999
deallocate(nzbr,idisp,stat=info) deallocate(nzbr,idisp,stat=info)
if (info /= psb_success_) then if (info /= psb_success_) then
info = psb_err_alloc_dealloc_ info = psb_err_alloc_dealloc_

1
mlprec/mld_daggrmat_nosmth_asb.F90
Unescape Escape View File

@ -183,6 +183,7 @@ subroutine mld_daggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,p,info)
if (p%iprcparm(mld_coarse_mat_) == mld_repl_mat_) then if (p%iprcparm(mld_coarse_mat_) == mld_repl_mat_) then
call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.) call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.)
if (info == psb_success_) call psb_cdasb(p%desc_ac,info)
if(info /= psb_success_) then if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_cdall') call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_cdall')
goto 9999 goto 9999

8
mlprec/mld_daggrmat_smth_asb.F90
Unescape Escape View File

@ -572,10 +572,11 @@ subroutine mld_daggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,p,info)
! !
! !
call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.) call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.)
call psb_gather(p%ac,b,p%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.) if (info == psb_success_) call psb_cdasb(p%desc_ac,info)
if (info == psb_success_) &
& call psb_gather(p%ac,b,p%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.)
if(info /= psb_success_) goto 9999 if (info /= psb_success_) goto 9999
if(info /= psb_success_) goto 9999
deallocate(nzbr,idisp,stat=info) deallocate(nzbr,idisp,stat=info)
if (info /= psb_success_) then if (info /= psb_success_) then
@ -609,6 +610,7 @@ subroutine mld_daggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,p,info)
! !
! !
call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.) call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.)
if (info == psb_success_) call psb_cdasb(p%desc_ac,info)
if(info /= psb_success_) then if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_cdall') call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_cdall')
goto 9999 goto 9999

1
mlprec/mld_saggrmat_nosmth_asb.F90
Unescape Escape View File

@ -183,6 +183,7 @@ subroutine mld_saggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,p,info)
if (p%iprcparm(mld_coarse_mat_) == mld_repl_mat_) then if (p%iprcparm(mld_coarse_mat_) == mld_repl_mat_) then
call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.) call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.)
if (info == psb_success_) call psb_cdasb(p%desc_ac,info)
if(info /= psb_success_) then if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_cdall') call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_cdall')
goto 9999 goto 9999

10
mlprec/mld_saggrmat_smth_asb.F90
Unescape Escape View File

@ -572,9 +572,11 @@ subroutine mld_saggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,p,info)
! !
! !
call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.) call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.)
call psb_gather(p%ac,b,p%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.) if (info == psb_success_) call psb_cdasb(p%desc_ac,info)
if(info /= psb_success_) goto 9999 if (info == psb_success_) &
if(info /= psb_success_) goto 9999 & call psb_gather(p%ac,b,p%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.)
if (info /= psb_success_) goto 9999
deallocate(nzbr,idisp,stat=info) deallocate(nzbr,idisp,stat=info)
if (info /= psb_success_) then if (info /= psb_success_) then
@ -608,6 +610,7 @@ subroutine mld_saggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,p,info)
! !
! !
call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.) call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.)
if (info == psb_success_) call psb_cdasb(p%desc_ac,info)
if(info /= psb_success_) then if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_cdall') call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_cdall')
goto 9999 goto 9999
@ -615,7 +618,6 @@ subroutine mld_saggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,p,info)
call psb_gather(p%ac,b,p%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.) call psb_gather(p%ac,b,p%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.)
if(info /= psb_success_) goto 9999 if(info /= psb_success_) goto 9999
deallocate(nzbr,idisp,stat=info) deallocate(nzbr,idisp,stat=info)
if (info /= psb_success_) then if (info /= psb_success_) then
info = psb_err_alloc_dealloc_ info = psb_err_alloc_dealloc_

1
mlprec/mld_zaggrmat_nosmth_asb.F90
Unescape Escape View File

@ -183,6 +183,7 @@ subroutine mld_zaggrmat_nosmth_asb(a,desc_a,ilaggr,nlaggr,p,info)
if (p%iprcparm(mld_coarse_mat_) == mld_repl_mat_) then if (p%iprcparm(mld_coarse_mat_) == mld_repl_mat_) then
call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.) call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.)
if (info == psb_success_) call psb_cdasb(p%desc_ac,info)
if(info /= psb_success_) then if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_cdall') call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_cdall')
goto 9999 goto 9999

8
mlprec/mld_zaggrmat_smth_asb.F90
Unescape Escape View File

@ -572,10 +572,11 @@ subroutine mld_zaggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,p,info)
! !
! !
call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.) call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.)
call psb_gather(p%ac,b,p%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.) if (info == psb_success_) call psb_cdasb(p%desc_ac,info)
if (info == psb_success_) &
& call psb_gather(p%ac,b,p%desc_ac,info,dupl=psb_dupl_add_,keeploc=.false.)
if(info /= psb_success_) goto 9999 if (info /= psb_success_) goto 9999
if(info /= psb_success_) goto 9999
deallocate(nzbr,idisp,stat=info) deallocate(nzbr,idisp,stat=info)
if (info /= psb_success_) then if (info /= psb_success_) then
@ -609,6 +610,7 @@ subroutine mld_zaggrmat_smth_asb(a,desc_a,ilaggr,nlaggr,p,info)
! !
! !
call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.) call psb_cdall(ictxt,p%desc_ac,info,mg=ntaggr,repl=.true.)
if (info == psb_success_) call psb_cdasb(p%desc_ac,info)
if(info /= psb_success_) then if(info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_cdall') call psb_errpush(psb_err_from_subroutine_,name,a_err='psb_cdall')
goto 9999 goto 9999

116
tests/fileread/cf_sample.f90
Unescape Escape View File

@ -52,11 +52,13 @@ program cf_sample
character(len=20) :: descr ! verbose description of the prec character(len=20) :: descr ! verbose description of the prec
character(len=10) :: prec ! overall prectype character(len=10) :: prec ! overall prectype
integer :: novr ! number of overlap layers integer :: novr ! number of overlap layers
integer :: jsweeps ! Jacobi/smoother sweeps
character(len=16) :: restr ! restriction over application of AS character(len=16) :: restr ! restriction over application of AS
character(len=16) :: prol ! prolongation over application of AS character(len=16) :: prol ! prolongation over application of AS
character(len=16) :: solve ! factorization type: ILU, SuperLU, UMFPACK character(len=16) :: solve ! factorization type: ILU, SuperLU, UMFPACK
integer :: fill ! fillin for factorization integer :: fill ! fillin for factorization
real(psb_spk_) :: thr ! threshold for fact. ILU(T) real(psb_spk_) :: thr ! threshold for fact. ILU(T)
character(len=16) :: smther ! Smoother
integer :: nlev ! number of levels in multilevel prec. integer :: nlev ! number of levels in multilevel prec.
character(len=16) :: aggrkind ! smoothed, raw aggregation character(len=16) :: aggrkind ! smoothed, raw aggregation
character(len=16) :: aggr_alg ! aggregation algorithm (currently only decoupled) character(len=16) :: aggr_alg ! aggregation algorithm (currently only decoupled)
@ -153,24 +155,24 @@ program cf_sample
case default case default
info = -1 info = -1
write(0,*) 'Wrong choice for fileformat ', filefmt write(psb_err_unit,*) 'Wrong choice for fileformat ', filefmt
end select end select
if (info /= psb_success_) then if (info /= psb_success_) then
write(0,*) 'Error while reading input matrix ' write(psb_err_unit,*) 'Error while reading input matrix '
call psb_abort(ictxt) call psb_abort(ictxt)
end if end if
m_problem = aux_a%m m_problem = aux_a%get_nrows()
call psb_bcast(ictxt,m_problem) call psb_bcast(ictxt,m_problem)
! At this point aux_b may still be unallocated ! At this point aux_b may still be unallocated
if (psb_size(aux_b,dim=1) == m_problem) then if (psb_size(aux_b,dim=1) == m_problem) then
! if any rhs were present, broadcast the first one ! if any rhs were present, broadcast the first one
write(0,'("Ok, got an rhs ")') write(psb_err_unit,'("Ok, got an rhs ")')
b_col_glob =>aux_b(:,1) b_col_glob =>aux_b(:,1)
else else
write(*,'("Generating an rhs...")') write(psb_out_unit,'("Generating an rhs...")')
write(*,'(" ")') write(psb_out_unit,'(" ")')
call psb_realloc(m_problem,1,aux_b,ircode) call psb_realloc(m_problem,1,aux_b,ircode)
if (ircode /= 0) then if (ircode /= 0) then
call psb_errpush(psb_err_alloc_dealloc_,name) call psb_errpush(psb_err_alloc_dealloc_,name)
@ -197,7 +199,7 @@ program cf_sample
! switch over different partition types ! switch over different partition types
if (ipart == 0) then if (ipart == 0) then
call psb_barrier(ictxt) call psb_barrier(ictxt)
if (iam == psb_root_) write(*,'("Partition type: block")') if (iam == psb_root_) write(psb_out_unit,'("Partition type: block")')
allocate(ivg(m_problem),ipv(np)) allocate(ivg(m_problem),ipv(np))
do i=1,m_problem do i=1,m_problem
call part_block(i,m_problem,np,ipv,nv) call part_block(i,m_problem,np,ipv,nv)
@ -207,19 +209,19 @@ program cf_sample
& desc_a,b_col_glob,b_col,info,fmt=afmt,v=ivg) & desc_a,b_col_glob,b_col,info,fmt=afmt,v=ivg)
else if (ipart == 2) then else if (ipart == 2) then
if (iam == psb_root_) then if (iam == psb_root_) then
write(*,'("Partition type: graph")') write(psb_out_unit,'("Partition type: graph")')
write(*,'(" ")') write(psb_out_unit,'(" ")')
! write(0,'("Build type: graph")') ! write(psb_err_unit,'("Build type: graph")')
call build_mtpart(aux_a%m,aux_a%fida,aux_a%ia1,aux_a%ia2,np) call build_mtpart(aux_a,np)
endif endif
call psb_barrier(ictxt) !!$ call psb_barrier(ictxt)
call distr_mtpart(psb_root_,ictxt) call distr_mtpart(psb_root_,ictxt)
call getv_mtpart(ivg) call getv_mtpart(ivg)
call psb_matdist(aux_a, a, ictxt, & call psb_matdist(aux_a, a, ictxt, &
& desc_a,b_col_glob,b_col,info,fmt=afmt,v=ivg) & desc_a,b_col_glob,b_col,info,fmt=afmt,v=ivg)
else else
if (iam == psb_root_) write(*,'("Partition type: block")') if (iam == psb_root_) write(psb_out_unit,'("Partition type: block")')
call psb_matdist(aux_a, a, ictxt, & call psb_matdist(aux_a, a, ictxt, &
& desc_a,b_col_glob,b_col,info,fmt=afmt,parts=part_block) & desc_a,b_col_glob,b_col,info,fmt=afmt,parts=part_block)
end if end if
@ -235,27 +237,25 @@ program cf_sample
call psb_amx(ictxt, t2) call psb_amx(ictxt, t2)
if (iam == psb_root_) then if (iam == psb_root_) then
write(*,'(" ")') write(psb_out_unit,'(" ")')
write(*,'("Time to read and partition matrix : ",es12.5)')t2 write(psb_out_unit,'("Time to read and partition matrix : ",es12.5)')t2
write(*,'(" ")') write(psb_out_unit,'(" ")')
write(*,*) 'Preconditioner: ',prec_choice%descr write(psb_out_unit,*) 'Preconditioner: ',prec_choice%descr
end if end if
! !
if (psb_toupper(prec_choice%prec) == 'ML') then if (psb_toupper(prec_choice%prec) == 'ML') then
nlv = prec_choice%nlev nlv = prec_choice%nlev
else call mld_precinit(prec,prec_choice%prec,info,nlev=nlv)
nlv = 1 call mld_precset(prec,mld_smoother_type_, prec_choice%smther, info)
end if call mld_precset(prec,mld_smoother_sweeps_, prec_choice%jsweeps, info)
call mld_precinit(prec,prec_choice%prec,info,nlev=nlv) call mld_precset(prec,mld_sub_ovr_, prec_choice%novr, info)
call mld_precset(prec,mld_sub_ovr_, prec_choice%novr, info) call mld_precset(prec,mld_sub_restr_, prec_choice%restr, info)
call mld_precset(prec,mld_sub_restr_, prec_choice%restr,info) call mld_precset(prec,mld_sub_prol_, prec_choice%prol, info)
call mld_precset(prec,mld_sub_prol_, prec_choice%prol, info) call mld_precset(prec,mld_sub_solve_, prec_choice%solve, info)
call mld_precset(prec,mld_sub_solve_, prec_choice%solve,info) call mld_precset(prec,mld_sub_fillin_, prec_choice%fill, info)
call mld_precset(prec,mld_sub_fillin_,prec_choice%fill,info) call mld_precset(prec,mld_sub_iluthrs_, prec_choice%thr, info)
call mld_precset(prec,mld_sub_iluthrs_, prec_choice%thr, info)
if (psb_toupper(prec_choice%prec) == 'ML') then
call mld_precset(prec,mld_aggr_kind_, prec_choice%aggrkind,info) call mld_precset(prec,mld_aggr_kind_, prec_choice%aggrkind,info)
call mld_precset(prec,mld_aggr_alg_, prec_choice%aggr_alg,info) call mld_precset(prec,mld_aggr_alg_, prec_choice%aggr_alg,info)
call mld_precset(prec,mld_ml_type_, prec_choice%mltype, info) call mld_precset(prec,mld_ml_type_, prec_choice%mltype, info)
@ -267,6 +267,16 @@ program cf_sample
call mld_precset(prec,mld_coarse_fillin_, prec_choice%cfill, info) call mld_precset(prec,mld_coarse_fillin_, prec_choice%cfill, info)
call mld_precset(prec,mld_coarse_iluthrs_, prec_choice%cthres, info) call mld_precset(prec,mld_coarse_iluthrs_, prec_choice%cthres, info)
call mld_precset(prec,mld_coarse_sweeps_, prec_choice%cjswp, info) call mld_precset(prec,mld_coarse_sweeps_, prec_choice%cjswp, info)
else
nlv = 1
call mld_precinit(prec,prec_choice%prec,info)
call mld_precset(prec,mld_smoother_sweeps_, prec_choice%jsweeps, info)
call mld_precset(prec,mld_sub_ovr_, prec_choice%novr, info)
call mld_precset(prec,mld_sub_restr_, prec_choice%restr, info)
call mld_precset(prec,mld_sub_prol_, prec_choice%prol, info)
call mld_precset(prec,mld_sub_solve_, prec_choice%solve, info)
call mld_precset(prec,mld_sub_fillin_, prec_choice%fill, info)
call mld_precset(prec,mld_sub_iluthrs_, prec_choice%thr, info)
end if end if
! building the preconditioner ! building the preconditioner
@ -278,12 +288,11 @@ program cf_sample
goto 9999 goto 9999
end if end if
call psb_amx(ictxt, tprec) call psb_amx(ictxt, tprec)
if(iam == psb_root_) then if(iam == psb_root_) then
write(*,'("Preconditioner time: ",es12.5)')tprec write(psb_out_unit,'("Preconditioner time: ",es12.5)')tprec
write(*,'(" ")') write(psb_out_unit,'(" ")')
end if end if
iparm = 0 iparm = 0
@ -308,36 +317,37 @@ program cf_sample
call psb_sum(ictxt,precsize) call psb_sum(ictxt,precsize)
if (iam == psb_root_) then if (iam == psb_root_) then
call mld_precdescr(prec,info) call mld_precdescr(prec,info)
write(*,'("Matrix: ",a)')mtrx_file write(psb_out_unit,'("Matrix: ",a)')mtrx_file
write(*,'("Computed solution on ",i8," processors")')np write(psb_out_unit,'("Computed solution on ",i8," processors")')np
write(*,'("Iterations to convergence : ",i6)')iter write(psb_out_unit,'("Iterations to convergence : ",i6)')iter
write(*,'("Error estimate on exit : ",es12.5)')err write(psb_out_unit,'("Error estimate on exit : ",es12.5)')err
write(*,'("Time to buil prec. : ",es12.5)')tprec write(psb_out_unit,'("Time to buil prec. : ",es12.5)')tprec
write(*,'("Time to solve matrix : ",es12.5)')t2 write(psb_out_unit,'("Time to solve matrix : ",es12.5)')t2
write(*,'("Time per iteration : ",es12.5)')t2/(iter) write(psb_out_unit,'("Time per iteration : ",es12.5)')t2/(iter)
write(*,'("Total time : ",es12.5)')t2+tprec write(psb_out_unit,'("Total time : ",es12.5)')t2+tprec
write(*,'("Residual norm 2 : ",es12.5)')resmx write(psb_out_unit,'("Residual norm 2 : ",es12.5)')resmx
write(*,'("Residual norm inf : ",es12.5)')resmxp write(psb_out_unit,'("Residual norm inf : ",es12.5)')resmxp
write(*,'("Total memory occupation for A : ",i12)')amatsize write(psb_out_unit,'("Total memory occupation for A : ",i12)')amatsize
write(*,'("Total memory occupation for DESC_A : ",i12)')descsize write(psb_out_unit,'("Total memory occupation for DESC_A : ",i12)')descsize
write(*,'("Total memory occupation for PREC : ",i12)')precsize write(psb_out_unit,'("Total memory occupation for PREC : ",i12)')precsize
end if end if
allocate(x_col_glob(m_problem),r_col_glob(m_problem),stat=ierr) allocate(x_col_glob(m_problem),r_col_glob(m_problem),stat=ierr)
if (ierr /= 0) then if (ierr /= 0) then
write(0,*) 'allocation error: no data collection' write(psb_err_unit,*) 'allocation error: no data collection'
else else
call psb_gather(x_col_glob,x_col,desc_a,info,root=psb_root_) call psb_gather(x_col_glob,x_col,desc_a,info,root=psb_root_)
call psb_gather(r_col_glob,r_col,desc_a,info,root=psb_root_) call psb_gather(r_col_glob,r_col,desc_a,info,root=psb_root_)
if (iam == psb_root_) then if (iam == psb_root_) then
write(0,'(" ")') write(psb_err_unit,'(" ")')
write(0,'("Saving x on file")') write(psb_err_unit,'("Saving x on file")')
write(20,*) 'matrix: ',mtrx_file write(20,*) 'matrix: ',mtrx_file
write(20,*) 'computed solution on ',np,' processors.' write(20,*) 'computed solution on ',np,' processors.'
write(20,*) 'iterations to convergence: ',iter write(20,*) 'iterations to convergence: ',iter
write(20,*) 'error estimate (infinity norm) on exit:', & write(20,*) 'error estimate (infinity norm) on exit:', &
& ' ||r||/(||a||||x||+||b||) = ',err & ' ||r||/(||a||||x||+||b||) = ',err
write(20,*) 'max residual = ',resmx, resmxp write(20,'("Residual norm 2 : ",es12.5)')resmx
write(20,'("Residual norm inf : ",es12.5)')resmxp
write(20,'(a8,4(2x,a20))') 'I','X(I)','R(I)','B(I)' write(20,'(a8,4(2x,a20))') 'I','X(I)','R(I)','B(I)'
do i=1,m_problem do i=1,m_problem
write(20,998) i,x_col_glob(i),r_col_glob(i),b_col_glob(i) write(20,998) i,x_col_glob(i),r_col_glob(i),b_col_glob(i)
@ -401,7 +411,9 @@ contains
call read_data(prec%solve,5) ! Factorization type: ILU, SuperLU, UMFPACK. call read_data(prec%solve,5) ! Factorization type: ILU, SuperLU, UMFPACK.
call read_data(prec%fill,5) ! Fill-in for factorization call read_data(prec%fill,5) ! Fill-in for factorization
call read_data(prec%thr,5) ! Threshold for fact. ILU(T) call read_data(prec%thr,5) ! Threshold for fact. ILU(T)
call read_data(prec%jsweeps,5) ! Jacobi sweeps for PJAC
if (psb_toupper(prec%prec) == 'ML') then if (psb_toupper(prec%prec) == 'ML') then
call read_data(prec%smther,5) ! Smoother type.
call read_data(prec%nlev,5) ! Number of levels in multilevel prec. call read_data(prec%nlev,5) ! Number of levels in multilevel prec.
call read_data(prec%aggrkind,5) ! smoothed/raw aggregatin call read_data(prec%aggrkind,5) ! smoothed/raw aggregatin
call read_data(prec%aggr_alg,5) ! local or global aggregation call read_data(prec%aggr_alg,5) ! local or global aggregation
@ -437,7 +449,9 @@ contains
call psb_bcast(icontxt,prec%solve) ! Factorization type: ILU, SuperLU, UMFPACK. call psb_bcast(icontxt,prec%solve) ! Factorization type: ILU, SuperLU, UMFPACK.
call psb_bcast(icontxt,prec%fill) ! Fill-in for factorization call psb_bcast(icontxt,prec%fill) ! Fill-in for factorization
call psb_bcast(icontxt,prec%thr) ! Threshold for fact. ILU(T) call psb_bcast(icontxt,prec%thr) ! Threshold for fact. ILU(T)
call psb_bcast(icontxt,prec%jsweeps) ! Jacobi sweeps
if (psb_toupper(prec%prec) == 'ML') then if (psb_toupper(prec%prec) == 'ML') then
call psb_bcast(icontxt,prec%smther) ! Smoother type.
call psb_bcast(icontxt,prec%nlev) ! Number of levels in multilevel prec. call psb_bcast(icontxt,prec%nlev) ! Number of levels in multilevel prec.
call psb_bcast(icontxt,prec%aggrkind) ! smoothed/raw aggregatin call psb_bcast(icontxt,prec%aggrkind) ! smoothed/raw aggregatin
call psb_bcast(icontxt,prec%aggr_alg) ! local or global aggregation call psb_bcast(icontxt,prec%aggr_alg) ! local or global aggregation

73
tests/fileread/df_sample.f90
Unescape Escape View File

@ -105,8 +105,8 @@ program df_sample
! other variables ! other variables
integer :: i,info,j,m_problem integer :: i,info,j,m_problem
integer :: internal, m,ii,nnzero integer :: internal, m,ii,nnzero
real(psb_dpk_) :: t1, t2, tprec, r_amax, b_amax,& real(psb_dpk_) :: t1, t2, tprec
&scale,resmx,resmxp real(psb_dpk_) :: r_amax, b_amax, scale,resmx,resmxp
integer :: nrhs, nrow, n_row, dim, nv, ne integer :: nrhs, nrow, n_row, dim, nv, ne
integer, allocatable :: ivg(:), ipv(:) integer, allocatable :: ivg(:), ipv(:)
@ -155,10 +155,10 @@ program df_sample
case default case default
info = -1 info = -1
write(0,*) 'Wrong choice for fileformat ', filefmt write(psb_err_unit,*) 'Wrong choice for fileformat ', filefmt
end select end select
if (info /= psb_success_) then if (info /= psb_success_) then
write(0,*) 'Error while reading input matrix ' write(psb_err_unit,*) 'Error while reading input matrix '
call psb_abort(ictxt) call psb_abort(ictxt)
end if end if
@ -168,11 +168,11 @@ program df_sample
! At this point aux_b may still be unallocated ! At this point aux_b may still be unallocated
if (psb_size(aux_b,dim=1) == m_problem) then if (psb_size(aux_b,dim=1) == m_problem) then
! if any rhs were present, broadcast the first one ! if any rhs were present, broadcast the first one
write(0,'("Ok, got an rhs ")') write(psb_err_unit,'("Ok, got an rhs ")')
b_col_glob =>aux_b(:,1) b_col_glob =>aux_b(:,1)
else else
write(*,'("Generating an rhs...")') write(psb_out_unit,'("Generating an rhs...")')
write(*,'(" ")') write(psb_out_unit,'(" ")')
call psb_realloc(m_problem,1,aux_b,ircode) call psb_realloc(m_problem,1,aux_b,ircode)
if (ircode /= 0) then if (ircode /= 0) then
call psb_errpush(psb_err_alloc_dealloc_,name) call psb_errpush(psb_err_alloc_dealloc_,name)
@ -199,7 +199,7 @@ program df_sample
! switch over different partition types ! switch over different partition types
if (ipart == 0) then if (ipart == 0) then
call psb_barrier(ictxt) call psb_barrier(ictxt)
if (iam == psb_root_) write(*,'("Partition type: block")') if (iam == psb_root_) write(psb_out_unit,'("Partition type: block")')
allocate(ivg(m_problem),ipv(np)) allocate(ivg(m_problem),ipv(np))
do i=1,m_problem do i=1,m_problem
call part_block(i,m_problem,np,ipv,nv) call part_block(i,m_problem,np,ipv,nv)
@ -209,9 +209,9 @@ program df_sample
& desc_a,b_col_glob,b_col,info,fmt=afmt,v=ivg) & desc_a,b_col_glob,b_col,info,fmt=afmt,v=ivg)
else if (ipart == 2) then else if (ipart == 2) then
if (iam == psb_root_) then if (iam == psb_root_) then
write(*,'("Partition type: graph")') write(psb_out_unit,'("Partition type: graph")')
write(*,'("")') write(psb_out_unit,'(" ")')
! write(0,'("Build type: graph")') ! write(psb_err_unit,'("Build type: graph")')
call build_mtpart(aux_a,np) call build_mtpart(aux_a,np)
endif endif
!!$ call psb_barrier(ictxt) !!$ call psb_barrier(ictxt)
@ -220,8 +220,8 @@ program df_sample
call psb_matdist(aux_a, a, ictxt, & call psb_matdist(aux_a, a, ictxt, &
& desc_a,b_col_glob,b_col,info,fmt=afmt,v=ivg) & desc_a,b_col_glob,b_col,info,fmt=afmt,v=ivg)
else else
if (iam == psb_root_) write(*,'("Partition type: block")') if (iam == psb_root_) write(psb_out_unit,'("Partition type: block")')
call psb_matdist(aux_a, a, ictxt, & call psb_matdist(aux_a, a, ictxt, &
& desc_a,b_col_glob,b_col,info,fmt=afmt,parts=part_block) & desc_a,b_col_glob,b_col,info,fmt=afmt,parts=part_block)
end if end if
@ -237,10 +237,10 @@ program df_sample
call psb_amx(ictxt, t2) call psb_amx(ictxt, t2)
if (iam == psb_root_) then if (iam == psb_root_) then
write(*,'(" ")') write(psb_out_unit,'(" ")')
write(*,'("Time to read and partition matrix : ",es12.5)')t2 write(psb_out_unit,'("Time to read and partition matrix : ",es12.5)')t2
write(*,'(" ")') write(psb_out_unit,'(" ")')
write(*,*) 'Preconditioner: ',prec_choice%descr write(psb_out_unit,*) 'Preconditioner: ',prec_choice%descr
end if end if
! !
@ -291,8 +291,8 @@ program df_sample
call psb_amx(ictxt, tprec) call psb_amx(ictxt, tprec)
if(iam == psb_root_) then if(iam == psb_root_) then
write(*,'("Preconditioner time: ",es12.5)')tprec write(psb_out_unit,'("Preconditioner time: ",es12.5)')tprec
write(*,'(" ")') write(psb_out_unit,'(" ")')
end if end if
iparm = 0 iparm = 0
@ -317,36 +317,37 @@ program df_sample
call psb_sum(ictxt,precsize) call psb_sum(ictxt,precsize)
if (iam == psb_root_) then if (iam == psb_root_) then
call mld_precdescr(prec,info) call mld_precdescr(prec,info)
write(*,'("Matrix: ",a)')mtrx_file write(psb_out_unit,'("Matrix: ",a)')mtrx_file
write(*,'("Computed solution on ",i8," processors")')np write(psb_out_unit,'("Computed solution on ",i8," processors")')np
write(*,'("Iterations to convergence : ",i6)')iter write(psb_out_unit,'("Iterations to convergence : ",i6)')iter
write(*,'("Error estimate on exit : ",es12.5)')err write(psb_out_unit,'("Error estimate on exit : ",es12.5)')err
write(*,'("Time to buil prec. : ",es12.5)')tprec write(psb_out_unit,'("Time to buil prec. : ",es12.5)')tprec
write(*,'("Time to solve matrix : ",es12.5)')t2 write(psb_out_unit,'("Time to solve matrix : ",es12.5)')t2
write(*,'("Time per iteration : ",es12.5)')t2/(iter) write(psb_out_unit,'("Time per iteration : ",es12.5)')t2/(iter)
write(*,'("Total time : ",es12.5)')t2+tprec write(psb_out_unit,'("Total time : ",es12.5)')t2+tprec
write(*,'("Residual norm 2 : ",es12.5)')resmx write(psb_out_unit,'("Residual norm 2 : ",es12.5)')resmx
write(*,'("Residual norm inf : ",es12.5)')resmxp write(psb_out_unit,'("Residual norm inf : ",es12.5)')resmxp
write(*,'("Total memory occupation for A : ",i12)')amatsize write(psb_out_unit,'("Total memory occupation for A : ",i12)')amatsize
write(*,'("Total memory occupation for DESC_A : ",i12)')descsize write(psb_out_unit,'("Total memory occupation for DESC_A : ",i12)')descsize
write(*,'("Total memory occupation for PREC : ",i12)')precsize write(psb_out_unit,'("Total memory occupation for PREC : ",i12)')precsize
end if end if
allocate(x_col_glob(m_problem),r_col_glob(m_problem),stat=ierr) allocate(x_col_glob(m_problem),r_col_glob(m_problem),stat=ierr)
if (ierr /= 0) then if (ierr /= 0) then
write(0,*) 'allocation error: no data collection' write(psb_err_unit,*) 'allocation error: no data collection'
else else
call psb_gather(x_col_glob,x_col,desc_a,info,root=psb_root_) call psb_gather(x_col_glob,x_col,desc_a,info,root=psb_root_)
call psb_gather(r_col_glob,r_col,desc_a,info,root=psb_root_) call psb_gather(r_col_glob,r_col,desc_a,info,root=psb_root_)
if (iam == psb_root_) then if (iam == psb_root_) then
write(0,'(" ")') write(psb_err_unit,'(" ")')
write(0,'("Saving x on file")') write(psb_err_unit,'("Saving x on file")')
write(20,*) 'matrix: ',mtrx_file write(20,*) 'matrix: ',mtrx_file
write(20,*) 'computed solution on ',np,' processors.' write(20,*) 'computed solution on ',np,' processors.'
write(20,*) 'iterations to convergence: ',iter write(20,*) 'iterations to convergence: ',iter
write(20,*) 'error estimate (infinity norm) on exit:', & write(20,*) 'error estimate (infinity norm) on exit:', &
& ' ||r||/(||a||||x||+||b||) = ',err & ' ||r||/(||a||||x||+||b||) = ',err
write(20,*) 'max residual = ',resmx, resmxp write(20,'("Residual norm 2 : ",es12.5)')resmx
write(20,'("Residual norm inf : ",es12.5)')resmxp
write(20,'(a8,4(2x,a20))') 'I','X(I)','R(I)','B(I)' write(20,'(a8,4(2x,a20))') 'I','X(I)','R(I)','B(I)'
do i=1,m_problem do i=1,m_problem
write(20,998) i,x_col_glob(i),r_col_glob(i),b_col_glob(i) write(20,998) i,x_col_glob(i),r_col_glob(i),b_col_glob(i)

116
tests/fileread/sf_sample.f90
Unescape Escape View File

@ -52,11 +52,13 @@ program sf_sample
character(len=20) :: descr ! verbose description of the prec character(len=20) :: descr ! verbose description of the prec
character(len=10) :: prec ! overall prectype character(len=10) :: prec ! overall prectype
integer :: novr ! number of overlap layers integer :: novr ! number of overlap layers
integer :: jsweeps ! Jacobi/smoother sweeps
character(len=16) :: restr ! restriction over application of AS character(len=16) :: restr ! restriction over application of AS
character(len=16) :: prol ! prolongation over application of AS character(len=16) :: prol ! prolongation over application of AS
character(len=16) :: solve ! factorization type: ILU, SuperLU, UMFPACK character(len=16) :: solve ! factorization type: ILU, SuperLU, UMFPACK
integer :: fill ! fillin for factorization integer :: fill ! fillin for factorization
real(psb_spk_) :: thr ! threshold for fact. ILU(T) real(psb_spk_) :: thr ! threshold for fact. ILU(T)
character(len=16) :: smther ! Smoother
integer :: nlev ! number of levels in multilevel prec. integer :: nlev ! number of levels in multilevel prec.
character(len=16) :: aggrkind ! smoothed, raw aggregation character(len=16) :: aggrkind ! smoothed, raw aggregation
character(len=16) :: aggr_alg ! aggregation algorithm (currently only decoupled) character(len=16) :: aggr_alg ! aggregation algorithm (currently only decoupled)
@ -153,24 +155,24 @@ program sf_sample
case default case default
info = -1 info = -1
write(0,*) 'Wrong choice for fileformat ', filefmt write(psb_err_unit,*) 'Wrong choice for fileformat ', filefmt
end select end select
if (info /= psb_success_) then if (info /= psb_success_) then
write(0,*) 'Error while reading input matrix ' write(psb_err_unit,*) 'Error while reading input matrix '
call psb_abort(ictxt) call psb_abort(ictxt)
end if end if
m_problem = aux_a%m m_problem = aux_a%get_nrows()
call psb_bcast(ictxt,m_problem) call psb_bcast(ictxt,m_problem)
! At this point aux_b may still be unallocated ! At this point aux_b may still be unallocated
if (psb_size(aux_b,dim=1) == m_problem) then if (psb_size(aux_b,dim=1) == m_problem) then
! if any rhs were present, broadcast the first one ! if any rhs were present, broadcast the first one
write(0,'("Ok, got an rhs ")') write(psb_err_unit,'("Ok, got an rhs ")')
b_col_glob =>aux_b(:,1) b_col_glob =>aux_b(:,1)
else else
write(*,'("Generating an rhs...")') write(psb_out_unit,'("Generating an rhs...")')
write(*,'(" ")') write(psb_out_unit,'(" ")')
call psb_realloc(m_problem,1,aux_b,ircode) call psb_realloc(m_problem,1,aux_b,ircode)
if (ircode /= 0) then if (ircode /= 0) then
call psb_errpush(psb_err_alloc_dealloc_,name) call psb_errpush(psb_err_alloc_dealloc_,name)
@ -197,7 +199,7 @@ program sf_sample
! switch over different partition types ! switch over different partition types
if (ipart == 0) then if (ipart == 0) then
call psb_barrier(ictxt) call psb_barrier(ictxt)
if (iam == psb_root_) write(*,'("Partition type: block")') if (iam == psb_root_) write(psb_out_unit,'("Partition type: block")')
allocate(ivg(m_problem),ipv(np)) allocate(ivg(m_problem),ipv(np))
do i=1,m_problem do i=1,m_problem
call part_block(i,m_problem,np,ipv,nv) call part_block(i,m_problem,np,ipv,nv)
@ -207,19 +209,19 @@ program sf_sample
& desc_a,b_col_glob,b_col,info,fmt=afmt,v=ivg) & desc_a,b_col_glob,b_col,info,fmt=afmt,v=ivg)
else if (ipart == 2) then else if (ipart == 2) then
if (iam == psb_root_) then if (iam == psb_root_) then
write(*,'("Partition type: graph")') write(psb_out_unit,'("Partition type: graph")')
write(*,'(" ")') write(psb_out_unit,'(" ")')
! write(0,'("Build type: graph")') ! write(psb_err_unit,'("Build type: graph")')
call build_mtpart(aux_a%m,aux_a%fida,aux_a%ia1,aux_a%ia2,np) call build_mtpart(aux_a,np)
endif endif
call psb_barrier(ictxt) !!$ call psb_barrier(ictxt)
call distr_mtpart(psb_root_,ictxt) call distr_mtpart(psb_root_,ictxt)
call getv_mtpart(ivg) call getv_mtpart(ivg)
call psb_matdist(aux_a, a, ictxt, & call psb_matdist(aux_a, a, ictxt, &
& desc_a,b_col_glob,b_col,info,fmt=afmt,v=ivg) & desc_a,b_col_glob,b_col,info,fmt=afmt,v=ivg)
else else
if (iam == psb_root_) write(*,'("Partition type: block")') if (iam == psb_root_) write(psb_out_unit,'("Partition type: block")')
call psb_matdist(aux_a, a, ictxt, & call psb_matdist(aux_a, a, ictxt, &
& desc_a,b_col_glob,b_col,info,fmt=afmt,parts=part_block) & desc_a,b_col_glob,b_col,info,fmt=afmt,parts=part_block)
end if end if
@ -235,27 +237,25 @@ program sf_sample
call psb_amx(ictxt, t2) call psb_amx(ictxt, t2)
if (iam == psb_root_) then if (iam == psb_root_) then
write(*,'(" ")') write(psb_out_unit,'(" ")')
write(*,'("Time to read and partition matrix : ",es12.5)')t2 write(psb_out_unit,'("Time to read and partition matrix : ",es12.5)')t2
write(*,'(" ")') write(psb_out_unit,'(" ")')
write(*,*) 'Preconditioner: ',prec_choice%descr write(psb_out_unit,*) 'Preconditioner: ',prec_choice%descr
end if end if
! !
if (psb_toupper(prec_choice%prec) == 'ML') then if (psb_toupper(prec_choice%prec) == 'ML') then
nlv = prec_choice%nlev nlv = prec_choice%nlev
else call mld_precinit(prec,prec_choice%prec,info,nlev=nlv)
nlv = 1 call mld_precset(prec,mld_smoother_type_, prec_choice%smther, info)
end if call mld_precset(prec,mld_smoother_sweeps_, prec_choice%jsweeps, info)
call mld_precinit(prec,prec_choice%prec,info,nlev=nlv) call mld_precset(prec,mld_sub_ovr_, prec_choice%novr, info)
call mld_precset(prec,mld_sub_ovr_, prec_choice%novr, info) call mld_precset(prec,mld_sub_restr_, prec_choice%restr, info)
call mld_precset(prec,mld_sub_restr_, prec_choice%restr,info) call mld_precset(prec,mld_sub_prol_, prec_choice%prol, info)
call mld_precset(prec,mld_sub_prol_, prec_choice%prol, info) call mld_precset(prec,mld_sub_solve_, prec_choice%solve, info)
call mld_precset(prec,mld_sub_solve_, prec_choice%solve,info) call mld_precset(prec,mld_sub_fillin_, prec_choice%fill, info)
call mld_precset(prec,mld_sub_fillin_,prec_choice%fill,info) call mld_precset(prec,mld_sub_iluthrs_, prec_choice%thr, info)
call mld_precset(prec,mld_sub_iluthrs_, prec_choice%thr, info)
if (psb_toupper(prec_choice%prec) == 'ML') then
call mld_precset(prec,mld_aggr_kind_, prec_choice%aggrkind,info) call mld_precset(prec,mld_aggr_kind_, prec_choice%aggrkind,info)
call mld_precset(prec,mld_aggr_alg_, prec_choice%aggr_alg,info) call mld_precset(prec,mld_aggr_alg_, prec_choice%aggr_alg,info)
call mld_precset(prec,mld_ml_type_, prec_choice%mltype, info) call mld_precset(prec,mld_ml_type_, prec_choice%mltype, info)
@ -267,6 +267,16 @@ program sf_sample
call mld_precset(prec,mld_coarse_fillin_, prec_choice%cfill, info) call mld_precset(prec,mld_coarse_fillin_, prec_choice%cfill, info)
call mld_precset(prec,mld_coarse_iluthrs_, prec_choice%cthres, info) call mld_precset(prec,mld_coarse_iluthrs_, prec_choice%cthres, info)
call mld_precset(prec,mld_coarse_sweeps_, prec_choice%cjswp, info) call mld_precset(prec,mld_coarse_sweeps_, prec_choice%cjswp, info)
else
nlv = 1
call mld_precinit(prec,prec_choice%prec,info)
call mld_precset(prec,mld_smoother_sweeps_, prec_choice%jsweeps, info)
call mld_precset(prec,mld_sub_ovr_, prec_choice%novr, info)
call mld_precset(prec,mld_sub_restr_, prec_choice%restr, info)
call mld_precset(prec,mld_sub_prol_, prec_choice%prol, info)
call mld_precset(prec,mld_sub_solve_, prec_choice%solve, info)
call mld_precset(prec,mld_sub_fillin_, prec_choice%fill, info)
call mld_precset(prec,mld_sub_iluthrs_, prec_choice%thr, info)
end if end if
! building the preconditioner ! building the preconditioner
@ -278,12 +288,11 @@ program sf_sample
goto 9999 goto 9999
end if end if
call psb_amx(ictxt, tprec) call psb_amx(ictxt, tprec)
if(iam == psb_root_) then if(iam == psb_root_) then
write(*,'("Preconditioner time: ",es12.5)')tprec write(psb_out_unit,'("Preconditioner time: ",es12.5)')tprec
write(*,'(" ")') write(psb_out_unit,'(" ")')
end if end if
iparm = 0 iparm = 0
@ -308,36 +317,37 @@ program sf_sample
call psb_sum(ictxt,precsize) call psb_sum(ictxt,precsize)
if (iam == psb_root_) then if (iam == psb_root_) then
call mld_precdescr(prec,info) call mld_precdescr(prec,info)
write(*,'("Matrix: ",a)')mtrx_file write(psb_out_unit,'("Matrix: ",a)')mtrx_file
write(*,'("Computed solution on ",i8," processors")')np write(psb_out_unit,'("Computed solution on ",i8," processors")')np
write(*,'("Iterations to convergence : ",i6)')iter write(psb_out_unit,'("Iterations to convergence : ",i6)')iter
write(*,'("Error estimate on exit : ",es12.5)')err write(psb_out_unit,'("Error estimate on exit : ",es12.5)')err
write(*,'("Time to buil prec. : ",es12.5)')tprec write(psb_out_unit,'("Time to buil prec. : ",es12.5)')tprec
write(*,'("Time to solve matrix : ",es12.5)')t2 write(psb_out_unit,'("Time to solve matrix : ",es12.5)')t2
write(*,'("Time per iteration : ",es12.5)')t2/(iter) write(psb_out_unit,'("Time per iteration : ",es12.5)')t2/(iter)
write(*,'("Total time : ",es12.5)')t2+tprec write(psb_out_unit,'("Total time : ",es12.5)')t2+tprec
write(*,'("Residual norm 2 : ",es12.5)')resmx write(psb_out_unit,'("Residual norm 2 : ",es12.5)')resmx
write(*,'("Residual norm inf : ",es12.5)')resmxp write(psb_out_unit,'("Residual norm inf : ",es12.5)')resmxp
write(*,'("Total memory occupation for A : ",i12)')amatsize write(psb_out_unit,'("Total memory occupation for A : ",i12)')amatsize
write(*,'("Total memory occupation for DESC_A : ",i12)')descsize write(psb_out_unit,'("Total memory occupation for DESC_A : ",i12)')descsize
write(*,'("Total memory occupation for PREC : ",i12)')precsize write(psb_out_unit,'("Total memory occupation for PREC : ",i12)')precsize
end if end if
allocate(x_col_glob(m_problem),r_col_glob(m_problem),stat=ierr) allocate(x_col_glob(m_problem),r_col_glob(m_problem),stat=ierr)
if (ierr /= 0) then if (ierr /= 0) then
write(0,*) 'allocation error: no data collection' write(psb_err_unit,*) 'allocation error: no data collection'
else else
call psb_gather(x_col_glob,x_col,desc_a,info,root=psb_root_) call psb_gather(x_col_glob,x_col,desc_a,info,root=psb_root_)
call psb_gather(r_col_glob,r_col,desc_a,info,root=psb_root_) call psb_gather(r_col_glob,r_col,desc_a,info,root=psb_root_)
if (iam == psb_root_) then if (iam == psb_root_) then
write(0,'(" ")') write(psb_err_unit,'(" ")')
write(0,'("Saving x on file")') write(psb_err_unit,'("Saving x on file")')
write(20,*) 'matrix: ',mtrx_file write(20,*) 'matrix: ',mtrx_file
write(20,*) 'computed solution on ',np,' processors.' write(20,*) 'computed solution on ',np,' processors.'
write(20,*) 'iterations to convergence: ',iter write(20,*) 'iterations to convergence: ',iter
write(20,*) 'error estimate (infinity norm) on exit:', & write(20,*) 'error estimate (infinity norm) on exit:', &
& ' ||r||/(||a||||x||+||b||) = ',err & ' ||r||/(||a||||x||+||b||) = ',err
write(20,*) 'max residual = ',resmx, resmxp write(20,'("Residual norm 2 : ",es12.5)')resmx
write(20,'("Residual norm inf : ",es12.5)')resmxp
write(20,'(a8,4(2x,a20))') 'I','X(I)','R(I)','B(I)' write(20,'(a8,4(2x,a20))') 'I','X(I)','R(I)','B(I)'
do i=1,m_problem do i=1,m_problem
write(20,998) i,x_col_glob(i),r_col_glob(i),b_col_glob(i) write(20,998) i,x_col_glob(i),r_col_glob(i),b_col_glob(i)
@ -401,7 +411,9 @@ contains
call read_data(prec%solve,5) ! Factorization type: ILU, SuperLU, UMFPACK. call read_data(prec%solve,5) ! Factorization type: ILU, SuperLU, UMFPACK.
call read_data(prec%fill,5) ! Fill-in for factorization call read_data(prec%fill,5) ! Fill-in for factorization
call read_data(prec%thr,5) ! Threshold for fact. ILU(T) call read_data(prec%thr,5) ! Threshold for fact. ILU(T)
call read_data(prec%jsweeps,5) ! Jacobi sweeps for PJAC
if (psb_toupper(prec%prec) == 'ML') then if (psb_toupper(prec%prec) == 'ML') then
call read_data(prec%smther,5) ! Smoother type.
call read_data(prec%nlev,5) ! Number of levels in multilevel prec. call read_data(prec%nlev,5) ! Number of levels in multilevel prec.
call read_data(prec%aggrkind,5) ! smoothed/raw aggregatin call read_data(prec%aggrkind,5) ! smoothed/raw aggregatin
call read_data(prec%aggr_alg,5) ! local or global aggregation call read_data(prec%aggr_alg,5) ! local or global aggregation
@ -437,7 +449,9 @@ contains
call psb_bcast(icontxt,prec%solve) ! Factorization type: ILU, SuperLU, UMFPACK. call psb_bcast(icontxt,prec%solve) ! Factorization type: ILU, SuperLU, UMFPACK.
call psb_bcast(icontxt,prec%fill) ! Fill-in for factorization call psb_bcast(icontxt,prec%fill) ! Fill-in for factorization
call psb_bcast(icontxt,prec%thr) ! Threshold for fact. ILU(T) call psb_bcast(icontxt,prec%thr) ! Threshold for fact. ILU(T)
call psb_bcast(icontxt,prec%jsweeps) ! Jacobi sweeps
if (psb_toupper(prec%prec) == 'ML') then if (psb_toupper(prec%prec) == 'ML') then
call psb_bcast(icontxt,prec%smther) ! Smoother type.
call psb_bcast(icontxt,prec%nlev) ! Number of levels in multilevel prec. call psb_bcast(icontxt,prec%nlev) ! Number of levels in multilevel prec.
call psb_bcast(icontxt,prec%aggrkind) ! smoothed/raw aggregatin call psb_bcast(icontxt,prec%aggrkind) ! smoothed/raw aggregatin
call psb_bcast(icontxt,prec%aggr_alg) ! local or global aggregation call psb_bcast(icontxt,prec%aggr_alg) ! local or global aggregation

124
tests/fileread/zf_sample.f90
Unescape Escape View File

@ -52,11 +52,13 @@ program zf_sample
character(len=20) :: descr ! verbose description of the prec character(len=20) :: descr ! verbose description of the prec
character(len=10) :: prec ! overall prectype character(len=10) :: prec ! overall prectype
integer :: novr ! number of overlap layers integer :: novr ! number of overlap layers
integer :: jsweeps ! Jacobi/smoother sweeps
character(len=16) :: restr ! restriction over application of AS character(len=16) :: restr ! restriction over application of AS
character(len=16) :: prol ! prolongation over application of AS character(len=16) :: prol ! prolongation over application of AS
character(len=16) :: solve ! factorization type: ILU, SuperLU, UMFPACK character(len=16) :: solve ! factorization type: ILU, SuperLU, UMFPACK
integer :: fill ! fillin for factorization integer :: fill ! fillin for factorization
real(psb_dpk_) :: thr ! threshold for fact. ILU(T) real(psb_dpk_) :: thr ! threshold for fact. ILU(T)
character(len=16) :: smther ! Smoother
integer :: nlev ! number of levels in multilevel prec. integer :: nlev ! number of levels in multilevel prec.
character(len=16) :: aggrkind ! smoothed, raw aggregation character(len=16) :: aggrkind ! smoothed, raw aggregation
character(len=16) :: aggr_alg ! aggregation algorithm (currently only decoupled) character(len=16) :: aggr_alg ! aggregation algorithm (currently only decoupled)
@ -103,8 +105,8 @@ program zf_sample
! other variables ! other variables
integer :: i,info,j,m_problem integer :: i,info,j,m_problem
integer :: internal, m,ii,nnzero integer :: internal, m,ii,nnzero
real(psb_dpk_) :: t1, t2, tprec, r_amax, b_amax,& real(psb_dpk_) :: t1, t2, tprec
&scale,resmx,resmxp real(psb_dpk_) :: r_amax, b_amax, scale,resmx,resmxp
integer :: nrhs, nrow, n_row, dim, nv, ne integer :: nrhs, nrow, n_row, dim, nv, ne
integer, allocatable :: ivg(:), ipv(:) integer, allocatable :: ivg(:), ipv(:)
@ -145,32 +147,32 @@ program zf_sample
call mm_vet_read(aux_b,info,iunit=iunit,filename=rhs_file) call mm_vet_read(aux_b,info,iunit=iunit,filename=rhs_file)
end if end if
end if end if
case ('HB') case ('HB')
! For Harwell-Boeing we have a single file which may or may not ! For Harwell-Boeing we have a single file which may or may not
! contain an RHS. ! contain an RHS.
call hb_read(aux_a,info,iunit=iunit,b=aux_b,filename=mtrx_file) call hb_read(aux_a,info,iunit=iunit,b=aux_b,filename=mtrx_file)
case default case default
info = -1 info = -1
write(0,*) 'Wrong choice for fileformat ', filefmt write(psb_err_unit,*) 'Wrong choice for fileformat ', filefmt
end select end select
if (info /= psb_success_) then if (info /= psb_success_) then
write(0,*) 'Error while reading input matrix ' write(psb_err_unit,*) 'Error while reading input matrix '
call psb_abort(ictxt) call psb_abort(ictxt)
end if end if
m_problem = aux_a%m m_problem = aux_a%get_nrows()
call psb_bcast(ictxt,m_problem) call psb_bcast(ictxt,m_problem)
! At this point aux_b may still be unallocated ! At this point aux_b may still be unallocated
if (psb_size(aux_b,dim=1) == m_problem) then if (psb_size(aux_b,dim=1) == m_problem) then
! if any rhs were present, broadcast the first one ! if any rhs were present, broadcast the first one
write(0,'("Ok, got an rhs ")') write(psb_err_unit,'("Ok, got an rhs ")')
b_col_glob =>aux_b(:,1) b_col_glob =>aux_b(:,1)
else else
write(*,'("Generating an rhs...")') write(psb_out_unit,'("Generating an rhs...")')
write(*,'(" ")') write(psb_out_unit,'(" ")')
call psb_realloc(m_problem,1,aux_b,ircode) call psb_realloc(m_problem,1,aux_b,ircode)
if (ircode /= 0) then if (ircode /= 0) then
call psb_errpush(psb_err_alloc_dealloc_,name) call psb_errpush(psb_err_alloc_dealloc_,name)
@ -197,7 +199,7 @@ program zf_sample
! switch over different partition types ! switch over different partition types
if (ipart == 0) then if (ipart == 0) then
call psb_barrier(ictxt) call psb_barrier(ictxt)
if (iam == psb_root_) write(*,'("Partition type: block")') if (iam == psb_root_) write(psb_out_unit,'("Partition type: block")')
allocate(ivg(m_problem),ipv(np)) allocate(ivg(m_problem),ipv(np))
do i=1,m_problem do i=1,m_problem
call part_block(i,m_problem,np,ipv,nv) call part_block(i,m_problem,np,ipv,nv)
@ -207,19 +209,19 @@ program zf_sample
& desc_a,b_col_glob,b_col,info,fmt=afmt,v=ivg) & desc_a,b_col_glob,b_col,info,fmt=afmt,v=ivg)
else if (ipart == 2) then else if (ipart == 2) then
if (iam == psb_root_) then if (iam == psb_root_) then
write(*,'("Partition type: graph")') write(psb_out_unit,'("Partition type: graph")')
write(*,'(" ")') write(psb_out_unit,'(" ")')
! write(0,'("Build type: graph")') ! write(psb_err_unit,'("Build type: graph")')
call build_mtpart(aux_a%m,aux_a%fida,aux_a%ia1,aux_a%ia2,np) call build_mtpart(aux_a,np)
endif endif
call psb_barrier(ictxt) !!$ call psb_barrier(ictxt)
call distr_mtpart(psb_root_,ictxt) call distr_mtpart(psb_root_,ictxt)
call getv_mtpart(ivg) call getv_mtpart(ivg)
call psb_matdist(aux_a, a, ictxt, & call psb_matdist(aux_a, a, ictxt, &
& desc_a,b_col_glob,b_col,info,fmt=afmt,v=ivg) & desc_a,b_col_glob,b_col,info,fmt=afmt,v=ivg)
else else
if (iam == psb_root_) write(*,'("Partition type: block")') if (iam == psb_root_) write(psb_out_unit,'("Partition type: block")')
call psb_matdist(aux_a, a, ictxt, & call psb_matdist(aux_a, a, ictxt, &
& desc_a,b_col_glob,b_col,info,fmt=afmt,parts=part_block) & desc_a,b_col_glob,b_col,info,fmt=afmt,parts=part_block)
end if end if
@ -235,27 +237,25 @@ program zf_sample
call psb_amx(ictxt, t2) call psb_amx(ictxt, t2)
if (iam == psb_root_) then if (iam == psb_root_) then
write(*,'(" ")') write(psb_out_unit,'(" ")')
write(*,'("Time to read and partition matrix : ",es12.5)')t2 write(psb_out_unit,'("Time to read and partition matrix : ",es12.5)')t2
write(*,'(" ")') write(psb_out_unit,'(" ")')
write(*,*) 'Preconditioner: ',prec_choice%descr write(psb_out_unit,*) 'Preconditioner: ',prec_choice%descr
end if end if
! !
if (psb_toupper(prec_choice%prec) == 'ML') then if (psb_toupper(prec_choice%prec) == 'ML') then
nlv = prec_choice%nlev nlv = prec_choice%nlev
else call mld_precinit(prec,prec_choice%prec,info,nlev=nlv)
nlv = 1 call mld_precset(prec,mld_smoother_type_, prec_choice%smther, info)
end if call mld_precset(prec,mld_smoother_sweeps_, prec_choice%jsweeps, info)
call mld_precinit(prec,prec_choice%prec,info,nlev=nlv) call mld_precset(prec,mld_sub_ovr_, prec_choice%novr, info)
call mld_precset(prec,mld_sub_ovr_, prec_choice%novr, info) call mld_precset(prec,mld_sub_restr_, prec_choice%restr, info)
call mld_precset(prec,mld_sub_restr_, prec_choice%restr,info) call mld_precset(prec,mld_sub_prol_, prec_choice%prol, info)
call mld_precset(prec,mld_sub_prol_, prec_choice%prol, info) call mld_precset(prec,mld_sub_solve_, prec_choice%solve, info)
call mld_precset(prec,mld_sub_solve_, prec_choice%solve,info) call mld_precset(prec,mld_sub_fillin_, prec_choice%fill, info)
call mld_precset(prec,mld_sub_fillin_,prec_choice%fill,info) call mld_precset(prec,mld_sub_iluthrs_, prec_choice%thr, info)
call mld_precset(prec,mld_sub_iluthrs_, prec_choice%thr, info)
if (psb_toupper(prec_choice%prec) == 'ML') then
call mld_precset(prec,mld_aggr_kind_, prec_choice%aggrkind,info) call mld_precset(prec,mld_aggr_kind_, prec_choice%aggrkind,info)
call mld_precset(prec,mld_aggr_alg_, prec_choice%aggr_alg,info) call mld_precset(prec,mld_aggr_alg_, prec_choice%aggr_alg,info)
call mld_precset(prec,mld_ml_type_, prec_choice%mltype, info) call mld_precset(prec,mld_ml_type_, prec_choice%mltype, info)
@ -267,6 +267,16 @@ program zf_sample
call mld_precset(prec,mld_coarse_fillin_, prec_choice%cfill, info) call mld_precset(prec,mld_coarse_fillin_, prec_choice%cfill, info)
call mld_precset(prec,mld_coarse_iluthrs_, prec_choice%cthres, info) call mld_precset(prec,mld_coarse_iluthrs_, prec_choice%cthres, info)
call mld_precset(prec,mld_coarse_sweeps_, prec_choice%cjswp, info) call mld_precset(prec,mld_coarse_sweeps_, prec_choice%cjswp, info)
else
nlv = 1
call mld_precinit(prec,prec_choice%prec,info)
call mld_precset(prec,mld_smoother_sweeps_, prec_choice%jsweeps, info)
call mld_precset(prec,mld_sub_ovr_, prec_choice%novr, info)
call mld_precset(prec,mld_sub_restr_, prec_choice%restr, info)
call mld_precset(prec,mld_sub_prol_, prec_choice%prol, info)
call mld_precset(prec,mld_sub_solve_, prec_choice%solve, info)
call mld_precset(prec,mld_sub_fillin_, prec_choice%fill, info)
call mld_precset(prec,mld_sub_iluthrs_, prec_choice%thr, info)
end if end if
! building the preconditioner ! building the preconditioner
@ -278,12 +288,11 @@ program zf_sample
goto 9999 goto 9999
end if end if
call psb_amx(ictxt, tprec) call psb_amx(ictxt, tprec)
if(iam == psb_root_) then if(iam == psb_root_) then
write(*,'("Preconditioner time: ",es12.5)')tprec write(psb_out_unit,'("Preconditioner time: ",es12.5)')tprec
write(*,'(" ")') write(psb_out_unit,'(" ")')
end if end if
iparm = 0 iparm = 0
@ -308,36 +317,37 @@ program zf_sample
call psb_sum(ictxt,precsize) call psb_sum(ictxt,precsize)
if (iam == psb_root_) then if (iam == psb_root_) then
call mld_precdescr(prec,info) call mld_precdescr(prec,info)
write(*,'("Matrix: ",a)')mtrx_file write(psb_out_unit,'("Matrix: ",a)')mtrx_file
write(*,'("Computed solution on ",i8," processors")')np write(psb_out_unit,'("Computed solution on ",i8," processors")')np
write(*,'("Iterations to convergence : ",i6)')iter write(psb_out_unit,'("Iterations to convergence : ",i6)')iter
write(*,'("Error estimate on exit : ",es12.5)')err write(psb_out_unit,'("Error estimate on exit : ",es12.5)')err
write(*,'("Time to buil prec. : ",es12.5)')tprec write(psb_out_unit,'("Time to buil prec. : ",es12.5)')tprec
write(*,'("Time to solve matrix : ",es12.5)')t2 write(psb_out_unit,'("Time to solve matrix : ",es12.5)')t2
write(*,'("Time per iteration : ",es12.5)')t2/(iter) write(psb_out_unit,'("Time per iteration : ",es12.5)')t2/(iter)
write(*,'("Total time : ",es12.5)')t2+tprec write(psb_out_unit,'("Total time : ",es12.5)')t2+tprec
write(*,'("Residual norm 2 : ",es12.5)')resmx write(psb_out_unit,'("Residual norm 2 : ",es12.5)')resmx
write(*,'("Residual norm inf : ",es12.5)')resmxp write(psb_out_unit,'("Residual norm inf : ",es12.5)')resmxp
write(*,'("Total memory occupation for A : ",i12)')amatsize write(psb_out_unit,'("Total memory occupation for A : ",i12)')amatsize
write(*,'("Total memory occupation for DESC_A : ",i12)')descsize write(psb_out_unit,'("Total memory occupation for DESC_A : ",i12)')descsize
write(*,'("Total memory occupation for PREC : ",i12)')precsize write(psb_out_unit,'("Total memory occupation for PREC : ",i12)')precsize
end if end if
allocate(x_col_glob(m_problem),r_col_glob(m_problem),stat=ierr) allocate(x_col_glob(m_problem),r_col_glob(m_problem),stat=ierr)
if (ierr /= 0) then if (ierr /= 0) then
write(0,*) 'allocation error: no data collection' write(psb_err_unit,*) 'allocation error: no data collection'
else else
call psb_gather(x_col_glob,x_col,desc_a,info,root=psb_root_) call psb_gather(x_col_glob,x_col,desc_a,info,root=psb_root_)
call psb_gather(r_col_glob,r_col,desc_a,info,root=psb_root_) call psb_gather(r_col_glob,r_col,desc_a,info,root=psb_root_)
if (iam == psb_root_) then if (iam == psb_root_) then
write(0,'(" ")') write(psb_err_unit,'(" ")')
write(0,'("Saving x on file")') write(psb_err_unit,'("Saving x on file")')
write(20,*) 'matrix: ',mtrx_file write(20,*) 'matrix: ',mtrx_file
write(20,*) 'computed solution on ',np,' processors.' write(20,*) 'computed solution on ',np,' processors.'
write(20,*) 'iterations to convergence: ',iter write(20,*) 'iterations to convergence: ',iter
write(20,*) 'error estimate (infinity norm) on exit:', & write(20,*) 'error estimate (infinity norm) on exit:', &
& ' ||r||/(||a||||x||+||b||) = ',err & ' ||r||/(||a||||x||+||b||) = ',err
write(20,*) 'max residual = ',resmx, resmxp write(20,'("Residual norm 2 : ",es12.5)')resmx
write(20,'("Residual norm inf : ",es12.5)')resmxp
write(20,'(a8,4(2x,a20))') 'I','X(I)','R(I)','B(I)' write(20,'(a8,4(2x,a20))') 'I','X(I)','R(I)','B(I)'
do i=1,m_problem do i=1,m_problem
write(20,998) i,x_col_glob(i),r_col_glob(i),b_col_glob(i) write(20,998) i,x_col_glob(i),r_col_glob(i),b_col_glob(i)
@ -401,7 +411,9 @@ contains
call read_data(prec%solve,5) ! Factorization type: ILU, SuperLU, UMFPACK. call read_data(prec%solve,5) ! Factorization type: ILU, SuperLU, UMFPACK.
call read_data(prec%fill,5) ! Fill-in for factorization call read_data(prec%fill,5) ! Fill-in for factorization
call read_data(prec%thr,5) ! Threshold for fact. ILU(T) call read_data(prec%thr,5) ! Threshold for fact. ILU(T)
call read_data(prec%jsweeps,5) ! Jacobi sweeps for PJAC
if (psb_toupper(prec%prec) == 'ML') then if (psb_toupper(prec%prec) == 'ML') then
call read_data(prec%smther,5) ! Smoother type.
call read_data(prec%nlev,5) ! Number of levels in multilevel prec. call read_data(prec%nlev,5) ! Number of levels in multilevel prec.
call read_data(prec%aggrkind,5) ! smoothed/raw aggregatin call read_data(prec%aggrkind,5) ! smoothed/raw aggregatin
call read_data(prec%aggr_alg,5) ! local or global aggregation call read_data(prec%aggr_alg,5) ! local or global aggregation
@ -437,7 +449,9 @@ contains
call psb_bcast(icontxt,prec%solve) ! Factorization type: ILU, SuperLU, UMFPACK. call psb_bcast(icontxt,prec%solve) ! Factorization type: ILU, SuperLU, UMFPACK.
call psb_bcast(icontxt,prec%fill) ! Fill-in for factorization call psb_bcast(icontxt,prec%fill) ! Fill-in for factorization
call psb_bcast(icontxt,prec%thr) ! Threshold for fact. ILU(T) call psb_bcast(icontxt,prec%thr) ! Threshold for fact. ILU(T)
call psb_bcast(icontxt,prec%jsweeps) ! Jacobi sweeps
if (psb_toupper(prec%prec) == 'ML') then if (psb_toupper(prec%prec) == 'ML') then
call psb_bcast(icontxt,prec%smther) ! Smoother type.
call psb_bcast(icontxt,prec%nlev) ! Number of levels in multilevel prec. call psb_bcast(icontxt,prec%nlev) ! Number of levels in multilevel prec.
call psb_bcast(icontxt,prec%aggrkind) ! smoothed/raw aggregatin call psb_bcast(icontxt,prec%aggrkind) ! smoothed/raw aggregatin
call psb_bcast(icontxt,prec%aggr_alg) ! local or global aggregation call psb_bcast(icontxt,prec%aggr_alg) ! local or global aggregation

2
tests/pdegen/runs/ppde.inp
Unescape Escape View File

@ -21,7 +21,7 @@ SMOOTHED ! Kind of aggregation: SMOOTHED, NONSMOOTHED, MINENE
DEC ! Type of aggregation DEC SYMDEC GLB DEC ! Type of aggregation DEC SYMDEC GLB
MULT ! Type of multilevel correction: ADD MULT MULT ! Type of multilevel correction: ADD MULT
POST ! Side of multiplicative correction PRE POST TWOSIDE (ignored for ADD) POST ! Side of multiplicative correction PRE POST TWOSIDE (ignored for ADD)
DIST ! Coarse level: matrix distribution DIST REPL REPL ! Coarse level: matrix distribution DIST REPL
BJAC ! Coarse level: solver JACOBI BJAC UMF SLU SLUDIST BJAC ! Coarse level: solver JACOBI BJAC UMF SLU SLUDIST
ILU ! Coarse level: subsolver DSCALE ILU UMF SLU SLUDIST ILU ! Coarse level: subsolver DSCALE ILU UMF SLU SLUDIST
1 ! Coarse level: Level-set N for ILU(N) 1 ! Coarse level: Level-set N for ILU(N)

Write Preview
Loading…
Cancel
Save