Merged changes from development.

Everything compiles.
psb_fnd_owner does not work right now.
ILmat
Salvatore Filippone 8 years ago
parent f409859a7d
commit 52c95853d6

@ -3,7 +3,7 @@ include ../../Make.inc
FOBJS = psi_compute_size.o psi_crea_bnd_elem.o psi_crea_index.o \
psi_crea_ovr_elem.o psi_bld_tmpovrl.o psi_dl_check.o \
psi_bld_tmphalo.o psi_sort_dl.o \
psi_desc_impl.o psi_exist_ovr_elem.o psi_list_search.o psi_srtlist.o
psi_desc_impl.o psi_list_search.o psi_srtlist.o
MPFOBJS = psi_desc_index.o psi_extrct_dl.o \
psi_fnd_owner.o psb_indx_map_fnd_owner.o

@ -66,15 +66,15 @@ subroutine psb_indx_map_fnd_owner(idx,iprc,idxmap,info)
integer(psb_ipk_), intent(out) :: info
integer(psb_lpk_), allocatable :: hproc(:)
integer(psb_lpk_), allocatable :: answers(:,:), idxsrch(:,:), hproc(:)
integer(psb_ipk_), allocatable :: helem(:),&
& answers(:,:),idxsrch(:,:), hhidx(:)
& hhidx(:)
integer(psb_mpk_), allocatable :: hsz(:),hidx(:), &
& sdsz(:),sdidx(:), rvsz(:), rvidx(:)
integer(psb_mpk_) :: icomm, minfo, iictxt
integer(psb_ipk_) :: i,n_row,n_col,err_act,ih,hsize,ip,isz,k,j,&
integer(psb_ipk_) :: i,n_row,n_col,err_act,hsize,ip,isz,j, k,&
& last_ih, last_j, nv
integer(psb_lpk_) :: mglob
integer(psb_lpk_) :: mglob, ih
integer(psb_ipk_) :: ictxt,np,me, nresp
logical, parameter :: gettime=.false.
real(psb_dpk_) :: t0, t1, t2, t3, t4, tamx, tidx
@ -230,8 +230,8 @@ subroutine psb_indx_map_fnd_owner(idx,iprc,idxmap,info)
rvidx(ip) = j
j = j + rvsz(ip)
end do
call mpi_alltoallv(hproc,sdsz,sdidx,psb_mpi_ipk_,&
& answers(:,1),rvsz,rvidx,psb_mpi_ipk_,&
call mpi_alltoallv(hproc,sdsz,sdidx,psb_mpi_lpk_,&
& answers(:,1),rvsz,rvidx,psb_mpi_lpk_,&
& icomm,minfo)
if (gettime) then
tamx = psb_wtime() - t3 + tamx

@ -65,7 +65,8 @@ subroutine psi_i_bld_tmpovrl(iv,desc,info)
!locals
integer(psb_ipk_) :: counter,i,j,np,me,loc_row,err,loc_col,nprocs,&
& l_ov_ix,l_ov_el,idx, err_act, itmpov, k, glx, icomm
& l_ov_ix,l_ov_el, err_act, itmpov, k, glx, icomm
integer(psb_ipk_) :: idx
integer(psb_ipk_), allocatable :: ov_idx(:),ov_el(:,:)
integer(psb_ipk_) :: ictxt,n_row, debug_unit, debug_level

@ -52,8 +52,7 @@
! nrcv - integer Total receive buffer size on the calling process
!
!
subroutine psi_i_crea_index(desc_a,index_in,index_out,glob_idx,nxch,nsnd,nrcv,info)
subroutine psi_i_crea_index(desc_a,index_in,index_out,nxch,nsnd,nrcv,info)
use psb_realloc_mod
use psb_desc_mod
use psb_error_mod
@ -63,9 +62,8 @@ subroutine psi_i_crea_index(desc_a,index_in,index_out,glob_idx,nxch,nsnd,nrcv,in
type(psb_desc_type), intent(in) :: desc_a
integer(psb_ipk_), intent(out) :: info,nxch,nsnd,nrcv
integer(psb_lpk_), intent(in) :: index_in(:)
integer(psb_ipk_), intent(in) :: index_in(:)
integer(psb_ipk_), allocatable, intent(inout) :: index_out(:)
logical :: glob_idx
! ....local scalars...
integer(psb_ipk_) :: ictxt, me, np, mode, err_act, dl_lda
@ -135,7 +133,7 @@ subroutine psi_i_crea_index(desc_a,index_in,index_out,glob_idx,nxch,nsnd,nrcv,in
& write(debug_unit,*) me,' ',trim(name),': calling psi_desc_index'
! Do the actual format conversion.
call psi_desc_index(desc_a,index_in,dep_list(1:,me),&
& length_dl(me),nsnd,nrcv, index_out,glob_idx,info)
& length_dl(me),nsnd,nrcv, index_out,info)
if(debug_level >= psb_debug_inner_) &
& write(debug_unit,*) me,' ',trim(name),': out of psi_desc_index',&
& size(index_out)

@ -61,10 +61,6 @@ subroutine psi_i_crea_ovr_elem(me,desc_overlap,ovr_elem,info)
integer(psb_ipk_) :: dim_ovr_elem
integer(psb_ipk_) :: pairtree(2)
! ...external function...
integer(psb_ipk_) :: psi_exist_ovr_elem
external :: psi_exist_ovr_elem
integer(psb_ipk_) :: nel, ip, ix, iel, insize, err_act, iproc
integer(psb_ipk_), allocatable :: telem(:,:)

@ -67,9 +67,9 @@ subroutine psi_i_cnv_dsc(halo_in,ovrlap_in,ext_in,cdesc, info, mold)
implicit none
! ....scalars parameters....
integer(psb_ipk_), intent(in) :: halo_in(:), ovrlap_in(:),ext_in(:)
integer(psb_ipk_), intent(in) :: halo_in(:), ovrlap_in(:),ext_in(:)
type(psb_desc_type), intent(inout) :: cdesc
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), intent(out) :: info
class(psb_i_base_vect_type), optional, intent(in) :: mold
! ....local scalars....
@ -102,7 +102,7 @@ subroutine psi_i_cnv_dsc(halo_in,ovrlap_in,ext_in,cdesc, info, mold)
! first the halo index
if (debug_level>0) write(debug_unit,*) me,'Calling crea_index on halo',&
& size(halo_in)
call psi_crea_index(cdesc,halo_in, idx_out,.false.,nxch,nsnd,nrcv,info)
call psi_crea_index(cdesc,halo_in, idx_out,nxch,nsnd,nrcv,info)
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='psi_crea_index')
goto 9999
@ -115,7 +115,7 @@ subroutine psi_i_cnv_dsc(halo_in,ovrlap_in,ext_in,cdesc, info, mold)
! then ext index
if (debug_level>0) write(debug_unit,*) me,'Calling crea_index on ext'
call psi_crea_index(cdesc,ext_in, idx_out,.false.,nxch,nsnd,nrcv,info)
call psi_crea_index(cdesc,ext_in, idx_out,nxch,nsnd,nrcv,info)
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='psi_crea_index')
goto 9999
@ -126,7 +126,7 @@ subroutine psi_i_cnv_dsc(halo_in,ovrlap_in,ext_in,cdesc, info, mold)
if (debug_level>0) write(debug_unit,*) me,'Calling crea_index on ovrlap'
! then the overlap index
call psi_crea_index(cdesc,ovrlap_in, idx_out,.true.,nxch,nsnd,nrcv,info)
call psi_crea_index(cdesc,ovrlap_in, idx_out,nxch,nsnd,nrcv,info)
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='psi_crea_index')
goto 9999
@ -150,7 +150,7 @@ subroutine psi_i_cnv_dsc(halo_in,ovrlap_in,ext_in,cdesc, info, mold)
if (debug_level>0) write(debug_unit,*) me,'Calling bld_ovr_mst'
call psi_bld_ovr_mst(me,cdesc%ovrlap_elem,tmp_mst_idx,info)
if (info == psb_success_) call psi_crea_index(cdesc,&
& tmp_mst_idx,idx_out,.false.,nxch,nsnd,nrcv,info)
& tmp_mst_idx,idx_out,nxch,nsnd,nrcv,info)
if (debug_level>0) write(debug_unit,*) me,'Done crea_indx'
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='psi_bld_ovr_mst')

