legacy_real_2stage_c_version.c 7.38 KB
Newer Older
Andreas Marek's avatar
Andreas Marek committed
1
2
3
4
5
/*     This file is part of ELPA. */
/*  */
/*     The ELPA library was originally created by the ELPA consortium, */
/*     consisting of the following organizations: */
/*  */
6
7
/*     - Max Planck Computing and Data Facility (MPCDF), formerly known as */
/*       Rechenzentrum Garching der Max-Planck-Gesellschaft (RZG), */
Andreas Marek's avatar
Andreas Marek committed
8
9
10
11
12
13
14
15
16
17
18
19
/*     - Bergische Universität Wuppertal, Lehrstuhl für angewandte */
/*       Informatik, */
/*     - Technische Universität München, Lehrstuhl für Informatik mit */
/*       Schwerpunkt Wissenschaftliches Rechnen , */
/*     - Fritz-Haber-Institut, Berlin, Abt. Theorie, */
/*     - Max-Plack-Institut für Mathematik in den Naturwissenschaften, */
/*       Leipzig, Abt. Komplexe Strukutren in Biologie und Kognition, */
/*       and */
/*     - IBM Deutschland GmbH */
/*  */
/*  */
/*     More information can be found here: */
20
/*     http://elpa.mpcdf.mpg.de/ */
Andreas Marek's avatar
Andreas Marek committed
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
/*  */
/*     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. */
/*  */
/*  */

#include "config-f90.h"

#include <stdio.h>
#include <stdlib.h>
48
#ifdef WITH_MPI
Andreas Marek's avatar
Andreas Marek committed
49
#include <mpi.h>
50
#endif
Andreas Marek's avatar
Andreas Marek committed
51
52
#include <math.h>

53
#include <elpa/elpa_legacy.h>
54
#include <test/shared/generated.h>
Andreas Marek's avatar
Andreas Marek committed
55

56
57
#define DOUBLE_PRECISION_REAL 1

58
int main(int argc, char** argv) {
Andreas Marek's avatar
Andreas Marek committed
59
60
   int myid;
   int nprocs;
61
62
63
#ifndef WITH_MPI
   int MPI_COMM_WORLD;
#endif
Andreas Marek's avatar
Andreas Marek committed
64
65
66
67
68
69
   int na, nev, nblk;

   int status;

   int np_cols, np_rows, np_colsStart;

70
   int my_blacs_ctxt, my_prow, my_pcol;
Andreas Marek's avatar
Andreas Marek committed
71
72
73
74
75
76
77
78
79
80

   int mpierr;

   int my_mpi_comm_world;
   int mpi_comm_rows, mpi_comm_cols;

   int info, *sc_desc;

   int na_rows, na_cols;
   double startVal;
81
#ifdef DOUBLE_PRECISION_REAL
82
   double *a, *z, *as, *ev;
83
#else
84
   float *a, *z, *as, *ev;
85
#endif
Andreas Marek's avatar
Andreas Marek committed
86
87
88

   int success;

89
   int useQr, THIS_REAL_ELPA_KERNEL_API, useGPU;
90

91
#ifdef WITH_MPI
Andreas Marek's avatar
Andreas Marek committed
92
93
94
   MPI_Init(&argc, &argv);
   MPI_Comm_size(MPI_COMM_WORLD, &nprocs);
   MPI_Comm_rank(MPI_COMM_WORLD, &myid);
95
96
97
98
99
#else
   nprocs = 1;
   myid=0;
   MPI_COMM_WORLD=1;
#endif
Andreas Marek's avatar
Andreas Marek committed
100
101
102
103
104
105
106
   na = 1000;
   nev = 500;
   nblk = 16;

   if (myid == 0) {
     printf("This is the c version of an ELPA test-programm\n");
     printf("\n");
107
108
     printf("It will call the 2stage ELPA real solver for an\n");
     printf("matrix of size %d. It will compute %d eigenvalues\n",na,nev);
Andreas Marek's avatar
Andreas Marek committed
109
110
111
112
     printf("and uses a blocksize of %d\n",nblk);
     printf("\n");
     printf("This is an example program with much less functionality\n");
     printf("as it's Fortran counterpart. It's only purpose is to show how \n");
113
     printf("to evoke ELPA2 from a c programm\n");
Andreas Marek's avatar
Andreas Marek committed
114
     printf("\n");
115
116
117
118
119
#ifdef DOUBLE_PRECISION_REAL
    printf(" Double precision version of ELPA2 is used. \n");
#else
    printf(" Single precision version of ELPA2 is used. \n");
#endif
Andreas Marek's avatar
Andreas Marek committed
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
   }

   status = 0;

   startVal = sqrt((double) nprocs);
   np_colsStart = (int) round(startVal);
   for (np_cols=np_colsStart;np_cols>1;np_cols--){
     if (nprocs %np_cols ==0){
     break;
     }
   }

   np_rows = nprocs/np_cols;

   if (myid == 0) {
     printf("\n");
     printf("Number of processor rows %d, cols %d, total %d \n",np_rows,np_cols,nprocs);
   }

   /* set up blacs */
   /* convert communicators before */
141
#ifdef WITH_MPI
Andreas Marek's avatar
Andreas Marek committed
142
   my_mpi_comm_world = MPI_Comm_c2f(MPI_COMM_WORLD);
143
144
145
#else
  my_mpi_comm_world = 1;
#endif
146
   set_up_blacsgrid_f(my_mpi_comm_world, np_rows, np_cols, 'C', &my_blacs_ctxt, &my_prow, &my_pcol);
Andreas Marek's avatar
Andreas Marek committed
147
148
149
150
151
152
153
154
155

   if (myid == 0) {
     printf("\n");
     printf("Past BLACS_Gridinfo...\n");
     printf("\n");
   }

   /* get the ELPA row and col communicators. */
   /* These are NOT usable in C without calling the MPI_Comm_f2c function on them !! */
156
#ifdef WITH_MPI
Andreas Marek's avatar
Andreas Marek committed
157
   my_mpi_comm_world = MPI_Comm_c2f(MPI_COMM_WORLD);
158
#endif
159
   mpierr = elpa_get_communicators(my_mpi_comm_world, my_prow, my_pcol, &mpi_comm_rows, &mpi_comm_cols);
Andreas Marek's avatar
Andreas Marek committed
160
161
162
163
164
165
166
167
168

   if (myid == 0) {
     printf("\n");
     printf("Past split communicator setup for rows and columns...\n");
     printf("\n");
   }

   sc_desc = malloc(9*sizeof(int));

169
   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);
