elpa1_auxiliary.F90 26 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
!    This file is part of ELPA.
!
!    The ELPA library was originally created by the ELPA consortium,
!    consisting of the following organizations:
!
!    - Max Planck Computing and Data Facility (MPCDF), formerly known as
!      Rechenzentrum Garching der Max-Planck-Gesellschaft (RZG),
!    - Bergische Universität Wuppertal, Lehrstuhl für angewandte
!      Informatik,
!    - Technische Universität München, Lehrstuhl für Informatik mit
!      Schwerpunkt Wissenschaftliches Rechnen ,
!    - Fritz-Haber-Institut, Berlin, Abt. Theorie,
13
!    - Max-Plack-Institut für Mathematik in den Naturwissenschaften,
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
!      Leipzig, Abt. Komplexe Strukutren in Biologie und Kognition,
!      and
!    - IBM Deutschland GmbH
!
!    This particular source code file contains additions, changes and
!    enhancements authored by Intel Corporation which is not part of
!    the ELPA consortium.
!
!    More information can be found here:
!    http://elpa.mpcdf.mpg.de/
!
!    ELPA is free software: you can redistribute it and/or modify
!    it under the terms of the version 3 of the license of the
!    GNU Lesser General Public License as published by the Free
!    Software Foundation.
!
!    ELPA is distributed in the hope that it will be useful,
!    but WITHOUT ANY WARRANTY; without even the implied warranty of
!    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
!    GNU Lesser General Public License for more details.
!
!    You should have received a copy of the GNU Lesser General Public License
!    along with ELPA.  If not, see <http://www.gnu.org/licenses/>
!
!    ELPA reflects a substantial effort on the part of the original
!    ELPA consortium, and we ask you to respect the spirit of the
!    license that we chose: i.e., please contribute any changes you
!    may have back to the original ELPA library distribution, and keep
!    any derivatives of ELPA under the same license that we chose for
!    the original distribution, the GNU Lesser General Public License.
!
!
! ELPA1 -- Faster replacements for ScaLAPACK symmetric eigenvalue routines
!
! Copyright of the original code rests with the authors inside the ELPA
! consortium. The copyright of any additional modifications shall rest
! with their original authors, but shall adhere to the licensing terms
! distributed along with the original code in the file "COPYING".
Andreas Marek's avatar
Andreas Marek committed
52
53
!
! This file has been rewritten by A. Marek, MPCDF
54
55
#include "config-f90.h"

Andreas Marek's avatar
Andreas Marek committed
56
!> \brief Fortran module which provides helper routines for matrix calculations
57
module ELPA1_AUXILIARY
58
  use elpa_utilities
Andreas Marek's avatar
Andreas Marek committed
59

60
61
  implicit none

62
  public :: elpa_mult_at_b_real_double      !< Multiply double-precision real matrices A**T * B
63

64
  public :: elpa_mult_ah_b_complex_double   !< Multiply double-precision complex matrices A**H * B
65

66
  public :: elpa_invert_trm_real_double     !< Invert double-precision real triangular matrix
67

68
  public :: elpa_invert_trm_complex_double  !< Invert double-precision complex triangular matrix
69

70
  public :: elpa_cholesky_real_double       !< Cholesky factorization of a double-precision real matrix
71

72
  public :: elpa_cholesky_complex_double    !< Cholesky factorization of a double-precision complex matrix
73

74
  public :: elpa_solve_tridi_double         !< Solve tridiagonal eigensystem for a double-precision matrix with divide and conquer method
75
76

#ifdef WANT_SINGLE_PRECISION_REAL
77
78
79
80
  public :: elpa_cholesky_real_single       !< Cholesky factorization of a single-precision real matrix
  public :: elpa_invert_trm_real_single     !< Invert single-precision real triangular matrix
  public :: elpa_mult_at_b_real_single      !< Multiply single-precision real matrices A**T * B
  public :: elpa_solve_tridi_single         !< Solve tridiagonal eigensystem for a single-precision matrix with divide and conquer method
