diff --git a/mlprec/mld_base_prec_type.F90 b/mlprec/mld_base_prec_type.F90 index 47b699e7..10c014dc 100644 --- a/mlprec/mld_base_prec_type.F90 +++ b/mlprec/mld_base_prec_type.F90 @@ -132,7 +132,8 @@ module mld_base_prec_type integer, parameter :: mld_sub_prol_ = 4 integer, parameter :: mld_sub_ren_ = 5 integer, parameter :: mld_sub_ovr_ = 6 - integer, parameter :: mld_sub_fillin_ = 8 + integer, parameter :: mld_sub_fillin_ = 7 + integer, parameter :: mld_ilu_scale_ = 8 !! 2 ints for 64 bit versions integer, parameter :: mld_slu_ptr_ = 10 integer, parameter :: mld_umf_symptr_ = 12 @@ -173,14 +174,14 @@ module mld_base_prec_type ! ! Legal values for entry: mld_sub_solve_ ! - integer, parameter :: mld_slv_delta_ = mld_max_prec_+1 - integer, parameter :: mld_f_none_ = mld_slv_delta_+0 + integer, parameter :: mld_slv_delta_ = mld_max_prec_+1 + integer, parameter :: mld_f_none_ = mld_slv_delta_+0 integer, parameter :: mld_diag_scale_ = mld_slv_delta_+1 - integer, parameter :: mld_ilu_n_ = mld_slv_delta_+2 + integer, parameter :: mld_ilu_n_ = mld_slv_delta_+2 integer, parameter :: mld_milu_n_ = mld_slv_delta_+3 - integer, parameter :: mld_ilu_t_ = mld_slv_delta_+4 + integer, parameter :: mld_ilu_t_ = mld_slv_delta_+4 integer, parameter :: mld_slu_ = mld_slv_delta_+5 - integer, parameter :: mld_umf_ = mld_slv_delta_+6 + integer, parameter :: mld_umf_ = mld_slv_delta_+6 integer, parameter :: mld_sludist_ = mld_slv_delta_+7 integer, parameter :: mld_max_sub_solve_= mld_slv_delta_+7 integer, parameter :: mld_min_sub_solve_= mld_diag_scale_ @@ -193,6 +194,14 @@ module mld_base_prec_type ! For the time being we are disabling GPS renumbering. integer, parameter :: mld_max_renum_=1 ! + ! Legal values for entry: mld_ilu_scale_ + ! + integer, parameter :: mld_ilu_scale_none_ = 0 + integer, parameter :: mld_ilu_scale_maxval_ = 1 + integer, parameter :: mld_ilu_scale_diag_ = 2 + ! For the time being enable only maxval scale + integer, parameter :: mld_max_ilu_scale_ = 1 + ! ! Legal values for entry: mld_ml_type_ ! integer, parameter :: mld_no_ml_ = 0 @@ -615,6 +624,13 @@ contains is_legal_renum = ((ip >= 0).and.(ip <= mld_max_renum_)) return end function is_legal_renum + function is_legal_ilu_scale(ip) + implicit none + integer, intent(in) :: ip + logical :: is_legal_ilu_scale + is_legal_ilu_scale = ((ip >= mld_ilu_scale_none_).and.(ip <= mld_max_ilu_scale_)) + return + end function is_legal_ilu_scale function is_legal_jac_sweeps(ip) implicit none integer, intent(in) :: ip diff --git a/mlprec/mld_c_id_solver.f90 b/mlprec/mld_c_id_solver.f90 index a19c3846..9d85e34b 100644 --- a/mlprec/mld_c_id_solver.f90 +++ b/mlprec/mld_c_id_solver.f90 @@ -56,14 +56,13 @@ module mld_c_id_solver procedure, pass(sv) :: setc => c_id_solver_setc procedure, pass(sv) :: setr => c_id_solver_setr procedure, pass(sv) :: descr => c_id_solver_descr - procedure, pass(sv) :: sizeof => c_id_solver_sizeof end type mld_c_id_solver_type private :: c_id_solver_bld, c_id_solver_apply, & & c_id_solver_free, c_id_solver_seti, & & c_id_solver_setc, c_id_solver_setr,& - & c_id_solver_descr, c_id_solver_sizeof + & c_id_solver_descr contains @@ -265,17 +264,4 @@ contains end subroutine c_id_solver_descr - function c_id_solver_sizeof(sv) result(val) - use psb_base_mod - implicit none - ! Arguments - class(mld_c_id_solver_type), intent(in) :: sv - integer(psb_long_int_k_) :: val - integer :: i - - val = 0 - - return - end function c_id_solver_sizeof - end module mld_c_id_solver diff --git a/mlprec/mld_d_as_smoother.f90 b/mlprec/mld_d_as_smoother.f90 index caecf0ab..cc383db3 100644 --- a/mlprec/mld_d_as_smoother.f90 +++ b/mlprec/mld_d_as_smoother.f90 @@ -45,7 +45,7 @@ module mld_d_as_smoother use mld_d_prec_type - + type, extends(mld_d_base_smoother_type) :: mld_d_as_smoother_type ! The local solver component is inherited from the ! parent type. @@ -67,6 +67,7 @@ module mld_d_as_smoother procedure, pass(sm) :: descr => d_as_smoother_descr procedure, pass(sm) :: sizeof => d_as_smoother_sizeof procedure, pass(sm) :: default => d_as_smoother_default + procedure, pass(sm) :: get_nzeros => d_as_smoother_get_nzeros end type mld_d_as_smoother_type @@ -75,8 +76,9 @@ module mld_d_as_smoother & d_as_smoother_setc, d_as_smoother_setr,& & d_as_smoother_descr, d_as_smoother_sizeof, & & d_as_smoother_check, d_as_smoother_default,& - & d_as_smoother_dmp, d_as_smoother_apply_vect - + & d_as_smoother_dmp, d_as_smoother_apply_vect,& + & d_as_smoother_get_nzeros + character(len=6), parameter, private :: & & restrict_names(0:4)=(/'none ','halo ',' ',' ',' '/) character(len=12), parameter, private :: & @@ -94,12 +96,12 @@ contains ! Arguments class(mld_d_as_smoother_type), intent(inout) :: sm - + sm%restr = psb_halo_ sm%prol = psb_none_ sm%novr = 1 - + if (allocated(sm%sv)) then call sm%sv%default() end if @@ -129,7 +131,7 @@ contains call mld_check_def(sm%novr,& & 'Overlap layers ',0,is_legal_n_ovr) - + if (allocated(sm%sv)) then call sm%sv%check(info) else @@ -139,7 +141,7 @@ contains end if if (info /= psb_success_) goto 9999 - + call psb_erractionrestore(err_act) return @@ -251,7 +253,7 @@ contains vx = x%getCopy() - + call psb_geall(vtx,sm%desc_data,info) call psb_geasb(vtx,sm%desc_data,info,mold=x%v) call psb_geall(vty,sm%desc_data,info) @@ -601,7 +603,7 @@ contains call vww%free(info) call vtx%free(info) call vty%free(info) - + call psb_erractionrestore(err_act) return @@ -1123,7 +1125,7 @@ contains & write(debug_unit,*) me,' ',trim(name),& & ' From cdbldext _:',sm%desc_data%get_local_rows(),& & sm%desc_data%get_local_cols() - + if (info /= psb_success_) then info=psb_err_from_subroutine_ ch_err='psb_cdbldext' @@ -1141,14 +1143,14 @@ contains ! matrix data_ = psb_comm_ext_ Call psb_sphalo(a,sm%desc_data,blck,info,data=data_,rowscale=.true.) - + if (info /= psb_success_) then info=psb_err_from_subroutine_ ch_err='psb_sphalo' call psb_errpush(info,name,a_err=ch_err) goto 9999 end if - + if (debug_level >=psb_debug_outer_) & & write(debug_unit,*) me,' ',trim(name),& & 'After psb_sphalo ',& @@ -1162,13 +1164,13 @@ contains nrow_a = a%get_nrows() n_row = sm%desc_data%get_local_rows() n_col = sm%desc_data%get_local_cols() - + if (info == psb_success_) call a%csclip(sm%nd,info,& & jmin=nrow_a+1,rscale=.false.,cscale=.false.) if (info == psb_success_) call blck%csclip(atmp,info,& & jmin=nrow_a+1,rscale=.false.,cscale=.false.) if (info == psb_success_) call psb_rwextd(n_row,sm%nd,info,b=atmp) - + if (info == psb_success_) then if (present(amold)) then call sm%nd%cscnv(info,& @@ -1432,6 +1434,18 @@ contains return end function d_as_smoother_sizeof + function d_as_smoother_get_nzeros(sm) result(val) + implicit none + class(mld_d_as_smoother_type), intent(in) :: sm + integer(psb_long_int_k_) :: val + integer :: i + val = 0 + if (allocated(sm%sv)) & + & val = sm%sv%get_nzeros() + val = val + sm%nd%get_nzeros() + + end function d_as_smoother_get_nzeros + subroutine d_as_smoother_dmp(sm,ictxt,level,info,prefix,head,smoother,solver) use psb_base_mod implicit none diff --git a/mlprec/mld_d_diag_solver.f90 b/mlprec/mld_d_diag_solver.f90 index 5c17657a..fc0a4df6 100644 --- a/mlprec/mld_d_diag_solver.f90 +++ b/mlprec/mld_d_diag_solver.f90 @@ -60,13 +60,15 @@ module mld_d_diag_solver procedure, pass(sv) :: setr => d_diag_solver_setr procedure, pass(sv) :: descr => d_diag_solver_descr procedure, pass(sv) :: sizeof => d_diag_solver_sizeof + procedure, pass(sv) :: get_nzeros => d_diag_solver_get_nzeros end type mld_d_diag_solver_type private :: d_diag_solver_bld, d_diag_solver_apply, & & d_diag_solver_free, d_diag_solver_seti, & & d_diag_solver_setc, d_diag_solver_setr,& - & d_diag_solver_descr, d_diag_solver_sizeof + & d_diag_solver_descr, d_diag_solver_sizeof,& + & d_diag_solver_get_nzeros contains @@ -487,7 +489,9 @@ contains call psb_erractionsave(err_act) info = psb_success_ - + + call sv%dv%free(info) + if (allocated(sv%d)) then deallocate(sv%d,stat=info) if (info /= psb_success_) then @@ -549,9 +553,23 @@ contains integer :: i val = 0 - if (allocated(sv%d)) val = val + psb_sizeof_dp * size(sv%d) + if (allocated(sv%dv)) val = val + sv%dv%sizeof() return end function d_diag_solver_sizeof + function d_diag_solver_get_nzeros(sv) result(val) + use psb_base_mod + implicit none + ! Arguments + class(mld_d_diag_solver_type), intent(in) :: sv + integer(psb_long_int_k_) :: val + integer :: i + + val = 0 + if (allocated(sv%dv)) val = val + sv%dv%get_nrows() + + return + end function d_diag_solver_get_nzeros + end module mld_d_diag_solver diff --git a/mlprec/mld_d_id_solver.f90 b/mlprec/mld_d_id_solver.f90 index 68493ee1..7c8e50a6 100644 --- a/mlprec/mld_d_id_solver.f90 +++ b/mlprec/mld_d_id_solver.f90 @@ -57,14 +57,13 @@ module mld_d_id_solver procedure, pass(sv) :: setc => d_id_solver_setc procedure, pass(sv) :: setr => d_id_solver_setr procedure, pass(sv) :: descr => d_id_solver_descr - procedure, pass(sv) :: sizeof => d_id_solver_sizeof end type mld_d_id_solver_type private :: d_id_solver_bld, d_id_solver_apply, & & d_id_solver_free, d_id_solver_seti, & & d_id_solver_setc, d_id_solver_setr,& - & d_id_solver_descr, d_id_solver_sizeof + & d_id_solver_descr contains @@ -316,17 +315,4 @@ contains end subroutine d_id_solver_descr - function d_id_solver_sizeof(sv) result(val) - use psb_base_mod - implicit none - ! Arguments - class(mld_d_id_solver_type), intent(in) :: sv - integer(psb_long_int_k_) :: val - integer :: i - - val = 0 - - return - end function d_id_solver_sizeof - end module mld_d_id_solver diff --git a/mlprec/mld_d_ilu_fact_mod.f90 b/mlprec/mld_d_ilu_fact_mod.f90 index dbec54e3..624ea16e 100644 --- a/mlprec/mld_d_ilu_fact_mod.f90 +++ b/mlprec/mld_d_ilu_fact_mod.f90 @@ -28,15 +28,16 @@ module mld_d_ilu_fact_mod end interface interface mld_ilut_fact - subroutine mld_dilut_fact(fill_in,thres,a,l,u,d,info,blck) + subroutine mld_dilut_fact(fill_in,thres,a,l,u,d,info,blck,iscale) import psb_dspmat_type, psb_dpk_ - integer, intent(in) :: fill_in - real(psb_dpk_), intent(in) :: thres - integer, intent(out) :: info + integer, intent(in) :: fill_in + real(psb_dpk_), intent(in) :: thres + integer, intent(out) :: info type(psb_dspmat_type),intent(in) :: a type(psb_dspmat_type),intent(inout) :: l,u + real(psb_dpk_), intent(inout) :: d(:) type(psb_dspmat_type),intent(in), optional, target :: blck - real(psb_dpk_), intent(inout) :: d(:) + integer, intent(in), optional :: iscale end subroutine mld_dilut_fact end interface diff --git a/mlprec/mld_d_ilu_solver.f90 b/mlprec/mld_d_ilu_solver.f90 index e88f7c9b..6ea57a5d 100644 --- a/mlprec/mld_d_ilu_solver.f90 +++ b/mlprec/mld_d_ilu_solver.f90 @@ -51,7 +51,7 @@ module mld_d_ilu_solver type, extends(mld_d_base_solver_type) :: mld_d_ilu_solver_type type(psb_dspmat_type) :: l, u real(psb_dpk_), allocatable :: d(:) - type(psb_d_vect_type), allocatable :: dv + type(psb_d_vect_type) :: dv integer :: fact_type, fill_in real(psb_dpk_) :: thresh contains @@ -64,9 +64,9 @@ module mld_d_ilu_solver procedure, pass(sv) :: setc => d_ilu_solver_setc procedure, pass(sv) :: setr => d_ilu_solver_setr procedure, pass(sv) :: descr => d_ilu_solver_descr + procedure, pass(sv) :: default => d_ilu_solver_default procedure, pass(sv) :: sizeof => d_ilu_solver_sizeof procedure, pass(sv) :: get_nzeros => d_ilu_solver_get_nzeros - procedure, pass(sv) :: default => d_ilu_solver_default end type mld_d_ilu_solver_type @@ -193,11 +193,6 @@ contains call psb_errpush(info,name,i_err=(/3,n_row,0,0,0/)) goto 9999 end if - if (.not.allocated(sv%dv)) then - info = 1124 - call psb_errpush(info,name,a_err="preconditioner: DV") - goto 9999 - end if if (sv%dv%get_nrows() < n_row) then info = 1124 call psb_errpush(info,name,a_err="preconditioner: DV") @@ -449,19 +444,8 @@ contains deallocate(sv%d) endif endif - if (.not.allocated(sv%d)) then - allocate(sv%d(n_row),stat=info) - endif - if (info == psb_success_) then - allocate(sv%dv, stat=info) - if (info == 0) then - if (present(vmold)) then - allocate(sv%dv%v,mold=vmold,stat=info) - else - allocate(psb_d_base_vect_type :: sv%dv%v,stat=info) - end if - end if - end if + if (.not.allocated(sv%d)) allocate(sv%d(n_row),stat=info) + if (info /= psb_success_) then call psb_errpush(psb_err_from_subroutine_,name,a_err='Allocate') goto 9999 @@ -562,7 +546,7 @@ contains call sv%l%trim() call sv%u%set_asb() call sv%u%trim() - call sv%dv%bld(sv%d) + call sv%dv%bld(sv%d,mold=vmold) if (present(amold)) then call sv%l%cscnv(info,mold=amold) @@ -725,6 +709,7 @@ contains end if call sv%l%free() call sv%u%free() + call sv%dv%free(info) call psb_erractionrestore(err_act) return @@ -795,7 +780,7 @@ contains integer :: i val = 0 - if (allocated(sv%dv)) val = val + sv%dv%get_nrows() + val = val + sv%dv%get_nrows() val = val + sv%l%get_nzeros() val = val + sv%u%get_nzeros() @@ -811,9 +796,9 @@ contains integer :: i val = 2*psb_sizeof_int + psb_sizeof_dp - if (allocated(sv%d)) val = val + psb_sizeof_dp * size(sv%d) - val = val + psb_sizeof(sv%l) - val = val + psb_sizeof(sv%u) + val = val + sv%dv%sizeof() + val = val + sv%l%sizeof() + val = val + sv%u%sizeof() return end function d_ilu_solver_sizeof @@ -868,5 +853,4 @@ contains end subroutine d_ilu_solver_dmp - end module mld_d_ilu_solver diff --git a/mlprec/mld_d_jac_smoother.f90 b/mlprec/mld_d_jac_smoother.f90 index 8aed6d2b..8592ba2e 100644 --- a/mlprec/mld_d_jac_smoother.f90 +++ b/mlprec/mld_d_jac_smoother.f90 @@ -63,6 +63,7 @@ module mld_d_jac_smoother procedure, pass(sm) :: setr => d_jac_smoother_setr procedure, pass(sm) :: descr => d_jac_smoother_descr procedure, pass(sm) :: sizeof => d_jac_smoother_sizeof + procedure, pass(sm) :: get_nzeros => d_jac_smoother_get_nzeros end type mld_d_jac_smoother_type @@ -70,7 +71,7 @@ module mld_d_jac_smoother & d_jac_smoother_free, d_jac_smoother_seti, & & d_jac_smoother_setc, d_jac_smoother_setr,& & d_jac_smoother_descr, d_jac_smoother_sizeof, & - & d_jac_smoother_apply_vect + & d_jac_smoother_apply_vect, d_jac_smoother_get_nzeros @@ -687,4 +688,19 @@ contains return end function d_jac_smoother_sizeof + function d_jac_smoother_get_nzeros(sm) result(val) + use psb_base_mod + implicit none + ! Arguments + class(mld_d_jac_smoother_type), intent(in) :: sm + integer(psb_long_int_k_) :: val + integer :: i + + val = 0 + if (allocated(sm%sv)) val = val + sm%sv%get_nzeros() + val = val + sm%nd%get_nzeros() + + return + end function d_jac_smoother_get_nzeros + end module mld_d_jac_smoother diff --git a/mlprec/mld_dilut_fact.f90 b/mlprec/mld_dilut_fact.f90 index c1b6a7d3..6a9a3bf5 100644 --- a/mlprec/mld_dilut_fact.f90 +++ b/mlprec/mld_dilut_fact.f90 @@ -92,7 +92,7 @@ ! greater than 0. If the overlap is 0 or the matrix has been reordered ! (see mld_fact_bld), then blck does not contain any row. ! -subroutine mld_dilut_fact(fill_in,thres,a,l,u,d,info,blck) +subroutine mld_dilut_fact(fill_in,thres,a,l,u,d,info,blck,iscale) use psb_base_mod use mld_d_ilu_fact_mod, mld_protect_name => mld_dilut_fact @@ -105,13 +105,15 @@ subroutine mld_dilut_fact(fill_in,thres,a,l,u,d,info,blck) integer, intent(out) :: info type(psb_dspmat_type),intent(in) :: a type(psb_dspmat_type),intent(inout) :: l,u + real(psb_dpk_), intent(inout) :: d(:) type(psb_dspmat_type),intent(in), optional, target :: blck - real(psb_dpk_), intent(inout) :: d(:) + integer, intent(in), optional :: iscale ! Local Variables - integer :: l1, l2, m, err_act + integer :: l1, l2, m, err_act, iscale_ type(psb_dspmat_type), pointer :: blck_ type(psb_d_csr_sparse_mat) :: ll, uu + real(psb_dpk_) :: scale character(len=20) :: name, ch_err name='mld_dilut_fact' @@ -138,7 +140,24 @@ subroutine mld_dilut_fact(fill_in,thres,a,l,u,d,info,blck) goto 9999 end if endif - + if (present(iscale)) then + iscale_ = iscale + else + iscale_ = mld_ilu_scale_none_ + end if + + select case(iscale_) + case(mld_ilu_scale_none_) + scale = done + case(mld_ilu_scale_maxval_) + scale = max(a%maxval(),blck_%maxval()) + scale = done/scale + case default + info=psb_err_input_asize_invalid_i_ + call psb_errpush(info,name,i_err=(/9,iscale_,0,0,0/)) + goto 9999 + end select + m = a%get_nrows() + blck_%get_nrows() if ((m /= l%get_nrows()).or.(m /= u%get_nrows()).or.& & (m > size(d)) ) then @@ -155,7 +174,7 @@ subroutine mld_dilut_fact(fill_in,thres,a,l,u,d,info,blck) ! Compute the ILU(k,t) factorization ! call mld_dilut_factint(fill_in,thres,a,blck_,& - & d,ll%val,ll%ja,ll%irp,uu%val,uu%ja,uu%irp,l1,l2,info) + & d,ll%val,ll%ja,ll%irp,uu%val,uu%ja,uu%irp,l1,l2,info,scale) if (info /= psb_success_) then info=psb_err_from_subroutine_ ch_err='mld_dilut_factint' @@ -270,7 +289,7 @@ contains ! Error code. ! subroutine mld_dilut_factint(fill_in,thres,a,b,& - & d,lval,lja,lirp,uval,uja,uirp,l1,l2,info) + & d,lval,lja,lirp,uval,uja,uirp,l1,l2,info,scale) use psb_base_mod @@ -279,21 +298,22 @@ contains ! Arguments integer, intent(in) :: fill_in real(psb_dpk_), intent(in) :: thres - type(psb_dspmat_type),intent(in) :: a,b + type(psb_dspmat_type),intent(in) :: a,b integer,intent(inout) :: l1,l2,info integer, allocatable, intent(inout) :: lja(:),lirp(:),uja(:),uirp(:) real(psb_dpk_), allocatable, intent(inout) :: lval(:),uval(:) real(psb_dpk_), intent(inout) :: d(:) + real(psb_dpk_), intent(in), optional :: scale ! Local Variables integer :: i, ktrw,err_act,nidx,nlw,nup,jmaxup, ma, mb, m - real(psb_dpk_) :: nrmi - integer, allocatable :: idxs(:) - real(psb_dpk_), allocatable :: row(:) + real(psb_dpk_) :: nrmi, weight + integer, allocatable :: idxs(:) + real(psb_dpk_), allocatable :: row(:) type(psb_int_heap) :: heap - type(psb_d_coo_sparse_mat) :: trw - character(len=20), parameter :: name='mld_dilut_factint' - character(len=20) :: ch_err + type(psb_d_coo_sparse_mat) :: trw + character(len=20), parameter :: name='mld_dilut_factint' + character(len=20) :: ch_err if (psb_get_errstatus() /= 0) return info = psb_success_ @@ -333,7 +353,8 @@ contains end if row(:) = dzero - + weight = done + if (present(scale)) weight = abs(scale) ! ! Cycle over the matrix rows ! @@ -349,9 +370,11 @@ contains ! d(i) = dzero if (i<=ma) then - call ilut_copyin(i,ma,a,i,1,m,nlw,nup,jmaxup,nrmi,row,heap,ktrw,trw,info) + call ilut_copyin(i,ma,a,i,1,m,nlw,nup,jmaxup,nrmi,weight,& + & row,heap,ktrw,trw,info) else - call ilut_copyin(i-ma,mb,b,i,1,m,nlw,nup,jmaxup,nrmi,row,heap,ktrw,trw,info) + call ilut_copyin(i-ma,mb,b,i,1,m,nlw,nup,jmaxup,nrmi,weight,& + & row,heap,ktrw,trw,info) endif ! @@ -362,7 +385,8 @@ contains ! ! Copy the row into lval/d(i)/uval ! - if (info == psb_success_) call ilut_copyout(fill_in,thres,i,m,nlw,nup,jmaxup,nrmi,row,nidx,idxs,& + if (info == psb_success_) call ilut_copyout(fill_in,thres,i,m,& + & nlw,nup,jmaxup,nrmi,row,nidx,idxs,& & l1,l2,lja,lirp,lval,d,uja,uirp,uval,info) if (info /= psb_success_) then @@ -372,6 +396,12 @@ contains end if end do + ! + ! Adjust diagonal accounting for scale factor + ! + if (weight /= done) then + d(1:m) = d(1:m)*weight + end if ! ! And we're done, so deallocate the memory @@ -483,14 +513,16 @@ contains ! until we empty the buffer. Thus we will make a call to psb_sp_getblk ! every nrb calls to copyin. If A is in CSR format it is unused. ! - subroutine ilut_copyin(i,m,a,jd,jmin,jmax,nlw,nup,jmaxup,nrmi,row,heap,ktrw,trw,info) + subroutine ilut_copyin(i,m,a,jd,jmin,jmax,nlw,nup,jmaxup,& + & nrmi,weight,row,heap,ktrw,trw,info) use psb_base_mod implicit none - type(psb_dspmat_type), intent(in) :: a + type(psb_dspmat_type), intent(in) :: a type(psb_d_coo_sparse_mat), intent(inout) :: trw integer, intent(in) :: i, m,jmin,jmax,jd integer, intent(inout) :: ktrw,nlw,nup,jmaxup,info real(psb_dpk_), intent(inout) :: nrmi,row(:) + real(psb_dpk_), intent(in) :: weight type(psb_int_heap), intent(inout) :: heap integer :: k,j,irb,kin,nz @@ -532,7 +564,7 @@ contains do j = aa%irp(i), aa%irp(i+1) - 1 k = aa%ja(j) if ((jmin<=k).and.(k<=jmax)) then - row(k) = aa%val(j) + row(k) = aa%val(j)*weight call psb_insert_heap(k,heap,info) if (info /= psb_success_) exit end if @@ -552,7 +584,7 @@ contains end if nz = aa%irp(i+1) - aa%irp(i) - nrmi = dnrm2(nz,aa%val(aa%irp(i)),ione) + nrmi = weight*dnrm2(nz,aa%val(aa%irp(i)),ione) class default @@ -583,7 +615,7 @@ contains if (trw%ia(ktrw) > i) exit k = trw%ja(ktrw) if ((jmin<=k).and.(k<=jmax)) then - row(k) = trw%val(ktrw) + row(k) = trw%val(ktrw)*weight call psb_insert_heap(k,heap,info) if (info /= psb_success_) exit @@ -599,7 +631,7 @@ contains ktrw = ktrw + 1 enddo nz = ktrw - kin - nrmi = dnrm2(nz,trw%val(kin),ione) + nrmi = weight*dnrm2(nz,trw%val(kin),ione) end select call psb_erractionrestore(err_act) @@ -860,8 +892,9 @@ contains ! The array where the entries of the row corresponding to the ! U factor are copied. ! - subroutine ilut_copyout(fill_in,thres,i,m,nlw,nup,jmaxup,nrmi,row, & - & nidx,idxs,l1,l2,lja,lirp,lval,d,uja,uirp,uval,info) + subroutine ilut_copyout(fill_in,thres,i,m,nlw,nup,jmaxup,& + & nrmi,row, nidx,idxs,l1,l2,lja,lirp,lval,& + & d,uja,uirp,uval,info) use psb_base_mod @@ -877,9 +910,9 @@ contains real(psb_dpk_), intent(inout) :: row(:), d(:) ! Local variables - real(psb_dpk_),allocatable :: xw(:) + real(psb_dpk_),allocatable :: xw(:) integer, allocatable :: xwid(:), indx(:) - real(psb_dpk_) :: witem + real(psb_dpk_) :: witem integer :: widx integer :: k,isz,err_act,int_err(5),idxp, nz type(psb_double_idx_heap) :: heap @@ -1134,7 +1167,6 @@ contains uja(l2) = xwid(k) uval(l2) = d(i)*xw(indx(k)) end do - ! ! Set row to zero ! diff --git a/mlprec/mld_s_id_solver.f90 b/mlprec/mld_s_id_solver.f90 index 3fcb40bc..fc4b3b12 100644 --- a/mlprec/mld_s_id_solver.f90 +++ b/mlprec/mld_s_id_solver.f90 @@ -56,14 +56,13 @@ module mld_s_id_solver procedure, pass(sv) :: setc => s_id_solver_setc procedure, pass(sv) :: setr => s_id_solver_setr procedure, pass(sv) :: descr => s_id_solver_descr - procedure, pass(sv) :: sizeof => s_id_solver_sizeof end type mld_s_id_solver_type private :: s_id_solver_bld, s_id_solver_apply, & & s_id_solver_free, s_id_solver_seti, & & s_id_solver_setc, s_id_solver_setr,& - & s_id_solver_descr, s_id_solver_sizeof + & s_id_solver_descr contains @@ -265,17 +264,4 @@ contains end subroutine s_id_solver_descr - function s_id_solver_sizeof(sv) result(val) - use psb_base_mod - implicit none - ! Arguments - class(mld_s_id_solver_type), intent(in) :: sv - integer(psb_long_int_k_) :: val - integer :: i - - val = 0 - - return - end function s_id_solver_sizeof - end module mld_s_id_solver diff --git a/mlprec/mld_z_id_solver.f90 b/mlprec/mld_z_id_solver.f90 index 3b8a0d8f..924dcf2c 100644 --- a/mlprec/mld_z_id_solver.f90 +++ b/mlprec/mld_z_id_solver.f90 @@ -56,14 +56,13 @@ module mld_z_id_solver procedure, pass(sv) :: setc => z_id_solver_setc procedure, pass(sv) :: setr => z_id_solver_setr procedure, pass(sv) :: descr => z_id_solver_descr - procedure, pass(sv) :: sizeof => z_id_solver_sizeof end type mld_z_id_solver_type private :: z_id_solver_bld, z_id_solver_apply, & & z_id_solver_free, z_id_solver_seti, & & z_id_solver_setc, z_id_solver_setr,& - & z_id_solver_descr, z_id_solver_sizeof + & z_id_solver_descr contains @@ -265,17 +264,4 @@ contains end subroutine z_id_solver_descr - function z_id_solver_sizeof(sv) result(val) - use psb_base_mod - implicit none - ! Arguments - class(mld_z_id_solver_type), intent(in) :: sv - integer(psb_long_int_k_) :: val - integer :: i - - val = 0 - - return - end function z_id_solver_sizeof - end module mld_z_id_solver