mod_compute_hh_trafo_real.F90 54.3 KB
Newer Older
1
2
3
4
5
!    This file is part of ELPA.
!
!    The ELPA library was originally created by the ELPA consortium,
!    consisting of the following organizations:
!
6
7
!    - Max Planck Computing and Data Facility (MPCDF), formerly known as
!      Rechenzentrum Garching der Max-Planck-Gesellschaft (RZG),
8
9
10
11
12
!    - Bergische Universität Wuppertal, Lehrstuhl für angewandte
!      Informatik,
!    - Technische Universität München, Lehrstuhl für Informatik mit
!      Schwerpunkt Wissenschaftliches Rechnen ,
!    - Fritz-Haber-Institut, Berlin, Abt. Theorie,
13
!    - Max-Plack-Institut für Mathematik in den Naturwissenschaften,
14
15
16
17
18
19
!      Leipzig, Abt. Komplexe Strukutren in Biologie und Kognition,
!      and
!    - IBM Deutschland GmbH
!
!
!    More information can be found here:
20
!    http://elpa.mpcdf.mpg.de/
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
!
!    ELPA is free software: you can redistribute it and/or modify
!    it under the terms of the version 3 of the license of the
!    GNU Lesser General Public License as published by the Free
!    Software Foundation.
!
!    ELPA is distributed in the hope that it will be useful,
!    but WITHOUT ANY WARRANTY; without even the implied warranty of
!    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
!    GNU Lesser General Public License for more details.
!
!    You should have received a copy of the GNU Lesser General Public License
!    along with ELPA.  If not, see <http://www.gnu.org/licenses/>
!
!    ELPA reflects a substantial effort on the part of the original
!    ELPA consortium, and we ask you to respect the spirit of the
!    license that we chose: i.e., please contribute any changes you
!    may have back to the original ELPA library distribution, and keep
!    any derivatives of ELPA under the same license that we chose for
!    the original distribution, the GNU Lesser General Public License.
!
! This file was written by A. Marek, MPCDF

44
45
module compute_hh_trafo_real
#include "config-f90.h"
46
  use elpa_mpi
47
48
49
  implicit none

#ifdef WITH_OPENMP
50
51
52
53
54
55
56
57
  public compute_hh_trafo_real_cpu_openmp_double
#else
  public compute_hh_trafo_real_cpu_double
#endif

#ifdef WANT_SINGLE_PRECISION_REAL
#ifdef WITH_OPENMP
  public compute_hh_trafo_real_cpu_openmp_single
58
#else
59
60
61
62
  public compute_hh_trafo_real_cpu_single
#endif


63
64
65
66
67
#endif

  contains

#ifdef WITH_OPENMP
68
       subroutine compute_hh_trafo_real_cpu_openmp_double(a, a_dev, stripe_width, a_dim2, stripe_count, max_threads, l_nev, &
69
                                                   a_off, nbw, max_blk_size, bcast_buffer, bcast_buffer_dev, hh_dot_dev,    &
70
                                                   hh_tau_dev, kernel_flops, kernel_time, off, ncols, istripe,              &
71
                                                   my_thread, thread_width, THIS_REAL_ELPA_KERNEL)
72
#else
73
       subroutine compute_hh_trafo_real_cpu_double       (a, a_dev, stripe_width, a_dim2, stripe_count,                       &
74
75
76
                                                   a_off, nbw, max_blk_size, bcast_buffer,  bcast_buffer_dev, hh_dot_dev,     &
                                                   hh_tau_dev, kernel_flops, kernel_time, off, ncols, istripe,                &
                                                   last_stripe_width, THIS_REAL_ELPA_KERNEL)
77
78
79
80
#endif


         use precision
81
         use iso_c_binding
82
83
         use elpa2_utilities
         use single_hh_trafo_real
84
85
86
         use cuda_c_kernel
         use cuda_functions

87
#if defined(WITH_REAL_GENERIC_SIMPLE_KERNEL) && !(defined(USE_ASSUMED_SIZE))
88
         use real_generic_simple_kernel !, only : double_hh_trafo_generic_simple
89
90
#endif

91
#if defined(WITH_REAL_GENERIC_KERNEL) && !(defined(USE_ASSUMED_SIZE))
92
        use real_generic_kernel !, only : double_hh_trafo_generic
93
#endif
94
95

#if defined(WITH_REAL_BGP_KERNEL)
96
         use real_bgp_kernel !, only : double_hh_trafo_bgp
97
98
99
#endif

#if defined(WITH_REAL_BGQ_KERNEL)
100
         use real_bgq_kernel !, only : double_hh_trafo_bgq
101
102
103
104
#endif
#ifdef HAVE_DETAILED_TIMINGS
         use timings
#endif
105

106
#if defined(HAVE_AVX) || defined(HAVE_AVX2) || defined(HAVE_SSE_INTRINSICS) || defined(HAVE_SSE_ASSEMBLY)
107
108
         use kernel_interfaces
#endif
109
         implicit none
110
         real(kind=c_double), intent(inout) :: kernel_time  ! MPI_WTIME always needs double
111
112
         integer(kind=lik)            :: kernel_flops
         integer(kind=ik), intent(in) :: nbw, max_blk_size
113
         real(kind=rk8)                :: bcast_buffer(nbw,max_blk_size)
114
115
116
117
118
119
         integer(kind=ik), intent(in) :: a_off

         integer(kind=ik), intent(in) :: stripe_width,a_dim2,stripe_count

#ifndef WITH_OPENMP
         integer(kind=ik), intent(in) :: last_stripe_width
120
!         real(kind=rk8)                :: a(stripe_width,a_dim2,stripe_count)
121
         real(kind=rk8), pointer      :: a(:,:,:)
122
#else
123
         integer(kind=ik), intent(in) :: max_threads, l_nev, thread_width
124
!         real(kind=rk8)                :: a(stripe_width,a_dim2,stripe_count,max_threads)
125
         real(kind=rk8), pointer      :: a(:,:,:,:)
126
127
128
#endif
         integer(kind=ik), intent(in) :: THIS_REAL_ELPA_KERNEL

129
130
131
132
133
         integer(kind=c_intptr_t)     :: a_dev
         integer(kind=c_intptr_t)     :: bcast_buffer_dev
         integer(kind=c_size_t)       :: dev_offset
         integer(kind=c_intptr_t)     :: hh_dot_dev
         integer(kind=c_intptr_t)     :: hh_tau_dev
