Line data Source code
1 : ! This file is part of ELPA.
2 : !
3 : ! The ELPA library was originally created by the ELPA consortium,
4 : ! consisting of the following organizations:
5 : !
6 : ! - Max Planck Computing and Data Facility (MPCDF), formerly known as
7 : ! Rechenzentrum Garching der Max-Planck-Gesellschaft (RZG),
8 : ! - Bergische Universität Wuppertal, Lehrstuhl für angewandte
9 : ! Informatik,
10 : ! - Technische Universität München, Lehrstuhl für Informatik mit
11 : ! Schwerpunkt Wissenschaftliches Rechnen ,
12 : ! - Fritz-Haber-Institut, Berlin, Abt. Theorie,
13 : ! - Max-Plack-Institut für Mathematik in den Naturwissenschaften,
14 : ! Leipzig, Abt. Komplexe Strukutren in Biologie und Kognition,
15 : ! and
16 : ! - IBM Deutschland GmbH
17 : !
18 : ! This particular source code file contains additions, changes and
19 : ! enhancements authored by Intel Corporation which is not part of
20 : ! the ELPA consortium.
21 : !
22 : ! More information can be found here:
23 : ! http://elpa.mpcdf.mpg.de/
24 : !
25 : ! ELPA is free software: you can redistribute it and/or modify
26 : ! it under the terms of the version 3 of the license of the
27 : ! GNU Lesser General Public License as published by the Free
28 : ! Software Foundation.
29 : !
30 : ! ELPA is distributed in the hope that it will be useful,
31 : ! but WITHOUT ANY WARRANTY; without even the implied warranty of
32 : ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
33 : ! GNU Lesser General Public License for more details.
34 : !
35 : ! You should have received a copy of the GNU Lesser General Public License
36 : ! along with ELPA. If not, see <http://www.gnu.org/licenses/>
37 : !
38 : ! ELPA reflects a substantial effort on the part of the original
39 : ! ELPA consortium, and we ask you to respect the spirit of the
40 : ! license that we chose: i.e., please contribute any changes you
41 : ! may have back to the original ELPA library distribution, and keep
42 : ! any derivatives of ELPA under the same license that we chose for
43 : ! the original distribution, the GNU Lesser General Public License.
44 : !
45 : !
46 : ! ELPA1 -- Faster replacements for ScaLAPACK symmetric eigenvalue routines
47 : !
48 : ! Copyright of the original code rests with the authors inside the ELPA
49 : ! consortium. The copyright of any additional modifications shall rest
50 : ! with their original authors, but shall adhere to the licensing terms
51 : ! distributed along with the original code in the file "COPYING".
52 : !
53 : ! This file has been rewritten by A. Marek, MPCDF
54 : #include "config-f90.h"
55 :
56 : !> \brief Fortran module which provides helper routines for matrix calculations
57 : module elpa1_auxiliary_impl
58 : use elpa_utilities
59 :
60 : implicit none
61 :
62 : public :: elpa_mult_at_b_real_double_impl !< Multiply double-precision real matrices A**T * B
63 :
64 : public :: elpa_mult_ah_b_complex_double_impl !< Multiply double-precision complex matrices A**H * B
65 :
66 : public :: elpa_invert_trm_real_double_impl !< Invert double-precision real triangular matrix
67 :
68 : public :: elpa_invert_trm_complex_double_impl !< Invert double-precision complex triangular matrix
69 :
70 : public :: elpa_cholesky_real_double_impl !< Cholesky factorization of a double-precision real matrix
71 :
72 : public :: elpa_cholesky_complex_double_impl !< Cholesky factorization of a double-precision complex matrix
73 :
74 : public :: elpa_solve_tridi_double_impl !< Solve tridiagonal eigensystem for a double-precision matrix with divide and conquer method
75 :
76 : #ifdef WANT_SINGLE_PRECISION_REAL
77 : public :: elpa_cholesky_real_single_impl !< Cholesky factorization of a single-precision real matrix
78 : public :: elpa_invert_trm_real_single_impl !< Invert single-precision real triangular matrix
79 : public :: elpa_mult_at_b_real_single_impl !< Multiply single-precision real matrices A**T * B
80 : public :: elpa_solve_tridi_single_impl !< Solve tridiagonal eigensystem for a single-precision matrix with divide and conquer method
81 : #endif
82 :
83 : #ifdef WANT_SINGLE_PRECISION_COMPLEX
84 : public :: elpa_cholesky_complex_single_impl !< Cholesky factorization of a single-precision complex matrix
85 : public :: elpa_invert_trm_complex_single_impl !< Invert single-precision complex triangular matrix
86 : public :: elpa_mult_ah_b_complex_single_impl !< Multiply single-precision complex matrices A**H * B
87 : #endif
88 :
89 : contains
90 :
91 : #define REALCASE 1
92 : #define DOUBLE_PRECISION
93 : #include "../general/precision_macros.h"
94 :
95 864 : function elpa_cholesky_real_double_impl (obj, a) result(success)
96 : #include "elpa_cholesky_template.F90"
97 :
98 864 : end function elpa_cholesky_real_double_impl
99 : #undef DOUBLE_PRECISION
100 : #undef REALCASE
101 :
102 : #ifdef WANT_SINGLE_PRECISION_REAL
103 : #define REALCASE 1
104 : #define SINGLE_PRECISION
105 : #include "../general/precision_macros.h"
106 :
107 432 : function elpa_cholesky_real_single_impl(obj, a) result(success)
108 : #include "elpa_cholesky_template.F90"
109 :
110 432 : end function elpa_cholesky_real_single_impl
111 : #undef SINGLE_PRECISION
112 : #undef REALCASE
113 :
114 : #endif /* WANT_SINGLE_PRECSION_REAL */
115 :
116 : #define REALCASE 1
117 : #define DOUBLE_PRECISION
118 : #include "../general/precision_macros.h"
119 : !> \brief elpa_invert_trm_real_double: Inverts a double-precision real upper triangular matrix
120 : !> \details
121 : !> \param obj elpa_t object contains:
122 : !> \param - obj%na Order of matrix
123 : !> \param - obj%local_nrows Leading dimension of a
124 : !> \param - obj%local_ncols local columns of matrix a
125 : !> \param - obj%nblk blocksize of cyclic distribution, must be the same in both directions!
126 : !> \param - obj%mpi_comm_rows MPI communicator for rows
127 : !> \param - obj%mpi_comm_cols MPI communicator for columns
128 : !> \param - obj%wantDebug logical, more debug information on failure
129 : !> \param a(lda,matrixCols) Distributed matrix which should be inverted
130 : !> Distribution is like in Scalapack.
131 : !> Only upper triangle needs to be set.
132 : !> The lower triangle is not referenced.
133 : !> \result succes logical, reports success or failure
134 288 : function elpa_invert_trm_real_double_impl(obj, a) result(success)
135 : #include "elpa_invert_trm.F90"
136 288 : end function elpa_invert_trm_real_double_impl
137 : #undef DOUBLE_PRECISION
138 : #undef REALCASE
139 :
140 : #if WANT_SINGLE_PRECISION_REAL
141 : #define REALCASE 1
142 : #define SINGLE_PRECISION
143 : #include "../general/precision_macros.h"
144 :
145 : !> \brief elpa_invert_trm_real_single_impl: Inverts a single-precision real upper triangular matrix
146 : !> \details
147 : !> \param obj elpa_t object contains:
148 : !> \param - obj%na Order of matrix
149 : !> \param - obj%local_nrows Leading dimension of a
150 : !> \param - obj%local_ncols local columns of matrix a
151 : !> \param - obj%nblk blocksize of cyclic distribution, must be the same in both directions!
152 : !> \param - obj%mpi_comm_rows MPI communicator for rows
153 : !> \param - obj%mpi_comm_cols MPI communicator for columns
154 : !> \param - obj%wantDebug logical, more debug information on failure
155 : !> \param a(lda,matrixCols) Distributed matrix which should be inverted
156 : !> Distribution is like in Scalapack.
157 : !> Only upper triangle needs to be set.
158 : !> The lower triangle is not referenced.
159 : !> \result succes logical, reports success or failure
160 :
161 144 : function elpa_invert_trm_real_single_impl(obj, a) result(success)
162 : #include "elpa_invert_trm.F90"
163 144 : end function elpa_invert_trm_real_single_impl
164 : #undef SINGLE_PRECISION
165 : #undef REALCASE
166 :
167 : #endif /* WANT_SINGLE_PRECISION_REAL */
168 :
169 :
170 : #define COMPLEXCASE 1
171 : #define DOUBLE_PRECISION
172 : #include "../general/precision_macros.h"
173 :
174 : !> \brief elpa_cholesky_complex_double_impl: Cholesky factorization of a double-precision complex hermitian matrix
175 : !> \details
176 : !> \param obj elpa_t object contains:
177 : !> \param - obj%na Order of matrix
178 : !> \param - obj%local_nrows Leading dimension of a
179 : !> \param - obj%local_ncols local columns of matrix a
180 : !> \param - obj%nblk blocksize of cyclic distribution, must be the same in both directions!
181 : !> \param - obj%mpi_comm_rows MPI communicator for rows
182 : !> \param - obj%mpi_comm_cols MPI communicator for columns
183 : !> \param - obj%wantDebug logical, more debug information on failure
184 : !> \param a(lda,matrixCols) Distributed matrix which should be inverted
185 : !> Distribution is like in Scalapack.
186 : !> Only upper triangle needs to be set.
187 : !> The lower triangle is not referenced.
188 : !> \result succes logical, reports success or failure
189 864 : function elpa_cholesky_complex_double_impl(obj, a) result(success)
190 :
191 : #include "elpa_cholesky_template.F90"
192 :
193 864 : end function elpa_cholesky_complex_double_impl
194 : #undef DOUBLE_PRECISION
195 : #undef COMPLEXCASE
196 :
197 : #ifdef WANT_SINGLE_PRECISION_COMPLEX
198 : #define COMPLEXCASE 1
199 : #define SINGLE_PRECISION
200 : #include "../general/precision_macros.h"
201 :
202 : !> \brief elpa_cholesky_complex_single_impl: Cholesky factorization of a single-precision complex hermitian matrix
203 : !> \details
204 : !> \param obj elpa_t object contains:
205 : !> \param - obj%na Order of matrix
206 : !> \param - obj%local_nrows Leading dimension of a
207 : !> \param - obj%local_ncols local columns of matrix a
208 : !> \param - obj%nblk blocksize of cyclic distribution, must be the same in both directions!
209 : !> \param - obj%mpi_comm_rows MPI communicator for rows
210 : !> \param - obj%mpi_comm_cols MPI communicator for columns
211 : !> \param - obj%wantDebug logical, more debug information on failure
212 : !> \param a(lda,matrixCols) Distributed matrix which should be inverted
213 : !> Distribution is like in Scalapack.
214 : !> Only upper triangle needs to be set.
215 : !> The lower triangle is not referenced.
216 : !> \result succes logical, reports success or failure
217 432 : function elpa_cholesky_complex_single_impl(obj, a) result(success)
218 :
219 : #include "elpa_cholesky_template.F90"
220 :
221 432 : end function elpa_cholesky_complex_single_impl
222 : #undef SINGLE_PRECISION
223 : #undef COMPLEXCASE
224 :
225 : #endif /* WANT_SINGLE_PRECISION_COMPLEX */
226 :
227 : #define COMPLEXCASE 1
228 : #define DOUBLE_PRECISION
229 : #include "../general/precision_macros.h"
230 :
231 : !> \brief elpa_invert_trm_complex_double_impl: Inverts a double-precision complex upper triangular matrix
232 : !> \details
233 : !> \param obj elpa_t object contains:
234 : !> \param - obj%na Order of matrix
235 : !> \param - obj%local_nrows Leading dimension of a
236 : !> \param - obj%local_ncols local columns of matrix a
237 : !> \param - obj%nblk blocksize of cyclic distribution, must be the same in both directions!
238 : !> \param - obj%mpi_comm_rows MPI communicator for rows
239 : !> \param - obj%mpi_comm_cols MPI communicator for columns
240 : !> \param - obj%wantDebug logical, more debug information on failure
241 : !> \param a(lda,matrixCols) Distributed matrix which should be inverted
242 : !> Distribution is like in Scalapack.
243 : !> Only upper triangle needs to be set.
244 : !> The lower triangle is not referenced.
245 : !> \result succes logical, reports success or failure
246 288 : function elpa_invert_trm_complex_double_impl(obj, a) result(success)
247 : #include "elpa_invert_trm.F90"
248 288 : end function elpa_invert_trm_complex_double_impl
249 : #undef DOUBLE_PRECISION
250 : #undef COMPLEXCASE
251 :
252 : #ifdef WANT_SINGLE_PRECISION_COMPLEX
253 : #define COMPLEXCASE 1
254 : #define SINGLE_PRECISION
255 : #include "../general/precision_macros.h"
256 :
257 : !> \brief elpa_invert_trm_complex_single_impl: Inverts a single-precision complex upper triangular matrix
258 : !> \details
259 : !> \param obj elpa_t object contains:
260 : !> \param - obj%na Order of matrix
261 : !> \param - obj%local_nrows Leading dimension of a
262 : !> \param - obj%local_ncols local columns of matrix a
263 : !> \param - obj%nblk blocksize of cyclic distribution, must be the same in both directions!
264 : !> \param - obj%mpi_comm_rows MPI communicator for rows
265 : !> \param - obj%mpi_comm_cols MPI communicator for columns
266 : !> \param - obj%wantDebug logical, more debug information on failure
267 : !> \param a(lda,matrixCols) Distributed matrix which should be inverted
268 : !> Distribution is like in Scalapack.
269 : !> Only upper triangle needs to be set.
270 : !> The lower triangle is not referenced.
271 : !> \result succes logical, reports success or failure
272 144 : function elpa_invert_trm_complex_single_impl(obj, a) result(success)
273 : #include "elpa_invert_trm.F90"
274 144 : end function elpa_invert_trm_complex_single_impl
275 : #undef SINGLE_PRECISION
276 : #undef COMPLEXCASE
277 :
278 : #endif /* WANT_SINGE_PRECISION_COMPLEX */
279 :
280 : #define REALCASE 1
281 : #define DOUBLE_PRECISION
282 : #include "../general/precision_macros.h"
283 864 : function elpa_mult_at_b_real_double_impl(obj, uplo_a, uplo_c, ncb, a, b, ldb, ldbCols, &
284 864 : c, ldc, ldcCols) result(success)
285 : #include "elpa_multiply_a_b.F90"
286 864 : end function elpa_mult_at_b_real_double_impl
287 : #undef DOUBLE_PRECISION
288 : #undef REALCASE
289 :
290 : #if WANT_SINGLE_PRECISION_REAL
291 : #define REALCASE 1
292 : #define SINGLE_PRECISION
293 : #include "../general/precision_macros.h"
294 :
295 : !> \brief elpa_mult_at_b_real_single_impl: Performs C : = A**T * B
296 : !> where A is a square matrix (obj%na,obj%na) which is optionally upper or lower triangular
297 : !> B is a (obj%na,ncb) matrix
298 : !> C is a (obj%na,ncb) matrix where optionally only the upper or lower
299 : !> triangle may be computed
300 : !> \details
301 :
302 : !> \param uplo_a 'U' if A is upper triangular
303 : !> 'L' if A is lower triangular
304 : !> anything else if A is a full matrix
305 : !> Please note: This pertains to the original A (as set in the calling program)
306 : !> whereas the transpose of A is used for calculations
307 : !> If uplo_a is 'U' or 'L', the other triangle is not used at all,
308 : !> i.e. it may contain arbitrary numbers
309 : !> \param uplo_c 'U' if only the upper diagonal part of C is needed
310 : !> 'L' if only the upper diagonal part of C is needed
311 : !> anything else if the full matrix C is needed
312 : !> Please note: Even when uplo_c is 'U' or 'L', the other triangle may be
313 : !> written to a certain extent, i.e. one shouldn't rely on the content there!
314 : !> \param na Number of rows/columns of A, number of rows of B and C
315 : !> \param ncb Number of columns of B and C
316 : !> \param a matrix a
317 : !> \param obj%local_nrows leading dimension of matrix a, set with class method obj%set("local_nrows",value)
318 : !> \param b matrix b
319 : !> \param ldb leading dimension of matrix b
320 : !> \param nblk blocksize of cyclic distribution, must be the same in both directions!
321 : !> \param mpi_comm_rows MPI communicator for rows
322 : !> \param mpi_comm_cols MPI communicator for columns
323 : !> \param c matrix c
324 : !> \param ldc leading dimension of matrix c
325 : !> \result success
326 :
327 288 : function elpa_mult_at_b_real_single_impl(obj, uplo_a, uplo_c, ncb, a, b, ldb, ldbCols, &
328 288 : c, ldc, ldcCols) result(success)
329 :
330 : #include "elpa_multiply_a_b.F90"
331 :
332 288 : end function elpa_mult_at_b_real_single_impl
333 : #undef SINGLE_PRECISION
334 : #undef REALCASE
335 : #endif /* WANT_SINGLE_PRECISION_REAL */
336 :
337 :
338 : #define COMPLEXCASE 1
339 : #define DOUBLE_PRECISION
340 : #include "../general/precision_macros.h"
341 :
342 : !> \brief elpa_mult_ah_b_complex_double_impl: Performs C : = A**H * B
343 : !> where A is a square matrix (obj%na,obj%na) which is optionally upper or lower triangular
344 : !> B is a (obj%na,ncb) matrix
345 : !> C is a (obj%na,ncb) matrix where optionally only the upper or lower
346 : !> triangle may be computed
347 : !> \details
348 : !>
349 : !> \param uplo_a 'U' if A is upper triangular
350 : !> 'L' if A is lower triangular
351 : !> anything else if A is a full matrix
352 : !> Please note: This pertains to the original A (as set in the calling program)
353 : !> whereas the transpose of A is used for calculations
354 : !> If uplo_a is 'U' or 'L', the other triangle is not used at all,
355 : !> i.e. it may contain arbitrary numbers
356 : !> \param uplo_c 'U' if only the upper diagonal part of C is needed
357 : !> 'L' if only the upper diagonal part of C is needed
358 : !> anything else if the full matrix C is needed
359 : !> Please note: Even when uplo_c is 'U' or 'L', the other triangle may be
360 : !> written to a certain extent, i.e. one shouldn't rely on the content there!
361 : !> \param na Number of rows/columns of A, number of rows of B and C
362 : !> \param ncb Number of columns of B and C
363 : !> \param a matrix a
364 : !> \param obj%local_ncols leading dimension of matrix a, set with class method obj%set("local_nrows",value)
365 : !> \param ldaCols columns of matrix a
366 : !> \param b matrix b
367 : !> \param ldb leading dimension of matrix b
368 : !> \param ldbCols columns of matrix b
369 : !> \param nblk blocksize of cyclic distribution, must be the same in both directions!
370 : !> \param mpi_comm_rows MPI communicator for rows
371 : !> \param mpi_comm_cols MPI communicator for columns
372 : !> \param c matrix c
373 : !> \param ldc leading dimension of matrix c
374 : !> \result success
375 :
376 576 : function elpa_mult_ah_b_complex_double_impl(obj, uplo_a, uplo_c, ncb, a, b, ldb, ldbCols, &
377 576 : c, ldc, ldcCols) result(success)
378 : #include "elpa_multiply_a_b.F90"
379 :
380 576 : end function elpa_mult_ah_b_complex_double_impl
381 : #undef DOUBLE_PRECISION
382 : #undef COMPLEXCASE
383 :
384 : #ifdef WANT_SINGLE_PRECISION_COMPLEX
385 : #define COMPLEXCASE 1
386 : #define SINGLE_PRECISION
387 : #include "../general/precision_macros.h"
388 :
389 : !> \brief elpa_mult_ah_b_complex_single_impl: Performs C : = A**H * B
390 : !> where A is a square matrix (obj%na,obj%na) which is optionally upper or lower triangular
391 : !> B is a (obj%na,ncb) matrix
392 : !> C is a (obj%na,ncb) matrix where optionally only the upper or lower
393 : !> triangle may be computed
394 : !> \details
395 : !>
396 : !> \param uplo_a 'U' if A is upper triangular
397 : !> 'L' if A is lower triangular
398 : !> anything else if A is a full matrix
399 : !> Please note: This pertains to the original A (as set in the calling program)
400 : !> whereas the transpose of A is used for calculations
401 : !> If uplo_a is 'U' or 'L', the other triangle is not used at all,
402 : !> i.e. it may contain arbitrary numbers
403 : !> \param uplo_c 'U' if only the upper diagonal part of C is needed
404 : !> 'L' if only the upper diagonal part of C is needed
405 : !> anything else if the full matrix C is needed
406 : !> Please note: Even when uplo_c is 'U' or 'L', the other triangle may be
407 : !> written to a certain extent, i.e. one shouldn't rely on the content there!
408 : !> \param na Number of rows/columns of A, number of rows of B and C
409 : !> \param ncb Number of columns of B and C
410 : !> \param a matrix a
411 : !> \param lda leading dimension of matrix a
412 : !> \param ldaCols columns of matrix a
413 : !> \param b matrix b
414 : !> \param ldb leading dimension of matrix b
415 : !> \param ldbCols columns of matrix b
416 : !> \param nblk blocksize of cyclic distribution, must be the same in both directions!
417 : !> \param mpi_comm_rows MPI communicator for rows
418 : !> \param mpi_comm_cols MPI communicator for columns
419 : !> \param c matrix c
420 : !> \param ldc leading dimension of matrix c
421 : !> \result success
422 :
423 288 : function elpa_mult_ah_b_complex_single_impl(obj, uplo_a, uplo_c, ncb, a, b, ldb, ldbCols, &
424 288 : c, ldc, ldcCols) result(success)
425 :
426 : #include "elpa_multiply_a_b.F90"
427 :
428 288 : end function elpa_mult_ah_b_complex_single_impl
429 : #undef SINGLE_PRECISION
430 : #undef COMPLEXCASE
431 : #endif /* WANT_SINGLE_PRECISION_COMPLEX */
432 :
433 : #define REALCASE 1
434 : #define DOUBLE_PRECISION
435 : #include "../general/precision_macros.h"
436 :
437 : !> \brief elpa_solve_tridi_double_impl: Solve tridiagonal eigensystem for a double-precision matrix with divide and conquer method
438 : !> \details
439 : !> \param obj elpa_t object contains:
440 : !> \param - obj%na Order of matrix
441 : !> \param - obj%nev number of eigenvalues/vectors to be computed
442 : !> \param - obj%local_nrows Leading dimension of q
443 : !> \param - obj%local_ncols local columns of matrix q
444 : !> \param - obj%nblk blocksize of cyclic distribution, must be the same in both directions!
445 : !> \param - obj%mpi_comm_rows MPI communicator for rows
446 : !> \param - obj%mpi_comm_cols MPI communicator for columns
447 : !> \param - obj%wantDebug logical, more debug information on failure
448 : !> \param d array d(na) on input diagonal elements of tridiagonal matrix, on
449 : !> output the eigenvalues in ascending order
450 : !> \param e array e(na) on input subdiagonal elements of matrix, on exit destroyed
451 : !> \param q on exit : matrix q(ldq,matrixCols) contains the eigenvectors
452 : !> \result succes logical, reports success or failure
453 576 : function elpa_solve_tridi_double_impl(obj, d, e, q) result(success)
454 :
455 : #include "elpa_solve_tridi_impl_public.F90"
456 :
457 576 : end function
458 : #undef DOUBLE_PRECISION
459 : #undef REALCASE
460 :
461 : #ifdef WANT_SINGLE_PRECISION_REAL
462 : #define REALCASE 1
463 : #define SINGLE_PRECISION
464 : #include "../general/precision_macros.h"
465 :
466 : !> \brief elpa_solve_tridi_single_impl: Solve tridiagonal eigensystem for a single-precision matrix with divide and conquer method
467 : !> \details
468 : !> \param obj elpa_t object contains:
469 : !> \param - obj%na Order of matrix
470 : !> \param - obj%nev number of eigenvalues/vectors to be computed
471 : !> \param - obj%local_nrows Leading dimension of q
472 : !> \param - obj%local_ncols local columns of matrix q
473 : !> \param - obj%nblk blocksize of cyclic distribution, must be the same in both directions!
474 : !> \param - obj%mpi_comm_rows MPI communicator for rows
475 : !> \param - obj%mpi_comm_cols MPI communicator for columns
476 : !> \param - obj%wantDebug logical, more debug information on failure
477 : !> \param d array d(na) on input diagonal elements of tridiagonal matrix, on
478 : !> output the eigenvalues in ascending order
479 : !> \param e array e(na) on input subdiagonal elements of matrix, on exit destroyed
480 : !> \param q on exit : matrix q(ldq,matrixCols) contains the eigenvectors
481 : !> \result succes logical, reports success or failure
482 144 : function elpa_solve_tridi_single_impl(obj, d, e, q) result(success)
483 :
484 : #include "elpa_solve_tridi_impl_public.F90"
485 :
486 144 : end function
487 : #undef SINGLE_PRECISION
488 : #undef REALCASE
489 : #endif /* WANT_SINGLE_PRECISION_REAL */
490 :
491 :
492 :
493 :
494 : end module elpa1_auxiliary_impl
495 :
|