/[escript]/trunk/escript/src/Utils.cpp
ViewVC logotype

Contents of /trunk/escript/src/Utils.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2642 - (show annotations)
Tue Sep 1 04:15:50 2009 UTC (10 years, 7 months ago) by jfenwick
File size: 13945 byte(s)
pre cleanup
1
2 /*******************************************************
3 *
4 * Copyright (c) 2003-2009 by University of Queensland
5 * Earth Systems Science Computational Center (ESSCC)
6 * http://www.uq.edu.au/esscc
7 *
8 * Primary Business: Queensland, Australia
9 * Licensed under the Open Software License version 3.0
10 * http://www.opensource.org/licenses/osl-3.0.php
11 *
12 *******************************************************/
13
14
15 #include <fstream>
16 #include <string.h>
17
18 // added for saveCSV
19 #include <boost/python.hpp>
20 #include <boost/scoped_ptr.hpp>
21 #include "Data.h"
22
23 #include "Utils.h"
24 #include "DataVector.h"
25
26 #ifdef _OPENMP
27 #include <omp.h>
28 #endif
29
30 #ifdef PASO_MPI
31 #include <mpi.h>
32 #endif
33
34 #ifdef _WIN32
35 #include <WinSock2.h>
36 #else
37 #include <unistd.h>
38 #endif
39
40 namespace escript {
41
42 int getSvnVersion()
43 {
44 #ifdef SVN_VERSION
45 return SVN_VERSION;
46 #else
47 return 0;
48 #endif
49 }
50
51 /* This is probably not very robust, but it works on Savanna today and is useful for performance analysis */
52 int get_core_id() {
53 int processor_num=-1;
54 #ifdef CORE_ID1
55 FILE *fp;
56 int i, count_spaces=0;
57 char fname[100];
58 char buf[1000];
59
60 sprintf(fname, "/proc/%d/stat", getpid());
61 fp = fopen(fname, "r");
62 if (fp == NULL) return(-1);
63 fgets(buf, 1000, fp);
64 fclose(fp);
65
66 for (i=strlen(buf)-1; i>=0; i--) {
67 if (buf[i] == ' ') count_spaces++;
68 if (count_spaces == 4) break;
69 }
70 processor_num = atoi(&buf[i+1]);
71 #endif
72 return(processor_num);
73 }
74
75
76 void printParallelThreadCnt()
77 {
78 int mpi_iam=0, mpi_num=1;
79 char hname[64];
80
81 #ifdef HAVE_GETHOSTNAME
82 gethostname(hname, 64);
83 hname[63] = '\0';
84 #else
85 strcpy(hname, "unknown host");
86 #endif
87
88 #ifdef PASO_MPI
89 MPI_Comm_rank(MPI_COMM_WORLD, &mpi_iam);
90 MPI_Comm_size(MPI_COMM_WORLD, &mpi_num);
91 #endif
92
93 #pragma omp parallel
94 {
95 int omp_iam=0, omp_num=1;
96 #ifdef _OPENMP
97 omp_iam = omp_get_thread_num(); /* Call in a parallel region */
98 omp_num = omp_get_num_threads();
99 #endif
100 #pragma omp critical (printthrdcount)
101 printf("printParallelThreadCounts: MPI=%03d/%03d OpenMP=%03d/%03d running on %s core %d\n",
102 mpi_iam, mpi_num, omp_iam, omp_num, hname, get_core_id());
103 }
104 }
105
106 void setNumberOfThreads(const int num_threads)
107 {
108
109 #ifdef _OPENMP
110 omp_set_num_threads(num_threads);
111 #endif
112
113 }
114
115 int getNumberOfThreads()
116 {
117 #ifdef _OPENMP
118 return omp_get_max_threads();
119 #else
120 return 1;
121 #endif
122
123 }
124
125 ESCRIPT_DLL_API int getMPISizeWorld() {
126 int mpi_num = 1;
127 #ifdef PASO_MPI
128 MPI_Comm_size(MPI_COMM_WORLD, &mpi_num);
129 #endif
130 return mpi_num;
131 }
132
133 ESCRIPT_DLL_API int getMPIRankWorld() {
134 int mpi_iam = 0;
135 #ifdef PASO_MPI
136 MPI_Comm_rank(MPI_COMM_WORLD, &mpi_iam);
137 #endif
138 return mpi_iam;
139 }
140
141 ESCRIPT_DLL_API int getMPIWorldMax(const int val) {
142 #ifdef PASO_MPI
143 int val2 = val;
144 int out = val;
145 MPI_Allreduce( &val2, &out, 1, MPI_INT, MPI_MAX, MPI_COMM_WORLD );
146 #else
147 int out = val;
148 #endif
149 return out;
150 }
151
152 ESCRIPT_DLL_API int getMPIWorldSum(const int val) {
153 #ifdef PASO_MPI
154 int val2 = val;
155 int out = 0;
156 MPI_Allreduce( &val2, &out, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD );
157 #else
158 int out = val;
159 #endif
160 return out;
161 }
162
163 ESCRIPT_DLL_API double getMachinePrecision() {
164 return DBL_EPSILON;
165 }
166 ESCRIPT_DLL_API double getMaxFloat() {
167 return DBL_MAX;
168 }
169 ESCRIPT_DLL_API void MPIBarrierWorld() {
170 #ifdef PASO_MPI
171 MPI_Barrier(MPI_COMM_WORLD );
172 #endif
173 }
174
175
176 /*
177 ESCRIPT_DLL_API
178 void
179 saveDataCSV(const std::string& filename, boost::python::dict arg, const std::string& sep, const std::string& csep,
180 bool append)
181 {
182 boost::python::list keys=arg.keys();
183 int numdata = boost::python::extract<int>(arg.attr("__len__")());
184 bool hasmask=arg.has_key("mask");
185 Data mask;
186 if (hasmask)
187 {
188 mask=boost::python::extract<escript::Data>(arg["mask"]);
189 keys.remove("mask");
190 if (mask.getDataPointRank()!=0)
191 {
192 throw DataException("saveDataCSVcpp: masks must be scalar.");
193 }
194 }
195 if (numdata<1)
196 {
197 throw DataException("saveDataCSVcpp: no data to save specified.");
198 }
199 std::vector<std::string> names(numdata);
200 std::vector<Data> data(numdata);
201 std::vector<int> fstypes(numdata); // FunctionSpace types for each data
202
203 if (hasmask)
204 {
205 numdata--;
206 }
207
208 // We need to interpret the samples correctly even if they are different types
209 // for this reason, we should interate over samples
210 for (int i=0;i<numdata;++i)
211 {
212 names[i]=boost::python::extract<std::string>(keys[i]);
213 data[i]=boost::python::extract<escript::Data>(arg[keys[i]]);
214 fstypes[i]=data[i].getFunctionSpace().getTypeCode();
215 if (i>0)
216 {
217 if (data[i].getDomain()!=data[i-1].getDomain())
218 {
219 throw DataException("saveDataCSVcpp: all data must be on the same domain.");
220 }
221 }
222 }
223 const_Domain_ptr dom0=data[0].getDomain();
224 if (hasmask)
225 {
226 if (dom0!=mask.getDomain())
227 {
228 throw DataException("saveDataCSVcpp: mask be on the same FunctionSpace as data.");
229 }
230 fstypes[numdata]=mask.getFunctionSpace().getTypeCode();
231 names[numdata]="mask";
232 data[numdata]=mask;
233 }
234 int bestfnspace=0;
235 if (!dom0->commonFunctionSpace(fstypes, bestfnspace))
236 {
237 throw DataException("saveDataCSVcpp: FunctionSpaces of data are incompatible");
238 }
239 // now we interpolate all data to the same type
240 FunctionSpace best(dom0,bestfnspace);
241 for (int i=0;i<data.size();++i)
242 {
243 data[i]=data[i].interpolate(best);
244 }
245 dom0->saveDataCSV(filename, data, names, sep, csep, append, hasmask);
246 }
247 */
248
249
250 ESCRIPT_DLL_API
251 void
252 saveDataCSV(const std::string& filename, boost::python::dict arg, const std::string& sep, const std::string& csep,
253 bool append)
254 {
255 using std::cout;
256 using std::endl;
257 boost::python::list keys=arg.keys();
258 int numdata = boost::python::extract<int>(arg.attr("__len__")());
259 bool hasmask=arg.has_key("mask");
260 Data mask;
261 if (hasmask)
262 {
263 mask=boost::python::extract<escript::Data>(arg["mask"]);
264 keys.remove("mask");
265 numdata--;
266 if (mask.getDataPointRank()!=0)
267 {
268 throw DataException("saveDataCSVcpp: masks must be scalar.");
269 }
270 }
271 if (numdata<1)
272 {
273 throw DataException("saveDataCSVcpp: no data to save specified.");
274 }
275 std::vector<int> step(numdata);
276 std::vector<std::string> names(numdata);
277 std::vector<Data> data(numdata);
278 std::vector<const DataAbstract::ValueType::value_type*> samples(numdata);
279 std::vector<int> offset(numdata);
280 std::vector<int> fstypes(numdata); // FunctionSpace types for each data
281
282 // We need to interpret the samples correctly even if they are different types
283 // for this reason, we should interate over samples
284 for (int i=0;i<numdata;++i)
285 {
286 names[i]=boost::python::extract<std::string>(keys[i]);
287 data[i]=boost::python::extract<escript::Data>(arg[keys[i]]);
288 step[i]=(data[i].actsExpanded()?DataTypes::noValues(data[i].getDataPointShape()):0);
289 fstypes[i]=data[i].getFunctionSpace().getTypeCode();
290 if (i>0)
291 {
292 if (data[i].getDomain()!=data[i-1].getDomain())
293 {
294 throw DataException("saveDataCSVcpp: all data must be on the same domain.");
295 }
296 }
297 }
298 int bestfnspace=0;
299 if (!data[0].getDomain()->commonFunctionSpace(fstypes, bestfnspace))
300 {
301 throw DataException("saveDataCSVcpp: FunctionSpaces of data are incompatible");
302 }
303 // now we interpolate all data to the same type
304 FunctionSpace best(data[0].getDomain(),bestfnspace);
305 for (int i=0;i<numdata;++i)
306 {
307 data[i]=data[i].interpolate(best);
308 }
309 int numsamples=data[0].getNumSamples(); // these must be the same for all data
310 int dpps=data[0].getNumDataPointsPerSample();
311
312
313 std::ostringstream os;
314
315
316 bool first=true;
317
318 if (data[0].getDomain()->getMPIRank()==0)
319 {
320 for (int i=0;i<numdata;++i)
321 {
322 const DataTypes::ShapeType& s=data[i].getDataPointShape();
323 switch (data[i].getDataPointRank())
324 {
325 case 0: if (!first)
326 {
327 os << sep;
328 }
329 else
330 {
331 first=false;
332 }
333 os << names[i]; break;
334 case 1: for (int j=0;j<s[0];++j)
335 {
336 if (!first)
337 {
338 os << sep;
339 }
340 else
341 {
342 first=false;
343 }
344 os << names[i] << csep << j;
345 }
346 break;
347 case 2: for (int j=0;j<s[0];++j)
348 {
349 for (int k=0;k<s[1];++k)
350 {
351 if (!first)
352 {
353 os << sep;
354 }
355 else
356 {
357 first=false;
358 }
359 os << names[i] << csep << k << csep << j;
360 }
361 }
362 break;
363 case 3: for (int j=0;j<s[0];++j)
364 {
365 for (int k=0;k<s[1];++k)
366 {
367 for (int l=0;l<s[2];++l)
368 {
369 if (!first)
370 {
371 os << sep;
372 }
373 else
374 {
375 first=false;
376 }
377 os << names[i] << csep << k << csep << j << csep << l;
378 }
379 }
380 }
381 break;
382 case 4: for (int j=0;j<s[0];++j)
383 {
384 for (int k=0;k<s[1];++k)
385 {
386 for (int l=0;l<s[2];++l)
387 {
388 for (int m=0;m<s[3];++m)
389 {
390 if (!first)
391 {
392 os << sep;
393 }
394 else
395 {
396 first=false;
397 }
398 os << names[i] << csep << k << csep << j << csep << l << csep << m;
399 }
400 }
401 }
402 }
403 break;
404 default:
405 throw DataException("saveDataCSV: Illegal rank");
406 }
407 }
408 os << endl;
409 }
410 boost::scoped_ptr<BufferGroup> maskbuffer; // sample buffer for the mask [if we have one]
411 const double* masksample=0;
412 int maskoffset=0;
413 //the use of shared_ptr here is just to ensure the buffer group is freed
414 //I would have used scoped_ptr but they don't work in vectors
415 std::vector<boost::shared_ptr<BufferGroup> > bg(numdata);
416 for (int d=0;d<numdata;++d)
417 {
418 bg[d].reset(data[d].allocSampleBuffer());
419 }
420
421 bool expandedmask=false; // does the mask act expanded. Are there mask value for each point in the sample
422 bool wantrow=true; // do we output this row?
423 if (hasmask)
424 {
425 maskbuffer.reset(mask.allocSampleBuffer());
426 if (mask.actsExpanded())
427 {
428 maskoffset=DataTypes::noValues(mask.getDataPointShape());
429 expandedmask=true;
430 }
431 }
432 os.setf(std::ios_base::scientific, std::ios_base::floatfield);
433 os.precision(15);
434
435 // errors prior to this point will occur on all processes anyway
436 // so there is no need to explicitly notify other ranks
437 int error=0;
438 try{
439 for (int i=0;i<numsamples;++i)
440 {
441 if (!best.ownSample(i))
442 {
443 continue;
444 }
445 wantrow=true;
446 for (int d=0;d<numdata;++d)
447 {
448 samples[d]=data[d].getSampleDataRO(i,bg[d].get());
449 }
450 if (hasmask)
451 {
452 masksample=mask.getSampleDataRO(i, maskbuffer.get());
453 if (!expandedmask) // mask controls whole sample
454 {
455 if (masksample[0]<=0) // masks are scalar
456 {
457 wantrow=false;
458 }
459 }
460 }
461 for (int j=0;j<dpps;++j)
462 {
463 // now we need to check if this point is masked off
464 if (expandedmask)
465 {
466 wantrow=(masksample[j]>0); // masks are scalar to the relevant value is at [j]
467 }
468 if (wantrow)
469 {
470 bool needsep=false;
471 for (int d=0;d<numdata;++d)
472 {
473 DataTypes::pointToStream(os, samples[d], data[d].getDataPointShape(), offset[d], needsep, sep);
474 needsep=true;
475 offset[d]+=step[d];
476 }
477 os << endl;
478 }
479 }
480 for (int d=0;d<numdata;++d)
481 {
482 offset[d]=0;
483 }
484 }
485 } catch (...)
486 {
487 error=1;
488 #ifndef PASO_MPI
489 throw;
490 #endif
491 }
492 #ifdef PASO_MPI
493 MPI_Comm com=data[0].getDomain()->getMPIComm();
494 int rerror=0;
495 MPI_Allreduce( &error, &rerror, 1, MPI_INT, MPI_MAX, com );
496 error=rerror;
497 if (error)
498 {
499 throw DataException("saveDataCSVcpp: error building output");
500 }
501 #endif
502
503 // at this point os will contain the text to be written
504 #ifndef PASO_MPI
505
506 std::ofstream ofs;
507 if (append)
508 {
509 ofs.open(filename.c_str(), std::ios_base::app);
510 }
511 else
512 {
513 ofs.open(filename.c_str());
514 }
515 if (!ofs.is_open())
516 {
517 throw DataException("saveDataCSVcpp: unable to open file for writing");
518 }
519 ofs << os.str();
520 ofs.close();
521
522 #else
523 // here we have MPI
524 const char* mpistr=0;
525 MPI_File mpi_fileHandle_p;
526 MPI_Status mpi_status;
527 MPI_Info mpi_info = MPI_INFO_NULL;
528 char* fname_c=new char[filename.size()+1];
529 strcpy(fname_c,filename.c_str());
530 boost::scoped_ptr<char> fname_p(fname_c);
531
532 int amode = MPI_MODE_CREATE|MPI_MODE_WRONLY|MPI_MODE_UNIQUE_OPEN;
533 if (append)
534 {
535 amode |= MPI_MODE_APPEND;
536 }
537 else
538 {
539 if (data[0].getDomain()->getMPIRank()==0)
540 {
541 std::ifstream ifs(fname_p.get()); // if file exists, remove it
542 if (ifs.is_open())
543 {
544 ifs.close();
545 if (remove(fname_p.get()))
546 {
547 error=1;
548 }
549 }
550 }
551 data[0].getDomain()->MPIBarrier();
552 int rerror=0;
553 MPI_Allreduce( &error, &rerror, 1, MPI_INT, MPI_MAX, com );
554 if (rerror!=0)
555 {
556 std::ostringstream oss;
557 oss << "saveDataCSVcpp: File " << filename << " already exists and could not be removed in preparation for new output.";
558 throw DataException(oss.str());
559 }
560 }
561 int ierr;
562 ierr = MPI_File_open(com, fname_p.get(), amode, mpi_info, &mpi_fileHandle_p);
563 if (ierr != MPI_SUCCESS)
564 {
565 std::ostringstream oss;
566 oss << "saveDataCSVcpp: File " << filename << " could not be opened for writing in parallel";
567 // file is not open so we can throw
568 throw DataException(oss.str());
569 }
570 else
571 {
572 ierr=MPI_File_set_view(mpi_fileHandle_p,MPI_DISPLACEMENT_CURRENT,
573 MPI_CHAR, MPI_CHAR, "native", mpi_info);
574 // here we are assuming that std::string holds the same type of char as MPI_CHAR
575 }
576
577 std::string contents=os.str();
578 char* con=new char[contents.size()+1];
579 strcpy(con, contents.c_str());
580 boost::scoped_ptr<char> buff(con);
581 ierr=MPI_File_write_ordered(mpi_fileHandle_p, buff.get(), contents.size(), MPI_CHAR, &mpi_status);
582 if (ierr != MPI_SUCCESS)
583 {
584 error=1;
585 }
586
587 if (MPI_File_close(&mpi_fileHandle_p)!= MPI_SUCCESS)
588 {
589 error=1;
590 }
591 data[0].getDomain()->MPIBarrier();
592 if (error) // any errors at this stage are from collective routines
593 { // so there is no need to reduce_all
594 throw DataException("saveDataCSVcpp: Error writing and closing file");
595 }
596
597 #endif
598 }
599
600
601 } // end of namespace

  ViewVC Help
Powered by ViewVC 1.1.26