134
135
136
137
138
139
         ! Private variables in OMP regions (my_thread) should better be in the argument list!
         integer(kind=ik)             :: off, ncols, istripe
#ifdef WITH_OPENMP
         integer(kind=ik)             :: my_thread, noff
#endif
         integer(kind=ik)             :: j, nl, jj, jjj
140
         real(kind=rk8)                :: w(nbw,6)
141
         real(kind=c_double)          :: ttt ! MPI_WTIME always needs double
142

143
144
145
146
147
         if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_GPU) then
           ! ncols - indicates the number of HH reflectors to apply; at least 1 must be available
           if (ncols < 1) return
         endif

148
149
#ifdef HAVE_DETAILED_TIMINGS
#ifdef WITH_OPENMP
150
         call timer%start("compute_hh_trafo_real_cpu_openmp_double")
151
#else
152
         call timer%start("compute_hh_trafo_real_cpu_double")
153
154
#endif
#endif
Andreas Marek's avatar
Andreas Marek committed
155
156
157
158
159
160
161
162
163

#ifdef WITH_OPENMP
         if (my_thread==1) then
#endif
           ttt = mpi_wtime()
#ifdef WITH_OPENMP
         endif
#endif

164
165
166
167

#ifndef WITH_OPENMP
         nl = merge(stripe_width, last_stripe_width, istripe<stripe_count)
#else
168
169
170
171
         if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_GPU) then
           print *,"compute_hh_trafo_real GPU OPENMP: not yet implemented"
           stop
         endif
172
173
174
175
176
177
178
179
180

         if (istripe<stripe_count) then
           nl = stripe_width
         else
           noff = (my_thread-1)*thread_width + (istripe-1)*stripe_width
           nl = min(my_thread*thread_width-noff, l_nev-noff)
           if (nl<=0) then
#ifdef HAVE_DETAILED_TIMINGS
#ifdef WITH_OPENMP
181
             call timer%stop("compute_hh_trafo_real_cpu_openmp_double")
182
#else
183
             call timer%stop("compute_hh_trafo_real_cpu_double")
184
185
186
187
188
#endif
#endif
             return
           endif
         endif
189
190
191
#endif /* not WITH_OPENMP */

         if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_GPU) then
192
193
           dev_offset = (0 + (a_off * stripe_width) + ( (istripe - 1) * stripe_width *a_dim2 )) *size_of_double_real_datatype
           call launch_compute_hh_trafo_c_kernel_real_double(a_dev + dev_offset, bcast_buffer_dev, hh_dot_dev, &
194
195
                                                      hh_tau_dev, nl, nbw, stripe_width, off, ncols)
         else ! not CUDA kernel
196
197
198

#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
         if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_AVX_BLOCK2 .or. &
199
             THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_AVX2_BLOCK2 .or. &
200
             THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_SSE_BLOCK2 .or. &
201
202
203
204
205
206
207
             THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_GENERIC    .or. &
             THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_GENERIC_SIMPLE .or. &
             THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_SSE .or.        &
             THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_BGP .or.        &
             THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_BGQ) then
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */

208
             !FORTRAN CODE / X86 INRINISIC CODE / BG ASSEMBLER USING 2 HOUSEHOLDER VECTORS
209
210
#if defined(WITH_REAL_GENERIC_KERNEL)
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
211
             if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_GENERIC) then
212
213
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */

214
215
216
               do j = ncols, 2, -2
                 w(:,1) = bcast_buffer(1:nbw,j+off)
                 w(:,2) = bcast_buffer(1:nbw,j+off-1)
217
218

#ifdef WITH_OPENMP
219
220

#ifdef USE_ASSUMED_SIZE
221
                 call double_hh_trafo_generic_double(a(1,j+off+a_off-1,istripe,my_thread), w, &
222
223
                                            nbw, nl, stripe_width, nbw)

224
#else
225
                 call double_hh_trafo_generic_double(a(1:stripe_width,j+off+a_off-1:j+off+a_off+nbw-1, &
226
                                              istripe,my_thread), w(1:nbw,1:6), &
227
                                              nbw, nl, stripe_width, nbw)
228
229
230
231
#endif

#else /* WITH_OPENMP */

232
#ifdef USE_ASSUMED_SIZE
233
                 call double_hh_trafo_generic_double(a(1,j+off+a_off-1,istripe),w, &
234
235
236
                                            nbw, nl, stripe_width, nbw)

#else
237
                 call double_hh_trafo_generic_double(a(1:stripe_width,j+off+a_off-1:j+off+a_off+nbw-1,istripe),w(1:nbw,1:6), &
238
                                            nbw, nl, stripe_width, nbw)
239
#endif
240
241
#endif /* WITH_OPENMP */

242
               enddo
243
244

#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
245
             endif
246
247
248
249
250
251
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
#endif /* WITH_REAL_GENERIC_KERNEL */


#if defined(WITH_REAL_GENERIC_SIMPLE_KERNEL)
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
252
             if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_GENERIC_SIMPLE) then
253
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
254
255
256
               do j = ncols, 2, -2
                 w(:,1) = bcast_buffer(1:nbw,j+off)
                 w(:,2) = bcast_buffer(1:nbw,j+off-1)
257
#ifdef WITH_OPENMP
258
259

#ifdef USE_ASSUMED_SIZE
260
                 call double_hh_trafo_generic_simple_double(a(1,j+off+a_off-1,istripe,my_thread), &
261
                                                       w, nbw, nl, stripe_width, nbw)
262
#else
263
                 call double_hh_trafo_generic_simple_double(a(1:stripe_width,j+off+a_off-1:j+off+a_off-1+nbw,istripe,my_thread), &
264
265
266
267
268
                                                     w, nbw, nl, stripe_width, nbw)

#endif

#else /* WITH_OPENMP */
269
270

#ifdef USE_ASSUMED_SIZE
271
                 call double_hh_trafo_generic_simple_double(a(1,j+off+a_off-1,istripe), &
272
                                                       w, nbw, nl, stripe_width, nbw)
273
#else
274
                 call double_hh_trafo_generic_simple_double(a(1:stripe_width,j+off+a_off-1:j+off+a_off-1+nbw,istripe), &
275
276
                                                     w, nbw, nl, stripe_width, nbw)

277
#endif
278
279
280

#endif /* WITH_OPENMP */

281
               enddo
282
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
283
             endif
284
285
286
287
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
#endif /* WITH_REAL_GENERIC_SIMPLE_KERNEL */


