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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2642 - (hide annotations)
Tue Sep 1 04:15:50 2009 UTC (10 years, 7 months ago) by jfenwick
File size: 13945 byte(s)
pre cleanup
1 gross 391
2 ksteube 1312 /*******************************************************
3 ksteube 1811 *
4 jfenwick 2548 * Copyright (c) 2003-2009 by University of Queensland
5 ksteube 1811 * 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 ksteube 1312
14 ksteube 1811
15 jfenwick 2641 #include <fstream>
16 ksteube 1806 #include <string.h>
17    
18 jfenwick 2635 // added for saveCSV
19     #include <boost/python.hpp>
20 jfenwick 2641 #include <boost/scoped_ptr.hpp>
21 jfenwick 2635 #include "Data.h"
22    
23 jgs 474 #include "Utils.h"
24 gross 797 #include "DataVector.h"
25 gross 391
26 jgs 478 #ifdef _OPENMP
27     #include <omp.h>
28     #endif
29    
30 ksteube 1561 #ifdef PASO_MPI
31     #include <mpi.h>
32     #endif
33    
34 phornby 1628 #ifdef _WIN32
35     #include <WinSock2.h>
36 phornby 1835 #else
37     #include <unistd.h>
38 phornby 1628 #endif
39    
40 gross 391 namespace escript {
41    
42 ksteube 1247 int getSvnVersion()
43     {
44     #ifdef SVN_VERSION
45     return SVN_VERSION;
46     #else
47     return 0;
48     #endif
49     }
50    
51 ksteube 1620 /* 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 ksteube 1561 void printParallelThreadCnt()
77     {
78     int mpi_iam=0, mpi_num=1;
79 ksteube 1568 char hname[64];
80 ksteube 1561
81 ksteube 1705 #ifdef HAVE_GETHOSTNAME
82 ksteube 1568 gethostname(hname, 64);
83 ksteube 1806 hname[63] = '\0';
84 ksteube 1705 #else
85     strcpy(hname, "unknown host");
86     #endif
87 ksteube 1567
88 ksteube 1561 #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 jfenwick 2607 #pragma omp critical (printthrdcount)
101 ksteube 1620 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 ksteube 1561 }
104     }
105    
106 gross 391 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 gross 2313 ESCRIPT_DLL_API int getMPISizeWorld() {
126 ksteube 1805 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 gross 2313 ESCRIPT_DLL_API int getMPIRankWorld() {
134 ksteube 1805 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 gross 2308 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 jfenwick 2607 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 gross 2308
163 gross 2313 ESCRIPT_DLL_API double getMachinePrecision() {
164 gross 2100 return DBL_EPSILON;
165     }
166 gross 2313 ESCRIPT_DLL_API double getMaxFloat() {
167 gross 2100 return DBL_MAX;
168     }
169 gross 2313 ESCRIPT_DLL_API void MPIBarrierWorld() {
170     #ifdef PASO_MPI
171     MPI_Barrier(MPI_COMM_WORLD );
172     #endif
173     }
174 gross 2100
175 jfenwick 2642
176     /*
177 jfenwick 2635 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 jfenwick 2642 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 jfenwick 2635 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 jfenwick 2637 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 jfenwick 2635 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 gross 2100
282 jfenwick 2635 // 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 jfenwick 2637 step[i]=(data[i].actsExpanded()?DataTypes::noValues(data[i].getDataPointShape()):0);
289 jfenwick 2635 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 jfenwick 2641 std::ostringstream os;
314 jfenwick 2635
315 jfenwick 2641
316 jfenwick 2635 bool first=true;
317 jfenwick 2641
318     if (data[0].getDomain()->getMPIRank()==0)
319 jfenwick 2635 {
320 jfenwick 2641 for (int i=0;i<numdata;++i)
321     {
322 jfenwick 2635 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 jfenwick 2641 }
408     os << endl;
409 jfenwick 2635 }
410 jfenwick 2637 boost::scoped_ptr<BufferGroup> maskbuffer; // sample buffer for the mask [if we have one]
411     const double* masksample=0;
412     int maskoffset=0;
413 jfenwick 2635 //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 jfenwick 2637 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 jfenwick 2641 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 jfenwick 2635 try{
439     for (int i=0;i<numsamples;++i)
440     {
441 jfenwick 2642 if (!best.ownSample(i))
442     {
443     continue;
444     }
445 jfenwick 2637 wantrow=true;
446 jfenwick 2635 for (int d=0;d<numdata;++d)
447     {
448     samples[d]=data[d].getSampleDataRO(i,bg[d].get());
449     }
450 jfenwick 2637 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 jfenwick 2635 for (int j=0;j<dpps;++j)
462     {
463 jfenwick 2637 // now we need to check if this point is masked off
464     if (expandedmask)
465 jfenwick 2635 {
466 jfenwick 2637 wantrow=(masksample[j]>0); // masks are scalar to the relevant value is at [j]
467 jfenwick 2635 }
468 jfenwick 2637 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 jfenwick 2635 }
480     for (int d=0;d<numdata;++d)
481     {
482     offset[d]=0;
483     }
484     }
485     } catch (...)
486     {
487 jfenwick 2641 error=1;
488     #ifndef PASO_MPI
489 jfenwick 2635 throw;
490 jfenwick 2641 #endif
491 jfenwick 2635 }
492 jfenwick 2641 #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 jfenwick 2637
503 jfenwick 2641 // at this point os will contain the text to be written
504     #ifndef PASO_MPI
505 jfenwick 2637
506 jfenwick 2641 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 jfenwick 2642 if (remove(fname_p.get()))
546 jfenwick 2641 {
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 jfenwick 2642
577 jfenwick 2641 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 jfenwick 2635 }
599    
600 jfenwick 2642
601 gross 391 } // end of namespace

  ViewVC Help
Powered by ViewVC 1.1.26