aboutsummaryrefslogtreecommitdiff
path: root/Carpet/CarpetIOHDF5/src/GetAllActive.cc
blob: bb92f4a7ba68217f43ec8ef43aeeba0f2d3264b3 (plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
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
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
 /*@@
   @file      GetAllActive.cc
   @date      Sun Sep 19 20:34:27 EDT 2010
   @author    Erik Schnetter, Roland Haas
   @desc
              a piece of dh::regrid to compute the active region
   @enddesc
 @@*/

#include <cassert>

#include <vector>

#include "cctk.h"

#include "CarpetIOHDF5.hh"

namespace CarpetIOHDF5 {

  /********************************************************************
   *********************     Internal Routines     ********************
   ********************************************************************/

  static
  void
  assert_error (char const * restrict const checkstring,
                char const * restrict const file, int const line,
                int const ml, int const rl,
                char const * restrict const message);
  static
  void
  assert_error (char const * restrict const checkstring,
                char const * restrict const file, int const line,
                int const ml, int const rl, int const c,
                char const * restrict const message);
  static
  void
  assert_error (char const * restrict const checkstring,
                char const * restrict const file, int const line,
                int const ml, int const rl, int const c, int const cc,
                char const * restrict const message);

  /********************************************************************
   *********************     Local Data   *****************************
   ********************************************************************/

  // accumulates any consistency errors
  static bool there_was_an_error = false;


   /*@@
     @routine    
     @date       Sun Sep 19 20:46:35 EDT 2010
     @author     Erik Schnetter, Roland Haas
     @desc
                 Output warning for consistency checks for a refinement level.
                
     @enddesc
     @history
     @endhistory
     @var        checkstring
     @vdesc      the check that failed
     @vtype      char const * restrict const
     @vio        in
     @endvar
     @var        file
     @vdesc      name of source file where error occured
     @vtype      char const * restrict const
     @vio        in
     @endvar
     @var        line
     @vdesc      source line in file where error occured
     @vtype      int const
     @vio        in
     @endvar
     @var        ml
     @vdesc      meta level where error occured
     @vtype      int const
     @vio        in
     @endvar
     @var        rl
     @vdesc      refinement level where error occured
     @vtype      int const
     @vio        in
     @endvar
     @var        message
     @vdesc      error message to display
     @vtype      char const * restrict const
     @vio        in
     @endvar
     @returntype void
  @@*/
  static
  void
  assert_error (char const * restrict const checkstring,
                char const * restrict const file, int const line,
                int const ml, int const rl,
                char const * restrict const message)
  {
    CCTK_VWarn (CCTK_WARN_ALERT, __LINE__, __FILE__, CCTK_THORNSTRING,
                "\n%s:%d:\n   [ml=%d rl=%d] The following grid structure consistency check failed:\n   %s\n   %s",
                file, line, ml, rl, message, checkstring);
    there_was_an_error = true;
  }

   /*@@
     @routine    
     @date       Sun Sep 19 20:46:35 EDT 2010
     @author     Erik Schnetter, Roland Haas
     @desc
                 Output warning for consistency checks for a single component.
                
     @enddesc
     @history
     @endhistory
     @var        checkstring
     @vdesc      the check that failed
     @vtype      char const * restrict const
     @vio        in
     @endvar
     @var        file
     @vdesc      name of source file where error occured
     @vtype      char const * restrict const
     @vio        in
     @endvar
     @var        line
     @vdesc      source line in file where error occured
     @vtype      int const
     @vio        in
     @endvar
     @var        ml
     @vdesc      meta level where error occured
     @vtype      int const
     @vio        in
     @endvar
     @var        rl
     @vdesc      refinement level where error occured
     @vtype      int const
     @vio        in
     @endvar
     @var        c
     @vdesc      component where error occured
     @vtype      int const
     @vio        in
     @endvar
     @var        message
     @vdesc      error message to display
     @vtype      char const * restrict const
     @vio        in
     @endvar
     @returntype void
  @@*/
  static
  void
  assert_error (char const * restrict const checkstring,
                char const * restrict const file, int const line,
                int const ml, int const rl, int const c,
                char const * restrict const message)
  {
    CCTK_VWarn (CCTK_WARN_ALERT, __LINE__, __FILE__, CCTK_THORNSTRING,
                "\n%s:%d:\n   [ml=%d rl=%d c=%d] The following grid structure consistency check failed:\n   %s\n   %s",
                file, line, ml, rl, c, message, checkstring);
    there_was_an_error = true;
  }

   /*@@
     @routine    
     @date       Sun Sep 19 20:46:35 EDT 2010
     @author     Erik Schnetter, Roland Haas
     @desc
                 Output warning for consistency checks between two components.
                
     @enddesc
     @history
     @endhistory
     @var        checkstring
     @vdesc      the check that failed
     @vtype      char const * restrict const
     @vio        in
     @endvar
     @var        file
     @vdesc      name of source file where error occured
     @vtype      char const * restrict const
     @vio        in
     @endvar
     @var        line
     @vdesc      source line in file where error occured
     @vtype      int const
     @vio        in
     @endvar
     @var        ml
     @vdesc      meta level where error occured
     @vtype      int const
     @vio        in
     @endvar
     @var        rl
     @vdesc      refinement level where error occured
     @vtype      int const
     @vio        in
     @endvar
     @var        c
     @vdesc      first component where error occured
     @vtype      int const
     @vio        in
     @endvar
     @var        cc
     @vdesc      second component where error occured
     @vtype      int const
     @vio        in
     @endvar
     @var        message
     @vdesc      error message to display
     @vtype      char const * restrict const
     @vio        in
     @endvar
     @returntype void
  @@*/
  static
  void
  assert_error (char const * restrict const checkstring,
                char const * restrict const file, int const line,
                int const ml, int const rl, int const c, int const cc,
                char const * restrict const message)
  {
    CCTK_VWarn (CCTK_WARN_ALERT, __LINE__, __FILE__, CCTK_THORNSTRING,
                "\n%s:%d:\n   [ml=%d rl=%d c=%d cc=%d] The following grid structure consistency check failed:\n   %s\n   %s",
                file, line, ml, rl, c, cc, message, checkstring);
    there_was_an_error = true;
  }

#ifdef CARPET_OPTIMISE

  // For highest efficiency, omit all self-checks
#define ASSERT_rl(check, message)
#define ASSERT_c(check, message)
#define ASSERT_cc(check, message)

#else

#define ASSERT_rl(check, message)                                       \
    do {                                                                  \
      if (not (check)) {                                                  \
        assert_error (#check, __FILE__, __LINE__, ml, rl, message);       \
      }                                                                   \
    } while (false)

#define ASSERT_c(check, message)                                        \
    do {                                                                  \
      if (not (check)) {                                                  \
        assert_error (#check, __FILE__, __LINE__, ml, rl, c, message);    \
      }                                                                   \
    } while (false)

#define ASSERT_cc(check, message)                                       \
    do {                                                                  \
      if (not (check)) {                                                  \
        assert_error (#check, __FILE__, __LINE__, ml, rl, c, cc, message); \
      }                                                                   \
    } while (false)

#endif

   /*@@
     @routine    
     @date       Sun Sep 19 20:46:35 EDT 2010
     @author     Erik Schnetter, Roland Haas
     @desc
                 Recompute the union of all active regions. This routine is a
                 stripped down copy of dh::regrid.
     @enddesc
     @history
     @endhistory
     @var        dh
     @vdesc      data hirarchy for this group
     @vtype      const dh *
     @vio        in
     @endvar
     @var        hh
     @vdesc      grid hirarchy
     @vtype      const gh *
     @vio        in
     @endvar
     @var        ml
     @vdesc      meta level to compute active region on
     @vtype      int
     @vio        in
     @endvar
     @var        rl
     @vdesc      refinement level to compute active region on
     @vtype      int
     @vio        in
     @endvar
     @var        allactive
     @vdesc      receives the list of active boxes
     @vtype      ibset &
     @vio        out
     @endvar
     @returntype void
  @@*/
  void
  GetAllActive(const dh *dd, const gh *hh, int ml, int rl, ibset &allactive)
  {
    DECLARE_CCTK_PARAMETERS;
      
    // All owned regions
    ibset allowned;
    
    // Domain:
    
    ibbox const & domain_exterior = hh->baseextent(ml,rl);
    // Variables may have size zero
    // ASSERT_rl (not domain_exterior.empty(),
    //            "The exterior of the domain must not be empty");
    
    i2vect const & buffer_width = dd->buffer_widths.AT(rl);
    i2vect const & boundary_width = hh->boundary_width;
    ASSERT_rl (all (all (boundary_width >= 0)),
               "The gh boundary widths must not be negative");
    
    ibbox const domain_active = domain_exterior.expand (- boundary_width);
    // Variables may have size zero
    // ASSERT_rl (not domain_active.empty(),
    //            "The active part of the domain must not be empty");
    ASSERT_rl (domain_active <= domain_exterior,
               "The active part of the domain must be contained in the exterior part of the domain");
    
    for (int c = 0; c < hh->components(rl); ++ c) {
      
      // Interior:
      
      ibbox intr;
      intr = ibbox::poison();
      
      // The interior of the grid has the extent as specified by the
      // regridding thorn
      intr = hh->extent (ml,rl,c);
      
      // (The interior must not be empty)
      // Variables may have size zero
      // ASSERT_c (not intr.empty(),
      //           "The interior must not be empty");
      
      // The interior must be contained in the domain
      ASSERT_c (intr <= hh->baseextent(ml,rl),
                "The interior must be contained in the domain");
      
      
      
      // Outer boundary faces:
      
      b2vect is_outer_boundary;
      
      // The outer boundary faces are where the interior extends up
      // to the outer boundary of the domain.  It is not possible to
      // check whether it extends past the active part of the
      // domain, since this would be wrong when the outer boundary
      // width is zero.
      is_outer_boundary[0] = intr.lower() == domain_exterior.lower(); 
      is_outer_boundary[1] = intr.upper() == domain_exterior.upper(); 
      
      
      
      // Owned region:
      
      ibbox owned;
      owned = ibbox::poison();
      
      owned = intr.expand (i2vect (is_outer_boundary) * (- boundary_width));
      
      // (The owned region must not be empty)
      // Variables may have size zero
      // ASSERT_c (not owned.empty(),
      //           "The owned region must not be empty");
      
      // The owned region must be contained in the active part of
      // the domain
      ASSERT_c (owned <= domain_active,
                "The owned region must be contained in the active part of the domain");
      
       allowned |= owned;
       ASSERT_rl (allowned <= domain_active,
                "The owned regions must be contained in the active part of the domain");
       
     } // for c
     
     // Enlarge active part of domain
     i2vect const safedist = i2vect (0);
     ibbox const domain_enlarged = domain_active.expand (safedist);
     
     // All not-owned regions
     ibset const notowned = domain_enlarged - allowned;
     
     // All not-active points
     ibset const notactive = notowned.expand (buffer_width);
     
     // All active points
     allactive = allowned - notactive;

     // Fail if any of the consistency checks failed
     assert(there_was_an_error == false);
  }

} // namespace CarpetIOHDF5