288
#if defined(WITH_REAL_SSE_ASSEMBLY_KERNEL)
289
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
290
             if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_SSE) then
291
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
292
293
294
               do j = ncols, 2, -2
                 w(:,1) = bcast_buffer(1:nbw,j+off)
                 w(:,2) = bcast_buffer(1:nbw,j+off-1)
295
#ifdef WITH_OPENMP
296
                 call double_hh_trafo_double(c_loc(a(1,j+off+a_off-1,istripe,my_thread)), w, nbw, nl, &
297
298
                                      stripe_width, nbw)
#else
299
                 call double_hh_trafo_double(c_loc(a(1,j+off+a_off-1,istripe)), w, nbw, nl, &
300
301
                                      stripe_width, nbw)
#endif
302
               enddo
303
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
304
             endif
305
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
306
#endif /* WITH_REAL_SSE_ASSEMBLY_KERNEL */
307

308
309
310
311
#if defined(WITH_REAL_SSE_BLOCK2_KERNEL)
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
           if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_SSE_BLOCK2) then
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
312
313

#if defined(WITH_NO_SPECIFIC_REAL_KERNEL) || (defined(WITH_ONE_SPECIFIC_REAL_KERNEL) && !defined(WITH_REAL_SSE_BLOCK6_KERNEL) && !defined(WITH_REAL_SSE_BLOCK4_KERNEL))
314
315
316
317
             do j = ncols, 2, -2
               w(:,1) = bcast_buffer(1:nbw,j+off)
               w(:,2) = bcast_buffer(1:nbw,j+off-1)
#ifdef WITH_OPENMP
318
               call double_hh_trafo_real_sse_2hv_double(c_loc(a(1,j+off+a_off-1,istripe,my_thread)), &
319
320
                                                       w, nbw, nl, stripe_width, nbw)
#else
321
               call double_hh_trafo_real_sse_2hv_double(c_loc(a(1,j+off+a_off-1,istripe)), &
322
323
324
                                                       w, nbw, nl, stripe_width, nbw)
#endif
             enddo
325
326
#endif /* defined(WITH_NO_SPECIFIC_REAL_KERNEL) || (defined(WITH_ONE_SPECIFIC_REAL_KERNEL) && !defined(WITH_REAL_SSE_BLOCK6_KERNEL) && !defined(WITH_REAL_SSE_BLOCK4_KERNEL)) */

327
328
329
330
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
           endif
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
#endif /* WITH_REAL_SSE_BLOCK2_KERNEL */
331

332
#if defined(WITH_REAL_AVX_BLOCK2_KERNEL) || defined(WITH_REAL_AVX2_BLOCK2_KERNEL)
333
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
334

335
336
           if ((THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_AVX_BLOCK2) .or. &
               (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_AVX2_BLOCK2))  then
337
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
338
339

#if defined(WITH_NO_SPECIFIC_REAL_KERNEL) || (defined(WITH_ONE_SPECIFIC_REAL_KERNEL) && !defined(WITH_REAL_AVX_BLOCK6_KERNEL) && !defined(WITH_REAL_AVX_BLOCK4_KERNEL) && !defined(WITH_REAL_AVX2_BLOCK6_KERNEL) && !defined(WITH_REAL_AVX2_BLOCK4_KERNEL))
340
341
342
               do j = ncols, 2, -2
                 w(:,1) = bcast_buffer(1:nbw,j+off)
                 w(:,2) = bcast_buffer(1:nbw,j+off-1)
343
#ifdef WITH_OPENMP
344

345
               call double_hh_trafo_real_avx_avx2_2hv_double(c_loc(a(1,j+off+a_off-1,istripe,my_thread)), &
346
347
                                                       w, nbw, nl, stripe_width, nbw)
#else
348
               call double_hh_trafo_real_avx_avx2_2hv_double(c_loc(a(1,j+off+a_off-1,istripe)), &
349
350
351
                                                       w, nbw, nl, stripe_width, nbw)
#endif
               enddo
352
353
#endif /* defined(WITH_NO_SPECIFIC_REAL_KERNEL) || (defined(WITH_ONE_SPECIFIC_REAL_KERNEL) ... */

354
355
356
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
             endif
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
357
#endif /* WITH_REAL_AVX_BLOCK2_KERNEL || WITH_REAL_AVX2_BLOCK2_KERNEL */
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400

#if defined(WITH_REAL_BGP_KERNEL)
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
             if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_BGP) then
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
               do j = ncols, 2, -2
                 w(:,1) = bcast_buffer(1:nbw,j+off)
                 w(:,2) = bcast_buffer(1:nbw,j+off-1)
#ifdef WITH_OPENMP
                 call double_hh_trafo_bgp_double(a(1,j+off+a_off-1,istripe,my_thread), w, nbw, nl, &
                                          stripe_width, nbw)
#else
                 call double_hh_trafo_bgp_double(a(1,j+off+a_off-1,istripe), w, nbw, nl, &
                                          stripe_width, nbw)
#endif
               enddo
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
             endif
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
#endif /* WITH_REAL_BGP_KERNEL */

#if defined(WITH_REAL_BGQ_KERNEL)
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
             if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_BGQ) then
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
               do j = ncols, 2, -2
                 w(:,1) = bcast_buffer(1:nbw,j+off)
                 w(:,2) = bcast_buffer(1:nbw,j+off-1)
#ifdef WITH_OPENMP
                 call double_hh_trafo_bgq_double(a(1,j+off+a_off-1,istripe,my_thread), w, nbw, nl, &
                                          stripe_width, nbw)
#else
                 call double_hh_trafo_bgq_double(a(1,j+off+a_off-1,istripe), w, nbw, nl, &
                                          stripe_width, nbw)
#endif
               enddo
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
             endif
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
#endif /* WITH_REAL_BGQ_KERNEL */


!#if defined(WITH_AVX_SANDYBRIDGE)
401
!              call double_hh_trafo_real_avx_avx2_2hv(a(1,j+off+a_off-1,istripe), w, nbw, nl, stripe_width, nbw)
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
!#endif