@ -42,12 +42,9 @@
! mapping parts are used.
! index_in(:) - integer The index list, build format
! index_out(:) - integer(psb_ipk_), allocatable The index list, assembled format
! glob_idx - logical Whether the input indices are in local or global
! numbering; the global numbering is used when
! converting the overlap exchange lists.
! nxch - integer The number of data exchanges on the calling process
! nsnd - integer Total send buffer size on the calling process
! nrcv - integer Total receive buffer size on the calling process
! nxch - integer The number of data exchanges on the calling process
! nsnd - integer Total send buffer size on the calling process
! nrcv - integer Total receive buffer size on the calling process
!
! The format of the index lists. Copied from base/modules/psb_desc_type
!
@ -99,7 +96,7 @@
!
!
subroutine psi_i_desc_index(desc,index_in,dep_list,&
& length_dl,nsnd,nrcv,desc_index,isglob_in,info)
& length_dl,nsnd,nrcv,desc_index,info)
use psb_desc_mod
use psb_realloc_mod
use psb_error_mod
@ -116,11 +113,10 @@ subroutine psi_i_desc_index(desc,index_in,dep_list,&
! ...array parameters.....
type(psb_desc_type) :: desc
integer(psb_lpk_) :: index_in(:)
integer(psb_ipk_) :: index_in(:)
integer(psb_ipk_) :: dep_list(:)
integer(psb_ipk_),allocatable :: desc_index(:)
integer(psb_ipk_) :: length_dl,nsnd,nrcv,info
logical :: isglob_in
! ....local scalars...
integer(psb_ipk_) :: j,me,np,i,proc
! ...parameters...
@ -256,24 +252,15 @@ subroutine psi_i_desc_index(desc,index_in,dep_list,&
!
! note that here bsdinx is zero-based, hence the following loop
!
if (isglob_in) then
do j=1, nerv
sndbuf(bsdindx(proc+1)+j) = (index_in(i+j))
end do
else
sndbuf(bsdindx(proc+1)+1:bsdindx(proc+1)+nerv) = index_in(i+1:i+nerv)
call desc%indxmap%l2gip(sndbuf(bsdindx(proc+1)+1:bsdindx(proc+1)+nerv),&
& info)
!!$ call desc%indxmap%l2g(index_in(i+1:i+nerv),&
!!$ & sndbuf(bsdindx(proc+1)+1:bsdindx(proc+1)+nerv),&
!!$ & info)
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='l2g')
goto 9999
end if
endif
call desc%indxmap%l2g(index_in(i+1:i+nerv),&
& sndbuf(bsdindx(proc+1)+1:bsdindx(proc+1)+nerv),&
& info)
if (info /= psb_success_) then
call psb_errpush(psb_err_from_subroutine_,name,a_err='l2g')
goto 9999
end if
bsdindx(proc+1) = bsdindx(proc+1) + nerv
i = i + nerv + 1
end do

@ -1,73 +0,0 @@
!
! Parallel Sparse BLAS version 3.5
! (C) Copyright 2006, 2010, 2015, 2017
! Salvatore Filippone
! Alfredo Buttari CNRS-IRIT, Toulouse
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
!
integer function psi_exist_ovr_elem(ovr_elem, dim_list,elem_searched)
use psb_const_mod
! PURPOSE:
! == = ====
!
! If ELEM_SEARCHED exist in the list OVR_ELEM returns its position in
! the list, else returns -1
!
!
! INPUT
! == = ===
! OVRLAP_ELEMENT_D.: Contains for all overlap points belonging to
! the current process:
! 1. overlap point index
! 2. Number of domains sharing that overlap point
! the end is marked by a -1...............................
!
! DIM_LIST..........: Dimension of list OVRLAP_ELEMENT_D
!
! ELEM_SEARCHED.....:point's Local index identifier to be searched.
implicit none
! ....Scalars parameters....
integer(psb_ipk_) :: dim_list,elem_searched
! ...array parameters....
integer(psb_ipk_) :: ovr_elem(dim_list,*)
! ...local scalars....
integer(psb_ipk_) :: i
i=1
do while ((i.le.dim_list).and.(ovr_elem(i,1).ne.elem_searched))
i=i+1
enddo
if ((i.le.dim_list).and.(ovr_elem(i,1).eq.elem_searched)) then
psi_exist_ovr_elem=i
else
psi_exist_ovr_elem=-1
endif
end function psi_exist_ovr_elem

@ -136,7 +136,7 @@ subroutine psi_i_extract_dep_list(ictxt,is_bld,is_upd,desc_str,dep_list,&
integer(psb_ipk_) :: np,dl_lda,mode, info
! ....array parameters....
integer(psb_lpk_) :: desc_str(*)
integer(psb_ipk_) :: desc_str(*)
integer(psb_ipk_) :: dep_list(dl_lda,0:np),length_dl(0:np)
integer(psb_ipk_), allocatable :: itmp(:)
! .....local arrays....

@ -142,9 +142,9 @@ module psb_desc_mod
! psb_ovrl subroutine.
!
! 8. When the descriptor is in the BLD state the INDEX vectors contains only
! the indices to be received, organized as a sequence
! of entries of the form (proc,N,(lx1,lx2,...,lxn)) with owning process,
! number of indices (most often but not necessarily N=1), list of local indices.
! the indices to be received, organized as a sequence of entries of
! the form (proc,N,(lx1,lx2,...,lxn)) with owning process, number of indices
! (most often but not necessarily N=1), list of local indices.
! This is because we only know the list of halo indices to be received
! as we go about building the sparse matrix pattern, and we want the build
! phase to be loosely synchronized. Thus we record the indices we have to ask
@ -211,7 +211,6 @@ module psb_desc_mod
type(psb_i_vect_type) :: v_ovrlap_index
type(psb_i_vect_type) :: v_ovr_mst_idx
integer(psb_lpk_), allocatable :: tmp_ovrlap_index(:)
integer(psb_ipk_), allocatable :: ovrlap_elem(:,:)
integer(psb_ipk_), allocatable :: bnd_elem(:)
integer(psb_ipk_), allocatable :: lprm(:)
@ -283,6 +282,7 @@ module psb_desc_mod
module procedure psb_cdfree
end interface psb_free
private :: nullify_desc, cd_get_fmt,&
& cd_l2gs1, cd_l2gs2, cd_l2gv1, cd_l2gv2, cd_g2ls1,&
& cd_g2ls2, cd_g2lv1, cd_g2lv2, cd_g2ls1_ins,&
@ -1073,8 +1073,7 @@ contains
end subroutine psb_cd_clone
Subroutine psb_cd_get_recv_idx_loc(tmp,desc,data,info)
Subroutine psb_cd_get_recv_idx(tmp,desc,data,info)
use psb_error_mod
use psb_penv_mod
use psb_realloc_mod
@ -1157,7 +1156,7 @@ contains
return
end Subroutine psb_cd_get_recv_idx_loc
end Subroutine psb_cd_get_recv_idx
Subroutine psb_cd_get_recv_idx_glob(tmp,desc,data,info)

@ -93,7 +93,7 @@ module psb_const_mod
integer, parameter :: psb_ipk_ = psb_mpk_
integer, parameter :: psb_lpk_ = psb_epk_
#elif defined(INT_I8_L8)
integer, parameter :: psb_ipk_ = psb_mpk_
integer, parameter :: psb_ipk_ = psb_epk_
integer, parameter :: psb_lpk_ = psb_epk_
#else
! Unsupported combination, compilation will stop later on
@ -101,13 +101,13 @@ module psb_const_mod
integer, parameter :: psb_lpk_ = -1
#endif
integer(psb_ipk_), save :: psb_sizeof_sp
integer(psb_ipk_), save :: psb_sizeof_dp
integer(psb_ipk_), save :: psb_sizeof_i2p
integer(psb_ipk_), save :: psb_sizeof_mp
integer(psb_ipk_), save :: psb_sizeof_ep
integer(psb_ipk_), save :: psb_sizeof_ip
integer(psb_ipk_), save :: psb_sizeof_lp
integer(psb_ipk_), save :: psb_sizeof_sp
integer(psb_ipk_), save :: psb_sizeof_dp
integer(psb_ipk_), save :: psb_sizeof_i2p
integer(psb_ipk_), save :: psb_sizeof_mp
integer(psb_ipk_), save :: psb_sizeof_ep
integer(psb_ipk_), save :: psb_sizeof_ip
integer(psb_ipk_), save :: psb_sizeof_lp
!
! Integer type identifiers for MPI operations.
!

@ -30,12 +30,12 @@
!
!
module psi_c_mod
use psb_desc_mod, only : psb_desc_type, psb_ipk_, psb_spk_, psb_i_base_vect_type
use psi_c_comm_a_mod
use psb_c_base_vect_mod, only : psb_c_base_vect_type
use psb_c_base_multivect_mod, only : psb_c_base_multivect_type
use psi_c_comm_v_mod
end module psi_c_mod

@ -30,12 +30,12 @@
!
!
module psi_d_mod
use psb_desc_mod, only : psb_desc_type, psb_ipk_, psb_dpk_, psb_i_base_vect_type
use psi_d_comm_a_mod
use psb_d_base_vect_mod, only : psb_d_base_vect_type
use psb_d_base_multivect_mod, only : psb_d_base_multivect_type
use psi_d_comm_v_mod
end module psi_d_mod

@ -30,13 +30,13 @@
!
!
module psi_i_mod
use psb_desc_mod, only : psb_desc_type, psb_ipk_, psb_mpk_, psb_epk_, psb_lpk_
use psi_m_comm_a_mod
use psi_e_comm_a_mod
use psb_i_base_vect_mod, only : psb_i_base_vect_type
use psb_i_base_multivect_mod, only : psb_i_base_multivect_type
use psi_i_comm_v_mod
interface psi_compute_size
subroutine psi_i_compute_size(desc_data,&
@ -58,13 +58,12 @@ module psi_i_mod
end interface
interface psi_crea_index
subroutine psi_i_crea_index(desc_a,index_in,index_out,glob_idx,nxch,nsnd,nrcv,info)
subroutine psi_i_crea_index(desc_a,index_in,index_out,nxch,nsnd,nrcv,info)
import
type(psb_desc_type), intent(in) :: desc_a
integer(psb_ipk_), intent(out) :: nxch,nsnd,nrcv
integer(psb_ipk_), intent(in) :: index_in(:)
integer(psb_ipk_), allocatable, intent(inout) :: index_out(:)
logical :: glob_idx
integer(psb_ipk_), intent(out) :: info
end subroutine psi_i_crea_index
end interface
@ -80,14 +79,12 @@ module psi_i_mod
interface psi_desc_index
subroutine psi_i_desc_index(desc,index_in,dep_list,&
& length_dl,nsnd,nrcv,desc_index,isglob_in,info)
& length_dl,nsnd,nrcv,desc_index,info)
import
type(psb_desc_type) :: desc
integer(psb_lpk_) :: index_in(:)
integer(psb_ipk_) :: dep_list(:)
integer(psb_ipk_) :: index_in(:),dep_list(:)
integer(psb_ipk_),allocatable :: desc_index(:)
integer(psb_ipk_) :: length_dl,nsnd,nrcv
logical :: isglob_in
integer(psb_ipk_) :: info
end subroutine psi_i_desc_index
end interface
@ -116,8 +113,7 @@ module psi_i_mod
logical :: is_bld, is_upd
integer(psb_ipk_) :: ictxt
integer(psb_ipk_) :: dl_lda,mode
integer(psb_lpk_) :: desc_str(*)
integer(psb_ipk_) :: dep_list(dl_lda,0:np),length_dl(0:np)
integer(psb_ipk_) :: desc_str(*),dep_list(dl_lda,0:np),length_dl(0:np)
integer(psb_mpk_) :: np
integer(psb_ipk_) :: info
end subroutine psi_i_extract_dep_list
@ -146,7 +142,7 @@ module psi_i_mod
interface psi_bld_tmpovrl
subroutine psi_i_bld_tmpovrl(iv,desc,info)
import
integer(psb_lpk_), intent(in) :: iv(:)
integer(psb_ipk_), intent(in) :: iv(:)
type(psb_desc_type), intent(inout) :: desc
integer(psb_ipk_), intent(out) :: info
end subroutine psi_i_bld_tmpovrl
@ -206,5 +202,6 @@ module psi_i_mod
integer(psb_ipk_), intent(out) :: info
end subroutine psi_i_bld_ovr_mst
end interface
end module psi_i_mod

@ -30,13 +30,13 @@
!
!
module psi_l_mod
use psb_desc_mod, only : psb_desc_type, psb_ipk_, psb_mpk_, psb_epk_, psb_lpk_
use psi_m_comm_a_mod
use psi_e_comm_a_mod
use psb_l_base_vect_mod, only : psb_l_base_vect_type
use psb_l_base_multivect_mod, only : psb_l_base_multivect_type
use psi_l_comm_v_mod
end module psi_l_mod

@ -30,12 +30,12 @@
!
!
module psi_s_mod
use psb_desc_mod, only : psb_desc_type, psb_ipk_, psb_spk_, psb_i_base_vect_type
use psi_s_comm_a_mod
use psb_s_base_vect_mod, only : psb_s_base_vect_type
use psb_s_base_multivect_mod, only : psb_s_base_multivect_type
use psi_s_comm_v_mod
end module psi_s_mod

@ -30,12 +30,12 @@
!
!
module psi_z_mod
use psb_desc_mod, only : psb_desc_type, psb_ipk_, psb_dpk_, psb_i_base_vect_type
use psi_z_comm_a_mod
use psb_z_base_vect_mod, only : psb_z_base_vect_type
use psb_z_base_multivect_mod, only : psb_z_base_multivect_type
use psi_z_comm_v_mod
end module psi_z_mod

@ -258,25 +258,19 @@ Subroutine psb_ccdbldext(a,desc_a,novr,desc_ov,info, extype)
Do j=0,n_elem_recv-1
idx = ovrlap(counter+psb_elem_recv_+j)
call desc_ov%indxmap%l2g(idx,gidx,info)
If (gidx < 0) then
info=-3
call psb_errpush(info,name)
goto 9999
endif
call psb_ensure_size((cntov_o+3),orig_ovr,info,pad=-ione)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='psb_ensure_size')
goto 9999
end if
orig_ovr(cntov_o)=proc
orig_ovr(cntov_o+1)=1
orig_ovr(cntov_o+2)=gidx
orig_ovr(cntov_o+3)=-1
orig_ovr(cntov_o) = proc
orig_ovr(cntov_o+1) = 1
orig_ovr(cntov_o+2) = idx
orig_ovr(cntov_o+3) = -1
cntov_o=cntov_o+3
end Do
counter=counter+n_elem_recv+n_elem_send+3
counter = counter+n_elem_recv+n_elem_send+3
end Do
@ -327,16 +321,16 @@ Subroutine psb_ccdbldext(a,desc_a,novr,desc_ov,info, extype)
n_col_prev = desc_ov%get_local_cols()
Do While (halo(counter) /= -1)
tot_elem=0
proc=halo(counter+psb_proc_id_)
n_elem_recv=halo(counter+psb_n_elem_recv_)
n_elem_send=halo(counter+n_elem_recv+psb_n_elem_send_)
tot_elem = 0
proc = halo(counter+psb_proc_id_)
n_elem_recv = halo(counter+psb_n_elem_recv_)
n_elem_send = halo(counter+n_elem_recv+psb_n_elem_send_)
If ((counter+n_elem_recv+n_elem_send) > Size(halo)) then
info = -1
call psb_errpush(info,name)
goto 9999
end If
tot_recv=tot_recv+n_elem_recv
tot_recv = tot_recv+n_elem_recv
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),&
& ': tot_recv:',proc,n_elem_recv,tot_recv
@ -359,12 +353,6 @@ Subroutine psb_ccdbldext(a,desc_a,novr,desc_ov,info, extype)
end If
idx = halo(counter+psb_elem_recv_+j)
call desc_ov%l2g(idx,gidx,info)
If (gidx < 0) then
info=-3
call psb_errpush(info,name)
goto 9999
endif
call psb_ensure_size((counter_o+3),tmp_ovr_idx,info,pad=-ione)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
@ -374,7 +362,7 @@ Subroutine psb_ccdbldext(a,desc_a,novr,desc_ov,info, extype)
tmp_ovr_idx(counter_o) = proc
tmp_ovr_idx(counter_o+1) = 1
tmp_ovr_idx(counter_o+2) = gidx
tmp_ovr_idx(counter_o+2) = idx
tmp_ovr_idx(counter_o+3) = -1
counter_o=counter_o+3
call psb_ensure_size((counter_h+3),tmp_halo,info,pad=-ione)
@ -403,12 +391,6 @@ Subroutine psb_ccdbldext(a,desc_a,novr,desc_ov,info, extype)
Do j=0,n_elem_send-1
idx = halo(counter+psb_elem_send_+j)
call desc_ov%l2g(idx,gidx,info)
If (gidx < 0) then
info=-3
call psb_errpush(info,name)
goto 9999
endif
call psb_ensure_size((counter_o+3),tmp_ovr_idx,info,pad=-ione)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
@ -418,7 +400,7 @@ Subroutine psb_ccdbldext(a,desc_a,novr,desc_ov,info, extype)
tmp_ovr_idx(counter_o) = proc
tmp_ovr_idx(counter_o+1) = 1
tmp_ovr_idx(counter_o+2) = gidx
tmp_ovr_idx(counter_o+2) = idx
tmp_ovr_idx(counter_o+3) = -1
counter_o=counter_o+3
@ -427,6 +409,7 @@ Subroutine psb_ccdbldext(a,desc_a,novr,desc_ov,info, extype)
!
If (i_ovr <= (novr)) Then
call a%csget(idx,idx,n_elem,irow,icol,info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='csget')
@ -507,8 +490,8 @@ Subroutine psb_ccdbldext(a,desc_a,novr,desc_ov,info, extype)
lworkr = max(iszr,1)
end if
call mpi_alltoallv(works,sdsz,bsdindx,psb_mpi_ipk_,&
& workr,rvsz,brvindx,psb_mpi_ipk_,icomm,minfo)
call mpi_alltoallv(works,sdsz,bsdindx,psb_mpi_lpk_,&
& workr,rvsz,brvindx,psb_mpi_lpk_,icomm,minfo)
if (minfo /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='mpi_alltoallv')
@ -534,7 +517,7 @@ Subroutine psb_ccdbldext(a,desc_a,novr,desc_ov,info, extype)
j = 0
do i=1,iszr
if (maskr(i) < 0) then
j=j+1
j = j+1
works(j) = workr(i)
end if
end do
@ -557,8 +540,8 @@ Subroutine psb_ccdbldext(a,desc_a,novr,desc_ov,info, extype)
& ': Done fnd_owner', desc_ov%indxmap%get_state()
do i=1,iszs
gidx = works(i)
n_col = desc_ov%get_local_cols()
gidx = works(i)
n_col = desc_ov%get_local_cols()
call desc_ov%indxmap%g2l_ins(gidx,lidx,info)
if (desc_ov%get_local_cols() > n_col ) then
!
@ -603,7 +586,7 @@ Subroutine psb_ccdbldext(a,desc_a,novr,desc_ov,info, extype)
write(debug_unit,*) me,' ',trim(name),':Calling Crea_index'
end if
call psi_crea_index(desc_ov,t_halo_in,t_halo_out,.false.,&
call psi_crea_index(desc_ov,t_halo_in,t_halo_out,&
& nxch,nsnd,nrcv,info)
if (debug_level >= psb_debug_outer_) then
@ -661,7 +644,7 @@ Subroutine psb_ccdbldext(a,desc_a,novr,desc_ov,info, extype)
goto 9999
end if
orig_ovr(cntov_o:cntov_o+counter_o-1) = tmp_ovr_idx(1:counter_o)
cntov_o = cntov_o+counter_o-1
cntov_o = cntov_o+counter_o-1
orig_ovr(cntov_o:) = -1
call psb_move_alloc(orig_ovr,desc_ov%ovrlap_index,info)
deallocate(tmp_ovr_idx,stat=info)

@ -64,8 +64,8 @@ subroutine psb_cd_inloc(v, ictxt, desc, info, globalcheck,idx)
integer(psb_lpk_) :: m, n, nrt, il
integer(psb_ipk_) :: int_err(5),exch(3)
integer(psb_ipk_), allocatable :: tmpgidx(:,:), &
& nov(:), ov_idx(:,:)
integer(psb_lpk_), allocatable :: vl(:), ix(:), temp_ovrlap(:)
& nov(:), ov_idx(:,:), temp_ovrlap(:)
integer(psb_lpk_), allocatable :: vl(:), ix(:), l_temp_ovrlap(:)
integer(psb_ipk_) :: debug_level, debug_unit
integer(psb_mpk_) :: iictxt
logical :: check_, islarge
@ -304,7 +304,7 @@ subroutine psb_cd_inloc(v, ictxt, desc, info, globalcheck,idx)
end if
! allocate work vector
allocate(temp_ovrlap(max(1,2*loc_row)),desc%lprm(1),&
allocate(l_temp_ovrlap(max(1,2*loc_row)),desc%lprm(1),&
& stat=info)
if (info == psb_success_) then
desc%lprm(1) = 0
@ -316,8 +316,7 @@ subroutine psb_cd_inloc(v, ictxt, desc, info, globalcheck,idx)
goto 9999
endif
temp_ovrlap(:) = -1
l_temp_ovrlap(:) = -1
if (debug_level >= psb_debug_ext_) &
& write(debug_unit,*) me,' ',trim(name),': starting main loop' ,info
@ -338,17 +337,17 @@ subroutine psb_cd_inloc(v, ictxt, desc, info, globalcheck,idx)
if (ov_idx(j,1) == i) exit
j = j + 1
end do
call psb_ensure_size((itmpov+3+nprocs),temp_ovrlap,info,pad=-1_psb_lpk_)
call psb_ensure_size((itmpov+3+nprocs),l_temp_ovrlap,info,pad=-1_psb_lpk_)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='psb_ensure_size')
goto 9999
end if
itmpov = itmpov + 1
temp_ovrlap(itmpov) = il
l_temp_ovrlap(itmpov) = il
itmpov = itmpov + 1
temp_ovrlap(itmpov) = nprocs
temp_ovrlap(itmpov+1:itmpov+nprocs) = ov_idx(j:j+nprocs-1,2)
l_temp_ovrlap(itmpov) = nprocs
l_temp_ovrlap(itmpov+1:itmpov+nprocs) = ov_idx(j:j+nprocs-1,2)
itmpov = itmpov + nprocs
end if
@ -371,7 +370,29 @@ subroutine psb_cd_inloc(v, ictxt, desc, info, globalcheck,idx)
call aa%init(iictxt,vl(1:nlu),info)
end select
call psi_bld_tmpovrl(temp_ovrlap,desc,info)
!
! Now that we have initialized indxmap we can convert the
! indices to local numbering.
!
block
integer(psb_ipk_) :: i,nprocs
allocate(temp_ovrlap(size(l_temp_ovrlap)),stat=info)
if (info == psb_success_) then
temp_ovrlap = -1
i = 1
do while (l_temp_ovrlap(i) /= -1)
call desc%indxmap%g2l(l_temp_ovrlap(i),temp_ovrlap(i),info)
i = i + 1
temp_ovrlap(i) = l_temp_ovrlap(i)
nprocs = temp_ovrlap(i)
temp_ovrlap(i+1:i+nprocs) = l_temp_ovrlap(i+1:i+nprocs)
i = i + 1
i = i + nprocs
enddo
end if
end block
if (info == psb_success_) call psi_bld_tmpovrl(temp_ovrlap,desc,info)
if (info == psb_success_) deallocate(temp_ovrlap,vl,ix,stat=info)
if ((info == psb_success_).and.(allocated(tmpgidx)))&

@ -45,8 +45,7 @@ Subroutine psb_cd_reinit(desc,info)
! .. Local Scalars ..
integer(psb_ipk_) :: np, me, ictxt
integer(psb_ipk_), allocatable :: tmp_halo(:),tmp_ext(:)
integer(psb_lpk_), allocatable :: tmp_ovr(:)
integer(psb_ipk_), allocatable :: tmp_halo(:),tmp_ext(:), tmp_ovr(:)
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name, ch_err
@ -62,11 +61,11 @@ Subroutine psb_cd_reinit(desc,info)
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),': start'
if (desc%is_asb()) then
call psb_cd_get_recv_idx_glob(tmp_ovr,desc,psb_comm_ovr_,info)
call psb_cd_get_recv_idx_loc(tmp_halo,desc,psb_comm_halo_,info)
call psb_cd_get_recv_idx_loc(tmp_ext,desc,psb_comm_ext_,info)
call psb_cd_get_recv_idx(tmp_ovr,desc,psb_comm_ovr_,info)
call psb_cd_get_recv_idx(tmp_halo,desc,psb_comm_halo_,info)
call psb_cd_get_recv_idx(tmp_ext,desc,psb_comm_ext_,info)
call psb_move_alloc(tmp_ovr,desc%tmp_ovrlap_index,info)
call psb_move_alloc(tmp_ovr,desc%ovrlap_index,info)
call psb_move_alloc(tmp_halo,desc%halo_index,info)
call psb_move_alloc(tmp_ext,desc%ext_index,info)
call desc%indxmap%reinit(info)

