elpa.F90 16.1 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
!
!    Copyright 2017, L. Hüdepohl and A. Marek, MPCDF
!
!    This file is part of ELPA.
!
!    The ELPA library was originally created by the ELPA consortium,
!    consisting of the following organizations:
!
!    - Max Planck Computing and Data Facility (MPCDF), formerly known as
!      Rechenzentrum Garching der Max-Planck-Gesellschaft (RZG),
!    - 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,
!    - Max-Plack-Institut für Mathematik in den Naturwissenschaften,
!      Leipzig, Abt. Komplexe Strukutren in Biologie und Kognition,
!      and
!    - IBM Deutschland GmbH
!
!    This particular source code file contains additions, changes and
!    enhancements authored by Intel Corporation which is not part of
!    the ELPA consortium.
!
!    More information can be found here:
!    http://elpa.mpcdf.mpg.de/
!
!    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.
!

49
50
! The ELPA public API

51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73

!> \mainpage
!> Eigenvalue SoLvers for Petaflop-Applications (ELPA)
!> \par
!> http://elpa.mpcdf.mpg.de
!>
!> \par
!>    The ELPA library was originally created by the ELPA consortium,
!>    consisting of the following organizations:
!>
!>    - Max Planck Computing and Data Facility (MPCDF) formerly known as
!>      Rechenzentrum Garching der Max-Planck-Gesellschaft (RZG),
!>    - 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,
!>    - Max-Plack-Institut für Mathematik in den Naturwissenschaften,
!>      Leipzig, Abt. Komplexe Strukutren in Biologie und Kognition,
!>      and
!>    - IBM Deutschland GmbH
!>
!>   Some parts and enhancements of ELPA have been contributed and authored
Andreas Marek's avatar
Andreas Marek committed
74
!>   by the Intel Corporation and Nvidia Corporation, which are not part of
75
76
!>   the ELPA consortium.
!>
Andreas Marek's avatar
Andreas Marek committed
77
78
79
80
!>   Maintainance and development of the ELPA library is done by the
!>   Max Planck Computing and Data Facility (MPCDF)
!>
!>
Andreas Marek's avatar
Andreas Marek committed
81
!>   Futher support of the ELPA library is done by the ELPA-AEO consortium,
Andreas Marek's avatar
Andreas Marek committed
82
83
84
85
86
87
88
89
90
91
92
93
!>   consisting of the following organizations:
!>
!>    - Max Planck Computing and Data Facility (MPCDF) formerly known as
!>      Rechenzentrum Garching der Max-Planck-Gesellschaft (RZG),
!>    - Bergische Universität Wuppertal, Lehrstuhl für angewandte
!>      Informatik,
!>    - Technische Universität München, Lehrstuhl für Informatik mit
!>      Schwerpunkt Wissenschaftliches Rechnen ,
!>    - Technische Universität München, Lehrstuhl für theoretische Chemie,
!>    - Fritz-Haber-Institut, Berlin, Abt. Theorie
!>
!>
94
95
96
97
98
!>   Contributions to the ELPA source have been authored by (in alphabetical order):
!>
!> \author T. Auckenthaler, Volker Blum, A. Heinecke, L. Huedepohl, R. Johanni, Werner Jürgens, Pavel Kus, and A. Marek
!>
!> All the important information is in the \ref elpa_api::elpa_t derived type
Andreas Marek's avatar
Andreas Marek committed
99
100
101
102
103
104
105
106
107
108
!>
!> \brief Abstract definition of the elpa_t type
!>
!>
!> A typical usage of ELPA might look like this:
!>
!> Fortran synopsis
!>
!> \code{.f90}
!>  use elpa
109
!>  class(elpa_t), pointer :: elpaInstance
110
!>  integer :: success
Andreas Marek's avatar
Andreas Marek committed
111
!>
Andreas Marek's avatar
Andreas Marek committed
112
113
!>  ! We urge the user to always check the error code of all ELPA functions
!>
114
!>  if (elpa_init(20200417) /= ELPA_OK) then
Andreas Marek's avatar
Andreas Marek committed
115
116
117
!>     print *, "ELPA API version not supported"
!>     stop
!>   endif
Andreas Marek's avatar
Andreas Marek committed
118
!>   elpa => elpa_allocate(success)
119
!>   if (success /= ELPA_OK) then
Andreas Marek's avatar
Andreas Marek committed
120
121
!>     print *,"Could not allocate ELPA"
!>   endif
Andreas Marek's avatar
Andreas Marek committed
122
123
!>
!>   ! set parameters decribing the matrix and it's MPI distribution
124
!>   call elpaIstance%set("na", na, success, success)
125
!>   if (success /= ELPA_OK) then
Andreas Marek's avatar
Andreas Marek committed
126
127
!>     print *,"Could not set entry"
!>   endif
128
!>   call elpaInstance%set("nev", nev, success, success)
Andreas Marek's avatar
Andreas Marek committed
129
130
!>   ! check success code ...
!>
131
!>   call elpaInstance%set("local_nrows", na_rows, success)
Andreas Marek's avatar
Andreas Marek committed
132
133
!>   ! check success code ...
!>
134
135
136
137
138
!>   call elpaInstance%set("local_ncols", na_cols, success)
!>   call elpaInstance%set("nblk", nblk, success)
!>   call elpaInstance%set("mpi_comm_parent", MPI_COMM_WORLD, success)
!>   call elpaInstance%set("process_row", my_prow, success)
!>   call elpaInstance%set("process_col", my_pcol, success)
Andreas Marek's avatar
Andreas Marek committed
139
!>
140
!>   ! set up the elpa object
141
!>   success = elpaInstance%setup()
142
!>   if (succes /= ELPA_OK) then
Andreas Marek's avatar
Andreas Marek committed
143
144
!>     print *,"Could not setup ELPA object"
!>   endif
Andreas Marek's avatar
Andreas Marek committed
145
!>
146
147
148
149
150
151
152
153
154
155
!>   ! settings for GPU
!>   call elpaInstance%set("gpu", 1, success) ! 1=on, 2=off
!>   ! in case of GPU usage you have the choice whether ELPA
!>   ! should automatically assign each MPI task to a certain GPU
!>   ! (this is default) or whether you want to set this assignment
!>   ! for _each_ task yourself
!>   ! set assignment your self (only using one task here and assigning it 
!>   ! to GPU id 1)
!>   if (my_rank .eq. 0) call elpaInstance%set("use_gpu_id", 1, success)
!>
Andreas Marek's avatar
Andreas Marek committed
156
!>   ! if desired, set tunable run-time options
157
!>   ! here we want to use the 2-stage solver
158
!>   call elpaInstance%set("solver", ELPA_SOLVER_2STAGE, success)
159
160
!>
!>   ! and set a specific kernel (must be supported on the machine)
161
!>   call elpaInstance%set("real_kernel", ELPA_2STAGE_REAL_AVX_BLOCK2)
Andreas Marek's avatar
Andreas Marek committed
162
163
164
165
!> \endcode
!>   ... set and get all other options that are desired
!> \code{.f90}
!>
166
!>   ! if wanted you can store the settings and load them in another program
Andreas Marek's avatar
Andreas Marek committed
167
!>   call elpa%store_settings("save_to_disk.txt", success)
168
!>
169
170
!>   ! use method solve to solve the eigenvalue problem to obtain eigenvalues
!>   ! and eigenvectors
Andreas Marek's avatar
Andreas Marek committed
171
!>   ! other possible methods are desribed in \ref elpa_api::elpa_t derived type
172
!>   call elpaInstance%eigenvectors(a, ev, z, success)
Andreas Marek's avatar
Andreas Marek committed
173
174
!>
!>   ! cleanup
175
!>   call elpa_deallocate(e, success)
Andreas Marek's avatar
Andreas Marek committed
176
177
178
!>
!>   call elpa_uninit()
!> \endcode
179
180
181
182
183
184
!>
!>
!> C synopsis
!>
!>  \code{.c}
!>   #include <elpa/elpa.h>
185
!>
186
!>   elpa_t handle;
187
!>   int error;
188
!>
Andreas Marek's avatar
Andreas Marek committed
189
190
!>   /*  We urge the user to always check the error code of all ELPA functions */
!>
191
!>   if (elpa_init(20200417) != ELPA_OK) {
192
193
194
195
!>     fprintf(stderr, "Error: ELPA API version not supported");
!>     exit(1);
!>   }
!>
196
!>   
197
!>   handle = elpa_allocate(&error);
Andreas Marek's avatar
Andreas Marek committed
198
199
200
!>   if (error != ELPA_OK) {
!>   /* do sth. */
!>   }
201
202
203
204
205
206
207
208
209
210
211
212
!>
!>   /* Set parameters the matrix and it's MPI distribution */
!>   elpa_set(handle, "na", na, &error);
!>   elpa_set(handle, "nev", nev, &error);
!>   elpa_set(handle, "local_nrows", na_rows, &error);
!>   elpa_set(handle, "local_ncols", na_cols, &error);
!>   elpa_set(handle, "nblk", nblk, &error);
!>   elpa_set(handle, "mpi_comm_parent", MPI_Comm_c2f(MPI_COMM_WORLD), &error);
!>   elpa_set(handle, "process_row", my_prow, &error);
!>   elpa_set(handle, "process_col", my_pcol, &error);
!>
!>   /* Setup */
Andreas Marek's avatar
Andreas Marek committed
213
!>   error = elpa_setup(handle);
214
215
!>
!>   /* if desired, set tunable run-time options */
Andreas Marek's avatar
Andreas Marek committed
216
!>   /* here we want to use the 2-stage solver */
217
!>   elpa_set(handle, "solver", ELPA_SOLVER_2STAGE, &error);
218
!>
219
220
221
222
223
224
225
226
227
228
!>   /* settings for GPU */
!>   elpa_set(handle, "gpu", 1, &error);  /* 1=on, 2=off */
!>   /* in case of GPU usage you have the choice whether ELPA
!>      should automatically assign each MPI task to a certain GPU
!>      (this is default) or whether you want to set this assignment
!>      for _each_ task yourself
!>      set assignment your self (only using one task here and assigning it 
!>      to GPU id 1) */
!>   if (my_rank == 0) elpa_set(handle, "use_gpu_id", 1, &error);
!>
229
!>   elpa_set(handle,"real_kernel", ELPA_2STAGE_REAL_AVX_BLOCK2, &error);
230
231
232
233
!>  \endcode
!>   ... set and get all other options that are desired
!>  \code{.c}
!>
Andreas Marek's avatar
Andreas Marek committed
234
!>   /* if you want you can store the settings and load them in another program */
Andreas Marek's avatar
Andreas Marek committed
235
!>   elpa_store_settings(handle, "save_to_disk.txt");
236
!>
237
238
!>   /* use method solve to solve the eigenvalue problem */
!>   /* other possible methods are desribed in \ref elpa_api::elpa_t derived type */
239
!>   elpa_eigenvectors(handle, a, ev, z, &error);
240
241
!>
!>   /* cleanup */
242
!>   elpa_deallocate(handle, &error);
243
244
245
!>   elpa_uninit();
!> \endcode
!>
246
247
248
249
250
251
252
253
254
255
!> the autotuning could be used like this:
!>
!> Fortran synopsis
!>
!> \code{.f90}
!>  use elpa
!>  class(elpa_t), pointer :: elpa
!>  class(elpa_autotune_t), pointer :: tune_state
!>  integer :: success
!>
256
!>  if (elpa_init(20200417) /= ELPA_OK) then
257
258
259
!>     print *, "ELPA API version not supported"
!>     stop
!>   endif
Andreas Marek's avatar
Andreas Marek committed
260
!>   elpa => elpa_allocate(success)
261
262
263
264
265
266
267
268
269
270
271
272
!>
!>   ! set parameters decribing the matrix and it's MPI distribution
!>   call elpa%set("na", na, success)
!>   call elpa%set("nev", nev, success)
!>   call elpa%set("local_nrows", na_rows, success)
!>   call elpa%set("local_ncols", na_cols, success)
!>   call elpa%set("nblk", nblk, success)
!>   call elpa%set("mpi_comm_parent", MPI_COMM_WORLD, success)
!>   call elpa%set("process_row", my_prow, success)
!>   call elpa%set("process_col", my_pcol, success)
!>
!>   ! set up the elpa object
Andreas Marek's avatar
Andreas Marek committed
273
!>   success = elpa%setup()
274
275
!>
!>   ! create autotune object
Andreas Marek's avatar
Andreas Marek committed
276
!>   tune_state => elpa%autotune_setup(ELPA_AUTOTUNE_FAST, ELPA_AUTOTUNE_DOMAIN_REAL, success)
277
278
279
280
281
282
283
!>
!>   ! you can set some options, these will be then FIXED for the autotuning step
!>   ! if desired, set tunable run-time options
!>   ! here we want to use the 2-stage solver
!>   call e%set("solver", ELPA_SOLVER_2STAGE, success)
!>
!>   ! and set a specific kernel (must be supported on the machine)
Andreas Marek's avatar
Andreas Marek committed
284
!>   call e%set("real_kernel", ELPA_2STAGE_REAL_AVX_BLOCK2, success)
285
286
287
288
289
!> \endcode
!>   ... set and get all other options that are desired
!> \code{.f90}
!>
!>   iter = 0
Andreas Marek's avatar
Andreas Marek committed
290
!>   do while (elpa%autotune_step(tune_state, success))
291
292
293
294
295
296
!>     iter = iter + 1
!>     call e%eigenvectors(a, ev, z, success)
!>
!>     ! if needed you can save the autotune state at any point
!>     ! and resume it
!>     if (iter > MAX_ITER) then
Andreas Marek's avatar
Andreas Marek committed
297
!>       call elpa%autotune_save_state(tune_state,"autotune_checkpoint.txt", success)
298
299
300
301
302
!>       exit
!>     endif
!>   enddo
!>
!>   !set and print the finished autotuning
Andreas Marek's avatar
Andreas Marek committed
303
!>   call elpa%autotune_set_best(tune_state, success)
304
305
!>   
!>   ! store _TUNED_ ELPA object, if needed
Andreas Marek's avatar
Andreas Marek committed
306
!>   call elpa%store("autotuned_object.txt", success)
307
308
!>
!>   !deallocate autotune object
309
!>   call elpa_autotune_deallocate(tune_state, success)
310
311
!>
!>   ! cleanup
312
!>   call elpa_deallocate(e, success)
313
314
315
316
317
318
319
!>
!>   call elpa_uninit()
!> \endcode
!>
!> More examples can be found in the folder "test", where Fortran and C example programs
!> are stored