81
82
83
#endif

#ifdef WANT_SINGLE_PRECISION_COMPLEX
84
85
86
  public :: elpa_cholesky_complex_single    !< Cholesky factorization of a single-precision complex matrix
  public :: elpa_invert_trm_complex_single  !< Invert single-precision complex triangular matrix
  public :: elpa_mult_ah_b_complex_single   !< Multiply single-precision complex matrices A**H * B
87
88
#endif

89
90
  contains

91

Andreas Marek's avatar
Andreas Marek committed
92
93
#define REALCASE 1
#define DOUBLE_PRECISION
94
#include "../../general/precision_macros.h"
Andreas Marek's avatar
Andreas Marek committed
95

96

97
98
    function elpa_cholesky_real_double(na, a, lda, nblk, matrixCols, mpi_comm_rows, mpi_comm_cols, &
                                             wantDebug) result(success)
99
#include "./elpa_cholesky_template.F90"
100
    end function
101

102
#ifdef WANT_SINGLE_PRECISION_REAL
Andreas Marek's avatar
Andreas Marek committed
103
104
#define REALCASE 1
#define SINGLE_PRECISION
105
#include "../../general/precision_macros.h"
106

107
!> \brief  cholesky_real_single: Cholesky factorization of a single-precision real symmetric matrix
108
!> \details
109
!>
110
!> \param  na                   Order of matrix
111
!> \param  a(lda,matrixCols)    Distributed matrix which should be factorized.
112
!>                              Distribution is like in Scalapack.
Andreas Marek's avatar
Andreas Marek committed
113
!>                              Only upper triangle needs to be set.
114
115
!>                              On return, the upper triangle contains the Cholesky factor
!>                              and the lower triangle is set to 0.
116
117
!> \param  lda                  Leading dimension of a
!> \param  nblk                 blocksize of cyclic distribution, must be the same in both directions!
Andreas Marek's avatar
Andreas Marek committed
118
!> \param  matrixCols           local columns of matrix a
119
120
121
122
123
!> \param  mpi_comm_rows        MPI communicator for rows
!> \param  mpi_comm_cols        MPI communicator for columns
!> \param wantDebug             logical, more debug information on failure
!> \param succes                logical, reports success or failure

124
    function elpa_cholesky_real_single(na, a, lda, nblk, matrixCols, mpi_comm_rows, mpi_comm_cols, &
125
                                            wantDebug) result(success)
126
#include "./elpa_cholesky_template.F90"
127
    end function
128

129
#endif /* WANT_SINGLE_PRECSION_REAL */
130

Andreas Marek's avatar
Andreas Marek committed
131
132
#define REALCASE 1
#define DOUBLE_PRECISION
133
#include "../../general/precision_macros.h"
134
!> \brief  elpa_invert_trm_real_double: Inverts a double-precision real upper triangular matrix
135
!> \details
136
137
138
!> \param  na                   Order of matrix
!> \param  a(lda,matrixCols)    Distributed matrix which should be inverted
!>                              Distribution is like in Scalapack.
Andreas Marek's avatar
Andreas Marek committed
139
!>                              Only upper triangle needs to be set.
140
141
142
!>                              The lower triangle is not referenced.
!> \param  lda                  Leading dimension of a
!> \param  nblk                 blocksize of cyclic distribution, must be the same in both directions!
Andreas Marek's avatar
Andreas Marek committed
143
!> \param  matrixCols           local columns of matrix a
144
145
146
!> \param  mpi_comm_rows        MPI communicator for rows
!> \param  mpi_comm_cols        MPI communicator for columns
!> \param wantDebug             logical, more debug information on failure
147
!> \result succes               logical, reports success or failure
148
     function elpa_invert_trm_real_double(na, a, lda, nblk, matrixCols, mpi_comm_rows, mpi_comm_cols, wantDebug) result(success)
