Skip to content
GitLab
Projects
Groups
Snippets
/
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in
Toggle navigation
Menu
Open sidebar
elpa
elpa
Commits
eda9513f
Commit
eda9513f
authored
Apr 08, 2017
by
Andreas Marek
Browse files
New interface: map to use new routines for old API
parent
9ec10a93
Changes
8
Pipelines
1
Hide whitespace changes
Inline
Side-by-side
src/elpa1/elpa1_template.X90
View file @
eda9513f
...
...
@@ -61,8 +61,6 @@ function elpa_solve_evp_&
& (na, nev, a, lda, ev, q, ldq, nblk, matrixCols, mpi_comm_rows, mpi_comm_cols, mpi_comm_all, &
useGPU) result(success)
use precision
use cuda_functions
use mod_check_for_gpu
#ifdef HAVE_DETAILED_TIMINGS
use timings
#else
...
...
@@ -70,196 +68,105 @@ function elpa_solve_evp_&
#endif
use iso_c_binding
use elpa_mpi
use elpa
1_comput
e
use elpa
_typ
e
implicit none
integer(kind=c_int), intent(in) :: na, nev, lda, ldq, nblk, matrixCols, mpi_comm_rows, &
mpi_comm_cols, mpi_comm_all
real(kind=REAL_DATATYPE), intent(out) :: ev(na)
integer(kind=c_int) :: my_pe, n_pes, my_prow, my_pcol, mpierr
#if REALCASE == 1
#ifdef USE_ASSUMED_SIZE
real(kind=C_DATATYPE_KIND), intent(inout) :: a(lda,*)
real(kind=C_DATATYPE_KIND), intent(out) :: q(ldq,*)
#else
!
#ifdef USE_ASSUMED_SIZE
!
real(kind=C_DATATYPE_KIND), intent(inout) :: a(lda,*)
!
real(kind=C_DATATYPE_KIND), intent(out) :: q(ldq,*)
!
#else
real(kind=C_DATATYPE_KIND), intent(inout) :: a(lda,matrixCols)
real(kind=C_DATATYPE_KIND), intent(out) :: q(ldq,matrixCols)
#endif
real(kind=C_DATATYPE_KIND), allocatable :: tau(:)
!#endif
#endif /* REALCASE */
#if COMPLEXCASE == 1
#ifdef USE_ASSUMED_SIZE
complex(kind=C_DATATYPE_KIND), intent(inout) :: a(lda,*)
complex(kind=C_DATATYPE_KIND), intent(out) :: q(ldq,*)
#else
!
#ifdef USE_ASSUMED_SIZE
!
complex(kind=C_DATATYPE_KIND), intent(inout) :: a(lda,*)
!
complex(kind=C_DATATYPE_KIND), intent(out) :: q(ldq,*)
!
#else
complex(kind=C_DATATYPE_KIND), intent(inout) :: a(lda,matrixCols)
complex(kind=C_DATATYPE_KIND), intent(out) :: q(ldq,matrixCols)
#endif
!
#endif
real(kind=REAL_DATATYPE), allocatable :: q_real(:,:)
complex(kind=C_DATATYPE_KIND), allocatable :: tau(:)
integer(kind=c_int) :: l_cols, l_rows, l_cols_nev, np_rows, np_cols
#endif /* COMPLEXCASE */
logical, intent(in), optional :: useGPU
logical :: success
logical :: do_useGPU
integer(kind=ik) :: numberOfGPUDevices
integer(kind=c_int) :: my_pe, n_pes, my_prow, my_pcol, mpierr
real(kind=C_DATATYPE_KIND), allocatable :: e(:)
real(kind=c_double) :: ttt0, ttt1 ! MPI_WTIME always needs double
logical, save :: firstCall = .true.
logical :: wantDebug
integer(kind=c_int) :: istat
character(200) :: errorMessage
integer(kind=c_int) :: successInternal
type(elpa_t) :: elpa1stage
call timer%start("elpa_solve_evp_&
&MATH_DATATYPE&
&_1stage_&
&PRECISION&
&")
call timer%start("mpi_communication")
call mpi_comm_rank(mpi_comm_all,my_pe,mpierr)
call mpi_comm_size(mpi_comm_all,n_pes,mpierr)
call mpi_comm_rank(mpi_comm_rows,my_prow,mpierr)
call mpi_comm_rank(mpi_comm_cols,my_pcol,mpierr)
#if COMPLEXCASE == 1
call mpi_comm_size(mpi_comm_rows,np_rows,mpierr)
call mpi_comm_size(mpi_comm_cols,np_cols,mpierr)
#endif
call timer%stop("mpi_communication")
success = .true.
if (elpa_init(20170403) /= ELPA_OK) then
success = .false.
error stop "ELPA API version not supported"
endif
wantDebug = .false.
if (firstCall) then
! are debug messages desired?
wantDebug = debug_messages_via_environment_variable()
firstCall = .false.
elpa1stage = elpa_create(na, nev, lda, matrixCols, nblk, mpi_comm_all, &
my_prow, my_pcol, successInternal)
if (successInternal .ne. ELPA_OK) then
print *, "Cannot create elpa object"
success = .false.
stop
return
endif
do_useGPU = .false.
call elpa1stage%set("solver", ELPA_SOLVER_1STAGE, successInternal)
if (successInternal .ne. ELPA_OK) then
print *, "Cannot set ELPA 1stage solver"
stop
success = .false.
return
endif
if (present(useGPU)) then
! user defined GPU usage via the optional argument in the API call
if (useGPU) then
if (check_for_gpu(my_pe,numberOfGPUDevices, wantDebug=wantDebug)) then
do_useGPU = .true.
! set the neccessary parameters
cudaMemcpyHostToDevice = cuda_memcpyHostToDevice()
cudaMemcpyDeviceToHost = cuda_memcpyDeviceToHost()
cudaMemcpyDeviceToDevice = cuda_memcpyDeviceToDevice()
cudaHostRegisterPortable = cuda_hostRegisterPortable()
cudaHostRegisterMapped = cuda_hostRegisterMapped()
else
print *,"GPUs are requested but not detected! Aborting..."
call elpa1stage%set("gpu", 1, successInternal)
if (successInternal .ne. ELPA_OK) then
print *, "Cannot set gpu"
stop
success = .false.
return
endif
else
call elpa1stage%set("gpu", 0, successInternal)
if (successInternal .ne. ELPA_OK) then
print *, "Cannot set gpu"
stop
success = .false.
return
endif
endif
else
! check whether set by environment variable
do_useGPU = gpu_usage_via_environment_variable()
endif
#if COMPLEXCASE == 1
l_rows = local_index(na, my_prow, np_rows, nblk, -1) ! Local rows of a and q
l_cols = local_index(na, my_pcol, np_cols, nblk, -1) ! Local columns of q
l_cols_nev = local_index(nev, my_pcol, np_cols, nblk, -1) ! Local columns corresponding to nev
allocate(q_real(l_rows,l_cols), stat=istat, errmsg=errorMessage)
if (istat .ne. 0) then
print *,"solve_evp_&
&MATH_DATATYPE&
&_1stage_&
&PRECISION&
&" // ": error when allocating q_real "//errorMessage
stop 1
endif
#endif
allocate(e(na), tau(na), stat=istat, errmsg=errorMessage)
if (istat .ne. 0) then
print *,"solve_evp_&
&MATH_DATATYPE&
&_1stage_&
&PRECISION&
&" // ": error when allocating e, tau "//errorMessage
stop 1
call elpa1stage%solve(a(1:lda,1:matrixCols), ev, q(1:ldq,1:matrixCols), successInternal)
if (successInternal .ne. ELPA_OK) then
print *, "Cannot solve with ELPA 1stage"
stop
success = .false.
return
endif
ttt0 = MPI_Wtime()
call tridiag_&
&MATH_DATATYPE&
&_&
&PRECISION&
& (na, a, lda, nblk, matrixCols, mpi_comm_rows, mpi_comm_cols, ev, e, tau, do_useGPU)
ttt1 = MPI_Wtime()
if(my_prow==0 .and. my_pcol==0 .and. elpa_print_times) write(error_unit,*) 'Time tridiag_real :',ttt1-ttt0
time_evp_fwd = ttt1-ttt0
ttt0 = MPI_Wtime()
call solve_tridi_&
&PRECISION&
& (na, nev, ev, e, &
#if REALCASE == 1
q, ldq, &
#endif
#if COMPLEXCASE == 1
q_real, l_rows, &
#endif
nblk, matrixCols, mpi_comm_rows, mpi_comm_cols, wantDebug, success)
if (.not.(success)) return
ttt1 = MPI_Wtime()
if(my_prow==0 .and. my_pcol==0 .and. elpa_print_times) write(error_unit,*) 'Time solve_tridi :',ttt1-ttt0
time_evp_solve = ttt1-ttt0
ttt0 = MPI_Wtime()
#if COMPLEXCASE == 1
q(1:l_rows,1:l_cols_nev) = q_real(1:l_rows,1:l_cols_nev)
#endif
call trans_ev_&
&MATH_DATATYPE&
&_&
&PRECISION&
& (na, nev, a, lda, tau, q, ldq, nblk, matrixCols, mpi_comm_rows, mpi_comm_cols, do_useGPU)
ttt1 = MPI_Wtime()
if(my_prow==0 .and. my_pcol==0 .and. elpa_print_times) write(error_unit,*) 'Time trans_ev_real:',ttt1-ttt0
time_evp_back = ttt1-ttt0
#if COMPLEXCASE == 1
deallocate(q_real, stat=istat, errmsg=errorMessage)
if (istat .ne. 0) then
print *,"solve_evp_&
&MATH_DATATYPE&
&_1stage_&
&PRECISION&
&" // ": error when deallocating q_real "//errorMessage
stop 1
endif
#endif
call elpa1stage%destroy()
deallocate(e, tau, stat=istat, errmsg=errorMessage)
if (istat .ne. 0) then
print *,"solve_evp_&
&MATH_DATATYPE&
&_1stage_&
&PRECISION&
&" // ": error when deallocating e, tau "//errorMessage
stop 1
endif
call elpa_uninit()
call timer%stop("elpa_solve_evp_&
&MATH_DATATYPE&
...
...
src/elpa1/elpa1_template_new_interface.X90
View file @
eda9513f
...
...
@@ -119,7 +119,7 @@ function elpa_solve_evp_&
&MATH_DATATYPE&
&_1stage_&
&PRECISION&
&")
&
_new
")
call timer%start("mpi_communication")
...
...
@@ -279,7 +279,7 @@ function elpa_solve_evp_&
&MATH_DATATYPE&
&_1stage_&
&PRECISION&
&")
&
_new
")
end function
...
...
src/elpa2/elpa2_compute.F90
View file @
eda9513f
...
...
@@ -60,7 +60,7 @@ module ELPA2_compute
use
ELPA_utilities
USE
ELPA1_compute
use
elpa1
,
only
:
elpa_print_times
,
time_evp_back
,
time_evp_fwd
,
time_evp_solve
use
elpa1
_new
,
only
:
elpa_print_times
,
time_evp_back
,
time_evp_fwd
,
time_evp_solve
use
elpa2_utilities
use
elpa_pdgeqrf
use
precision
...
...
src/elpa2/elpa2_new_interface.F90
View file @
eda9513f
...
...
@@ -59,7 +59,7 @@ module ELPA2_new
! Version 1.1.2, 2011-02-21
use
elpa_utilities
use
elpa1
,
only
:
elpa_print_times
,
time_evp_back
,
time_evp_fwd
,
time_evp_solve
use
elpa1
_new
,
only
:
elpa_print_times
,
time_evp_back
,
time_evp_fwd
,
time_evp_solve
use
elpa2_utilities
implicit
none
...
...
src/elpa2/qr/elpa_pdgeqrf.F90
View file @
eda9513f
...
...
@@ -71,7 +71,7 @@ module elpa_pdgeqrf
work
,
workLength
,
lwork
,
m
,
n
,
mb
,
nb
,
rowidx
,
colidx
,
&
rev
,
trans
,
PQRPARAM
,
mpicomm_rows
,
mpicomm_cols
,
blockheuristic
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -388,7 +388,7 @@ module elpa_pdgeqrf
subroutine
qr_pdgeqrf_1dcomm_double
(
a
,
lda
,
v
,
ldv
,
tau
,
t
,
ldt
,
work
,
lwork
,
m
,
n
,
mb
,
baseidx
,
rowidx
,
rev
,
trans
,
&
PQRPARAM
,
mpicomm
,
blockheuristic
)
use
precision
use
ELPA1
use
elpa1_new
#ifdef HAVE_DETAILED_TIMINGS
use
timings
#endif
...
...
@@ -528,7 +528,7 @@ module elpa_pdgeqrf
subroutine
qr_pdgeqr2_1dcomm_double
(
a
,
lda
,
v
,
ldv
,
tau
,
t
,
ldt
,
work
,
lwork
,
m
,
n
,
mb
,
baseidx
,
rowidx
,
rev
,
&
trans
,
PQRPARAM
,
mpicomm
,
blockheuristic
)
use
precision
use
ELPA1
use
elpa1_new
#ifdef HAVE_DETAILED_TIMINGS
use
timings
#endif
...
...
@@ -721,7 +721,7 @@ module elpa_pdgeqrf
subroutine
qr_pdlarfg_1dcomm_double
(
x
,
incx
,
tau
,
work
,
lwork
,
n
,
idx
,
nb
,
hgmode
,
rev
,
communicator
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -1153,7 +1153,7 @@ module elpa_pdgeqrf
subroutine
qr_pdlarfg2_1dcomm_seed_double
(
a
,
lda
,
work
,
lwork
,
seed
,
n
,
nb
,
idx
,
rev
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -1353,7 +1353,7 @@ module elpa_pdgeqrf
! id=1: second vector
subroutine
qr_pdlarfg2_1dcomm_vector_double
(
x
,
incx
,
tau
,
seed
,
n
,
nb
,
idx
,
id
,
rev
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -1449,7 +1449,7 @@ module elpa_pdgeqrf
subroutine
qr_pdlarfg2_1dcomm_update_double
(
v
,
incv
,
baseidx
,
a
,
lda
,
seed
,
n
,
idx
,
nb
,
rev
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -1688,7 +1688,7 @@ module elpa_pdgeqrf
subroutine
qr_pdlarfgk_1dcomm_seed_double
(
a
,
lda
,
baseidx
,
work
,
lwork
,
seedC
,
seedD
,
m
,
k
,
mb
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -2202,7 +2202,7 @@ module elpa_pdgeqrf
!rank: actual rank (k >= rank)
subroutine
qr_pdlarfgk_1dcomm_vector_double
(
x
,
incx
,
baseidx
,
tau
,
seedC
,
seedD
,
k
,
sidx
,
n
,
nb
,
rev
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -2311,7 +2311,7 @@ module elpa_pdgeqrf
! computed
subroutine
qr_pdlarfgk_1dcomm_update_double
(
a
,
lda
,
baseidx
,
work
,
lwork
,
seedC
,
seedD
,
k
,
rank
,
sidx
,
tau
,
n
,
nb
,
rev
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -2548,7 +2548,7 @@ module elpa_pdgeqrf
!direction=1: unpack from work buffer
subroutine
qr_pdgeqrf_pack_unpack_double
(
v
,
ldv
,
work
,
lwork
,
m
,
n
,
mb
,
baseidx
,
rowidx
,
rev
,
direction
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -2614,7 +2614,7 @@ module elpa_pdgeqrf
!direction=1: unpack from work buffer
subroutine
qr_pdgeqrf_pack_unpack_tmatrix_double
(
tau
,
t
,
ldt
,
work
,
lwork
,
n
,
direction
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -2756,7 +2756,7 @@ module elpa_pdgeqrf
subroutine
qr_pdlarfg_copy_1dcomm_double
(
x
,
incx
,
v
,
incv
,
n
,
baseidx
,
idx
,
nb
,
rev
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -2818,7 +2818,7 @@ module elpa_pdgeqrf
work
,
workLength
,
lwork
,
m
,
n
,
mb
,
nb
,
rowidx
,
colidx
,
&
rev
,
trans
,
PQRPARAM
,
mpicomm_rows
,
mpicomm_cols
,
blockheuristic
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -3135,7 +3135,7 @@ module elpa_pdgeqrf
subroutine
qr_pdgeqrf_1dcomm_single
(
a
,
lda
,
v
,
ldv
,
tau
,
t
,
ldt
,
work
,
lwork
,
m
,
n
,
mb
,
baseidx
,
&
rowidx
,
rev
,
trans
,
PQRPARAM
,
mpicomm
,
blockheuristic
)
use
precision
use
ELPA1
use
elpa1_new
#ifdef HAVE_DETAILED_TIMINGS
use
timings
#endif
...
...
@@ -3275,7 +3275,7 @@ module elpa_pdgeqrf
subroutine
qr_pdgeqr2_1dcomm_single
(
a
,
lda
,
v
,
ldv
,
tau
,
t
,
ldt
,
work
,
lwork
,
m
,
n
,
mb
,
baseidx
,
&
rowidx
,
rev
,
trans
,
PQRPARAM
,
mpicomm
,
blockheuristic
)
use
precision
use
ELPA1
use
elpa1_new
#ifdef HAVE_DETAILED_TIMINGS
use
timings
#endif
...
...
@@ -3467,7 +3467,7 @@ module elpa_pdgeqrf
subroutine
qr_pdlarfg_1dcomm_single
(
x
,
incx
,
tau
,
work
,
lwork
,
n
,
idx
,
nb
,
hgmode
,
rev
,
communicator
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -3905,7 +3905,7 @@ module elpa_pdgeqrf
subroutine
qr_pdlarfg2_1dcomm_seed_single
(
a
,
lda
,
work
,
lwork
,
seed
,
n
,
nb
,
idx
,
rev
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -4108,7 +4108,7 @@ module elpa_pdgeqrf
! id=1: second vector
subroutine
qr_pdlarfg2_1dcomm_vector_single
(
x
,
incx
,
tau
,
seed
,
n
,
nb
,
idx
,
id
,
rev
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -4204,7 +4204,7 @@ module elpa_pdgeqrf
subroutine
qr_pdlarfg2_1dcomm_update_single
(
v
,
incv
,
baseidx
,
a
,
lda
,
seed
,
n
,
idx
,
nb
,
rev
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -4444,7 +4444,7 @@ module elpa_pdgeqrf
subroutine
qr_pdlarfgk_1dcomm_seed_single
(
a
,
lda
,
baseidx
,
work
,
lwork
,
seedC
,
seedD
,
m
,
k
,
mb
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -4958,7 +4958,7 @@ module elpa_pdgeqrf
!rank: actual rank (k >= rank)
subroutine
qr_pdlarfgk_1dcomm_vector_single
(
x
,
incx
,
baseidx
,
tau
,
seedC
,
seedD
,
k
,
sidx
,
n
,
nb
,
rev
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -5066,7 +5066,7 @@ module elpa_pdgeqrf
! computed
subroutine
qr_pdlarfgk_1dcomm_update_single
(
a
,
lda
,
baseidx
,
work
,
lwork
,
seedC
,
seedD
,
k
,
rank
,
sidx
,
tau
,
n
,
nb
,
rev
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -5303,7 +5303,7 @@ module elpa_pdgeqrf
!direction=1: unpack from work buffer
subroutine
qr_pdgeqrf_pack_unpack_single
(
v
,
ldv
,
work
,
lwork
,
m
,
n
,
mb
,
baseidx
,
rowidx
,
rev
,
direction
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -5369,7 +5369,7 @@ module elpa_pdgeqrf
!direction=1: unpack from work buffer
subroutine
qr_pdgeqrf_pack_unpack_tmatrix_single
(
tau
,
t
,
ldt
,
work
,
lwork
,
n
,
direction
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
@@ -5471,7 +5471,7 @@ module elpa_pdgeqrf
subroutine
qr_pdlarfg_copy_1dcomm_single
(
x
,
incx
,
v
,
incv
,
n
,
baseidx
,
idx
,
nb
,
rev
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
#ifdef HAVE_DETAILED_TIMINGS
use
timings
...
...
src/elpa2/qr/elpa_pdlarfb.F90
View file @
eda9513f
...
...
@@ -416,7 +416,7 @@ end subroutine qr_pdlarft_tree_merge_1dcomm_double
! - assume right positions for v
subroutine
qr_pdlarfl_1dcomm_double
(
v
,
incv
,
baseidx
,
a
,
lda
,
tau
,
work
,
lwork
,
m
,
n
,
idx
,
mb
,
rev
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
implicit
none
...
...
@@ -500,7 +500,7 @@ end subroutine qr_pdlarfl_1dcomm_double
subroutine
qr_pdlarfl2_tmatrix_1dcomm_double
(
v
,
ldv
,
baseidx
,
a
,
lda
,
t
,
ldt
,
work
,
lwork
,
m
,
n
,
idx
,
mb
,
rev
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
implicit
none
...
...
@@ -1171,7 +1171,7 @@ end subroutine qr_pdlarft_tree_merge_1dcomm_single
! - assume right positions for v
subroutine
qr_pdlarfl_1dcomm_single
(
v
,
incv
,
baseidx
,
a
,
lda
,
tau
,
work
,
lwork
,
m
,
n
,
idx
,
mb
,
rev
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
implicit
none
...
...
@@ -1255,7 +1255,7 @@ end subroutine qr_pdlarfl_1dcomm_single
subroutine
qr_pdlarfl2_tmatrix_1dcomm_single
(
v
,
ldv
,
baseidx
,
a
,
lda
,
t
,
ldt
,
work
,
lwork
,
m
,
n
,
idx
,
mb
,
rev
,
mpicomm
)
use
precision
use
ELPA1
use
elpa1_new
use
qr_utils_mod
implicit
none
...
...
src/elpa_c_interface.F90
View file @
eda9513f
...
...
@@ -56,14 +56,14 @@
mpi_comm_rows
,
mpi_comm_cols
)
&
result
(
mpierr
)
bind
(
C
,
name
=
"get_elpa_row_col_comms"
)
use
,
intrinsic
::
iso_c_binding
use
elpa1
,
only
:
get_
elpa_
row_col_comm
s
use
elpa1
,
only
:
elpa_
get_communicator
s
implicit
none
integer
(
kind
=
c_int
)
::
mpierr
integer
(
kind
=
c_int
),
value
::
mpi_comm_world
,
my_prow
,
my_pcol
integer
(
kind
=
c_int
)
::
mpi_comm_rows
,
mpi_comm_cols
mpierr
=
get_
elpa_
row_col_comm
s
(
mpi_comm_world
,
my_prow
,
my_pcol
,
&
mpierr
=
elpa_
get_communicator
s
(
mpi_comm_world
,
my_prow
,
my_pcol
,
&
mpi_comm_rows
,
mpi_comm_cols
)
end
function
...
...
@@ -82,14 +82,14 @@
mpi_comm_rows
,
mpi_comm_cols
)
&
result
(
mpierr
)
bind
(
C
,
name
=
"get_elpa_communicators"
)
use
,
intrinsic
::
iso_c_binding
use
elpa1
,
only
:
get_
elpa_communicators
use
elpa1
,
only
:
elpa
_get
_communicators
implicit
none
integer
(
kind
=
c_int
)
::
mpierr
integer
(
kind
=
c_int
),
value
::
mpi_comm_world
,
my_prow
,
my_pcol
integer
(
kind
=
c_int
)
::
mpi_comm_rows
,
mpi_comm_cols
mpierr
=
get_
elpa_communicators
(
mpi_comm_world
,
my_prow
,
my_pcol
,
&
mpierr
=
elpa
_get
_communicators
(
mpi_comm_world
,
my_prow
,
my_pcol
,
&
mpi_comm_rows
,
mpi_comm_cols
)
end
function
...
...
src/elpa_t.F90
View file @
eda9513f
...
...
@@ -49,8 +49,8 @@ module elpa_type
type
::
elpa_t
private
type
(
c_ptr
)
::
options
=
C_NULL_PTR
integer
::
mpi_comm_parent
=
0
type
(
c_ptr
)
::
options
=
C_NULL_PTR
integer
::
mpi_comm_parent
=
0
integer
(
kind
=
c_int
)
::
mpi_comm_rows
=
0
integer
(
kind
=
c_int
)
::
mpi_comm_cols
=
0
integer
(
kind
=
c_int
)
::
na
=
0
...
...
@@ -105,10 +105,10 @@ module elpa_type
interface
function
get_int_option
(
options
,
name
,
success
)
result
(
value
)
bind
(
C
,
name
=
"get_int_option"
)
import
c_ptr
,
c_char
,
c_int
type
(
c_ptr
),
value
::
options
type
(
c_ptr
),
value
::
options
character
(
kind
=
c_char
),
intent
(
in
)
::
name
(
*
)
integer
(
kind
=
c_int
)
::
value
integer
(
kind
=
c_int
),
optional
::
success
integer
(
kind
=
c_int
)
::
value
integer
(
kind
=
c_int
),
optional
::
success
end
function
end
interface
...
...
@@ -116,10 +116,10 @@ module elpa_type
interface
function
set_int_option
(
options
,
name
,
value
)
result
(
success
)
bind
(
C
,
name
=
"set_int_option"
)
import
c_ptr
,
c_char
,
c_int
type
(
c_ptr
),
value
::
options
character
(
kind
=
c_char
),
intent
(
in
)
::
name
(
*
)
type
(
c_ptr
),
value
::
options
character
(
kind
=
c_char
),
intent
(
in
)
::
name
(
*
)
integer
(
kind
=
c_int
),
intent
(
in
),
value
::
value
integer
(
kind
=
c_int
)
::
success
integer
(
kind
=
c_int
)