320
!> \brief Fortran module to use the ELPA library. No other module shoule be used
Andreas Marek's avatar
Andreas Marek committed
321
322

#include "config-f90.h"
323
324
module elpa
  use elpa_constants
325
  use elpa_api
326
327
328
329
330
331

  implicit none
  public

  contains

332
333
334
    !> \brief function to allocate an ELPA instance
    !> Parameters
    !> \details
335
    !> \params  error      integer, optional : error code
336
    !> \result  obj        class(elpa_t), pointer : pointer to allocated object
337
    function elpa_allocate(error) result(obj)
338
      use elpa_impl
339
340
341
342
343
344
345
346
347
      class(elpa_t), pointer         :: obj
#ifdef USE_FORTRAN2008
      integer, optional, intent(out) :: error
#else
      integer, intent(out)           :: error
#endif
      integer                        :: error2

      obj => elpa_impl_allocate(error2)
Andreas Marek's avatar
Andreas Marek committed
348
349

#ifdef USE_FORTRAN2008
350
      if (present(error)) then
Andreas Marek's avatar
Andreas Marek committed
351
#endif
352
353
354
355
356
357
        error = error2
        if (error .ne. ELPA_OK) then
          write(*,*) "Cannot allocate the ELPA object!"
          write(*,*) "This is a critical error!"
          write(*,*) "ELPA not usable with this error"
        endif