149
#include "./elpa_invert_trm.F90"
150
     end function
151
152

#if WANT_SINGLE_PRECISION_REAL
Andreas Marek's avatar
Andreas Marek committed
153
154
#define REALCASE 1
#define SINGLE_PRECISION
155
#include "../../general/precision_macros.h"
Andreas Marek's avatar
Andreas Marek committed
156

157
!> \brief  elpa_invert_trm_real_single: Inverts a single-precision real upper triangular matrix
158
159
160
161
!> \details
!> \param  na                   Order of matrix
!> \param  a(lda,matrixCols)    Distributed matrix which should be inverted
!>                              Distribution is like in Scalapack.
162
!>                              Only upper triangle needs to be set.
163
164
165
166
167
168
169
!>                              The lower triangle is not referenced.
!> \param  lda                  Leading dimension of a
!> \param                       matrixCols  local columns of matrix a
!> \param  nblk                 blocksize of cyclic distribution, must be the same in both directions!
!> \param  mpi_comm_rows        MPI communicator for rows
!> \param  mpi_comm_cols        MPI communicator for columns
!> \param wantDebug             logical, more debug information on failure
170
171
!> \result succes               logical, reports success or failure
    function elpa_invert_trm_real_single(na, a, lda, nblk, matrixCols, mpi_comm_rows, mpi_comm_cols, wantDebug) result(success)
172
#include "./elpa_invert_trm.F90"
173
    end function
174
175
176

#endif /* WANT_SINGLE_PRECISION_REAL */

177

Andreas Marek's avatar
Andreas Marek committed
178
179
#define COMPLEXCASE 1
#define DOUBLE_PRECISION
180
#include "../../general/precision_macros.h"
Andreas Marek's avatar
Andreas Marek committed
181

182
!> \brief  elpa_cholesky_complex_double: Cholesky factorization of a double-precision complex hermitian matrix
183
184
185
186
!> \details
!> \param  na                   Order of matrix
!> \param  a(lda,matrixCols)    Distributed matrix which should be factorized.
!>                              Distribution is like in Scalapack.
187
!>                              Only upper triangle needs to be set.
188
189
190
191
192
193
194
195
!>                              On return, the upper triangle contains the Cholesky factor
!>                              and the lower triangle is set to 0.
!> \param  lda                  Leading dimension of a
!> \param                       matrixCols  local columns of matrix a
!> \param  nblk                 blocksize of cyclic distribution, must be the same in both directions!
!> \param  mpi_comm_rows        MPI communicator for rows
!> \param  mpi_comm_cols        MPI communicator for columns
!> \param wantDebug             logical, more debug information on failure
196
197
!> \result succes               logical, reports success or failure
    function elpa_cholesky_complex_double(na, a, lda, nblk, matrixCols, mpi_comm_rows, mpi_comm_cols, wantDebug) result(success)
198
#include "./elpa_cholesky_template.F90"
199
    end function
200

Andreas Marek's avatar
Andreas Marek committed
201

202
#ifdef WANT_SINGLE_PRECISION_COMPLEX
Andreas Marek's avatar
Andreas Marek committed
203
204
#define COMPLEXCASE 1
#define SINGLE_PRECISION
205
#include "../../general/precision_macros.h"
Andreas Marek's avatar
Andreas Marek committed
206

207
!> \brief  elpa_cholesky_complex_single: Cholesky factorization of a single-precision complex hermitian matrix
208
209
210
211
!> \details
!> \param  na                   Order of matrix
!> \param  a(lda,matrixCols)    Distributed matrix which should be factorized.
!>                              Distribution is like in Scalapack.
212
!>                              Only upper triangle needs to be set.
213
214
215
216
217
218
219
220
!>                              On return, the upper triangle contains the Cholesky factor
!>                              and the lower triangle is set to 0.
!> \param  lda                  Leading dimension of a
!> \param                       matrixCols  local columns of matrix a
!> \param  nblk                 blocksize of cyclic distribution, must be the same in both directions!
!> \param  mpi_comm_rows        MPI communicator for rows
!> \param  mpi_comm_cols        MPI communicator for columns
!> \param wantDebug             logical, more debug information on failure
221
222
!> \result succes               logical, reports success or failure
    function elpa_cholesky_complex_single(na, a, lda, nblk, matrixCols, mpi_comm_rows, mpi_comm_cols, wantDebug) result(success)
