psblas3-integer8:

base/modules/psi_penv_mod.F90
 base/modules/psi_reduce_mod.F90

Support versinos of amx/amn for integer4 in 8-bytes compilation mode.
psblas3-type-indexed
Salvatore Filippone 13 years ago
parent 04d97d04bb
commit a58bb22ece

@ -56,6 +56,7 @@ module psi_penv_mod
#else
integer(psb_mpik_), save :: mpi_iamx_op, mpi_iamn_op
integer(psb_mpik_), save :: mpi_i4amx_op, mpi_i4amn_op
integer(psb_mpik_), save :: mpi_i8amx_op, mpi_i8amn_op
integer(psb_mpik_), save :: mpi_samx_op, mpi_samn_op
integer(psb_mpik_), save :: mpi_damx_op, mpi_damn_op
@ -69,6 +70,7 @@ module psi_penv_mod
private :: psi_get_sizes, psi_register_mpi_extras
private :: psi_iamx_op, psi_iamn_op
private :: psi_i4amx_op, psi_i4amn_op
private :: psi_i8amx_op, psi_i8amn_op
private :: psi_samx_op, psi_samn_op
private :: psi_damx_op, psi_damn_op
@ -121,6 +123,8 @@ contains
#else
if (info == 0) call mpi_op_create(psi_iamx_op,.true.,mpi_iamx_op,info)
if (info == 0) call mpi_op_create(psi_iamn_op,.true.,mpi_iamn_op,info)
if (info == 0) call mpi_op_create(psi_i4amx_op,.true.,mpi_i4amx_op,info)
if (info == 0) call mpi_op_create(psi_i4amn_op,.true.,mpi_i4amn_op,info)
if (info == 0) call mpi_op_create(psi_i8amx_op,.true.,mpi_i8amx_op,info)
if (info == 0) call mpi_op_create(psi_i8amn_op,.true.,mpi_i8amn_op,info)
if (info == 0) call mpi_op_create(psi_samx_op,.true.,mpi_samx_op,info)
@ -520,6 +524,26 @@ contains
end do
end subroutine psi_iamn_op
subroutine psi_i4amx_op(inv, outv,len,type)
integer(psb_mpik_) :: inv(*),outv(*)
integer(psb_mpik_) :: len,type
integer(psb_mpik_) :: i
do i=1, len
if (abs(inv(i)) > abs(outv(i))) outv(i) = inv(i)
end do
end subroutine psi_i4amx_op
subroutine psi_i4amn_op(inv, outv,len,type)
integer(psb_mpik_) :: inv(*),outv(*)
integer(psb_mpik_) :: len,type
integer(psb_mpik_) :: i
do i=1, len
if (abs(inv(i)) < abs(outv(i))) outv(i) = inv(i)
end do
end subroutine psi_i4amn_op
subroutine psi_i8amx_op(inv, outv,len,type)
integer(psb_long_int_k_) :: inv(*),outv(*)
integer(psb_mpik_) :: len,type
@ -531,19 +555,10 @@ contains
end subroutine psi_i8amx_op
subroutine psi_i8amn_op(inv, outv,len,type)
#ifdef MPI_MOD
use mpi
#endif
implicit none
#ifdef MPI_H
include 'mpif.h'
#endif
integer(psb_long_int_k_) :: inv(*),outv(*)
integer(psb_mpik_) :: len,type
integer(psb_mpik_) :: i
if (type /= psb_mpi_lng_integer) then
write(psb_err_unit,*) 'Invalid type !!!'
end if
do i=1, len
if (abs(inv(i)) < abs(outv(i))) outv(i) = inv(i)
end do

