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 :
19 : More information can be found here:
20 : http://elpa.mpcdf.mpg.de/
21 :
22 : ELPA is free software: you can redistribute it and/or modify
23 : it under the terms of the version 3 of the license of the
24 : GNU Lesser General Public License as published by the Free
25 : Software Foundation.
26 :
27 : ELPA is distributed in the hope that it will be useful,
28 : but WITHOUT ANY WARRANTY; without even the implied warranty of
29 : MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30 : GNU Lesser General Public License for more details.
31 :
32 : You should have received a copy of the GNU Lesser General Public License
33 : along with ELPA. If not, see <http://www.gnu.org/licenses/>
34 :
35 : ELPA reflects a substantial effort on the part of the original
36 : ELPA consortium, and we ask you to respect the spirit of the
37 : license that we chose: i.e., please contribute any changes you
38 : may have back to the original ELPA library distribution, and keep
39 : any derivatives of ELPA under the same license that we chose for
40 : the original distribution, the GNU Lesser General Public License.
41 : */
42 :
43 : #include "config.h"
44 :
45 : #include <stdio.h>
46 : #include <stdlib.h>
47 : #ifdef WITH_MPI
48 : #include <mpi.h>
49 : #endif
50 : #include <math.h>
51 :
52 : #include <elpa/elpa.h>
53 : #include <assert.h>
54 :
55 : #include "test/shared/generated.h"
56 :
57 : #if !(defined(TEST_REAL) ^ defined(TEST_COMPLEX))
58 : #error "define exactly one of TEST_REAL or TEST_COMPLEX"
59 : #endif
60 :
61 : #if !(defined(TEST_SINGLE) ^ defined(TEST_DOUBLE))
62 : #error "define exactly one of TEST_SINGLE or TEST_DOUBLE"
63 : #endif
64 :
65 : #if !(defined(TEST_SOLVER_1STAGE) ^ defined(TEST_SOLVER_2STAGE))
66 : #error "define exactly one of TEST_SOLVER_1STAGE or TEST_SOLVER_2STAGE"
67 : #endif
68 :
69 : #ifdef TEST_SINGLE
70 : # define EV_TYPE float
71 : # ifdef TEST_REAL
72 : # define MATRIX_TYPE float
73 : # else
74 : # define MATRIX_TYPE complex float
75 : # endif
76 : #else
77 : # define EV_TYPE double
78 : # ifdef TEST_REAL
79 : # define MATRIX_TYPE double
80 : # else
81 : # define MATRIX_TYPE complex double
82 : # endif
83 : #endif
84 :
85 : #define assert_elpa_ok(x) assert(x == ELPA_OK)
86 :
87 :
88 1008 : int main(int argc, char** argv) {
89 : /* matrix dimensions */
90 : int na, nev, nblk;
91 :
92 : /* mpi */
93 : int myid, nprocs;
94 : int na_cols, na_rows;
95 : int np_cols, np_rows;
96 : int my_prow, my_pcol;
97 : int mpi_comm;
98 : int provided_mpi_thread_level;
99 :
100 : /* blacs */
101 : int my_blacs_ctxt, sc_desc[9], info;
102 :
103 : /* The Matrix */
104 : MATRIX_TYPE *a, *as, *z;
105 : EV_TYPE *ev;
106 :
107 : int error, status;
108 :
109 : elpa_t handle;
110 :
111 : int value;
112 : #ifdef WITH_MPI
113 : #ifndef WITH_OPENMP
114 336 : MPI_Init(&argc, &argv);
115 : #else
116 336 : MPI_Init_thread(&argc, &argv, MPI_THREAD_MULTIPLE, &provided_mpi_thread_level);
117 :
118 336 : if (provided_mpi_thread_level != MPI_THREAD_MULTIPLE) {
119 0 : fprintf(stderr, "MPI ERROR: MPI_THREAD_MULTIPLE is not provided on this system\n");
120 0 : MPI_Finalize();
121 0 : exit(77);
122 : }
123 : #endif
124 :
125 672 : MPI_Comm_size(MPI_COMM_WORLD, &nprocs);
126 672 : MPI_Comm_rank(MPI_COMM_WORLD, &myid);
127 :
128 : #else
129 336 : nprocs = 1;
130 336 : myid = 0;
131 : #endif
132 :
133 1008 : if (argc == 4) {
134 1008 : na = atoi(argv[1]);
135 1008 : nev = atoi(argv[2]);
136 1008 : nblk = atoi(argv[3]);
137 : } else {
138 0 : na = 500;
139 0 : nev = 250;
140 0 : nblk = 16;
141 : }
142 :
143 1008 : for (np_cols = (int) sqrt((double) nprocs); np_cols > 1; np_cols--) {
144 0 : if (nprocs % np_cols == 0) {
145 0 : break;
146 : }
147 : }
148 :
149 1008 : np_rows = nprocs/np_cols;
150 :
151 : /* set up blacs */
152 : /* convert communicators before */
153 : #ifdef WITH_MPI
154 672 : mpi_comm = MPI_Comm_c2f(MPI_COMM_WORLD);
155 : #else
156 336 : mpi_comm = 0;
157 : #endif
158 1008 : set_up_blacsgrid_f(mpi_comm, np_rows, np_cols, 'C', &my_blacs_ctxt, &my_prow, &my_pcol);
159 1008 : set_up_blacs_descriptor_f(na, nblk, my_prow, my_pcol, np_rows, np_cols, &na_rows, &na_cols, sc_desc, my_blacs_ctxt, &info);
160 :
161 : /* allocate the matrices needed for elpa */
162 1008 : a = calloc(na_rows*na_cols, sizeof(MATRIX_TYPE));
163 1008 : z = calloc(na_rows*na_cols, sizeof(MATRIX_TYPE));
164 1008 : as = calloc(na_rows*na_cols, sizeof(MATRIX_TYPE));
165 1008 : ev = calloc(na, sizeof(EV_TYPE));
166 :
167 : #ifdef TEST_REAL
168 : #ifdef TEST_DOUBLE
169 384 : prepare_matrix_random_real_double_f(na, myid, na_rows, na_cols, sc_desc, a, z, as);
170 : #else
171 192 : prepare_matrix_random_real_single_f(na, myid, na_rows, na_cols, sc_desc, a, z, as);
172 : #endif
173 : #else
174 : #ifdef TEST_DOUBLE
175 288 : prepare_matrix_random_complex_double_f(na, myid, na_rows, na_cols, sc_desc, a, z, as);
176 : #else
177 144 : prepare_matrix_random_complex_single_f(na, myid, na_rows, na_cols, sc_desc, a, z, as);
178 : #endif
179 : #endif
180 :
181 1008 : if (elpa_init(CURRENT_API_VERSION) != ELPA_OK) {
182 0 : fprintf(stderr, "Error: ELPA API version not supported");
183 0 : exit(1);
184 : }
185 :
186 1008 : handle = elpa_allocate(&error);
187 1008 : assert_elpa_ok(error);
188 :
189 : /* Set parameters */
190 1008 : elpa_set(handle, "na", na, &error);
191 1008 : assert_elpa_ok(error);
192 :
193 1008 : elpa_set(handle, "nev", nev, &error);
194 1008 : assert_elpa_ok(error);
195 :
196 1008 : if (myid == 0) {
197 672 : printf("Setting the matrix parameters na=%d, nev=%d \n",na,nev);
198 : }
199 1008 : elpa_set(handle, "local_nrows", na_rows, &error);
200 1008 : assert_elpa_ok(error);
201 :
202 1008 : elpa_set(handle, "local_ncols", na_cols, &error);
203 1008 : assert_elpa_ok(error);
204 :
205 1008 : elpa_set(handle, "nblk", nblk, &error);
206 1008 : assert_elpa_ok(error);
207 :
208 : #ifdef WITH_MPI
209 672 : elpa_set(handle, "mpi_comm_parent", MPI_Comm_c2f(MPI_COMM_WORLD), &error);
210 672 : assert_elpa_ok(error);
211 :
212 672 : elpa_set(handle, "process_row", my_prow, &error);
213 672 : assert_elpa_ok(error);
214 :
215 672 : elpa_set(handle, "process_col", my_pcol, &error);
216 672 : assert_elpa_ok(error);
217 : #endif
218 :
219 : /* Setup */
220 1008 : assert_elpa_ok(elpa_setup(handle));
221 :
222 : /* Set tunables */
223 : #ifdef TEST_SOLVER_1STAGE
224 576 : elpa_set(handle, "solver", ELPA_SOLVER_1STAGE, &error);
225 : #else
226 432 : elpa_set(handle, "solver", ELPA_SOLVER_2STAGE, &error);
227 : #endif
228 1008 : assert_elpa_ok(error);
229 :
230 1008 : elpa_set(handle, "gpu", TEST_GPU, &error);
231 1008 : assert_elpa_ok(error);
232 :
233 : #if defined(TEST_SOLVE_2STAGE) && defined(TEST_KERNEL)
234 : # ifdef TEST_COMPLEX
235 : elpa_set(handle, "complex_kernel", TEST_KERNEL, &error);
236 : # else
237 : elpa_set(handle, "real_kernel", TEST_KERNEL, &error);
238 : # endif
239 : assert_elpa_ok(error);
240 : #endif
241 :
242 1008 : elpa_get(handle, "solver", &value, &error);
243 1008 : if (myid == 0) {
244 672 : printf("Solver is set to %d \n", value);
245 : }
246 : /* Solve EV problem */
247 1008 : elpa_eigenvectors(handle, a, ev, z, &error);
248 1008 : assert_elpa_ok(error);
249 :
250 1008 : elpa_deallocate(handle);
251 1008 : elpa_uninit();
252 :
253 :
254 : /* check the results */
255 : #ifdef TEST_REAL
256 : #ifdef TEST_DOUBLE
257 384 : status = check_correctness_evp_numeric_residuals_real_double_f(na, nev, na_rows, na_cols, as, z, ev, sc_desc, myid);
258 : #else
259 192 : status = check_correctness_evp_numeric_residuals_real_single_f(na, nev, na_rows, na_cols, as, z, ev, sc_desc, myid);
260 : #endif
261 : #else
262 : #ifdef TEST_DOUBLE
263 288 : status = check_correctness_evp_numeric_residuals_complex_double_f(na, nev, na_rows, na_cols, as, z, ev, sc_desc, myid);
264 : #else
265 144 : status = check_correctness_evp_numeric_residuals_complex_single_f(na, nev, na_rows, na_cols, as, z, ev, sc_desc, myid);
266 : #endif
267 : #endif
268 :
269 1008 : if (status !=0){
270 0 : printf("The computed EVs are not correct !\n");
271 : }
272 1008 : if (status ==0){
273 1008 : printf("All ok!\n");
274 : }
275 :
276 1008 : free(a);
277 1008 : free(z);
278 1008 : free(as);
279 1008 : free(ev);
280 :
281 : #ifdef WITH_MPI
282 672 : MPI_Finalize();
283 : #endif
284 :
285 1008 : return !!status;
286 : }
|