223
#include "./elpa_cholesky_template.F90"
224
    end function
225
226

#endif /* WANT_SINGLE_PRECISION_COMPLEX */
227

Andreas Marek's avatar
Andreas Marek committed
228
229
#define COMPLEXCASE 1
#define DOUBLE_PRECISION
230
#include "../../general/precision_macros.h"
231

232
!> \brief  elpa_invert_trm_complex_double: Inverts a double-precision complex upper triangular matrix
233
234
235
236
!> \details
!> \param  na                   Order of matrix
!> \param  a(lda,matrixCols)    Distributed matrix which should be inverted
!>                              Distribution is like in Scalapack.
237
!>                              Only upper triangle needs to be set.
238
239
240
241
242
243
244
!>                              The lower triangle is not referenced.
!> \param  lda                  Leading dimension of a
!> \param                       matrixCols  local columns of matrix a
!> \param  nblk                 blocksize of cyclic distribution, must be the same in both directions!
!> \param  mpi_comm_rows        MPI communicator for rows
!> \param  mpi_comm_cols        MPI communicator for columns
!> \param wantDebug             logical, more debug information on failure
245
246
!> \result succes               logical, reports success or failure

247
    function elpa_invert_trm_complex_double(na, a, lda, nblk, matrixCols, mpi_comm_rows, mpi_comm_cols, wantDebug) result(success)
248
#include "./elpa_invert_trm.F90"
249
    end function
250
251

#ifdef WANT_SINGLE_PRECISION_COMPLEX
Andreas Marek's avatar
Andreas Marek committed
252
253
#define COMPLEXCASE 1
#define SINGLE_PRECISION
254
#include "../../general/precision_macros.h"
Andreas Marek's avatar
Andreas Marek committed
255

256
!> \brief  elpa_invert_trm_complex_single: Inverts a single-precision complex upper triangular matrix
257
258
259
260
!> \details
!> \param  na                   Order of matrix
!> \param  a(lda,matrixCols)    Distributed matrix which should be inverted
!>                              Distribution is like in Scalapack.
261
!>                              Only upper triangle needs to be set.
262
263
264
265
266
267
268
!>                              The lower triangle is not referenced.
!> \param  lda                  Leading dimension of a
!> \param                       matrixCols  local columns of matrix a
!> \param  nblk                 blocksize of cyclic distribution, must be the same in both directions!
!> \param  mpi_comm_rows        MPI communicator for rows
!> \param  mpi_comm_cols        MPI communicator for columns
!> \param wantDebug             logical, more debug information on failure
269
270
271
!> \result succes               logical, reports success or failure

    function elpa_invert_trm_complex_single(na, a, lda, nblk, matrixCols, mpi_comm_rows, mpi_comm_cols, wantDebug) result(success)
272
#include "./elpa_invert_trm.F90"
273
    end function
274

275
#endif /* WANT_SINGE_PRECISION_COMPLEX */
276

Andreas Marek's avatar
Andreas Marek committed
277
278
#define REALCASE 1
#define DOUBLE_PRECISION
279
#include "../../general/precision_macros.h"
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
!> \brief  mult_at_b_real_double: Performs C : = A**T * B
!>         where   A is a square matrix (na,na) which is optionally upper or lower triangular
!>                 B is a (na,ncb) matrix
!>                 C is a (na,ncb) matrix where optionally only the upper or lower
!>                   triangle may be computed
!> \details