@ -54,7 +54,7 @@ subroutine psb_cdall(ictxt, desc, info,mg,ng,parts,vg,vl,flag,nl,repl, globalche
character(len=20) :: name
integer(psb_ipk_) :: err_act, flag_, i, me, np, nlp, nnv, lr
integer(psb_lpk_) :: n_
integer(psb_lpk_), allocatable :: itmpsz(:)
integer(psb_ipk_), allocatable :: itmpsz(:)
integer(psb_mpk_) :: iictxt

@ -62,8 +62,8 @@ subroutine psb_cdals(m, n, parts, ictxt, desc, info)
& l_ov_ix,l_ov_el,idx, err_act, itmpov, k, glx, nlx
integer(psb_lpk_) :: iglob
integer(psb_ipk_) :: int_err(5),exch(3)
integer(psb_lpk_), allocatable :: loc_idx(:)
integer(psb_lpk_), allocatable :: temp_ovrlap(:)
integer(psb_ipk_), allocatable :: temp_ovrlap(:)
integer(psb_lpk_), allocatable :: l_temp_ovrlap(:), loc_idx(:)
integer(psb_ipk_), allocatable :: prc_v(:)
integer(psb_ipk_) :: debug_level, debug_unit
integer(psb_ipk_) :: me, np, nprocs
@ -125,7 +125,7 @@ subroutine psb_cdals(m, n, parts, ictxt, desc, info)
! count local rows number
loc_row = max(1,(m+np-1)/np)
! allocate work vector
allocate(temp_ovrlap(max(1,2*loc_row)), prc_v(np),stat=info)
allocate(l_temp_ovrlap(max(1,2*loc_row)), prc_v(np),stat=info)
if (info /= psb_success_) then
info=psb_err_alloc_request_
@ -139,7 +139,7 @@ subroutine psb_cdals(m, n, parts, ictxt, desc, info)
& write(debug_unit,*) me,' ',trim(name),': starting main loop' ,info
counter = 0
itmpov = 0
temp_ovrlap(:) = -1
l_temp_ovrlap(:) = -1
!
! We have to decide whether we have a "large" index space.
!
@ -230,17 +230,17 @@ subroutine psb_cdals(m, n, parts, ictxt, desc, info)
loc_idx(k) = iglob
if (nprocs > 1) then
call psb_ensure_size((itmpov+3+nprocs),temp_ovrlap,info,pad=-1_psb_lpk_)
call psb_ensure_size((itmpov+3+nprocs),l_temp_ovrlap,info,pad=-1_psb_lpk_)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='psb_ensure_size')
goto 9999
end if
itmpov = itmpov + 1
temp_ovrlap(itmpov) = iglob
l_temp_ovrlap(itmpov) = iglob
itmpov = itmpov + 1
temp_ovrlap(itmpov) = nprocs
temp_ovrlap(itmpov+1:itmpov+nprocs) = prc_v(1:nprocs)
l_temp_ovrlap(itmpov) = nprocs
l_temp_ovrlap(itmpov+1:itmpov+nprocs) = prc_v(1:nprocs)
itmpov = itmpov + nprocs
endif
end if
@ -269,9 +269,28 @@ subroutine psb_cdals(m, n, parts, ictxt, desc, info)
if (debug_level >= psb_debug_ext_) &
& write(debug_unit,*) me,' ',trim(name),': error check:' ,err
call psi_bld_tmpovrl(temp_ovrlap,desc,info)
!
! Now that we have initialized indxmap we can convert the
! indices to local numbering.
!
block
integer(psb_ipk_) :: i,nprocs
allocate(temp_ovrlap(size(l_temp_ovrlap)),stat=info)
if (info == psb_success_) then
temp_ovrlap = -1
i = 1
do while (l_temp_ovrlap(i) /= -1)
call desc%indxmap%g2l(l_temp_ovrlap(i),temp_ovrlap(i),info)
i = i + 1
temp_ovrlap(i) = l_temp_ovrlap(i)
nprocs = temp_ovrlap(i)
temp_ovrlap(i+1:i+nprocs) = l_temp_ovrlap(i+1:i+nprocs)
i = i + 1
i = i + nprocs
enddo
end if
end block
if (info == psb_success_) call psi_bld_tmpovrl(temp_ovrlap,desc,info)
if (info == psb_success_) deallocate(prc_v,temp_ovrlap,stat=info)
if (info /= psb_no_err_) then