#ifdef WITH_OPENMP
             if (j==1) call single_hh_trafo_real_cpu_openmp_double(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_real_cpu_double(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 !
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */

420
421
422
423
#if defined(WITH_REAL_SSE_BLOCK4_KERNEL)
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
           if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_SSE_BLOCK4) then
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
424
425

#if defined(WITH_NO_SPECIFIC_REAL_KERNEL) || (defined(WITH_ONE_SPECIFIC_REAL_KERNEL) && !defined(WITH_REAL_SSE_BLOCK6_KERNEL))
426
427
428
429
430
431
432
             ! X86 INTRINSIC CODE, USING 4 HOUSEHOLDER VECTORS
             do j = ncols, 4, -4
               w(:,1) = bcast_buffer(1:nbw,j+off)
               w(:,2) = bcast_buffer(1:nbw,j+off-1)
               w(:,3) = bcast_buffer(1:nbw,j+off-2)
               w(:,4) = bcast_buffer(1:nbw,j+off-3)
#ifdef WITH_OPENMP
433
               call quad_hh_trafo_real_sse_4hv_double(c_loc(a(1,j+off+a_off-3,istripe,my_thread)), w, &
434
435
                                                  nbw, nl, stripe_width, nbw)
#else
436
               call quad_hh_trafo_real_sse_4hv_double(c_loc(a(1,j+off+a_off-3,istripe)), w, &
437
438
439
440
441
442
443
                                                  nbw, nl, stripe_width, nbw)
#endif
             enddo
             do jj = j, 2, -2
               w(:,1) = bcast_buffer(1:nbw,jj+off)
               w(:,2) = bcast_buffer(1:nbw,jj+off-1)
#ifdef WITH_OPENMP
444
               call double_hh_trafo_real_sse_2hv_double(c_loc(a(1,jj+off+a_off-1,istripe,my_thread)), &
445
446
                                                    w, nbw, nl, stripe_width, nbw)
#else
447
               call double_hh_trafo_real_sse_2hv_double(c_loc(a(1,jj+off+a_off-1,istripe)), &
448
449
450
451
                                                    w, nbw, nl, stripe_width, nbw)
#endif
             enddo
#ifdef WITH_OPENMP
452
453
             if (jj==1) call single_hh_trafo_real_cpu_openmp_double(a(1:stripe_width,1+off+a_off:1+off+a_off+nbw-1, &
                                                                     istripe,my_thread), &
454
455
456
457
458
                                          bcast_buffer(1:nbw,off+1), nbw, nl, stripe_width)
#else
             if (jj==1) call single_hh_trafo_real_cpu_double(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
459
460
461

#endif /* defined(WITH_NO_SPECIFIC_REAL_KERNEL) || (defined(WITH_ONE_SPECIFIC_REAL_KERNEL) && !defined(WITH_REAL_SSE_BLOCK6_KERNEL)) */

462
463
464
465
466
467
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
           endif
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
#endif /* WITH_REAL_SSE_BLOCK4_KERNEL */


468
#if defined(WITH_REAL_AVX_BLOCK4_KERNEL) || defined(WITH_REAL_AVX2_BLOCK4_KERNEL)
469
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
470
471
           if ((THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_AVX_BLOCK4) .or. &
               (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_AVX2_BLOCK4) ) then
472
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
473
474

#if defined(WITH_NO_SPECIFIC_REAL_KERNEL) || (defined(WITH_ONE_SPECIFIC_REAL_KERNEL) && !defined(WITH_REAL_AVX_BLOCK6_KERNEL) && !defined(WITH_REAL_AVX2_BLOCK6_KERNEL))
475
476
477
478
479
480
481
             ! X86 INTRINSIC CODE, USING 4 HOUSEHOLDER VECTORS
             do j = ncols, 4, -4
               w(:,1) = bcast_buffer(1:nbw,j+off)
               w(:,2) = bcast_buffer(1:nbw,j+off-1)
               w(:,3) = bcast_buffer(1:nbw,j+off-2)
               w(:,4) = bcast_buffer(1:nbw,j+off-3)
#ifdef WITH_OPENMP
482
               call quad_hh_trafo_real_avx_avx2_4hv_double(c_loc(a(1,j+off+a_off-3,istripe,my_thread)), w, &
483
484
                                                  nbw, nl, stripe_width, nbw)
#else
485
               call quad_hh_trafo_real_avx_avx2_4hv_double(c_loc(a(1,j+off+a_off-3,istripe)), w, &
486
487
488
489
490
491
492
                                                  nbw, nl, stripe_width, nbw)
#endif
             enddo
             do jj = j, 2, -2
               w(:,1) = bcast_buffer(1:nbw,jj+off)
               w(:,2) = bcast_buffer(1:nbw,jj+off-1)
#ifdef WITH_OPENMP
493
               call double_hh_trafo_real_avx_avx2_2hv_double(c_loc(a(1,jj+off+a_off-1,istripe,my_thread)), &
494
495
                                                    w, nbw, nl, stripe_width, nbw)
#else
496
               call double_hh_trafo_real_avx_avx2_2hv_double(c_loc(a(1,jj+off+a_off-1,istripe)), &
497
498
499
500
                                                    w, nbw, nl, stripe_width, nbw)
#endif
             enddo
#ifdef WITH_OPENMP
501
502
             if (jj==1) call single_hh_trafo_real_cpu_openmp_double(a(1:stripe_width,1+off+a_off:1+off+a_off+nbw-1, &
                                                                      istripe,my_thread), &
503
504
505
506
507
                                          bcast_buffer(1:nbw,off+1), nbw, nl, stripe_width)
#else
             if (jj==1) call single_hh_trafo_real_cpu_double(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
508
509
510

#endif /* defined(WITH_NO_SPECIFIC_REAL_KERNEL) || (defined(WITH_ONE_SPECIFIC_REAL_KERNEL) && !defined(WITH_REAL_AVX_BLOCK6_KERNEL) && !defined(WITH_REAL_AVX2_BLOCK6_KERNEL)) */

511
512
513
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
           endif
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
514
#endif /* WITH_REAL_AVX_BLOCK4_KERNEL || WITH_REAL_AVX2_BLOCK4_KERNEL */
515

516
517
518
519
520
521
522
523
524
525
526
527
528
#if defined(WITH_REAL_SSE_BLOCK6_KERNEL)
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
           if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_SSE_BLOCK6) then
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
             ! X86 INTRINSIC CODE, USING 6 HOUSEHOLDER VECTORS
             do j = ncols, 6, -6
               w(:,1) = bcast_buffer(1:nbw,j+off)
               w(:,2) = bcast_buffer(1:nbw,j+off-1)
               w(:,3) = bcast_buffer(1:nbw,j+off-2)
               w(:,4) = bcast_buffer(1:nbw,j+off-3)
               w(:,5) = bcast_buffer(1:nbw,j+off-4)
               w(:,6) = bcast_buffer(1:nbw,j+off-5)
#ifdef WITH_OPENMP
529
               call hexa_hh_trafo_real_sse_6hv_double(c_loc(a(1,j+off+a_off-5,istripe,my_thread)), w, &
530
531
                                                  nbw, nl, stripe_width, nbw)
#else
532
               call hexa_hh_trafo_real_sse_6hv_double(c_loc(a(1,j+off+a_off-5,istripe)), w, &
533
534
535
536
537
538
539
540
541
                                                  nbw, nl, stripe_width, nbw)
#endif
             enddo
             do jj = j, 4, -4
               w(:,1) = bcast_buffer(1:nbw,jj+off)
               w(:,2) = bcast_buffer(1:nbw,jj+off-1)
               w(:,3) = bcast_buffer(1:nbw,jj+off-2)
               w(:,4) = bcast_buffer(1:nbw,jj+off-3)
#ifdef WITH_OPENMP
542
               call quad_hh_trafo_real_sse_4hv_double(c_loc(a(1,jj+off+a_off-3,istripe,my_thread)), w, &
543
544
                                                  nbw, nl, stripe_width, nbw)
#else
545
               call quad_hh_trafo_real_sse_4hv_double(c_loc(a(1,jj+off+a_off-3,istripe)), w, &
546
547
548
549
550
551
552
                                                  nbw, nl, stripe_width, nbw)
#endif
             enddo
             do jjj = jj, 2, -2
               w(:,1) = bcast_buffer(1:nbw,jjj+off)
               w(:,2) = bcast_buffer(1:nbw,jjj+off-1)
#ifdef WITH_OPENMP
553
               call double_hh_trafo_real_sse_2hv_double(c_loc(a(1,jjj+off+a_off-1,istripe,my_thread)), &
554
555
                                                    w, nbw, nl, stripe_width, nbw)
#else
556
               call double_hh_trafo_real_sse_2hv_double(c_loc(a(1,jjj+off+a_off-1,istripe)), &
557
558
559
560
                                                    w, nbw, nl, stripe_width, nbw)
#endif
             enddo
#ifdef WITH_OPENMP
561
562
             if (jjj==1) call single_hh_trafo_real_cpu_openmp_double(a(1:stripe_width,1+off+a_off:1+off+a_off+nbw-1, &
                                                                        istripe,my_thread), &
563
564
565
566
567
568
569
570
571
572
                                           bcast_buffer(1:nbw,off+1), nbw, nl, stripe_width)
#else
             if (jjj==1) call single_hh_trafo_real_cpu_double(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
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
#endif /* WITH_REAL_SSE_BLOCK4_KERNEL */

573

574
#if defined(WITH_REAL_AVX_BLOCK6_KERNEL) || defined(WITH_REAL_AVX2_BLOCK6_KERNEL)
575
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
576
577
           if ((THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_AVX_BLOCK6) .or. &
               (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_AVX2_BLOCK6)) then
578
579
580
581
582
583
584
585
586
587
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
             ! X86 INTRINSIC CODE, USING 6 HOUSEHOLDER VECTORS
             do j = ncols, 6, -6
               w(:,1) = bcast_buffer(1:nbw,j+off)
               w(:,2) = bcast_buffer(1:nbw,j+off-1)
               w(:,3) = bcast_buffer(1:nbw,j+off-2)
               w(:,4) = bcast_buffer(1:nbw,j+off-3)
               w(:,5) = bcast_buffer(1:nbw,j+off-4)
               w(:,6) = bcast_buffer(1:nbw,j+off-5)
#ifdef WITH_OPENMP
588
               call hexa_hh_trafo_real_avx_avx2_6hv_double(c_loc(a(1,j+off+a_off-5,istripe,my_thread)), w, &
589
590
                                                  nbw, nl, stripe_width, nbw)
#else
591
               call hexa_hh_trafo_real_avx_avx2_6hv_double(c_loc(a(1,j+off+a_off-5,istripe)), w, &
592
593
594
595
596
597
598
599
600
                                                  nbw, nl, stripe_width, nbw)
#endif
             enddo
             do jj = j, 4, -4
               w(:,1) = bcast_buffer(1:nbw,jj+off)
               w(:,2) = bcast_buffer(1:nbw,jj+off-1)
               w(:,3) = bcast_buffer(1:nbw,jj+off-2)
               w(:,4) = bcast_buffer(1:nbw,jj+off-3)
#ifdef WITH_OPENMP
601
               call quad_hh_trafo_real_avx_avx2_4hv_double(c_loc(a(1,jj+off+a_off-3,istripe,my_thread)), w, &
602
603
                                                  nbw, nl, stripe_width, nbw)
#else
604
               call quad_hh_trafo_real_avx_avx2_4hv_double(c_loc(a(1,jj+off+a_off-3,istripe)), w, &
605
606
607
608
609
610
611
                                                  nbw, nl, stripe_width, nbw)
#endif
             enddo
             do jjj = jj, 2, -2
               w(:,1) = bcast_buffer(1:nbw,jjj+off)
               w(:,2) = bcast_buffer(1:nbw,jjj+off-1)
#ifdef WITH_OPENMP
612
               call double_hh_trafo_real_avx_avx2_2hv_double(c_loc(a(1,jjj+off+a_off-1,istripe,my_thread)), &
613
614
                                                    w, nbw, nl, stripe_width, nbw)
#else
615
               call double_hh_trafo_real_avx_avx2_2hv_double(c_loc(a(1,jjj+off+a_off-1,istripe)), &
616
617
618
619
                                                    w, nbw, nl, stripe_width, nbw)
#endif
             enddo
#ifdef WITH_OPENMP
620
621
             if (jjj==1) call single_hh_trafo_real_cpu_openmp_double(a(1:stripe_width,1+off+a_off:1+off+a_off+nbw-1, &
                                                                       istripe,my_thread), &
622
623
624
625
626
627
628
629
                                           bcast_buffer(1:nbw,off+1), nbw, nl, stripe_width)
#else
             if (jjj==1) call single_hh_trafo_real_cpu_double(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
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
630
#endif /* WITH_REAL_AVX_BLOCK6_KERNEL || WITH_REAL_AVX2_BLOCK6_KERNEL */
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678

         endif ! GPU_KERNEL

#ifdef WITH_OPENMP
         if (my_thread==1) then
#endif
           kernel_flops = kernel_flops + 4*int(nl,8)*int(ncols,8)*int(nbw,8)
           kernel_time = kernel_time + mpi_wtime()-ttt
#ifdef WITH_OPENMP
         endif
#endif
#ifdef HAVE_DETAILED_TIMINGS
#ifdef WITH_OPENMP
         call timer%stop("compute_hh_trafo_real_cpu_openmp_double")
#else
         call timer%stop("compute_hh_trafo_real_cpu_double")
#endif
#endif

#ifdef WITH_OPENMP
       end subroutine compute_hh_trafo_real_cpu_openmp_double
#else
       end subroutine compute_hh_trafo_real_cpu_double
#endif

#ifdef WANT_SINGLE_PRECISION_REAL
! SINGLE PRECISION IMPLEMENTATION AT THE MOMENT DUPLICATED !!!

#ifdef WITH_OPENMP
       subroutine compute_hh_trafo_real_cpu_openmp_single(a, a_dev, stripe_width, a_dim2, stripe_count, max_threads, l_nev, &
                                                   a_off, nbw, max_blk_size, bcast_buffer, bcast_buffer_dev, hh_dot_dev,    &
                                                   hh_tau_dev, kernel_flops, kernel_time, off, ncols, istripe,              &
                                                   my_thread, thread_width, THIS_REAL_ELPA_KERNEL)
#else
       subroutine compute_hh_trafo_real_cpu_single       (a, a_dev, stripe_width, a_dim2, stripe_count,                      &
                                                   a_off, nbw, max_blk_size, bcast_buffer,  bcast_buffer_dev, hh_dot_dev,     &
                                                   hh_tau_dev, kernel_flops, kernel_time, off, ncols, istripe,                &
                                                   last_stripe_width, THIS_REAL_ELPA_KERNEL)
#endif


         use precision
         use iso_c_binding
         use elpa2_utilities
         use single_hh_trafo_real
         use cuda_c_kernel
         use cuda_functions

679
#if defined(WITH_REAL_GENERIC_SIMPLE_KERNEL) && !(defined(USE_ASSUMED_SIZE))
680
681
682
         use real_generic_simple_kernel !, only : double_hh_trafo_generic_simple
#endif

683
#if defined(WITH_REAL_GENERIC_KERNEL) && !(defined(USE_ASSUMED_SIZE))
684
685
686
687
688
689
690
691
692
693
694
695
696
        use real_generic_kernel !, only : double_hh_trafo_generic
#endif

#if defined(WITH_REAL_BGP_KERNEL)
         use real_bgp_kernel !, only : double_hh_trafo_bgp
#endif

#if defined(WITH_REAL_BGQ_KERNEL)
         use real_bgq_kernel !, only : double_hh_trafo_bgq
#endif
#ifdef HAVE_DETAILED_TIMINGS
         use timings
#endif
697

698
#if defined(HAVE_AVX) || defined(HAVE_AVX2) || defined(HAVE_SSE_INTRINSICS) || defined(HAVE_SSE_ASSEMBLY)
699
700
         use kernel_interfaces
#endif
701
702
703
704
705
706
707
708
709
710
711
712
         implicit none
         real(kind=c_double), intent(inout) :: kernel_time  ! MPI_WTIME always needs double
         integer(kind=lik)            :: kernel_flops
         integer(kind=ik), intent(in) :: nbw, max_blk_size
         real(kind=rk4)                :: bcast_buffer(nbw,max_blk_size)
         integer(kind=ik), intent(in) :: a_off

         integer(kind=ik), intent(in) :: stripe_width,a_dim2,stripe_count

#ifndef WITH_OPENMP
         integer(kind=ik), intent(in) :: last_stripe_width
!         real(kind=rk4)                :: a(stripe_width,a_dim2,stripe_count)
713
         real(kind=rk4), pointer      :: a(:,:,:)
714
715
716
#else
         integer(kind=ik), intent(in) :: max_threads, l_nev, thread_width
!         real(kind=rk8)                :: a(stripe_width,a_dim2,stripe_count,max_threads)
717
         real(kind=rk4), pointer      :: a(:,:,:,:)
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
#endif
         integer(kind=ik), intent(in) :: THIS_REAL_ELPA_KERNEL

         integer(kind=c_intptr_t)     :: a_dev
         integer(kind=c_intptr_t)     :: bcast_buffer_dev
         integer(kind=c_size_t)       :: dev_offset
         integer(kind=c_intptr_t)     :: hh_dot_dev
         integer(kind=c_intptr_t)     :: hh_tau_dev
         ! Private variables in OMP regions (my_thread) should better be in the argument list!
         integer(kind=ik)             :: off, ncols, istripe
#ifdef WITH_OPENMP
         integer(kind=ik)             :: my_thread, noff
#endif
         integer(kind=ik)             :: j, nl, jj, jjj
         real(kind=rk4)                :: w(nbw,6)
         real(kind=c_double)          :: ttt ! MPI_WTIME always needs double

         if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_GPU) then
           ! ncols - indicates the number of HH reflectors to apply; at least 1 must be available
           if (ncols < 1) return
         endif

#ifdef HAVE_DETAILED_TIMINGS
#ifdef WITH_OPENMP
         call timer%start("compute_hh_trafo_real_cpu_openmp_single")
#else
         call timer%start("compute_hh_trafo_real_cpu_single")
#endif
#endif
         ttt = mpi_wtime()

#ifndef WITH_OPENMP
         nl = merge(stripe_width, last_stripe_width, istripe<stripe_count)
#else
         if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_GPU) then
           print *,"compute_hh_trafo_real GPU OPENMP: not yet implemented"
           stop
         endif

         if (istripe<stripe_count) then
           nl = stripe_width
         else
           noff = (my_thread-1)*thread_width + (istripe-1)*stripe_width
           nl = min(my_thread*thread_width-noff, l_nev-noff)
           if (nl<=0) then
#ifdef HAVE_DETAILED_TIMINGS
#ifdef WITH_OPENMP
             call timer%stop("compute_hh_trafo_real_cpu_openmp_single")
#else
             call timer%stop("compute_hh_trafo_real_cpu_single")
#endif
#endif
             return
           endif
         endif
#endif /* not WITH_OPENMP */

         if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_GPU) then
           dev_offset = (0 + (a_off * stripe_width) + ( (istripe - 1) * stripe_width *a_dim2 )) *size_of_single_real_datatype
           call launch_compute_hh_trafo_c_kernel_real_single(a_dev + dev_offset, bcast_buffer_dev, hh_dot_dev, &
                                                      hh_tau_dev, nl, nbw, stripe_width, off, ncols)
         else ! not CUDA kernel

#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
           if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_AVX_BLOCK2 .or. &
783
               THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_SSE_BLOCK2 .or. &
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
               THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_GENERIC    .or. &
               THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_GENERIC_SIMPLE .or. &
               THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_SSE .or.        &
               THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_BGP .or.        &
               THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_BGQ) then
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */

             !FORTRAN CODE / X86 INRINISIC CODE / BG ASSEMBLER USING 2 HOUSEHOLDER VECTORS