!> \param  uplo_a               'U' if A is upper triangular
!>                              'L' if A is lower triangular
!>                              anything else if A is a full matrix
!>                              Please note: This pertains to the original A (as set in the calling program)
!>                                           whereas the transpose of A is used for calculations
!>                              If uplo_a is 'U' or 'L', the other triangle is not used at all,
!>                              i.e. it may contain arbitrary numbers
!> \param uplo_c                'U' if only the upper diagonal part of C is needed
!>                              'L' if only the upper diagonal part of C is needed
!>                              anything else if the full matrix C is needed
!>                              Please note: Even when uplo_c is 'U' or 'L', the other triangle may be
!>                                            written to a certain extent, i.e. one shouldn't rely on the content there!
!> \param na                    Number of rows/columns of A, number of rows of B and C
!> \param ncb                   Number of columns  of B and C
!> \param a                     matrix a
!> \param lda                   leading dimension of matrix a
!> \param b                     matrix b
!> \param ldb                   leading dimension of matrix b
!> \param nblk                  blocksize of cyclic distribution, must be the same in both directions!
!> \param  mpi_comm_rows        MPI communicator for rows
!> \param  mpi_comm_cols        MPI communicator for columns
!> \param c                     matrix c
!> \param ldc                   leading dimension of matrix c
310
!> \result success
311

312
313
    function elpa_mult_at_b_real_double(uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, nblk, &
                              mpi_comm_rows, mpi_comm_cols, c, ldc, ldcCols) result(success)
314
#include "./elpa_multiply_a_b.F90"
315
    end function elpa_mult_at_b_real_double
316
317

#if WANT_SINGLE_PRECISION_REAL
Andreas Marek's avatar
Andreas Marek committed
318
319
#define REALCASE 1
#define SINGLE_PRECISION
320
#include "../../general/precision_macros.h"
Andreas Marek's avatar
Andreas Marek committed
321

322
!> \brief  elpa_mult_at_b_real_single: Performs C : = A**T * B
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
!>         where   A is a square matrix (na,na) which is optionally upper or lower triangular
!>                 B is a (na,ncb) matrix
!>                 C is a (na,ncb) matrix where optionally only the upper or lower
!>                   triangle may be computed
!> \details

!> \param  uplo_a               'U' if A is upper triangular
!>                              'L' if A is lower triangular
!>                              anything else if A is a full matrix
!>                              Please note: This pertains to the original A (as set in the calling program)
!>                                           whereas the transpose of A is used for calculations
!>                              If uplo_a is 'U' or 'L', the other triangle is not used at all,
!>                              i.e. it may contain arbitrary numbers
!> \param uplo_c                'U' if only the upper diagonal part of C is needed
!>                              'L' if only the upper diagonal part of C is needed
!>                              anything else if the full matrix C is needed
!>                              Please note: Even when uplo_c is 'U' or 'L', the other triangle may be
!>                                            written to a certain extent, i.e. one shouldn't rely on the content there!
!> \param na                    Number of rows/columns of A, number of rows of B and C
!> \param ncb                   Number of columns  of B and C
!> \param a                     matrix a
!> \param lda                   leading dimension of matrix a
!> \param b                     matrix b
!> \param ldb                   leading dimension of matrix b
!> \param nblk                  blocksize of cyclic distribution, must be the same in both directions!
!> \param  mpi_comm_rows        MPI communicator for rows
!> \param  mpi_comm_cols        MPI communicator for columns
!> \param c                     matrix c
!> \param ldc                   leading dimension of matrix c
352
!> \result success
353

354
355
    function elpa_mult_at_b_real_single(uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, nblk, &
                              mpi_comm_rows, mpi_comm_cols, c, ldc, ldcCols) result(success)
Andreas Marek's avatar
Andreas Marek committed
356

357
#include "./elpa_multiply_a_b.F90"
358

359
    end function elpa_mult_at_b_real_single
360

361
362
#endif /* WANT_SINGLE_PRECISION_REAL */

363