@ -62,7 +62,7 @@ subroutine psb_cdalv(v, ictxt, desc, info, flag)
& l_ov_ix,l_ov_el,idx, flag_, err_act
integer(psb_lpk_) :: m,n,i
integer(psb_ipk_) :: int_err(5),exch(3)
integer(psb_lpk_), allocatable :: temp_ovrlap(:)
integer(psb_ipk_), allocatable :: temp_ovrlap(:)
integer(psb_ipk_) :: debug_level, debug_unit
integer(psb_mpk_) :: iictxt
character(len=20) :: name

@ -116,7 +116,7 @@ subroutine psb_cdrep(m, ictxt, desc, info)
integer(psb_ipk_) :: i,np,me,err,err_act
integer(psb_lpk_) :: n
integer(psb_ipk_) :: int_err(5),exch(2)
integer(psb_lpk_) :: thalo(1), tovr(1), text(1)
integer(psb_ipk_) :: thalo(1), tovr(1), text(1)
integer(psb_ipk_) :: debug_level, debug_unit
integer(psb_mpk_) :: iictxt
character(len=20) :: name
@ -178,25 +178,6 @@ subroutine psb_cdrep(m, ictxt, desc, info)
!count local rows number
! allocate work vector
!!$ allocate(desc%matrix_data(psb_mdata_size_),&
!!$ & desc%ovrlap_elem(0,3),stat=info)
!!$ if (info /= psb_success_) then
!!$ info=psb_err_alloc_request_
!!$ int_err(1)=2*m+psb_mdata_size_+1
!!$ call psb_errpush(info,name,i_err=int_err,a_err='integer')
!!$ goto 9999
!!$ endif
!!$ ! If the index space is replicated there's no point in not having
!!$ ! the full map on the current process.
!!$
!!$ desc%matrix_data(psb_m_) = m
!!$ desc%matrix_data(psb_n_) = n
!!$ desc%matrix_data(psb_n_row_) = m
!!$ desc%matrix_data(psb_n_col_) = n
!!$ desc%matrix_data(psb_ctxt_) = ictxt
!!$ call psb_get_mpicomm(ictxt,desc%matrix_data(psb_mpi_c_))
!!$ desc%matrix_data(psb_dec_type_) = psb_desc_bld_
allocate(psb_repl_map :: desc%indxmap, stat=info)
select type(aa => desc%indxmap)