#if defined(WITH_REAL_GENERIC_KERNEL)
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
             if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_GENERIC) then
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */

               do j = ncols, 2, -2
                 w(:,1) = bcast_buffer(1:nbw,j+off)
                 w(:,2) = bcast_buffer(1:nbw,j+off-1)

#ifdef WITH_OPENMP
802
#ifdef USE_ASSUMED_SIZE
803
804
805
806
807
808
809
810
811
812
813
                 call double_hh_trafo_generic_single(a(1,j+off+a_off-1,istripe,my_thread), w, &
                                            nbw, nl, stripe_width, nbw)

#else
                 call double_hh_trafo_generic_single(a(1:stripe_width,j+off+a_off-1:j+off+a_off+nbw-1, &
                                              istripe,my_thread), w(1:nbw,1:6), &
                                              nbw, nl, stripe_width, nbw)
#endif

#else /* WITH_OPENMP */

814
#ifdef USE_ASSUMED_SIZE
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
                 call double_hh_trafo_generic_single(a(1,j+off+a_off-1,istripe),w, &
                                            nbw, nl, stripe_width, nbw)

#else
                 call double_hh_trafo_generic_single(a(1:stripe_width,j+off+a_off-1:j+off+a_off+nbw-1,istripe),w(1:nbw,1:6), &
                                            nbw, nl, stripe_width, nbw)