Andreas Marek's avatar
Andreas Marek committed
358
#ifdef USE_FORTRAN2008
359
360
361
362
363
364
365
366
      else
        if (error2 .ne. ELPA_OK) then
          write(*,*) "Cannot allocate the ELPA object!"
          write(*,*) "This is a critical error, but you do not check the error codes!"
          write(*,*) "ELPA not usable with this error"
          stop
        endif
      endif
Andreas Marek's avatar
Andreas Marek committed
367
368
#endif

369
370
    end function

371

372
373
374
    !> \brief function to deallocate an ELPA instance
    !> Parameters
    !> \details
375
    !> \param  obj        class(elpa_t), pointer : pointer to the ELPA object to be destroyed and deallocated
376
377
378
379
380
381
382
383
384
385
386
    !> \param  error      integer, optional : error code
    subroutine elpa_deallocate(obj, error)
      class(elpa_t), pointer         :: obj
#ifdef USE_FORTRAN2008
      integer, optional, intent(out) :: error
#else
      integer, intent(out)           :: error
#endif
      integer                        :: error2
        
      call obj%destroy(error2)
Andreas Marek's avatar
Andreas Marek committed
387
#ifdef USE_FORTRAN2008
388
      if (present(error)) then
Andreas Marek's avatar
Andreas Marek committed
389
#endif
390
391
392
393
394
395
396
397
        error = error2
        if (error .ne. ELPA_OK) then
          write(*,*) "Cannot destroy the ELPA object!"  
          write(*,*) "This is a critical error!"  
          write(*,*) "This might lead to a memory leak in your application!"
          error = ELPA_ERROR_CRITICAL
          return
        endif