Andreas Marek's avatar
Andreas Marek committed
364
365
#define COMPLEXCASE 1
#define DOUBLE_PRECISION
366
#include "../../general/precision_macros.h"
Andreas Marek's avatar
Andreas Marek committed
367

368
!> \brief  elpa_mult_ah_b_complex_double: Performs C : = A**H * B
369
370
371
372
!>         where   A is a square matrix (na,na) which is optionally upper or lower triangular
!>                 B is a (na,ncb) matrix
!>                 C is a (na,ncb) matrix where optionally only the upper or lower
!>                   triangle may be computed
373
!> \details
374
!>
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
!> \param  uplo_a               'U' if A is upper triangular
!>                              'L' if A is lower triangular
!>                              anything else if A is a full matrix
!>                              Please note: This pertains to the original A (as set in the calling program)
!>                                           whereas the transpose of A is used for calculations
!>                              If uplo_a is 'U' or 'L', the other triangle is not used at all,
!>                              i.e. it may contain arbitrary numbers
!> \param uplo_c                'U' if only the upper diagonal part of C is needed
!>                              'L' if only the upper diagonal part of C is needed
!>                              anything else if the full matrix C is needed
!>                              Please note: Even when uplo_c is 'U' or 'L', the other triangle may be
!>                                            written to a certain extent, i.e. one shouldn't rely on the content there!
!> \param na                    Number of rows/columns of A, number of rows of B and C
!> \param ncb                   Number of columns  of B and C
!> \param a                     matrix a
!> \param lda                   leading dimension of matrix a
391
!> \param ldaCols               columns of matrix a
392
393
!> \param b                     matrix b
!> \param ldb                   leading dimension of matrix b
394
!> \param ldbCols               columns of matrix b
395
396
397
398
399
!> \param nblk                  blocksize of cyclic distribution, must be the same in both directions!
!> \param  mpi_comm_rows        MPI communicator for rows
!> \param  mpi_comm_cols        MPI communicator for columns
!> \param c                     matrix c
!> \param ldc                   leading dimension of matrix c
400
!> \result success
401

402
403
    function elpa_mult_ah_b_complex_double(uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, nblk, &
                                 mpi_comm_rows, mpi_comm_cols, c, ldc, ldcCols) result(success)
404
#include "./elpa_multiply_a_b.F90"
405

406
    end function elpa_mult_ah_b_complex_double
407

408
#ifdef WANT_SINGLE_PRECISION_COMPLEX
Andreas Marek's avatar
Andreas Marek committed
409
410
#define COMPLEXCASE 1
#define SINGLE_PRECISION
411
#include "../../general/precision_macros.h"
412

413
!> \brief  elpa_mult_ah_b_complex_single: Performs C : = A**H * B
414
415
416
417
!>         where   A is a square matrix (na,na) which is optionally upper or lower triangular
!>                 B is a (na,ncb) matrix
!>                 C is a (na,ncb) matrix where optionally only the upper or lower
!>                   triangle may be computed
418
!> \details
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
!>
!> \param  uplo_a               'U' if A is upper triangular
!>                              'L' if A is lower triangular
!>                              anything else if A is a full matrix
!>                              Please note: This pertains to the original A (as set in the calling program)
!>                                           whereas the transpose of A is used for calculations
!>                              If uplo_a is 'U' or 'L', the other triangle is not used at all,
!>                              i.e. it may contain arbitrary numbers
!> \param uplo_c                'U' if only the upper diagonal part of C is needed
!>                              'L' if only the upper diagonal part of C is needed
!>                              anything else if the full matrix C is needed
!>                              Please note: Even when uplo_c is 'U' or 'L', the other triangle may be
!>                                            written to a certain extent, i.e. one shouldn't rely on the content there!
!> \param na                    Number of rows/columns of A, number of rows of B and C
!> \param ncb                   Number of columns  of B and C
!> \param a                     matrix a
!> \param lda                   leading dimension of matrix a
436
!> \param ldaCols               columns of matrix a
437
438
!> \param b                     matrix b
!> \param ldb                   leading dimension of matrix b
439
!> \param ldbCols               columns of matrix b
440
441
442
443
444
!> \param nblk                  blocksize of cyclic distribution, must be the same in both directions!
!> \param  mpi_comm_rows        MPI communicator for rows
!> \param  mpi_comm_cols        MPI communicator for columns
!> \param c                     matrix c
!> \param ldc                   leading dimension of matrix c
445
!> \result success
446