Andreas Marek's avatar
Andreas Marek committed
170
171
172
173
174
175
176
177
178
179
180
181
182

   if (myid == 0) {
     printf("\n");
     printf("Past scalapack descriptor setup...\n");
     printf("\n");
   }

   /* allocate the matrices needed for elpa */
   if (myid == 0) {
     printf("\n");
     printf("Allocating matrices with na_rows=%d and na_cols=%d\n",na_rows, na_cols);
     printf("\n");
   }
183
#ifdef DOUBLE_PRECISION_REAL
Andreas Marek's avatar
Andreas Marek committed
184
185
186
187
   a  = malloc(na_rows*na_cols*sizeof(double));
   z  = malloc(na_rows*na_cols*sizeof(double));
   as = malloc(na_rows*na_cols*sizeof(double));
   ev = malloc(na*sizeof(double));
188
189
190
191
192
193
194
#else
   a  = malloc(na_rows*na_cols*sizeof(float));
   z  = malloc(na_rows*na_cols*sizeof(float));
   as = malloc(na_rows*na_cols*sizeof(float));
   ev = malloc(na*sizeof(float));
#endif
#ifdef DOUBLE_PRECISION_REAL
195
   prepare_matrix_real_double_f(na, myid, na_rows, na_cols, sc_desc, a, z, as);
196
#else
197
   prepare_matrix_real_single_f(na, myid, na_rows, na_cols, sc_desc, a, z, as);
198
#endif
Andreas Marek's avatar
Andreas Marek committed
199
200
201
202
203
   if (myid == 0) {
     printf("\n");
     printf("Entering ELPA 2stage real solver\n");
     printf("\n");
   }
204
#ifdef WITH_MPI
Andreas Marek's avatar
Andreas Marek committed
205
   mpierr = MPI_Barrier(MPI_COMM_WORLD);
206
#endif
207
   useGPU =0 ;
Andreas Marek's avatar
Andreas Marek committed
208
   useQr = 0;
209
   THIS_REAL_ELPA_KERNEL_API = ELPA_2STAGE_REAL_DEFAULT;
210
#ifdef DOUBLE_PRECISION_REAL
211
   success = elpa_solve_evp_real_2stage_double_precision(na, nev, a, na_rows, ev, z, na_rows, nblk, na_cols, mpi_comm_rows, mpi_comm_cols, my_mpi_comm_world, THIS_REAL_ELPA_KERNEL_API, useQr, useGPU);
212
#else
213
   success = elpa_solve_evp_real_2stage_single_precision(na, nev, a, na_rows, ev, z, na_rows, nblk, na_cols, mpi_comm_rows, mpi_comm_cols, my_mpi_comm_world, THIS_REAL_ELPA_KERNEL_API, useQr, useGPU);
214
#endif
Andreas Marek's avatar
Andreas Marek committed
215
216
   if (success != 1) {
     printf("error in ELPA solve \n");
217
#ifdef WITH_MPI
Andreas Marek's avatar
Andreas Marek committed
218
     mpierr = MPI_Abort(MPI_COMM_WORLD, 99);
219
220
#else
     exit(99);
221
#endif
Andreas Marek's avatar
Andreas Marek committed
222
223
224
225
226
227
228
229
230
231
   }


   if (myid == 0) {
     printf("\n");
     printf("2stage ELPA real solver complete\n");
     printf("\n");
   }

   /* check the results */
232
#ifdef DOUBLE_PRECISION_REAL
233
   status = check_correctness_real_double_f(na, nev, na_rows, na_cols, as, z, ev, sc_desc, myid);
234
#else
235
   status = check_correctness_real_single_f(na, nev, na_rows, na_cols, as, z, ev, sc_desc, myid);
236
#endif
Andreas Marek's avatar
Andreas Marek committed
237
238
239
240
241
242
243
244
245
246
247
248
249
250

   if (status !=0){
     printf("The computed EVs are not correct !\n");
   }
   if (status ==0){
     if (myid ==0) {
       printf("All ok!\n");
     }
   }

   free(sc_desc);
   free(a);
   free(z);
   free(as);
251
   free(ev);
252

253
#ifdef WITH_MPI
Andreas Marek's avatar
Andreas Marek committed
254
   MPI_Finalize();
255
#endif
Andreas Marek's avatar
Andreas Marek committed
256
257
   return 0;
}