Andreas Marek's avatar
Andreas Marek committed
398
#ifdef USE_FORTRAN2008
399
400
401
402
403
404
405
406
407
      else
        if (error2 .ne. ELPA_OK) then
          write(*,*) "Cannot destroy the ELPA object!"
          write(*,*) "This is a critical error!"
          write(*,*) "This might lead to a memory leak in your application!"
          write(*,*) "But you do not check the error codes!"
          return
        endif
      endif
Andreas Marek's avatar
Andreas Marek committed
408
#endif
409
410
411
412
413
      deallocate(obj, stat=error2)
      if (error2 .ne. 0) then
        write(*,*) "Cannot deallocate the ELPA object!"  
        write(*,*) "This is a critical error!"  
        write(*,*) "This might lead to a memory leak in your application!"
Andreas Marek's avatar
Andreas Marek committed
414
#ifdef USE_FORTRAN2008
415
416
417
418
        if (present(error)) then
          error = ELPA_ERROR_CRITICAL
          return
        endif
Andreas Marek's avatar
Andreas Marek committed
419
420
421
422
#else
        error = ELPA_ERROR_CRITICAL
        return
#endif
423
      endif
424
425
    end subroutine

Andreas Marek's avatar
Andreas Marek committed
426
#ifdef ENABLE_AUTOTUNING
427
428
429
    !> \brief function to deallocate an ELPA autotune instance
    !> Parameters
    !> \details
