Put routine single_hh_trafo in a module

This routine has been contained in a subroutine.
It has been moved to a module and and renamed to
"single_hh_trafo_real" to make it's intention more clear
parent 5c284727
......@@ -15,6 +15,7 @@ libelpa@SUFFIX@_la_SOURCES = src/mod_precision.f90 \
src/elpa1.F90 \
src/elpa2_utilities.F90 \
src/mod_pack_unpack_real.F90 \
src/elpa2_kernels/mod_single_hh_trafo_real.F90 \
src/mod_pack_unpack_complex.F90 \
src/elpa2_compute.F90 \
src/elpa2.F90 \
......
......@@ -2463,6 +2463,7 @@ module ELPA2_compute
#endif
use precision
use single_hh_trafo_real
#if defined(WITH_REAL_GENERIC_SIMPLE_KERNEL)
use real_generic_simple_kernel, only : double_hh_trafo_generic_simple
#endif
......@@ -2652,12 +2653,12 @@ module ELPA2_compute
!#endif
#ifdef WITH_OPENMP
if (j==1) call single_hh_trafo(a(1,1+off+a_off,istripe,my_thread), &
bcast_buffer(1,off+1), nbw, nl, &
if (j==1) call single_hh_trafo_real_cpu_openmp(a(1:stripe_width,1+off+a_off:1+off_a_off+nbw-1,istripe,my_thread), &
bcast_buffer(1:nbw,off+1), nbw, nl, &
stripe_width)
#else
if (j==1) call single_hh_trafo(a(1,1+off+a_off,istripe), &
bcast_buffer(1,off+1), nbw, nl, &
if (j==1) call single_hh_trafo_real_cpu(a(1:stripe_width,1+off+a_off:1+off+a_off+nbw-1,istripe), &
bcast_buffer(1:nbw,off+1), nbw, nl, &
stripe_width)
#endif
......@@ -2698,11 +2699,11 @@ module ELPA2_compute
#endif
enddo
#ifdef WITH_OPENMP
if (jj==1) call single_hh_trafo(a(1,1+off+a_off,istripe,my_thread), &
bcast_buffer(1,off+1), nbw, nl, stripe_width)
if (jj==1) call single_hh_trafo_real_cpu_openmp(a(1:stripe_width,1+off+a_off:1+off+a_off+nbw-1,istripe,my_thread), &
bcast_buffer(1:nbw,off+1), nbw, nl, stripe_width)
#else
if (jj==1) call single_hh_trafo(a(1,1+off+a_off,istripe), &
bcast_buffer(1,off+1), nbw, nl, stripe_width)
if (jj==1) call single_hh_trafo_real_cpu(a(1:stripe_width,1+off+a_off:1+off+a_off+nbw-1,istripe), &
bcast_buffer(1:nbw,off+1), nbw, nl, stripe_width)
#endif
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
endif
......@@ -2755,11 +2756,11 @@ module ELPA2_compute
#endif
enddo
#ifdef WITH_OPENMP
if (jjj==1) call single_hh_trafo(a(1,1+off+a_off,istripe,my_thread), &
bcast_buffer(1,off+1), nbw, nl, stripe_width)
if (jjj==1) call single_hh_trafo_real_cpu_openmp(a(1:stripe_width,1+off+a_off:1+off+a_off+nbw-1,istripe,my_thread), &
bcast_buffer(1:nbw,off+1), nbw, nl, stripe_width)
#else
if (jjj==1) call single_hh_trafo(a(1,1+off+a_off,istripe), &
bcast_buffer(1,off+1), nbw, nl, stripe_width)
if (jjj==1) call single_hh_trafo_real_cpu(a(1:stripe_width,1+off+a_off:1+off+a_off+nbw-1,istripe), &
bcast_buffer(1:nbw,off+1), nbw, nl, stripe_width)
#endif
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
endif
......
module single_hh_trafo_real
implicit none
#include "config-f90.h"
#ifdef WITH_OPENMP
public single_hh_trafo_real_cpu_openmp
#else
public single_hh_trafo_real_cpu
#endif
contains
#ifdef WITH_OPENMP
subroutine single_hh_trafo_real_cpu_openmp(q, hh, nb, nq, ldq)
#else
subroutine single_hh_trafo_real_cpu(q, hh, nb, nq, ldq)
#endif
#ifdef HAVE_DETAILED_TIMINGS
use timings
#endif
use precision
! Perform single real Householder transformation.
! This routine is not performance critical and thus it is coded here in Fortran
implicit none
integer(kind=ik), intent(in) :: nb, nq, ldq
! real(kind=rk), intent(inout) :: q(ldq, *)
! real(kind=rk), intent(in) :: hh(*)
real(kind=rk), intent(inout) :: q(1:ldq, 1:nb)
real(kind=rk), intent(in) :: hh(1:nb)
integer(kind=ik) :: i
real(kind=rk) :: v(nq)
#ifdef HAVE_DETAILED_TIMINGS
#ifdef WITH_OPENMP
call timer%start("single_hh_trafo_real_cpu_openmp")
#else
call timer%start("single_hh_trafo_real_cpu")
#endif
#endif
! v = q * hh
v(:) = q(1:nq,1)
do i=2,nb
v(:) = v(:) + q(1:nq,i) * hh(i)
enddo
! v = v * tau
v(:) = v(:) * hh(1)
! q = q - v * hh**T
q(1:nq,1) = q(1:nq,1) - v(:)
do i=2,nb
q(1:nq,i) = q(1:nq,i) - v(:) * hh(i)
enddo
#ifdef HAVE_DETAILED_TIMINGS
#ifdef WITH_OPENMP
call timer%stop("single_hh_trafo_real_cpu_openmp")
#else
call timer%stop("single_hh_trafo_real_cpu")
#endif
#endif
end subroutine
end module
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment