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

Merge branch 'new-context' into dev-openmp

Browse Source
new-context
Salvatore Filippone 4 years ago
parent 926cb1e5ad c58e523918
commit 327003ad06
38 changed files with 5499 additions and 1518 deletions
Show Stats Download Patch File Download Diff File

2
base/modules/Makefile
Unescape Escape View File

@ -257,7 +257,7 @@ serial/psb_c_csc_mat_mod.o serial/psb_c_csr_mat_mod.o serial/psb_lc_csr_mat_mod.
serial/psb_z_csc_mat_mod.o serial/psb_z_csr_mat_mod.o serial/psb_lz_csr_mat_mod.o: serial/psb_z_base_mat_mod.o
serial/psb_mat_mod.o: serial/psb_vect_mod.o serial/psb_s_mat_mod.o serial/psb_d_mat_mod.o serial/psb_c_mat_mod.o serial/psb_z_mat_mod.o
serial/psb_serial_mod.o: serial/psb_s_serial_mod.o serial/psb_d_serial_mod.o serial/psb_c_serial_mod.o serial/psb_z_serial_mod.o
serial/psb_serial_mod.o: serial/psb_s_serial_mod.o serial/psb_d_serial_mod.o serial/psb_c_serial_mod.o serial/psb_z_serial_mod.o auxil/psi_serial_mod.o
serial/psb_i_vect_mod.o: serial/psb_i_base_vect_mod.o
serial/psb_l_vect_mod.o: serial/psb_l_base_vect_mod.o serial/psb_i_vect_mod.o
serial/psb_s_vect_mod.o: serial/psb_s_base_vect_mod.o serial/psb_i_vect_mod.o

34
base/modules/auxil/psi_c_serial_mod.f90
Unescape Escape View File

