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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2637 - (hide annotations)
Fri Aug 28 04:12:24 2009 UTC (10 years, 7 months ago) by jfenwick
File size: 8891 byte(s)
saveDataCSV looks good but does not do MPI yet.
Also not unit tested.

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 ksteube 1806 #include <string.h>
16    
17 jfenwick 2635 // added for saveCSV
18     #include <boost/python.hpp>
19     #include "Data.h"
20    
21 jgs 474 #include "Utils.h"
22 gross 797 #include "DataVector.h"
23 gross 391
24 jgs 478 #ifdef _OPENMP
25     #include <omp.h>
26     #endif
27    
28 ksteube 1561 #ifdef PASO_MPI
29     #include <mpi.h>
30     #endif
31    
32 phornby 1628 #ifdef _WIN32
33     #include <WinSock2.h>
34 phornby 1835 #else
35     #include <unistd.h>
36 phornby 1628 #endif
37    
38 gross 391 namespace escript {
39    
40 ksteube 1247 int getSvnVersion()
41     {
42     #ifdef SVN_VERSION
43     return SVN_VERSION;
44     #else
45     return 0;
46     #endif
47     }
48    
49 ksteube 1620 /* This is probably not very robust, but it works on Savanna today and is useful for performance analysis */
50     int get_core_id() {
51     int processor_num=-1;
52     #ifdef CORE_ID1
53     FILE *fp;
54     int i, count_spaces=0;
55     char fname[100];
56     char buf[1000];
57    
58     sprintf(fname, "/proc/%d/stat", getpid());
59     fp = fopen(fname, "r");
60     if (fp == NULL) return(-1);
61     fgets(buf, 1000, fp);
62     fclose(fp);
63    
64     for (i=strlen(buf)-1; i>=0; i--) {
65     if (buf[i] == ' ') count_spaces++;
66     if (count_spaces == 4) break;
67     }
68     processor_num = atoi(&buf[i+1]);
69     #endif
70     return(processor_num);
71     }
72    
73    
74 ksteube 1561 void printParallelThreadCnt()
75     {
76     int mpi_iam=0, mpi_num=1;
77 ksteube 1568 char hname[64];
78 ksteube 1561
79 ksteube 1705 #ifdef HAVE_GETHOSTNAME
80 ksteube 1568 gethostname(hname, 64);
81 ksteube 1806 hname[63] = '\0';
82 ksteube 1705 #else
83     strcpy(hname, "unknown host");
84     #endif
85 ksteube 1567
86 ksteube 1561 #ifdef PASO_MPI
87     MPI_Comm_rank(MPI_COMM_WORLD, &mpi_iam);
88     MPI_Comm_size(MPI_COMM_WORLD, &mpi_num);
89     #endif
90    
91     #pragma omp parallel
92     {
93     int omp_iam=0, omp_num=1;
94     #ifdef _OPENMP
95     omp_iam = omp_get_thread_num(); /* Call in a parallel region */
96     omp_num = omp_get_num_threads();
97     #endif
98 jfenwick 2607 #pragma omp critical (printthrdcount)
99 ksteube 1620 printf("printParallelThreadCounts: MPI=%03d/%03d OpenMP=%03d/%03d running on %s core %d\n",
100     mpi_iam, mpi_num, omp_iam, omp_num, hname, get_core_id());
101 ksteube 1561 }
102     }
103    
104 gross 391 void setNumberOfThreads(const int num_threads)
105     {
106    
107     #ifdef _OPENMP
108     omp_set_num_threads(num_threads);
109     #endif
110    
111     }
112    
113     int getNumberOfThreads()
114     {
115     #ifdef _OPENMP
116     return omp_get_max_threads();
117     #else
118     return 1;
119     #endif
120    
121     }
122    
123 gross 2313 ESCRIPT_DLL_API int getMPISizeWorld() {
124 ksteube 1805 int mpi_num = 1;
125     #ifdef PASO_MPI
126     MPI_Comm_size(MPI_COMM_WORLD, &mpi_num);
127     #endif
128     return mpi_num;
129     }
130    
131 gross 2313 ESCRIPT_DLL_API int getMPIRankWorld() {
132 ksteube 1805 int mpi_iam = 0;
133     #ifdef PASO_MPI
134     MPI_Comm_rank(MPI_COMM_WORLD, &mpi_iam);
135     #endif
136     return mpi_iam;
137     }
138    
139 gross 2308 ESCRIPT_DLL_API int getMPIWorldMax(const int val) {
140     #ifdef PASO_MPI
141     int val2 = val;
142     int out = val;
143     MPI_Allreduce( &val2, &out, 1, MPI_INT, MPI_MAX, MPI_COMM_WORLD );
144     #else
145     int out = val;
146     #endif
147     return out;
148     }
149    
150 jfenwick 2607 ESCRIPT_DLL_API int getMPIWorldSum(const int val) {
151     #ifdef PASO_MPI
152     int val2 = val;
153     int out = 0;
154     MPI_Allreduce( &val2, &out, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD );
155     #else
156     int out = val;
157     #endif
158     return out;
159     }
160 gross 2308
161 gross 2313 ESCRIPT_DLL_API double getMachinePrecision() {
162 gross 2100 return DBL_EPSILON;
163     }
164 gross 2313 ESCRIPT_DLL_API double getMaxFloat() {
165 gross 2100 return DBL_MAX;
166     }
167 gross 2313 ESCRIPT_DLL_API void MPIBarrierWorld() {
168     #ifdef PASO_MPI
169     MPI_Barrier(MPI_COMM_WORLD );
170     #endif
171     }
172 gross 2100
173 jfenwick 2635 ESCRIPT_DLL_API
174     void
175     saveDataCSV(const std::string& filename, boost::python::dict arg, const std::string& sep, const std::string& csep,
176     bool append)
177     {
178     using std::cout;
179     using std::endl;
180     boost::python::list keys=arg.keys();
181     int numdata = boost::python::extract<int>(arg.attr("__len__")());
182 jfenwick 2637 bool hasmask=arg.has_key("mask");
183     Data mask;
184     if (hasmask)
185     {
186     mask=boost::python::extract<escript::Data>(arg["mask"]);
187     keys.remove("mask");
188     numdata--;
189     if (mask.getDataPointRank()!=0)
190     {
191     throw DataException("saveDataCSVcpp: masks must be scalar.");
192     }
193     }
194 jfenwick 2635 if (numdata<1)
195     {
196     throw DataException("saveDataCSVcpp: no data to save specified.");
197     }
198     std::vector<int> step(numdata);
199     std::vector<std::string> names(numdata);
200     std::vector<Data> data(numdata);
201     std::vector<const DataAbstract::ValueType::value_type*> samples(numdata);
202     std::vector<int> offset(numdata);
203     std::vector<int> fstypes(numdata); // FunctionSpace types for each data
204 gross 2100
205 jfenwick 2635 // We need to interpret the samples correctly even if they are different types
206     // for this reason, we should interate over samples
207     for (int i=0;i<numdata;++i)
208     {
209     names[i]=boost::python::extract<std::string>(keys[i]);
210     data[i]=boost::python::extract<escript::Data>(arg[keys[i]]);
211 jfenwick 2637 step[i]=(data[i].actsExpanded()?DataTypes::noValues(data[i].getDataPointShape()):0);
212 jfenwick 2635 fstypes[i]=data[i].getFunctionSpace().getTypeCode();
213     if (i>0)
214     {
215     if (data[i].getDomain()!=data[i-1].getDomain())
216     {
217     throw DataException("saveDataCSVcpp: all data must be on the same domain.");
218     }
219     }
220     }
221     int bestfnspace=0;
222     if (!data[0].getDomain()->commonFunctionSpace(fstypes, bestfnspace))
223     {
224     throw DataException("saveDataCSVcpp: FunctionSpaces of data are incompatible");
225     }
226     // now we interpolate all data to the same type
227     FunctionSpace best(data[0].getDomain(),bestfnspace);
228     for (int i=0;i<numdata;++i)
229     {
230     data[i]=data[i].interpolate(best);
231     }
232     int numsamples=data[0].getNumSamples(); // these must be the same for all data
233     int dpps=data[0].getNumDataPointsPerSample();
234    
235    
236     std::ofstream os;
237     if (append)
238     {
239     os.open(filename.c_str(), std::ios_base::app);
240     }
241     else
242     {
243     os.open(filename.c_str());
244     }
245     if (!os.is_open())
246     {
247     throw DataException("saveDataCSVcpp: unable to open file for writing");
248     }
249    
250     bool first=true;
251     for (int i=0;i<numdata;++i)
252     {
253     const DataTypes::ShapeType& s=data[i].getDataPointShape();
254     switch (data[i].getDataPointRank())
255     {
256     case 0: if (!first)
257     {
258     os << sep;
259     }
260     else
261     {
262     first=false;
263     }
264     os << names[i]; break;
265     case 1: for (int j=0;j<s[0];++j)
266     {
267     if (!first)
268     {
269     os << sep;
270     }
271     else
272     {
273     first=false;
274     }
275     os << names[i] << csep << j;
276     }
277     break;
278     case 2: for (int j=0;j<s[0];++j)
279     {
280     for (int k=0;k<s[1];++k)
281     {
282     if (!first)
283     {
284     os << sep;
285     }
286     else
287     {
288     first=false;
289     }
290     os << names[i] << csep << k << csep << j;
291     }
292     }
293     break;
294     case 3: for (int j=0;j<s[0];++j)
295     {
296     for (int k=0;k<s[1];++k)
297     {
298     for (int l=0;l<s[2];++l)
299     {
300     if (!first)
301     {
302     os << sep;
303     }
304     else
305     {
306     first=false;
307     }
308     os << names[i] << csep << k << csep << j << csep << l;
309     }
310     }
311     }
312     break;
313     case 4: for (int j=0;j<s[0];++j)
314     {
315     for (int k=0;k<s[1];++k)
316     {
317     for (int l=0;l<s[2];++l)
318     {
319     for (int m=0;m<s[3];++m)
320     {
321     if (!first)
322     {
323     os << sep;
324     }
325     else
326     {
327     first=false;
328     }
329     os << names[i] << csep << k << csep << j << csep << l << csep << m;
330     }
331     }
332     }
333     }
334     break;
335     default:
336     throw DataException("saveDataCSV: Illegal rank");
337     }
338     }
339     os << endl;
340    
341 jfenwick 2637 boost::scoped_ptr<BufferGroup> maskbuffer; // sample buffer for the mask [if we have one]
342     const double* masksample=0;
343     int maskoffset=0;
344 jfenwick 2635 //the use of shared_ptr here is just to ensure the buffer group is freed
345     //I would have used scoped_ptr but they don't work in vectors
346     std::vector<boost::shared_ptr<BufferGroup> > bg(numdata);
347     for (int d=0;d<numdata;++d)
348     {
349     bg[d].reset(data[d].allocSampleBuffer());
350     }
351    
352 jfenwick 2637 bool expandedmask=false; // does the mask act expanded. Are there mask value for each point in the sample
353     bool wantrow=true; // do we output this row?
354     if (hasmask)
355     {
356     maskbuffer.reset(mask.allocSampleBuffer());
357     if (mask.actsExpanded())
358     {
359     maskoffset=DataTypes::noValues(mask.getDataPointShape());
360     expandedmask=true;
361     }
362     }
363 jfenwick 2635 try{
364     for (int i=0;i<numsamples;++i)
365     {
366 jfenwick 2637 wantrow=true;
367 jfenwick 2635 for (int d=0;d<numdata;++d)
368     {
369     samples[d]=data[d].getSampleDataRO(i,bg[d].get());
370     }
371 jfenwick 2637 if (hasmask)
372     {
373     masksample=mask.getSampleDataRO(i, maskbuffer.get());
374     if (!expandedmask) // mask controls whole sample
375     {
376     if (masksample[0]<=0) // masks are scalar
377     {
378     wantrow=false;
379     }
380     }
381     }
382 jfenwick 2635 for (int j=0;j<dpps;++j)
383     {
384 jfenwick 2637 // now we need to check if this point is masked off
385     if (expandedmask)
386 jfenwick 2635 {
387 jfenwick 2637 wantrow=(masksample[j]>0); // masks are scalar to the relevant value is at [j]
388 jfenwick 2635 }
389 jfenwick 2637 if (wantrow)
390     {
391     bool needsep=false;
392     for (int d=0;d<numdata;++d)
393     {
394     DataTypes::pointToStream(os, samples[d], data[d].getDataPointShape(), offset[d], needsep, sep);
395     needsep=true;
396     offset[d]+=step[d];
397     }
398     os << endl;
399     }
400 jfenwick 2635 }
401     for (int d=0;d<numdata;++d)
402     {
403     offset[d]=0;
404     }
405     }
406     } catch (...)
407     {
408     os.close();
409     throw;
410     }
411     os.close();
412 jfenwick 2637
413     cout << "This method is not MPI safe" << endl;
414    
415 jfenwick 2635 }
416    
417 gross 391 } // end of namespace

  ViewVC Help
Powered by ViewVC 1.1.26