430
431
432
    !> \param  obj        class(elpa_autotune_t), pointer : pointer to the autotune object to be destroyed and deallocated
    !> \param  error      integer, optional : error code
    subroutine elpa_autotune_deallocate(obj, error)
433
      class(elpa_autotune_t), pointer :: obj
434
435
436
437
438
439
440
#ifdef USE_FORTRAN2008
      integer, optional, intent(out)  :: error
#else
      integer, intent(out)            :: error
#endif
      integer                         :: error2
      call obj%destroy(error2)
Andreas Marek's avatar
Andreas Marek committed
441
#ifdef USE_FORTRAN2008
442
      if (present(error)) then
Andreas Marek's avatar
Andreas Marek committed
443
#endif
444
445
446
447
448
449
450
451
        error = error2
        if (error2 .ne. ELPA_OK) then
          write(*,*) "Cannot destroy the ELPA autotuning object!"
          write(*,*) "This is a critical error!"
          write(*,*) "This might lead to a memory leak in your application!"
          error = ELPA_ERROR_CRITICAL
          return
        endif
Andreas Marek's avatar
Andreas Marek committed
452
#ifdef USE_FORTRAN2008
453
454
455
456
457
458
459
460
461
      else
        if (error2 .ne. ELPA_OK) then
          write(*,*) "Cannot destroy the ELPA autotuning object!"
          write(*,*) "This is a critical error!"
          write(*,*) "This might lead to a memory leak in your application!"
          write(*,*) "But you do not check the error codes"
          return
        endif
      endif
Andreas Marek's avatar
Andreas Marek committed
462
#endif
463
464
465
466
467
      deallocate(obj, stat=error2)
      if (error2 .ne. 0) then
        write(*,*) "Cannot deallocate the ELPA autotuning object!"  
        write(*,*) "This is a critical error!"  
        write(*,*) "This might lead to a memory leak in your application!"
Andreas Marek's avatar
Andreas Marek committed
468
#ifdef USE_FORTRAN2008
469
470
471
472
        if (present(error)) then
          error = ELPA_ERROR_CRITICAL
          return
        endif
Andreas Marek's avatar
Andreas Marek committed
473
474
475
476
#else
        error = ELPA_ERROR_CRITICAL
        return
#endif
477
478
      endif

479
    end subroutine
Andreas Marek's avatar
Andreas Marek committed
480
#endif
481

482
end module