@ -32,24 +32,40 @@
module psi_c_serial_mod
use psb_const_mod, only : psb_ipk_, psb_lpk_, psb_mpk_, psb_epk_, psb_spk_
interface psb_gelp
interface psb_gelp
! 2-D version
subroutine psb_cgelp(trans,iperm,x,info)
import :: psb_ipk_, psb_spk_
subroutine psb_m_cgelp(trans,iperm,x,info)
import :: psb_ipk_, psb_mpk_, psb_spk_
implicit none
complex(psb_spk_), intent(inout) :: x(:,:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_cgelp
subroutine psb_cgelpv(trans,iperm,x,info)
import :: psb_ipk_, psb_spk_
end subroutine psb_m_cgelp
subroutine psb_m_cgelpv(trans,iperm,x,info)
import :: psb_ipk_, psb_mpk_,psb_spk_
implicit none
complex(psb_spk_), intent(inout) :: x(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_m_cgelpv
subroutine psb_e_cgelp(trans,iperm,x,info)
import :: psb_ipk_, psb_epk_, psb_spk_
implicit none
complex(psb_spk_), intent(inout) :: x(:,:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_e_cgelp
subroutine psb_e_cgelpv(trans,iperm,x,info)
import :: psb_ipk_, psb_epk_, psb_spk_
implicit none
complex(psb_spk_), intent(inout) :: x(:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_cgelpv
end subroutine psb_e_cgelpv
end interface psb_gelp
interface psb_geaxpby

34
base/modules/auxil/psi_d_serial_mod.f90
Unescape Escape View File

@ -32,24 +32,40 @@
module psi_d_serial_mod
use psb_const_mod, only : psb_ipk_, psb_lpk_, psb_mpk_, psb_epk_, psb_dpk_
interface psb_gelp
interface psb_gelp
! 2-D version
subroutine psb_dgelp(trans,iperm,x,info)
import :: psb_ipk_, psb_dpk_
subroutine psb_m_dgelp(trans,iperm,x,info)
import :: psb_ipk_, psb_mpk_, psb_dpk_
implicit none
real(psb_dpk_), intent(inout) :: x(:,:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_dgelp
subroutine psb_dgelpv(trans,iperm,x,info)
import :: psb_ipk_, psb_dpk_
end subroutine psb_m_dgelp
subroutine psb_m_dgelpv(trans,iperm,x,info)
import :: psb_ipk_, psb_mpk_,psb_dpk_
implicit none
real(psb_dpk_), intent(inout) :: x(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_m_dgelpv
subroutine psb_e_dgelp(trans,iperm,x,info)
import :: psb_ipk_, psb_epk_, psb_dpk_
implicit none
real(psb_dpk_), intent(inout) :: x(:,:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_e_dgelp
subroutine psb_e_dgelpv(trans,iperm,x,info)
import :: psb_ipk_, psb_epk_, psb_dpk_
implicit none
real(psb_dpk_), intent(inout) :: x(:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_dgelpv
end subroutine psb_e_dgelpv
end interface psb_gelp
interface psb_geaxpby

30
base/modules/auxil/psi_e_serial_mod.f90
Unescape Escape View File

@ -32,24 +32,40 @@
module psi_e_serial_mod
use psb_const_mod, only : psb_ipk_, psb_lpk_, psb_mpk_, psb_epk_
interface psb_gelp
interface psb_gelp
! 2-D version
subroutine psb_egelp(trans,iperm,x,info)
subroutine psb_m_egelp(trans,iperm,x,info)
import :: psb_ipk_, psb_lpk_,psb_mpk_, psb_epk_
implicit none
integer(psb_epk_), intent(inout) :: x(:,:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_egelp
subroutine psb_egelpv(trans,iperm,x,info)
end subroutine psb_m_egelp
subroutine psb_m_egelpv(trans,iperm,x,info)
import :: psb_ipk_, psb_lpk_,psb_mpk_, psb_epk_
implicit none
integer(psb_epk_), intent(inout) :: x(:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_egelpv
end subroutine psb_m_egelpv
subroutine psb_e_egelp(trans,iperm,x,info)
import :: psb_ipk_, psb_lpk_,psb_mpk_, psb_epk_
implicit none
integer(psb_epk_), intent(inout) :: x(:,:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_e_egelp
subroutine psb_e_egelpv(trans,iperm,x,info)
import :: psb_ipk_, psb_lpk_,psb_mpk_, psb_epk_
implicit none
integer(psb_epk_), intent(inout) :: x(:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_e_egelpv
end interface psb_gelp
interface psb_geaxpby

30
base/modules/auxil/psi_i2_serial_mod.f90
Unescape Escape View File

@ -32,24 +32,40 @@
module psi_i2_serial_mod
use psb_const_mod, only : psb_ipk_, psb_lpk_, psb_mpk_, psb_epk_
interface psb_gelp
interface psb_gelp
! 2-D version
subroutine psb_i2gelp(trans,iperm,x,info)
subroutine psb_m_i2gelp(trans,iperm,x,info)
import :: psb_ipk_, psb_lpk_,psb_mpk_, psb_epk_
implicit none
integer(psb_i2pk_), intent(inout) :: x(:,:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_i2gelp
subroutine psb_i2gelpv(trans,iperm,x,info)
end subroutine psb_m_i2gelp
subroutine psb_m_i2gelpv(trans,iperm,x,info)
import :: psb_ipk_, psb_lpk_,psb_mpk_, psb_epk_
implicit none
integer(psb_i2pk_), intent(inout) :: x(:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_i2gelpv
end subroutine psb_m_i2gelpv
subroutine psb_e_i2gelp(trans,iperm,x,info)
import :: psb_ipk_, psb_lpk_,psb_mpk_, psb_epk_
implicit none
integer(psb_i2pk_), intent(inout) :: x(:,:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_e_i2gelp
subroutine psb_e_i2gelpv(trans,iperm,x,info)
import :: psb_ipk_, psb_lpk_,psb_mpk_, psb_epk_
implicit none
integer(psb_i2pk_), intent(inout) :: x(:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_e_i2gelpv
end interface psb_gelp
interface psb_geaxpby

30
base/modules/auxil/psi_m_serial_mod.f90
Unescape Escape View File

@ -32,24 +32,40 @@
module psi_m_serial_mod
use psb_const_mod, only : psb_ipk_, psb_lpk_, psb_mpk_, psb_epk_
interface psb_gelp
interface psb_gelp
! 2-D version
subroutine psb_mgelp(trans,iperm,x,info)
subroutine psb_m_mgelp(trans,iperm,x,info)
import :: psb_ipk_, psb_lpk_,psb_mpk_, psb_epk_
implicit none
integer(psb_mpk_), intent(inout) :: x(:,:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_mgelp
subroutine psb_mgelpv(trans,iperm,x,info)
end subroutine psb_m_mgelp
subroutine psb_m_mgelpv(trans,iperm,x,info)
import :: psb_ipk_, psb_lpk_,psb_mpk_, psb_epk_
implicit none
integer(psb_mpk_), intent(inout) :: x(:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_mgelpv
end subroutine psb_m_mgelpv
subroutine psb_e_mgelp(trans,iperm,x,info)
import :: psb_ipk_, psb_lpk_,psb_mpk_, psb_epk_
implicit none
integer(psb_mpk_), intent(inout) :: x(:,:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_e_mgelp
subroutine psb_e_mgelpv(trans,iperm,x,info)
import :: psb_ipk_, psb_lpk_,psb_mpk_, psb_epk_
implicit none
integer(psb_mpk_), intent(inout) :: x(:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_e_mgelpv
end interface psb_gelp
interface psb_geaxpby

34
base/modules/auxil/psi_s_serial_mod.f90
Unescape Escape View File

@ -32,24 +32,40 @@
module psi_s_serial_mod
use psb_const_mod, only : psb_ipk_, psb_lpk_, psb_mpk_, psb_epk_, psb_spk_
interface psb_gelp
interface psb_gelp
! 2-D version
subroutine psb_sgelp(trans,iperm,x,info)
import :: psb_ipk_, psb_spk_
subroutine psb_m_sgelp(trans,iperm,x,info)
import :: psb_ipk_, psb_mpk_, psb_spk_
implicit none
real(psb_spk_), intent(inout) :: x(:,:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_sgelp
subroutine psb_sgelpv(trans,iperm,x,info)
import :: psb_ipk_, psb_spk_
end subroutine psb_m_sgelp
subroutine psb_m_sgelpv(trans,iperm,x,info)
import :: psb_ipk_, psb_mpk_,psb_spk_
implicit none
real(psb_spk_), intent(inout) :: x(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_m_sgelpv
subroutine psb_e_sgelp(trans,iperm,x,info)
import :: psb_ipk_, psb_epk_, psb_spk_
implicit none
real(psb_spk_), intent(inout) :: x(:,:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_e_sgelp
subroutine psb_e_sgelpv(trans,iperm,x,info)
import :: psb_ipk_, psb_epk_, psb_spk_
implicit none
real(psb_spk_), intent(inout) :: x(:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_sgelpv
end subroutine psb_e_sgelpv
end interface psb_gelp
interface psb_geaxpby

34
base/modules/auxil/psi_z_serial_mod.f90
Unescape Escape View File

@ -32,24 +32,40 @@
module psi_z_serial_mod
use psb_const_mod, only : psb_ipk_, psb_lpk_, psb_mpk_, psb_epk_, psb_dpk_
interface psb_gelp
interface psb_gelp
! 2-D version
subroutine psb_zgelp(trans,iperm,x,info)
import :: psb_ipk_, psb_dpk_
subroutine psb_m_zgelp(trans,iperm,x,info)
import :: psb_ipk_, psb_mpk_, psb_dpk_
implicit none
complex(psb_dpk_), intent(inout) :: x(:,:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_zgelp
subroutine psb_zgelpv(trans,iperm,x,info)
import :: psb_ipk_, psb_dpk_
end subroutine psb_m_zgelp
subroutine psb_m_zgelpv(trans,iperm,x,info)
import :: psb_ipk_, psb_mpk_,psb_dpk_
implicit none
complex(psb_dpk_), intent(inout) :: x(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_m_zgelpv
subroutine psb_e_zgelp(trans,iperm,x,info)
import :: psb_ipk_, psb_epk_, psb_dpk_
implicit none
complex(psb_dpk_), intent(inout) :: x(:,:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_e_zgelp
subroutine psb_e_zgelpv(trans,iperm,x,info)
import :: psb_ipk_, psb_epk_, psb_dpk_
implicit none
complex(psb_dpk_), intent(inout) :: x(:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
end subroutine psb_zgelpv
end subroutine psb_e_zgelpv
end interface psb_gelp
interface psb_geaxpby

2
base/serial/Makefile
Unescape Escape View File

@ -11,10 +11,10 @@ FOBJS = psb_lsame.o psi_m_serial_impl.o psi_e_serial_impl.o \
smmp.o lsmmp.o \
psb_sgeprt.o psb_dgeprt.o psb_cgeprt.o psb_zgeprt.o\
psb_spdot_srtd.o psb_aspxpby.o psb_spge_dot.o\
psb_sgelp.o psb_dgelp.o psb_cgelp.o psb_zgelp.o \
psb_samax_s.o psb_damax_s.o psb_camax_s.o psb_zamax_s.o \
psb_sasum_s.o psb_dasum_s.o psb_casum_s.o psb_zasum_s.o
LIBDIR=..
INCDIR=..
MODDIR=../modules

246
base/serial/psb_cgelp.f90
Unescape Escape View File

@ -1,246 +0,0 @@
!
! Parallel Sparse BLAS version 3.5
! (C) Copyright 2006-2018
! Salvatore Filippone
! Alfredo Buttari University of Rome Tor Vergata
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
!
! File: psb_cgelp.f90
!
!
! Subroutine: psb_cgelp
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:,:).
! info - integer. Return code.
subroutine psb_cgelp(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_cgelp
use psb_const_mod
use psb_error_mod
implicit none
complex(psb_spk_), intent(inout) :: x(:,:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
complex(psb_spk_),allocatable :: temp(:)
integer(psb_ipk_) :: int_err(5), i1sz, i2sz, err_act,i,j
integer(psb_ipk_), allocatable :: itemp(:)
complex(psb_spk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name, ch_err
name = 'psb_cgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = size(x,dim=1)
i2sz = size(x,dim=2)
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz,i2sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do j=1,i2sz
do i=1,i1sz
temp(i) = x(itemp(i),j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case('T')
do j=1,i2sz
do i=1,i1sz
temp(itemp(i)) = x(i,j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case default
info=psb_err_from_subroutine_
ch_err='dgelp'
call psb_errpush(info,name,a_err=ch_err)
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_cgelp
!!$
!!$ Parallel Sparse BLAS version 3.5
!!$ (C) Copyright 2006-2018
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
!
!
! Subroutine: psb_cgelpv
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:).
! info - integer. Return code.
subroutine psb_cgelpv(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_cgelpv
use psb_const_mod
use psb_error_mod
implicit none
complex(psb_spk_), intent(inout) :: x(:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_ipk_) :: int_err(5), i1sz, err_act, i
complex(psb_spk_),allocatable :: temp(:)
integer(psb_ipk_), allocatable :: itemp(:)
complex(psb_spk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name, ch_err
name = 'psb_cgelpv'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = min(size(x),size(iperm))
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do i=1,i1sz
temp(i) = x(itemp(i))
end do
do i=1,i1sz
x(i) = temp(i)
end do
case('T')
do i=1,i1sz
temp(itemp(i)) = x(i)
end do
do i=1,i1sz
x(i) = temp(i)
end do
case default
info=psb_err_from_subroutine_
ch_err='dgelp'
call psb_errpush(info,name,a_err=ch_err)
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_cgelpv

246
base/serial/psb_dgelp.f90
Unescape Escape View File

@ -1,246 +0,0 @@
!
! Parallel Sparse BLAS version 3.5
! (C) Copyright 2006-2018
! Salvatore Filippone
! Alfredo Buttari University of Rome Tor Vergata
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
!
! File: psb_dgelp.f90
!
!
! Subroutine: psb_dgelp
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:,:).
! info - integer. Return code.
subroutine psb_dgelp(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_dgelp
use psb_const_mod
use psb_error_mod
implicit none
real(psb_dpk_), intent(inout) :: x(:,:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
real(psb_dpk_),allocatable :: temp(:)
integer(psb_ipk_) :: int_err(5), i1sz, i2sz, err_act,i,j
integer(psb_ipk_), allocatable :: itemp(:)
real(psb_dpk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name, ch_err
name = 'psb_dgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = size(x,dim=1)
i2sz = size(x,dim=2)
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz,i2sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do j=1,i2sz
do i=1,i1sz
temp(i) = x(itemp(i),j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case('T')
do j=1,i2sz
do i=1,i1sz
temp(itemp(i)) = x(i,j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case default
info=psb_err_from_subroutine_
ch_err='dgelp'
call psb_errpush(info,name,a_err=ch_err)
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_dgelp
!!$
!!$ Parallel Sparse BLAS version 3.5
!!$ (C) Copyright 2006-2018
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
!
!
! Subroutine: psb_dgelpv
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:).
! info - integer. Return code.
subroutine psb_dgelpv(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_dgelpv
use psb_const_mod
use psb_error_mod
implicit none
real(psb_dpk_), intent(inout) :: x(:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_ipk_) :: int_err(5), i1sz, err_act, i
real(psb_dpk_),allocatable :: temp(:)
integer(psb_ipk_), allocatable :: itemp(:)
real(psb_dpk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name, ch_err
name = 'psb_dgelpv'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = min(size(x),size(iperm))
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do i=1,i1sz
temp(i) = x(itemp(i))
end do
do i=1,i1sz
x(i) = temp(i)
end do
case('T')
do i=1,i1sz
temp(itemp(i)) = x(i)
end do
do i=1,i1sz
x(i) = temp(i)
end do
case default
info=psb_err_from_subroutine_
ch_err='dgelp'
call psb_errpush(info,name,a_err=ch_err)
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_dgelpv

246
base/serial/psb_sgelp.f90
Unescape Escape View File

@ -1,246 +0,0 @@
!
! Parallel Sparse BLAS version 3.5
! (C) Copyright 2006-2018
! Salvatore Filippone
! Alfredo Buttari University of Rome Tor Vergata
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
!
! File: psb_sgelp.f90
!
!
! Subroutine: psb_sgelp
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:,:).
! info - integer. Return code.
subroutine psb_sgelp(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_sgelp
use psb_const_mod
use psb_error_mod
implicit none
real(psb_spk_), intent(inout) :: x(:,:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
real(psb_spk_),allocatable :: temp(:)
integer(psb_ipk_) :: int_err(5), i1sz, i2sz, err_act,i,j
integer(psb_ipk_), allocatable :: itemp(:)
real(psb_spk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name, ch_err
name = 'psb_sgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = size(x,dim=1)
i2sz = size(x,dim=2)
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz,i2sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do j=1,i2sz
do i=1,i1sz
temp(i) = x(itemp(i),j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case('T')
do j=1,i2sz
do i=1,i1sz
temp(itemp(i)) = x(i,j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case default
info=psb_err_from_subroutine_
ch_err='dgelp'
call psb_errpush(info,name,a_err=ch_err)
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_sgelp
!!$
!!$ Parallel Sparse BLAS version 3.5
!!$ (C) Copyright 2006-2018
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
!
!
! Subroutine: psb_sgelpv
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:).
! info - integer. Return code.
subroutine psb_sgelpv(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_sgelpv
use psb_const_mod
use psb_error_mod
implicit none
real(psb_spk_), intent(inout) :: x(:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_ipk_) :: int_err(5), i1sz, err_act, i
real(psb_spk_),allocatable :: temp(:)
integer(psb_ipk_), allocatable :: itemp(:)
real(psb_spk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name, ch_err
name = 'psb_sgelpv'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = min(size(x),size(iperm))
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do i=1,i1sz
temp(i) = x(itemp(i))
end do
do i=1,i1sz
x(i) = temp(i)
end do
case('T')
do i=1,i1sz
temp(itemp(i)) = x(i)
end do
do i=1,i1sz
x(i) = temp(i)
end do
case default
info=psb_err_from_subroutine_
ch_err='dgelp'
call psb_errpush(info,name,a_err=ch_err)
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_sgelpv

246
base/serial/psb_zgelp.f90
Unescape Escape View File

@ -1,246 +0,0 @@
!
! Parallel Sparse BLAS version 3.5
! (C) Copyright 2006-2018
! Salvatore Filippone
! Alfredo Buttari University of Rome Tor Vergata
!
! Redistribution and use in source and binary forms, with or without
! modification, are permitted provided that the following conditions
! are met:
! 1. Redistributions of source code must retain the above copyright
! notice, this list of conditions and the following disclaimer.
! 2. Redistributions in binary form must reproduce the above copyright
! notice, this list of conditions, and the following disclaimer in the
! documentation and/or other materials provided with the distribution.
! 3. The name of the PSBLAS group or the names of its contributors may
! not be used to endorse or promote products derived from this
! software without specific written permission.
!
! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
! ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
! TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
! PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
! BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
! CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
! SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
! CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
! ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
! POSSIBILITY OF SUCH DAMAGE.
!
!
! File: psb_zgelp.f90
!
!
! Subroutine: psb_zgelp
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:,:).
! info - integer. Return code.
subroutine psb_zgelp(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_zgelp
use psb_const_mod
use psb_error_mod
implicit none
complex(psb_dpk_), intent(inout) :: x(:,:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
complex(psb_dpk_),allocatable :: temp(:)
integer(psb_ipk_) :: int_err(5), i1sz, i2sz, err_act,i,j
integer(psb_ipk_), allocatable :: itemp(:)
complex(psb_dpk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name, ch_err
name = 'psb_zgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = size(x,dim=1)
i2sz = size(x,dim=2)
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz,i2sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do j=1,i2sz
do i=1,i1sz
temp(i) = x(itemp(i),j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case('T')
do j=1,i2sz
do i=1,i1sz
temp(itemp(i)) = x(i,j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case default
info=psb_err_from_subroutine_
ch_err='dgelp'
call psb_errpush(info,name,a_err=ch_err)
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_zgelp
!!$
!!$ Parallel Sparse BLAS version 3.5
!!$ (C) Copyright 2006-2018
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
!
!
! Subroutine: psb_zgelpv
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:).
! info - integer. Return code.
subroutine psb_zgelpv(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_zgelpv
use psb_const_mod
use psb_error_mod
implicit none
complex(psb_dpk_), intent(inout) :: x(:)
integer(psb_ipk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_ipk_) :: int_err(5), i1sz, err_act, i
complex(psb_dpk_),allocatable :: temp(:)
integer(psb_ipk_), allocatable :: itemp(:)
complex(psb_dpk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name, ch_err
name = 'psb_zgelpv'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = min(size(x),size(iperm))
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do i=1,i1sz
temp(i) = x(itemp(i))
end do
do i=1,i1sz
x(i) = temp(i)
end do
case('T')
do i=1,i1sz
temp(itemp(i)) = x(i)
end do
do i=1,i1sz
x(i) = temp(i)
end do
case default
info=psb_err_from_subroutine_
ch_err='dgelp'
call psb_errpush(info,name,a_err=ch_err)
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_zgelpv

401
base/serial/psi_c_serial_impl.F90
Unescape Escape View File

@ -29,6 +29,407 @@
! POSSIBILITY OF SUCH DAMAGE.
!
!
subroutine psb_m_cgelp(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_m_cgelp
use psb_const_mod
use psb_error_mod
implicit none
complex(psb_spk_), intent(inout) :: x(:,:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
complex(psb_spk_),allocatable :: temp(:)
integer(psb_ipk_) :: int_err(5), i1sz, i2sz, err_act,i,j
integer(psb_ipk_), allocatable :: itemp(:)
complex(psb_spk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_cgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = size(x,dim=1)
i2sz = size(x,dim=2)
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz,i2sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do j=1,i2sz
do i=1,i1sz
temp(i) = x(itemp(i),j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case('T')
do j=1,i2sz
do i=1,i1sz
temp(itemp(i)) = x(i,j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='cgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_m_cgelp
!!$
!!$ Parallel Sparse BLAS version 3.5
!!$ (C) Copyright 2006-2018
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
!
!
! Subroutine: psb_cgelpv
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:).
! info - integer. Return code.
subroutine psb_m_cgelpv(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_m_cgelpv
use psb_const_mod
use psb_error_mod
implicit none
complex(psb_spk_), intent(inout) :: x(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_ipk_) :: int_err(5), i1sz, err_act, i
complex(psb_spk_),allocatable :: temp(:)
integer(psb_ipk_), allocatable :: itemp(:)
complex(psb_spk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_cgelpv'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = min(size(x),size(iperm))
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do i=1,i1sz
temp(i) = x(itemp(i))
end do
do i=1,i1sz
x(i) = temp(i)
end do
case('T')
do i=1,i1sz
temp(itemp(i)) = x(i)
end do
do i=1,i1sz
x(i) = temp(i)
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='cgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_m_cgelpv
subroutine psb_e_cgelp(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_e_cgelp
use psb_const_mod
use psb_error_mod
implicit none
complex(psb_spk_), intent(inout) :: x(:,:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
complex(psb_spk_),allocatable :: temp(:)
integer(psb_ipk_) :: int_err(5), err_act
integer(psb_epk_) :: i1sz, i2sz, i, j
integer(psb_epk_), allocatable :: itemp(:)
complex(psb_spk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_cgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = size(x,dim=1)
i2sz = size(x,dim=2)
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz,i2sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do j=1,i2sz
do i=1,i1sz
temp(i) = x(itemp(i),j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case('T')
do j=1,i2sz
do i=1,i1sz
temp(itemp(i)) = x(i,j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='cgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_e_cgelp
!!$
!!$ Parallel Sparse BLAS version 3.5
!!$ (C) Copyright 2006-2018
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
!
!
! Subroutine: psb_cgelpv
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:).
! info - integer. Return code.
subroutine psb_e_cgelpv(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_e_cgelpv
use psb_const_mod
use psb_error_mod
implicit none
complex(psb_spk_), intent(inout) :: x(:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_ipk_) :: int_err(5), err_act
complex(psb_spk_),allocatable :: temp(:)
integer(psb_epk_) :: i1sz, i
integer(psb_epk_), allocatable :: itemp(:)
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_cgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = min(size(x),size(iperm))
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do i=1,i1sz
temp(i) = x(itemp(i))
end do
do i=1,i1sz
x(i) = temp(i)
end do
case('T')
do i=1,i1sz
temp(itemp(i)) = x(i)
end do
do i=1,i1sz
x(i) = temp(i)
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='cgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_e_cgelpv
subroutine psi_caxpby(m,n,alpha, x, beta, y, info)
use psb_const_mod

401
base/serial/psi_d_serial_impl.F90
Unescape Escape View File

@ -29,6 +29,407 @@
! POSSIBILITY OF SUCH DAMAGE.
!
!
subroutine psb_m_dgelp(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_m_dgelp
use psb_const_mod
use psb_error_mod
implicit none
real(psb_dpk_), intent(inout) :: x(:,:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
real(psb_dpk_),allocatable :: temp(:)
integer(psb_ipk_) :: int_err(5), i1sz, i2sz, err_act,i,j
integer(psb_ipk_), allocatable :: itemp(:)
real(psb_dpk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_dgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = size(x,dim=1)
i2sz = size(x,dim=2)
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz,i2sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do j=1,i2sz
do i=1,i1sz
temp(i) = x(itemp(i),j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case('T')
do j=1,i2sz
do i=1,i1sz
temp(itemp(i)) = x(i,j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='dgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_m_dgelp
!!$
!!$ Parallel Sparse BLAS version 3.5
!!$ (C) Copyright 2006-2018
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
!
!
! Subroutine: psb_dgelpv
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:).
! info - integer. Return code.
subroutine psb_m_dgelpv(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_m_dgelpv
use psb_const_mod
use psb_error_mod
implicit none
real(psb_dpk_), intent(inout) :: x(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_ipk_) :: int_err(5), i1sz, err_act, i
real(psb_dpk_),allocatable :: temp(:)
integer(psb_ipk_), allocatable :: itemp(:)
real(psb_dpk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_dgelpv'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = min(size(x),size(iperm))
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do i=1,i1sz
temp(i) = x(itemp(i))
end do
do i=1,i1sz
x(i) = temp(i)
end do
case('T')
do i=1,i1sz
temp(itemp(i)) = x(i)
end do
do i=1,i1sz
x(i) = temp(i)
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='dgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_m_dgelpv
subroutine psb_e_dgelp(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_e_dgelp
use psb_const_mod
use psb_error_mod
implicit none
real(psb_dpk_), intent(inout) :: x(:,:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
real(psb_dpk_),allocatable :: temp(:)
integer(psb_ipk_) :: int_err(5), err_act
integer(psb_epk_) :: i1sz, i2sz, i, j
integer(psb_epk_), allocatable :: itemp(:)
real(psb_dpk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_dgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = size(x,dim=1)
i2sz = size(x,dim=2)
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz,i2sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do j=1,i2sz
do i=1,i1sz
temp(i) = x(itemp(i),j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case('T')
do j=1,i2sz
do i=1,i1sz
temp(itemp(i)) = x(i,j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='dgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_e_dgelp
!!$
!!$ Parallel Sparse BLAS version 3.5
!!$ (C) Copyright 2006-2018
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
!
!
! Subroutine: psb_dgelpv
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:).
! info - integer. Return code.
subroutine psb_e_dgelpv(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_e_dgelpv
use psb_const_mod
use psb_error_mod
implicit none
real(psb_dpk_), intent(inout) :: x(:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_ipk_) :: int_err(5), err_act
real(psb_dpk_),allocatable :: temp(:)
integer(psb_epk_) :: i1sz, i
integer(psb_epk_), allocatable :: itemp(:)
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_dgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = min(size(x),size(iperm))
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do i=1,i1sz
temp(i) = x(itemp(i))
end do
do i=1,i1sz
x(i) = temp(i)
end do
case('T')
do i=1,i1sz
temp(itemp(i)) = x(i)
end do
do i=1,i1sz
x(i) = temp(i)
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='dgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_e_dgelpv
subroutine psi_daxpby(m,n,alpha, x, beta, y, info)
use psb_const_mod

401
base/serial/psi_e_serial_impl.F90
Unescape Escape View File

@ -29,6 +29,407 @@
! POSSIBILITY OF SUCH DAMAGE.
!
!
subroutine psb_m_egelp(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_m_egelp
use psb_const_mod
use psb_error_mod
implicit none
integer(psb_epk_), intent(inout) :: x(:,:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_epk_),allocatable :: temp(:)
integer(psb_ipk_) :: int_err(5), i1sz, i2sz, err_act,i,j
integer(psb_ipk_), allocatable :: itemp(:)
integer(psb_epk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_egelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = size(x,dim=1)
i2sz = size(x,dim=2)
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz,i2sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do j=1,i2sz
do i=1,i1sz
temp(i) = x(itemp(i),j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case('T')
do j=1,i2sz
do i=1,i1sz
temp(itemp(i)) = x(i,j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='egelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_m_egelp
!!$
!!$ Parallel Sparse BLAS version 3.5
!!$ (C) Copyright 2006-2018
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
!
!
! Subroutine: psb_egelpv
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:).
! info - integer. Return code.
subroutine psb_m_egelpv(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_m_egelpv
use psb_const_mod
use psb_error_mod
implicit none
integer(psb_epk_), intent(inout) :: x(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_ipk_) :: int_err(5), i1sz, err_act, i
integer(psb_epk_),allocatable :: temp(:)
integer(psb_ipk_), allocatable :: itemp(:)
integer(psb_epk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_egelpv'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = min(size(x),size(iperm))
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do i=1,i1sz
temp(i) = x(itemp(i))
end do
do i=1,i1sz
x(i) = temp(i)
end do
case('T')
do i=1,i1sz
temp(itemp(i)) = x(i)
end do
do i=1,i1sz
x(i) = temp(i)
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='egelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_m_egelpv
subroutine psb_e_egelp(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_e_egelp
use psb_const_mod
use psb_error_mod
implicit none
integer(psb_epk_), intent(inout) :: x(:,:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_epk_),allocatable :: temp(:)
integer(psb_ipk_) :: int_err(5), err_act
integer(psb_epk_) :: i1sz, i2sz, i, j
integer(psb_epk_), allocatable :: itemp(:)
integer(psb_epk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_egelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = size(x,dim=1)
i2sz = size(x,dim=2)
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz,i2sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do j=1,i2sz
do i=1,i1sz
temp(i) = x(itemp(i),j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case('T')
do j=1,i2sz
do i=1,i1sz
temp(itemp(i)) = x(i,j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='egelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_e_egelp
!!$
!!$ Parallel Sparse BLAS version 3.5
!!$ (C) Copyright 2006-2018
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
!
!
! Subroutine: psb_egelpv
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:).
! info - integer. Return code.
subroutine psb_e_egelpv(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_e_egelpv
use psb_const_mod
use psb_error_mod
implicit none
integer(psb_epk_), intent(inout) :: x(:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_ipk_) :: int_err(5), err_act
integer(psb_epk_),allocatable :: temp(:)
integer(psb_epk_) :: i1sz, i
integer(psb_epk_), allocatable :: itemp(:)
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_egelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = min(size(x),size(iperm))
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do i=1,i1sz
temp(i) = x(itemp(i))
end do
do i=1,i1sz
x(i) = temp(i)
end do
case('T')
do i=1,i1sz
temp(itemp(i)) = x(i)
end do
do i=1,i1sz
x(i) = temp(i)
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='egelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_e_egelpv
subroutine psi_eaxpby(m,n,alpha, x, beta, y, info)
use psb_const_mod

401
base/serial/psi_i2_serial_impl.F90
Unescape Escape View File

@ -29,6 +29,407 @@
! POSSIBILITY OF SUCH DAMAGE.
!
!
subroutine psb_m_i2gelp(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_m_i2gelp
use psb_const_mod
use psb_error_mod
implicit none
integer(psb_i2pk_), intent(inout) :: x(:,:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_i2pk_),allocatable :: temp(:)
integer(psb_ipk_) :: int_err(5), i1sz, i2sz, err_act,i,j
integer(psb_ipk_), allocatable :: itemp(:)
integer(psb_i2pk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_i2gelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = size(x,dim=1)
i2sz = size(x,dim=2)
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz,i2sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do j=1,i2sz
do i=1,i1sz
temp(i) = x(itemp(i),j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case('T')
do j=1,i2sz
do i=1,i1sz
temp(itemp(i)) = x(i,j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='i2gelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_m_i2gelp
!!$
!!$ Parallel Sparse BLAS version 3.5
!!$ (C) Copyright 2006-2018
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
!
!
! Subroutine: psb_i2gelpv
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:).
! info - integer. Return code.
subroutine psb_m_i2gelpv(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_m_i2gelpv
use psb_const_mod
use psb_error_mod
implicit none
integer(psb_i2pk_), intent(inout) :: x(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_ipk_) :: int_err(5), i1sz, err_act, i
integer(psb_i2pk_),allocatable :: temp(:)
integer(psb_ipk_), allocatable :: itemp(:)
integer(psb_i2pk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_i2gelpv'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = min(size(x),size(iperm))
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do i=1,i1sz
temp(i) = x(itemp(i))
end do
do i=1,i1sz
x(i) = temp(i)
end do
case('T')
do i=1,i1sz
temp(itemp(i)) = x(i)
end do
do i=1,i1sz
x(i) = temp(i)
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='i2gelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_m_i2gelpv
subroutine psb_e_i2gelp(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_e_i2gelp
use psb_const_mod
use psb_error_mod
implicit none
integer(psb_i2pk_), intent(inout) :: x(:,:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_i2pk_),allocatable :: temp(:)
integer(psb_ipk_) :: int_err(5), err_act
integer(psb_epk_) :: i1sz, i2sz, i, j
integer(psb_epk_), allocatable :: itemp(:)
integer(psb_i2pk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_i2gelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = size(x,dim=1)
i2sz = size(x,dim=2)
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz,i2sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do j=1,i2sz
do i=1,i1sz
temp(i) = x(itemp(i),j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case('T')
do j=1,i2sz
do i=1,i1sz
temp(itemp(i)) = x(i,j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='i2gelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_e_i2gelp
!!$
!!$ Parallel Sparse BLAS version 3.5
!!$ (C) Copyright 2006-2018
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
!
!
! Subroutine: psb_i2gelpv
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:).
! info - integer. Return code.
subroutine psb_e_i2gelpv(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_e_i2gelpv
use psb_const_mod
use psb_error_mod
implicit none
integer(psb_i2pk_), intent(inout) :: x(:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_ipk_) :: int_err(5), err_act
integer(psb_i2pk_),allocatable :: temp(:)
integer(psb_epk_) :: i1sz, i
integer(psb_epk_), allocatable :: itemp(:)
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_i2gelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = min(size(x),size(iperm))
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do i=1,i1sz
temp(i) = x(itemp(i))
end do
do i=1,i1sz
x(i) = temp(i)
end do
case('T')
do i=1,i1sz
temp(itemp(i)) = x(i)
end do
do i=1,i1sz
x(i) = temp(i)
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='i2gelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_e_i2gelpv
subroutine psi_i2axpby(m,n,alpha, x, beta, y, info)
use psb_const_mod

401
base/serial/psi_m_serial_impl.F90
Unescape Escape View File

@ -29,6 +29,407 @@
! POSSIBILITY OF SUCH DAMAGE.
!
!
subroutine psb_m_mgelp(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_m_mgelp
use psb_const_mod
use psb_error_mod
implicit none
integer(psb_mpk_), intent(inout) :: x(:,:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_mpk_),allocatable :: temp(:)
integer(psb_ipk_) :: int_err(5), i1sz, i2sz, err_act,i,j
integer(psb_ipk_), allocatable :: itemp(:)
integer(psb_mpk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_mgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = size(x,dim=1)
i2sz = size(x,dim=2)
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz,i2sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do j=1,i2sz
do i=1,i1sz
temp(i) = x(itemp(i),j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case('T')
do j=1,i2sz
do i=1,i1sz
temp(itemp(i)) = x(i,j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='mgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_m_mgelp
!!$
!!$ Parallel Sparse BLAS version 3.5
!!$ (C) Copyright 2006-2018
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
!
!
! Subroutine: psb_mgelpv
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:).
! info - integer. Return code.
subroutine psb_m_mgelpv(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_m_mgelpv
use psb_const_mod
use psb_error_mod
implicit none
integer(psb_mpk_), intent(inout) :: x(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_ipk_) :: int_err(5), i1sz, err_act, i
integer(psb_mpk_),allocatable :: temp(:)
integer(psb_ipk_), allocatable :: itemp(:)
integer(psb_mpk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_mgelpv'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = min(size(x),size(iperm))
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do i=1,i1sz
temp(i) = x(itemp(i))
end do
do i=1,i1sz
x(i) = temp(i)
end do
case('T')
do i=1,i1sz
temp(itemp(i)) = x(i)
end do
do i=1,i1sz
x(i) = temp(i)
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='mgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_m_mgelpv
subroutine psb_e_mgelp(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_e_mgelp
use psb_const_mod
use psb_error_mod
implicit none
integer(psb_mpk_), intent(inout) :: x(:,:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_mpk_),allocatable :: temp(:)
integer(psb_ipk_) :: int_err(5), err_act
integer(psb_epk_) :: i1sz, i2sz, i, j
integer(psb_epk_), allocatable :: itemp(:)
integer(psb_mpk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_mgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = size(x,dim=1)
i2sz = size(x,dim=2)
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz,i2sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do j=1,i2sz
do i=1,i1sz
temp(i) = x(itemp(i),j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case('T')
do j=1,i2sz
do i=1,i1sz
temp(itemp(i)) = x(i,j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='mgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_e_mgelp
!!$
!!$ Parallel Sparse BLAS version 3.5
!!$ (C) Copyright 2006-2018
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
!
!
! Subroutine: psb_mgelpv
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:).
! info - integer. Return code.
subroutine psb_e_mgelpv(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_e_mgelpv
use psb_const_mod
use psb_error_mod
implicit none
integer(psb_mpk_), intent(inout) :: x(:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_ipk_) :: int_err(5), err_act
integer(psb_mpk_),allocatable :: temp(:)
integer(psb_epk_) :: i1sz, i
integer(psb_epk_), allocatable :: itemp(:)
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_mgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = min(size(x),size(iperm))
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do i=1,i1sz
temp(i) = x(itemp(i))
end do
do i=1,i1sz
x(i) = temp(i)
end do
case('T')
do i=1,i1sz
temp(itemp(i)) = x(i)
end do
do i=1,i1sz
x(i) = temp(i)
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='mgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_e_mgelpv
subroutine psi_maxpby(m,n,alpha, x, beta, y, info)
use psb_const_mod

401
base/serial/psi_s_serial_impl.F90
Unescape Escape View File

@ -29,6 +29,407 @@
! POSSIBILITY OF SUCH DAMAGE.
!
!
subroutine psb_m_sgelp(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_m_sgelp
use psb_const_mod
use psb_error_mod
implicit none
real(psb_spk_), intent(inout) :: x(:,:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
real(psb_spk_),allocatable :: temp(:)
integer(psb_ipk_) :: int_err(5), i1sz, i2sz, err_act,i,j
integer(psb_ipk_), allocatable :: itemp(:)
real(psb_spk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_sgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = size(x,dim=1)
i2sz = size(x,dim=2)
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz,i2sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do j=1,i2sz
do i=1,i1sz
temp(i) = x(itemp(i),j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case('T')
do j=1,i2sz
do i=1,i1sz
temp(itemp(i)) = x(i,j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='sgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_m_sgelp
!!$
!!$ Parallel Sparse BLAS version 3.5
!!$ (C) Copyright 2006-2018
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
!
!
! Subroutine: psb_sgelpv
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:).
! info - integer. Return code.
subroutine psb_m_sgelpv(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_m_sgelpv
use psb_const_mod
use psb_error_mod
implicit none
real(psb_spk_), intent(inout) :: x(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_ipk_) :: int_err(5), i1sz, err_act, i
real(psb_spk_),allocatable :: temp(:)
integer(psb_ipk_), allocatable :: itemp(:)
real(psb_spk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_sgelpv'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = min(size(x),size(iperm))
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do i=1,i1sz
temp(i) = x(itemp(i))
end do
do i=1,i1sz
x(i) = temp(i)
end do
case('T')
do i=1,i1sz
temp(itemp(i)) = x(i)
end do
do i=1,i1sz
x(i) = temp(i)
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='sgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_m_sgelpv
subroutine psb_e_sgelp(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_e_sgelp
use psb_const_mod
use psb_error_mod
implicit none
real(psb_spk_), intent(inout) :: x(:,:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
real(psb_spk_),allocatable :: temp(:)
integer(psb_ipk_) :: int_err(5), err_act
integer(psb_epk_) :: i1sz, i2sz, i, j
integer(psb_epk_), allocatable :: itemp(:)
real(psb_spk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_sgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = size(x,dim=1)
i2sz = size(x,dim=2)
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz,i2sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do j=1,i2sz
do i=1,i1sz
temp(i) = x(itemp(i),j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case('T')
do j=1,i2sz
do i=1,i1sz
temp(itemp(i)) = x(i,j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='sgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_e_sgelp
!!$
!!$ Parallel Sparse BLAS version 3.5
!!$ (C) Copyright 2006-2018
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
!
!
! Subroutine: psb_sgelpv
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:).
! info - integer. Return code.
subroutine psb_e_sgelpv(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_e_sgelpv
use psb_const_mod
use psb_error_mod
implicit none
real(psb_spk_), intent(inout) :: x(:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_ipk_) :: int_err(5), err_act
real(psb_spk_),allocatable :: temp(:)
integer(psb_epk_) :: i1sz, i
integer(psb_epk_), allocatable :: itemp(:)
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_sgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = min(size(x),size(iperm))
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do i=1,i1sz
temp(i) = x(itemp(i))
end do
do i=1,i1sz
x(i) = temp(i)
end do
case('T')
do i=1,i1sz
temp(itemp(i)) = x(i)
end do
do i=1,i1sz
x(i) = temp(i)
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='sgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_e_sgelpv
subroutine psi_saxpby(m,n,alpha, x, beta, y, info)
use psb_const_mod

401
base/serial/psi_z_serial_impl.F90
Unescape Escape View File

@ -29,6 +29,407 @@
! POSSIBILITY OF SUCH DAMAGE.
!
!
subroutine psb_m_zgelp(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_m_zgelp
use psb_const_mod
use psb_error_mod
implicit none
complex(psb_dpk_), intent(inout) :: x(:,:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
complex(psb_dpk_),allocatable :: temp(:)
integer(psb_ipk_) :: int_err(5), i1sz, i2sz, err_act,i,j
integer(psb_ipk_), allocatable :: itemp(:)
complex(psb_dpk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_zgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = size(x,dim=1)
i2sz = size(x,dim=2)
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz,i2sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do j=1,i2sz
do i=1,i1sz
temp(i) = x(itemp(i),j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case('T')
do j=1,i2sz
do i=1,i1sz
temp(itemp(i)) = x(i,j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='zgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_m_zgelp
!!$
!!$ Parallel Sparse BLAS version 3.5
!!$ (C) Copyright 2006-2018
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
!
!
! Subroutine: psb_zgelpv
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:).
! info - integer. Return code.
subroutine psb_m_zgelpv(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_m_zgelpv
use psb_const_mod
use psb_error_mod
implicit none
complex(psb_dpk_), intent(inout) :: x(:)
integer(psb_mpk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_ipk_) :: int_err(5), i1sz, err_act, i
complex(psb_dpk_),allocatable :: temp(:)
integer(psb_ipk_), allocatable :: itemp(:)
complex(psb_dpk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_zgelpv'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = min(size(x),size(iperm))
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do i=1,i1sz
temp(i) = x(itemp(i))
end do
do i=1,i1sz
x(i) = temp(i)
end do
case('T')
do i=1,i1sz
temp(itemp(i)) = x(i)
end do
do i=1,i1sz
x(i) = temp(i)
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='zgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_m_zgelpv
subroutine psb_e_zgelp(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_e_zgelp
use psb_const_mod
use psb_error_mod
implicit none
complex(psb_dpk_), intent(inout) :: x(:,:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
complex(psb_dpk_),allocatable :: temp(:)
integer(psb_ipk_) :: int_err(5), err_act
integer(psb_epk_) :: i1sz, i2sz, i, j
integer(psb_epk_), allocatable :: itemp(:)
complex(psb_dpk_),parameter :: one=1
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_zgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = size(x,dim=1)
i2sz = size(x,dim=2)
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz,i2sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do j=1,i2sz
do i=1,i1sz
temp(i) = x(itemp(i),j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case('T')
do j=1,i2sz
do i=1,i1sz
temp(itemp(i)) = x(i,j)
end do
do i=1,i1sz
x(i,j) = temp(i)
end do
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='zgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_e_zgelp
!!$
!!$ Parallel Sparse BLAS version 3.5
!!$ (C) Copyright 2006-2018
!!$ Salvatore Filippone University of Rome Tor Vergata
!!$ Alfredo Buttari University of Rome Tor Vergata
!!$
!!$ Redistribution and use in source and binary forms, with or without
!!$ modification, are permitted provided that the following conditions
!!$ are met:
!!$ 1. Redistributions of source code must retain the above copyright
!!$ notice, this list of conditions and the following disclaimer.
!!$ 2. Redistributions in binary form must reproduce the above copyright
!!$ notice, this list of conditions, and the following disclaimer in the
!!$ documentation and/or other materials provided with the distribution.
!!$ 3. The name of the PSBLAS group or the names of its contributors may
!!$ not be used to endorse or promote products derived from this
!!$ software without specific written permission.
!!$
!!$ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
!!$ ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
!!$ TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
!!$ PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE PSBLAS GROUP OR ITS CONTRIBUTORS
!!$ BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
!!$ CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
!!$ SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
!!$ INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!!$ CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
!!$ ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
!!$ POSSIBILITY OF SUCH DAMAGE.
!!$
!!$
!
!
! Subroutine: psb_zgelpv
! Apply a left permutation to a dense matrix
!
! Arguments:
! trans - character.
! iperm - integer.
! x - real, dimension(:).
! info - integer. Return code.
subroutine psb_e_zgelpv(trans,iperm,x,info)
use psb_serial_mod, psb_protect_name => psb_e_zgelpv
use psb_const_mod
use psb_error_mod
implicit none
complex(psb_dpk_), intent(inout) :: x(:)
integer(psb_epk_), intent(in) :: iperm(:)
integer(psb_ipk_), intent(out) :: info
character, intent(in) :: trans
! local variables
integer(psb_ipk_) :: int_err(5), err_act
complex(psb_dpk_),allocatable :: temp(:)
integer(psb_epk_) :: i1sz, i
integer(psb_epk_), allocatable :: itemp(:)
integer(psb_ipk_) :: debug_level, debug_unit
character(len=20) :: name
name = 'psb_zgelp'
if(psb_get_errstatus() /= 0) return
info=psb_success_
call psb_erractionsave(err_act)
debug_unit = psb_get_debug_unit()
debug_level = psb_get_debug_level()
i1sz = min(size(x),size(iperm))
if (debug_level >= psb_debug_serial_)&
& write(debug_unit,*) trim(name),': size',i1sz
allocate(temp(i1sz),itemp(size(iperm)),stat=info)
if (info /= psb_success_) then
info=2040
call psb_errpush(info,name)
goto 9999
end if
itemp(:) = iperm(:)
if (.not.psb_isaperm(i1sz,itemp)) then
info=psb_err_iarg_invalid_value_
int_err(1) = 1
call psb_errpush(info,name,i_err=int_err)
goto 9999
endif
select case( psb_toupper(trans))
case('N')
do i=1,i1sz
temp(i) = x(itemp(i))
end do
do i=1,i1sz
x(i) = temp(i)
end do
case('T')
do i=1,i1sz
temp(itemp(i)) = x(i)
end do
do i=1,i1sz
x(i) = temp(i)
end do
case default
info=psb_err_from_subroutine_
call psb_errpush(info,name,a_err='zgelp')
end select
deallocate(temp,itemp)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_e_zgelpv
subroutine psi_zaxpby(m,n,alpha, x, beta, y, info)
use psb_const_mod

4
docs/src/datastruct.tex
Unescape Escape View File

@ -64,7 +64,7 @@ First we describe the \verb|psb_indx_map| type. This is a data
structure that keeps track of a certain number of basic issues such
as:
\begin{itemize}
\item The value of the communication/MPI context;
\item The value of the communication context;
\item The number of indices in the index space, i.e. global number of
rows and columns of a sparse matrix;
\item The local set of indices, including:
@ -309,7 +309,7 @@ Type: {\bf optional}; default: \verb|.true.|.\\
\subsubsection{get\_context --- Get communication context}
\begin{verbatim}
ictxt = desc%get_context()
ctxt = desc%get_context()
\end{verbatim}
\begin{description}

74
docs/src/penv.tex
Unescape Escape View File

@ -7,7 +7,7 @@
environment}
\begin{verbatim}
call psb_init(icontxt, np, basectxt, ids)
call psb_init(ctxt, np, basectxt, ids)
\end{verbatim}
This subroutine initializes the PSBLAS parallel environment, defining
@ -41,8 +41,8 @@ Default: use the indices $(0\dots np-1)$.
\begin{description}
\item[\bf On Return]
\item[icontxt] the communication context identifying the virtual
parallel machine. Note that this is always a duplicate of
\item[ctxt] the communication context identifying the virtual
parallel machine, type \verb|psb_ctxt_type|. Note that this is always a duplicate of
\verb|basectxt|, so that library communications are completely
separated from other communication operations.\\
Scope: {\bf global}.\\
@ -65,7 +65,7 @@ Specified as: an integer variable.
environment}
\begin{verbatim}
call psb_info(icontxt, iam, np)
call psb_info(ctxt, iam, np)
\end{verbatim}
This subroutine returns information about the PSBLAS parallel environment, defining
@ -73,7 +73,7 @@ a virtual parallel machine.
\begin{description}
\item[Type:] Asynchronous.
\item[\bf On Entry ]
\item[icontxt] the communication context identifying the virtual
\item[ctxt] the communication context identifying the virtual
parallel machine.\\
Scope: {\bf global}.\\
Type: {\bf required}.\\
@ -103,7 +103,7 @@ Specified as: an integer variable. \
\item If the user has requested on \verb|psb_init| a number of
processes less than the total available in the parallel execution
environment, the remaining processes will have on return $iam=-1$;
the only call involving \verb|icontxt| that any such process may
the only call involving \verb|ctxt| that any such process may
execute is to \verb|psb_exit|.
\end{enumerate}
@ -112,22 +112,22 @@ Specified as: an integer variable. \
environment}
\begin{verbatim}
call psb_exit(icontxt)
call psb_exit(icontxt,close)
call psb_exit(ctxt)
call psb_exit(ctxt,close)
\end{verbatim}
This subroutine exits from the PSBLAS parallel virtual machine.
\begin{description}
\item[Type:] Synchronous.
\item[\bf On Entry ]
\item[icontxt] the communication context identifying the virtual
\item[ctxt] the communication context identifying the virtual
parallel machine.\\
Scope: {\bf global}.\\
Type: {\bf required}.\\
Intent: {\bf in}.\\
Specified as: an integer variable.
\item[close] Whether to close all data structures related to the
virtual parallel machine, besides those associated with icontxt.\\
virtual parallel machine, besides those associated with ctxt.\\
Scope: {\bf global}.\\
Type: {\bf optional}.\\
Intent: {\bf in}.\\
@ -138,7 +138,7 @@ Specified as: a logical variable, default value: true.
\begin{enumerate}
\item This routine may be called even if a previous call to
\verb|psb_info| has returned with $iam=-1$; indeed, it it is the only
routine that may be called with argument \verb|icontxt| in this
routine that may be called with argument \verb|ctxt| in this
situation.
\item A call to this routine with \verb|close=.true.| implies a call
to \verb|MPI_Finalize|, after which no parallel routine may be called.
@ -154,14 +154,14 @@ Specified as: a logical variable, default value: true.
\clearpage\subsection{psb\_get\_mpi\_comm --- Get the MPI communicator}
\begin{verbatim}
icomm = psb_get_mpi_comm(icontxt)
icomm = psb_get_mpi_comm(ctxt)
\end{verbatim}
This function returns the MPI communicator associated with a PSBLAS context
\begin{description}
\item[Type:] Asynchronous.
\item[\bf On Entry ]
\item[icontxt] the communication context identifying the virtual
\item[ctxt] the communication context identifying the virtual
parallel machine.\\
Scope: {\bf global}.\\
Type: {\bf required}.\\
@ -184,14 +184,14 @@ is deprecated.
\clearpage\subsection{psb\_get\_mpi\_rank --- Get the MPI rank}
\begin{verbatim}
rank = psb_get_mpi_rank(icontxt, id)
rank = psb_get_mpi_rank(ctxt, id)
\end{verbatim}
This function returns the MPI rank of the PSBLAS process $id$
\begin{description}
\item[Type:] Asynchronous.
\item[\bf On Entry ]
\item[icontxt] the communication context identifying the virtual
\item[ctxt] the communication context identifying the virtual
parallel machine.\\
Scope: {\bf global}.\\
Type: {\bf required}.\\
@ -238,7 +238,7 @@ Returned as: a \verb|real(psb_dpk_)| variable.
environment}
\begin{verbatim}
call psb_barrier(icontxt)
call psb_barrier(ctxt)
\end{verbatim}
This subroutine acts as an explicit synchronization point for the PSBLAS
@ -246,7 +246,7 @@ parallel virtual machine.
\begin{description}
\item[Type:] Synchronous.
\item[\bf On Entry ]
\item[icontxt] the communication context identifying the virtual
\item[ctxt] the communication context identifying the virtual
parallel machine.\\
Scope: {\bf global}.\\
Type: {\bf required}.\\
@ -258,14 +258,14 @@ Specified as: an integer variable.
\clearpage\subsection{psb\_abort --- Abort a computation}
\begin{verbatim}
call psb_abort(icontxt)
call psb_abort(ctxt)
\end{verbatim}
This subroutine aborts computation on the parallel virtual machine.
\begin{description}
\item[Type:] Asynchronous.
\item[\bf On Entry ]
\item[icontxt] the communication context identifying the virtual
\item[ctxt] the communication context identifying the virtual
parallel machine.\\
Scope: {\bf global}.\\
Type: {\bf required}.\\
@ -280,7 +280,7 @@ Specified as: an integer variable.
\clearpage\subsection{psb\_bcast --- Broadcast data}
\begin{verbatim}
call psb_bcast(icontxt, dat, root)
call psb_bcast(ctxt, dat, root)
\end{verbatim}
This subroutine implements a broadcast operation based on the
@ -288,7 +288,7 @@ underlying communication library.
\begin{description}
\item[Type:] Synchronous.
\item[\bf On Entry ]
\item[icontxt] the communication context identifying the virtual
\item[ctxt] the communication context identifying the virtual
parallel machine.\\
Scope: {\bf global}.\\
Type: {\bf required}.\\
@ -325,7 +325,7 @@ Type, kind, rank and size must agree on all processes.
\clearpage\subsection{psb\_sum --- Global sum}
\begin{verbatim}
call psb_sum(icontxt, dat, root)
call psb_sum(ctxt, dat, root)
\end{verbatim}
This subroutine implements a sum reduction operation based on the
@ -333,7 +333,7 @@ underlying communication library.
\begin{description}
\item[Type:] Synchronous.
\item[\bf On Entry ]
\item[icontxt] the communication context identifying the virtual
\item[ctxt] the communication context identifying the virtual
parallel machine.\\
Scope: {\bf global}.\\
Type: {\bf required}.\\
@ -379,7 +379,7 @@ Type, kind, rank and size must agree on all processes.
\clearpage\subsection{psb\_max --- Global maximum}
\begin{verbatim}
call psb_max(icontxt, dat, root)
call psb_max(ctxt, dat, root)
\end{verbatim}
This subroutine implements a maximum valuereduction
@ -387,7 +387,7 @@ operation based on the underlying communication library.
\begin{description}
\item[Type:] Synchronous.
\item[\bf On Entry ]
\item[icontxt] the communication context identifying the virtual
\item[ctxt] the communication context identifying the virtual
parallel machine.\\
Scope: {\bf global}.\\
Type: {\bf required}.\\
@ -432,7 +432,7 @@ Type, kind, rank and size must agree on all processes.
\clearpage\subsection{psb\_min --- Global minimum}
\begin{verbatim}
call psb_min(icontxt, dat, root)
call psb_min(ctxt, dat, root)
\end{verbatim}
This subroutine implements a minimum value reduction
@ -440,7 +440,7 @@ operation based on the underlying communication library.
\begin{description}
\item[Type:] Synchronous.
\item[\bf On Entry ]
\item[icontxt] the communication context identifying the virtual
\item[ctxt] the communication context identifying the virtual
parallel machine.\\
Scope: {\bf global}.\\
Type: {\bf required}.\\
@ -485,7 +485,7 @@ Type, kind, rank and size must agree on all processes.
\clearpage\subsection{psb\_amx --- Global maximum absolute value}
\begin{verbatim}
call psb_amx(icontxt, dat, root)
call psb_amx(ctxt, dat, root)
\end{verbatim}
This subroutine implements a maximum absolute value reduction
@ -493,7 +493,7 @@ operation based on the underlying communication library.
\begin{description}
\item[Type:] Synchronous.
\item[\bf On Entry ]
\item[icontxt] the communication context identifying the virtual
\item[ctxt] the communication context identifying the virtual
parallel machine.\\
Scope: {\bf global}.\\
Type: {\bf required}.\\
@ -538,7 +538,7 @@ Type, kind, rank and size must agree on all processes.
\clearpage\subsection{psb\_amn --- Global minimum absolute value}
\begin{verbatim}
call psb_amn(icontxt, dat, root)
call psb_amn(ctxt, dat, root)
\end{verbatim}
This subroutine implements a minimum absolute value reduction
@ -546,7 +546,7 @@ operation based on the underlying communication library.
\begin{description}
\item[Type:] Synchronous.
\item[\bf On Entry ]
\item[icontxt] the communication context identifying the virtual
\item[ctxt] the communication context identifying the virtual
parallel machine.\\
Scope: {\bf global}.\\
Type: {\bf required}.\\
@ -591,7 +591,7 @@ Type, kind, rank and size must agree on all processes.
\clearpage\subsection{psb\_nrm2 --- Global 2-norm reduction}
\begin{verbatim}
call psb_nrm2(icontxt, dat, root)
call psb_nrm2(ctxt, dat, root)
\end{verbatim}
This subroutine implements a 2-norm value reduction
@ -599,7 +599,7 @@ operation based on the underlying communication library.
\begin{description}
\item[Type:] Synchronous.
\item[\bf On Entry ]
\item[icontxt] the communication context identifying the virtual
\item[ctxt] the communication context identifying the virtual
parallel machine.\\
Scope: {\bf global}.\\
Type: {\bf required}.\\
@ -651,14 +651,14 @@ Kind, rank and size must agree on all processes.
\clearpage\subsection{psb\_snd --- Send data}
\begin{verbatim}
call psb_snd(icontxt, dat, dst, m)
call psb_snd(ctxt, dat, dst, m)
\end{verbatim}
This subroutine sends a packet of data to a destination.
\begin{description}
\item[Type:] Synchronous: see usage notes.
\item[\bf On Entry ]
\item[icontxt] the communication context identifying the virtual
\item[ctxt] the communication context identifying the virtual
parallel machine.\\
Scope: {\bf global}.\\
Type: {\bf required}.\\
@ -702,14 +702,14 @@ same value on sending and receiving processes.
\clearpage\subsection{psb\_rcv --- Receive data}
\begin{verbatim}
call psb_rcv(icontxt, dat, src, m)
call psb_rcv(ctxt, dat, src, m)
\end{verbatim}
This subroutine receives a packet of data to a destination.
\begin{description}
\item[Type:] Synchronous: see usage notes.
\item[\bf On Entry ]
\item[icontxt] the communication context identifying the virtual
\item[ctxt] the communication context identifying the virtual
parallel machine.\\
Scope: {\bf global}.\\
Type: {\bf required}.\\

8
docs/src/psbrout.tex
Unescape Escape View File

@ -215,7 +215,7 @@ An integer value; 0 means no error has been detected.
vres(1) = psb_gedot(x1,y1,desc_a,info,global=.false.)
vres(2) = psb_gedot(x2,y2,desc_a,info,global=.false.)
vres(3) = psb_gedot(x3,y3,desc_a,info,global=.false.)
call psb_sum(ictxt,vres(1:3))
call psb_sum(ctxt,vres(1:3))
\end{lstlisting}
In this way the global communication, which for small sizes is a
latency-bound operation, is invoked only once.
@ -391,7 +391,7 @@ An integer value; 0 means no error has been detected.
vres(1) = psb_geamax(x1,desc_a,info,global=.false.)
vres(2) = psb_geamax(x2,desc_a,info,global=.false.)
vres(3) = psb_geamax(x3,desc_a,info,global=.false.)
call psb_amx(ictxt,vres(1:3))
call psb_amx(ctxt,vres(1:3))
\end{lstlisting}
In this way the global communication, which for small sizes is a
latency-bound operation, is invoked only once.
@ -544,7 +544,7 @@ An integer value; 0 means no error has been detected.
vres(1) = psb_geasum(x1,desc_a,info,global=.false.)
vres(2) = psb_geasum(x2,desc_a,info,global=.false.)
vres(3) = psb_geasum(x3,desc_a,info,global=.false.)
call psb_sum(ictxt,vres(1:3))
call psb_sum(ctxt,vres(1:3))
\end{lstlisting}
In this way the global communication, which for small sizes is a
latency-bound operation, is invoked only once.
@ -714,7 +714,7 @@ An integer value; 0 means no error has been detected.
vres(1) = psb_genrm2(x1,desc_a,info,global=.false.)
vres(2) = psb_genrm2(x2,desc_a,info,global=.false.)
vres(3) = psb_genrm2(x3,desc_a,info,global=.false.)
call psb_nrm2(ictxt,vres(1:3))
call psb_nrm2(ctxt,vres(1:3))
\end{lstlisting}
In this way the global communication, which for small sizes is a
latency-bound operation, is invoked only once.

78
prec/impl/psb_c_bjacprec_impl.f90
Unescape Escape View File

@ -566,38 +566,52 @@ subroutine psb_c_bjac_precbld(a,desc_a,prec,info,amold,vmold,imold)
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
! Check if the variant for the AINV is known to the library
if( (prec%iprcparm(psb_f_type_) == psb_f_ainv_).and.(&
& (iinvalg == psb_ainv_llk_).or.(iinvalg == psb_ainv_s_llk_).or. &
& (iinvalg == psb_ainv_s_ft_llk_).or.(iinvalg == psb_ainv_llk_noth_).or.&
& (iinvalg == psb_ainv_mlk_).or.(iinvalg == psb_ainv_lmx_ ) ) ) then
! Do nothing, these are okay
else
info=psb_err_from_subroutine_
ch_err='psb_ainv_alg_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
! Check if the drop-tolerance make sense, if fact_eps > 1 and we are requiring
! either ILUT, or INVT we give an error.
if( (fact_eps > 1).and.( &
& (prec%iprcparm(psb_f_type_) == psb_f_ilu_t_).or.&
& (prec%iprcparm(psb_f_type_) == psb_f_invt_) )) then
info=psb_err_from_subroutine_
ch_err='psb_fact_eps_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
! It the drop-tolerance for the inverse factors, inv_thresh > 1 and we are
! requiring AINV or, or INVT we give an error
if( (inv_thresh > 1).and.( &
& (prec%iprcparm(psb_f_type_) == psb_f_ainv_).or.&
& (prec%iprcparm(psb_f_type_) == psb_f_invt_) )) then
info=psb_err_from_subroutine_
ch_err='psb_inv_thresh_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
select case (prec%iprcparm(psb_f_type_))
case (psb_f_ainv_)
! Check if the variant for the AINV is known to the library
select case (iinvalg)
case(psb_ainv_llk_,psb_ainv_s_llk_,psb_ainv_s_ft_llk_,psb_ainv_llk_noth_,&
& psb_ainv_mlk_)
! Do nothing these are okay
case default
info=psb_err_from_subroutine_
ch_err='psb_ainv_alg_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end select ! AINV Variant
! Check if the drop-tolerance make sense
if( inv_thresh > 1) then
info=psb_err_from_subroutine_
ch_err='psb_inv_thresh_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
case (psb_f_ilu_t_)
if (fact_eps > 1) then
! Check if the drop-tolerance make sense
info=psb_err_from_subroutine_
ch_err='psb_fact_eps_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
case (psb_f_invt_)
! Check both tolerances
if (fact_eps > 1) then
info=psb_err_from_subroutine_
ch_err='psb_fact_eps_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if( inv_thresh > 1) then
info=psb_err_from_subroutine_
ch_err='psb_inv_thresh_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
case default
end select
! Checks relative to the fill-in parameters
if (prec%iprcparm(psb_f_type_) == psb_f_ilu_n_) then
if(fill_in < 0) then

78
prec/impl/psb_d_bjacprec_impl.f90
Unescape Escape View File

@ -566,38 +566,52 @@ subroutine psb_d_bjac_precbld(a,desc_a,prec,info,amold,vmold,imold)
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
! Check if the variant for the AINV is known to the library
if( (prec%iprcparm(psb_f_type_) == psb_f_ainv_).and.(&
& (iinvalg == psb_ainv_llk_).or.(iinvalg == psb_ainv_s_llk_).or. &
& (iinvalg == psb_ainv_s_ft_llk_).or.(iinvalg == psb_ainv_llk_noth_).or.&
& (iinvalg == psb_ainv_mlk_).or.(iinvalg == psb_ainv_lmx_ ) ) ) then
! Do nothing, these are okay
else
info=psb_err_from_subroutine_
ch_err='psb_ainv_alg_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
! Check if the drop-tolerance make sense, if fact_eps > 1 and we are requiring
! either ILUT, or INVT we give an error.
if( (fact_eps > 1).and.( &
& (prec%iprcparm(psb_f_type_) == psb_f_ilu_t_).or.&
& (prec%iprcparm(psb_f_type_) == psb_f_invt_) )) then
info=psb_err_from_subroutine_
ch_err='psb_fact_eps_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
! It the drop-tolerance for the inverse factors, inv_thresh > 1 and we are
! requiring AINV or, or INVT we give an error
if( (inv_thresh > 1).and.( &
& (prec%iprcparm(psb_f_type_) == psb_f_ainv_).or.&
& (prec%iprcparm(psb_f_type_) == psb_f_invt_) )) then
info=psb_err_from_subroutine_
ch_err='psb_inv_thresh_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
select case (prec%iprcparm(psb_f_type_))
case (psb_f_ainv_)
! Check if the variant for the AINV is known to the library
select case (iinvalg)
case(psb_ainv_llk_,psb_ainv_s_llk_,psb_ainv_s_ft_llk_,psb_ainv_llk_noth_,&
& psb_ainv_mlk_)
! Do nothing these are okay
case default
info=psb_err_from_subroutine_
ch_err='psb_ainv_alg_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end select ! AINV Variant
! Check if the drop-tolerance make sense
if( inv_thresh > 1) then
info=psb_err_from_subroutine_
ch_err='psb_inv_thresh_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
case (psb_f_ilu_t_)
if (fact_eps > 1) then
! Check if the drop-tolerance make sense
info=psb_err_from_subroutine_
ch_err='psb_fact_eps_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
case (psb_f_invt_)
! Check both tolerances
if (fact_eps > 1) then
info=psb_err_from_subroutine_
ch_err='psb_fact_eps_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if( inv_thresh > 1) then
info=psb_err_from_subroutine_
ch_err='psb_inv_thresh_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
case default
end select
! Checks relative to the fill-in parameters
if (prec%iprcparm(psb_f_type_) == psb_f_ilu_n_) then
if(fill_in < 0) then

78
prec/impl/psb_s_bjacprec_impl.f90
Unescape Escape View File

@ -566,38 +566,52 @@ subroutine psb_s_bjac_precbld(a,desc_a,prec,info,amold,vmold,imold)
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
! Check if the variant for the AINV is known to the library
if( (prec%iprcparm(psb_f_type_) == psb_f_ainv_).and.(&
& (iinvalg == psb_ainv_llk_).or.(iinvalg == psb_ainv_s_llk_).or. &
& (iinvalg == psb_ainv_s_ft_llk_).or.(iinvalg == psb_ainv_llk_noth_).or.&
& (iinvalg == psb_ainv_mlk_).or.(iinvalg == psb_ainv_lmx_ ) ) ) then
! Do nothing, these are okay
else
info=psb_err_from_subroutine_
ch_err='psb_ainv_alg_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
! Check if the drop-tolerance make sense, if fact_eps > 1 and we are requiring
! either ILUT, or INVT we give an error.
if( (fact_eps > 1).and.( &
& (prec%iprcparm(psb_f_type_) == psb_f_ilu_t_).or.&
& (prec%iprcparm(psb_f_type_) == psb_f_invt_) )) then
info=psb_err_from_subroutine_
ch_err='psb_fact_eps_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
! It the drop-tolerance for the inverse factors, inv_thresh > 1 and we are
! requiring AINV or, or INVT we give an error
if( (inv_thresh > 1).and.( &
& (prec%iprcparm(psb_f_type_) == psb_f_ainv_).or.&
& (prec%iprcparm(psb_f_type_) == psb_f_invt_) )) then
info=psb_err_from_subroutine_
ch_err='psb_inv_thresh_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
select case (prec%iprcparm(psb_f_type_))
case (psb_f_ainv_)
! Check if the variant for the AINV is known to the library
select case (iinvalg)
case(psb_ainv_llk_,psb_ainv_s_llk_,psb_ainv_s_ft_llk_,psb_ainv_llk_noth_,&
& psb_ainv_mlk_)
! Do nothing these are okay
case default
info=psb_err_from_subroutine_
ch_err='psb_ainv_alg_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end select ! AINV Variant
! Check if the drop-tolerance make sense
if( inv_thresh > 1) then
info=psb_err_from_subroutine_
ch_err='psb_inv_thresh_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
case (psb_f_ilu_t_)
if (fact_eps > 1) then
! Check if the drop-tolerance make sense
info=psb_err_from_subroutine_
ch_err='psb_fact_eps_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
case (psb_f_invt_)
! Check both tolerances
if (fact_eps > 1) then
info=psb_err_from_subroutine_
ch_err='psb_fact_eps_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if( inv_thresh > 1) then
info=psb_err_from_subroutine_
ch_err='psb_inv_thresh_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
case default
end select
! Checks relative to the fill-in parameters
if (prec%iprcparm(psb_f_type_) == psb_f_ilu_n_) then
if(fill_in < 0) then

78
prec/impl/psb_z_bjacprec_impl.f90
Unescape Escape View File

@ -566,38 +566,52 @@ subroutine psb_z_bjac_precbld(a,desc_a,prec,info,amold,vmold,imold)
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
! Check if the variant for the AINV is known to the library
if( (prec%iprcparm(psb_f_type_) == psb_f_ainv_).and.(&
& (iinvalg == psb_ainv_llk_).or.(iinvalg == psb_ainv_s_llk_).or. &
& (iinvalg == psb_ainv_s_ft_llk_).or.(iinvalg == psb_ainv_llk_noth_).or.&
& (iinvalg == psb_ainv_mlk_).or.(iinvalg == psb_ainv_lmx_ ) ) ) then
! Do nothing, these are okay
else
info=psb_err_from_subroutine_
ch_err='psb_ainv_alg_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
! Check if the drop-tolerance make sense, if fact_eps > 1 and we are requiring
! either ILUT, or INVT we give an error.
if( (fact_eps > 1).and.( &
& (prec%iprcparm(psb_f_type_) == psb_f_ilu_t_).or.&
& (prec%iprcparm(psb_f_type_) == psb_f_invt_) )) then
info=psb_err_from_subroutine_
ch_err='psb_fact_eps_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
! It the drop-tolerance for the inverse factors, inv_thresh > 1 and we are
! requiring AINV or, or INVT we give an error
if( (inv_thresh > 1).and.( &
& (prec%iprcparm(psb_f_type_) == psb_f_ainv_).or.&
& (prec%iprcparm(psb_f_type_) == psb_f_invt_) )) then
info=psb_err_from_subroutine_
ch_err='psb_inv_thresh_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
select case (prec%iprcparm(psb_f_type_))
case (psb_f_ainv_)
! Check if the variant for the AINV is known to the library
select case (iinvalg)
case(psb_ainv_llk_,psb_ainv_s_llk_,psb_ainv_s_ft_llk_,psb_ainv_llk_noth_,&
& psb_ainv_mlk_)
! Do nothing these are okay
case default
info=psb_err_from_subroutine_
ch_err='psb_ainv_alg_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end select ! AINV Variant
! Check if the drop-tolerance make sense
if( inv_thresh > 1) then
info=psb_err_from_subroutine_
ch_err='psb_inv_thresh_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
case (psb_f_ilu_t_)
if (fact_eps > 1) then
! Check if the drop-tolerance make sense
info=psb_err_from_subroutine_
ch_err='psb_fact_eps_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
case (psb_f_invt_)
! Check both tolerances
if (fact_eps > 1) then
info=psb_err_from_subroutine_
ch_err='psb_fact_eps_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
if( inv_thresh > 1) then
info=psb_err_from_subroutine_
ch_err='psb_inv_thresh_'
call psb_errpush(info,name,a_err=ch_err)
goto 9999
end if
case default
end select
! Checks relative to the fill-in parameters
if (prec%iprcparm(psb_f_type_) == psb_f_ilu_n_) then
if(fill_in < 0) then

4
util/Makefile
Unescape Escape View File

@ -10,7 +10,8 @@ HERE=.
BASEOBJS= psb_blockpart_mod.o psb_metispart_mod.o psb_partidx_mod.o \
psb_hbio_mod.o psb_mmio_mod.o psb_mat_dist_mod.o \
psb_s_mat_dist_mod.o psb_d_mat_dist_mod.o psb_c_mat_dist_mod.o psb_z_mat_dist_mod.o \
psb_renum_mod.o psb_gps_mod.o
psb_renum_mod.o psb_gps_mod.o \
psb_s_renum_mod.o psb_d_renum_mod.o psb_c_renum_mod.o psb_z_renum_mod.o
IMPLOBJS= psb_s_hbio_impl.o psb_d_hbio_impl.o \
psb_c_hbio_impl.o psb_z_hbio_impl.o \
psb_s_mmio_impl.o psb_d_mmio_impl.o \
@ -40,6 +41,7 @@ $(OBJS): $(MODDIR)/$(BASEMODNAME)$(.mod)
psb_util_mod.o: $(BASEOBJS)
psb_metispart_mod.o: psb_metis_int.o
psb_mat_dist_mod.o: psb_s_mat_dist_mod.o psb_d_mat_dist_mod.o psb_c_mat_dist_mod.o psb_z_mat_dist_mod.o
psb_renum_mod.o: psb_s_renum_mod.o psb_d_renum_mod.o psb_c_renum_mod.o psb_z_renum_mod.o
$(IMPLOBJS): $(BASEOBJS)

370
util/psb_c_renum_impl.F90
Unescape Escape View File

@ -29,37 +29,61 @@
! POSSIBILITY OF SUCH DAMAGE.
!
!
subroutine psb_c_mat_renums(alg,mat,info,perm)
subroutine psb_c_mat_renum(alg,mat,info,perm)
use psb_base_mod
use psb_renum_mod, psb_protect_name => psb_c_mat_renums
use psb_renum_mod, psb_protect_name => psb_c_mat_renum
implicit none
character(len=*), intent(in) :: alg
type(psb_cspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
integer(psb_ipk_) :: err_act, nr, nc, ialg
integer(psb_ipk_) :: err_act, nr, nc, i, ierr(5)
character(len=20) :: name
info = psb_success_
name = 'mat_renum'
call psb_erractionsave(err_act)
nr = mat%get_nrows()
nc = mat%get_ncols()
if (nr /= nc) then
info = psb_err_rectangular_mat_unsupported_
ierr(1) = nr; ierr(2) = nc;
call psb_errpush(info,name,i_err=ierr)
goto 9999
end if
info = psb_success_
select case (psb_toupper(alg))
case ('GPS')
ialg = psb_mat_renum_gps_
case ('AMD')
ialg = psb_mat_renum_amd_
call psb_mat_renum_gps(mat,info,perm)
case('AMD')
call psb_mat_renum_amd(mat,info,perm)
case ('NONE', 'ID')
ialg = psb_mat_renum_identity_
nr = mat%get_nrows()
allocate(perm(nr),stat=info)
if (info == 0) then
do i=1,nr
perm(i) = i
end do
else
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
endif
case default
write(0,*) 'Unknown algorithm "',psb_toupper(alg),'"'
ialg = -1
info = psb_err_input_value_invalid_i_
ierr(1) = 1;
call psb_errpush(info,name,i_err=ierr)
goto 9999
end select
call psb_mat_renum(ialg,mat,info,perm)
if (info /= psb_success_) then
info = psb_err_from_subroutine_non_
call psb_errpush(info,name)
@ -71,26 +95,219 @@ subroutine psb_c_mat_renums(alg,mat,info,perm)
9999 call psb_error_handler(err_act)
return
end subroutine psb_c_mat_renums
subroutine psb_c_mat_renum(alg,mat,info,perm)
contains
subroutine psb_mat_renum_gps(a,info,operm)
use psb_base_mod
use psb_gps_mod
implicit none
type(psb_cspmat_type), intent(inout) :: a
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: operm(:)
!
class(psb_c_base_sparse_mat), allocatable :: aa
type(psb_c_csr_sparse_mat) :: acsr
type(psb_c_coo_sparse_mat) :: acoo
integer(psb_ipk_) :: err_act
character(len=20) :: name
integer(psb_ipk_), allocatable :: ndstk(:,:), iold(:), ndeg(:), perm(:)
integer(psb_ipk_) :: i, j, k, ideg, nr, ibw, ipf, idpth
info = psb_success_
name = 'mat_renum_gps'
call psb_erractionsave(err_act)
info = psb_success_
call a%mold(aa)
call a%mv_to(aa)
call aa%mv_to_fmt(acsr,info)
! Insert call to gps_reduce
nr = acsr%get_nrows()
ideg = 0
do i=1, nr
ideg = max(ideg,acsr%irp(i+1)-acsr%irp(i))
end do
allocate(ndstk(nr,ideg), iold(nr), perm(nr+1), ndeg(nr),stat=info)
if (info /= 0) then
info = psb_err_alloc_dealloc_
call psb_errpush(info, name)
goto 9999
end if
do i=1, nr
iold(i) = i
ndstk(i,:) = 0
k = 0
do j=acsr%irp(i),acsr%irp(i+1)-1
k = k + 1
ndstk(i,k) = acsr%ja(j)
end do
end do
perm = 0
call psb_gps_reduce(ndstk,nr,ideg,iold,perm,ndeg,ibw,ipf,idpth)
if (.not.psb_isaperm(nr,perm)) then
write(0,*) 'Something wrong: bad perm from gps_reduce'
info = psb_err_from_subroutine_
call psb_errpush(info,name)
goto 9999
end if
! Move to coordinate to apply renumbering
call acsr%mv_to_coo(acoo,info)
do i=1, acoo%get_nzeros()
acoo%ia(i) = perm(acoo%ia(i))
acoo%ja(i) = perm(acoo%ja(i))
end do
call acoo%fix(info)
! Get back to where we started from
call aa%mv_from_coo(acoo,info)
call a%mv_from(aa)
if (present(operm)) then
call psb_realloc(nr,operm,info)
if (info /= psb_success_) then
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
end if
operm(1:nr) = perm(1:nr)
end if
deallocate(aa)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_mat_renum_gps
subroutine psb_mat_renum_amd(a,info,operm)
#if defined(HAVE_AMD)
use iso_c_binding
#endif
use psb_base_mod
implicit none
type(psb_cspmat_type), intent(inout) :: a
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: operm(:)
!
#if defined(HAVE_AMD)
interface
function psb_amd_order(n,ap,ai,p)&
& result(res) bind(c,name='psb_amd_order')
use iso_c_binding
integer(c_int) :: res
integer(c_int), value :: n
integer(c_int) :: ap(*), ai(*), p(*)
end function psb_amd_order
end interface
#endif
type(psb_c_csc_sparse_mat) :: acsc
class(psb_c_base_sparse_mat), allocatable :: aa
type(psb_c_coo_sparse_mat) :: acoo
integer(psb_ipk_), allocatable :: perm(:)
integer(psb_ipk_) :: err_act
character(len=20) :: name
integer(psb_ipk_) :: i, j, k, ideg, nr, ibw, ipf, idpth, nz
info = psb_success_
name = 'mat_renum_amd'
call psb_erractionsave(err_act)
#if defined(HAVE_AMD) && defined(IPK4)
info = psb_success_
nr = a%get_nrows()
nz = a%get_nzeros()
allocate(perm(nr),stat=info)
if (info /= 0) then
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
end if
allocate(aa, mold=a%a)
call a%mv_to(acsc)
acsc%ia(:) = acsc%ia(:) - 1
acsc%icp(:) = acsc%icp(:) - 1
info = psb_amd_order(nr,acsc%icp,acsc%ia,perm)
if (info /= psb_success_) then
info = psb_err_from_subroutine_
call psb_errpush(info,name,a_err='psb_amd_order')
goto 9999
end if
perm(:) = perm(:) + 1
acsc%ia(:) = acsc%ia(:) + 1
acsc%icp(:) = acsc%icp(:) + 1
call acsc%mv_to_coo(acoo,info)
do i=1, acoo%get_nzeros()
acoo%ia(i) = perm(acoo%ia(i))
acoo%ja(i) = perm(acoo%ja(i))
end do
call acoo%fix(info)
! Get back to where we started from
call aa%mv_from_coo(acoo,info)
call a%mv_from(aa)
if (present(operm)) then
call psb_realloc(nr,operm,info)
if (info /= psb_success_) then
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
end if
operm(1:nr) = perm(1:nr)
end if
deallocate(aa,perm)
#else
info = psb_err_missing_aux_lib_
call psb_errpush(info, name)
goto 9999
#endif
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_mat_renum_amd
end subroutine psb_c_mat_renum
subroutine psb_lc_mat_renum(alg,mat,info,perm)
use psb_base_mod
use psb_renum_mod, psb_protect_name => psb_c_mat_renum
use psb_renum_mod, psb_protect_name => psb_lc_mat_renum
implicit none
integer(psb_ipk_), intent(in) :: alg
type(psb_cspmat_type), intent(inout) :: mat
character(len=*), intent(in) :: alg
type(psb_lcspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
integer(psb_lpk_), allocatable, optional, intent(out) :: perm(:)
integer(psb_ipk_) :: err_act, nr, nc, i, ierr(5)
integer(psb_lpk_) :: nr, nc, i
integer(psb_ipk_) :: err_act, ierr(5)
character(len=20) :: name
info = psb_success_
name = 'mat_renum'
call psb_erractionsave(err_act)
info = psb_success_
nr = mat%get_nrows()
nc = mat%get_ncols()
if (nr /= nc) then
@ -100,16 +317,17 @@ subroutine psb_c_mat_renum(alg,mat,info,perm)
goto 9999
end if
select case (alg)
case(psb_mat_renum_gps_)
info = psb_success_
select case (psb_toupper(alg))
case ('GPS')
call psb_mat_renum_gps(mat,info,perm)
call psb_lmat_renum_gps(mat,info,perm)
case(psb_mat_renum_amd_)
case('AMD')
call psb_mat_renum_amd(mat,info,perm)
call psb_lmat_renum_amd(mat,info,perm)
case(psb_mat_renum_identity_)
case ('NONE', 'ID')
nr = mat%get_nrows()
allocate(perm(nr),stat=info)
if (info == 0) then
@ -122,8 +340,9 @@ subroutine psb_c_mat_renum(alg,mat,info,perm)
goto 9999
endif
case default
write(0,*) 'Unknown algorithm "',psb_toupper(alg),'"'
info = psb_err_input_value_invalid_i_
ierr(1) = 1; ierr(2) = alg;
ierr(1) = 1;
call psb_errpush(info,name,i_err=ierr)
goto 9999
end select
@ -142,26 +361,26 @@ subroutine psb_c_mat_renum(alg,mat,info,perm)
contains
subroutine psb_mat_renum_gps(a,info,operm)
subroutine psb_lmat_renum_gps(a,info,operm)
use psb_base_mod
use psb_gps_mod
use psb_lgps_mod
implicit none
type(psb_cspmat_type), intent(inout) :: a
type(psb_lcspmat_type), intent(inout) :: a
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: operm(:)
integer(psb_lpk_), allocatable, optional, intent(out) :: operm(:)
!
class(psb_c_base_sparse_mat), allocatable :: aa
type(psb_c_csr_sparse_mat) :: acsr
type(psb_c_coo_sparse_mat) :: acoo
class(psb_lc_base_sparse_mat), allocatable :: aa
type(psb_lc_csr_sparse_mat) :: acsr
type(psb_lc_coo_sparse_mat) :: acoo
integer(psb_ipk_) :: err_act
character(len=20) :: name
integer(psb_ipk_), allocatable :: ndstk(:,:), iold(:), ndeg(:), perm(:)
integer(psb_ipk_) :: i, j, k, ideg, nr, ibw, ipf, idpth
integer(psb_lpk_), allocatable :: ndstk(:,:), iold(:), ndeg(:), perm(:)
integer(psb_lpk_) :: i, j, k, ideg, nr, ibw, ipf, idpth
info = psb_success_
name = 'mat_renum'
name = 'mat_renum_gps'
call psb_erractionsave(err_act)
info = psb_success_
@ -192,7 +411,7 @@ contains
end do
perm = 0
call psb_gps_reduce(ndstk,nr,ideg,iold,perm,ndeg,ibw,ipf,idpth)
call psb_lgps_reduce(ndstk,nr,ideg,iold,perm,ndeg,ibw,ipf,idpth)
if (.not.psb_isaperm(nr,perm)) then
write(0,*) 'Something wrong: bad perm from gps_reduce'
@ -228,18 +447,18 @@ contains
9999 call psb_error_handler(err_act)
return
end subroutine psb_mat_renum_gps
end subroutine psb_lmat_renum_gps
subroutine psb_mat_renum_amd(a,info,operm)
subroutine psb_lmat_renum_amd(a,info,operm)
#if defined(HAVE_AMD)
use iso_c_binding
#endif
use psb_base_mod
implicit none
type(psb_cspmat_type), intent(inout) :: a
type(psb_lcspmat_type), intent(inout) :: a
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: operm(:)
integer(psb_lpk_), allocatable, optional, intent(out) :: operm(:)
!
#if defined(HAVE_AMD)
@ -254,20 +473,20 @@ contains
end interface
#endif
type(psb_c_csc_sparse_mat) :: acsc
class(psb_c_base_sparse_mat), allocatable :: aa
type(psb_c_coo_sparse_mat) :: acoo
type(psb_lc_csc_sparse_mat) :: acsc
class(psb_lc_base_sparse_mat), allocatable :: aa
type(psb_lc_coo_sparse_mat) :: acoo
integer(psb_ipk_), allocatable :: perm(:)
integer(psb_ipk_) :: err_act
character(len=20) :: name
integer(psb_ipk_) :: i, j, k, ideg, nr, ibw, ipf, idpth, nz
integer(psb_lpk_) :: i, j, k, ideg, nr, ibw, ipf, idpth, nz
info = psb_success_
name = 'mat_renum_amd'
call psb_erractionsave(err_act)
#if defined(HAVE_AMD) && defined(IPK4)
#if defined(HAVE_AMD) && defined(LPK4)
info = psb_success_
nr = a%get_nrows()
@ -331,11 +550,9 @@ contains
9999 call psb_error_handler(err_act)
return
end subroutine psb_lmat_renum_amd
end subroutine psb_mat_renum_amd
end subroutine psb_c_mat_renum
end subroutine psb_lc_mat_renum
subroutine psb_c_cmp_bwpf(mat,bwl,bwu,prf,info)
use psb_base_mod
@ -386,3 +603,52 @@ subroutine psb_c_cmp_bwpf(mat,bwl,bwu,prf,info)
end select
end subroutine psb_c_cmp_bwpf
subroutine psb_lc_cmp_bwpf(mat,bwl,bwu,prf,info)
use psb_base_mod
use psb_renum_mod, psb_protect_name => psb_lc_cmp_bwpf
implicit none
type(psb_lcspmat_type), intent(in) :: mat
integer(psb_lpk_), intent(out) :: bwl, bwu
integer(psb_lpk_), intent(out) :: prf
integer(psb_ipk_), intent(out) :: info
!
integer(psb_lpk_), allocatable :: irow(:), icol(:)
complex(psb_spk_), allocatable :: val(:)
integer(psb_lpk_) :: nz, i, j, lrbu, lrbl
info = psb_success_
bwl = 0
bwu = 0
prf = 0
select type (aa=>mat%a)
class is (psb_lc_csr_sparse_mat)
do i=1, aa%get_nrows()
lrbl = 0
lrbu = 0
do j = aa%irp(i), aa%irp(i+1) - 1
lrbl = max(lrbl,i-aa%ja(j))
lrbu = max(lrbu,aa%ja(j)-i)
end do
prf = prf + lrbl+lrbu
bwu = max(bwu,lrbu)
bwl = max(bwl,lrbu)
end do
class default
do i=1, aa%get_nrows()
lrbl = 0
lrbu = 0
call aa%csget(i,i,nz,irow,icol,val,info)
if (info /= psb_success_) return
do j=1, nz
lrbl = max(lrbl,i-icol(j))
lrbu = max(lrbu,icol(j)-i)
end do
prf = prf + lrbl+lrbu
bwu = max(bwu,lrbu)
bwl = max(bwl,lrbu)
end do
end select
end subroutine psb_lc_cmp_bwpf

69
util/psb_c_renum_mod.f90
Unescape Escape View File

@ -0,0 +1,69 @@
!
! Parallel Sparse BLAS version 3.5
! (C) Copyright 2006-2018
! Salvatore Filippone
! Alfredo Buttari
!
! 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.
!
!
module psb_c_renum_mod
use psb_base_mod
interface psb_mat_renum
subroutine psb_c_mat_renum(alg,mat,info,perm)
import :: psb_ipk_, psb_cspmat_type
character(len=*), intent(in) :: alg
type(psb_cspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
end subroutine psb_c_mat_renum
subroutine psb_lc_mat_renum(alg,mat,info,perm)
import :: psb_ipk_, psb_lpk_, psb_lcspmat_type
character(len=*), intent(in) :: alg
type(psb_lcspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_lpk_), allocatable, optional, intent(out) :: perm(:)
end subroutine psb_lc_mat_renum
end interface psb_mat_renum
interface psb_cmp_bwpf
subroutine psb_c_cmp_bwpf(mat,bwl,bwu,prf,info)
import :: psb_ipk_, psb_cspmat_type
type(psb_cspmat_type), intent(in) :: mat
integer(psb_ipk_), intent(out) :: bwl, bwu
integer(psb_ipk_), intent(out) :: prf
integer(psb_ipk_), intent(out) :: info
end subroutine psb_c_cmp_bwpf
subroutine psb_lc_cmp_bwpf(mat,bwl,bwu,prf,info)
import :: psb_ipk_, psb_lpk_, psb_lcspmat_type
type(psb_lcspmat_type), intent(in) :: mat
integer(psb_lpk_), intent(out) :: bwl, bwu
integer(psb_lpk_), intent(out) :: prf
integer(psb_ipk_), intent(out) :: info
end subroutine psb_lc_cmp_bwpf
end interface psb_cmp_bwpf
end module psb_c_renum_mod

367
util/psb_d_renum_impl.F90
Unescape Escape View File

@ -29,37 +29,61 @@
! POSSIBILITY OF SUCH DAMAGE.
!
!
subroutine psb_d_mat_renums(alg,mat,info,perm)
subroutine psb_d_mat_renum(alg,mat,info,perm)
use psb_base_mod
use psb_renum_mod, psb_protect_name => psb_d_mat_renums
use psb_renum_mod, psb_protect_name => psb_d_mat_renum
implicit none
character(len=*), intent(in) :: alg
type(psb_dspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
integer(psb_ipk_) :: err_act, nr, nc, ialg
integer(psb_ipk_) :: err_act, nr, nc, i, ierr(5)
character(len=20) :: name
info = psb_success_
name = 'mat_renum'
call psb_erractionsave(err_act)
nr = mat%get_nrows()
nc = mat%get_ncols()
if (nr /= nc) then
info = psb_err_rectangular_mat_unsupported_
ierr(1) = nr; ierr(2) = nc;
call psb_errpush(info,name,i_err=ierr)
goto 9999
end if
info = psb_success_
select case (psb_toupper(alg))
case ('GPS')
ialg = psb_mat_renum_gps_
case ('AMD')
ialg = psb_mat_renum_amd_
call psb_mat_renum_gps(mat,info,perm)
case('AMD')
call psb_mat_renum_amd(mat,info,perm)
case ('NONE', 'ID')
ialg = psb_mat_renum_identity_
nr = mat%get_nrows()
allocate(perm(nr),stat=info)
if (info == 0) then
do i=1,nr
perm(i) = i
end do
else
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
endif
case default
write(0,*) 'Unknown algorithm "',psb_toupper(alg),'"'
ialg = -1
info = psb_err_input_value_invalid_i_
ierr(1) = 1;
call psb_errpush(info,name,i_err=ierr)
goto 9999
end select
call psb_mat_renum(ialg,mat,info,perm)
if (info /= psb_success_) then
info = psb_err_from_subroutine_non_
call psb_errpush(info,name)
@ -71,26 +95,219 @@ subroutine psb_d_mat_renums(alg,mat,info,perm)
9999 call psb_error_handler(err_act)
return
end subroutine psb_d_mat_renums
subroutine psb_d_mat_renum(alg,mat,info,perm)
contains
subroutine psb_mat_renum_gps(a,info,operm)
use psb_base_mod
use psb_gps_mod
implicit none
type(psb_dspmat_type), intent(inout) :: a
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: operm(:)
!
class(psb_d_base_sparse_mat), allocatable :: aa
type(psb_d_csr_sparse_mat) :: acsr
type(psb_d_coo_sparse_mat) :: acoo
integer(psb_ipk_) :: err_act
character(len=20) :: name
integer(psb_ipk_), allocatable :: ndstk(:,:), iold(:), ndeg(:), perm(:)
integer(psb_ipk_) :: i, j, k, ideg, nr, ibw, ipf, idpth
info = psb_success_
name = 'mat_renum_gps'
call psb_erractionsave(err_act)
info = psb_success_
call a%mold(aa)
call a%mv_to(aa)
call aa%mv_to_fmt(acsr,info)
! Insert call to gps_reduce
nr = acsr%get_nrows()
ideg = 0
do i=1, nr
ideg = max(ideg,acsr%irp(i+1)-acsr%irp(i))
end do
allocate(ndstk(nr,ideg), iold(nr), perm(nr+1), ndeg(nr),stat=info)
if (info /= 0) then
info = psb_err_alloc_dealloc_
call psb_errpush(info, name)
goto 9999
end if
do i=1, nr
iold(i) = i
ndstk(i,:) = 0
k = 0
do j=acsr%irp(i),acsr%irp(i+1)-1
k = k + 1
ndstk(i,k) = acsr%ja(j)
end do
end do
perm = 0
call psb_gps_reduce(ndstk,nr,ideg,iold,perm,ndeg,ibw,ipf,idpth)
if (.not.psb_isaperm(nr,perm)) then
write(0,*) 'Something wrong: bad perm from gps_reduce'
info = psb_err_from_subroutine_
call psb_errpush(info,name)
goto 9999
end if
! Move to coordinate to apply renumbering
call acsr%mv_to_coo(acoo,info)
do i=1, acoo%get_nzeros()
acoo%ia(i) = perm(acoo%ia(i))
acoo%ja(i) = perm(acoo%ja(i))
end do
call acoo%fix(info)
! Get back to where we started from
call aa%mv_from_coo(acoo,info)
call a%mv_from(aa)
if (present(operm)) then
call psb_realloc(nr,operm,info)
if (info /= psb_success_) then
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
end if
operm(1:nr) = perm(1:nr)
end if
deallocate(aa)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_mat_renum_gps
subroutine psb_mat_renum_amd(a,info,operm)
#if defined(HAVE_AMD)
use iso_c_binding
#endif
use psb_base_mod
implicit none
type(psb_dspmat_type), intent(inout) :: a
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: operm(:)
!
#if defined(HAVE_AMD)
interface
function psb_amd_order(n,ap,ai,p)&
& result(res) bind(c,name='psb_amd_order')
use iso_c_binding
integer(c_int) :: res
integer(c_int), value :: n
integer(c_int) :: ap(*), ai(*), p(*)
end function psb_amd_order
end interface
#endif
type(psb_d_csc_sparse_mat) :: acsc
class(psb_d_base_sparse_mat), allocatable :: aa
type(psb_d_coo_sparse_mat) :: acoo
integer(psb_ipk_), allocatable :: perm(:)
integer(psb_ipk_) :: err_act
character(len=20) :: name
integer(psb_ipk_) :: i, j, k, ideg, nr, ibw, ipf, idpth, nz
info = psb_success_
name = 'mat_renum_amd'
call psb_erractionsave(err_act)
#if defined(HAVE_AMD) && defined(IPK4)
info = psb_success_
nr = a%get_nrows()
nz = a%get_nzeros()
allocate(perm(nr),stat=info)
if (info /= 0) then
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
end if
allocate(aa, mold=a%a)
call a%mv_to(acsc)
acsc%ia(:) = acsc%ia(:) - 1
acsc%icp(:) = acsc%icp(:) - 1
info = psb_amd_order(nr,acsc%icp,acsc%ia,perm)
if (info /= psb_success_) then
info = psb_err_from_subroutine_
call psb_errpush(info,name,a_err='psb_amd_order')
goto 9999
end if
perm(:) = perm(:) + 1
acsc%ia(:) = acsc%ia(:) + 1
acsc%icp(:) = acsc%icp(:) + 1
call acsc%mv_to_coo(acoo,info)
do i=1, acoo%get_nzeros()
acoo%ia(i) = perm(acoo%ia(i))
acoo%ja(i) = perm(acoo%ja(i))
end do
call acoo%fix(info)
! Get back to where we started from
call aa%mv_from_coo(acoo,info)
call a%mv_from(aa)
if (present(operm)) then
call psb_realloc(nr,operm,info)
if (info /= psb_success_) then
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
end if
operm(1:nr) = perm(1:nr)
end if
deallocate(aa,perm)
#else
info = psb_err_missing_aux_lib_
call psb_errpush(info, name)
goto 9999
#endif
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_mat_renum_amd
end subroutine psb_d_mat_renum
subroutine psb_ld_mat_renum(alg,mat,info,perm)
use psb_base_mod
use psb_renum_mod, psb_protect_name => psb_d_mat_renum
use psb_renum_mod, psb_protect_name => psb_ld_mat_renum
implicit none
integer(psb_ipk_), intent(in) :: alg
type(psb_dspmat_type), intent(inout) :: mat
character(len=*), intent(in) :: alg
type(psb_ldspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
integer(psb_lpk_), allocatable, optional, intent(out) :: perm(:)
integer(psb_ipk_) :: err_act, nr, nc, i, ierr(5)
integer(psb_lpk_) :: nr, nc, i
integer(psb_ipk_) :: err_act, ierr(5)
character(len=20) :: name
info = psb_success_
name = 'mat_renum'
call psb_erractionsave(err_act)
info = psb_success_
nr = mat%get_nrows()
nc = mat%get_ncols()
if (nr /= nc) then
@ -100,16 +317,17 @@ subroutine psb_d_mat_renum(alg,mat,info,perm)
goto 9999
end if
select case (alg)
case(psb_mat_renum_gps_)
info = psb_success_
select case (psb_toupper(alg))
case ('GPS')
call psb_mat_renum_gps(mat,info,perm)
call psb_lmat_renum_gps(mat,info,perm)
case(psb_mat_renum_amd_)
case('AMD')
call psb_mat_renum_amd(mat,info,perm)
call psb_lmat_renum_amd(mat,info,perm)
case(psb_mat_renum_identity_)
case ('NONE', 'ID')
nr = mat%get_nrows()
allocate(perm(nr),stat=info)
if (info == 0) then
@ -122,8 +340,9 @@ subroutine psb_d_mat_renum(alg,mat,info,perm)
goto 9999
endif
case default
write(0,*) 'Unknown algorithm "',psb_toupper(alg),'"'
info = psb_err_input_value_invalid_i_
ierr(1) = 1; ierr(2) = alg;
ierr(1) = 1;
call psb_errpush(info,name,i_err=ierr)
goto 9999
end select
@ -142,23 +361,23 @@ subroutine psb_d_mat_renum(alg,mat,info,perm)
contains
subroutine psb_mat_renum_gps(a,info,operm)
subroutine psb_lmat_renum_gps(a,info,operm)
use psb_base_mod
use psb_gps_mod
use psb_lgps_mod
implicit none
type(psb_dspmat_type), intent(inout) :: a
type(psb_ldspmat_type), intent(inout) :: a
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: operm(:)
integer(psb_lpk_), allocatable, optional, intent(out) :: operm(:)
!
class(psb_d_base_sparse_mat), allocatable :: aa
type(psb_d_csr_sparse_mat) :: acsr
type(psb_d_coo_sparse_mat) :: acoo
class(psb_ld_base_sparse_mat), allocatable :: aa
type(psb_ld_csr_sparse_mat) :: acsr
type(psb_ld_coo_sparse_mat) :: acoo
integer(psb_ipk_) :: err_act
character(len=20) :: name
integer(psb_ipk_), allocatable :: ndstk(:,:), iold(:), ndeg(:), perm(:)
integer(psb_ipk_) :: i, j, k, ideg, nr, ibw, ipf, idpth
integer(psb_lpk_), allocatable :: ndstk(:,:), iold(:), ndeg(:), perm(:)
integer(psb_lpk_) :: i, j, k, ideg, nr, ibw, ipf, idpth
info = psb_success_
name = 'mat_renum_gps'
@ -192,7 +411,7 @@ contains
end do
perm = 0
call psb_gps_reduce(ndstk,nr,ideg,iold,perm,ndeg,ibw,ipf,idpth)
call psb_lgps_reduce(ndstk,nr,ideg,iold,perm,ndeg,ibw,ipf,idpth)
if (.not.psb_isaperm(nr,perm)) then
write(0,*) 'Something wrong: bad perm from gps_reduce'
@ -228,18 +447,18 @@ contains
9999 call psb_error_handler(err_act)
return
end subroutine psb_mat_renum_gps
end subroutine psb_lmat_renum_gps
subroutine psb_mat_renum_amd(a,info,operm)
subroutine psb_lmat_renum_amd(a,info,operm)
#if defined(HAVE_AMD)
use iso_c_binding
#endif
use psb_base_mod
implicit none
type(psb_dspmat_type), intent(inout) :: a
type(psb_ldspmat_type), intent(inout) :: a
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: operm(:)
integer(psb_lpk_), allocatable, optional, intent(out) :: operm(:)
!
#if defined(HAVE_AMD)
@ -254,20 +473,20 @@ contains
end interface
#endif
type(psb_d_csc_sparse_mat) :: acsc
class(psb_d_base_sparse_mat), allocatable :: aa
type(psb_d_coo_sparse_mat) :: acoo
type(psb_ld_csc_sparse_mat) :: acsc
class(psb_ld_base_sparse_mat), allocatable :: aa
type(psb_ld_coo_sparse_mat) :: acoo
integer(psb_ipk_), allocatable :: perm(:)
integer(psb_ipk_) :: err_act
character(len=20) :: name
integer(psb_ipk_) :: i, j, k, ideg, nr, ibw, ipf, idpth, nz
integer(psb_lpk_) :: i, j, k, ideg, nr, ibw, ipf, idpth, nz
info = psb_success_
name = 'mat_renum_amd'
call psb_erractionsave(err_act)
#if defined(HAVE_AMD) && defined(IPK4)
#if defined(HAVE_AMD) && defined(LPK4)
info = psb_success_
nr = a%get_nrows()
@ -331,10 +550,9 @@ contains
9999 call psb_error_handler(err_act)
return
end subroutine psb_mat_renum_amd
end subroutine psb_d_mat_renum
end subroutine psb_lmat_renum_amd
end subroutine psb_ld_mat_renum
subroutine psb_d_cmp_bwpf(mat,bwl,bwu,prf,info)
use psb_base_mod
@ -385,3 +603,52 @@ subroutine psb_d_cmp_bwpf(mat,bwl,bwu,prf,info)
end select
end subroutine psb_d_cmp_bwpf
subroutine psb_ld_cmp_bwpf(mat,bwl,bwu,prf,info)
use psb_base_mod
use psb_renum_mod, psb_protect_name => psb_ld_cmp_bwpf
implicit none
type(psb_ldspmat_type), intent(in) :: mat
integer(psb_lpk_), intent(out) :: bwl, bwu
integer(psb_lpk_), intent(out) :: prf
integer(psb_ipk_), intent(out) :: info
!
integer(psb_lpk_), allocatable :: irow(:), icol(:)
real(psb_dpk_), allocatable :: val(:)
integer(psb_lpk_) :: nz, i, j, lrbu, lrbl
info = psb_success_
bwl = 0
bwu = 0
prf = 0
select type (aa=>mat%a)
class is (psb_ld_csr_sparse_mat)
do i=1, aa%get_nrows()
lrbl = 0
lrbu = 0
do j = aa%irp(i), aa%irp(i+1) - 1
lrbl = max(lrbl,i-aa%ja(j))
lrbu = max(lrbu,aa%ja(j)-i)
end do
prf = prf + lrbl+lrbu
bwu = max(bwu,lrbu)
bwl = max(bwl,lrbu)
end do
class default
do i=1, aa%get_nrows()
lrbl = 0
lrbu = 0
call aa%csget(i,i,nz,irow,icol,val,info)
if (info /= psb_success_) return
do j=1, nz
lrbl = max(lrbl,i-icol(j))
lrbu = max(lrbu,icol(j)-i)
end do
prf = prf + lrbl+lrbu
bwu = max(bwu,lrbu)
bwl = max(bwl,lrbu)
end do
end select
end subroutine psb_ld_cmp_bwpf

69
util/psb_d_renum_mod.f90
Unescape Escape View File

@ -0,0 +1,69 @@
!
! Parallel Sparse BLAS version 3.5
! (C) Copyright 2006-2018
! Salvatore Filippone
! Alfredo Buttari
!
! 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.
!
!
module psb_d_renum_mod
use psb_base_mod
interface psb_mat_renum
subroutine psb_d_mat_renum(alg,mat,info,perm)
import :: psb_ipk_, psb_dspmat_type
character(len=*), intent(in) :: alg
type(psb_dspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
end subroutine psb_d_mat_renum
subroutine psb_ld_mat_renum(alg,mat,info,perm)
import :: psb_ipk_, psb_lpk_, psb_ldspmat_type
character(len=*), intent(in) :: alg
type(psb_ldspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_lpk_), allocatable, optional, intent(out) :: perm(:)
end subroutine psb_ld_mat_renum
end interface psb_mat_renum
interface psb_cmp_bwpf
subroutine psb_d_cmp_bwpf(mat,bwl,bwu,prf,info)
import :: psb_ipk_, psb_dspmat_type
type(psb_dspmat_type), intent(in) :: mat
integer(psb_ipk_), intent(out) :: bwl, bwu
integer(psb_ipk_), intent(out) :: prf
integer(psb_ipk_), intent(out) :: info
end subroutine psb_d_cmp_bwpf
subroutine psb_ld_cmp_bwpf(mat,bwl,bwu,prf,info)
import :: psb_ipk_, psb_lpk_, psb_ldspmat_type
type(psb_ldspmat_type), intent(in) :: mat
integer(psb_lpk_), intent(out) :: bwl, bwu
integer(psb_lpk_), intent(out) :: prf
integer(psb_ipk_), intent(out) :: info
end subroutine psb_ld_cmp_bwpf
end interface psb_cmp_bwpf
end module psb_d_renum_mod

742
util/psb_gps_mod.f90
Unescape Escape View File

@ -774,4 +774,744 @@ CONTAINS
RETURN
END SUBROUTINE REALLOC
!
END MODULE psb_gps_mod
end module psb_gps_mod
module psb_lgps_mod
!
use psb_base_mod, only : psb_lpk_
public psb_lgps_reduce
!
! COMMON /GRA/ N, IDPTH, IDEG
!
private
! common /CC/ XCC,SIZEG,STPT
integer(psb_lpk_), save :: xcc,n,idpth,ideg
integer(psb_lpk_),allocatable,SAVE :: SIZEG(:),STPT(:)
!
! COMMON /LVLW/ NHIGH,NLOW,NACUM
integer(psb_lpk_),allocatable,target,save :: NHIGH(:),NLOW(:),NACUM(:),AUX(:)
integer(psb_lpk_),PARAMETER :: INIT=500
!
CONTAINS
!
!!$ SUBROUTINE psb_gps_reduce(NDSTK, NR, IOLD, RENUM, NDEG, LVL, LVLS1, LVLS2,&
!!$ & CCSTOR, IBW2, IPF2,NE,IDPTHE,IDEGE)
SUBROUTINE psb_lgps_reduce(NDSTK, NR, IDEGE, IOLD, RENUM, NDEG,ibw2,ipf2,IDPTHE)
! SUBROUTINE REDUCE DETERMINES A ROW AND COLUMN PERMUTATION WHICH,
! WHEN APPLIED TO A GIVEN SPARSE MATRIX, PRODUCES A PERMUTED
! MATRIX WITH A SMALLER BANDWIDTH AND PROFILE.
! THE INPUT ARRAY IS A CONNECTION TABLE WHICH REPRESENTS THE
! INDICES OF THE NONZERO ELEMENTS OF THE MATRIX, A. THE ALGO-
! RITHM IS DESCRIBED IN TERMS OF THE ADJACENCY GRAPH WHICH
! HAS THE CHARACTERISTIC THAT THERE IS AN EDGE (CONNECTION)
! BETWEEN NODES I AND J IF A(I,J) /= 0 AND I /= J.
! DIMENSIONING INFORMATION--THE FOLLOWING INTEGER ARRAYS MUST BE
! DIMENSIONED IN THE CALLING ROUTINE.
! NDSTK(NR,D1) D1 IS >= MAXIMUM DEGREE OF ALL NODES.
! IOLD(D2) D2 AND NR ARE >= THE TOTAL NUMBER OF
! RENUM(D2+1) NODES IN THE GRAPH.
! NDEG(D2) STORAGE REQUIREMENTS CAN BE SIGNIFICANTLY
! LVL(D2) DECREASED FOR IBM 360 AND 370 COMPUTERS
! LVLS1(D2) BY REPLACING INTEGER NDSTK BY
! LVLS2(D2) INTEGER*2 NDSTK IN SUBROUTINES REDUCE,
! CCSTOR(D2) DGREE, FNDIAM, TREE AND NUMBER.
! COMMON INFORMATION--THE FOLLOWING COMMON BLOCK MUST BE IN THE
! CALLING ROUTINE.
! COMMON/GRA/N,IDPTH,IDEG
! EXPLANATION OF INPUT VARIABLES--
! NDSTK- CONNECTION TABLE REPRESENTING GRAPH.
! NDSTK(I,J)=NODE NUMBER OF JTH CONNECTION TO NODE
! NUMBER I. A CONNECTION OF A NODE TO ITSELF IS NOT
! LISTED. EXTRA POSITIONS MUST HAVE ZERO FILL.
! NR- ROW DIMENSION ASSIGNED NDSTK IN CALLING PROGRAM.
! IOLD(I)- NUMBERING OF ITH NODE UPON INPUT.
! IF NO NUMBERING EXISTS THEN IOLD(I)=I.
! N- NUMBER OF NODES IN GRAPH (EQUAL TO ORDER OF MATRIX).
! IDEG- MAXIMUM DEGREE OF ANY NODE IN THE GRAPH.
! EXPLANATION OF OUTPUT VARIABLES--
! RENUM(I)- THE NEW NUMBER FOR THE ITH NODE.
! NDEG(I)- THE DEGREE OF THE ITH NODE.
! IBW2- THE BANDWIDTH AFTER RENUMBERING.
! IPF2- THE PROFILE AFTER RENUMBERING.
! IDPTH- NUMBER OF LEVELS IN REDUCE LEVEL STRUCTURE.
! THE FOLLOWING ONLY HAVE MEANING IF THE GRAPH WAS CONNECTED--
! LVL(I)- INDEX INTO LVLS1 TO THE FIRST NODE IN LEVEL I.
! LVL(I+1)-LVL(I)= NUMBER OF NODES IN ITH LEVEL
! LVLS1- NODE NUMBERS LISTED BY LEVEL.
! LVLS2(I)- THE LEVEL ASSIGNED TO NODE I BY REDUCE.
! WORKING STORAGE VARIABLE--
! CCSTOR
! LOCAL STORAGE--
! COMMON/CC/-SUBROUTINES REDUCE, SORT2 AND PIKLVL ASSUME THAT
! THE GRAPH HAS AT MOST 50 CONNECTED COMPONENTS.
! SUBROUTINE FNDIAM ASSUMES THAT THERE ARE AT MOST
! 100 NODES IN THE LAST LEVEL.
! COMMON/LVLW/-SUBROUTINES SETUP AND PIKLVL ASSUME THAT THERE
! ARE AT MOST 100 LEVELS.
! USE INTEGER*2 NDSTK WITH AN IBM 360 OR 370.
use psb_base_mod
implicit none
INTEGER(psb_lpk_) :: NR, IDEGE, IBW2, IPF2, IDPTHE
! COMMON /GRA/ N, IDPTH, IDEG
! IT IS ASSUMED THAT THE GRAPH HAS AT MOST 50 CONNECTED COMPONENTS.
! COMMON /CC/ XCC, SIZEG(50), STPT(50)
! COMMON /LVLW/ NHIGH(100), NLOW(100), NACUM(100)
integer(psb_lpk_) :: stnode, rvnode, stnum, sbnum
integer(psb_lpk_) :: ndstk(nr,idege), iold(nr), renum(nr+1), ndeg(nr)
integer(psb_lpk_) :: lvl(nr), lvls1(nr), lvls2(nr), ccstor(nr)
integer(psb_lpk_) :: nflg, info, i, ibw1, ipf1, idflt, isdir, lroot, lowdg
integer(psb_lpk_) :: lvlbot, lvln, lvlwth, maxlw, num
n = nr
ideg = idege
idpth = 0
ALLOCATE(SIZEG(NR),STPT(NR), STAT=INFO)
IF(INFO /= psb_success_) THEN
write(psb_out_unit,*) 'ERROR! MEMORY ALLOCATION # 1 FAILED IN GPS'
STOP
END IF
!
ALLOCATE(NHIGH(INIT), NLOW(INIT), NACUM(INIT), AUX(INIT), STAT=INFO)
IF(INFO /= psb_success_) THEN
write(psb_out_unit,*) 'ERROR! MEMORY ALLOCATION # 2 FAILED IN GPS'
STOP
END IF
!
IBW2 = 0
IPF2 = 0
! SET RENUM(I)=0 FOR ALL I TO INDICATE NODE I IS UNNUMBERED
DO I=1,N
RENUM(I) = 0
END DO
! COMPUTE DEGREE OF EACH NODE AND ORIGINAL BANDWIDTH AND PROFILE
CALL DGREE(NDSTK, NR, NDEG, IOLD, IBW1, IPF1)
! SBNUM= LOW END OF AVAILABLE NUMBERS FOR RENUMBERING
! STNUM= HIGH END OF AVAILABLE NUMBERS FOR RENUMBERING
SBNUM = 1
STNUM = N
! NUMBER THE NODES OF DEGREE ZERO
DO I=1,N
IF (NDEG(I) > 0) CYCLE
RENUM(I) = STNUM
STNUM = STNUM - 1
END DO
! FIND AN UNNUMBERED NODE OF MIN DEGREE TO START ON
do
LOWDG = IDEG + 1
NFLG = 1
ISDIR = 1
DO I=1,N
IF (NDEG(I) >= LOWDG) CYCLE
IF (RENUM(I) > 0) CYCLE
LOWDG = NDEG(I)
STNODE = I
END DO
! FIND PSEUDO-DIAMETER AND ASSOCIATED LEVEL STRUCTURES.
! STNODE AND RVNODE ARE THE ENDS OF THE DIAM AND LVLS1 AND LVLS2
! ARE THE RESPECTIVE LEVEL STRUCTURES.
CALL FNDIAM(STNODE, RVNODE, NDSTK, NR, NDEG, LVL, LVLS1,LVLS2, CCSTOR, IDFLT)
IF (.not.(ndeg(stnode) <= ndeg(rvnode))) then
! NFLG INDICATES THE END TO BEGIN NUMBERING ON
NFLG = -1
STNODE = RVNODE
endif
CALL SETUP(LVL, LVLS1, LVLS2)
! FIND ALL THE CONNECTED COMPONENTS (XCC COUNTS THEM)
XCC = 0
LROOT = 1
LVLN = 1
DO I=1,N
IF (LVL(I) /= 0) CYCLE
XCC = XCC + 1
STPT(XCC) = LROOT
CALL TREE(I, NDSTK, NR, LVL, CCSTOR, NDEG, LVLWTH, LVLBOT,LVLN, MAXLW, N)
SIZEG(XCC) = LVLBOT + LVLWTH - LROOT
LROOT = LVLBOT + LVLWTH
LVLN = LROOT
END DO
if (sort2() /= 0) then
CALL PIKLVL(LVLS1, LVLS2, CCSTOR, IDFLT, ISDIR)
endif
! ON RETURN FROM PIKLVL, ISDIR INDICATES THE DIRECTION THE LARGEST
! COMPONENT FELL. ISDIR IS MODIFIED NOW TO INDICATE THE NUMBERING
! DIRECTION. NUM IS SET TO THE PROPER VALUE FOR THIS DIRECTION.
ISDIR = ISDIR*NFLG
NUM = SBNUM
IF (ISDIR < 0) NUM = STNUM
CALL NUMBER(STNODE, NUM, NDSTK, LVLS2, NDEG, RENUM, LVLS1,LVL,&
& NR, NFLG, IBW2, IPF2, CCSTOR, ISDIR)
! UPDATE STNUM OR SBNUM AFTER NUMBERING
IF (ISDIR < 0) STNUM = NUM
IF (ISDIR > 0) SBNUM = NUM
IF (.not.(sbnum <= stnum)) exit
end do
IF (IBW2 > IBW1) then
! IF ORIGINAL NUMBERING IS BETTER THAN NEW ONE, SET UP TO RETURN IT
DO I=1,N
RENUM(I) = IOLD(I)
END DO
IBW2 = IBW1
IPF2 = IPF1
!
endif
DEALLOCATE(SIZEG,STPT,NHIGH,NLOW,AUX,NACUM)
idpthe = idpth
RETURN
end subroutine psb_lgps_reduce
!
SUBROUTINE DGREE(NDSTK, NR, NDEG, IOLD, IBW1, IPF1)
! DGREE COMPUTES THE DEGREE OF EACH NODE IN NDSTK AND STORES
! IT IN THE ARRAY NDEG. THE BANDWIDTH AND PROFILE FOR THE ORIGINAL
! OR INPUT RENUMBERING OF THE GRAPH IS COMPUTED ALSO.
! USE INTEGER*2 NDSTK WITH AN IBM 360 OR 370.
implicit none
INTEGER(psb_lpk_) :: NR, IBW1, IPF1, NDSTK(NR,IDEG), NDEG(N), IOLD(N)
! COMMON /GRA/ N, IDPTH, IDEG
integer(psb_lpk_) :: i, itst, j, idif, irw
IBW1 = 0
IPF1 = 0
DO I=1,N
NDEG(I) = 0
IRW = 0
DO J=1,IDEG
ITST = NDSTK(I,J)
IF(ITST <= 0) EXIT
NDEG(I) = NDEG(I) + 1
IDIF = IOLD(I) - IOLD(ITST)
IF (IRW < IDIF) IRW = IDIF
END DO
IPF1 = IPF1 + IRW
IF (IRW > IBW1) IBW1 = IRW
END DO
RETURN
END SUBROUTINE DGREE
!
SUBROUTINE FNDIAM(SND1, SND2, NDSTK, NR, NDEG, LVL, LVLS1,LVLS2, IWK, IDFLT)
! FNDIAM IS THE CONTROL PROCEDURE FOR FINDING THE PSEUDO-DIAMETER OF
! NDSTK AS WELL AS THE LEVEL STRUCTURE FROM EACH END
! SND1- ON INPUT THIS IS THE NODE NUMBER OF THE FIRST
! ATTEMPT AT FINDING A DIAMETER. ON OUTPUT IT
! CONTAINS THE ACTUAL NUMBER USED.
! SND2- ON OUTPUT CONTAINS OTHER END OF DIAMETER
! LVLS1- ARRAY CONTAINING LEVEL STRUCTURE WITH SND1 AS ROOT
! LVLS2- ARRAY CONTAINING LEVEL STRUCTURE WITH SND2 AS ROOT
! IDFLT- FLAG USED IN PICKING FINAL LEVEL STRUCTURE, SET
! =1 IF WIDTH OF LVLS1 <= WIDTH OF LVLS2, OTHERWISE =2
! LVL,IWK- WORKING STORAGE
! USE INTEGER*2 NDSTK WITH AN IBM 360 OR 370.
implicit none
INTEGER(psb_lpk_) :: FLAG, SND, SND1, SND2, NR, idflt
! COMMON /GRA/ N, IDPTH, IDEG
! IT IS ASSUMED THAT THE LAST LEVEL HAS AT MOST 100 NODES.
! COMMON /CC/ NDLST(100)
integer(psb_lpk_),POINTER :: NDLST(:)
integer(psb_lpk_) :: NDSTK(NR,IDEG), NDEG(1), LVL(N), LVLS1(N), LVLS2(N),IWK(N)
integer(psb_lpk_) :: i, j, mtw2, ndxn, ndxl, inow, lvlbot,lvln,lvlwth
integer(psb_lpk_) :: k,mtw1, maxlw
!
NDLST => AUX
!
FLAG = 0
MTW2 = N
SND = SND1
! ZERO LVL TO INDICATE ALL NODES ARE AVAILABLE TO TREE
10 DO 20 I=1,N
LVL(I) = 0
20 END DO
LVLN = 1
! DROP A TREE FROM SND
CALL TREE(SND, NDSTK, NR, LVL, IWK, NDEG, LVLWTH, LVLBOT,LVLN, MAXLW, MTW2)
IF (FLAG >= 1) GO TO 50
FLAG = 1
30 IDPTH = LVLN - 1
MTW1 = MAXLW
! COPY LEVEL STRUCTURE INTO LVLS1
DO 40 I=1,N
LVLS1(I) = LVL(I)
40 END DO
NDXN = 1
NDXL = 0
MTW2 = N
! SORT LAST LEVEL BY DEGREE AND STORE IN NDLST
CALL SORTDG(NDLST, IWK(LVLBOT), NDXL, LVLWTH, NDEG)
SND = NDLST(1)
GO TO 10
50 IF (IDPTH >= LVLN-1) GO TO 60
! START AGAIN WITH NEW STARTING NODE
SND1 = SND
GO TO 30
60 IF (MAXLW >= MTW2) GO TO 80
MTW2 = MAXLW
SND2 = SND
! STORE NARROWEST REVERSE LEVEL STRUCTURE IN LVLS2
DO 70 I=1,N
LVLS2(I) = LVL(I)
70 END DO
80 IF (NDXN == NDXL) GO TO 90
! TRY NEXT NODE IN NDLST
NDXN = NDXN + 1
SND = NDLST(NDXN)
GO TO 10
90 IDFLT = 1
IF (MTW2 <= MTW1) IDFLT = 2
NULLIFY(NDLST)
RETURN
END SUBROUTINE FNDIAM
!
SUBROUTINE TREE(IROOT, NDSTK, NR, LVL, IWK, NDEG, LVLWTH, LVLBOT, LVLN, MAXLW, IBORT)
! TREE DROPS A TREE IN NDSTK FROM IROOT
! LVL- ARRAY INDICATING AVAILABLE NODES IN NDSTK WITH ZERO
! ENTRIES. TREE ENTERS LEVEL NUMBERS ASSIGNED
! DURING EXECUTION OF THIS PROCEDURE
! IWK- ON OUTPUT CONTAINS NODE NUMBERS USED IN TREE
! ARRANGED BY LEVELS (IWK(LVLN) CONTAINS IROOT
! AND IWK(LVLBOT+LVLWTH-1) CONTAINS LAST NODE ENTERED)
! LVLWTH- ON OUTPUT CONTAINS WIDTH OF LAST LEVEL
! LVLBOT- ON OUTPUT CONTAINS INDEX INTO IWK OF FIRST
! NODE IN LAST LEVEL
! MAXLW- ON OUTPUT CONTAINS THE MAXIMUM LEVEL WIDTH
! LVLN- ON INPUT THE FIRST AVAILABLE LOCATION IN IWK
! USUALLY ONE BUT IF IWK IS USED TO STORE PREVIOUS
! CONNECTED COMPONENTS, LVLN IS NEXT AVAILABLE LOCATION.
! ON OUTPUT THE TOTAL NUMBER OF LEVELS + 1
! IBORT- INPUT PARAM WHICH TRIGGERS EARLY RETURN IF
! MAXLW BECOMES >= IBORT
! USE INTEGER*2 NDSTK WITH AN IBM 360 OR 370.
implicit none
integer(psb_lpk_) :: IROOT, NR, NDSTK(NR,*), LVL(*), IWK(*), NDEG(*)
integer(psb_lpk_) :: LVLWTH, LVLBOT, LVLN, MAXLW, IBORT
integer(psb_lpk_) :: itest, iwknow, itop, lvltop,j , inow, ndrow
MAXLW = 0
ITOP = LVLN
INOW = LVLN
LVLBOT = LVLN
LVLTOP = LVLN + 1
LVLN = 1
LVL(IROOT) = 1
IWK(ITOP) = IROOT
10 LVLN = LVLN + 1
20 IWKNOW = IWK(INOW)
NDROW = NDEG(IWKNOW)
DO 30 J=1,NDROW
ITEST = NDSTK(IWKNOW,J)
IF (LVL(ITEST) /= 0) CYCLE
LVL(ITEST) = LVLN
ITOP = ITOP + 1
IWK(ITOP) = ITEST
30 END DO
INOW = INOW + 1
IF (INOW < LVLTOP) GO TO 20
LVLWTH = LVLTOP - LVLBOT
IF (MAXLW < LVLWTH) MAXLW = LVLWTH
IF (MAXLW >= IBORT) RETURN
IF (ITOP < LVLTOP) RETURN
LVLBOT = INOW
LVLTOP = ITOP + 1
GO TO 10
END SUBROUTINE TREE
!
SUBROUTINE SORTDG(STK1, STK2, X1, X2, NDEG)
! SORTDG SORTS STK2 BY DEGREE OF THE NODE AND ADDS IT TO THE END
! OF STK1 IN ORDER OF LOWEST TO HIGHEST DEGREE. X1 AND X2 ARE THE
! NUMBER OF NODES IN STK1 AND STK2 RESPECTIVELY.
implicit none
INTEGER(psb_lpk_) :: X1, X2, STK1, STK2, TEMP,NDEG
! COMMON /GRA/ N, IDPTH, IDEG
DIMENSION NDEG(N), STK1(X1+X2), STK2(X2)
integer(psb_lpk_) :: ind,itest,i,j,istk2,jstk2
IND = X2
10 ITEST = 0
IND = IND - 1
IF (IND < 1) GO TO 30
DO 20 I=1,IND
J = I + 1
ISTK2 = STK2(I)
JSTK2 = STK2(J)
IF (NDEG(ISTK2) <= NDEG(JSTK2)) CYCLE
ITEST = 1
TEMP = STK2(I)
STK2(I) = STK2(J)
STK2(J) = TEMP
20 END DO
IF (ITEST == 1) GO TO 10
30 DO 40 I=1,X2
X1 = X1 + 1
STK1(X1) = STK2(I)
40 END DO
RETURN
END SUBROUTINE SORTDG
!
SUBROUTINE SETUP(LVL, LVLS1, LVLS2)
! SETUP COMPUTES THE REVERSE LEVELING INFO FROM LVLS2 AND STORES
! IT INTO LVLS2. NACUM(I) IS INITIALIZED TO NODES/ITH LEVEL FOR NODES
! ON THE PSEUDO-DIAMETER OF THE GRAPH. LVL IS INITIALIZED TO NON-
! ZERO FOR NODES ON THE PSEUDO-DIAM AND NODES IN A DIFFERENT
! COMPONENT OF THE GRAPH.
! COMMON /GRA/ N, IDPTH, IDEG
! IT IS ASSUMED THAT THERE ARE AT MOST 100 LEVELS.
! COMMON /LVLW/ NHIGH(100), NLOW(100), NACUM(100)
use psb_base_mod
implicit none
integer(psb_lpk_) :: LVL(N), LVLS1(N), LVLS2(N)
integer(psb_lpk_) :: SZ,i,itemp
!-----------------------------------------------------
SZ=SIZE(NACUM)
IF(SZ < IDPTH) THEN
write(psb_out_unit,*) 'GPS_SETUP: on fly reallocation of NACUM'
CALL REALLOC(NACUM,SZ,IDPTH)
END IF
!-----------------------------------------------------
DO 10 I=1,IDPTH
NACUM(I) = 0
10 END DO
DO 30 I=1,N
LVL(I) = 1
LVLS2(I) = IDPTH + 1 - LVLS2(I)
ITEMP = LVLS2(I)
IF (ITEMP > IDPTH) CYCLE
IF (ITEMP /= LVLS1(I)) GO TO 20
NACUM(ITEMP) = NACUM(ITEMP) + 1
CYCLE
20 LVL(I) = 0
30 END DO
RETURN
END SUBROUTINE SETUP
!
FUNCTION SORT2() result(val)
implicit none
INTEGER(psb_lpk_) :: val
! SORT2 SORTS SIZEG AND STPT INTO DESCENDING ORDER ACCORDING TO
! VALUES OF SIZEG. XCC=NUMBER OF ENTRIES IN EACH ARRAY
INTEGER(psb_lpk_) :: TEMP,ind,itest,i,j
! IT IS ASSUMED THAT THE GRAPH HAS AT MOST 50 CONNECTED COMPONENTS.
!COMMON /CC/ XCC, SIZEG(50), STPT(50)
VAL = 0
IF (XCC == 0) RETURN
VAL = 1
IND = XCC
10 ITEST = 0
IND = IND - 1
IF (IND < 1) RETURN
DO 20 I=1,IND
J = I + 1
IF (SIZEG(I) >= SIZEG(J)) CYCLE
ITEST = 1
TEMP = SIZEG(I)
SIZEG(I) = SIZEG(J)
SIZEG(J) = TEMP
TEMP = STPT(I)
STPT(I) = STPT(J)
STPT(J) = TEMP
20 END DO
IF (ITEST == 1) GO TO 10
RETURN
END FUNCTION SORT2
!
SUBROUTINE PIKLVL(LVLS1, LVLS2, CCSTOR, IDFLT, ISDIR)
use psb_base_mod
implicit none
! PIKLVL CHOOSES THE LEVEL STRUCTURE USED IN NUMBERING GRAPH
! LVLS1- ON INPUT CONTAINS FORWARD LEVELING INFO
! LVLS2- ON INPUT CONTAINS REVERSE LEVELING INFO
! ON OUTPUT THE FINAL LEVEL STRUCTURE CHOSEN
! CCSTOR- ON INPUT CONTAINS CONNECTED COMPONENT INFO
! IDFLT- ON INPUT =1 IF WDTH LVLS1 <= WDTH LVLS2, =2 OTHERWISE
! NHIGH KEEPS TRACK OF LEVEL WIDTHS FOR HIGH NUMBERING
! NLOW- KEEPS TRACK OF LEVEL WIDTHS FOR LOW NUMBERING
! NACUM- KEEPS TRACK OF LEVEL WIDTHS FOR CHOSEN LEVEL STRUCTURE
! XCC- NUMBER OF CONNECTED COMPONENTS
! SIZEG(I)- SIZE OF ITH CONNECTED COMPONENT
! STPT(I)- INDEX INTO CCSTORE OF 1ST NODE IN ITH CON COMPT
! ISDIR- FLAG WHICH INDICATES WHICH WAY THE LARGEST CONNECTED
! COMPONENT FELL. =+1 IF LOW AND -1 IF HIGH
! COMMON /GRA/ N, IDPTH, IDEG
! IT IS ASSUMED THAT THE GRAPH HAS AT MOST 50 COMPONENTS AND
! THAT THERE ARE AT MOST 100 LEVELS.
! COMMON /LVLW/ NHIGH(100), NLOW(100), NACUM(100)
! COMMON /CC/ XCC, SIZEG(50), STPT(50)
INTEGER(psb_lpk_) :: LVLS1(N), LVLS2(N), CCSTOR(N)
integer(psb_lpk_) :: SZ, ENDC,i,j,max1,max2,inode
integer(psb_lpk_) :: lvlnh, it, k, lvlnl,idflt,isdir
! FOR EACH CONNECTED COMPONENT DO
DO 80 I=1,XCC
J = STPT(I)
ENDC= SIZEG(I) + J - 1
! SET NHIGH AND NLOW EQUAL TO NACUM
!-----------------------------------------------------
SZ=SIZE(NHIGH)
IF(SZ < IDPTH) THEN
write(psb_out_unit,*) 'GPS_PIKLVL: on fly reallocation of NHIGH'
CALL REALLOC(NHIGH,SZ,IDPTH)
END IF
SZ=SIZE(NLOW)
IF(SZ < IDPTH) THEN
write(psb_out_unit,*) 'GPS_PIKLVL: on fly reallocation of NLOW'
CALL REALLOC(NLOW,SZ,IDPTH)
END IF
!-----------------------------------------------------
DO 10 K=1,IDPTH
NHIGH(K) = NACUM(K)
NLOW(K) = NACUM(K)
10 END DO
! UPDATE NHIGH AND NLOW FOR EACH NODE IN CONNECTED COMPONENT
DO 20 K=J,ENDC
INODE = CCSTOR(K)
LVLNH = LVLS1(INODE)
NHIGH(LVLNH) = NHIGH(LVLNH) + 1
LVLNL = LVLS2(INODE)
NLOW(LVLNL) = NLOW(LVLNL) + 1
20 END DO
MAX1 = 0
MAX2 = 0
! SET MAX1=LARGEST NEW NUMBER IN NHIGH
! SET MAX2=LARGEST NEW NUMBER IN NLOW
DO 30 K=1,IDPTH
IF (2*NACUM(K) == NLOW(K)+NHIGH(K)) CYCLE
IF (NHIGH(K) > MAX1) MAX1 = NHIGH(K)
IF (NLOW(K) > MAX2) MAX2 = NLOW(K)
30 END DO
! SET IT= NUMBER OF LEVEL STRUCTURE TO BE USED
IT = 1
IF (MAX1 > MAX2) IT = 2
IF (MAX1 == MAX2) IT = IDFLT
IF (IT == 2) GO TO 60
IF (I == 1) ISDIR = -1
! COPY LVLS1 INTO LVLS2 FOR EACH NODE IN CONNECTED COMPONENT
DO 40 K=J,ENDC
INODE = CCSTOR(K)
LVLS2(INODE) = LVLS1(INODE)
40 END DO
! UPDATE NACUM TO BE THE SAME AS NHIGH
DO 50 K=1,IDPTH
NACUM(K) = NHIGH(K)
50 END DO
CYCLE
! UPDATE NACUM TO BE THE SAME AS NLOW
60 DO 70 K=1,IDPTH
NACUM(K) = NLOW(K)
70 END DO
80 END DO
RETURN
END SUBROUTINE PIKLVL
!
SUBROUTINE NUMBER(SND, NUM, NDSTK, LVLS2, NDEG, RENUM, LVLST,LSTPT,&
& NR, NFLG, IBW2, IPF2, IPFA, ISDIR)
use psb_base_mod
implicit none
! NUMBER PRODUCES THE NUMBERING OF THE GRAPH FOR MIN BANDWIDTH
! SND- ON INPUT THE NODE TO BEGIN NUMBERING ON
! NUM- ON INPUT AND OUTPUT, THE NEXT AVAILABLE NUMBER
! LVLS2- THE LEVEL STRUCTURE TO BE USED IN NUMBERING
! RENUM- THE ARRAY USED TO STORE THE NEW NUMBERING
! LVLST- ON OUTPUT CONTAINS LEVEL STRUCTURE
! LSTPT(I)- ON OUTPUT, INDEX INTO LVLST TO FIRST NODE IN ITH LVL
! LSTPT(I+1) - LSTPT(I) = NUMBER OF NODES IN ITH LVL
! NFLG- =+1 IF SND IS FORWARD END OF PSEUDO-DIAM
! =-1 IF SND IS REVERSE END OF PSEUDO-DIAM
! IBW2- BANDWIDTH OF NEW NUMBERING COMPUTED BY NUMBER
! IPF2- PROFILE OF NEW NUMBERING COMPUTED BY NUMBER
! IPFA- WORKING STORAGE USED TO COMPUTE PROFILE AND BANDWIDTH
! ISDIR- INDICATES STEP DIRECTION USED IN NUMBERING(+1 OR -1)
! USE INTEGER*2 NDSTK WITH AN IBM 360 OR 370.
INTEGER(psb_lpk_) :: SND, NUM, XA, XB, XC, XD, CX, ENDC, TEST, NR, ISDIR
! COMMON /GRA/ N, IDPTH, IDEG
! THE STORAGE IN COMMON BLOCKS CC AND LVLW IS NOW FREE AND CAN
! BE USED FOR STACKS.
!COMMON /LVLW/ STKA(100), STKB(100), STKC(100)
!COMMON /CC/ STKD(100)
INTEGER(psb_lpk_) :: IPFA(N), NDSTK(NR,IDEG), LVLS2(N),&
& NDEG(N), RENUM(N+1), LVLST(N),LSTPT(N),ipf2,ibw2,nflg, nbw
integer(psb_lpk_),POINTER :: STKA(:),STKB(:),STKC(:),STKD(:)
integer(psb_lpk_) :: SZ1,SZ2,i,j,nstpt,lvln, lst, lnd, inx, max, ipro,&
& lvlnl, k, it
!
STKA => NHIGH
STKB => NLOW
STKC => NACUM
STKD => AUX
!
! SET UP LVLST AND LSTPT FROM LVLS2
DO 10 I=1,N
IPFA(I) = 0
10 END DO
NSTPT = 1
DO 30 I=1,IDPTH
LSTPT(I) = NSTPT
DO 20 J=1,N
IF (LVLS2(J) /= I) CYCLE
LVLST(NSTPT) = J
NSTPT = NSTPT + 1
20 END DO
30 END DO
LSTPT(IDPTH+1) = NSTPT
! STKA, STKB, STKC AND STKD ARE STACKS WITH POINTERS
! XA,XB,XC, AND XD. CX IS A SPECIAL POINTER INTO STKC WHICH
! INDICATES THE PARTICULAR NODE BEING PROCESSED.
! LVLN KEEPS TRACK OF THE LEVEL WE ARE WORKING AT.
! INITIALLY STKC CONTAINS ONLY THE INITIAL NODE, SND.
LVLN = 0
IF (NFLG < 0) LVLN = IDPTH + 1
XC = 1
STKC(XC) = SND
40 CX = 1
XD = 0
LVLN = LVLN + NFLG
LST = LSTPT(LVLN)
LND = LSTPT(LVLN+1) - 1
! BEGIN PROCESSING NODE STKC(CX)
50 IPRO = STKC(CX)
RENUM(IPRO) = NUM
NUM = NUM + ISDIR
ENDC = NDEG(IPRO)
XA = 0
XB = 0
! CHECK ALL ADJACENT NODES
DO 80 I=1,ENDC
TEST = NDSTK(IPRO,I)
INX = RENUM(TEST)
! ONLY NODES NOT NUMBERED OR ALREADY ON A STACK ARE ADDED
IF (INX == 0) GO TO 60
IF (INX < 0) CYCLE
! DO PRELIMINARY BANDWIDTH AND PROFILE CALCULATIONS
NBW = (RENUM(IPRO)-INX)*ISDIR
IF (ISDIR > 0) INX = RENUM(IPRO)
IF (IPFA(INX) < NBW) IPFA(INX) = NBW
CYCLE
60 RENUM(TEST) = -1
! PUT NODES ON SAME LEVEL ON STKA, ALL OTHERS ON STKB
IF (LVLS2(TEST) == LVLS2(IPRO)) GO TO 70
XB = XB + 1
STKB(XB) = TEST
CYCLE
70 XA = XA + 1
STKA(XA) = TEST
80 END DO
! SORT STKA AND STKB INTO INCREASING DEGREE AND ADD STKA TO STKC
! AND STKB TO STKD
IF (XA == 0) GO TO 100
IF (XA == 1) GO TO 90
!-----------------------------------------------------------------
SZ1=SIZE(STKC)
SZ2=XC+XA
IF(SZ1 < SZ2) THEN
write(psb_out_unit,*) 'GPS_NUMBER - Check #1: on fly reallocation of STKC'
CALL REALLOC(NACUM,SZ1,SZ2)
STKC => NACUM
END IF
!-----------------------------------------------------------------
CALL SORTDG(STKC, STKA, XC, XA, NDEG)
GO TO 100
90 XC = XC + 1
!-----------------------------------------------------------------
SZ1=SIZE(STKC)
SZ2=XC
IF(SZ1 < SZ2) THEN
write(psb_out_unit,*) 'GPS_NUMBER - Check #2: on fly reallocation of STKC'
SZ2=SZ2+INIT
CALL REALLOC(NACUM,SZ1,SZ2)
STKC => NACUM
END IF
!-----------------------------------------------------------------
STKC(XC) = STKA(XA)
100 IF (XB == 0) GO TO 120
IF (XB == 1) GO TO 110
!-----------------------------------------------------------------
SZ1=SIZE(STKD)
SZ2=XD+XB
IF(SZ1 < SZ2) THEN
write(psb_out_unit,*) 'GPS_NUMBER - Check #3: on fly reallocation of STKD'
CALL REALLOC(AUX,SZ1,SZ2)
STKD => AUX
END IF
!-----------------------------------------------------------------
CALL SORTDG(STKD, STKB, XD, XB, NDEG)
GO TO 120
110 XD = XD + 1
!-----------------------------------------------------------------
SZ1=SIZE(STKD)
SZ2=XD
IF(SZ1 < SZ2) THEN
write(psb_out_unit,*) 'GPS_NUMBER - Check #4: on fly reallocation of STKD'
SZ2=SZ2+INIT
CALL REALLOC(AUX,SZ1,SZ2)
STKD => AUX
END IF
!-----------------------------------------------------------------
STKD(XD) = STKB(XB)
! BE SURE TO PROCESS ALL NODES IN STKC
120 CX = CX + 1
IF (XC >= CX) GO TO 50
! WHEN STKC IS EXHAUSTED LOOK FOR MIN DEGREE NODE IN SAME LEVEL
! WHICH HAS NOT BEEN PROCESSED
MAX = IDEG + 1
SND = N + 1
DO 130 I=LST,LND
TEST = LVLST(I)
IF (RENUM(TEST) /= 0) CYCLE
IF (NDEG(TEST) >= MAX) CYCLE
RENUM(SND) = 0
RENUM(TEST) = -1
MAX = NDEG(TEST)
SND = TEST
130 END DO
IF (SND == N+1) GO TO 140
XC = XC + 1
!-----------------------------------------------------------------
SZ1=SIZE(STKC)
SZ2=XC
IF(SZ1 < SZ2) THEN
write(psb_out_unit,*) 'GPS_NUMBER - Check #5: on fly reallocation of STKC'
SZ2=SZ2+INIT
CALL REALLOC(NACUM,SZ1,SZ2)
STKC => NACUM
END IF
!-----------------------------------------------------------------
STKC(XC) = SND
GO TO 50
! IF STKD IS EMPTY WE ARE DONE, OTHERWISE COPY STKD ONTO STKC
! AND BEGIN PROCESSING NEW STKC
140 IF (XD == 0) GO TO 160
!-----------------------------------------------------------------
SZ1=SIZE(STKC)
SZ2=XD
IF(SZ1 < SZ2) THEN
write(psb_out_unit,*) 'GPS_NUMBER - Check #6: on fly reallocation of STKC'
SZ2=SZ2+INIT
CALL REALLOC(NACUM,SZ1,SZ2)
STKC => NACUM
END IF
!-----------------------------------------------------------------
DO 150 I=1,XD
STKC(I) = STKD(I)
150 END DO
XC = XD
GO TO 40
! DO FINAL BANDWIDTH AND PROFILE CALCULATIONS
160 DO 170 I=1,N
IF (IPFA(I) > IBW2) IBW2 = IPFA(I)
IPF2 = IPF2 + IPFA(I)
170 END DO
!
RETURN
END SUBROUTINE NUMBER
!
! ---------------------------------------------------------
SUBROUTINE REALLOC(VET,SZ1,SZ2)
use psb_base_mod
! PERFORM ON FLY REALLOCATION OF POINTER VET INCREASING
! ITS SIZE FROM SZ1 TO SZ2
IMPLICIT NONE
integer(psb_lpk_),allocatable :: VET(:)
integer(psb_lpk_) :: SZ1,SZ2
integer(psb_ipk_) :: info
call psb_realloc(sz2,vet,info)
IF(INFO /= psb_success_) THEN
write(psb_out_unit,*) 'Error! Memory allocation failure in REALLOC'
STOP
END IF
RETURN
END SUBROUTINE REALLOC
!
end module psb_lgps_mod

101
util/psb_renum_mod.f90
Unescape Escape View File

@ -32,101 +32,8 @@
module psb_renum_mod
use psb_base_mod
integer(psb_ipk_), parameter :: psb_mat_renum_identity_ = 0
integer(psb_ipk_), parameter :: psb_mat_renum_gps_ = 456
integer(psb_ipk_), parameter :: psb_mat_renum_amd_ = psb_mat_renum_gps_ + 1
interface psb_mat_renum
subroutine psb_d_mat_renums(alg,mat,info,perm)
import :: psb_ipk_, psb_dspmat_type
character(len=*), intent(in) :: alg
type(psb_dspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
end subroutine psb_d_mat_renums
subroutine psb_d_mat_renum(alg,mat,info,perm)
import :: psb_ipk_, psb_dspmat_type
integer(psb_ipk_), intent(in) :: alg
type(psb_dspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
end subroutine psb_d_mat_renum
subroutine psb_s_mat_renums(alg,mat,info,perm)
import :: psb_ipk_, psb_sspmat_type
character(len=*), intent(in) :: alg
type(psb_sspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
end subroutine psb_s_mat_renums
subroutine psb_s_mat_renum(alg,mat,info,perm)
import :: psb_ipk_, psb_sspmat_type
integer(psb_ipk_), intent(in) :: alg
type(psb_sspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
end subroutine psb_s_mat_renum
subroutine psb_z_mat_renums(alg,mat,info,perm)
import :: psb_ipk_, psb_zspmat_type
character(len=*), intent(in) :: alg
type(psb_zspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
end subroutine psb_z_mat_renums
subroutine psb_z_mat_renum(alg,mat,info,perm)
import :: psb_ipk_, psb_zspmat_type
integer(psb_ipk_), intent(in) :: alg
type(psb_zspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
end subroutine psb_z_mat_renum
subroutine psb_c_mat_renums(alg,mat,info,perm)
import :: psb_ipk_, psb_cspmat_type
character(len=*), intent(in) :: alg
type(psb_cspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
end subroutine psb_c_mat_renums
subroutine psb_c_mat_renum(alg,mat,info,perm)
import :: psb_ipk_, psb_cspmat_type
integer(psb_ipk_), intent(in) :: alg
type(psb_cspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
end subroutine psb_c_mat_renum
end interface psb_mat_renum
interface psb_cmp_bwpf
subroutine psb_s_cmp_bwpf(mat,bwl,bwu,prf,info)
import :: psb_ipk_, psb_sspmat_type
type(psb_sspmat_type), intent(in) :: mat
integer(psb_ipk_), intent(out) :: bwl, bwu
integer(psb_ipk_), intent(out) :: prf
integer(psb_ipk_), intent(out) :: info
end subroutine psb_s_cmp_bwpf
subroutine psb_d_cmp_bwpf(mat,bwl,bwu,prf,info)
import :: psb_ipk_, psb_dspmat_type
type(psb_dspmat_type), intent(in) :: mat
integer(psb_ipk_), intent(out) :: bwl, bwu
integer(psb_ipk_), intent(out) :: prf
integer(psb_ipk_), intent(out) :: info
end subroutine psb_d_cmp_bwpf
subroutine psb_c_cmp_bwpf(mat,bwl,bwu,prf,info)
import :: psb_ipk_, psb_cspmat_type
type(psb_cspmat_type), intent(in) :: mat
integer(psb_ipk_), intent(out) :: bwl, bwu
integer(psb_ipk_), intent(out) :: prf
integer(psb_ipk_), intent(out) :: info
end subroutine psb_c_cmp_bwpf
subroutine psb_z_cmp_bwpf(mat,bwl,bwu,prf,info)
import :: psb_ipk_, psb_zspmat_type
type(psb_zspmat_type), intent(in) :: mat
integer(psb_ipk_), intent(out) :: bwl, bwu
integer(psb_ipk_), intent(out) :: prf
integer(psb_ipk_), intent(out) :: info
end subroutine psb_z_cmp_bwpf
end interface psb_cmp_bwpf
use psb_s_renum_mod
use psb_c_renum_mod
use psb_d_renum_mod
use psb_z_renum_mod
end module psb_renum_mod

369
util/psb_s_renum_impl.F90
Unescape Escape View File

@ -29,37 +29,61 @@
! POSSIBILITY OF SUCH DAMAGE.
!
!
subroutine psb_s_mat_renums(alg,mat,info,perm)
subroutine psb_s_mat_renum(alg,mat,info,perm)
use psb_base_mod
use psb_renum_mod, psb_protect_name => psb_s_mat_renums
use psb_renum_mod, psb_protect_name => psb_s_mat_renum
implicit none
character(len=*), intent(in) :: alg
type(psb_sspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
integer(psb_ipk_) :: err_act, nr, nc, ialg
integer(psb_ipk_) :: err_act, nr, nc, i, ierr(5)
character(len=20) :: name
info = psb_success_
name = 'mat_renum'
call psb_erractionsave(err_act)
nr = mat%get_nrows()
nc = mat%get_ncols()
if (nr /= nc) then
info = psb_err_rectangular_mat_unsupported_
ierr(1) = nr; ierr(2) = nc;
call psb_errpush(info,name,i_err=ierr)
goto 9999
end if
info = psb_success_
select case (psb_toupper(alg))
case ('GPS')
ialg = psb_mat_renum_gps_
case ('AMD')
ialg = psb_mat_renum_amd_
call psb_mat_renum_gps(mat,info,perm)
case('AMD')
call psb_mat_renum_amd(mat,info,perm)
case ('NONE', 'ID')
ialg = psb_mat_renum_identity_
nr = mat%get_nrows()
allocate(perm(nr),stat=info)
if (info == 0) then
do i=1,nr
perm(i) = i
end do
else
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
endif
case default
write(0,*) 'Unknown algorithm "',psb_toupper(alg),'"'
ialg = -1
info = psb_err_input_value_invalid_i_
ierr(1) = 1;
call psb_errpush(info,name,i_err=ierr)
goto 9999
end select
call psb_mat_renum(ialg,mat,info,perm)
if (info /= psb_success_) then
info = psb_err_from_subroutine_non_
call psb_errpush(info,name)
@ -72,26 +96,218 @@ subroutine psb_s_mat_renums(alg,mat,info,perm)
9999 call psb_error_handler(err_act)
return
end subroutine psb_s_mat_renums
subroutine psb_s_mat_renum(alg,mat,info,perm)
contains
subroutine psb_mat_renum_gps(a,info,operm)
use psb_base_mod
use psb_gps_mod
implicit none
type(psb_sspmat_type), intent(inout) :: a
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: operm(:)
!
class(psb_s_base_sparse_mat), allocatable :: aa
type(psb_s_csr_sparse_mat) :: acsr
type(psb_s_coo_sparse_mat) :: acoo
integer(psb_ipk_) :: err_act
character(len=20) :: name
integer(psb_ipk_), allocatable :: ndstk(:,:), iold(:), ndeg(:), perm(:)
integer(psb_ipk_) :: i, j, k, ideg, nr, ibw, ipf, idpth
info = psb_success_
name = 'mat_renum_gps'
call psb_erractionsave(err_act)
info = psb_success_
call a%mold(aa)
call a%mv_to(aa)
call aa%mv_to_fmt(acsr,info)
! Insert call to gps_reduce
nr = acsr%get_nrows()
ideg = 0
do i=1, nr
ideg = max(ideg,acsr%irp(i+1)-acsr%irp(i))
end do
allocate(ndstk(nr,ideg), iold(nr), perm(nr+1), ndeg(nr),stat=info)
if (info /= 0) then
info = psb_err_alloc_dealloc_
call psb_errpush(info, name)
goto 9999
end if
do i=1, nr
iold(i) = i
ndstk(i,:) = 0
k = 0
do j=acsr%irp(i),acsr%irp(i+1)-1
k = k + 1
ndstk(i,k) = acsr%ja(j)
end do
end do
perm = 0
call psb_gps_reduce(ndstk,nr,ideg,iold,perm,ndeg,ibw,ipf,idpth)
if (.not.psb_isaperm(nr,perm)) then
write(0,*) 'Something wrong: bad perm from gps_reduce'
info = psb_err_from_subroutine_
call psb_errpush(info,name)
goto 9999
end if
! Move to coordinate to apply renumbering
call acsr%mv_to_coo(acoo,info)
do i=1, acoo%get_nzeros()
acoo%ia(i) = perm(acoo%ia(i))
acoo%ja(i) = perm(acoo%ja(i))
end do
call acoo%fix(info)
! Get back to where we started from
call aa%mv_from_coo(acoo,info)
call a%mv_from(aa)
if (present(operm)) then
call psb_realloc(nr,operm,info)
if (info /= psb_success_) then
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
end if
operm(1:nr) = perm(1:nr)
end if
deallocate(aa)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_mat_renum_gps
subroutine psb_mat_renum_amd(a,info,operm)
#if defined(HAVE_AMD)
use iso_c_binding
#endif
use psb_base_mod
implicit none
type(psb_sspmat_type), intent(inout) :: a
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: operm(:)
!
#if defined(HAVE_AMD)
interface
function psb_amd_order(n,ap,ai,p)&
& result(res) bind(c,name='psb_amd_order')
use iso_c_binding
integer(c_int) :: res
integer(c_int), value :: n
integer(c_int) :: ap(*), ai(*), p(*)
end function psb_amd_order
end interface
#endif
type(psb_s_csc_sparse_mat) :: acsc
class(psb_s_base_sparse_mat), allocatable :: aa
type(psb_s_coo_sparse_mat) :: acoo
integer(psb_ipk_), allocatable :: perm(:)
integer(psb_ipk_) :: err_act
character(len=20) :: name
integer(psb_ipk_) :: i, j, k, ideg, nr, ibw, ipf, idpth, nz
info = psb_success_
name = 'mat_renum_amd'
call psb_erractionsave(err_act)
#if defined(HAVE_AMD) && defined(IPK4)
info = psb_success_
nr = a%get_nrows()
nz = a%get_nzeros()
allocate(perm(nr),stat=info)
if (info /= 0) then
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
end if
allocate(aa, mold=a%a)
call a%mv_to(acsc)
acsc%ia(:) = acsc%ia(:) - 1
acsc%icp(:) = acsc%icp(:) - 1
info = psb_amd_order(nr,acsc%icp,acsc%ia,perm)
if (info /= psb_success_) then
info = psb_err_from_subroutine_
call psb_errpush(info,name,a_err='psb_amd_order')
goto 9999
end if
perm(:) = perm(:) + 1
acsc%ia(:) = acsc%ia(:) + 1
acsc%icp(:) = acsc%icp(:) + 1
call acsc%mv_to_coo(acoo,info)
do i=1, acoo%get_nzeros()
acoo%ia(i) = perm(acoo%ia(i))
acoo%ja(i) = perm(acoo%ja(i))
end do
call acoo%fix(info)
! Get back to where we started from
call aa%mv_from_coo(acoo,info)
call a%mv_from(aa)
if (present(operm)) then
call psb_realloc(nr,operm,info)
if (info /= psb_success_) then
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
end if
operm(1:nr) = perm(1:nr)
end if
deallocate(aa,perm)
#else
info = psb_err_missing_aux_lib_
call psb_errpush(info, name)
goto 9999
#endif
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_mat_renum_amd
end subroutine psb_s_mat_renum
subroutine psb_ls_mat_renum(alg,mat,info,perm)
use psb_base_mod
use psb_renum_mod, psb_protect_name => psb_s_mat_renum
use psb_renum_mod, psb_protect_name => psb_ls_mat_renum
implicit none
integer(psb_ipk_), intent(in) :: alg
type(psb_sspmat_type), intent(inout) :: mat
character(len=*), intent(in) :: alg
type(psb_lsspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
integer(psb_lpk_), allocatable, optional, intent(out) :: perm(:)
integer(psb_ipk_) :: err_act, nr, nc, i, ierr(5)
integer(psb_lpk_) :: nr, nc, i
integer(psb_ipk_) :: err_act, ierr(5)
character(len=20) :: name
info = psb_success_
name = 'mat_renum'
call psb_erractionsave(err_act)
info = psb_success_
nr = mat%get_nrows()
nc = mat%get_ncols()
if (nr /= nc) then
@ -101,16 +317,17 @@ subroutine psb_s_mat_renum(alg,mat,info,perm)
goto 9999
end if
select case (alg)
case(psb_mat_renum_gps_)
info = psb_success_
select case (psb_toupper(alg))
case ('GPS')
call psb_mat_renum_gps(mat,info,perm)
call psb_lmat_renum_gps(mat,info,perm)
case(psb_mat_renum_amd_)
case('AMD')
call psb_mat_renum_amd(mat,info,perm)
call psb_lmat_renum_amd(mat,info,perm)
case(psb_mat_renum_identity_)
case ('NONE', 'ID')
nr = mat%get_nrows()
allocate(perm(nr),stat=info)
if (info == 0) then
@ -123,8 +340,9 @@ subroutine psb_s_mat_renum(alg,mat,info,perm)
goto 9999
endif
case default
write(0,*) 'Unknown algorithm "',psb_toupper(alg),'"'
info = psb_err_input_value_invalid_i_
ierr(1) = 1; ierr(2) = alg;
ierr(1) = 1;
call psb_errpush(info,name,i_err=ierr)
goto 9999
end select
@ -143,26 +361,26 @@ subroutine psb_s_mat_renum(alg,mat,info,perm)
contains
subroutine psb_mat_renum_gps(a,info,operm)
subroutine psb_lmat_renum_gps(a,info,operm)
use psb_base_mod
use psb_gps_mod
use psb_lgps_mod
implicit none
type(psb_sspmat_type), intent(inout) :: a
type(psb_lsspmat_type), intent(inout) :: a
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: operm(:)
integer(psb_lpk_), allocatable, optional, intent(out) :: operm(:)
!
class(psb_s_base_sparse_mat), allocatable :: aa
type(psb_s_csr_sparse_mat) :: acsr
type(psb_s_coo_sparse_mat) :: acoo
class(psb_ls_base_sparse_mat), allocatable :: aa
type(psb_ls_csr_sparse_mat) :: acsr
type(psb_ls_coo_sparse_mat) :: acoo
integer(psb_ipk_) :: err_act
character(len=20) :: name
integer(psb_ipk_), allocatable :: ndstk(:,:), iold(:), ndeg(:), perm(:)
integer(psb_ipk_) :: i, j, k, ideg, nr, ibw, ipf, idpth
integer(psb_lpk_), allocatable :: ndstk(:,:), iold(:), ndeg(:), perm(:)
integer(psb_lpk_) :: i, j, k, ideg, nr, ibw, ipf, idpth
info = psb_success_
name = 'mat_renum'
name = 'mat_renum_gps'
call psb_erractionsave(err_act)
info = psb_success_
@ -193,7 +411,7 @@ contains
end do
perm = 0
call psb_gps_reduce(ndstk,nr,ideg,iold,perm,ndeg,ibw,ipf,idpth)
call psb_lgps_reduce(ndstk,nr,ideg,iold,perm,ndeg,ibw,ipf,idpth)
if (.not.psb_isaperm(nr,perm)) then
write(0,*) 'Something wrong: bad perm from gps_reduce'
@ -229,18 +447,18 @@ contains
9999 call psb_error_handler(err_act)
return
end subroutine psb_mat_renum_gps
end subroutine psb_lmat_renum_gps
subroutine psb_mat_renum_amd(a,info,operm)
subroutine psb_lmat_renum_amd(a,info,operm)
#if defined(HAVE_AMD)
use iso_c_binding
#endif
use psb_base_mod
implicit none
type(psb_sspmat_type), intent(inout) :: a
type(psb_lsspmat_type), intent(inout) :: a
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: operm(:)
integer(psb_lpk_), allocatable, optional, intent(out) :: operm(:)
!
#if defined(HAVE_AMD)
@ -255,20 +473,20 @@ contains
end interface
#endif
type(psb_s_csc_sparse_mat) :: acsc
class(psb_s_base_sparse_mat), allocatable :: aa
type(psb_s_coo_sparse_mat) :: acoo
type(psb_ls_csc_sparse_mat) :: acsc
class(psb_ls_base_sparse_mat), allocatable :: aa
type(psb_ls_coo_sparse_mat) :: acoo
integer(psb_ipk_), allocatable :: perm(:)
integer(psb_ipk_) :: err_act
character(len=20) :: name
integer(psb_ipk_) :: i, j, k, ideg, nr, ibw, ipf, idpth, nz
integer(psb_lpk_) :: i, j, k, ideg, nr, ibw, ipf, idpth, nz
info = psb_success_
name = 'mat_renum_amd'
call psb_erractionsave(err_act)
#if defined(HAVE_AMD) && defined(IPK4)
#if defined(HAVE_AMD) && defined(LPK4)
info = psb_success_
nr = a%get_nrows()
@ -332,11 +550,9 @@ contains
9999 call psb_error_handler(err_act)
return
end subroutine psb_lmat_renum_amd
end subroutine psb_mat_renum_amd
end subroutine psb_s_mat_renum
end subroutine psb_ls_mat_renum
subroutine psb_s_cmp_bwpf(mat,bwl,bwu,prf,info)
use psb_base_mod
@ -387,3 +603,52 @@ subroutine psb_s_cmp_bwpf(mat,bwl,bwu,prf,info)
end select
end subroutine psb_s_cmp_bwpf
subroutine psb_ls_cmp_bwpf(mat,bwl,bwu,prf,info)
use psb_base_mod
use psb_renum_mod, psb_protect_name => psb_ls_cmp_bwpf
implicit none
type(psb_lsspmat_type), intent(in) :: mat
integer(psb_lpk_), intent(out) :: bwl, bwu
integer(psb_lpk_), intent(out) :: prf
integer(psb_ipk_), intent(out) :: info
!
integer(psb_lpk_), allocatable :: irow(:), icol(:)
real(psb_spk_), allocatable :: val(:)
integer(psb_lpk_) :: nz, i, j, lrbu, lrbl
info = psb_success_
bwl = 0
bwu = 0
prf = 0
select type (aa=>mat%a)
class is (psb_ls_csr_sparse_mat)
do i=1, aa%get_nrows()
lrbl = 0
lrbu = 0
do j = aa%irp(i), aa%irp(i+1) - 1
lrbl = max(lrbl,i-aa%ja(j))
lrbu = max(lrbu,aa%ja(j)-i)
end do
prf = prf + lrbl+lrbu
bwu = max(bwu,lrbu)
bwl = max(bwl,lrbu)
end do
class default
do i=1, aa%get_nrows()
lrbl = 0
lrbu = 0
call aa%csget(i,i,nz,irow,icol,val,info)
if (info /= psb_success_) return
do j=1, nz
lrbl = max(lrbl,i-icol(j))
lrbu = max(lrbu,icol(j)-i)
end do
prf = prf + lrbl+lrbu
bwu = max(bwu,lrbu)
bwl = max(bwl,lrbu)
end do
end select
end subroutine psb_ls_cmp_bwpf

69
util/psb_s_renum_mod.f90
Unescape Escape View File

@ -0,0 +1,69 @@
!
! Parallel Sparse BLAS version 3.5
! (C) Copyright 2006-2018
! Salvatore Filippone
! Alfredo Buttari
!
! 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.
!
!
module psb_s_renum_mod
use psb_base_mod
interface psb_mat_renum
subroutine psb_s_mat_renum(alg,mat,info,perm)
import :: psb_ipk_, psb_sspmat_type
character(len=*), intent(in) :: alg
type(psb_sspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
end subroutine psb_s_mat_renum
subroutine psb_ls_mat_renum(alg,mat,info,perm)
import :: psb_ipk_, psb_lpk_, psb_lsspmat_type
character(len=*), intent(in) :: alg
type(psb_lsspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_lpk_), allocatable, optional, intent(out) :: perm(:)
end subroutine psb_ls_mat_renum
end interface psb_mat_renum
interface psb_cmp_bwpf
subroutine psb_s_cmp_bwpf(mat,bwl,bwu,prf,info)
import :: psb_ipk_, psb_sspmat_type
type(psb_sspmat_type), intent(in) :: mat
integer(psb_ipk_), intent(out) :: bwl, bwu
integer(psb_ipk_), intent(out) :: prf
integer(psb_ipk_), intent(out) :: info
end subroutine psb_s_cmp_bwpf
subroutine psb_ls_cmp_bwpf(mat,bwl,bwu,prf,info)
import :: psb_ipk_, psb_lpk_, psb_lsspmat_type
type(psb_lsspmat_type), intent(in) :: mat
integer(psb_lpk_), intent(out) :: bwl, bwu
integer(psb_lpk_), intent(out) :: prf
integer(psb_ipk_), intent(out) :: info
end subroutine psb_ls_cmp_bwpf
end interface psb_cmp_bwpf
end module psb_s_renum_mod

369
util/psb_z_renum_impl.F90
Unescape Escape View File

@ -29,37 +29,61 @@
! POSSIBILITY OF SUCH DAMAGE.
!
!
subroutine psb_z_mat_renums(alg,mat,info,perm)
subroutine psb_z_mat_renum(alg,mat,info,perm)
use psb_base_mod
use psb_renum_mod, psb_protect_name => psb_z_mat_renums
use psb_renum_mod, psb_protect_name => psb_z_mat_renum
implicit none
character(len=*), intent(in) :: alg
type(psb_zspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
integer(psb_ipk_) :: err_act, nr, nc, ialg
integer(psb_ipk_) :: err_act, nr, nc, i, ierr(5)
character(len=20) :: name
info = psb_success_
name = 'mat_renum'
call psb_erractionsave(err_act)
nr = mat%get_nrows()
nc = mat%get_ncols()
if (nr /= nc) then
info = psb_err_rectangular_mat_unsupported_
ierr(1) = nr; ierr(2) = nc;
call psb_errpush(info,name,i_err=ierr)
goto 9999
end if
info = psb_success_
select case (psb_toupper(alg))
case ('GPS')
ialg = psb_mat_renum_gps_
case ('AMD')
ialg = psb_mat_renum_amd_
call psb_mat_renum_gps(mat,info,perm)
case('AMD')
call psb_mat_renum_amd(mat,info,perm)
case ('NONE', 'ID')
ialg = psb_mat_renum_identity_
nr = mat%get_nrows()
allocate(perm(nr),stat=info)
if (info == 0) then
do i=1,nr
perm(i) = i
end do
else
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
endif
case default
write(0,*) 'Unknown algorithm "',psb_toupper(alg),'"'
ialg = -1
info = psb_err_input_value_invalid_i_
ierr(1) = 1;
call psb_errpush(info,name,i_err=ierr)
goto 9999
end select
call psb_mat_renum(ialg,mat,info,perm)
if (info /= psb_success_) then
info = psb_err_from_subroutine_non_
call psb_errpush(info,name)
@ -72,26 +96,218 @@ subroutine psb_z_mat_renums(alg,mat,info,perm)
9999 call psb_error_handler(err_act)
return
end subroutine psb_z_mat_renums
subroutine psb_z_mat_renum(alg,mat,info,perm)
contains
subroutine psb_mat_renum_gps(a,info,operm)
use psb_base_mod
use psb_gps_mod
implicit none
type(psb_zspmat_type), intent(inout) :: a
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: operm(:)
!
class(psb_z_base_sparse_mat), allocatable :: aa
type(psb_z_csr_sparse_mat) :: acsr
type(psb_z_coo_sparse_mat) :: acoo
integer(psb_ipk_) :: err_act
character(len=20) :: name
integer(psb_ipk_), allocatable :: ndstk(:,:), iold(:), ndeg(:), perm(:)
integer(psb_ipk_) :: i, j, k, ideg, nr, ibw, ipf, idpth
info = psb_success_
name = 'mat_renum_gps'
call psb_erractionsave(err_act)
info = psb_success_
call a%mold(aa)
call a%mv_to(aa)
call aa%mv_to_fmt(acsr,info)
! Insert call to gps_reduce
nr = acsr%get_nrows()
ideg = 0
do i=1, nr
ideg = max(ideg,acsr%irp(i+1)-acsr%irp(i))
end do
allocate(ndstk(nr,ideg), iold(nr), perm(nr+1), ndeg(nr),stat=info)
if (info /= 0) then
info = psb_err_alloc_dealloc_
call psb_errpush(info, name)
goto 9999
end if
do i=1, nr
iold(i) = i
ndstk(i,:) = 0
k = 0
do j=acsr%irp(i),acsr%irp(i+1)-1
k = k + 1
ndstk(i,k) = acsr%ja(j)
end do
end do
perm = 0
call psb_gps_reduce(ndstk,nr,ideg,iold,perm,ndeg,ibw,ipf,idpth)
if (.not.psb_isaperm(nr,perm)) then
write(0,*) 'Something wrong: bad perm from gps_reduce'
info = psb_err_from_subroutine_
call psb_errpush(info,name)
goto 9999
end if
! Move to coordinate to apply renumbering
call acsr%mv_to_coo(acoo,info)
do i=1, acoo%get_nzeros()
acoo%ia(i) = perm(acoo%ia(i))
acoo%ja(i) = perm(acoo%ja(i))
end do
call acoo%fix(info)
! Get back to where we started from
call aa%mv_from_coo(acoo,info)
call a%mv_from(aa)
if (present(operm)) then
call psb_realloc(nr,operm,info)
if (info /= psb_success_) then
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
end if
operm(1:nr) = perm(1:nr)
end if
deallocate(aa)
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_mat_renum_gps
subroutine psb_mat_renum_amd(a,info,operm)
#if defined(HAVE_AMD)
use iso_c_binding
#endif
use psb_base_mod
implicit none
type(psb_zspmat_type), intent(inout) :: a
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: operm(:)
!
#if defined(HAVE_AMD)
interface
function psb_amd_order(n,ap,ai,p)&
& result(res) bind(c,name='psb_amd_order')
use iso_c_binding
integer(c_int) :: res
integer(c_int), value :: n
integer(c_int) :: ap(*), ai(*), p(*)
end function psb_amd_order
end interface
#endif
type(psb_z_csc_sparse_mat) :: acsc
class(psb_z_base_sparse_mat), allocatable :: aa
type(psb_z_coo_sparse_mat) :: acoo
integer(psb_ipk_), allocatable :: perm(:)
integer(psb_ipk_) :: err_act
character(len=20) :: name
integer(psb_ipk_) :: i, j, k, ideg, nr, ibw, ipf, idpth, nz
info = psb_success_
name = 'mat_renum_amd'
call psb_erractionsave(err_act)
#if defined(HAVE_AMD) && defined(IPK4)
info = psb_success_
nr = a%get_nrows()
nz = a%get_nzeros()
allocate(perm(nr),stat=info)
if (info /= 0) then
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
end if
allocate(aa, mold=a%a)
call a%mv_to(acsc)
acsc%ia(:) = acsc%ia(:) - 1
acsc%icp(:) = acsc%icp(:) - 1
info = psb_amd_order(nr,acsc%icp,acsc%ia,perm)
if (info /= psb_success_) then
info = psb_err_from_subroutine_
call psb_errpush(info,name,a_err='psb_amd_order')
goto 9999
end if
perm(:) = perm(:) + 1
acsc%ia(:) = acsc%ia(:) + 1
acsc%icp(:) = acsc%icp(:) + 1
call acsc%mv_to_coo(acoo,info)
do i=1, acoo%get_nzeros()
acoo%ia(i) = perm(acoo%ia(i))
acoo%ja(i) = perm(acoo%ja(i))
end do
call acoo%fix(info)
! Get back to where we started from
call aa%mv_from_coo(acoo,info)
call a%mv_from(aa)
if (present(operm)) then
call psb_realloc(nr,operm,info)
if (info /= psb_success_) then
info = psb_err_alloc_dealloc_
call psb_errpush(info,name)
goto 9999
end if
operm(1:nr) = perm(1:nr)
end if
deallocate(aa,perm)
#else
info = psb_err_missing_aux_lib_
call psb_errpush(info, name)
goto 9999
#endif
call psb_erractionrestore(err_act)
return
9999 call psb_error_handler(err_act)
return
end subroutine psb_mat_renum_amd
end subroutine psb_z_mat_renum
subroutine psb_lz_mat_renum(alg,mat,info,perm)
use psb_base_mod
use psb_renum_mod, psb_protect_name => psb_z_mat_renum
use psb_renum_mod, psb_protect_name => psb_lz_mat_renum
implicit none
integer(psb_ipk_), intent(in) :: alg
type(psb_zspmat_type), intent(inout) :: mat
character(len=*), intent(in) :: alg
type(psb_lzspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
integer(psb_lpk_), allocatable, optional, intent(out) :: perm(:)
integer(psb_ipk_) :: err_act, nr, nc, i, ierr(5)
integer(psb_lpk_) :: nr, nc, i
integer(psb_ipk_) :: err_act, ierr(5)
character(len=20) :: name
info = psb_success_
name = 'mat_renum'
call psb_erractionsave(err_act)
info = psb_success_
nr = mat%get_nrows()
nc = mat%get_ncols()
if (nr /= nc) then
@ -101,16 +317,17 @@ subroutine psb_z_mat_renum(alg,mat,info,perm)
goto 9999
end if
select case (alg)
case(psb_mat_renum_gps_)
info = psb_success_
select case (psb_toupper(alg))
case ('GPS')
call psb_mat_renum_gps(mat,info,perm)
call psb_lmat_renum_gps(mat,info,perm)
case(psb_mat_renum_amd_)
case('AMD')
call psb_mat_renum_amd(mat,info,perm)
call psb_lmat_renum_amd(mat,info,perm)
case(psb_mat_renum_identity_)
case ('NONE', 'ID')
nr = mat%get_nrows()
allocate(perm(nr),stat=info)
if (info == 0) then
@ -123,8 +340,9 @@ subroutine psb_z_mat_renum(alg,mat,info,perm)
goto 9999
endif
case default
write(0,*) 'Unknown algorithm "',psb_toupper(alg),'"'
info = psb_err_input_value_invalid_i_
ierr(1) = 1; ierr(2) = alg;
ierr(1) = 1;
call psb_errpush(info,name,i_err=ierr)
goto 9999
end select
@ -143,26 +361,26 @@ subroutine psb_z_mat_renum(alg,mat,info,perm)
contains
subroutine psb_mat_renum_gps(a,info,operm)
subroutine psb_lmat_renum_gps(a,info,operm)
use psb_base_mod
use psb_gps_mod
use psb_lgps_mod
implicit none
type(psb_zspmat_type), intent(inout) :: a
type(psb_lzspmat_type), intent(inout) :: a
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: operm(:)
integer(psb_lpk_), allocatable, optional, intent(out) :: operm(:)
!
class(psb_z_base_sparse_mat), allocatable :: aa
type(psb_z_csr_sparse_mat) :: acsr
type(psb_z_coo_sparse_mat) :: acoo
class(psb_lz_base_sparse_mat), allocatable :: aa
type(psb_lz_csr_sparse_mat) :: acsr
type(psb_lz_coo_sparse_mat) :: acoo
integer(psb_ipk_) :: err_act
character(len=20) :: name
integer(psb_ipk_), allocatable :: ndstk(:,:), iold(:), ndeg(:), perm(:)
integer(psb_ipk_) :: i, j, k, ideg, nr, ibw, ipf, idpth
integer(psb_lpk_), allocatable :: ndstk(:,:), iold(:), ndeg(:), perm(:)
integer(psb_lpk_) :: i, j, k, ideg, nr, ibw, ipf, idpth
info = psb_success_
name = 'mat_renum'
name = 'mat_renum_gps'
call psb_erractionsave(err_act)
info = psb_success_
@ -193,7 +411,7 @@ contains
end do
perm = 0
call psb_gps_reduce(ndstk,nr,ideg,iold,perm,ndeg,ibw,ipf,idpth)
call psb_lgps_reduce(ndstk,nr,ideg,iold,perm,ndeg,ibw,ipf,idpth)
if (.not.psb_isaperm(nr,perm)) then
write(0,*) 'Something wrong: bad perm from gps_reduce'
@ -228,18 +446,19 @@ contains
9999 call psb_error_handler(err_act)
return
end subroutine psb_mat_renum_gps
end subroutine psb_lmat_renum_gps
subroutine psb_mat_renum_amd(a,info,operm)
subroutine psb_lmat_renum_amd(a,info,operm)
#if defined(HAVE_AMD)
use iso_c_binding
#endif
use psb_base_mod
implicit none
type(psb_zspmat_type), intent(inout) :: a
type(psb_lzspmat_type), intent(inout) :: a
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: operm(:)
integer(psb_lpk_), allocatable, optional, intent(out) :: operm(:)
!
#if defined(HAVE_AMD)
@ -254,20 +473,20 @@ contains
end interface
#endif
type(psb_z_csc_sparse_mat) :: acsc
class(psb_z_base_sparse_mat), allocatable :: aa
type(psb_z_coo_sparse_mat) :: acoo
type(psb_lz_csc_sparse_mat) :: acsc
class(psb_lz_base_sparse_mat), allocatable :: aa
type(psb_lz_coo_sparse_mat) :: acoo
integer(psb_ipk_), allocatable :: perm(:)
integer(psb_ipk_) :: err_act
character(len=20) :: name
integer(psb_ipk_) :: i, j, k, ideg, nr, ibw, ipf, idpth, nz
integer(psb_lpk_) :: i, j, k, ideg, nr, ibw, ipf, idpth, nz
info = psb_success_
name = 'mat_renum_amd'
call psb_erractionsave(err_act)
#if defined(HAVE_AMD) && defined(IPK4)
#if defined(HAVE_AMD) && defined(LPK4)
info = psb_success_
nr = a%get_nrows()
@ -331,10 +550,9 @@ contains
9999 call psb_error_handler(err_act)
return
end subroutine psb_mat_renum_amd
end subroutine psb_z_mat_renum
end subroutine psb_lmat_renum_amd
end subroutine psb_lz_mat_renum
subroutine psb_z_cmp_bwpf(mat,bwl,bwu,prf,info)
use psb_base_mod
@ -385,3 +603,52 @@ subroutine psb_z_cmp_bwpf(mat,bwl,bwu,prf,info)
end select
end subroutine psb_z_cmp_bwpf
subroutine psb_lz_cmp_bwpf(mat,bwl,bwu,prf,info)
use psb_base_mod
use psb_renum_mod, psb_protect_name => psb_lz_cmp_bwpf
implicit none
type(psb_lzspmat_type), intent(in) :: mat
integer(psb_lpk_), intent(out) :: bwl, bwu
integer(psb_lpk_), intent(out) :: prf
integer(psb_ipk_), intent(out) :: info
!
integer(psb_lpk_), allocatable :: irow(:), icol(:)
complex(psb_dpk_), allocatable :: val(:)
integer(psb_lpk_) :: nz, i, j, lrbu, lrbl
info = psb_success_
bwl = 0
bwu = 0
prf = 0
select type (aa=>mat%a)
class is (psb_lz_csr_sparse_mat)
do i=1, aa%get_nrows()
lrbl = 0
lrbu = 0
do j = aa%irp(i), aa%irp(i+1) - 1
lrbl = max(lrbl,i-aa%ja(j))
lrbu = max(lrbu,aa%ja(j)-i)
end do
prf = prf + lrbl+lrbu
bwu = max(bwu,lrbu)
bwl = max(bwl,lrbu)
end do
class default
do i=1, aa%get_nrows()
lrbl = 0
lrbu = 0
call aa%csget(i,i,nz,irow,icol,val,info)
if (info /= psb_success_) return
do j=1, nz
lrbl = max(lrbl,i-icol(j))
lrbu = max(lrbu,icol(j)-i)
end do
prf = prf + lrbl+lrbu
bwu = max(bwu,lrbu)
bwl = max(bwl,lrbu)
end do
end select
end subroutine psb_lz_cmp_bwpf

69
util/psb_z_renum_mod.f90
Unescape Escape View File

@ -0,0 +1,69 @@
!
! Parallel Sparse BLAS version 3.5
! (C) Copyright 2006-2018
! Salvatore Filippone
! Alfredo Buttari
!
! 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.
!
!
module psb_z_renum_mod
use psb_base_mod
interface psb_mat_renum
subroutine psb_z_mat_renum(alg,mat,info,perm)
import :: psb_ipk_, psb_zspmat_type
character(len=*), intent(in) :: alg
type(psb_zspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_ipk_), allocatable, optional, intent(out) :: perm(:)
end subroutine psb_z_mat_renum
subroutine psb_lz_mat_renum(alg,mat,info,perm)
import :: psb_ipk_, psb_lpk_, psb_lzspmat_type
character(len=*), intent(in) :: alg
type(psb_lzspmat_type), intent(inout) :: mat
integer(psb_ipk_), intent(out) :: info
integer(psb_lpk_), allocatable, optional, intent(out) :: perm(:)
end subroutine psb_lz_mat_renum
end interface psb_mat_renum
interface psb_cmp_bwpf
subroutine psb_z_cmp_bwpf(mat,bwl,bwu,prf,info)
import :: psb_ipk_, psb_zspmat_type
type(psb_zspmat_type), intent(in) :: mat
integer(psb_ipk_), intent(out) :: bwl, bwu
integer(psb_ipk_), intent(out) :: prf
integer(psb_ipk_), intent(out) :: info
end subroutine psb_z_cmp_bwpf
subroutine psb_lz_cmp_bwpf(mat,bwl,bwu,prf,info)
import :: psb_ipk_, psb_lpk_, psb_lzspmat_type
type(psb_lzspmat_type), intent(in) :: mat
integer(psb_lpk_), intent(out) :: bwl, bwu
integer(psb_lpk_), intent(out) :: prf
integer(psb_ipk_), intent(out) :: info
end subroutine psb_lz_cmp_bwpf
end interface psb_cmp_bwpf
end module psb_z_renum_mod
Write Preview
Loading…
Cancel
Save