@ -35,6 +35,11 @@ module psi_reduce_mod
module procedure psb_i8amxs, psb_i8amxv, psb_i8amxm
end interface
#endif
#if defined(LONG_INTEGERS)
interface psb_amx
module procedure psb_i4amxs, psb_i4amxv, psb_i4amxm
end interface
#endif
interface psb_amn
module procedure psb_iamns, psb_iamnv, psb_iamnm,&
@ -43,6 +48,11 @@ module psi_reduce_mod
& psb_damns, psb_damnv, psb_damnm,&
& psb_zamns, psb_zamnv, psb_zamnm
end interface
#if defined(LONG_INTEGERS)
interface psb_amn
module procedure psb_i4amns, psb_i4amnv, psb_i4amnm
end interface
#endif
#if !defined(LONG_INTEGERS)
interface psb_amn
module procedure psb_i8amns, psb_i8amnv, psb_i8amnm
@ -1284,6 +1294,134 @@ contains
end subroutine psb_iamxm
#if defined(LONG_INTEGERS)
subroutine psb_i4amxs(ictxt,dat,root)
#ifdef MPI_MOD
use mpi
#endif
implicit none
#ifdef MPI_H
include 'mpif.h'
#endif
integer(psb_mpik_), intent(in) :: ictxt
integer(psb_mpik_), intent(inout) :: dat
integer(psb_mpik_), intent(in), optional :: root
integer(psb_mpik_) :: root_
integer(psb_mpik_) :: dat_
integer(psb_mpik_) :: iam, np, info
integer(psb_mpik_) :: iinfo
#if !defined(SERIAL_MPI)
call psb_info(ictxt,iam,np)
if (present(root)) then
root_ = root
else
root_ = -1
endif
if (root_ == -1) then
call mpi_allreduce(dat,dat_,1,psb_mpi_def_integer,mpi_i4amx_op,ictxt,info)
dat = dat_
else
call mpi_reduce(dat,dat_,1,psb_mpi_def_integer,mpi_i4amx_op,root_,ictxt,info)
dat = dat_
endif
#endif
end subroutine psb_i4amxs
subroutine psb_i4amxv(ictxt,dat,root)
use psb_realloc_mod
#ifdef MPI_MOD
use mpi
#endif
implicit none
#ifdef MPI_H
include 'mpif.h'
#endif
integer(psb_mpik_), intent(in) :: ictxt
integer(psb_mpik_), intent(inout) :: dat(:)
integer(psb_mpik_), intent(in), optional :: root
integer(psb_mpik_) :: root_
integer(psb_mpik_), allocatable :: dat_(:)
integer(psb_mpik_) :: iam, np, info
integer(psb_mpik_) :: iinfo
#if !defined(SERIAL_MPI)
call psb_info(ictxt,iam,np)
if (present(root)) then
root_ = root
else
root_ = -1
endif
if (root_ == -1) then
call psb_realloc(size(dat),dat_,iinfo)
dat_=dat
if (iinfo == psb_success_) call mpi_allreduce(dat_,dat,size(dat),&
& psb_mpi_def_integer,mpi_i4amx_op,ictxt,info)
else
if (iam == root_) then
call psb_realloc(size(dat),dat_,iinfo)
dat_=dat
call mpi_reduce(dat_,dat,size(dat),psb_mpi_def_integer,mpi_i4amx_op,root_,ictxt,info)
else
call psb_realloc(1,dat_,iinfo)
call mpi_reduce(dat,dat_,size(dat),psb_mpi_def_integer,mpi_i4amx_op,root_,ictxt,info)
end if
endif
#endif
end subroutine psb_i4amxv
subroutine psb_i4amxm(ictxt,dat,root)
use psb_realloc_mod
#ifdef MPI_MOD
use mpi
#endif
implicit none
#ifdef MPI_H
include 'mpif.h'
#endif
integer(psb_mpik_), intent(in) :: ictxt
integer(psb_mpik_), intent(inout) :: dat(:,:)
integer(psb_mpik_), intent(in), optional :: root
integer(psb_mpik_) :: root_
integer(psb_mpik_), allocatable :: dat_(:,:)
integer(psb_mpik_) :: iam, np, info
integer(psb_mpik_) :: iinfo
#if !defined(SERIAL_MPI)
call psb_info(ictxt,iam,np)
if (present(root)) then
root_ = root
else
root_ = -1
endif
if (root_ == -1) then
call psb_realloc(size(dat,1),size(dat,2),dat_,iinfo)
dat_=dat
if (iinfo == psb_success_) call mpi_allreduce(dat_,dat,size(dat),&
& psb_mpi_def_integer,mpi_i4amx_op,ictxt,info)
else
if (iam == root_) then
call psb_realloc(size(dat,1),size(dat,2),dat_,iinfo)
dat_=dat
call mpi_reduce(dat_,dat,size(dat),psb_mpi_def_integer,mpi_i4amx_op,root_,ictxt,info)
else
call psb_realloc(1,1,dat_,iinfo)
call mpi_reduce(dat,dat_,size(dat),psb_mpi_def_integer,mpi_i4amx_op,root_,ictxt,info)
end if
endif
#endif
end subroutine psb_i4amxm
#endif
#if !defined(LONG_INTEGERS)
subroutine psb_i8amxs(ictxt,dat,root)
@ -2043,6 +2181,134 @@ contains
end subroutine psb_iamnm
#if defined(LONG_INTEGERS)
subroutine psb_i4amns(ictxt,dat,root)
#ifdef MPI_MOD
use mpi
#endif
implicit none
#ifdef MPI_H
include 'mpif.h'
#endif
integer(psb_mpik_), intent(in) :: ictxt
integer(psb_mpik_), intent(inout) :: dat
integer(psb_mpik_), intent(in), optional :: root
integer(psb_mpik_) :: root_
integer(psb_mpik_) :: dat_
integer(psb_mpik_) :: iam, np, info
integer(psb_mpik_) :: iinfo
#if !defined(SERIAL_MPI)
call psb_info(ictxt,iam,np)
if (present(root)) then
root_ = root
else
root_ = -1
endif
if (root_ == -1) then
call mpi_allreduce(dat,dat_,1,psb_mpi_def_integer,mpi_i4amn_op,ictxt,info)
dat = dat_
else
call mpi_reduce(dat,dat_,1,psb_mpi_def_integer,mpi_i4amn_op,root_,ictxt,info)
dat = dat_
endif
#endif
end subroutine psb_i4amns
subroutine psb_i4amnv(ictxt,dat,root)
use psb_realloc_mod
#ifdef MPI_MOD
use mpi
#endif
implicit none
#ifdef MPI_H
include 'mpif.h'
#endif
integer(psb_mpik_), intent(in) :: ictxt
integer(psb_mpik_), intent(inout) :: dat(:)
integer(psb_mpik_), intent(in), optional :: root
integer(psb_mpik_) :: root_
integer(psb_mpik_), allocatable :: dat_(:)
integer(psb_mpik_) :: iam, np, info
integer(psb_mpik_) :: iinfo
#if !defined(SERIAL_MPI)
call psb_info(ictxt,iam,np)
if (present(root)) then
root_ = root
else
root_ = -1
endif
if (root_ == -1) then
call psb_realloc(size(dat),dat_,iinfo)
dat_=dat
if (iinfo == psb_success_) call mpi_allreduce(dat_,dat,size(dat),&
& psb_mpi_def_integer,mpi_i4amn_op,ictxt,info)
else
if (iam == root_) then
call psb_realloc(size(dat),dat_,iinfo)
dat_=dat
call mpi_reduce(dat_,dat,size(dat),psb_mpi_def_integer,mpi_i4amn_op,root_,ictxt,info)
else
call psb_realloc(1,dat_,iinfo)
call mpi_reduce(dat,dat_,size(dat),psb_mpi_def_integer,mpi_i4amn_op,root_,ictxt,info)
end if
endif
#endif
end subroutine psb_i4amnv
subroutine psb_i4amnm(ictxt,dat,root)
use psb_realloc_mod
#ifdef MPI_MOD
use mpi
#endif
implicit none
#ifdef MPI_H
include 'mpif.h'
#endif
integer(psb_mpik_), intent(in) :: ictxt
integer(psb_mpik_), intent(inout) :: dat(:,:)
integer(psb_mpik_), intent(in), optional :: root
integer(psb_mpik_) :: root_
integer(psb_mpik_), allocatable :: dat_(:,:)
integer(psb_mpik_) :: iam, np, info
integer(psb_mpik_) :: iinfo
#if !defined(SERIAL_MPI)
call psb_info(ictxt,iam,np)
if (present(root)) then
root_ = root
else
root_ = -1
endif
if (root_ == -1) then
call psb_realloc(size(dat,1),size(dat,2),dat_,iinfo)
dat_=dat
if (iinfo == psb_success_) call mpi_allreduce(dat_,dat,size(dat),&
& psb_mpi_def_integer,mpi_i4amn_op,ictxt,info)
else
if (iam == root_) then
call psb_realloc(size(dat,1),size(dat,2),dat_,iinfo)
dat_=dat
call mpi_reduce(dat_,dat,size(dat),psb_mpi_def_integer,mpi_i4amn_op,root_,ictxt,info)
else
call psb_realloc(1,1,dat_,iinfo)
call mpi_reduce(dat,dat_,size(dat),psb_mpi_def_integer,mpi_i4amn_op,root_,ictxt,info)
end if
endif
#endif
end subroutine psb_i4amnm
#endif
#if !defined(LONG_INTEGERS)
subroutine psb_i8amns(ictxt,dat,root)

Loading…
Cancel
Save