#endif
#endif /* WITH_OPENMP */

               enddo

#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
             endif
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
#endif /* WITH_REAL_GENERIC_KERNEL */


#if defined(WITH_REAL_GENERIC_SIMPLE_KERNEL)
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
             if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_GENERIC_SIMPLE) then
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
               do j = ncols, 2, -2
                 w(:,1) = bcast_buffer(1:nbw,j+off)
                 w(:,2) = bcast_buffer(1:nbw,j+off-1)
#ifdef WITH_OPENMP
840
#ifdef USE_ASSUMED_SIZE
841
842
843
844
845
846
847
848
849
                 call double_hh_trafo_generic_simple_single(a(1,j+off+a_off-1,istripe,my_thread), &
                                                       w, nbw, nl, stripe_width, nbw)
#else
                 call double_hh_trafo_generic_simple_single(a(1:stripe_width,j+off+a_off-1:j+off+a_off-1+nbw,istripe,my_thread), &
                                                     w, nbw, nl, stripe_width, nbw)

#endif

#else /* WITH_OPENMP */
850
#ifdef USE_ASSUMED_SIZE
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
                 call double_hh_trafo_generic_simple_single(a(1,j+off+a_off-1,istripe), &
                                                       w, nbw, nl, stripe_width, nbw)
#else
                 call double_hh_trafo_generic_simple_single(a(1:stripe_width,j+off+a_off-1:j+off+a_off-1+nbw,istripe), &
                                                     w, nbw, nl, stripe_width, nbw)

