Line data Source code
1 : !
2 : ! Copyright 2017, L. Hüdepohl and A. Marek, MPCDF
3 : !
4 : ! This file is part of ELPA.
5 : !
6 : ! The ELPA library was originally created by the ELPA consortium,
7 : ! consisting of the following organizations:
8 : !
9 : ! - Max Planck Computing and Data Facility (MPCDF), formerly known as
10 : ! Rechenzentrum Garching der Max-Planck-Gesellschaft (RZG),
11 : ! - Bergische Universität Wuppertal, Lehrstuhl für angewandte
12 : ! Informatik,
13 : ! - Technische Universität München, Lehrstuhl für Informatik mit
14 : ! Schwerpunkt Wissenschaftliches Rechnen ,
15 : ! - Fritz-Haber-Institut, Berlin, Abt. Theorie,
16 : ! - Max-Plack-Institut für Mathematik in den Naturwissenschaften,
17 : ! Leipzig, Abt. Komplexe Strukutren in Biologie und Kognition,
18 : ! and
19 : ! - IBM Deutschland GmbH
20 : !
21 : ! This particular source code file contains additions, changes and
22 : ! enhancements authored by Intel Corporation which is not part of
23 : ! the ELPA consortium.
24 : !
25 : ! More information can be found here:
26 : ! http://elpa.mpcdf.mpg.de/
27 : !
28 : ! ELPA is free software: you can redistribute it and/or modify
29 : ! it under the terms of the version 3 of the license of the
30 : ! GNU Lesser General Public License as published by the Free
31 : ! Software Foundation.
32 : !
33 : ! ELPA is distributed in the hope that it will be useful,
34 : ! but WITHOUT ANY WARRANTY; without even the implied warranty of
35 : ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
36 : ! GNU Lesser General Public License for more details.
37 : !
38 : ! You should have received a copy of the GNU Lesser General Public License
39 : ! along with ELPA. If not, see <http://www.gnu.org/licenses/>
40 : !
41 : ! ELPA reflects a substantial effort on the part of the original
42 : ! ELPA consortium, and we ask you to respect the spirit of the
43 : ! license that we chose: i.e., please contribute any changes you
44 : ! may have back to the original ELPA library distribution, and keep
45 : ! any derivatives of ELPA under the same license that we chose for
46 : ! the original distribution, the GNU Lesser General Public License.
47 : !
48 : #include "config-f90.h"
49 :
50 : !> \brief Fortran module which provides the actual implementation of the API. Do not use directly! Use the module "elpa"
51 : module elpa_impl
52 : use precision
53 : use elpa2_impl
54 : use elpa1_impl
55 : use elpa1_auxiliary_impl
56 : #ifdef WITH_MPI
57 : use elpa_mpi
58 : #endif
59 : use elpa_generated_fortran_interfaces
60 : use elpa_utilities, only : error_unit
61 :
62 : use elpa_abstract_impl
63 : use elpa_autotune_impl
64 : use, intrinsic :: iso_c_binding
65 : implicit none
66 :
67 : private
68 : public :: elpa_impl_allocate
69 :
70 : !> \brief Definition of the extended elpa_impl_t type
71 : type, extends(elpa_abstract_impl_t) :: elpa_impl_t
72 : private
73 : integer :: communicators_owned
74 :
75 : !> \brief methods available with the elpa_impl_t type
76 : contains
77 : !> \brief the puplic methods
78 : ! con-/destructor
79 : procedure, public :: setup => elpa_setup !< a setup method: implemented in elpa_setup
80 : procedure, public :: destroy => elpa_destroy !< a destroy method: implemented in elpa_destroy
81 :
82 : ! KV store
83 : procedure, public :: is_set => elpa_is_set !< a method to check whether a key/value pair has been set : implemented
84 : !< in elpa_is_set
85 : procedure, public :: can_set => elpa_can_set !< a method to check whether a key/value pair can be set : implemented
86 : !< in elpa_can_set
87 :
88 :
89 : ! timer
90 : procedure, public :: get_time => elpa_get_time
91 : procedure, public :: print_times => elpa_print_times
92 : procedure, public :: timer_start => elpa_timer_start
93 : procedure, public :: timer_stop => elpa_timer_stop
94 :
95 :
96 : !> \brief the implemenation methods
97 :
98 : procedure, public :: elpa_eigenvectors_d !< public methods to implement the solve step for real/complex
99 : !< double/single matrices
100 : procedure, public :: elpa_eigenvectors_f
101 : procedure, public :: elpa_eigenvectors_dc
102 : procedure, public :: elpa_eigenvectors_fc
103 :
104 : procedure, public :: elpa_eigenvalues_d !< public methods to implement the solve step for real/complex
105 : !< double/single matrices; only the eigenvalues are computed
106 : procedure, public :: elpa_eigenvalues_f
107 : procedure, public :: elpa_eigenvalues_dc
108 : procedure, public :: elpa_eigenvalues_fc
109 :
110 : #if 0
111 : procedure, public :: elpa_generalized_eigenvectors_d !< public methods to implement the solve step for generalized
112 : !< eigenproblem and real/complex double/single matrices
113 : procedure, public :: elpa_generalized_eigenvectors_f
114 : procedure, public :: elpa_generalized_eigenvectors_dc
115 : procedure, public :: elpa_generalized_eigenvectors_fc
116 : #endif
117 :
118 : procedure, public :: elpa_hermitian_multiply_d !< public methods to implement a "hermitian" multiplication of matrices a and b
119 : procedure, public :: elpa_hermitian_multiply_f !< for real valued matrices: a**T * b
120 : procedure, public :: elpa_hermitian_multiply_dc !< for complex valued matrices: a**H * b
121 : procedure, public :: elpa_hermitian_multiply_fc
122 :
123 : procedure, public :: elpa_cholesky_d !< public methods to implement the cholesky factorisation of
124 : !< real/complex double/single matrices
125 : procedure, public :: elpa_cholesky_f
126 : procedure, public :: elpa_cholesky_dc
127 : procedure, public :: elpa_cholesky_fc
128 :
129 : procedure, public :: elpa_invert_trm_d !< public methods to implement the inversion of a triangular
130 : !< real/complex double/single matrix
131 : procedure, public :: elpa_invert_trm_f
132 : procedure, public :: elpa_invert_trm_dc
133 : procedure, public :: elpa_invert_trm_fc
134 :
135 : procedure, public :: elpa_solve_tridiagonal_d !< public methods to implement the solve step for a real valued
136 : procedure, public :: elpa_solve_tridiagonal_f !< double/single tridiagonal matrix
137 :
138 : procedure, public :: associate_int => elpa_associate_int !< public method to set some pointers
139 :
140 : #if 0
141 : procedure, private :: elpa_transform_generalized_d
142 : procedure, private :: elpa_transform_back_generalized_d
143 : procedure, private :: elpa_transform_generalized_dc
144 : procedure, private :: elpa_transform_back_generalized_dc
145 : #ifdef WANT_SINGLE_PRECISION_REAL
146 : procedure, private :: elpa_transform_generalized_f
147 : procedure, private :: elpa_transform_back_generalized_f
148 : #endif
149 : #ifdef WANT_SINGLE_PRECISION_COMPLEX
150 : procedure, private :: elpa_transform_generalized_fc
151 : procedure, private :: elpa_transform_back_generalized_fc
152 : #endif
153 : #endif
154 :
155 : procedure, public :: autotune_setup => elpa_autotune_setup
156 : procedure, public :: autotune_step => elpa_autotune_step
157 : procedure, public :: autotune_set_best => elpa_autotune_set_best
158 :
159 : end type elpa_impl_t
160 :
161 : !> \brief the implementation of the generic methods
162 : contains
163 :
164 :
165 : !> \brief function to allocate an ELPA object
166 : !> Parameters
167 : !> \param error integer, optional to get an error code
168 : !> \result obj class(elpa_impl_t) allocated ELPA object
169 21312 : function elpa_impl_allocate(error) result(obj)
170 : type(elpa_impl_t), pointer :: obj
171 : integer, optional :: error
172 :
173 21312 : allocate(obj)
174 :
175 : ! check whether init has ever been called
176 21312 : if ( elpa_initialized() .ne. ELPA_OK) then
177 0 : write(error_unit, *) "elpa_allocate(): you must call elpa_init() once before creating instances of ELPA"
178 0 : if(present(error)) then
179 0 : error = ELPA_ERROR
180 : endif
181 0 : return
182 : endif
183 :
184 21312 : obj%index = elpa_index_instance_c()
185 :
186 : ! Associate some important integer pointers for convenience
187 21312 : obj%na => obj%associate_int("na")
188 21312 : obj%nev => obj%associate_int("nev")
189 21312 : obj%local_nrows => obj%associate_int("local_nrows")
190 21312 : obj%local_ncols => obj%associate_int("local_ncols")
191 21312 : obj%nblk => obj%associate_int("nblk")
192 :
193 21312 : if(present(error)) then
194 2016 : error = ELPA_OK
195 : endif
196 21312 : end function
197 :
198 : !c> /*! \brief C interface for the implementation of the elpa_allocate method
199 : !c> *
200 : !c> * \param none
201 : !c> * \result elpa_t handle
202 : !c> */
203 : !c> elpa_t elpa_allocate();
204 2016 : function elpa_impl_allocate_c(error) result(ptr) bind(C, name="elpa_allocate")
205 : integer(kind=c_int) :: error
206 : type(c_ptr) :: ptr
207 : type(elpa_impl_t), pointer :: obj
208 :
209 2016 : obj => elpa_impl_allocate(error)
210 2016 : ptr = c_loc(obj)
211 2016 : end function
212 :
213 : !c> /*! \brief C interface for the implementation of the elpa_deallocate method
214 : !c> *
215 : !c> * \param elpa_t handle of ELPA object to be deallocated
216 : !c> * \result void
217 : !c> */
218 : !c> void elpa_deallocate(elpa_t handle);
219 2016 : subroutine elpa_impl_deallocate_c(handle) bind(C, name="elpa_deallocate")
220 : type(c_ptr), value :: handle
221 : type(elpa_impl_t), pointer :: self
222 :
223 2016 : call c_f_pointer(handle, self)
224 2016 : call self%destroy()
225 2016 : deallocate(self)
226 2016 : end subroutine
227 :
228 :
229 : !c> /*! \brief C interface for the implementation of the elpa_autotune_deallocate method
230 : !c> *
231 : !c> * \param elpa_autotune_impl_t handle of ELPA autotune object to be deallocated
232 : !c> * \result void
233 : !c> */
234 : !c> void elpa_autotune_deallocate(elpa_autotune_t handle);
235 0 : subroutine elpa_autotune_impl_deallocate_c( autotune_handle) bind(C, name="elpa_autotune_deallocate")
236 : type(c_ptr), value :: autotune_handle
237 :
238 : type(elpa_autotune_impl_t), pointer :: self
239 :
240 0 : call c_f_pointer(autotune_handle, self)
241 0 : call self%destroy()
242 0 : deallocate(self)
243 0 : end subroutine
244 :
245 :
246 : !> \brief function to setup an ELPA object and to store the MPI communicators internally
247 : !> Parameters
248 : !> \param self class(elpa_impl_t), the allocated ELPA object
249 : !> \result error integer, the error code
250 21312 : function elpa_setup(self) result(error)
251 : class(elpa_impl_t), intent(inout) :: self
252 : integer :: error, timings
253 :
254 : #ifdef WITH_MPI
255 : integer :: mpi_comm_parent, mpi_comm_rows, mpi_comm_cols, &
256 : mpierr, mpierr2, process_row, process_col, mpi_string_length
257 : character(len=MPI_MAX_ERROR_STRING) :: mpierr_string
258 : #endif
259 :
260 : #ifdef HAVE_DETAILED_TIMINGS
261 21312 : call self%get("timings",timings, error)
262 21312 : if (timings == 1) then
263 15840 : call self%timer%enable()
264 : endif
265 : #endif
266 :
267 21312 : error = ELPA_OK
268 :
269 : #ifdef WITH_MPI
270 : ! Create communicators ourselves
271 : if (self%is_set("mpi_comm_parent") == 1 .and. &
272 14336 : self%is_set("process_row") == 1 .and. &
273 : self%is_set("process_col") == 1) then
274 :
275 8320 : call self%get("mpi_comm_parent", mpi_comm_parent, error)
276 8320 : call self%get("process_row", process_row, error)
277 8320 : call self%get("process_col", process_col, error)
278 :
279 : ! mpi_comm_rows is used for communicating WITHIN rows, i.e. all processes
280 : ! having the same column coordinate share one mpi_comm_rows.
281 : ! So the "color" for splitting is process_col and the "key" is my row coordinate.
282 : ! Analogous for mpi_comm_cols
283 :
284 8320 : call mpi_comm_split(mpi_comm_parent,process_col,process_row,mpi_comm_rows,mpierr)
285 :
286 8320 : if (mpierr .ne. MPI_SUCCESS) then
287 0 : call MPI_ERROR_STRING(mpierr,mpierr_string, mpi_string_length, mpierr2)
288 0 : write(error_unit,*) "MPI ERROR occured during mpi_comm_split for row communicator: ", trim(mpierr_string)
289 0 : return
290 : endif
291 :
292 8320 : call mpi_comm_split(mpi_comm_parent,process_row,process_col,mpi_comm_cols, mpierr)
293 8320 : if (mpierr .ne. MPI_SUCCESS) then
294 0 : call MPI_ERROR_STRING(mpierr,mpierr_string, mpi_string_length, mpierr2)
295 0 : write(error_unit,*) "MPI ERROR occured during mpi_comm_split for col communicator: ", trim(mpierr_string)
296 0 : return
297 : endif
298 :
299 8320 : call self%set("mpi_comm_rows", mpi_comm_rows,error)
300 8320 : if (error .ne. ELPA_OK) then
301 0 : print *,"Problem setting option. Aborting..."
302 0 : stop
303 : endif
304 8320 : call self%set("mpi_comm_cols", mpi_comm_cols,error)
305 8320 : if (error .ne. ELPA_OK) then
306 0 : print *,"Problem setting option. Aborting..."
307 0 : stop
308 : endif
309 :
310 : ! remember that we created those communicators and we need to free them later
311 8320 : self%communicators_owned = 1
312 :
313 8320 : error = ELPA_OK
314 8320 : return
315 : endif
316 :
317 : ! Externally supplied communicators
318 6016 : if (self%is_set("mpi_comm_rows") == 1 .and. self%is_set("mpi_comm_cols") == 1) then
319 6016 : self%communicators_owned = 0
320 6016 : error = ELPA_OK
321 6016 : return
322 : endif
323 :
324 : ! Otherwise parameters are missing
325 0 : error = ELPA_ERROR
326 : #endif
327 :
328 6976 : end function
329 :
330 : !c> /*! \brief C interface for the implementation of the elpa_setup method
331 : !c> *
332 : !c> * \param elpa_t handle of the ELPA object which describes the problem to
333 : !c> * be set up
334 : !c> * \result int error code, which can be queried with elpa_strerr
335 : !c> */
336 : !c> int elpa_setup(elpa_t handle);
337 2016 : function elpa_setup_c(handle) result(error) bind(C, name="elpa_setup")
338 : type(c_ptr), intent(in), value :: handle
339 : type(elpa_impl_t), pointer :: self
340 : integer(kind=c_int) :: error
341 :
342 2016 : call c_f_pointer(handle, self)
343 2016 : error = self%setup()
344 2016 : end function
345 :
346 :
347 : !c> /*! \brief C interface for the implementation of the elpa_set_integer method
348 : !c> * This method is available to the user as C generic elpa_set method
349 : !c> *
350 : !c> * \param handle handle of the ELPA object for which a key/value pair should be set
351 : !c> * \param name the name of the key
352 : !c> * \param value the value to be set for the key
353 : !c> * \param error on return the error code, which can be queried with elpa_strerr()
354 : !c> * \result void
355 : !c> */
356 : !c> void elpa_set_integer(elpa_t handle, const char *name, int value, int *error);
357 18144 : subroutine elpa_set_integer_c(handle, name_p, value, error) bind(C, name="elpa_set_integer")
358 : type(c_ptr), intent(in), value :: handle
359 : type(elpa_impl_t), pointer :: self
360 : type(c_ptr), intent(in), value :: name_p
361 18144 : character(len=elpa_strlen_c(name_p)), pointer :: name
362 : integer(kind=c_int), intent(in), value :: value
363 :
364 : #ifdef USE_FORTRAN2008
365 : integer(kind=c_int) , intent(in), optional :: error
366 : #else
367 : integer(kind=c_int) , intent(in) :: error
368 : #endif
369 :
370 18144 : call c_f_pointer(handle, self)
371 18144 : call c_f_pointer(name_p, name)
372 18144 : call elpa_set_integer(self, name, value, error)
373 18144 : end subroutine
374 :
375 :
376 : !c> /*! \brief C interface for the implementation of the elpa_get_integer method
377 : !c> * This method is available to the user as C generic elpa_get method
378 : !c> *
379 : !c> * \param handle handle of the ELPA object for which a key/value pair should be queried
380 : !c> * \param name the name of the key
381 : !c> * \param value the value to be obtain for the key
382 : !c> * \param error on return the error code, which can be queried with elpa_strerr()
383 : !c> * \result void
384 : !c> */
385 : !c> void elpa_get_integer(elpa_t handle, const char *name, int *value, int *error);
386 2016 : subroutine elpa_get_integer_c(handle, name_p, value, error) bind(C, name="elpa_get_integer")
387 : type(c_ptr), intent(in), value :: handle
388 : type(elpa_impl_t), pointer :: self
389 : type(c_ptr), intent(in), value :: name_p
390 2016 : character(len=elpa_strlen_c(name_p)), pointer :: name
391 : integer(kind=c_int) :: value
392 : #ifdef ISE_FORTRAN2008
393 : integer(kind=c_int), intent(inout), optional :: error
394 : #else
395 : integer(kind=c_int), intent(inout) :: error
396 : #endif
397 2016 : call c_f_pointer(handle, self)
398 2016 : call c_f_pointer(name_p, name)
399 2016 : call elpa_get_integer(self, name, value, error)
400 2016 : end subroutine
401 :
402 :
403 : !> \brief function to check whether a key/value pair is set
404 : !> Parameters
405 : !> \param self class(elpa_impl_t) the allocated ELPA object
406 : !> \param name string, the key
407 : !> \result state integer, the state of the key/value pair
408 78992 : function elpa_is_set(self, name) result(state)
409 : class(elpa_impl_t) :: self
410 : character(*), intent(in) :: name
411 : integer :: state
412 :
413 78992 : state = elpa_index_value_is_set_c(self%index, name // c_null_char)
414 157984 : end function
415 :
416 : !> \brief function to check whether a key/value pair can be set
417 : !> Parameters
418 : !> \param self class(elpa_impl_t) the allocated ELPA object
419 : !> \param name string, the key
420 : !> \param value integer, value
421 : !> \result error integer, error code
422 0 : function elpa_can_set(self, name, value) result(error)
423 : class(elpa_impl_t) :: self
424 : character(*), intent(in) :: name
425 : integer(kind=c_int), intent(in) :: value
426 : integer :: error
427 :
428 0 : error = elpa_index_int_is_valid_c(self%index, name // c_null_char, value)
429 0 : end function
430 :
431 :
432 : !> \brief function to convert a value to an human readable string
433 : !> Parameters
434 : !> \param self class(elpa_impl_t) the allocated ELPA object
435 : !> \param option_name string: the name of the options, whose value should be converted
436 : !> \param error integer: errpr code
437 : !> \result string string: the humanreadable string
438 0 : function elpa_value_to_string(self, option_name, error) result(string)
439 : class(elpa_impl_t), intent(in) :: self
440 : character(kind=c_char, len=*), intent(in) :: option_name
441 : type(c_ptr) :: ptr
442 : integer, intent(out), optional :: error
443 : integer :: val, actual_error
444 : character(kind=c_char, len=elpa_index_int_value_to_strlen_c(self%index, option_name // C_NULL_CHAR)), pointer :: string
445 :
446 0 : nullify(string)
447 :
448 0 : call self%get(option_name, val, actual_error)
449 0 : if (actual_error /= ELPA_OK) then
450 0 : if (present(error)) then
451 0 : error = actual_error
452 : endif
453 0 : return
454 : endif
455 :
456 0 : actual_error = elpa_int_value_to_string_c(option_name // C_NULL_CHAR, val, ptr)
457 0 : if (c_associated(ptr)) then
458 0 : call c_f_pointer(ptr, string)
459 : endif
460 :
461 0 : if (present(error)) then
462 0 : error = actual_error
463 : endif
464 0 : end function
465 :
466 :
467 : !c> /*! \brief C interface for the implementation of the elpa_set_double method
468 : !c> * This method is available to the user as C generic elpa_set method
469 : !c> *
470 : !c> * \param handle handle of the ELPA object for which a key/value pair should be set
471 : !c> * \param name the name of the key
472 : !c> * \param value the value to be set for the key
473 : !c> * \param error on return the error code, which can be queried with elpa_strerr()
474 : !c> * \result void
475 : !c> */
476 : !c> void elpa_set_double(elpa_t handle, const char *name, double value, int *error);
477 0 : subroutine elpa_set_double_c(handle, name_p, value, error) bind(C, name="elpa_set_double")
478 : type(c_ptr), intent(in), value :: handle
479 : type(elpa_impl_t), pointer :: self
480 : type(c_ptr), intent(in), value :: name_p
481 0 : character(len=elpa_strlen_c(name_p)), pointer :: name
482 : real(kind=c_double), intent(in), value :: value
483 : #ifdef USE_FORTRAN2008
484 : integer(kind=c_int), intent(in), optional :: error
485 : #else
486 : integer(kind=c_int), intent(in) :: error
487 : #endif
488 0 : call c_f_pointer(handle, self)
489 0 : call c_f_pointer(name_p, name)
490 0 : call elpa_set_double(self, name, value, error)
491 0 : end subroutine
492 :
493 :
494 : !c> /*! \brief C interface for the implementation of the elpa_get_double method
495 : !c> * This method is available to the user as C generic elpa_get method
496 : !c> *
497 : !c> * \param handle handle of the ELPA object for which a key/value pair should be queried
498 : !c> * \param name the name of the key
499 : !c> * \param value the value to be obtain for the key
500 : !c> * \param error on return the error code, which can be queried with elpa_strerr()
501 : !c> * \result void
502 : !c> */
503 : !c> void elpa_get_double(elpa_t handle, const char *name, double *value, int *error);
504 0 : subroutine elpa_get_double_c(handle, name_p, value, error) bind(C, name="elpa_get_double")
505 : type(c_ptr), intent(in), value :: handle
506 : type(elpa_impl_t), pointer :: self
507 : type(c_ptr), intent(in), value :: name_p
508 0 : character(len=elpa_strlen_c(name_p)), pointer :: name
509 : real(kind=c_double) :: value
510 : #ifdef USE_FORTRAN2008
511 : integer(kind=c_int), intent(inout), optional :: error
512 : #else
513 : integer(kind=c_int), intent(inout) :: error
514 : #endif
515 0 : call c_f_pointer(handle, self)
516 0 : call c_f_pointer(name_p, name)
517 0 : call elpa_get_double(self, name, value, error)
518 0 : end subroutine
519 :
520 :
521 : !> \brief function to associate a pointer with an integer value
522 : !> Parameters
523 : !> \param self class(elpa_impl_t) the allocated ELPA object
524 : !> \param name string: the name of the entry
525 : !> \result value integer, pointer: the value for the entry
526 106560 : function elpa_associate_int(self, name) result(value)
527 : class(elpa_impl_t) :: self
528 : character(*), intent(in) :: name
529 : integer(kind=c_int), pointer :: value
530 :
531 : type(c_ptr) :: value_p
532 :
533 106560 : value_p = elpa_index_get_int_loc_c(self%index, name // c_null_char)
534 106560 : if (.not. c_associated(value_p)) then
535 0 : write(error_unit, '(a,a,a)') "ELPA: Warning, received NULL pointer for entry '", name, "'"
536 : endif
537 106560 : call c_f_pointer(value_p, value)
538 213120 : end function
539 :
540 :
541 : !> \brief function to querry the timing information at a certain level
542 : !> Parameters
543 : !> \param self class(elpa_impl_t) the allocated ELPA object
544 : !> \param name1 .. name6 string: the string identifier for the timer region.
545 : !> at the moment 6 nested levels can be queried
546 : !> \result s double: the timer metric for the region. Might be seconds,
547 : !> or any other supported metric
548 32256 : function elpa_get_time(self, name1, name2, name3, name4, name5, name6) result(s)
549 : class(elpa_impl_t), intent(in) :: self
550 : ! this is clunky, but what can you do..
551 : character(len=*), intent(in), optional :: name1, name2, name3, name4, name5, name6
552 : real(kind=c_double) :: s
553 :
554 : #ifdef HAVE_DETAILED_TIMINGS
555 32256 : s = self%timer%get(name1, name2, name3, name4, name5, name6)
556 : #else
557 : s = -1.0
558 : #endif
559 64512 : end function
560 :
561 :
562 : !> \brief function to print the timing tree below at a certain level
563 : !> Parameters
564 : !> \param self class(elpa_impl_t) the allocated ELPA object
565 : !> \param name1 .. name6 string: the string identifier for the timer region.
566 : !> at the moment 4 nested levels can be specified
567 20128 : subroutine elpa_print_times(self, name1, name2, name3, name4)
568 : class(elpa_impl_t), intent(in) :: self
569 : character(len=*), intent(in), optional :: name1, name2, name3, name4
570 : #ifdef HAVE_DETAILED_TIMINGS
571 20128 : call self%timer%print(name1, name2, name3, name4)
572 : #endif
573 40256 : end subroutine
574 :
575 :
576 : !> \brief function to start the timing of a code region
577 : !> Parameters
578 : !> \param self class(elpa_impl_t) the allocated ELPA object
579 : !> \param name string: a chosen identifier name for the code region
580 53184 : subroutine elpa_timer_start(self, name)
581 : class(elpa_impl_t), intent(inout) :: self
582 : character(len=*), intent(in) :: name
583 : #ifdef HAVE_DETAILED_TIMINGS
584 53184 : call self%timer%start(name)
585 : #endif
586 106368 : end subroutine
587 :
588 :
589 : !> \brief function to stop the timing of a code region
590 : !> Parameters
591 : !> \param self class(elpa_impl_t) the allocated ELPA object
592 : !> \param name string: identifier name for the code region to stop
593 53184 : subroutine elpa_timer_stop(self, name)
594 : class(elpa_impl_t), intent(inout) :: self
595 : character(len=*), intent(in) :: name
596 : #ifdef HAVE_DETAILED_TIMINGS
597 53184 : call self%timer%stop(name)
598 : #endif
599 106368 : end subroutine
600 :
601 :
602 : !> \brief elpa_eigenvectors_d: class method to solve the eigenvalue problem for double real matrices
603 : !>
604 : !> The dimensions of the matrix a (locally ditributed and global), the block-cyclic distribution
605 : !> blocksize, the number of eigenvectors
606 : !> to be computed and the MPI communicators are already known to the object and MUST be set BEFORE
607 : !> with the class method "setup"
608 : !>
609 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
610 : !> class method "set"
611 : !>
612 : !> Parameters
613 : !>
614 : !> \param a Distributed matrix for which eigenvalues are to be computed.
615 : !> Distribution is like in Scalapack.
616 : !> The full matrix must be set (not only one half like in scalapack).
617 : !> Destroyed on exit (upper and lower half).
618 : !>
619 : !> \param ev On output: eigenvalues of a, every processor gets the complete set
620 : !>
621 : !> \param q On output: Eigenvectors of a
622 : !> Distribution is like in Scalapack.
623 : !> Must be always dimensioned to the full size (corresponding to (na,na))
624 : !> even if only a part of the eigenvalues is needed.
625 : !>
626 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
627 15360 : subroutine elpa_eigenvectors_d(self, a, ev, q, error)
628 : class(elpa_impl_t) :: self
629 :
630 : #ifdef USE_ASSUMED_SIZE
631 : real(kind=c_double) :: a(self%local_nrows, *), q(self%local_nrows, *)
632 : #else
633 : real(kind=c_double) :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols)
634 : #endif
635 : real(kind=c_double) :: ev(self%na)
636 :
637 : #ifdef USE_FORTRAN2008
638 : integer, optional :: error
639 : #else
640 : integer :: error
641 : #endif
642 : integer :: error2
643 : integer(kind=c_int) :: solver
644 : logical :: success_l
645 :
646 :
647 15360 : call self%get("solver", solver,error2)
648 15360 : if (error2 .ne. ELPA_OK) then
649 0 : print *,"Problem setting option. Aborting..."
650 0 : stop
651 : endif
652 : #ifdef USE_FORTRAN2008
653 15360 : if (present(error)) then
654 15360 : error = error2
655 : endif
656 : #else
657 : error = error2
658 : #endif
659 15360 : if (solver .eq. ELPA_SOLVER_1STAGE) then
660 1920 : call self%autotune_timer%start("accumulator")
661 1920 : success_l = elpa_solve_evp_real_1stage_double_impl(self, a, ev, q)
662 1920 : call self%autotune_timer%stop("accumulator")
663 :
664 13440 : else if (solver .eq. ELPA_SOLVER_2STAGE) then
665 13440 : call self%autotune_timer%start("accumulator")
666 13440 : success_l = elpa_solve_evp_real_2stage_double_impl(self, a, ev, q)
667 13440 : call self%autotune_timer%stop("accumulator")
668 :
669 : else
670 0 : print *,"unknown solver"
671 0 : stop
672 : endif
673 :
674 : #ifdef USE_FORTRAN2008
675 15360 : if (present(error)) then
676 15360 : if (success_l) then
677 15360 : error = ELPA_OK
678 : else
679 0 : error = ELPA_ERROR
680 : endif
681 0 : else if (.not. success_l) then
682 0 : write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
683 : endif
684 : #else
685 : if (success_l) then
686 : error = ELPA_OK
687 : else
688 : error = ELPA_ERROR
689 : endif
690 : #endif
691 15360 : end subroutine
692 :
693 : !c> void elpa_eigenvectors_d(elpa_t handle, double *a, double *ev, double *q, int *error);
694 768 : subroutine elpa_eigenvectors_d_c(handle, a_p, ev_p, q_p, error) bind(C, name="elpa_eigenvectors_d")
695 : type(c_ptr), intent(in), value :: handle, a_p, ev_p, q_p
696 : #ifdef USE_FORTRAN2008
697 : integer(kind=c_int), optional, intent(in) :: error
698 : #else
699 : integer(kind=c_int), intent(in) :: error
700 : #endif
701 :
702 : real(kind=c_double), pointer :: a(:, :), q(:, :), ev(:)
703 : type(elpa_impl_t), pointer :: self
704 :
705 768 : call c_f_pointer(handle, self)
706 768 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
707 768 : call c_f_pointer(ev_p, ev, [self%na])
708 768 : call c_f_pointer(q_p, q, [self%local_nrows, self%local_ncols])
709 :
710 768 : call elpa_eigenvectors_d(self, a, ev, q, error)
711 768 : end subroutine
712 :
713 :
714 : !> \brief elpa_eigenvectors_f: class method to solve the eigenvalue problem for float real matrices
715 : !>
716 : !> The dimensions of the matrix a (locally ditributed and global), the block-cyclic distribution
717 : !> blocksize, the number of eigenvectors
718 : !> to be computed and the MPI communicators are already known to the object and MUST be set BEFORE
719 : !> with the class method "setup"
720 : !>
721 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
722 : !> class method "set"
723 : !>
724 : !> Parameters
725 : !>
726 : !> \param a Distributed matrix for which eigenvalues are to be computed.
727 : !> Distribution is like in Scalapack.
728 : !> The full matrix must be set (not only one half like in scalapack).
729 : !> Destroyed on exit (upper and lower half).
730 : !>
731 : !> \param ev On output: eigenvalues of a, every processor gets the complete set
732 : !>
733 : !> \param q On output: Eigenvectors of a
734 : !> Distribution is like in Scalapack.
735 : !> Must be always dimensioned to the full size (corresponding to (na,na))
736 : !> even if only a part of the eigenvalues is needed.
737 : !>
738 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
739 5904 : subroutine elpa_eigenvectors_f(self, a, ev, q, error)
740 : class(elpa_impl_t) :: self
741 : #ifdef USE_ASSUMED_SIZE
742 : real(kind=c_float) :: a(self%local_nrows, *), q(self%local_nrows, *)
743 : #else
744 : real(kind=c_float) :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols)
745 : #endif
746 : real(kind=c_float) :: ev(self%na)
747 :
748 : #ifdef USE_FORTRAN2008
749 : integer, optional :: error
750 : #else
751 : integer :: error
752 : #endif
753 : integer :: error2
754 : integer(kind=c_int) :: solver
755 : #ifdef WANT_SINGLE_PRECISION_REAL
756 : logical :: success_l
757 :
758 5904 : call self%get("solver",solver, error2)
759 5904 : if (error2 .ne. ELPA_OK) then
760 0 : print *,"Problem getting option. Aborting..."
761 0 : stop
762 : endif
763 : #if USE_FORTRAN2008
764 5904 : if (present(error)) then
765 5904 : error = error2
766 : endif
767 : #else
768 : error = error2
769 : #endif
770 5904 : if (solver .eq. ELPA_SOLVER_1STAGE) then
771 768 : call self%autotune_timer%start("accumulator")
772 768 : success_l = elpa_solve_evp_real_1stage_single_impl(self, a, ev, q)
773 768 : call self%autotune_timer%stop("accumulator")
774 :
775 5136 : else if (solver .eq. ELPA_SOLVER_2STAGE) then
776 5136 : call self%autotune_timer%start("accumulator")
777 5136 : success_l = elpa_solve_evp_real_2stage_single_impl(self, a, ev, q)
778 5136 : call self%autotune_timer%stop("accumulator")
779 :
780 : else
781 0 : print *,"unknown solver"
782 0 : stop
783 : endif
784 :
785 : #ifdef USE_FORTRAN2008
786 5904 : if (present(error)) then
787 5904 : if (success_l) then
788 5904 : error = ELPA_OK
789 : else
790 0 : error = ELPA_ERROR
791 : endif
792 0 : else if (.not. success_l) then
793 0 : write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
794 : endif
795 : #else
796 : if (success_l) then
797 : error = ELPA_OK
798 : else
799 : error = ELPA_ERROR
800 : endif
801 : #endif
802 :
803 : #else
804 0 : print *,"This installation of the ELPA library has not been build with single-precision support"
805 0 : error = ELPA_ERROR
806 : #endif
807 5904 : end subroutine
808 :
809 :
810 : !c> void elpa_eigenvectors_f(elpa_t handle, float *a, float *ev, float *q, int *error);
811 384 : subroutine elpa_eigenvectors_f_c(handle, a_p, ev_p, q_p, error) bind(C, name="elpa_eigenvectors_f")
812 : type(c_ptr), intent(in), value :: handle, a_p, ev_p, q_p
813 : #ifdef USE_FORTRAN2008
814 : integer(kind=c_int), optional, intent(in) :: error
815 : #else
816 : integer(kind=c_int), intent(in) :: error
817 : #endif
818 :
819 : real(kind=c_float), pointer :: a(:, :), q(:, :), ev(:)
820 : type(elpa_impl_t), pointer :: self
821 :
822 384 : call c_f_pointer(handle, self)
823 384 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
824 384 : call c_f_pointer(ev_p, ev, [self%na])
825 384 : call c_f_pointer(q_p, q, [self%local_nrows, self%local_ncols])
826 :
827 384 : call elpa_eigenvectors_f(self, a, ev, q, error)
828 384 : end subroutine
829 :
830 :
831 : !> \brief elpa_eigenvectors_dc: class method to solve the eigenvalue problem for double complex matrices
832 : !>
833 : !> The dimensions of the matrix a (locally ditributed and global), the block-cyclic distribution
834 : !> blocksize, the number of eigenvectors
835 : !> to be computed and the MPI communicators are already known to the object and MUST be set BEFORE
836 : !> with the class method "setup"
837 : !>
838 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
839 : !> class method "set"
840 : !>
841 : !> Parameters
842 : !>
843 : !> \param a Distributed matrix for which eigenvalues are to be computed.
844 : !> Distribution is like in Scalapack.
845 : !> The full matrix must be set (not only one half like in scalapack).
846 : !> Destroyed on exit (upper and lower half).
847 : !>
848 : !> \param ev On output: eigenvalues of a, every processor gets the complete set
849 : !>
850 : !> \param q On output: Eigenvectors of a
851 : !> Distribution is like in Scalapack.
852 : !> Must be always dimensioned to the full size (corresponding to (na,na))
853 : !> even if only a part of the eigenvalues is needed.
854 : !>
855 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
856 10176 : subroutine elpa_eigenvectors_dc(self, a, ev, q, error)
857 : class(elpa_impl_t) :: self
858 :
859 : #ifdef USE_ASSUMED_SIZE
860 : complex(kind=c_double_complex) :: a(self%local_nrows, *), q(self%local_nrows, *)
861 : #else
862 : complex(kind=c_double_complex) :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols)
863 : #endif
864 : real(kind=c_double) :: ev(self%na)
865 : #ifdef USE_FORTRAN2008
866 : integer, optional :: error
867 : #else
868 : integer :: error
869 : #endif
870 : integer :: error2
871 : integer(kind=c_int) :: solver
872 : logical :: success_l
873 :
874 10176 : call self%get("solver", solver,error2)
875 10176 : if (error2 .ne. ELPA_OK) then
876 0 : print *,"Problem getting option. Aborting..."
877 0 : stop
878 : endif
879 : #ifdef USE_FORTRAN2008
880 10176 : if (present(error)) then
881 10176 : error = error2
882 : endif
883 : #else
884 : error = error2
885 : #endif
886 :
887 10176 : if (solver .eq. ELPA_SOLVER_1STAGE) then
888 1920 : call self%autotune_timer%start("accumulator")
889 1920 : success_l = elpa_solve_evp_complex_1stage_double_impl(self, a, ev, q)
890 1920 : call self%autotune_timer%stop("accumulator")
891 :
892 8256 : else if (solver .eq. ELPA_SOLVER_2STAGE) then
893 8256 : call self%autotune_timer%start("accumulator")
894 8256 : success_l = elpa_solve_evp_complex_2stage_double_impl(self, a, ev, q)
895 8256 : call self%autotune_timer%stop("accumulator")
896 :
897 : else
898 0 : print *,"unknown solver"
899 0 : stop
900 : endif
901 :
902 : #ifdef USE_FORTRAN2008
903 10176 : if (present(error)) then
904 10176 : if (success_l) then
905 10176 : error = ELPA_OK
906 : else
907 0 : error = ELPA_ERROR
908 : endif
909 0 : else if (.not. success_l) then
910 0 : write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
911 : endif
912 : #else
913 : if (success_l) then
914 : error = ELPA_OK
915 : else
916 : error = ELPA_ERROR
917 : endif
918 : #endif
919 10176 : end subroutine
920 :
921 :
922 : !c> void elpa_eigenvectors_dc(elpa_t handle, double complex *a, double *ev, double complex *q, int *error);
923 576 : subroutine elpa_eigenvectors_dc_c(handle, a_p, ev_p, q_p, error) bind(C, name="elpa_eigenvectors_dc")
924 : type(c_ptr), intent(in), value :: handle, a_p, ev_p, q_p
925 : #ifdef USE_FORTRAN2008
926 : integer(kind=c_int), optional, intent(in) :: error
927 : #else
928 : integer(kind=c_int), intent(in) :: error
929 : #endif
930 :
931 : complex(kind=c_double_complex), pointer :: a(:, :), q(:, :)
932 : real(kind=c_double), pointer :: ev(:)
933 : type(elpa_impl_t), pointer :: self
934 :
935 576 : call c_f_pointer(handle, self)
936 576 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
937 576 : call c_f_pointer(ev_p, ev, [self%na])
938 576 : call c_f_pointer(q_p, q, [self%local_nrows, self%local_ncols])
939 :
940 576 : call elpa_eigenvectors_dc(self, a, ev, q, error)
941 576 : end subroutine
942 :
943 :
944 : !> \brief elpa_eigenvectors_fc: class method to solve the eigenvalue problem for float complex matrices
945 : !>
946 : !> The dimensions of the matrix a (locally ditributed and global), the block-cyclic distribution
947 : !> blocksize, the number of eigenvectors
948 : !> to be computed and the MPI communicators are already known to the object and MUST be set BEFORE
949 : !> with the class method "setup"
950 : !>
951 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
952 : !> class method "set"
953 : !>
954 : !> Parameters
955 : !>
956 : !> \param a Distributed matrix for which eigenvalues are to be computed.
957 : !> Distribution is like in Scalapack.
958 : !> The full matrix must be set (not only one half like in scalapack).
959 : !> Destroyed on exit (upper and lower half).
960 : !>
961 : !> \param ev On output: eigenvalues of a, every processor gets the complete set
962 : !>
963 : !> \param q On output: Eigenvectors of a
964 : !> Distribution is like in Scalapack.
965 : !> Must be always dimensioned to the full size (corresponding to (na,na))
966 : !> even if only a part of the eigenvalues is needed.
967 : !>
968 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
969 4800 : subroutine elpa_eigenvectors_fc(self, a, ev, q, error)
970 : class(elpa_impl_t) :: self
971 : #ifdef USE_ASSUMED_SIZE
972 : complex(kind=c_float_complex) :: a(self%local_nrows, *), q(self%local_nrows, *)
973 : #else
974 : complex(kind=c_float_complex) :: a(self%local_nrows, self%local_ncols), q(self%local_nrows, self%local_ncols)
975 : #endif
976 : real(kind=c_float) :: ev(self%na)
977 : #ifdef USE_FORTRAN2008
978 : integer, optional :: error
979 : #else
980 : integer :: error
981 : #endif
982 : integer :: error2
983 : integer(kind=c_int) :: solver
984 : #ifdef WANT_SINGLE_PRECISION_COMPLEX
985 : logical :: success_l
986 :
987 4800 : call self%get("solver", solver,error2)
988 4800 : if (error2 .ne. ELPA_OK) then
989 0 : print *,"Problem getting option. Aborting..."
990 0 : stop
991 : endif
992 : #ifdef USE_FORTRAN2008
993 4800 : if (present(error)) then
994 4800 : error = error2
995 : endif
996 : #else
997 : error = error2
998 : #endif
999 4800 : if (solver .eq. ELPA_SOLVER_1STAGE) then
1000 768 : call self%autotune_timer%start("accumulator")
1001 768 : success_l = elpa_solve_evp_complex_1stage_single_impl(self, a, ev, q)
1002 768 : call self%autotune_timer%stop("accumulator")
1003 :
1004 4032 : else if (solver .eq. ELPA_SOLVER_2STAGE) then
1005 4032 : call self%autotune_timer%start("accumulator")
1006 4032 : success_l = elpa_solve_evp_complex_2stage_single_impl(self, a, ev, q)
1007 4032 : call self%autotune_timer%stop("accumulator")
1008 :
1009 : else
1010 0 : print *,"unknown solver"
1011 0 : stop
1012 : endif
1013 : #ifdef USE_FORTRAN2008
1014 4800 : if (present(error)) then
1015 4800 : if (success_l) then
1016 4800 : error = ELPA_OK
1017 : else
1018 0 : error = ELPA_ERROR
1019 : endif
1020 0 : else if (.not. success_l) then
1021 0 : write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
1022 : endif
1023 : #else
1024 : if (success_l) then
1025 : error = ELPA_OK
1026 : else
1027 : error = ELPA_ERROR
1028 : endif
1029 : #endif
1030 :
1031 : #else
1032 0 : print *,"This installation of the ELPA library has not been build with single-precision support"
1033 0 : error = ELPA_ERROR
1034 : #endif
1035 4800 : end subroutine
1036 :
1037 :
1038 : !c> void elpa_eigenvectors_fc(elpa_t handle, float complex *a, float *ev, float complex *q, int *error);
1039 288 : subroutine elpa_eigenvectors_fc_c(handle, a_p, ev_p, q_p, error) bind(C, name="elpa_eigenvectors_fc")
1040 : type(c_ptr), intent(in), value :: handle, a_p, ev_p, q_p
1041 : #ifdef USE_FORTRAN2008
1042 : integer(kind=c_int), optional, intent(in) :: error
1043 : #else
1044 : integer(kind=c_int), intent(in) :: error
1045 : #endif
1046 : complex(kind=c_float_complex), pointer :: a(:, :), q(:, :)
1047 : real(kind=c_float), pointer :: ev(:)
1048 : type(elpa_impl_t), pointer :: self
1049 :
1050 288 : call c_f_pointer(handle, self)
1051 288 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
1052 288 : call c_f_pointer(ev_p, ev, [self%na])
1053 288 : call c_f_pointer(q_p, q, [self%local_nrows, self%local_ncols])
1054 :
1055 288 : call elpa_eigenvectors_fc(self, a, ev, q, error)
1056 288 : end subroutine
1057 :
1058 :
1059 :
1060 : !> \brief elpa_eigenvalues_d: class method to solve the eigenvalue problem for double real matrices
1061 : !>
1062 : !> The dimensions of the matrix a (locally ditributed and global), the block-cyclic distribution
1063 : !> blocksize, the number of eigenvectors
1064 : !> to be computed and the MPI communicators are already known to the object and MUST be set BEFORE
1065 : !> with the class method "setup"
1066 : !>
1067 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
1068 : !> class method "set"
1069 : !>
1070 : !> Parameters
1071 : !>
1072 : !> \param a Distributed matrix for which eigenvalues are to be computed.
1073 : !> Distribution is like in Scalapack.
1074 : !> The full matrix must be set (not only one half like in scalapack).
1075 : !> Destroyed on exit (upper and lower half).
1076 : !>
1077 : !> \param ev On output: eigenvalues of a, every processor gets the complete set
1078 : !>
1079 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
1080 768 : subroutine elpa_eigenvalues_d(self, a, ev, error)
1081 : class(elpa_impl_t) :: self
1082 : #ifdef USE_ASSUMED_SIZE
1083 : real(kind=c_double) :: a(self%local_nrows, *)
1084 : #else
1085 : real(kind=c_double) :: a(self%local_nrows, self%local_ncols)
1086 : #endif
1087 : real(kind=c_double) :: ev(self%na)
1088 : #ifdef USE_FORTRAN2008
1089 : integer, optional :: error
1090 : #else
1091 : integer :: error
1092 : #endif
1093 : integer :: error2
1094 : integer(kind=c_int) :: solver
1095 : logical :: success_l
1096 :
1097 :
1098 768 : call self%get("solver", solver,error2)
1099 768 : if (error2 .ne. ELPA_OK) then
1100 0 : print *,"Problem getting option. Aborting..."
1101 0 : stop
1102 : endif
1103 : #ifdef USE_FORTRAN2008
1104 768 : if (present(error)) then
1105 768 : error = error2
1106 : endif
1107 : #else
1108 : error = error2
1109 : #endif
1110 768 : if (solver .eq. ELPA_SOLVER_1STAGE) then
1111 384 : call self%autotune_timer%start("accumulator")
1112 384 : success_l = elpa_solve_evp_real_1stage_double_impl(self, a, ev)
1113 384 : call self%autotune_timer%stop("accumulator")
1114 :
1115 384 : else if (solver .eq. ELPA_SOLVER_2STAGE) then
1116 384 : call self%autotune_timer%start("accumulator")
1117 384 : success_l = elpa_solve_evp_real_2stage_double_impl(self, a, ev)
1118 384 : call self%autotune_timer%stop("accumulator")
1119 :
1120 : else
1121 0 : print *,"unknown solver"
1122 0 : stop
1123 : endif
1124 :
1125 : #ifdef USE_FORTRAN2008
1126 768 : if (present(error)) then
1127 768 : if (success_l) then
1128 768 : error = ELPA_OK
1129 : else
1130 0 : error = ELPA_ERROR
1131 : endif
1132 0 : else if (.not. success_l) then
1133 0 : write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
1134 : endif
1135 : #else
1136 : if (success_l) then
1137 : error = ELPA_OK
1138 : else
1139 : error = ELPA_ERROR
1140 : endif
1141 : #endif
1142 768 : end subroutine
1143 :
1144 : !c> void elpa_eigenvalues_d(elpa_t handle, double *a, double *ev, int *error);
1145 0 : subroutine elpa_eigenvalues_d_c(handle, a_p, ev_p, error) bind(C, name="elpa_eigenvalues_d")
1146 : type(c_ptr), intent(in), value :: handle, a_p, ev_p
1147 : integer(kind=c_int), intent(in) :: error
1148 :
1149 : real(kind=c_double), pointer :: a(:, :), ev(:)
1150 : type(elpa_impl_t), pointer :: self
1151 :
1152 0 : call c_f_pointer(handle, self)
1153 0 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
1154 0 : call c_f_pointer(ev_p, ev, [self%na])
1155 :
1156 0 : call elpa_eigenvalues_d(self, a, ev, error)
1157 0 : end subroutine
1158 :
1159 :
1160 : !> \brief elpa_eigenvectors_f: class method to solve the eigenvalue problem for float real matrices
1161 : !>
1162 : !> The dimensions of the matrix a (locally ditributed and global), the block-cyclic distribution
1163 : !> blocksize, the number of eigenvectors
1164 : !> to be computed and the MPI communicators are already known to the object and MUST be set BEFORE
1165 : !> with the class method "setup"
1166 : !>
1167 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
1168 : !> class method "set"
1169 : !>
1170 : !> Parameters
1171 : !>
1172 : !> \param a Distributed matrix for which eigenvalues are to be computed.
1173 : !> Distribution is like in Scalapack.
1174 : !> The full matrix must be set (not only one half like in scalapack).
1175 : !> Destroyed on exit (upper and lower half).
1176 : !>
1177 : !> \param ev On output: eigenvalues of a, every processor gets the complete set
1178 : !>
1179 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
1180 192 : subroutine elpa_eigenvalues_f(self, a, ev, error)
1181 : class(elpa_impl_t) :: self
1182 : #ifdef USE_ASSUMED_SIZE
1183 : real(kind=c_float) :: a(self%local_nrows, *)
1184 : #else
1185 : real(kind=c_float) :: a(self%local_nrows, self%local_ncols)
1186 : #endif
1187 : real(kind=c_float) :: ev(self%na)
1188 : #ifdef USE_FORTRAN2008
1189 : integer, optional :: error
1190 : #else
1191 : integer :: error
1192 : #endif
1193 : integer :: error2
1194 : integer(kind=c_int) :: solver
1195 : #ifdef WANT_SINGLE_PRECISION_REAL
1196 : logical :: success_l
1197 :
1198 192 : call self%get("solver",solver,error2)
1199 192 : if (error2 .ne. ELPA_OK) then
1200 0 : print *,"Problem getting option. Aborting..."
1201 0 : stop
1202 : endif
1203 : #ifdef USE_FORTRAN2008
1204 192 : if (present(error)) then
1205 192 : error = error2
1206 : endif
1207 : #else
1208 : error = error2
1209 : #endif
1210 192 : if (solver .eq. ELPA_SOLVER_1STAGE) then
1211 96 : call self%autotune_timer%start("accumulator")
1212 96 : success_l = elpa_solve_evp_real_1stage_single_impl(self, a, ev)
1213 96 : call self%autotune_timer%stop("accumulator")
1214 :
1215 96 : else if (solver .eq. ELPA_SOLVER_2STAGE) then
1216 96 : call self%autotune_timer%start("accumulator")
1217 96 : success_l = elpa_solve_evp_real_2stage_single_impl(self, a, ev)
1218 96 : call self%autotune_timer%stop("accumulator")
1219 :
1220 : else
1221 0 : print *,"unknown solver"
1222 0 : stop
1223 : endif
1224 :
1225 : #ifdef USE_FORTRAN2008
1226 192 : if (present(error)) then
1227 192 : if (success_l) then
1228 192 : error = ELPA_OK
1229 : else
1230 0 : error = ELPA_ERROR
1231 : endif
1232 0 : else if (.not. success_l) then
1233 0 : write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
1234 : endif
1235 : #else
1236 : if (success_l) then
1237 : error = ELPA_OK
1238 : else
1239 : error = ELPA_ERROR
1240 : endif
1241 : #endif
1242 :
1243 : #else
1244 0 : print *,"This installation of the ELPA library has not been build with single-precision support"
1245 0 : error = ELPA_ERROR
1246 : #endif
1247 192 : end subroutine
1248 :
1249 :
1250 : !c> void elpa_eigenvalues_f(elpa_t handle, float *a, float *ev, int *error);
1251 0 : subroutine elpa_eigenvalues_f_c(handle, a_p, ev_p, error) bind(C, name="elpa_eigenvalues_f")
1252 : type(c_ptr), intent(in), value :: handle, a_p, ev_p
1253 : #ifdef USE_FORTRAN2008
1254 : integer(kind=c_int), intent(in), optional :: error
1255 : #else
1256 : integer(kind=c_int), intent(in) :: error
1257 : #endif
1258 : real(kind=c_float), pointer :: a(:, :), ev(:)
1259 : type(elpa_impl_t), pointer :: self
1260 :
1261 0 : call c_f_pointer(handle, self)
1262 0 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
1263 0 : call c_f_pointer(ev_p, ev, [self%na])
1264 :
1265 0 : call elpa_eigenvalues_f(self, a, ev, error)
1266 0 : end subroutine
1267 :
1268 :
1269 : !> \brief elpa_eigenvalues_dc: class method to solve the eigenvalue problem for double complex matrices
1270 : !>
1271 : !> The dimensions of the matrix a (locally ditributed and global), the block-cyclic distribution
1272 : !> blocksize, the number of eigenvectors
1273 : !> to be computed and the MPI communicators are already known to the object and MUST be set BEFORE
1274 : !> with the class method "setup"
1275 : !>
1276 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
1277 : !> class method "set"
1278 : !>
1279 : !> Parameters
1280 : !>
1281 : !> \param a Distributed matrix for which eigenvalues are to be computed.
1282 : !> Distribution is like in Scalapack.
1283 : !> The full matrix must be set (not only one half like in scalapack).
1284 : !> Destroyed on exit (upper and lower half).
1285 : !>
1286 : !> \param ev On output: eigenvalues of a, every processor gets the complete set
1287 : !>
1288 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
1289 384 : subroutine elpa_eigenvalues_dc(self, a, ev, error)
1290 : class(elpa_impl_t) :: self
1291 : #ifdef USE_ASSUMED_SIZE
1292 : complex(kind=c_double_complex) :: a(self%local_nrows, *)
1293 : #else
1294 : complex(kind=c_double_complex) :: a(self%local_nrows, self%local_ncols)
1295 : #endif
1296 : real(kind=c_double) :: ev(self%na)
1297 :
1298 : #ifdef USE_FORTRAN2008
1299 : integer, optional :: error
1300 : #else
1301 : integer :: error
1302 : #endif
1303 : integer :: error2
1304 : integer(kind=c_int) :: solver
1305 : logical :: success_l
1306 :
1307 384 : call self%get("solver", solver,error2)
1308 384 : if (error2 .ne. ELPA_OK) then
1309 0 : print *,"Problem getting option. Aborting..."
1310 0 : stop
1311 : endif
1312 : #ifdef USE_FORTRAN2008
1313 384 : if (present(error)) then
1314 384 : error = error2
1315 : endif
1316 : #else
1317 : error = error2
1318 : #endif
1319 384 : if (solver .eq. ELPA_SOLVER_1STAGE) then
1320 192 : call self%autotune_timer%start("accumulator")
1321 192 : success_l = elpa_solve_evp_complex_1stage_double_impl(self, a, ev)
1322 192 : call self%autotune_timer%stop("accumulator")
1323 :
1324 192 : else if (solver .eq. ELPA_SOLVER_2STAGE) then
1325 192 : call self%autotune_timer%start("accumulator")
1326 192 : success_l = elpa_solve_evp_complex_2stage_double_impl(self, a, ev)
1327 192 : call self%autotune_timer%stop("accumulator")
1328 :
1329 : else
1330 0 : print *,"unknown solver"
1331 0 : stop
1332 : endif
1333 :
1334 : #ifdef USE_FORTRAN2008
1335 384 : if (present(error)) then
1336 384 : if (success_l) then
1337 384 : error = ELPA_OK
1338 : else
1339 0 : error = ELPA_ERROR
1340 : endif
1341 0 : else if (.not. success_l) then
1342 0 : write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
1343 : endif
1344 : #else
1345 : if (success_l) then
1346 : error = ELPA_OK
1347 : else
1348 : error = ELPA_ERROR
1349 : endif
1350 : #endif
1351 384 : end subroutine
1352 :
1353 :
1354 : !c> void elpa_eigenvalues_dc(elpa_t handle, double complex *a, double *ev, int *error);
1355 0 : subroutine elpa_eigenvalues_dc_c(handle, a_p, ev_p, error) bind(C, name="elpa_eigenvalues_dc")
1356 : type(c_ptr), intent(in), value :: handle, a_p, ev_p
1357 : #ifdef USE_FORTRAN2008
1358 : integer(kind=c_int), optional, intent(in) :: error
1359 : #else
1360 : integer(kind=c_int), intent(in) :: error
1361 : #endif
1362 : complex(kind=c_double_complex), pointer :: a(:, :)
1363 : real(kind=c_double), pointer :: ev(:)
1364 : type(elpa_impl_t), pointer :: self
1365 :
1366 0 : call c_f_pointer(handle, self)
1367 0 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
1368 0 : call c_f_pointer(ev_p, ev, [self%na])
1369 :
1370 0 : call elpa_eigenvalues_dc(self, a, ev, error)
1371 0 : end subroutine
1372 :
1373 :
1374 : !> \brief elpa_eigenvalues_fc: class method to solve the eigenvalue problem for float complex matrices
1375 : !>
1376 : !> The dimensions of the matrix a (locally ditributed and global), the block-cyclic distribution
1377 : !> blocksize, the number of eigenvectors
1378 : !> to be computed and the MPI communicators are already known to the object and MUST be set BEFORE
1379 : !> with the class method "setup"
1380 : !>
1381 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
1382 : !> class method "set"
1383 : !>
1384 : !> Parameters
1385 : !>
1386 : !> \param a Distributed matrix for which eigenvalues are to be computed.
1387 : !> Distribution is like in Scalapack.
1388 : !> The full matrix must be set (not only one half like in scalapack).
1389 : !> Destroyed on exit (upper and lower half).
1390 : !>
1391 : !> \param ev On output: eigenvalues of a, every processor gets the complete set
1392 : !>
1393 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
1394 192 : subroutine elpa_eigenvalues_fc(self, a, ev, error)
1395 : class(elpa_impl_t) :: self
1396 : #ifdef USE_ASSUMED_SIZE
1397 : complex(kind=c_float_complex) :: a(self%local_nrows, *)
1398 : #else
1399 : complex(kind=c_float_complex) :: a(self%local_nrows, self%local_ncols)
1400 : #endif
1401 : real(kind=c_float) :: ev(self%na)
1402 :
1403 : #ifdef USE_FORTRAN2008
1404 : integer, optional :: error
1405 : #else
1406 : integer :: error
1407 : #endif
1408 : integer :: error2
1409 : integer(kind=c_int) :: solver
1410 : #ifdef WANT_SINGLE_PRECISION_COMPLEX
1411 : logical :: success_l
1412 :
1413 192 : call self%get("solver", solver,error2)
1414 192 : if (error2 .ne. ELPA_OK) then
1415 0 : print *,"Problem getting option. Aborting..."
1416 0 : stop
1417 : endif
1418 : #ifdef USE_FORTRAN2008
1419 192 : if (present(error)) then
1420 192 : error = error2
1421 : endif
1422 : #else
1423 : error = error2
1424 : #endif
1425 192 : if (solver .eq. ELPA_SOLVER_1STAGE) then
1426 96 : call self%autotune_timer%start("accumulator")
1427 96 : success_l = elpa_solve_evp_complex_1stage_single_impl(self, a, ev)
1428 96 : call self%autotune_timer%stop("accumulator")
1429 :
1430 96 : else if (solver .eq. ELPA_SOLVER_2STAGE) then
1431 96 : call self%autotune_timer%start("accumulator")
1432 96 : success_l = elpa_solve_evp_complex_2stage_single_impl(self, a, ev)
1433 96 : call self%autotune_timer%stop("accumulator")
1434 :
1435 : else
1436 0 : print *,"unknown solver"
1437 0 : stop
1438 : endif
1439 :
1440 : #ifdef USE_FORTRAN2008
1441 192 : if (present(error)) then
1442 192 : if (success_l) then
1443 192 : error = ELPA_OK
1444 : else
1445 0 : error = ELPA_ERROR
1446 : endif
1447 0 : else if (.not. success_l) then
1448 0 : write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
1449 : endif
1450 : #else
1451 : if (success_l) then
1452 : error = ELPA_OK
1453 : else
1454 : error = ELPA_ERROR
1455 : endif
1456 : #endif
1457 :
1458 : #else
1459 0 : print *,"This installation of the ELPA library has not been build with single-precision support"
1460 0 : error = ELPA_ERROR
1461 : #endif
1462 192 : end subroutine
1463 :
1464 :
1465 : !c> void elpa_eigenvalues_fc(elpa_t handle, float complex *a, float *ev, int *error);
1466 0 : subroutine elpa_eigenvalues_fc_c(handle, a_p, ev_p, error) bind(C, name="elpa_eigenvalues_fc")
1467 : type(c_ptr), intent(in), value :: handle, a_p, ev_p
1468 : #ifdef USE_FORTRAN2008
1469 : integer(kind=c_int), optional, intent(in) :: error
1470 : #else
1471 : integer(kind=c_int), intent(in) :: error
1472 : #endif
1473 : complex(kind=c_float_complex), pointer :: a(:, :)
1474 : real(kind=c_float), pointer :: ev(:)
1475 : type(elpa_impl_t), pointer :: self
1476 :
1477 0 : call c_f_pointer(handle, self)
1478 0 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
1479 0 : call c_f_pointer(ev_p, ev, [self%na])
1480 :
1481 0 : call elpa_eigenvalues_fc(self, a, ev, error)
1482 0 : end subroutine
1483 :
1484 : #if 0
1485 : !********************************************************************************************************
1486 : ! GENERALIZED EIGENVECTOR PROBLEM
1487 : !********************************************************************************************************
1488 :
1489 : !> \brief elpa_generalized_eigenvectors_d: class method to solve the eigenvalue problem for double real matrices
1490 : !>
1491 : !> The dimensions of the matrix a (locally ditributed and global), the block-cyclic distribution
1492 : !> blocksize, the number of eigenvectors
1493 : !> to be computed and the MPI communicators are already known to the object and MUST be set BEFORE
1494 : !> with the class method "setup"
1495 : !>
1496 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
1497 : !> class method "set"
1498 : !>
1499 : !> Parameters
1500 : !>
1501 : !> \param a Distributed matrix for which eigenvalues are to be computed.
1502 : !> Distribution is like in Scalapack.
1503 : !> The full matrix must be set (not only one half like in scalapack).
1504 : !> Destroyed on exit (upper and lower half).
1505 : !>
1506 : !> \param ev On output: eigenvalues of a, every processor gets the complete set
1507 : !>
1508 : !> \param q On output: Eigenvectors of a
1509 : !> Distribution is like in Scalapack.
1510 : !> Must be always dimensioned to the full size (corresponding to (na,na))
1511 : !> even if only a part of the eigenvalues is needed.
1512 : !>
1513 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
1514 : subroutine elpa_generalized_eigenvectors_d(self, a, b, ev, q, sc_desc, error)
1515 : use elpa2_impl
1516 : use elpa1_impl
1517 : use elpa_utilities, only : error_unit
1518 : use iso_c_binding
1519 : class(elpa_impl_t) :: self
1520 :
1521 : #ifdef USE_ASSUMED_SIZE
1522 : real(kind=c_double) :: a(self%local_nrows, *), b(self%local_nrows, *), q(self%local_nrows, *)
1523 : #else
1524 : real(kind=c_double) :: a(self%local_nrows, self%local_ncols), b(self%local_nrows, self%local_ncols), &
1525 : q(self%local_nrows, self%local_ncols)
1526 : #endif
1527 : real(kind=c_double) :: ev(self%na)
1528 : integer :: sc_desc(SC_DESC_LEN)
1529 :
1530 : integer, optional :: error
1531 : integer :: error_l
1532 : integer(kind=c_int) :: solver
1533 : logical :: success_l
1534 :
1535 : call self%elpa_transform_generalized_d(a, b, sc_desc, error_l)
1536 : if (present(error)) then
1537 : error = error_l
1538 : else if (error_l .ne. ELPA_OK) then
1539 : write(error_unit,'(a)') "ELPA: Error in transform_generalized() and you did not check for errors!"
1540 : endif
1541 :
1542 : call self%get("solver", solver)
1543 : if (solver .eq. ELPA_SOLVER_1STAGE) then
1544 : success_l = elpa_solve_evp_real_1stage_double_impl(self, a, ev, q)
1545 :
1546 : else if (solver .eq. ELPA_SOLVER_2STAGE) then
1547 : success_l = elpa_solve_evp_real_2stage_double_impl(self, a, ev, q)
1548 : else
1549 : print *,"unknown solver"
1550 : stop
1551 : endif
1552 :
1553 : if (present(error)) then
1554 : if (success_l) then
1555 : error = ELPA_OK
1556 : else
1557 : error = ELPA_ERROR
1558 : endif
1559 : else if (.not. success_l) then
1560 : write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
1561 : endif
1562 :
1563 : call self%elpa_transform_back_generalized_d(b, q, sc_desc, error_l)
1564 : if (present(error)) then
1565 : error = error_l
1566 : else if (error_l .ne. ELPA_OK) then
1567 : write(error_unit,'(a)') "ELPA: Error in transform_back_generalized() and you did not check for errors!"
1568 : endif
1569 : end subroutine
1570 :
1571 : !c> void elpa_generalized_eigenvectors_d(elpa_t handle, double *a, double *ev, double *q, int *error);
1572 : subroutine elpa_generalized_eigenvectors_d_c(handle, a_p, b_p, ev_p, q_p, sc_desc_p, error) &
1573 : bind(C, name="elpa_generalized_eigenvectors_d")
1574 : type(c_ptr), intent(in), value :: handle, a_p, b_p, ev_p, q_p, sc_desc_p
1575 : integer(kind=c_int), optional, intent(in) :: error
1576 :
1577 : real(kind=c_double), pointer :: a(:, :), b(:, :), q(:, :), ev(:)
1578 : integer(kind=c_int), pointer :: sc_desc(:)
1579 : type(elpa_impl_t), pointer :: self
1580 :
1581 : call c_f_pointer(handle, self)
1582 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
1583 : call c_f_pointer(b_p, b, [self%local_nrows, self%local_ncols])
1584 : call c_f_pointer(ev_p, ev, [self%na])
1585 : call c_f_pointer(q_p, q, [self%local_nrows, self%local_ncols])
1586 : call c_f_pointer(sc_desc_p, sc_desc, [SC_DESC_LEN])
1587 :
1588 : call elpa_generalized_eigenvectors_d(self, a, b, ev, q, sc_desc, error)
1589 : end subroutine
1590 :
1591 :
1592 : !> \brief elpa_generalized_eigenvectors_f: class method to solve the eigenvalue problem for float real matrices
1593 : !>
1594 : !> The dimensions of the matrix a (locally ditributed and global), the block-cyclic distribution
1595 : !> blocksize, the number of eigenvectors
1596 : !> to be computed and the MPI communicators are already known to the object and MUST be set BEFORE
1597 : !> with the class method "setup"
1598 : !>
1599 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
1600 : !> class method "set"
1601 : !>
1602 : !> Parameters
1603 : !>
1604 : !> \param a Distributed matrix for which eigenvalues are to be computed.
1605 : !> Distribution is like in Scalapack.
1606 : !> The full matrix must be set (not only one half like in scalapack).
1607 : !> Destroyed on exit (upper and lower half).
1608 : !>
1609 : !> \param ev On output: eigenvalues of a, every processor gets the complete set
1610 : !>
1611 : !> \param q On output: Eigenvectors of a
1612 : !> Distribution is like in Scalapack.
1613 : !> Must be always dimensioned to the full size (corresponding to (na,na))
1614 : !> even if only a part of the eigenvalues is needed.
1615 : !>
1616 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
1617 : subroutine elpa_generalized_eigenvectors_f(self, a, b, ev, q, sc_desc, error)
1618 : use elpa2_impl
1619 : use elpa1_impl
1620 : use elpa_utilities, only : error_unit
1621 : use iso_c_binding
1622 : class(elpa_impl_t) :: self
1623 : #ifdef USE_ASSUMED_SIZE
1624 : real(kind=c_float) :: a(self%local_nrows, *), b(self%local_nrows, *), q(self%local_nrows, *)
1625 : #else
1626 : real(kind=c_float) :: a(self%local_nrows, self%local_ncols), b(self%local_nrows, self%local_ncols), &
1627 : q(self%local_nrows, self%local_ncols)
1628 : #endif
1629 : real(kind=c_float) :: ev(self%na)
1630 : integer :: sc_desc(SC_DESC_LEN)
1631 :
1632 : integer, optional :: error
1633 : integer :: error_l
1634 : integer(kind=c_int) :: solver
1635 : #ifdef WANT_SINGLE_PRECISION_REAL
1636 : logical :: success_l
1637 :
1638 : call self%elpa_transform_generalized_f(a, b, sc_desc, error_l)
1639 : if (present(error)) then
1640 : error = error_l
1641 : else if (error_l .ne. ELPA_OK) then
1642 : write(error_unit,'(a)') "ELPA: Error in transform_generalized() and you did not check for errors!"
1643 : endif
1644 :
1645 : call self%get("solver",solver)
1646 : if (solver .eq. ELPA_SOLVER_1STAGE) then
1647 : success_l = elpa_solve_evp_real_1stage_single_impl(self, a, ev, q)
1648 :
1649 : else if (solver .eq. ELPA_SOLVER_2STAGE) then
1650 : success_l = elpa_solve_evp_real_2stage_single_impl(self, a, ev, q)
1651 : else
1652 : print *,"unknown solver"
1653 : stop
1654 : endif
1655 :
1656 : if (present(error)) then
1657 : if (success_l) then
1658 : error = ELPA_OK
1659 : else
1660 : error = ELPA_ERROR
1661 : endif
1662 : else if (.not. success_l) then
1663 : write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
1664 : endif
1665 :
1666 : call self%elpa_transform_back_generalized_f(b, q, sc_desc, error_l)
1667 : if (present(error)) then
1668 : error = error_l
1669 : else if (error_l .ne. ELPA_OK) then
1670 : write(error_unit,'(a)') "ELPA: Error in transform_back_generalized() and you did not check for errors!"
1671 : endif
1672 : #else
1673 : print *,"This installation of the ELPA library has not been build with single-precision support"
1674 : error = ELPA_ERROR
1675 : #endif
1676 : end subroutine
1677 :
1678 :
1679 : !c> void elpa_generalized_eigenvectors_f(elpa_t handle, float *a, float *ev, float *q, int *error);
1680 : subroutine elpa_generalized_eigenvectors_f_c(handle, a_p, b_p, ev_p, q_p, sc_desc_p, error) &
1681 : bind(C, name="elpa_generalized_eigenvectors_f")
1682 : type(c_ptr), intent(in), value :: handle, a_p, b_p, ev_p, q_p, sc_desc_p
1683 : integer(kind=c_int), optional, intent(in) :: error
1684 :
1685 : real(kind=c_float), pointer :: a(:, :), b(:, :), q(:, :), ev(:)
1686 : integer(kind=c_int), pointer :: sc_desc(:)
1687 : type(elpa_impl_t), pointer :: self
1688 :
1689 : call c_f_pointer(handle, self)
1690 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
1691 : call c_f_pointer(b_p, b, [self%local_nrows, self%local_ncols])
1692 : call c_f_pointer(ev_p, ev, [self%na])
1693 : call c_f_pointer(q_p, q, [self%local_nrows, self%local_ncols])
1694 : call c_f_pointer(sc_desc_p, sc_desc, [SC_DESC_LEN])
1695 :
1696 : call elpa_generalized_eigenvectors_f(self, a, b, ev, q, sc_desc, error)
1697 : end subroutine
1698 :
1699 :
1700 : !> \brief elpa_generalized_eigenvectors_dc: class method to solve the eigenvalue problem for double complex matrices
1701 : !>
1702 : !> The dimensions of the matrix a (locally ditributed and global), the block-cyclic distribution
1703 : !> blocksize, the number of eigenvectors
1704 : !> to be computed and the MPI communicators are already known to the object and MUST be set BEFORE
1705 : !> with the class method "setup"
1706 : !>
1707 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
1708 : !> class method "set"
1709 : !>
1710 : !> Parameters
1711 : !>
1712 : !> \param a Distributed matrix for which eigenvalues are to be computed.
1713 : !> Distribution is like in Scalapack.
1714 : !> The full matrix must be set (not only one half like in scalapack).
1715 : !> Destroyed on exit (upper and lower half).
1716 : !>
1717 : !> \param ev On output: eigenvalues of a, every processor gets the complete set
1718 : !>
1719 : !> \param q On output: Eigenvectors of a
1720 : !> Distribution is like in Scalapack.
1721 : !> Must be always dimensioned to the full size (corresponding to (na,na))
1722 : !> even if only a part of the eigenvalues is needed.
1723 : !>
1724 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
1725 : subroutine elpa_generalized_eigenvectors_dc(self, a, b, ev, q, sc_desc, error)
1726 : use elpa2_impl
1727 : use elpa1_impl
1728 : use elpa_utilities, only : error_unit
1729 : use iso_c_binding
1730 : class(elpa_impl_t) :: self
1731 :
1732 : #ifdef USE_ASSUMED_SIZE
1733 : complex(kind=c_double_complex) :: a(self%local_nrows, *), b(self%local_nrows, *), q(self%local_nrows, *)
1734 : #else
1735 : complex(kind=c_double_complex) :: a(self%local_nrows, self%local_ncols), b(self%local_nrows, self%local_ncols), &
1736 : q(self%local_nrows, self%local_ncols)
1737 : #endif
1738 : real(kind=c_double) :: ev(self%na)
1739 : integer :: sc_desc(SC_DESC_LEN)
1740 :
1741 : integer, optional :: error
1742 : integer :: error_l
1743 : integer(kind=c_int) :: solver
1744 : logical :: success_l
1745 :
1746 : call self%elpa_transform_generalized_dc(a, b, sc_desc, error_l)
1747 : if (present(error)) then
1748 : error = error_l
1749 : else if (error_l .ne. ELPA_OK) then
1750 : write(error_unit,'(a)') "ELPA: Error in transform_generalized() and you did not check for errors!"
1751 : endif
1752 :
1753 : call self%get("solver", solver)
1754 : if (solver .eq. ELPA_SOLVER_1STAGE) then
1755 : success_l = elpa_solve_evp_complex_1stage_double_impl(self, a, ev, q)
1756 :
1757 : else if (solver .eq. ELPA_SOLVER_2STAGE) then
1758 : success_l = elpa_solve_evp_complex_2stage_double_impl(self, a, ev, q)
1759 : else
1760 : print *,"unknown solver"
1761 : stop
1762 : endif
1763 :
1764 : if (present(error)) then
1765 : if (success_l) then
1766 : error = ELPA_OK
1767 : else
1768 : error = ELPA_ERROR
1769 : endif
1770 : else if (.not. success_l) then
1771 : write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
1772 : endif
1773 :
1774 : call self%elpa_transform_back_generalized_dc(b, q, sc_desc, error_l)
1775 : if (present(error)) then
1776 : error = error_l
1777 : else if (error_l .ne. ELPA_OK) then
1778 : write(error_unit,'(a)') "ELPA: Error in transform_back_generalized() and you did not check for errors!"
1779 : endif
1780 : end subroutine
1781 :
1782 :
1783 : !c> void elpa_generalized_eigenvectors_dc(elpa_t handle, double complex *a, double *ev, double complex *q, int *error);
1784 : subroutine elpa_generalized_eigenvectors_dc_c(handle, a_p, b_p, ev_p, q_p, sc_desc_p, error) &
1785 : bind(C, name="elpa_generalized_eigenvectors_dc")
1786 : type(c_ptr), intent(in), value :: handle, a_p, b_p, ev_p, q_p, sc_desc_p
1787 : integer(kind=c_int), optional, intent(in) :: error
1788 :
1789 : complex(kind=c_double_complex), pointer :: a(:, :), b(:, :), q(:, :)
1790 : real(kind=c_double), pointer :: ev(:)
1791 : integer(kind=c_int), pointer :: sc_desc(:)
1792 : type(elpa_impl_t), pointer :: self
1793 :
1794 : call c_f_pointer(handle, self)
1795 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
1796 : call c_f_pointer(b_p, b, [self%local_nrows, self%local_ncols])
1797 : call c_f_pointer(ev_p, ev, [self%na])
1798 : call c_f_pointer(q_p, q, [self%local_nrows, self%local_ncols])
1799 : call c_f_pointer(sc_desc_p, sc_desc, [SC_DESC_LEN])
1800 :
1801 : call elpa_generalized_eigenvectors_dc(self, a, b, ev, q, sc_desc, error)
1802 : end subroutine
1803 :
1804 :
1805 : !> \brief elpa_generalized_eigenvectors_fc: class method to solve the eigenvalue problem for float complex matrices
1806 : !>
1807 : !> The dimensions of the matrix a (locally ditributed and global), the block-cyclic distribution
1808 : !> blocksize, the number of eigenvectors
1809 : !> to be computed and the MPI communicators are already known to the object and MUST be set BEFORE
1810 : !> with the class method "setup"
1811 : !>
1812 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
1813 : !> class method "set"
1814 : !>
1815 : !> Parameters
1816 : !>
1817 : !> \param a Distributed matrix for which eigenvalues are to be computed.
1818 : !> Distribution is like in Scalapack.
1819 : !> The full matrix must be set (not only one half like in scalapack).
1820 : !> Destroyed on exit (upper and lower half).
1821 : !>
1822 : !> \param ev On output: eigenvalues of a, every processor gets the complete set
1823 : !>
1824 : !> \param q On output: Eigenvectors of a
1825 : !> Distribution is like in Scalapack.
1826 : !> Must be always dimensioned to the full size (corresponding to (na,na))
1827 : !> even if only a part of the eigenvalues is needed.
1828 : !>
1829 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
1830 : subroutine elpa_generalized_eigenvectors_fc(self, a, b, ev, q, sc_desc, error)
1831 : use elpa2_impl
1832 : use elpa1_impl
1833 : use elpa_utilities, only : error_unit
1834 :
1835 : use iso_c_binding
1836 : class(elpa_impl_t) :: self
1837 : #ifdef USE_ASSUMED_SIZE
1838 : complex(kind=c_float_complex) :: a(self%local_nrows, *), b(self%local_nrows, *), q(self%local_nrows, *)
1839 : #else
1840 : complex(kind=c_float_complex) :: a(self%local_nrows, self%local_ncols), b(self%local_nrows, self%local_ncols), &
1841 : q(self%local_nrows, self%local_ncols)
1842 : #endif
1843 : real(kind=c_float) :: ev(self%na)
1844 : integer :: sc_desc(SC_DESC_LEN)
1845 :
1846 : integer, optional :: error
1847 : integer :: error_l
1848 : integer(kind=c_int) :: solver
1849 : #ifdef WANT_SINGLE_PRECISION_COMPLEX
1850 : logical :: success_l
1851 :
1852 : call self%elpa_transform_generalized_fc(a, b, sc_desc, error_l)
1853 : if (present(error)) then
1854 : error = error_l
1855 : else if (error_l .ne. ELPA_OK) then
1856 : write(error_unit,'(a)') "ELPA: Error in transform_generalized() and you did not check for errors!"
1857 : endif
1858 :
1859 : call self%get("solver", solver)
1860 : if (solver .eq. ELPA_SOLVER_1STAGE) then
1861 : success_l = elpa_solve_evp_complex_1stage_single_impl(self, a, ev, q)
1862 :
1863 : else if (solver .eq. ELPA_SOLVER_2STAGE) then
1864 : success_l = elpa_solve_evp_complex_2stage_single_impl(self, a, ev, q)
1865 : else
1866 : print *,"unknown solver"
1867 : stop
1868 : endif
1869 :
1870 : if (present(error)) then
1871 : if (success_l) then
1872 : error = ELPA_OK
1873 : else
1874 : error = ELPA_ERROR
1875 : endif
1876 : else if (.not. success_l) then
1877 : write(error_unit,'(a)') "ELPA: Error in solve() and you did not check for errors!"
1878 : endif
1879 :
1880 : call self%elpa_transform_back_generalized_fc(b, q, sc_desc, error_l)
1881 : if (present(error)) then
1882 : error = error_l
1883 : else if (error_l .ne. ELPA_OK) then
1884 : write(error_unit,'(a)') "ELPA: Error in transform_back_generalized() and you did not check for errors!"
1885 : endif
1886 : #else
1887 : print *,"This installation of the ELPA library has not been build with single-precision support"
1888 : error = ELPA_ERROR
1889 : #endif
1890 : end subroutine
1891 :
1892 :
1893 : !c> void elpa_generalized_eigenvectors_fc(elpa_t handle, float complex *a, float *ev, float complex *q, int *error);
1894 : subroutine elpa_generalized_eigenvectors_fc_c(handle, a_p, b_p, ev_p, q_p, sc_desc_p, error) &
1895 : bind(C, name="elpa_generalized_eigenvectors_fc")
1896 : type(c_ptr), intent(in), value :: handle, a_p, b_p, ev_p, q_p, sc_desc_p
1897 : integer(kind=c_int), optional, intent(in) :: error
1898 :
1899 : complex(kind=c_float_complex), pointer :: a(:, :), b(:, :), q(:, :)
1900 : real(kind=c_float), pointer :: ev(:)
1901 : integer(kind=c_int), pointer :: sc_desc(:)
1902 : type(elpa_impl_t), pointer :: self
1903 :
1904 : call c_f_pointer(handle, self)
1905 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
1906 : call c_f_pointer(b_p, b, [self%local_nrows, self%local_ncols])
1907 : call c_f_pointer(ev_p, ev, [self%na])
1908 : call c_f_pointer(q_p, q, [self%local_nrows, self%local_ncols])
1909 : call c_f_pointer(sc_desc_p, sc_desc, [SC_DESC_LEN])
1910 :
1911 : call elpa_generalized_eigenvectors_fc(self, a, b, ev, q, sc_desc, error)
1912 : end subroutine
1913 :
1914 : #endif
1915 :
1916 :
1917 : !********************************************************************************************************
1918 : ! HERMITIAN MULTIPLY
1919 : !********************************************************************************************************
1920 :
1921 : !> \brief elpa_hermitian_multiply_d: class method to perform C : = A**T * B for double real matrices
1922 : !> where A is a square matrix (self%na,self%na) which is optionally upper or lower triangular
1923 : !> B is a (self%na,ncb) matrix
1924 : !> C is a (self%na,ncb) matrix where optionally only the upper or lower
1925 : !> triangle may be computed
1926 : !>
1927 : !> the MPI commicators and the block-cyclic distribution block size are already known to the type.
1928 : !> Thus the class method "setup" must be called BEFORE this method is used
1929 : !>
1930 : !> \details
1931 : !>
1932 : !> \param self class(elpa_t), the ELPA object
1933 : !> \param uplo_a 'U' if A is upper triangular
1934 : !> 'L' if A is lower triangular
1935 : !> anything else if A is a full matrix
1936 : !> Please note: This pertains to the original A (as set in the calling program)
1937 : !> whereas the transpose of A is used for calculations
1938 : !> If uplo_a is 'U' or 'L', the other triangle is not used at all,
1939 : !> i.e. it may contain arbitrary numbers
1940 : !> \param uplo_c 'U' if only the upper diagonal part of C is needed
1941 : !> 'L' if only the upper diagonal part of C is needed
1942 : !> anything else if the full matrix C is needed
1943 : !> Please note: Even when uplo_c is 'U' or 'L', the other triangle may be
1944 : !> written to a certain extent, i.e. one shouldn't rely on the content there!
1945 : !> \param ncb Number of columns of global matrices B and C
1946 : !> \param a matrix a
1947 : !> \param local_nrows number of rows of local (sub) matrix a, set with class method set("local_nrows",value)
1948 : !> \param local_ncols number of columns of local (sub) matrix a, set with class method set("local_ncols",value)
1949 : !> \param b matrix b
1950 : !> \param nrows_b number of rows of local (sub) matrix b
1951 : !> \param ncols_b number of columns of local (sub) matrix b
1952 : !> \param c matrix c
1953 : !> \param nrows_c number of rows of local (sub) matrix c
1954 : !> \param ncols_c number of columns of local (sub) matrix c
1955 : !> \param error optional argument, error code which can be queried with elpa_strerr
1956 576 : subroutine elpa_hermitian_multiply_d (self, uplo_a, uplo_c, ncb, a, b, nrows_b, ncols_b, &
1957 576 : c, nrows_c, ncols_c, error)
1958 : class(elpa_impl_t) :: self
1959 : character*1 :: uplo_a, uplo_c
1960 : integer(kind=c_int), intent(in) :: nrows_b, ncols_b, nrows_c, ncols_c, ncb
1961 : #ifdef USE_ASSUMED_SIZE
1962 : real(kind=c_double) :: a(self%local_nrows,*), b(nrows_b,*), c(nrows_c,*)
1963 : #else
1964 : real(kind=c_double) :: a(self%local_nrows,self%local_ncols), b(nrows_b,ncols_b), c(nrows_c,ncols_c)
1965 : #endif
1966 : #ifdef USE_FORTRAN2008
1967 : integer, optional :: error
1968 : #else
1969 : integer :: error
1970 : #endif
1971 : logical :: success_l
1972 :
1973 : success_l = elpa_mult_at_b_real_double_impl(self, uplo_a, uplo_c, ncb, a, b, nrows_b, ncols_b, &
1974 576 : c, nrows_c, ncols_c)
1975 :
1976 : #ifdef USE_FORTRAN2008
1977 576 : if (present(error)) then
1978 576 : if (success_l) then
1979 576 : error = ELPA_OK
1980 : else
1981 0 : error = ELPA_ERROR
1982 : endif
1983 0 : else if (.not. success_l) then
1984 0 : write(error_unit,'(a)') "ELPA: Error in hermitian_multiply() and you did not check for errors!"
1985 : endif
1986 : #else
1987 : if (success_l) then
1988 : error = ELPA_OK
1989 : else
1990 : error = ELPA_ERROR
1991 : endif
1992 : #endif
1993 576 : end subroutine
1994 :
1995 : !c> void elpa_hermitian_multiply_d(elpa_t handle, char uplo_a, char uplo_c, int ncb, double *a, double *b, int nrows_b, int ncols_b, double *c, int nrows_c, int ncols_c, int *error);
1996 0 : subroutine elpa_hermitian_multiply_d_c(handle, uplo_a, uplo_c, ncb, a_p, b, nrows_b, &
1997 0 : ncols_b, c, nrows_c, ncols_c, error) &
1998 : bind(C, name="elpa_hermitian_multiply_d")
1999 : type(c_ptr), intent(in), value :: handle, a_p
2000 : character(1,C_CHAR), value :: uplo_a, uplo_c
2001 : integer(kind=c_int), value :: ncb, nrows_b, ncols_b, nrows_c, ncols_c
2002 : #ifdef USE_FORTRAN2008
2003 : integer(kind=c_int), optional, intent(in) :: error
2004 : #else
2005 : integer(kind=c_int), intent(in) :: error
2006 : #endif
2007 : real(kind=c_double), pointer :: a(:, :)
2008 : #ifdef USE_ASSUMED_SIZE
2009 : real(kind=c_double) :: b(nrows_b,*), c(nrows_c,*)
2010 : #else
2011 : real(kind=c_double) :: b(nrows_b,ncols_b), c(nrows_c,ncols_c)
2012 : #endif
2013 : type(elpa_impl_t), pointer :: self
2014 :
2015 0 : call c_f_pointer(handle, self)
2016 0 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
2017 :
2018 : call elpa_hermitian_multiply_d(self, uplo_a, uplo_c, ncb, a, b, nrows_b, &
2019 0 : ncols_b, c, nrows_c, ncols_c, error)
2020 0 : end subroutine
2021 :
2022 : !> \brief elpa_hermitian_multiply_f: class method to perform C : = A**T * B for float real matrices
2023 : !> where A is a square matrix (self%na,self%na) which is optionally upper or lower triangular
2024 : !> B is a (self%na,ncb) matrix
2025 : !> C is a (self%na,ncb) matrix where optionally only the upper or lower
2026 : !> triangle may be computed
2027 : !>
2028 : !> the MPI commicators and the block-cyclic distribution block size are already known to the type.
2029 : !> Thus the class method "setup" must be called BEFORE this method is used
2030 : !>
2031 : !> \details
2032 : !>
2033 : !> \param self class(elpa_t), the ELPA object
2034 : !> \param uplo_a 'U' if A is upper triangular
2035 : !> 'L' if A is lower triangular
2036 : !> anything else if A is a full matrix
2037 : !> Please note: This pertains to the original A (as set in the calling program)
2038 : !> whereas the transpose of A is used for calculations
2039 : !> If uplo_a is 'U' or 'L', the other triangle is not used at all,
2040 : !> i.e. it may contain arbitrary numbers
2041 : !> \param uplo_c 'U' if only the upper diagonal part of C is needed
2042 : !> 'L' if only the upper diagonal part of C is needed
2043 : !> anything else if the full matrix C is needed
2044 : !> Please note: Even when uplo_c is 'U' or 'L', the other triangle may be
2045 : !> written to a certain extent, i.e. one shouldn't rely on the content there!
2046 : !> \param ncb Number of columns of global matrices B and C
2047 : !> \param a matrix a
2048 : !> \param self%local_nrows number of rows of local (sub) matrix a, set with class method set("local_nrows",value)
2049 : !> \param self%local_ncols number of columns of local (sub) matrix a, set with class method set("local_ncols",value)
2050 : !> \param b matrix b
2051 : !> \param nrows_b number of rows of local (sub) matrix b
2052 : !> \param ncols_b number of columns of local (sub) matrix b
2053 : !> \param c matrix c
2054 : !> \param nrows_c number of rows of local (sub) matrix c
2055 : !> \param ncols_c number of columns of local (sub) matrix c
2056 : !> \param error optional argument, returns an error code, which can be queried with elpa_strerr
2057 192 : subroutine elpa_hermitian_multiply_f (self,uplo_a, uplo_c, ncb, a, b, nrows_b, ncols_b, &
2058 192 : c, nrows_c, ncols_c, error)
2059 : class(elpa_impl_t) :: self
2060 : character*1 :: uplo_a, uplo_c
2061 : integer(kind=c_int), intent(in) :: nrows_b, ncols_b, nrows_c, ncols_c, ncb
2062 : #ifdef USE_ASSUMED_SIZE
2063 : real(kind=c_float) :: a(self%local_nrows,*), b(self%local_nrows,*), c(nrows_c,*)
2064 : #else
2065 : real(kind=c_float) :: a(self%local_nrows,self%local_ncols), b(nrows_b,ncols_b), c(nrows_c,ncols_c)
2066 : #endif
2067 : #ifdef USE_FORTRAN2008
2068 : integer, optional :: error
2069 : #else
2070 : integer :: error
2071 : #endif
2072 :
2073 : #ifdef WANT_SINGLE_PRECISION_REAL
2074 : logical :: success_l
2075 :
2076 : success_l = elpa_mult_at_b_real_single_impl(self, uplo_a, uplo_c, ncb, a, b, nrows_b, ncols_b, &
2077 192 : c, nrows_c, ncols_c)
2078 : #ifdef USE_FORTRAN2008
2079 192 : if (present(error)) then
2080 192 : if (success_l) then
2081 192 : error = ELPA_OK
2082 : else
2083 0 : error = ELPA_ERROR
2084 : endif
2085 0 : else if (.not. success_l) then
2086 0 : write(error_unit,'(a)') "ELPA: Error in hermitian_multiply() and you did not check for errors!"
2087 : endif
2088 : #else
2089 : if (success_l) then
2090 : error = ELPA_OK
2091 : else
2092 : error = ELPA_ERROR
2093 : endif
2094 : #endif
2095 :
2096 : #else
2097 0 : print *,"This installation of the ELPA library has not been build with single-precision support"
2098 0 : error = ELPA_ERROR
2099 : #endif
2100 192 : end subroutine
2101 :
2102 : !c> void elpa_hermitian_multiply_f(elpa_t handle, char uplo_a, char uplo_c, int ncb, float *a, float *b, int nrows_b, int ncols_b, float *c, int nrows_c, int ncols_c, int *error);
2103 0 : subroutine elpa_hermitian_multiply_f_c(handle, uplo_a, uplo_c, ncb, a_p, b, nrows_b, &
2104 0 : ncols_b, c, nrows_c, ncols_c, error) &
2105 : bind(C, name="elpa_hermitian_multiply_f")
2106 : type(c_ptr), intent(in), value :: handle, a_p
2107 : character(1,C_CHAR), value :: uplo_a, uplo_c
2108 : integer(kind=c_int), value :: ncb, nrows_b, ncols_b, nrows_c, ncols_c
2109 : #ifdef USE_FORTRAN2008
2110 : integer(kind=c_int), optional, intent(in) :: error
2111 : #else
2112 : integer(kind=c_int), intent(in) :: error
2113 : #endif
2114 :
2115 : real(kind=c_float), pointer :: a(:, :)
2116 : #ifdef USE_ASSUMED_SIZE
2117 : real(kind=c_float) :: b(nrows_b,*), c(nrows_c,*)
2118 : #else
2119 : real(kind=c_float) :: b(nrows_b,ncols_b), c(nrows_c,ncols_c)
2120 : #endif
2121 : type(elpa_impl_t), pointer :: self
2122 :
2123 0 : call c_f_pointer(handle, self)
2124 0 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
2125 :
2126 : call elpa_hermitian_multiply_f(self, uplo_a, uplo_c, ncb, a, b, nrows_b, &
2127 0 : ncols_b, c, nrows_c, ncols_c, error)
2128 0 : end subroutine
2129 :
2130 :
2131 : !> \brief elpa_hermitian_multiply_dc: class method to perform C : = A**H * B for double complex matrices
2132 : !> where A is a square matrix (self%na,self%na) which is optionally upper or lower triangular
2133 : !> B is a (self%na,ncb) matrix
2134 : !> C is a (self%na,ncb) matrix where optionally only the upper or lower
2135 : !> triangle may be computed
2136 : !>
2137 : !> the MPI commicators and the block-cyclic distribution block size are already known to the type.
2138 : !> Thus the class method "setup" must be called BEFORE this method is used
2139 : !>
2140 : !> \details
2141 : !>
2142 : !> \param self class(elpa_t), the ELPA object
2143 : !> \param uplo_a 'U' if A is upper triangular
2144 : !> 'L' if A is lower triangular
2145 : !> anything else if A is a full matrix
2146 : !> Please note: This pertains to the original A (as set in the calling program)
2147 : !> whereas the transpose of A is used for calculations
2148 : !> If uplo_a is 'U' or 'L', the other triangle is not used at all,
2149 : !> i.e. it may contain arbitrary numbers
2150 : !> \param uplo_c 'U' if only the upper diagonal part of C is needed
2151 : !> 'L' if only the upper diagonal part of C is needed
2152 : !> anything else if the full matrix C is needed
2153 : !> Please note: Even when uplo_c is 'U' or 'L', the other triangle may be
2154 : !> written to a certain extent, i.e. one shouldn't rely on the content there!
2155 : !> \param ncb Number of columns of global matrices B and C
2156 : !> \param a matrix a
2157 : !> \param self%local_nrows number of rows of local (sub) matrix a, set with class method set("local_nows",value)
2158 : !> \param self%local_ncols number of columns of local (sub) matrix a, set with class method set("local_ncols",value)
2159 : !> \param b matrix b
2160 : !> \param nrows_b number of rows of local (sub) matrix b
2161 : !> \param ncols_b number of columns of local (sub) matrix b
2162 : !> \param c matrix c
2163 : !> \param nrows_c number of rows of local (sub) matrix c
2164 : !> \param ncols_c number of columns of local (sub) matrix c
2165 : !> \param error optional argument, returns an error code, which can be queried with elpa_strerr
2166 384 : subroutine elpa_hermitian_multiply_dc (self,uplo_a, uplo_c, ncb, a, b, nrows_b, ncols_b, &
2167 384 : c, nrows_c, ncols_c, error)
2168 : class(elpa_impl_t) :: self
2169 : character*1 :: uplo_a, uplo_c
2170 : integer(kind=c_int), intent(in) :: nrows_b, ncols_b, nrows_c, ncols_c, ncb
2171 : #ifdef USE_ASSUMED_SIZE
2172 : complex(kind=c_double_complex) :: a(self%local_nrows,*), b(nrows_b,*), c(nrows_c,*)
2173 : #else
2174 : complex(kind=c_double_complex) :: a(self%local_nrows,self%local_ncols), b(nrows_b,ncols_b), c(nrows_c,ncols_c)
2175 : #endif
2176 :
2177 : #ifdef USE_FORTRAN2008
2178 : integer, optional :: error
2179 : #else
2180 : integer :: error
2181 : #endif
2182 : logical :: success_l
2183 :
2184 : success_l = elpa_mult_ah_b_complex_double_impl(self, uplo_a, uplo_c, ncb, a, b, nrows_b, ncols_b, &
2185 384 : c, nrows_c, ncols_c)
2186 : #ifdef USE_FORTRAN2008
2187 384 : if (present(error)) then
2188 384 : if (success_l) then
2189 384 : error = ELPA_OK
2190 : else
2191 0 : error = ELPA_ERROR
2192 : endif
2193 0 : else if (.not. success_l) then
2194 0 : write(error_unit,'(a)') "ELPA: Error in hermitian_multiply() and you did not check for errors!"
2195 : endif
2196 : #else
2197 : if (success_l) then
2198 : error = ELPA_OK
2199 : else
2200 : error = ELPA_ERROR
2201 : endif
2202 : #endif
2203 384 : end subroutine
2204 :
2205 :
2206 : !c> void elpa_hermitian_multiply_dc(elpa_t handle, char uplo_a, char uplo_c, int ncb, double complex *a, double complex *b, int nrows_b, int ncols_b, double complex *c, int nrows_c, int ncols_c, int *error);
2207 0 : subroutine elpa_hermitian_multiply_dc_c(handle, uplo_a, uplo_c, ncb, a_p, b, nrows_b, &
2208 0 : ncols_b, c, nrows_c, ncols_c, error) &
2209 : bind(C, name="elpa_hermitian_multiply_dc")
2210 : type(c_ptr), intent(in), value :: handle, a_p
2211 : character(1,C_CHAR), value :: uplo_a, uplo_c
2212 : integer(kind=c_int), value :: ncb, nrows_b, ncols_b, nrows_c, ncols_c
2213 : #ifdef USE_FORTRAN2008
2214 : integer(kind=c_int), intent(in), optional :: error
2215 : #else
2216 : integer(kind=c_int), intent(in) :: error
2217 : #endif
2218 :
2219 : complex(kind=c_double_complex), pointer :: a(:, :)
2220 : #ifdef USE_ASSUMED_SIZE
2221 : complex(kind=c_double_complex) :: b(nrows_b,*), c(nrows_c,*)
2222 : #else
2223 : complex(kind=c_double_complex) :: b(nrows_b,ncols_b), c(nrows_c,ncols_c)
2224 : #endif
2225 : type(elpa_impl_t), pointer :: self
2226 :
2227 0 : call c_f_pointer(handle, self)
2228 0 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
2229 :
2230 : call elpa_hermitian_multiply_dc(self, uplo_a, uplo_c, ncb, a, b, nrows_b, &
2231 0 : ncols_b, c, nrows_c, ncols_c, error)
2232 0 : end subroutine
2233 :
2234 :
2235 : !> \brief elpa_hermitian_multiply_fc: class method to perform C : = A**H * B for float complex matrices
2236 : !> where A is a square matrix (self%na,self%na) which is optionally upper or lower triangular
2237 : !> B is a (self%na,ncb) matrix
2238 : !> C is a (self%na,ncb) matrix where optionally only the upper or lower
2239 : !> triangle may be computed
2240 : !>
2241 : !> the MPI commicators and the block-cyclic distribution block size are already known to the type.
2242 : !> Thus the class method "setup" must be called BEFORE this method is used
2243 : !>
2244 : !> \details
2245 : !>
2246 : !> \param self class(elpa_t), the ELPA object
2247 : !> \param uplo_a 'U' if A is upper triangular
2248 : !> 'L' if A is lower triangular
2249 : !> anything else if A is a full matrix
2250 : !> Please note: This pertains to the original A (as set in the calling program)
2251 : !> whereas the transpose of A is used for calculations
2252 : !> If uplo_a is 'U' or 'L', the other triangle is not used at all,
2253 : !> i.e. it may contain arbitrary numbers
2254 : !> \param uplo_c 'U' if only the upper diagonal part of C is needed
2255 : !> 'L' if only the upper diagonal part of C is needed
2256 : !> anything else if the full matrix C is needed
2257 : !> Please note: Even when uplo_c is 'U' or 'L', the other triangle may be
2258 : !> written to a certain extent, i.e. one shouldn't rely on the content there!
2259 : !> \param ncb Number of columns of global matrices B and C
2260 : !> \param a matrix a
2261 : !> \param self%local_nrows number of rows of local (sub) matrix a, set with class method set("local_nrows",value)
2262 : !> \param self%local_ncols number of columns of local (sub) matrix a, set with class method set("local_ncols",value)
2263 : !> \param b matrix b
2264 : !> \param nrows_b number of rows of local (sub) matrix b
2265 : !> \param ncols_b number of columns of local (sub) matrix b
2266 : !> \param c matrix c
2267 : !> \param nrows_c number of rows of local (sub) matrix c
2268 : !> \param ncols_c number of columns of local (sub) matrix c
2269 : !> \param error optional argument, returns an error code, which can be queried with elpa_strerr
2270 192 : subroutine elpa_hermitian_multiply_fc (self,uplo_a, uplo_c, ncb, a, b, nrows_b, ncols_b, &
2271 192 : c, nrows_c, ncols_c, error)
2272 : class(elpa_impl_t) :: self
2273 : character*1 :: uplo_a, uplo_c
2274 : integer(kind=c_int), intent(in) :: nrows_b, ncols_b, nrows_c, ncols_c, ncb
2275 : #ifdef USE_ASSUMED_SIZE
2276 : complex(kind=c_float_complex) :: a(self%local_nrows,*), b(nrows_b,*), c(nrows_c,*)
2277 : #else
2278 : complex(kind=c_float_complex) :: a(self%local_nrows,self%local_ncols), b(nrows_b,ncols_b), c(nrows_c,ncols_c)
2279 : #endif
2280 : #ifdef USE_FORTRAN2008
2281 : integer, optional :: error
2282 : #else
2283 : integer :: error
2284 : #endif
2285 : #ifdef WANT_SINGLE_PRECISION_COMPLEX
2286 : logical :: success_l
2287 :
2288 : success_l = elpa_mult_ah_b_complex_single_impl(self, uplo_a, uplo_c, ncb, a, b, nrows_b, ncols_b, &
2289 192 : c, nrows_c, ncols_c)
2290 : #ifdef USE_FORTRAN2008
2291 192 : if (present(error)) then
2292 192 : if (success_l) then
2293 192 : error = ELPA_OK
2294 : else
2295 0 : error = ELPA_ERROR
2296 : endif
2297 0 : else if (.not. success_l) then
2298 0 : write(error_unit,'(a)') "ELPA: Error in hermitian_multiply() and you did not check for errors!"
2299 : endif
2300 : #else
2301 : if (success_l) then
2302 : error = ELPA_OK
2303 : else
2304 : error = ELPA_ERROR
2305 : endif
2306 : #endif
2307 :
2308 : #else
2309 0 : print *,"This installation of the ELPA library has not been build with single-precision support"
2310 0 : error = ELPA_ERROR
2311 : #endif
2312 192 : end subroutine
2313 :
2314 :
2315 : !c> void elpa_hermitian_multiply_fc(elpa_t handle, char uplo_a, char uplo_c, int ncb, float complex *a, float complex *b, int nrows_b, int ncols_b, float complex *c, int nrows_c, int ncols_c, int *error);
2316 0 : subroutine elpa_hermitian_multiply_fc_c(handle, uplo_a, uplo_c, ncb, a_p, b, nrows_b, &
2317 0 : ncols_b, c, nrows_c, ncols_c, error) &
2318 : bind(C, name="elpa_hermitian_multiply_fc")
2319 : type(c_ptr), intent(in), value :: handle, a_p
2320 : character(1,C_CHAR), value :: uplo_a, uplo_c
2321 : integer(kind=c_int), value :: ncb, nrows_b, ncols_b, nrows_c, ncols_c
2322 : #ifdef USE_FORTRAN2008
2323 : integer(kind=c_int), optional, intent(in) :: error
2324 : #else
2325 : integer(kind=c_int), intent(in) :: error
2326 : #endif
2327 :
2328 : complex(kind=c_float_complex), pointer :: a(:, :)
2329 : #ifdef USE_ASSUMED_SIZE
2330 : complex(kind=c_float_complex) :: b(nrows_b,*), c(nrows_c,*)
2331 : #else
2332 : complex(kind=c_float_complex) :: b(nrows_b,ncols_b), c(nrows_c,ncols_c)
2333 : #endif
2334 : type(elpa_impl_t), pointer :: self
2335 :
2336 0 : call c_f_pointer(handle, self)
2337 0 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
2338 :
2339 : call elpa_hermitian_multiply_fc(self, uplo_a, uplo_c, ncb, a, b, nrows_b, &
2340 0 : ncols_b, c, nrows_c, ncols_c, error)
2341 0 : end subroutine
2342 :
2343 :
2344 : !> \brief elpa_choleksy_d: class method to do a cholesky factorization for a double real matrix
2345 : !>
2346 : !> The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
2347 : !> block size, and the MPI communicators are already known to the object and MUST be set BEFORE
2348 : !> with the class method "setup"
2349 : !>
2350 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
2351 : !> class method "set"
2352 : !>
2353 : !> Parameters
2354 : !>
2355 : !> \param a Distributed matrix for which eigenvalues are to be computed.
2356 : !> Distribution is like in Scalapack.
2357 : !> The full matrix must be set (not only one half like in scalapack).
2358 : !> Destroyed on exit (upper and lower half).
2359 : !>
2360 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
2361 576 : subroutine elpa_cholesky_d (self, a, error)
2362 : class(elpa_impl_t) :: self
2363 : #ifdef USE_ASSUMED_SIZE
2364 : real(kind=rk8) :: a(self%local_nrows,*)
2365 : #else
2366 : real(kind=rk8) :: a(self%local_nrows,self%local_ncols)
2367 : #endif
2368 : #ifdef USE_FORTRAN2008
2369 : integer, optional :: error
2370 : #else
2371 : integer :: error
2372 : #endif
2373 : logical :: success_l
2374 :
2375 576 : success_l = elpa_cholesky_real_double_impl (self, a)
2376 : #ifdef USE_FORTRAN2008
2377 576 : if (present(error)) then
2378 576 : if (success_l) then
2379 576 : error = ELPA_OK
2380 : else
2381 0 : error = ELPA_ERROR
2382 : endif
2383 0 : else if (.not. success_l) then
2384 0 : write(error_unit,'(a)') "ELPA: Error in cholesky() and you did not check for errors!"
2385 : endif
2386 : #else
2387 : if (success_l) then
2388 : error = ELPA_OK
2389 : else
2390 : error = ELPA_ERROR
2391 : endif
2392 : #endif
2393 576 : end subroutine
2394 :
2395 :
2396 : !c> void elpa_cholesky_d(elpa_t handle, double *a, int *error);
2397 0 : subroutine elpa_choleksy_d_c(handle, a_p, error) bind(C, name="elpa_cholesky_d")
2398 : type(c_ptr), intent(in), value :: handle, a_p
2399 : #ifdef USE_FORTRAN2008
2400 : integer(kind=c_int), optional, intent(in) :: error
2401 : #else
2402 : integer(kind=c_int), intent(in) :: error
2403 : #endif
2404 : real(kind=c_double), pointer :: a(:, :)
2405 : type(elpa_impl_t), pointer :: self
2406 :
2407 0 : call c_f_pointer(handle, self)
2408 0 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
2409 :
2410 0 : call elpa_cholesky_d(self, a, error)
2411 0 : end subroutine
2412 :
2413 : !> \brief elpa_choleksy_f: class method to do a cholesky factorization for a float real matrix
2414 : !>
2415 : !> The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
2416 : !> block size, and the MPI communicators are already known to the object and MUST be set BEFORE
2417 : !> with the class method "setup"
2418 : !>
2419 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
2420 : !> class method "set"
2421 : !>
2422 : !> Parameters
2423 : !>
2424 : !> \param a Distributed matrix for which eigenvalues are to be computed.
2425 : !> Distribution is like in Scalapack.
2426 : !> The full matrix must be set (not only one half like in scalapack).
2427 : !> Destroyed on exit (upper and lower half).
2428 : !>
2429 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
2430 288 : subroutine elpa_cholesky_f(self, a, error)
2431 : class(elpa_impl_t) :: self
2432 : #ifdef USE_ASSUMED_SIZE
2433 : real(kind=rk4) :: a(self%local_nrows,*)
2434 : #else
2435 : real(kind=rk4) :: a(self%local_nrows,self%local_ncols)
2436 : #endif
2437 : #ifdef USE_FORTRAN2008
2438 : integer, optional :: error
2439 : #else
2440 : integer :: error
2441 : #endif
2442 : #if WANT_SINGLE_PRECISION_REAL
2443 : logical :: success_l
2444 :
2445 288 : success_l = elpa_cholesky_real_single_impl (self, a)
2446 : #ifdef USE_FORTRAN2008
2447 288 : if (present(error)) then
2448 288 : if (success_l) then
2449 288 : error = ELPA_OK
2450 : else
2451 0 : error = ELPA_ERROR
2452 : endif
2453 0 : else if (.not. success_l) then
2454 0 : write(error_unit,'(a)') "ELPA: Error in cholesky() and you did not check for errors!"
2455 : endif
2456 : #else
2457 : if (success_l) then
2458 : error = ELPA_OK
2459 : else
2460 : error = ELPA_ERROR
2461 : endif
2462 : #endif
2463 : #else
2464 0 : print *,"This installation of the ELPA library has not been build with single-precision support"
2465 0 : error = ELPA_ERROR
2466 : #endif
2467 288 : end subroutine
2468 :
2469 :
2470 : !c> void elpa_cholesky_f(elpa_t handle, float *a, int *error);
2471 0 : subroutine elpa_choleksy_f_c(handle, a_p, error) bind(C, name="elpa_cholesky_f")
2472 : type(c_ptr), intent(in), value :: handle, a_p
2473 : #ifdef USE_FORTRAN2008
2474 : integer(kind=c_int), optional, intent(in) :: error
2475 : #else
2476 : integer(kind=c_int), intent(in) :: error
2477 : #endif
2478 : real(kind=c_float), pointer :: a(:, :)
2479 : type(elpa_impl_t), pointer :: self
2480 :
2481 0 : call c_f_pointer(handle, self)
2482 0 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
2483 :
2484 0 : call elpa_cholesky_f(self, a, error)
2485 0 : end subroutine
2486 :
2487 : !> \brief elpa_choleksy_d: class method to do a cholesky factorization for a double complex matrix
2488 : !>
2489 : !> The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
2490 : !> block size, and the MPI communicators are already known to the object and MUST be set BEFORE
2491 : !> with the class method "setup"
2492 : !>
2493 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
2494 : !> class method "set"
2495 : !>
2496 : !> Parameters
2497 : !>
2498 : !> \param a Distributed matrix for which eigenvalues are to be computed.
2499 : !> Distribution is like in Scalapack.
2500 : !> The full matrix must be set (not only one half like in scalapack).
2501 : !> Destroyed on exit (upper and lower half).
2502 : !>
2503 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
2504 576 : subroutine elpa_cholesky_dc (self, a, error)
2505 : class(elpa_impl_t) :: self
2506 : #ifdef USE_ASSUMED_SIZE
2507 : complex(kind=ck8) :: a(self%local_nrows,*)
2508 : #else
2509 : complex(kind=ck8) :: a(self%local_nrows,self%local_ncols)
2510 : #endif
2511 : #ifdef USE_FORTRAN2008
2512 : integer, optional :: error
2513 : #else
2514 : integer :: error
2515 : #endif
2516 : logical :: success_l
2517 :
2518 576 : success_l = elpa_cholesky_complex_double_impl (self, a)
2519 : #if USE_FORTRAN2008
2520 576 : if (present(error)) then
2521 576 : if (success_l) then
2522 576 : error = ELPA_OK
2523 : else
2524 0 : error = ELPA_ERROR
2525 : endif
2526 0 : else if (.not. success_l) then
2527 0 : write(error_unit,'(a)') "ELPA: Error in cholesky() and you did not check for errors!"
2528 : endif
2529 : #else
2530 : if (success_l) then
2531 : error = ELPA_OK
2532 : else
2533 : error = ELPA_ERROR
2534 : endif
2535 : #endif
2536 576 : end subroutine
2537 :
2538 :
2539 : !c> void elpa_cholesky_dc(elpa_t handle, double complex *a, int *error);
2540 0 : subroutine elpa_choleksy_dc_c(handle, a_p, error) bind(C, name="elpa_cholesky_dc")
2541 : type(c_ptr), intent(in), value :: handle, a_p
2542 : #ifdef USE_FORTRAN2008
2543 : integer(kind=c_int), optional, intent(in) :: error
2544 : #else
2545 : integer(kind=c_int), intent(in) :: error
2546 : #endif
2547 : complex(kind=c_double_complex), pointer :: a(:, :)
2548 : type(elpa_impl_t), pointer :: self
2549 :
2550 0 : call c_f_pointer(handle, self)
2551 0 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
2552 :
2553 0 : call elpa_cholesky_dc(self, a, error)
2554 0 : end subroutine
2555 :
2556 : !> \brief elpa_choleksy_fc: class method to do a cholesky factorization for a float complex matrix
2557 : !>
2558 : !> The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
2559 : !> block size, and the MPI communicators are already known to the object and MUST be set BEFORE
2560 : !> with the class method "setup"
2561 : !>
2562 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
2563 : !> class method "set"
2564 : !>
2565 : !> Parameters
2566 : !>
2567 : !> \param a Distributed matrix for which eigenvalues are to be computed.
2568 : !> Distribution is like in Scalapack.
2569 : !> The full matrix must be set (not only one half like in scalapack).
2570 : !> Destroyed on exit (upper and lower half).
2571 : !>
2572 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
2573 288 : subroutine elpa_cholesky_fc (self, a, error)
2574 : class(elpa_impl_t) :: self
2575 : #ifdef USE_ASSUMED_SIZE
2576 : complex(kind=c_float_complex) :: a(self%local_nrows,*)
2577 : #else
2578 : complex(kind=c_float_complex) :: a(self%local_nrows,self%local_ncols)
2579 : #endif
2580 : #ifdef USE_FORTRAN2008
2581 : integer, optional :: error
2582 : #else
2583 : integer :: error
2584 : #endif
2585 : #if WANT_SINGLE_PRECISION_COMPLEX
2586 : logical :: success_l
2587 :
2588 288 : success_l = elpa_cholesky_complex_single_impl (self, a)
2589 : #ifdef USE_FORTRAN2008
2590 288 : if (present(error)) then
2591 288 : if (success_l) then
2592 288 : error = ELPA_OK
2593 : else
2594 0 : error = ELPA_ERROR
2595 : endif
2596 0 : else if (.not. success_l) then
2597 0 : write(error_unit,'(a)') "ELPA: Error in cholesky() and you did not check for errors!"
2598 : endif
2599 : #else
2600 : if (success_l) then
2601 : error = ELPA_OK
2602 : else
2603 : error = ELPA_ERROR
2604 : endif
2605 : #endif
2606 : #else
2607 0 : print *,"This installation of the ELPA library has not been build with single-precision support"
2608 0 : error = ELPA_ERROR
2609 : #endif
2610 288 : end subroutine
2611 :
2612 :
2613 : !c> void elpa_cholesky_fc(elpa_t handle, float complex *a, int *error);
2614 0 : subroutine elpa_choleksy_fc_c(handle, a_p, error) bind(C, name="elpa_cholesky_fc")
2615 : type(c_ptr), intent(in), value :: handle, a_p
2616 : #ifdef USE_FORTRAN2008
2617 : integer(kind=c_int), optional, intent(in) :: error
2618 : #else
2619 : integer(kind=c_int), intent(in) :: error
2620 : #endif
2621 : complex(kind=c_float_complex), pointer :: a(:, :)
2622 : type(elpa_impl_t), pointer :: self
2623 :
2624 0 : call c_f_pointer(handle, self)
2625 0 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
2626 :
2627 0 : call elpa_cholesky_fc(self, a, error)
2628 0 : end subroutine
2629 :
2630 : !> \brief elpa_invert_trm_d: class method to invert a triangular double real matrix
2631 : !>
2632 : !> The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
2633 : !> block size, and the MPI communicators are already known to the object and MUST be set BEFORE
2634 : !> with the class method "setup"
2635 : !>
2636 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
2637 : !> class method "set"
2638 : !>
2639 : !> Parameters
2640 : !>
2641 : !> \param a Distributed matrix for which eigenvalues are to be computed.
2642 : !> Distribution is like in Scalapack.
2643 : !> The full matrix must be set (not only one half like in scalapack).
2644 : !> Destroyed on exit (upper and lower half).
2645 : !>
2646 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
2647 192 : subroutine elpa_invert_trm_d (self, a, error)
2648 : class(elpa_impl_t) :: self
2649 : #ifdef USE_ASSUMED_SIZE
2650 : real(kind=c_double) :: a(self%local_nrows,*)
2651 : #else
2652 : real(kind=c_double) :: a(self%local_nrows,self%local_ncols)
2653 : #endif
2654 : #ifdef USE_FORTRAN2008
2655 : integer, optional :: error
2656 : #else
2657 : integer :: error
2658 : #endif
2659 : logical :: success_l
2660 :
2661 192 : success_l = elpa_invert_trm_real_double_impl (self, a)
2662 : #ifdef USE_FORTRAN2008
2663 192 : if (present(error)) then
2664 192 : if (success_l) then
2665 192 : error = ELPA_OK
2666 : else
2667 0 : error = ELPA_ERROR
2668 : endif
2669 0 : else if (.not. success_l) then
2670 0 : write(error_unit,'(a)') "ELPA: Error in invert_trm() and you did not check for errors!"
2671 : endif
2672 : #else
2673 : if (success_l) then
2674 : error = ELPA_OK
2675 : else
2676 : error = ELPA_ERROR
2677 : endif
2678 : #endif
2679 192 : end subroutine
2680 :
2681 :
2682 : !c> void elpa_invert_trm_d(elpa_t handle, double *a, int *error);
2683 0 : subroutine elpa_invert_trm_d_c(handle, a_p, error) bind(C, name="elpa_invert_trm_d")
2684 : type(c_ptr), intent(in), value :: handle, a_p
2685 : #ifdef USE_FORTRAN2008
2686 : integer(kind=c_int), optional, intent(in) :: error
2687 : #else
2688 : integer(kind=c_int), intent(in) :: error
2689 : #endif
2690 : real(kind=c_double), pointer :: a(:, :)
2691 : type(elpa_impl_t), pointer :: self
2692 :
2693 0 : call c_f_pointer(handle, self)
2694 0 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
2695 :
2696 0 : call elpa_invert_trm_d(self, a, error)
2697 0 : end subroutine
2698 :
2699 : !> \brief elpa_invert_trm_f: class method to invert a triangular float real matrix
2700 : !>
2701 : !> The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
2702 : !> block size, and the MPI communicators are already known to the object and MUST be set BEFORE
2703 : !> with the class method "setup"
2704 : !>
2705 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
2706 : !> class method "set"
2707 : !>
2708 : !> Parameters
2709 : !>
2710 : !> \param a Distributed matrix for which eigenvalues are to be computed.
2711 : !> Distribution is like in Scalapack.
2712 : !> The full matrix must be set (not only one half like in scalapack).
2713 : !> Destroyed on exit (upper and lower half).
2714 : !>
2715 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
2716 96 : subroutine elpa_invert_trm_f (self, a, error)
2717 : class(elpa_impl_t) :: self
2718 : #ifdef USE_ASSUMED_SIZE
2719 : real(kind=c_float) :: a(self%local_nrows,*)
2720 : #else
2721 : real(kind=c_float) :: a(self%local_nrows,self%local_ncols)
2722 : #endif
2723 : #ifdef USE_FORTRAN2008
2724 : integer, optional :: error
2725 : #else
2726 : integer :: error
2727 : #endif
2728 : #if WANT_SINGLE_PRECISION_REAL
2729 : logical :: success_l
2730 :
2731 96 : success_l = elpa_invert_trm_real_single_impl (self, a)
2732 : #ifdef USE_FORTRAN2008
2733 96 : if (present(error)) then
2734 96 : if (success_l) then
2735 96 : error = ELPA_OK
2736 : else
2737 0 : error = ELPA_ERROR
2738 : endif
2739 0 : else if (.not. success_l) then
2740 0 : write(error_unit,'(a)') "ELPA: Error in invert_trm() and you did not check for errors!"
2741 : endif
2742 : #else
2743 : if (success_l) then
2744 : error = ELPA_OK
2745 : else
2746 : error = ELPA_ERROR
2747 : endif
2748 : #endif
2749 :
2750 : #else
2751 0 : print *,"This installation of the ELPA library has not been build with single-precision support"
2752 0 : error = ELPA_ERROR
2753 : #endif
2754 96 : end subroutine
2755 :
2756 :
2757 : !c> void elpa_invert_trm_f(elpa_t handle, float *a, int *error);
2758 0 : subroutine elpa_invert_trm_f_c(handle, a_p, error) bind(C, name="elpa_invert_trm_f")
2759 : type(c_ptr), intent(in), value :: handle, a_p
2760 : #ifdef USE_FORTRAN2008
2761 : integer(kind=c_int), optional, intent(in) :: error
2762 : #else
2763 : integer(kind=c_int), intent(in) :: error
2764 : #endif
2765 : real(kind=c_float), pointer :: a(:, :)
2766 : type(elpa_impl_t), pointer :: self
2767 :
2768 0 : call c_f_pointer(handle, self)
2769 0 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
2770 :
2771 0 : call elpa_invert_trm_f(self, a, error)
2772 0 : end subroutine
2773 :
2774 : !> \brief elpa_invert_trm_dc: class method to invert a triangular double complex matrix
2775 : !>
2776 : !> The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
2777 : !> block size, and the MPI communicators are already known to the object and MUST be set BEFORE
2778 : !> with the class method "setup"
2779 : !>
2780 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
2781 : !> class method "set"
2782 : !>
2783 : !> Parameters
2784 : !>
2785 : !> \param a Distributed matrix for which eigenvalues are to be computed.
2786 : !> Distribution is like in Scalapack.
2787 : !> The full matrix must be set (not only one half like in scalapack).
2788 : !> Destroyed on exit (upper and lower half).
2789 : !>
2790 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
2791 192 : subroutine elpa_invert_trm_dc (self, a, error)
2792 : class(elpa_impl_t) :: self
2793 : #ifdef USE_ASSUMED_SIZE
2794 : complex(kind=ck8) :: a(self%local_nrows,*)
2795 : #else
2796 : complex(kind=ck8) :: a(self%local_nrows,self%local_ncols)
2797 : #endif
2798 : #ifdef USE_FORTRAN2008
2799 : integer, optional :: error
2800 : #else
2801 : integer :: error
2802 : #endif
2803 : logical :: success_l
2804 :
2805 192 : success_l = elpa_invert_trm_complex_double_impl (self, a)
2806 : #ifdef USE_FORTRAN2008
2807 192 : if (present(error)) then
2808 192 : if (success_l) then
2809 192 : error = ELPA_OK
2810 : else
2811 0 : error = ELPA_ERROR
2812 : endif
2813 0 : else if (.not. success_l) then
2814 0 : write(error_unit,'(a)') "ELPA: Error in invert_trm() and you did not check for errors!"
2815 : endif
2816 : #else
2817 : if (success_l) then
2818 : error = ELPA_OK
2819 : else
2820 : error = ELPA_ERROR
2821 : endif
2822 : #endif
2823 192 : end subroutine
2824 :
2825 :
2826 : !c> void elpa_invert_trm_dc(elpa_t handle, double complex *a, int *error);
2827 0 : subroutine elpa_invert_trm_dc_c(handle, a_p, error) bind(C, name="elpa_invert_trm_dc")
2828 : type(c_ptr), intent(in), value :: handle, a_p
2829 : #ifdef USE_FORTRAN2008
2830 : integer(kind=c_int), optional, intent(in) :: error
2831 : #else
2832 : integer(kind=c_int), intent(in) :: error
2833 : #endif
2834 : complex(kind=c_double_complex), pointer :: a(:, :)
2835 : type(elpa_impl_t), pointer :: self
2836 :
2837 0 : call c_f_pointer(handle, self)
2838 0 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
2839 :
2840 0 : call elpa_invert_trm_dc(self, a, error)
2841 0 : end subroutine
2842 :
2843 : !> \brief elpa_invert_trm_fc: class method to invert a triangular float complex matrix
2844 : !>
2845 : !> The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
2846 : !> block size, and the MPI communicators are already known to the object and MUST be set BEFORE
2847 : !> with the class method "setup"
2848 : !>
2849 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
2850 : !> class method "set"
2851 : !>
2852 : !> Parameters
2853 : !>
2854 : !> \param a Distributed matrix for which eigenvalues are to be computed.
2855 : !> Distribution is like in Scalapack.
2856 : !> The full matrix must be set (not only one half like in scalapack).
2857 : !> Destroyed on exit (upper and lower half).
2858 : !>
2859 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
2860 96 : subroutine elpa_invert_trm_fc (self, a, error)
2861 : class(elpa_impl_t) :: self
2862 : #ifdef USE_ASSUMED_SIZE
2863 : complex(kind=c_float_complex) :: a(self%local_nrows,*)
2864 : #else
2865 : complex(kind=c_float_complex) :: a(self%local_nrows,self%local_ncols)
2866 : #endif
2867 : #ifdef USE_FORTRAN2008
2868 : integer, optional :: error
2869 : #else
2870 : integer :: error
2871 : #endif
2872 : #if WANT_SINGLE_PRECISION_COMPLEX
2873 : logical :: success_l
2874 :
2875 96 : success_l = elpa_invert_trm_complex_single_impl (self, a)
2876 : #ifdef USE_FORTRAN2008
2877 96 : if (present(error)) then
2878 96 : if (success_l) then
2879 96 : error = ELPA_OK
2880 : else
2881 0 : error = ELPA_ERROR
2882 : endif
2883 0 : else if (.not. success_l) then
2884 0 : write(error_unit,'(a)') "ELPA: Error in invert_trm() and you did not check for errors!"
2885 : endif
2886 : #else
2887 : if (success_l) then
2888 : error = ELPA_OK
2889 : else
2890 : error = ELPA_ERROR
2891 : endif
2892 : #endif
2893 : #else
2894 0 : print *,"This installation of the ELPA library has not been build with single-precision support"
2895 0 : error = ELPA_ERROR
2896 : #endif
2897 96 : end subroutine
2898 :
2899 :
2900 : !c> void elpa_invert_trm_fc(elpa_t handle, float complex *a, int *error);
2901 0 : subroutine elpa_invert_trm_fc_c(handle, a_p, error) bind(C, name="elpa_invert_trm_fc")
2902 : type(c_ptr), intent(in), value :: handle, a_p
2903 : #ifdef USE_FORTRAN2008
2904 : integer(kind=c_int), optional, intent(in) :: error
2905 : #else
2906 : integer(kind=c_int), intent(in) :: error
2907 : #endif
2908 : complex(kind=c_float_complex), pointer :: a(:, :)
2909 : type(elpa_impl_t), pointer :: self
2910 :
2911 0 : call c_f_pointer(handle, self)
2912 0 : call c_f_pointer(a_p, a, [self%local_nrows, self%local_ncols])
2913 :
2914 0 : call elpa_invert_trm_fc(self, a, error)
2915 0 : end subroutine
2916 :
2917 :
2918 : !> \brief elpa_solve_tridiagonal_d: class method to solve the eigenvalue problem for a double real tridiagonal matrix a
2919 : !>
2920 : !> The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
2921 : !> block size, and the MPI communicators are already known to the object and MUST be set BEFORE
2922 : !> with the class method "setup"
2923 : !>
2924 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
2925 : !> class method "set"
2926 : !>
2927 : !> Parameters
2928 : !>
2929 : !> \param d array d on input diagonal elements of tridiagonal matrix, on
2930 : !> output the eigenvalues in ascending order
2931 : !> \param e array e on input subdiagonal elements of matrix, on exit destroyed
2932 : !> \param q matrix on exit : contains the eigenvectors
2933 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
2934 : !> \todo e should have dimension (na - 1)
2935 384 : subroutine elpa_solve_tridiagonal_d (self, d, e, q, error)
2936 : class(elpa_impl_t) :: self
2937 : real(kind=rk8) :: d(self%na), e(self%na)
2938 : #ifdef USE_ASSUMED_SIZE
2939 : real(kind=rk8) :: q(self%local_nrows,*)
2940 : #else
2941 : real(kind=rk8) :: q(self%local_nrows,self%local_ncols)
2942 : #endif
2943 : #ifdef USE_FORTRAN2008
2944 : integer, optional :: error
2945 : #else
2946 : integer :: error
2947 : #endif
2948 : logical :: success_l
2949 :
2950 384 : success_l = elpa_solve_tridi_double_impl(self, d, e, q)
2951 : #ifdef USE_FORTRAN2008
2952 384 : if (present(error)) then
2953 384 : if (success_l) then
2954 384 : error = ELPA_OK
2955 : else
2956 0 : error = ELPA_ERROR
2957 : endif
2958 0 : else if (.not. success_l) then
2959 0 : write(error_unit,'(a)') "ELPA: Error in solve_tridiagonal() and you did not check for errors!"
2960 : endif
2961 : #else
2962 : if (success_l) then
2963 : error = ELPA_OK
2964 : else
2965 : error = ELPA_ERROR
2966 : endif
2967 : #endif
2968 384 : end subroutine
2969 :
2970 :
2971 : !> \brief elpa_solve_tridiagonal_f: class method to solve the eigenvalue problem for a float real tridiagonal matrix a
2972 : !>
2973 : !> The dimensions of the matrix a (locally ditributed and global), the block-cylic-distribution
2974 : !> block size, and the MPI communicators are already known to the object and MUST be set BEFORE
2975 : !> with the class method "setup"
2976 : !>
2977 : !> It is possible to change the behaviour of the method by setting tunable parameters with the
2978 : !> class method "set"
2979 : !>
2980 : !> Parameters
2981 : !>
2982 : !> \param d array d on input diagonal elements of tridiagonal matrix, on
2983 : !> output the eigenvalues in ascending order
2984 : !> \param e array e on input subdiagonal elements of matrix, on exit destroyed
2985 : !> \param q matrix on exit : contains the eigenvectors
2986 : !> \param error integer, optional: returns an error code, which can be queried with elpa_strerr
2987 : !> \todo e should have dimension (na - 1)
2988 96 : subroutine elpa_solve_tridiagonal_f (self, d, e, q, error)
2989 : class(elpa_impl_t) :: self
2990 : real(kind=rk4) :: d(self%na), e(self%na)
2991 : #ifdef USE_ASSUMED_SIZE
2992 : real(kind=rk4) :: q(self%local_nrows,*)
2993 : #else
2994 : real(kind=rk4) :: q(self%local_nrows,self%local_ncols)
2995 : #endif
2996 :
2997 : #ifdef USE_FORTRAN2008
2998 : integer, optional :: error
2999 : #else
3000 : integer :: error
3001 : #endif
3002 : #ifdef WANT_SINGLE_PRECISION_REAL
3003 : logical :: success_l
3004 :
3005 96 : success_l = elpa_solve_tridi_single_impl(self, d, e, q)
3006 : #ifdef USE_FORTRAN2008
3007 96 : if (present(error)) then
3008 96 : if (success_l) then
3009 96 : error = ELPA_OK
3010 : else
3011 0 : error = ELPA_ERROR
3012 : endif
3013 0 : else if (.not. success_l) then
3014 0 : write(error_unit,'(a)') "ELPA: Error in solve_tridiagonal() and you did not check for errors!"
3015 : endif
3016 : #else
3017 : if (success_l) then
3018 : error = ELPA_OK
3019 : else
3020 : error = ELPA_ERROR
3021 : endif
3022 : #endif
3023 : #else
3024 0 : print *,"This installation of the ELPA library has not been build with single-precision support"
3025 0 : error = ELPA_ERROR
3026 : #endif
3027 96 : end subroutine
3028 :
3029 :
3030 : !> \brief function to destroy an elpa object
3031 : !> Parameters
3032 : !> \param self class(elpa_impl_t) the allocated ELPA object
3033 21312 : subroutine elpa_destroy(self)
3034 : #ifdef WITH_MPI
3035 : integer :: mpi_comm_rows, mpi_comm_cols, mpierr, error
3036 : #endif
3037 : class(elpa_impl_t) :: self
3038 :
3039 : #ifdef WITH_MPI
3040 14336 : if (self%communicators_owned == 1) then
3041 8320 : call self%get("mpi_comm_rows", mpi_comm_rows,error)
3042 8320 : if (error .ne. ELPA_OK) then
3043 0 : print *,"Problem getting option. Aborting..."
3044 0 : stop
3045 : endif
3046 8320 : call self%get("mpi_comm_cols", mpi_comm_cols,error)
3047 8320 : if (error .ne. ELPA_OK) then
3048 0 : print *,"Problem getting option. Aborting..."
3049 0 : stop
3050 : endif
3051 :
3052 8320 : call mpi_comm_free(mpi_comm_rows, mpierr)
3053 8320 : call mpi_comm_free(mpi_comm_cols, mpierr)
3054 : endif
3055 : #endif
3056 :
3057 21312 : call timer_free(self%timer)
3058 21312 : call timer_free(self%autotune_timer)
3059 21312 : call elpa_index_free_c(self%index)
3060 :
3061 21312 : end subroutine
3062 :
3063 : #define REALCASE 1
3064 : #define DOUBLE_PRECISION 1
3065 : #include "general/precision_macros.h"
3066 : #include "elpa_impl_template.F90"
3067 : #undef REALCASE
3068 : #undef DOUBLE_PRECISION
3069 :
3070 : #ifdef WANT_SINGLE_PRECISION_REAL
3071 : #define REALCASE 1
3072 : #define SINGLE_PRECISION 1
3073 : #include "general/precision_macros.h"
3074 : #include "elpa_impl_template.F90"
3075 : #undef REALCASE
3076 : #undef SINGLE_PRECISION
3077 : #endif /* WANT_SINGLE_PRECISION_REAL */
3078 :
3079 : #define COMPLEXCASE 1
3080 : #define DOUBLE_PRECISION 1
3081 : #include "general/precision_macros.h"
3082 : #include "elpa_impl_template.F90"
3083 : #undef DOUBLE_PRECISION
3084 : #undef COMPLEXCASE
3085 :
3086 : #ifdef WANT_SINGLE_PRECISION_COMPLEX
3087 : #define COMPLEXCASE 1
3088 : #define SINGLE_PRECISION
3089 : #include "general/precision_macros.h"
3090 : #include "elpa_impl_template.F90"
3091 : #undef COMPLEXCASE
3092 : #undef SINGLE_PRECISION
3093 : #endif /* WANT_SINGLE_PRECISION_COMPLEX */
3094 :
3095 : !> \brief function to setup the ELPA autotuning and create the autotune object
3096 : !> Parameters
3097 : !> \param self the allocated ELPA object
3098 : !> \param level integer: the "thoroughness" of the planed autotuning
3099 : !> \param domain integer: the domain (real/complex) which should be tuned
3100 : !> \result tune_state the created autotuning object
3101 0 : function elpa_autotune_setup(self, level, domain, error) result(tune_state)
3102 : class(elpa_impl_t), intent(inout), target :: self
3103 : integer, intent(in) :: level, domain
3104 : type(elpa_autotune_impl_t), pointer :: ts_impl
3105 : class(elpa_autotune_t), pointer :: tune_state
3106 : #ifdef USE_FORTRAN2008
3107 : integer(kind=c_int), optional :: error
3108 : #else
3109 : integer(kind=c_int) :: error
3110 : #endif
3111 :
3112 : #ifdef USE_FORTRAN2008
3113 0 : if (present(error)) then
3114 0 : error = ELPA_OK
3115 : endif
3116 : #else
3117 : error = ELPA_OK
3118 : #endif
3119 0 : if (elpa_get_api_version() < EARLIEST_AUTOTUNE_VERSION) then
3120 0 : write(error_unit, "(a,i0,a)") "ELPA: Error API version: Autotuning does not support ", elpa_get_api_version()
3121 : #ifdef USE_FORTRAN2008
3122 0 : if (present(error)) then
3123 0 : error = ELPA_ERROR
3124 : endif
3125 : #else
3126 : error = ELPA_ERROR
3127 : #endif
3128 0 : return
3129 : endif
3130 :
3131 0 : allocate(ts_impl)
3132 0 : ts_impl%parent => self
3133 0 : ts_impl%level = level
3134 0 : ts_impl%domain = domain
3135 :
3136 0 : ts_impl%i = -1
3137 0 : ts_impl%min_loc = -1
3138 0 : ts_impl%N = elpa_index_autotune_cardinality_c(self%index, level, domain)
3139 :
3140 0 : tune_state => ts_impl
3141 :
3142 0 : call self%autotune_timer%enable()
3143 0 : end function
3144 :
3145 :
3146 :
3147 : !c> /*! \brief C interface for the implementation of the elpa_autotune_setup method
3148 : !c> *
3149 : !c> * \param elpa_t handle: of the ELPA object which should be tuned
3150 : !c> * \param int level: "thoroughness" of autotuning
3151 : !c> * \param int domain: real/complex autotuning
3152 : !c> * \result elpa_autotune_t handle: on the autotune object
3153 : !c> */
3154 : !c> elpa_autotune_t elpa_autotune_setup(elpa_t handle, int level, int domain, int *error);
3155 0 : function elpa_autotune_setup_c(handle ,level, domain, error) result(ptr) bind(C, name="elpa_autotune_setup")
3156 : type(c_ptr), intent(in), value :: handle
3157 : type(elpa_impl_t), pointer :: self
3158 : class(elpa_autotune_t), pointer :: tune_state
3159 : type(elpa_autotune_impl_t), pointer :: obj
3160 : integer(kind=c_int), intent(in), value :: level
3161 : integer(kind=c_int), intent(in), value :: domain
3162 : type(c_ptr) :: ptr
3163 : #ifdef USE_FORTRAN2008
3164 : integer(kind=c_int) , intent(in), optional :: error
3165 : #else
3166 : integer(kind=c_int) , intent(in) :: error
3167 : #endif
3168 :
3169 0 : call c_f_pointer(handle, self)
3170 :
3171 0 : tune_state => self%autotune_setup(level, domain, error)
3172 : select type(tune_state)
3173 : type is (elpa_autotune_impl_t)
3174 0 : obj => tune_state
3175 : class default
3176 0 : print *, "This should not happen"
3177 0 : stop
3178 : end select
3179 0 : ptr = c_loc(obj)
3180 :
3181 0 : end function
3182 :
3183 :
3184 : !> \brief function to do an autotunig step
3185 : !> Parameters
3186 : !> \param self class(elpa_impl_t) the allocated ELPA object
3187 : !> \param tune_state class(elpa_autotune_t): the autotuning object
3188 : !> \result unfinished logical: describes the state of the autotuning (completed/uncompleted)
3189 0 : function elpa_autotune_step(self, tune_state) result(unfinished)
3190 : implicit none
3191 : class(elpa_impl_t), intent(inout) :: self
3192 : class(elpa_autotune_t), intent(inout), target :: tune_state
3193 : type(elpa_autotune_impl_t), pointer :: ts_impl
3194 : logical :: unfinished
3195 : integer :: i
3196 : real(kind=C_DOUBLE) :: time_spent
3197 :
3198 : select type(tune_state)
3199 : type is (elpa_autotune_impl_t)
3200 0 : ts_impl => tune_state
3201 : class default
3202 0 : print *, "This should not happen"
3203 : end select
3204 :
3205 0 : unfinished = .false.
3206 :
3207 0 : if (ts_impl%i >= 0) then
3208 0 : time_spent = self%autotune_timer%get("accumulator")
3209 : !print *, time_spent
3210 0 : if (ts_impl%min_loc == -1 .or. (time_spent < ts_impl%min_val)) then
3211 0 : ts_impl%min_val = time_spent
3212 0 : ts_impl%min_loc = ts_impl%i
3213 : end if
3214 0 : call self%autotune_timer%free()
3215 : endif
3216 :
3217 0 : do while (ts_impl%i < ts_impl%N)
3218 0 : ts_impl%i = ts_impl%i + 1
3219 0 : if (elpa_index_set_autotune_parameters_c(self%index, ts_impl%level, ts_impl%domain, ts_impl%i) == 1) then
3220 0 : unfinished = .true.
3221 0 : return
3222 : end if
3223 : end do
3224 :
3225 0 : end function
3226 :
3227 :
3228 :
3229 : !c> /*! \brief C interface for the implementation of the elpa_autotune_step method
3230 : !c> *
3231 : !c> * \param elpa_t handle: of the ELPA object which should be tuned
3232 : !c> * \param elpa_autotune_t autotune_handle: the autotuning object
3233 : !c> * \result int unfinished: describes whether autotuning finished (0) or not (1)
3234 : !c> */
3235 : !c> int elpa_autotune_step(elpa_t handle, elpa_autotune_t autotune_handle);
3236 0 : function elpa_autotune_step_c(handle, autotune_handle) result(unfinished) bind(C, name="elpa_autotune_step")
3237 : type(c_ptr), intent(in), value :: handle
3238 : type(c_ptr), intent(in), value :: autotune_handle
3239 : type(elpa_impl_t), pointer :: self
3240 : type(elpa_autotune_impl_t), pointer :: tune_state
3241 : logical :: unfinished_f
3242 : integer(kind=c_int) :: unfinished
3243 :
3244 0 : call c_f_pointer(handle, self)
3245 0 : call c_f_pointer(autotune_handle, tune_state)
3246 :
3247 0 : unfinished_f = self%autotune_step(tune_state)
3248 0 : if (unfinished_f) then
3249 0 : unfinished = 1
3250 : else
3251 0 : unfinished = 0
3252 : endif
3253 :
3254 0 : end function
3255 :
3256 :
3257 :
3258 : !> \brief function to set the up-to-know best options of the autotuning
3259 : !> Parameters
3260 : !> \param self class(elpa_impl_t) the allocated ELPA object
3261 : !> \param tune_state class(elpa_autotune_t): the autotuning object
3262 0 : subroutine elpa_autotune_set_best(self, tune_state)
3263 : implicit none
3264 : class(elpa_impl_t), intent(inout) :: self
3265 : class(elpa_autotune_t), intent(in), target :: tune_state
3266 : type(elpa_autotune_impl_t), pointer :: ts_impl
3267 :
3268 : select type(tune_state)
3269 : type is (elpa_autotune_impl_t)
3270 0 : ts_impl => tune_state
3271 : class default
3272 0 : print *, "This should not happen"
3273 : end select
3274 :
3275 0 : print *, "set best, i = ", ts_impl%min_loc, "best time = ", ts_impl%min_val
3276 0 : if (elpa_index_set_autotune_parameters_c(self%index, ts_impl%level, ts_impl%domain, ts_impl%min_loc) /= 1) then
3277 0 : stop "This should not happen (in elpa_autotune_set_best())"
3278 : endif
3279 0 : end subroutine
3280 :
3281 :
3282 :
3283 : !c> /*! \brief C interface for the implementation of the elpa_autotune_set_best method
3284 : !c> *
3285 : !c> * \param elpa_t handle: of the ELPA object which should be tuned
3286 : !c> * \param elpa_autotune_t autotune_handle: the autotuning object
3287 : !c> * \result none
3288 : !c> */
3289 : !c> void elpa_autotune_set_best(elpa_t handle, elpa_autotune_t autotune_handle);
3290 0 : subroutine elpa_autotune_set_best_c(handle, autotune_handle) bind(C, name="elpa_autotune_set_best")
3291 : type(c_ptr), intent(in), value :: handle
3292 : type(c_ptr), intent(in), value :: autotune_handle
3293 : type(elpa_impl_t), pointer :: self
3294 : type(elpa_autotune_impl_t), pointer :: tune_state
3295 :
3296 0 : call c_f_pointer(handle, self)
3297 0 : call c_f_pointer(autotune_handle, tune_state)
3298 :
3299 0 : call self%autotune_set_best(tune_state)
3300 :
3301 0 : end subroutine
3302 :
3303 :
3304 :
3305 0 : end module
|