source: XIOS/dev/dev_ym/XIOS_COUPLING/src/distribution/reduce_transform_connector.hpp @ 2230

Last change on this file since 2230 was 2188, checked in by jderouillat, 3 years ago

Fix average reduction (update reference results for test_scalar_algo)
  • Property svn:executable set to *
File size: 11.8 KB
Line 
1#ifndef __REDUCE_TRANSFORM_CONNECTOR_HPP__
2#define __REDUCE_TRANSFORM_CONNECTOR_HPP__
3
4#include "xios_spl.hpp"
5#include "array_new.hpp"
6#include "local_view.hpp"
7#include "reduction_types.hpp"
8#include <algorithm>
9
10
11
12namespace xios
13{
14 
15  class CReduceTransformConnector
16  {
17   
18    private:
19      CLocalView* srcView_;
20      CLocalView* dstView_;
21
22      vector<int> connector_;  //  sizeof sum(nWeights_) 
23      vector<int> nSrc_ ;  //  sizeof dstSize_
24      int srcSize_ ;
25      int dstSize_ ;
26      bool detectMissingValue_ ;
27      EReduction type_ ;
28     
29      typedef void (CReduceTransformConnector::* transferPtr)(int, int, const CArray<double,1>& , CArray<double,1>&) ;
30      transferPtr transfer_ ;
31     
32     void computeConnector(unordered_map<int, std::vector<int>>& indexMap) ;
33     
34    public:
35
36    CReduceTransformConnector(CLocalView* srcView, CLocalView* dstView, EReduction op, unordered_map<int, std::vector<int>>& indexMap, 
37                              bool detectMissingValue=true) ;
38 
39    void transfer(const CArray<double,1>& dataIn, CArray<double,1>& dataOut)
40    {
41      transfer(1,1,dataIn, dataOut) ;
42    }
43        void transfer(int sizeT, const CArray<double,1>& dataIn, CArray<double,1>& dataOut)
44    { 
45      transfer(1,sizeT, dataIn, dataOut) ;
46    }
47
48    void transfer(int repeat, int sizeT, const CArray<double,1>& dataIn, CArray<double,1>& dataOut)
49    {
50      (this->*transfer_)(repeat, sizeT, dataIn, dataOut) ;
51    }
52   
53    private :
54   
55    void transferSum(int repeat, int sizeT, const CArray<double,1>& dataIn, CArray<double,1>& dataOut)
56    {
57      double defaultValue = std::numeric_limits<double>::quiet_NaN();
58
59      int dstSlice = dstSize_*sizeT ;
60      int srcSlice = srcSize_*sizeT ;
61      dataOut.resize(repeat* dstSlice) ;
62      dataOut=0 ;
63           
64      const double* input = dataIn.dataFirst()  ;
65      double* output = dataOut.dataFirst()  ;
66
67      if (detectMissingValue_)
68      {
69        vector<bool> touched(repeat* dstSlice, false) ;
70        int pos=0 ;
71        for(int r=0;r<repeat;r++, input+=srcSlice, output+=dstSlice, pos+=dstSlice)
72        {
73          const double* in = input;
74          double* out = output ;
75          int k=0 ;
76          for(int i=0; i<dstSize_; i++, out+=sizeT, pos+=sizeT)
77          {
78            if (nSrc_[i]==0) for(int l=0; l<sizeT; l++) out[l] = defaultValue ;
79            else 
80            {
81              for(int j=0 ; j<nSrc_[i] ; j++,k++)
82                for(int l=0; l<sizeT; l++) 
83                {
84                  if (! std::isnan(in[connector_[k]*sizeT+l]) ) 
85                  {
86                    if (touched[pos+l]) out[l] += in[connector_[k]*sizeT+l] ;
87                    else 
88                    {
89                      touched[pos+l] = true ;
90                      out[l] = in[connector_[k]*sizeT+l] ;
91                    }
92                  }
93                }
94            }
95          }
96        }
97
98        output = dataOut.dataFirst()  ;
99        pos=0 ;
100        for(int r=0;r<repeat;r++, input+=srcSlice, output+=dstSlice, pos+=dstSlice)
101        {
102          const double* in = input;
103          double* out = output ;
104          int k=0 ;
105          int pos=0 ;
106          for(int i=0; i<dstSize_; i++, out+=sizeT, pos+=sizeT)
107              for(int j=0 ; j<nSrc_[i] ; j++,k++)
108                for(int l=0; l<sizeT; l++) if (!touched[pos+l]) out[l] = defaultValue ;
109        }
110
111      }
112      else
113      {
114        for(int r=0;r<repeat;r++, input+=srcSlice, output+=dstSlice)
115        {
116          const double* in = input;
117          double* out = output ;
118          int k=0 ;
119          for(int i=0; i<dstSize_; i++, out+=sizeT)
120            if (nSrc_[i]==0)
121              for(int l=0; l<sizeT; l++) out[l] = defaultValue ;
122            else 
123            {
124              for(int j=0 ; j<nSrc_[i] ; j++,k++)
125                for(int l=0; l<sizeT; l++) out[l] += in[connector_[k]*sizeT+l] ;
126            }
127        }
128      }
129    }
130
131    void transferMin(int repeat, int sizeT, const CArray<double,1>& dataIn, CArray<double,1>& dataOut)
132    {
133      double defaultValue = std::numeric_limits<double>::quiet_NaN();
134
135      int dstSlice = dstSize_*sizeT ;
136      int srcSlice = srcSize_*sizeT ;
137      dataOut.resize(repeat* dstSlice) ;
138      dataOut=0 ;
139           
140      const double* input = dataIn.dataFirst()  ;
141      double* output = dataOut.dataFirst()  ;
142
143      if (detectMissingValue_)
144      {
145        vector<bool> touched(repeat* dstSlice, false) ;
146        int pos=0 ;
147        for(int r=0;r<repeat;r++, input+=srcSlice, output+=dstSlice, pos+=dstSlice)
148        {
149          const double* in = input;
150          double* out = output ;
151          int k=0 ;
152          for(int i=0; i<dstSize_; i++, out+=sizeT, pos+=sizeT)
153          {
154            if (nSrc_[i]==0) for(int l=0; l<sizeT; l++) out[l] = defaultValue ;
155            else 
156            {
157              for(int j=0 ; j<nSrc_[i] ; j++,k++)
158                for(int l=0; l<sizeT; l++) 
159                {
160                  if (! std::isnan(in[connector_[k]*sizeT+l]) ) 
161                  {
162                    if (touched[pos+l]) out[l]=min(out[l],in[connector_[k]*sizeT+l]) ;
163                    else 
164                    {
165                      touched[pos+l] = true ;
166                      out[l] = in[connector_[k]*sizeT+l] ;
167                    }
168                  }
169                }
170            }
171          }
172        }
173
174        output = dataOut.dataFirst()  ;
175        pos=0 ;
176        for(int r=0;r<repeat;r++, input+=srcSlice, output+=dstSlice, pos+=dstSlice)
177        {
178          const double* in = input;
179          double* out = output ;
180          int k=0 ;
181          int pos=0 ;
182          for(int i=0; i<dstSize_; i++, out+=sizeT, pos+=sizeT)
183              for(int j=0 ; j<nSrc_[i] ; j++,k++)
184                for(int l=0; l<sizeT; l++) if (!touched[pos+l]) out[l] = defaultValue ;
185        }
186
187      }
188      else
189      {
190        for(int r=0;r<repeat;r++, input+=srcSlice, output+=dstSlice)
191        {
192          const double* in = input;
193          double* out = output ;
194          int k=0 ;
195          for(int i=0; i<dstSize_; i++, out+=sizeT)
196            if (nSrc_[i]==0)
197              for(int l=0; l<sizeT; l++) out[l] = defaultValue ;
198            else 
199            {
200              for(int j=0 ; j<nSrc_[i] ; j++,k++)
201                for(int l=0; l<sizeT; l++) out[l]=min(out[l],in[connector_[k]*sizeT+l]) ;
202            }
203        }
204      }
205    }
206
207    void transferMax(int repeat, int sizeT, const CArray<double,1>& dataIn, CArray<double,1>& dataOut)
208    {
209      double defaultValue = std::numeric_limits<double>::quiet_NaN();
210
211      int dstSlice = dstSize_*sizeT ;
212      int srcSlice = srcSize_*sizeT ;
213      dataOut.resize(repeat* dstSlice) ;
214      dataOut=0 ;
215           
216      const double* input = dataIn.dataFirst()  ;
217      double* output = dataOut.dataFirst()  ;
218
219      if (detectMissingValue_)
220      {
221        vector<bool> touched(repeat* dstSlice, false) ;
222        int pos=0 ;
223        for(int r=0;r<repeat;r++, input+=srcSlice, output+=dstSlice, pos+=dstSlice)
224        {
225          const double* in = input;
226          double* out = output ;
227          int k=0 ;
228          for(int i=0; i<dstSize_; i++, out+=sizeT, pos+=sizeT)
229          {
230            if (nSrc_[i]==0) for(int l=0; l<sizeT; l++) out[l] = defaultValue ;
231            else 
232            {
233              for(int j=0 ; j<nSrc_[i] ; j++,k++)
234                for(int l=0; l<sizeT; l++) 
235                {
236                  if (! std::isnan(in[connector_[k]*sizeT+l]) ) 
237                  {
238                    if (touched[pos+l]) out[l]=max(out[l],in[connector_[k]*sizeT+l]) ;
239                    else 
240                    {
241                      touched[pos+l] = true ;
242                      out[l] = in[connector_[k]*sizeT+l] ;
243                    }
244                  }
245                }
246            }
247          }
248        }
249
250        output = dataOut.dataFirst()  ;
251        pos=0 ;
252        for(int r=0;r<repeat;r++, input+=srcSlice, output+=dstSlice, pos+=dstSlice)
253        {
254          const double* in = input;
255          double* out = output ;
256          int k=0 ;
257          int pos=0 ;
258          for(int i=0; i<dstSize_; i++, out+=sizeT, pos+=sizeT)
259              for(int j=0 ; j<nSrc_[i] ; j++,k++)
260                for(int l=0; l<sizeT; l++) if (!touched[pos+l]) out[l] = defaultValue ;
261        }
262
263      }
264      else
265      {
266        for(int r=0;r<repeat;r++, input+=srcSlice, output+=dstSlice)
267        {
268          const double* in = input;
269          double* out = output ;
270          int k=0 ;
271          for(int i=0; i<dstSize_; i++, out+=sizeT)
272            if (nSrc_[i]==0)
273              for(int l=0; l<sizeT; l++) out[l] = defaultValue ;
274            else 
275            {
276              for(int j=0 ; j<nSrc_[i] ; j++,k++)
277                for(int l=0; l<sizeT; l++) out[l]=max(out[l],in[connector_[k]*sizeT+l]) ;
278            }
279        }
280      }
281    }
282
283    void transferAverage(int repeat, int sizeT, const CArray<double,1>& dataIn, CArray<double,1>& dataOut)
284    {
285      double defaultValue = std::numeric_limits<double>::quiet_NaN();
286
287      int dstSlice = dstSize_*sizeT ;
288      int srcSlice = srcSize_*sizeT ;
289      dataOut.resize(repeat* dstSlice) ;
290      dataOut=0 ;
291           
292      const double* input = dataIn.dataFirst()  ;
293      double* output = dataOut.dataFirst()  ;
294
295      if (detectMissingValue_)
296      {
297        vector<int> touched(repeat* dstSlice, 0) ;
298        int pos=0 ;
299        for(int r=0;r<repeat;r++, input+=srcSlice, output+=dstSlice, pos+=dstSlice)
300        {
301          const double* in = input;
302          double* out = output ;
303          int k=0 ;
304          for(int i=0; i<dstSize_; i++, out+=sizeT, pos+=sizeT)
305          {
306            if (nSrc_[i]==0) for(int l=0; l<sizeT; l++) out[l] = defaultValue ;
307            else 
308            {
309              for(int j=0 ; j<nSrc_[i] ; j++,k++)
310                for(int l=0; l<sizeT; l++) 
311                {
312                  if (! std::isnan(in[connector_[k]*sizeT+l]) ) 
313                  {
314                    if (touched[pos+l]) out[l] += in[connector_[k]*sizeT+l] ;
315                    else 
316                    {
317                      out[l] = in[connector_[k]*sizeT+l] ;
318                    }
319                    touched[pos+l]++ ;
320
321                  }
322                }
323              for(int l=0; l<sizeT; l++) {
324                  if (touched[pos+l]!=0)
325                      out[l] /= touched[pos+l];
326              }
327
328            }
329          }
330        }
331
332        output = dataOut.dataFirst()  ;
333        pos=0 ;
334        for(int r=0;r<repeat;r++, input+=srcSlice, output+=dstSlice, pos+=dstSlice)
335        {
336          const double* in = input;
337          double* out = output ;
338          int k=0 ;
339          int pos=0 ;
340          for(int i=0; i<dstSize_; i++, out+=sizeT, pos+=sizeT)
341              for(int j=0 ; j<nSrc_[i] ; j++,k++)
342                for(int l=0; l<sizeT; l++) 
343                  if (touched[pos+l]==0) out[l] = defaultValue ;
344        }
345
346      }
347      else
348      {
349        for(int r=0;r<repeat;r++, input+=srcSlice, output+=dstSlice)
350        {
351          const double* in = input;
352          double* out = output ;
353          int k=0 ;
354          for(int i=0; i<dstSize_; i++, out+=sizeT)
355            if (nSrc_[i]==0)
356              for(int l=0; l<sizeT; l++) out[l] = defaultValue ;
357            else 
358            {
359              for(int j=0 ; j<nSrc_[i] ; j++,k++)
360                for(int l=0; l<sizeT; l++) out[l] += in[connector_[k]*sizeT+l] ;
361            }
362         
363          out = output ;
364          for(int i=0; i<dstSize_; i++, out+=sizeT)
365            if (nSrc_[i]!=0)
366              for(int l=0; l<sizeT; l++) out[l]/=nSrc_[i];
367        }
368      }
369    }
370 
371  };
372
373}
374
375#endif
Note: See TracBrowser for help on using the repository browser.