@ -258,25 +258,19 @@ Subroutine psb_dcdbldext(a,desc_a,novr,desc_ov,info, extype)
Do j=0,n_elem_recv-1
idx = ovrlap(counter+psb_elem_recv_+j)
call desc_ov%indxmap%l2g(idx,gidx,info)
If (gidx < 0) then
info=-3
call psb_errpush(info,name)
goto 9999
endif
call psb_ensure_size((cntov_o+3),orig_ovr,info,pad=-ione)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='psb_ensure_size')
goto 9999
end if
orig_ovr(cntov_o)=proc
orig_ovr(cntov_o+1)=1
orig_ovr(cntov_o+2)=gidx
orig_ovr(cntov_o+3)=-1
orig_ovr(cntov_o) = proc
orig_ovr(cntov_o+1) = 1
orig_ovr(cntov_o+2) = idx
orig_ovr(cntov_o+3) = -1
cntov_o=cntov_o+3
end Do
counter=counter+n_elem_recv+n_elem_send+3
counter = counter+n_elem_recv+n_elem_send+3
end Do
@ -327,16 +321,16 @@ Subroutine psb_dcdbldext(a,desc_a,novr,desc_ov,info, extype)
n_col_prev = desc_ov%get_local_cols()
Do While (halo(counter) /= -1)
tot_elem=0
proc=halo(counter+psb_proc_id_)
n_elem_recv=halo(counter+psb_n_elem_recv_)
n_elem_send=halo(counter+n_elem_recv+psb_n_elem_send_)
tot_elem = 0
proc = halo(counter+psb_proc_id_)
n_elem_recv = halo(counter+psb_n_elem_recv_)
n_elem_send = halo(counter+n_elem_recv+psb_n_elem_send_)
If ((counter+n_elem_recv+n_elem_send) > Size(halo)) then
info = -1
call psb_errpush(info,name)
goto 9999
end If
tot_recv=tot_recv+n_elem_recv
tot_recv = tot_recv+n_elem_recv
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),&
& ': tot_recv:',proc,n_elem_recv,tot_recv
@ -359,12 +353,6 @@ Subroutine psb_dcdbldext(a,desc_a,novr,desc_ov,info, extype)
end If
idx = halo(counter+psb_elem_recv_+j)
call desc_ov%l2g(idx,gidx,info)
If (gidx < 0) then
info=-3
call psb_errpush(info,name)
goto 9999
endif
call psb_ensure_size((counter_o+3),tmp_ovr_idx,info,pad=-ione)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
@ -374,7 +362,7 @@ Subroutine psb_dcdbldext(a,desc_a,novr,desc_ov,info, extype)
tmp_ovr_idx(counter_o) = proc
tmp_ovr_idx(counter_o+1) = 1
tmp_ovr_idx(counter_o+2) = gidx
tmp_ovr_idx(counter_o+2) = idx
tmp_ovr_idx(counter_o+3) = -1
counter_o=counter_o+3
call psb_ensure_size((counter_h+3),tmp_halo,info,pad=-ione)
@ -403,12 +391,6 @@ Subroutine psb_dcdbldext(a,desc_a,novr,desc_ov,info, extype)
Do j=0,n_elem_send-1
idx = halo(counter+psb_elem_send_+j)
call desc_ov%l2g(idx,gidx,info)
If (gidx < 0) then
info=-3
call psb_errpush(info,name)
goto 9999
endif
call psb_ensure_size((counter_o+3),tmp_ovr_idx,info,pad=-ione)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
@ -418,7 +400,7 @@ Subroutine psb_dcdbldext(a,desc_a,novr,desc_ov,info, extype)
tmp_ovr_idx(counter_o) = proc
tmp_ovr_idx(counter_o+1) = 1
tmp_ovr_idx(counter_o+2) = gidx
tmp_ovr_idx(counter_o+2) = idx
tmp_ovr_idx(counter_o+3) = -1
counter_o=counter_o+3
@ -427,6 +409,7 @@ Subroutine psb_dcdbldext(a,desc_a,novr,desc_ov,info, extype)
!
If (i_ovr <= (novr)) Then
call a%csget(idx,idx,n_elem,irow,icol,info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='csget')
@ -507,8 +490,8 @@ Subroutine psb_dcdbldext(a,desc_a,novr,desc_ov,info, extype)
lworkr = max(iszr,1)
end if
call mpi_alltoallv(works,sdsz,bsdindx,psb_mpi_ipk_,&
& workr,rvsz,brvindx,psb_mpi_ipk_,icomm,minfo)
call mpi_alltoallv(works,sdsz,bsdindx,psb_mpi_lpk_,&
& workr,rvsz,brvindx,psb_mpi_lpk_,icomm,minfo)
if (minfo /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='mpi_alltoallv')
@ -534,7 +517,7 @@ Subroutine psb_dcdbldext(a,desc_a,novr,desc_ov,info, extype)
j = 0
do i=1,iszr
if (maskr(i) < 0) then
j=j+1
j = j+1
works(j) = workr(i)
end if
end do
@ -557,8 +540,8 @@ Subroutine psb_dcdbldext(a,desc_a,novr,desc_ov,info, extype)
& ': Done fnd_owner', desc_ov%indxmap%get_state()
do i=1,iszs
gidx = works(i)
n_col = desc_ov%get_local_cols()
gidx = works(i)
n_col = desc_ov%get_local_cols()
call desc_ov%indxmap%g2l_ins(gidx,lidx,info)
if (desc_ov%get_local_cols() > n_col ) then
!
@ -603,7 +586,7 @@ Subroutine psb_dcdbldext(a,desc_a,novr,desc_ov,info, extype)
write(debug_unit,*) me,' ',trim(name),':Calling Crea_index'
end if
call psi_crea_index(desc_ov,t_halo_in,t_halo_out,.false.,&
call psi_crea_index(desc_ov,t_halo_in,t_halo_out,&
& nxch,nsnd,nrcv,info)
if (debug_level >= psb_debug_outer_) then
@ -661,7 +644,7 @@ Subroutine psb_dcdbldext(a,desc_a,novr,desc_ov,info, extype)
goto 9999
end if
orig_ovr(cntov_o:cntov_o+counter_o-1) = tmp_ovr_idx(1:counter_o)
cntov_o = cntov_o+counter_o-1
cntov_o = cntov_o+counter_o-1
orig_ovr(cntov_o:) = -1
call psb_move_alloc(orig_ovr,desc_ov%ovrlap_index,info)
deallocate(tmp_ovr_idx,stat=info)

@ -258,25 +258,19 @@ Subroutine psb_scdbldext(a,desc_a,novr,desc_ov,info, extype)
Do j=0,n_elem_recv-1
idx = ovrlap(counter+psb_elem_recv_+j)
call desc_ov%indxmap%l2g(idx,gidx,info)
If (gidx < 0) then
info=-3
call psb_errpush(info,name)
goto 9999
endif
call psb_ensure_size((cntov_o+3),orig_ovr,info,pad=-ione)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='psb_ensure_size')
goto 9999
end if
orig_ovr(cntov_o)=proc
orig_ovr(cntov_o+1)=1
orig_ovr(cntov_o+2)=gidx
orig_ovr(cntov_o+3)=-1
orig_ovr(cntov_o) = proc
orig_ovr(cntov_o+1) = 1
orig_ovr(cntov_o+2) = idx
orig_ovr(cntov_o+3) = -1
cntov_o=cntov_o+3
end Do
counter=counter+n_elem_recv+n_elem_send+3
counter = counter+n_elem_recv+n_elem_send+3
end Do
@ -327,16 +321,16 @@ Subroutine psb_scdbldext(a,desc_a,novr,desc_ov,info, extype)
n_col_prev = desc_ov%get_local_cols()
Do While (halo(counter) /= -1)
tot_elem=0
proc=halo(counter+psb_proc_id_)
n_elem_recv=halo(counter+psb_n_elem_recv_)
n_elem_send=halo(counter+n_elem_recv+psb_n_elem_send_)
tot_elem = 0
proc = halo(counter+psb_proc_id_)
n_elem_recv = halo(counter+psb_n_elem_recv_)
n_elem_send = halo(counter+n_elem_recv+psb_n_elem_send_)
If ((counter+n_elem_recv+n_elem_send) > Size(halo)) then
info = -1
call psb_errpush(info,name)
goto 9999
end If
tot_recv=tot_recv+n_elem_recv
tot_recv = tot_recv+n_elem_recv
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),&
& ': tot_recv:',proc,n_elem_recv,tot_recv
@ -359,12 +353,6 @@ Subroutine psb_scdbldext(a,desc_a,novr,desc_ov,info, extype)
end If
idx = halo(counter+psb_elem_recv_+j)
call desc_ov%l2g(idx,gidx,info)
If (gidx < 0) then
info=-3
call psb_errpush(info,name)
goto 9999
endif
call psb_ensure_size((counter_o+3),tmp_ovr_idx,info,pad=-ione)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
@ -374,7 +362,7 @@ Subroutine psb_scdbldext(a,desc_a,novr,desc_ov,info, extype)
tmp_ovr_idx(counter_o) = proc
tmp_ovr_idx(counter_o+1) = 1
tmp_ovr_idx(counter_o+2) = gidx
tmp_ovr_idx(counter_o+2) = idx
tmp_ovr_idx(counter_o+3) = -1
counter_o=counter_o+3
call psb_ensure_size((counter_h+3),tmp_halo,info,pad=-ione)
@ -403,12 +391,6 @@ Subroutine psb_scdbldext(a,desc_a,novr,desc_ov,info, extype)
Do j=0,n_elem_send-1
idx = halo(counter+psb_elem_send_+j)
call desc_ov%l2g(idx,gidx,info)
If (gidx < 0) then
info=-3
call psb_errpush(info,name)
goto 9999
endif
call psb_ensure_size((counter_o+3),tmp_ovr_idx,info,pad=-ione)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
@ -418,7 +400,7 @@ Subroutine psb_scdbldext(a,desc_a,novr,desc_ov,info, extype)
tmp_ovr_idx(counter_o) = proc
tmp_ovr_idx(counter_o+1) = 1
tmp_ovr_idx(counter_o+2) = gidx
tmp_ovr_idx(counter_o+2) = idx
tmp_ovr_idx(counter_o+3) = -1
counter_o=counter_o+3
@ -427,6 +409,7 @@ Subroutine psb_scdbldext(a,desc_a,novr,desc_ov,info, extype)
!
If (i_ovr <= (novr)) Then
call a%csget(idx,idx,n_elem,irow,icol,info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='csget')
@ -507,8 +490,8 @@ Subroutine psb_scdbldext(a,desc_a,novr,desc_ov,info, extype)
lworkr = max(iszr,1)
end if
call mpi_alltoallv(works,sdsz,bsdindx,psb_mpi_ipk_,&
& workr,rvsz,brvindx,psb_mpi_ipk_,icomm,minfo)
call mpi_alltoallv(works,sdsz,bsdindx,psb_mpi_lpk_,&
& workr,rvsz,brvindx,psb_mpi_lpk_,icomm,minfo)
if (minfo /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='mpi_alltoallv')
@ -534,7 +517,7 @@ Subroutine psb_scdbldext(a,desc_a,novr,desc_ov,info, extype)
j = 0
do i=1,iszr
if (maskr(i) < 0) then
j=j+1
j = j+1
works(j) = workr(i)
end if
end do
@ -557,8 +540,8 @@ Subroutine psb_scdbldext(a,desc_a,novr,desc_ov,info, extype)
& ': Done fnd_owner', desc_ov%indxmap%get_state()
do i=1,iszs
gidx = works(i)
n_col = desc_ov%get_local_cols()
gidx = works(i)
n_col = desc_ov%get_local_cols()
call desc_ov%indxmap%g2l_ins(gidx,lidx,info)
if (desc_ov%get_local_cols() > n_col ) then
!
@ -603,7 +586,7 @@ Subroutine psb_scdbldext(a,desc_a,novr,desc_ov,info, extype)
write(debug_unit,*) me,' ',trim(name),':Calling Crea_index'
end if
call psi_crea_index(desc_ov,t_halo_in,t_halo_out,.false.,&
call psi_crea_index(desc_ov,t_halo_in,t_halo_out,&
& nxch,nsnd,nrcv,info)
if (debug_level >= psb_debug_outer_) then
@ -661,7 +644,7 @@ Subroutine psb_scdbldext(a,desc_a,novr,desc_ov,info, extype)
goto 9999
end if
orig_ovr(cntov_o:cntov_o+counter_o-1) = tmp_ovr_idx(1:counter_o)
cntov_o = cntov_o+counter_o-1
cntov_o = cntov_o+counter_o-1
orig_ovr(cntov_o:) = -1
call psb_move_alloc(orig_ovr,desc_ov%ovrlap_index,info)
deallocate(tmp_ovr_idx,stat=info)

@ -258,25 +258,19 @@ Subroutine psb_zcdbldext(a,desc_a,novr,desc_ov,info, extype)
Do j=0,n_elem_recv-1
idx = ovrlap(counter+psb_elem_recv_+j)
call desc_ov%indxmap%l2g(idx,gidx,info)
If (gidx < 0) then
info=-3
call psb_errpush(info,name)
goto 9999
endif
call psb_ensure_size((cntov_o+3),orig_ovr,info,pad=-ione)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='psb_ensure_size')
goto 9999
end if
orig_ovr(cntov_o)=proc
orig_ovr(cntov_o+1)=1
orig_ovr(cntov_o+2)=gidx
orig_ovr(cntov_o+3)=-1
orig_ovr(cntov_o) = proc
orig_ovr(cntov_o+1) = 1
orig_ovr(cntov_o+2) = idx
orig_ovr(cntov_o+3) = -1
cntov_o=cntov_o+3
end Do
counter=counter+n_elem_recv+n_elem_send+3
counter = counter+n_elem_recv+n_elem_send+3
end Do
@ -327,16 +321,16 @@ Subroutine psb_zcdbldext(a,desc_a,novr,desc_ov,info, extype)
n_col_prev = desc_ov%get_local_cols()
Do While (halo(counter) /= -1)
tot_elem=0
proc=halo(counter+psb_proc_id_)
n_elem_recv=halo(counter+psb_n_elem_recv_)
n_elem_send=halo(counter+n_elem_recv+psb_n_elem_send_)
tot_elem = 0
proc = halo(counter+psb_proc_id_)
n_elem_recv = halo(counter+psb_n_elem_recv_)
n_elem_send = halo(counter+n_elem_recv+psb_n_elem_send_)
If ((counter+n_elem_recv+n_elem_send) > Size(halo)) then
info = -1
call psb_errpush(info,name)
goto 9999
end If
tot_recv=tot_recv+n_elem_recv
tot_recv = tot_recv+n_elem_recv
if (debug_level >= psb_debug_outer_) &
& write(debug_unit,*) me,' ',trim(name),&
& ': tot_recv:',proc,n_elem_recv,tot_recv
@ -359,12 +353,6 @@ Subroutine psb_zcdbldext(a,desc_a,novr,desc_ov,info, extype)
end If
idx = halo(counter+psb_elem_recv_+j)
call desc_ov%l2g(idx,gidx,info)
If (gidx < 0) then
info=-3
call psb_errpush(info,name)
goto 9999
endif
call psb_ensure_size((counter_o+3),tmp_ovr_idx,info,pad=-ione)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
@ -374,7 +362,7 @@ Subroutine psb_zcdbldext(a,desc_a,novr,desc_ov,info, extype)
tmp_ovr_idx(counter_o) = proc
tmp_ovr_idx(counter_o+1) = 1
tmp_ovr_idx(counter_o+2) = gidx
tmp_ovr_idx(counter_o+2) = idx
tmp_ovr_idx(counter_o+3) = -1
counter_o=counter_o+3
call psb_ensure_size((counter_h+3),tmp_halo,info,pad=-ione)
@ -403,12 +391,6 @@ Subroutine psb_zcdbldext(a,desc_a,novr,desc_ov,info, extype)
Do j=0,n_elem_send-1
idx = halo(counter+psb_elem_send_+j)
call desc_ov%l2g(idx,gidx,info)
If (gidx < 0) then
info=-3
call psb_errpush(info,name)
goto 9999
endif
call psb_ensure_size((counter_o+3),tmp_ovr_idx,info,pad=-ione)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
@ -418,7 +400,7 @@ Subroutine psb_zcdbldext(a,desc_a,novr,desc_ov,info, extype)
tmp_ovr_idx(counter_o) = proc
tmp_ovr_idx(counter_o+1) = 1
tmp_ovr_idx(counter_o+2) = gidx
tmp_ovr_idx(counter_o+2) = idx
tmp_ovr_idx(counter_o+3) = -1
counter_o=counter_o+3
@ -427,6 +409,7 @@ Subroutine psb_zcdbldext(a,desc_a,novr,desc_ov,info, extype)
!
If (i_ovr <= (novr)) Then
call a%csget(idx,idx,n_elem,irow,icol,info)
if (info /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='csget')
@ -507,8 +490,8 @@ Subroutine psb_zcdbldext(a,desc_a,novr,desc_ov,info, extype)
lworkr = max(iszr,1)
end if
call mpi_alltoallv(works,sdsz,bsdindx,psb_mpi_ipk_,&
& workr,rvsz,brvindx,psb_mpi_ipk_,icomm,minfo)
call mpi_alltoallv(works,sdsz,bsdindx,psb_mpi_lpk_,&
& workr,rvsz,brvindx,psb_mpi_lpk_,icomm,minfo)
if (minfo /= psb_success_) then
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='mpi_alltoallv')
@ -534,7 +517,7 @@ Subroutine psb_zcdbldext(a,desc_a,novr,desc_ov,info, extype)
j = 0
do i=1,iszr
if (maskr(i) < 0) then
j=j+1
j = j+1
works(j) = workr(i)
end if
end do
@ -557,8 +540,8 @@ Subroutine psb_zcdbldext(a,desc_a,novr,desc_ov,info, extype)
& ': Done fnd_owner', desc_ov%indxmap%get_state()
do i=1,iszs
gidx = works(i)
n_col = desc_ov%get_local_cols()
gidx = works(i)
n_col = desc_ov%get_local_cols()
call desc_ov%indxmap%g2l_ins(gidx,lidx,info)
if (desc_ov%get_local_cols() > n_col ) then
!
@ -603,7 +586,7 @@ Subroutine psb_zcdbldext(a,desc_a,novr,desc_ov,info, extype)
write(debug_unit,*) me,' ',trim(name),':Calling Crea_index'
end if
call psi_crea_index(desc_ov,t_halo_in,t_halo_out,.false.,&
call psi_crea_index(desc_ov,t_halo_in,t_halo_out,&
& nxch,nsnd,nrcv,info)
if (debug_level >= psb_debug_outer_) then
@ -661,7 +644,7 @@ Subroutine psb_zcdbldext(a,desc_a,novr,desc_ov,info, extype)
goto 9999
end if
orig_ovr(cntov_o:cntov_o+counter_o-1) = tmp_ovr_idx(1:counter_o)
cntov_o = cntov_o+counter_o-1
cntov_o = cntov_o+counter_o-1
orig_ovr(cntov_o:) = -1
call psb_move_alloc(orig_ovr,desc_ov%ovrlap_index,info)
deallocate(tmp_ovr_idx,stat=info)

@ -86,6 +86,56 @@ contains
end function d_null_func_2d
!
! functions parametrizing the differential equation
!
function b1(x,y)
use psb_base_mod, only : psb_dpk_, done, dzero
implicit none
real(psb_dpk_) :: b1
real(psb_dpk_), intent(in) :: x,y
b1=done/sqrt((2*done))
end function b1
function b2(x,y)
use psb_base_mod, only : psb_dpk_, done, dzero
implicit none
real(psb_dpk_) :: b2
real(psb_dpk_), intent(in) :: x,y
b2=done/sqrt((2*done))
end function b2
function c(x,y)
use psb_base_mod, only : psb_dpk_, done, dzero
implicit none
real(psb_dpk_) :: c
real(psb_dpk_), intent(in) :: x,y
c=0.d0
end function c
function a1(x,y)
use psb_base_mod, only : psb_dpk_, done, dzero
implicit none
real(psb_dpk_) :: a1
real(psb_dpk_), intent(in) :: x,y
a1=done/80
end function a1
function a2(x,y)
use psb_base_mod, only : psb_dpk_, done, dzero
implicit none
real(psb_dpk_) :: a2
real(psb_dpk_), intent(in) :: x,y
a2=done/80
end function a2
function g(x,y)
use psb_base_mod, only : psb_dpk_, done, dzero
implicit none
real(psb_dpk_) :: g
real(psb_dpk_), intent(in) :: x,y
g = dzero
if (x == done) then
g = done
else if (x == dzero) then
g = exp(-y**2)
end if
end function g
!
@ -465,51 +515,6 @@ contains
return
end subroutine psb_d_gen_pde2d
!
! functions parametrizing the differential equation
!
function b1(x,y)
use psb_base_mod, only : psb_dpk_
real(psb_dpk_) :: b1
real(psb_dpk_), intent(in) :: x,y
b1=done/sqrt((2*done))
end function b1
function b2(x,y)
use psb_base_mod, only : psb_dpk_
real(psb_dpk_) :: b2
real(psb_dpk_), intent(in) :: x,y
b2=done/sqrt((2*done))
end function b2
function c(x,y)
use psb_base_mod, only : psb_dpk_
real(psb_dpk_) :: c
real(psb_dpk_), intent(in) :: x,y
c=0.d0
end function c
function a1(x,y)
use psb_base_mod, only : psb_dpk_
real(psb_dpk_) :: a1
real(psb_dpk_), intent(in) :: x,y
a1=done/80
end function a1
function a2(x,y)
use psb_base_mod, only : psb_dpk_
real(psb_dpk_) :: a2
real(psb_dpk_), intent(in) :: x,y
a2=done/80
end function a2
function g(x,y)
use psb_base_mod, only : psb_dpk_, done, dzero
real(psb_dpk_) :: g
real(psb_dpk_), intent(in) :: x,y
g = dzero
if (x == done) then
g = done
else if (x == dzero) then
g = exp(-y**2)
end if
end function g
end module psb_d_pde2d_mod
program psb_d_pde2d
@ -787,6 +792,7 @@ contains
write(iout,*)' >= 1 do tracing every itrace'
write(iout,*)' iterations '
end subroutine pr_usage
end program psb_d_pde2d

@ -88,6 +88,70 @@ contains
val = dzero
end function d_null_func_3d
!
! functions parametrizing the differential equation
!
function b1(x,y,z)
use psb_base_mod, only : psb_dpk_, done, dzero
implicit none
real(psb_dpk_) :: b1
real(psb_dpk_), intent(in) :: x,y,z
b1=done/sqrt((3*done))
end function b1
function b2(x,y,z)
use psb_base_mod, only : psb_dpk_, done, dzero
implicit none
real(psb_dpk_) :: b2
real(psb_dpk_), intent(in) :: x,y,z
b2=done/sqrt((3*done))
end function b2
function b3(x,y,z)
use psb_base_mod, only : psb_dpk_, done, dzero
implicit none
real(psb_dpk_) :: b3
real(psb_dpk_), intent(in) :: x,y,z
b3=done/sqrt((3*done))
end function b3
function c(x,y,z)
use psb_base_mod, only : psb_dpk_, done, dzero
implicit none
real(psb_dpk_) :: c
real(psb_dpk_), intent(in) :: x,y,z
c=dzero
end function c
function a1(x,y,z)
use psb_base_mod, only : psb_dpk_, done, dzero
implicit none
real(psb_dpk_) :: a1
real(psb_dpk_), intent(in) :: x,y,z
a1=done/80
end function a1
function a2(x,y,z)
use psb_base_mod, only : psb_dpk_, done, dzero
implicit none
real(psb_dpk_) :: a2
real(psb_dpk_), intent(in) :: x,y,z
a2=done/80
end function a2
function a3(x,y,z)
use psb_base_mod, only : psb_dpk_, done, dzero
implicit none
real(psb_dpk_) :: a3
real(psb_dpk_), intent(in) :: x,y,z
a3=done/80
end function a3
function g(x,y,z)
use psb_base_mod, only : psb_dpk_, done, dzero
implicit none
real(psb_dpk_) :: g
real(psb_dpk_), intent(in) :: x,y,z
g = dzero
if (x == done) then
g = done
else if (x == dzero) then
g = exp(y**2-z**2)
end if
end function g
!
@ -114,7 +178,6 @@ contains
! Note that if b1=b2=b3=c=0., the PDE is the Laplace equation.
!
implicit none
! procedure(d_func_3d) :: b1,b2,b3,c,a1,a2,a3,g
integer(psb_ipk_) :: idim
type(psb_dspmat_type) :: a
type(psb_d_vect_type) :: xv,bv
@ -491,62 +554,7 @@ contains
return
end subroutine psb_d_gen_pde3d
!
! functions parametrizing the differential equation
!
function b1(x,y,z)
use psb_base_mod, only : psb_dpk_
real(psb_dpk_) :: b1
real(psb_dpk_), intent(in) :: x,y,z
b1=done/sqrt((3*done))
end function b1
function b2(x,y,z)
use psb_base_mod, only : psb_dpk_
real(psb_dpk_) :: b2
real(psb_dpk_), intent(in) :: x,y,z
b2=done/sqrt((3*done))
end function b2
function b3(x,y,z)
use psb_base_mod, only : psb_dpk_
real(psb_dpk_) :: b3
real(psb_dpk_), intent(in) :: x,y,z
b3=done/sqrt((3*done))
end function b3
function c(x,y,z)
use psb_base_mod, only : psb_dpk_
real(psb_dpk_) :: c
real(psb_dpk_), intent(in) :: x,y,z
c=dzero
end function c
function a1(x,y,z)
use psb_base_mod, only : psb_dpk_
real(psb_dpk_) :: a1
real(psb_dpk_), intent(in) :: x,y,z
a1=done/80
end function a1
function a2(x,y,z)
use psb_base_mod, only : psb_dpk_
real(psb_dpk_) :: a2
real(psb_dpk_), intent(in) :: x,y,z
a2=done/80
end function a2
function a3(x,y,z)
use psb_base_mod, only : psb_dpk_
real(psb_dpk_) :: a3
real(psb_dpk_), intent(in) :: x,y,z
a3=done/80
end function a3
function g(x,y,z)
use psb_base_mod, only : psb_dpk_, done, dzero
real(psb_dpk_) :: g
real(psb_dpk_), intent(in) :: x,y,z
g = dzero
if (x == done) then
g = done
else if (x == dzero) then
g = exp(y**2-z**2)
end if
end function g
end module psb_d_pde3d_mod

@ -86,6 +86,56 @@ contains
end function s_null_func_2d
!
! functions parametrizing the differential equation
!
function b1(x,y)
use psb_base_mod, only : psb_spk_, sone, szero
implicit none
real(psb_spk_) :: b1
real(psb_spk_), intent(in) :: x,y
b1=sone/sqrt((2*sone))
end function b1
function b2(x,y)
use psb_base_mod, only : psb_spk_, sone, szero
implicit none
real(psb_spk_) :: b2
real(psb_spk_), intent(in) :: x,y
b2=sone/sqrt((2*sone))
end function b2
function c(x,y)
use psb_base_mod, only : psb_spk_, sone, szero
implicit none
real(psb_spk_) :: c
real(psb_spk_), intent(in) :: x,y
c=0.d0
end function c
function a1(x,y)
use psb_base_mod, only : psb_spk_, sone, szero
implicit none
real(psb_spk_) :: a1
real(psb_spk_), intent(in) :: x,y
a1=sone/80
end function a1
function a2(x,y)
use psb_base_mod, only : psb_spk_, sone, szero
implicit none
real(psb_spk_) :: a2
real(psb_spk_), intent(in) :: x,y
a2=sone/80
end function a2
function g(x,y)
use psb_base_mod, only : psb_spk_, sone, szero
implicit none
real(psb_spk_) :: g
real(psb_spk_), intent(in) :: x,y
g = szero
if (x == sone) then
g = sone
else if (x == szero) then
g = exp(-y**2)
end if
end function g
!
@ -465,51 +515,6 @@ contains
return
end subroutine psb_s_gen_pde2d
!
! functions parametrizing the differential equation
!
function b1(x,y)
use psb_base_mod, only : psb_spk_
real(psb_spk_) :: b1
real(psb_spk_), intent(in) :: x,y
b1=sone/sqrt((2*sone))
end function b1
function b2(x,y)
use psb_base_mod, only : psb_spk_
real(psb_spk_) :: b2
real(psb_spk_), intent(in) :: x,y
b2=sone/sqrt((2*sone))
end function b2
function c(x,y)
use psb_base_mod, only : psb_spk_
real(psb_spk_) :: c
real(psb_spk_), intent(in) :: x,y
c=0.d0
end function c
function a1(x,y)
use psb_base_mod, only : psb_spk_
real(psb_spk_) :: a1
real(psb_spk_), intent(in) :: x,y
a1=sone/80
end function a1
function a2(x,y)
use psb_base_mod, only : psb_spk_
real(psb_spk_) :: a2
real(psb_spk_), intent(in) :: x,y
a2=sone/80
end function a2
function g(x,y)
use psb_base_mod, only : psb_spk_, sone, szero
real(psb_spk_) :: g
real(psb_spk_), intent(in) :: x,y
g = szero
if (x == sone) then
g = sone
else if (x == szero) then
g = exp(-y**2)
end if
end function g
end module psb_s_pde2d_mod
program psb_s_pde2d
@ -787,6 +792,7 @@ contains
write(iout,*)' >= 1 do tracing every itrace'
write(iout,*)' iterations '
end subroutine pr_usage
end program psb_s_pde2d

@ -88,6 +88,70 @@ contains
val = szero
end function s_null_func_3d
!
! functions parametrizing the differential equation
!
function b1(x,y,z)
use psb_base_mod, only : psb_spk_, sone, szero
implicit none
real(psb_spk_) :: b1
real(psb_spk_), intent(in) :: x,y,z
b1=sone/sqrt((3*sone))
end function b1
function b2(x,y,z)
use psb_base_mod, only : psb_spk_, sone, szero
implicit none
real(psb_spk_) :: b2
real(psb_spk_), intent(in) :: x,y,z
b2=sone/sqrt((3*sone))
end function b2
function b3(x,y,z)
use psb_base_mod, only : psb_spk_, sone, szero
implicit none
real(psb_spk_) :: b3
real(psb_spk_), intent(in) :: x,y,z
b3=sone/sqrt((3*sone))
end function b3
function c(x,y,z)
use psb_base_mod, only : psb_spk_, sone, szero
implicit none
real(psb_spk_) :: c
real(psb_spk_), intent(in) :: x,y,z
c=szero
end function c
function a1(x,y,z)
use psb_base_mod, only : psb_spk_, sone, szero
implicit none
real(psb_spk_) :: a1
real(psb_spk_), intent(in) :: x,y,z
a1=sone/80
end function a1
function a2(x,y,z)
use psb_base_mod, only : psb_spk_, sone, szero
implicit none
real(psb_spk_) :: a2
real(psb_spk_), intent(in) :: x,y,z
a2=sone/80
end function a2
function a3(x,y,z)
use psb_base_mod, only : psb_spk_, sone, szero
implicit none
real(psb_spk_) :: a3
real(psb_spk_), intent(in) :: x,y,z
a3=sone/80
end function a3
function g(x,y,z)
use psb_base_mod, only : psb_spk_, sone, szero
implicit none
real(psb_spk_) :: g
real(psb_spk_), intent(in) :: x,y,z
g = szero
if (x == sone) then
g = sone
else if (x == szero) then
g = exp(y**2-z**2)
end if
end function g
!
@ -114,7 +178,6 @@ contains
! Note that if b1=b2=b3=c=0., the PDE is the Laplace equation.
!
implicit none
! procedure(s_func_3d) :: b1,b2,b3,c,a1,a2,a3,g
integer(psb_ipk_) :: idim
type(psb_sspmat_type) :: a
type(psb_s_vect_type) :: xv,bv
@ -491,62 +554,7 @@ contains
return
end subroutine psb_s_gen_pde3d
!
! functions parametrizing the differential equation
!
function b1(x,y,z)
use psb_base_mod, only : psb_spk_
real(psb_spk_) :: b1
real(psb_spk_), intent(in) :: x,y,z
b1=sone/sqrt((3*sone))
end function b1
function b2(x,y,z)
use psb_base_mod, only : psb_spk_
real(psb_spk_) :: b2
real(psb_spk_), intent(in) :: x,y,z
b2=sone/sqrt((3*sone))
end function b2
function b3(x,y,z)
use psb_base_mod, only : psb_spk_
real(psb_spk_) :: b3
real(psb_spk_), intent(in) :: x,y,z
b3=sone/sqrt((3*sone))
end function b3
function c(x,y,z)
use psb_base_mod, only : psb_spk_
real(psb_spk_) :: c
real(psb_spk_), intent(in) :: x,y,z
c=szero
end function c
function a1(x,y,z)
use psb_base_mod, only : psb_spk_
real(psb_spk_) :: a1
real(psb_spk_), intent(in) :: x,y,z
a1=sone/80
end function a1
function a2(x,y,z)
use psb_base_mod, only : psb_spk_
real(psb_spk_) :: a2
real(psb_spk_), intent(in) :: x,y,z
a2=sone/80
end function a2
function a3(x,y,z)
use psb_base_mod, only : psb_spk_
real(psb_spk_) :: a3
real(psb_spk_), intent(in) :: x,y,z
a3=sone/80
end function a3
function g(x,y,z)
use psb_base_mod, only : psb_spk_, sone, szero
real(psb_spk_) :: g
real(psb_spk_), intent(in) :: x,y,z
g = szero
if (x == sone) then
g = sone
else if (x == szero) then
g = exp(y**2-z**2)
end if
end function g
end module psb_s_pde3d_mod

@ -2,7 +2,7 @@
BICGSTAB Iterative method BICGSTAB CGS BICG BICGSTABL RGMRES FCG CGR
BJAC Preconditioner NONE DIAG BJAC
CSR Storage format for matrix A: CSR COO
040 Domain size (acutal system is this**3 (pde3d) or **2 (pde2d) )
004 Domain size (acutal system is this**3 (pde3d) or **2 (pde2d) )
2 Stopping criterion 1 2
1000 MAXIT
-1 ITRACE

Loading…
Cancel
Save