#endif

#endif /* WITH_OPENMP */

               enddo
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
             endif
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
#endif /* WITH_REAL_GENERIC_SIMPLE_KERNEL */


868
#if defined(WITH_REAL_SSE_ASSEMBLY_KERNEL)
869
870
871
872
873
874
875
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
             if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_SSE) then
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
               do j = ncols, 2, -2
                 w(:,1) = bcast_buffer(1:nbw,j+off)
                 w(:,2) = bcast_buffer(1:nbw,j+off-1)
#ifdef WITH_OPENMP
876
                 call double_hh_trafo_single(c_loc(a(1,j+off+a_off-1,istripe,my_thread)), w, nbw, nl, &
877
                                      stripe_width, nbw)
878
#else
879
                 call double_hh_trafo_single(c_loc(a(1,j+off+a_off-1,istripe)), w, nbw, nl, &
880
                                      stripe_width, nbw)
881
882
883
884
885
#endif
               enddo
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
             endif
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
886
#endif /* WITH_REAL_SSE_ASSEMBLY_KERNEL */
887

888
889
890
891
#if defined(WITH_REAL_SSE_BLOCK2_KERNEL)
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
             if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_SSE_BLOCK2) then
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
892
893

#if defined(WITH_NO_SPECIFIC_REAL_KERNEL) || (defined(WITH_ONE_SPECIFIC_REAL_KERNEL) && !defined(WITH_REAL_SSE_BLOCK6_KERNEL) && !defined(WITH_REAL_SSE_BLOCK4_KERNEL)) 
894
895
896
897
               do j = ncols, 2, -2
                 w(:,1) = bcast_buffer(1:nbw,j+off)
                 w(:,2) = bcast_buffer(1:nbw,j+off-1)
#ifdef WITH_OPENMP
898
                 call double_hh_trafo_real_sse_2hv_single(c_loc(a(1,j+off+a_off-1,istripe,my_thread)), &
899
900
                                                       w, nbw, nl, stripe_width, nbw)
#else
901
                 call double_hh_trafo_real_sse_2hv_single(c_loc(a(1,j+off+a_off-1,istripe)), &
902
903
904
                                                       w, nbw, nl, stripe_width, nbw)
#endif
               enddo
905
906
#endif /* defined(WITH_NO_SPECIFIC_REAL_KERNEL) || (defined(WITH_ONE_SPECIFIC_REAL_KERNEL) && !defined(WITH_REAL_SSE_BLOCK6_KERNEL) && !defined(WITH_REAL_SSE_BLOCK4_KERNEL)) */

907
908
909
910
911
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
             endif
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
#endif /* WITH_REAL_SSE_BLOCK2_KERNEL */

Andreas Marek's avatar
Typo  
Andreas Marek committed
912
#if defined(WITH_REAL_AVX_BLOCK2_KERNEL) || defined(WITH_REAL_AVX2_BLOCK2_KERNEL)
913
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
Andreas Marek's avatar
Andreas Marek committed
914
915
             if ((THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_AVX_BLOCK2) .or.  &
                 (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_AVX2_BLOCK2))then
916
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
917
918