447
448
    function elpa_mult_ah_b_complex_single(uplo_a, uplo_c, na, ncb, a, lda, ldaCols, b, ldb, ldbCols, nblk, &
                                 mpi_comm_rows, mpi_comm_cols, c, ldc, ldcCols) result(success)
449

450
#include "./elpa_multiply_a_b.F90"
451

452
    end function elpa_mult_ah_b_complex_single
453
454
455

#endif /* WANT_SINGLE_PRECISION_COMPLEX */

456
457
#define REALCASE 1
#define DOUBLE_PRECISION
458
#include "../../general/precision_macros.h"
459
460

!> \brief  elpa_solve_tridi_double: Solve tridiagonal eigensystem for a double-precision matrix with divide and conquer method
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
!> \details
!>
!> \param na                    Matrix dimension
!> \param nev                   number of eigenvalues/vectors to be computed
!> \param d                     array d(na) on input diagonal elements of tridiagonal matrix, on
!>                              output the eigenvalues in ascending order
!> \param e                     array e(na) on input subdiagonal elements of matrix, on exit destroyed
!> \param q                     on exit : matrix q(ldq,matrixCols) contains the eigenvectors
!> \param ldq                   leading dimension of matrix q
!> \param nblk                  blocksize of cyclic distribution, must be the same in both directions!
!> \param matrixCols            columns of matrix q
!> \param mpi_comm_rows         MPI communicator for rows
!> \param mpi_comm_cols         MPI communicator for columns
!> \param wantDebug             logical, give more debug information if .true.
!> \result success              logical, .true. on success, else .false.

477
478
    function elpa_solve_tridi_double(na, nev, d, e, q, ldq, nblk, matrixCols, mpi_comm_rows, mpi_comm_cols, wantDebug) &
          result(success)
479

480
#include "./elpa_solve_tridi.F90"
481

482
    end function
483
484


485
#ifdef WANT_SINGLE_PRECISION_REAL
486
487
#define REALCASE 1
#define SINGLE_PRECISION
488
#include "../../general/precision_macros.h"
489

490
!> \brief  elpa_solve_tridi_single: Solve tridiagonal eigensystem for a single-precision matrix with divide and conquer method
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
!> \details
!>
!> \param na                    Matrix dimension
!> \param nev                   number of eigenvalues/vectors to be computed
!> \param d                     array d(na) on input diagonal elements of tridiagonal matrix, on
!>                              output the eigenvalues in ascending order
!> \param e                     array e(na) on input subdiagonal elements of matrix, on exit destroyed
!> \param q                     on exit : matrix q(ldq,matrixCols) contains the eigenvectors
!> \param ldq                   leading dimension of matrix q
!> \param nblk                  blocksize of cyclic distribution, must be the same in both directions!
!> \param matrixCols            columns of matrix q
!> \param mpi_comm_rows         MPI communicator for rows
!> \param mpi_comm_cols         MPI communicator for columns
!> \param wantDebug             logical, give more debug information if .true.
!> \result success              logical, .true. on success, else .false.

507
508
    function elpa_solve_tridi_single(na, nev, d, e, q, ldq, nblk, matrixCols, mpi_comm_rows, &
                                     mpi_comm_cols, wantDebug) result(success)
509

510
#include "./elpa_solve_tridi.F90"
511
512
513

    end function

514
515
#endif /* WANT_SINGLE_PRECISION_REAL */

516
517
518



519
end module elpa1_auxiliary