#if defined(WITH_NO_SPECIFIC_REAL_KERNEL) || (defined(WITH_ONE_SPECIFIC_REAL_KERNEL) && !defined(WITH_REAL_AVX_BLOCK6_KERNEL) && !defined(WITH_REAL_AVX2_BLOCK6_KERNEL) && !defined(WITH_REAL_AVX_BLOCK4_KERNEL) && !defined(WITH_REAL_AVX2_BLOCK4_KERNEL) )
919
920
921
922
               do j = ncols, 2, -2
                 w(:,1) = bcast_buffer(1:nbw,j+off)
                 w(:,2) = bcast_buffer(1:nbw,j+off-1)
#ifdef WITH_OPENMP
923
                 call double_hh_trafo_real_avx_avx2_2hv_single(c_loc(a(1,j+off+a_off-1,istripe,my_thread)), &
924
925
                                                       w, nbw, nl, stripe_width, nbw)
#else
926
                 call double_hh_trafo_real_avx_avx2_2hv_single(c_loc(a(1,j+off+a_off-1,istripe)), &
927
928
                                                       w, nbw, nl, stripe_width, nbw)
#endif
929
               enddo
930
931
932

#endif /* defined(WITH_NO_SPECIFIC_REAL_KERNEL) || (defined(WITH_ONE_SPECIFIC_REAL_KERNEL) && !defined(WITH_REAL_AVX_BLOCK6_KERNEL) && !defined(WITH_REAL_AVX2_BLOCK6_KERNEL) && !defined(WITH_REAL_AVX_BLOCK4_KERNEL) && !defined(WITH_REAL_AVX2_BLOCK4_KERNEL) ) */

933
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
934
             endif
935
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
Andreas Marek's avatar
Andreas Marek committed
936
#endif /* WITH_REAL_AVX_BLOCK2_KERNEL || WITH_REAL_AVX2_BLOCK2_KERNEL */
937
938
939

#if defined(WITH_REAL_BGP_KERNEL)
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
940
             if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_BGP) then
941
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
942
943
944
               do j = ncols, 2, -2
                 w(:,1) = bcast_buffer(1:nbw,j+off)
                 w(:,2) = bcast_buffer(1:nbw,j+off-1)
945
#ifdef WITH_OPENMP
946
                 call double_hh_trafo_bgp_single(a(1,j+off+a_off-1,istripe,my_thread), w, nbw, nl, &
947
948
                                          stripe_width, nbw)
#else
949
                 call double_hh_trafo_bgp_single(a(1,j+off+a_off-1,istripe), w, nbw, nl, &
950
951
                                          stripe_width, nbw)
#endif
952
               enddo
953
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
954
             endif
955
956
957
958
959
960
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
#endif /* WITH_REAL_BGP_KERNEL */


#if defined(WITH_REAL_BGQ_KERNEL)
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
961
             if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_BGQ) then
962
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
963
964
965
               do j = ncols, 2, -2
                 w(:,1) = bcast_buffer(1:nbw,j+off)
                 w(:,2) = bcast_buffer(1:nbw,j+off-1)
966
#ifdef WITH_OPENMP
967
                 call double_hh_trafo_bgq_single(a(1,j+off+a_off-1,istripe,my_thread), w, nbw, nl, &
968
969
                                          stripe_width, nbw)
#else
970
                 call double_hh_trafo_bgq_single(a(1,j+off+a_off-1,istripe), w, nbw, nl, &
971
972
                                          stripe_width, nbw)
#endif
973
               enddo
974
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
975
             endif
976
977
978
979
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
#endif /* WITH_REAL_BGQ_KERNEL */

#ifdef WITH_OPENMP
980
981
             if (j==1) call single_hh_trafo_real_cpu_openmp_single(a(1:stripe_width, &
                                   1+off+a_off:1+off+a_off+nbw-1,istripe,my_thread), &
982
983
984
                                      bcast_buffer(1:nbw,off+1), nbw, nl,     &
                                      stripe_width)
#else
985
986
             if (j==1) call single_hh_trafo_real_cpu_single(a(1:stripe_width, &
                                    1+off+a_off:1+off+a_off+nbw-1,istripe),           &
987
988
989
990
991
992
                                      bcast_buffer(1:nbw,off+1), nbw, nl,     &
                                      stripe_width)
#endif


#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
993
           endif !
994
995
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */

996
#if defined(WITH_REAL_SSE_BLOCK4_KERNEL)
997
#if defined(WITH_NO_SPECIFIC_REAL_KERNEL)
998
999
         if (THIS_REAL_ELPA_KERNEL .eq. REAL_ELPA_KERNEL_SSE_BLOCK4) then
#endif /* WITH_NO_SPECIFIC_REAL_KERNEL */
1000
1001
1002

#if defined(WITH_NO_SPECIFIC_REAL_KERNEL) || (defined(WITH_ONE_SPECIFIC_REAL_KERNEL) && !defined(WITH_REAL_SSE_BLOCK6_KERNEL) )

1003
1004
1005
1006
1007
1008
1009
           ! X86 INTRINSIC CODE, USING 4 HOUSEHOLDER VECTORS
           do j = ncols, 4, -4
             w(:,1) = bcast_buffer(1:nbw,j+off)
             w(:,2) = bcast_buffer(1:nbw,j+off-1)
             w(:,3) = bcast_buffer(1:nbw,j+off-2)
             w(:,4) = bcast_buffer(1:nbw,j+off-3)
#ifdef WITH_OPENMP
1010
             call quad_hh_trafo_real_sse_4hv_single(c_loc(a(1,j+off+a_off-3,istripe,my_thread)), w, &
1011
1012
                                                  nbw, nl, stripe_width, nbw)
#else
1013
             call quad_hh_trafo_real_sse_4hv_single(c_loc(a(1,j+off+a_off-3,istripe)), w, &
1014
1015
1016
1017
1018
1019
1020
                                                  nbw, nl, stripe_width, nbw)
#endif
           enddo
           do jj = j, 2, -2
             w(:,1) = bcast_buffer(1:nbw,jj+off)
             w(:,2) = bcast_buffer(1:nbw,jj+off-1)
#ifdef WITH_OPENMP
1021
             call double_hh_trafo_real_sse_2hv_single(c_loc(a(1,jj+off+a_off-1,istripe,my_thread)), &
1022
1023
                                                    w, nbw, nl, stripe_width, nbw)
#else
1024
             call double_hh_trafo_real_sse_2hv_single(c_loc(a(1,jj+off+a_off-1,istripe)), &
1025
1026
1027
1028
                                                    w, nbw, nl, stripe_width, n