boundarytypes.hh 12.8 KB
Newer Older
Andreas Lauser's avatar
Andreas Lauser committed
1
// $Id$
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
/*****************************************************************************
 *   Copyright (C) 2009 by Andreas Lauser                                    *
 *   Institute of Hydraulic Engineering                                      *
 *   University of Stuttgart, Germany                                        *
 *   email: <givenname>.<name>@iws.uni-stuttgart.de                          *
 *                                                                           *
 *   This program is free software; you can redistribute it and/or modify    *
 *   it under the terms of the GNU General Public License as published by    *
 *   the Free Software Foundation; either version 2 of the License, or       *
 *   (at your option) any later version, as long as this copyright notice    *
 *   is included in its original form.                                       *
 *                                                                           *
 *   This program is distributed WITHOUT ANY WARRANTY.                       *
 *****************************************************************************/
/*!
 * \file
 * \brief Class to specify the type of a boundary.
 */
#ifndef BOUNDARY_TYPES_HH
#define BOUNDARY_TYPES_HH



#include <dumux/common/valgrind.hh>

namespace Dumux
{

30
31
32
/*!
 * \brief Class to specify the type of a boundary.
 */
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
template <int numEq>
class BoundaryTypes
{
public:
    BoundaryTypes()
    { reset(); }

    /*!
     * \brief Reset the boundary types.
     *
     * After this method no equations will be disabled and neither
     * neumann nor dirichlet conditions will be evaluated. This
     * corrosponds to a neumann zero boundary.
     */
    void reset()
    {
        for (int i=0; i < numEq; ++i) {
            boundaryInfo_[i].visited = 0;

            boundaryInfo_[i].isDirichlet = 0;
Andreas Lauser's avatar
Andreas Lauser committed
53
            boundaryInfo_[i].isNeumann = 0;
54
            boundaryInfo_[i].isOutflow = 0;
55
56
            boundaryInfo_[i].isCouplingInflow = 0;
            boundaryInfo_[i].isCouplingOutflow = 0;
57
58
59
60
61
62
63

            eq2pvIdx_[i] = i;
            pv2eqIdx_[i] = i;
        };
    }

    /*!
64
65
66
67
     * \brief Returns true if the boundary types for a given equation
     *        has been specified.
     *
     * \param eqIdx The index of the equation
68
69
70
71
72
73
74
     */
    bool isSet(int eqIdx) const
    { return boundaryInfo_[eqIdx].visited; };

    /*!
     * \brief Make sure the boundary conditions are well-posed.
     *
75
76
     * If they are not, an assertation fails and the program aborts!
     * (if the NDEBUG macro is not defined)
77
78
79
     */
    void checkWellPosed() const
    {
80
81
#ifndef NDEBUG
        for (int i=0; i < numEq; ++i)
82
83
            // if this fails, at least one condition is missing.
            assert(boundaryInfo_[i].visited);
84
#endif
85
86
87
88
89
90
91
92
93
94
    };

    /*!
     * \brief Set all boundary conditions to neuman.
     */
    void setAllNeumann()
    {
        for (int eqIdx = 0; eqIdx < numEq; ++eqIdx) {
            boundaryInfo_[eqIdx].visited = 1;
            boundaryInfo_[eqIdx].isDirichlet = 0;
Andreas Lauser's avatar
Andreas Lauser committed
95
            boundaryInfo_[eqIdx].isNeumann = 1;
96
            boundaryInfo_[eqIdx].isOutflow = 0;
97
98
            boundaryInfo_[eqIdx].isCouplingInflow = 0;
            boundaryInfo_[eqIdx].isCouplingOutflow = 0;
99
100
101
102
103
104
105
106
107
108
109
110
111

            Valgrind::SetDefined(boundaryInfo_[eqIdx]);
        }
    }

    /*!
     * \brief Set all boundary conditions to dirichlet.
     */
    void setAllDirichlet()
    {
        for (int eqIdx = 0; eqIdx < numEq; ++ eqIdx) {
            boundaryInfo_[eqIdx].visited = 1;
            boundaryInfo_[eqIdx].isDirichlet = 1;
Andreas Lauser's avatar
Andreas Lauser committed
112
            boundaryInfo_[eqIdx].isNeumann = 0;
113
            boundaryInfo_[eqIdx].isOutflow = 0;
114
115
            boundaryInfo_[eqIdx].isCouplingInflow = 0;
            boundaryInfo_[eqIdx].isCouplingOutflow = 0;
116
117
118
119
120
121
122
123

            eq2pvIdx_[eqIdx] = eqIdx;
            pv2eqIdx_[eqIdx] = eqIdx;

            Valgrind::SetDefined(boundaryInfo_[eqIdx]);
        }
    }

124
125
126
127
128
129
130
131
132
133
    /*!
     * \brief Set all boundary conditions to neuman.
     */
    void setAllOutflow()
    {
        for (int eqIdx = 0; eqIdx < numEq; ++eqIdx) {
            boundaryInfo_[eqIdx].visited = 1;
            boundaryInfo_[eqIdx].isDirichlet = 0;
            boundaryInfo_[eqIdx].isNeumann = 0;
            boundaryInfo_[eqIdx].isOutflow = 1;
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
            boundaryInfo_[eqIdx].isCouplingInflow = 0;
            boundaryInfo_[eqIdx].isCouplingOutflow = 0;

            Valgrind::SetDefined(boundaryInfo_[eqIdx]);
        }
    }

    /*!
     * \brief Set all boundary conditions to coupling inflow.
     */
    void setAllCouplingInflow()
    {
        for (int eqIdx = 0; eqIdx < numEq; ++eqIdx) {
            boundaryInfo_[eqIdx].visited = 1;
            boundaryInfo_[eqIdx].isDirichlet = 0;
            boundaryInfo_[eqIdx].isNeumann = 0;
            boundaryInfo_[eqIdx].isOutflow = 0;
            boundaryInfo_[eqIdx].isCouplingInflow = 1;
            boundaryInfo_[eqIdx].isCouplingOutflow = 0;

            Valgrind::SetDefined(boundaryInfo_[eqIdx]);
        }
    }

    /*!
     * \brief Set all boundary conditions to coupling outflow.
     */
    void setAllCouplingOutflow()
    {
        for (int eqIdx = 0; eqIdx < numEq; ++eqIdx) {
            boundaryInfo_[eqIdx].visited = 1;
            boundaryInfo_[eqIdx].isDirichlet = 0;
            boundaryInfo_[eqIdx].isNeumann = 0;
167
            boundaryInfo_[eqIdx].isOutflow = 0;
168
169
            boundaryInfo_[eqIdx].isCouplingInflow = 0;
            boundaryInfo_[eqIdx].isCouplingOutflow = 1;
170
171
172
173
174

            Valgrind::SetDefined(boundaryInfo_[eqIdx]);
        }
    }

175
176
177
    /*!
     * \brief Set a neumann boundary condition for a single a single
     *        equation.
178
179
     *
     * \param eqIdx The index of the equation
180
181
182
183
     */
    void setNeumann(int eqIdx)
    {
        boundaryInfo_[eqIdx].visited = 1;
184
185
186
        boundaryInfo_[eqIdx].isDirichlet = 0;
        boundaryInfo_[eqIdx].isNeumann = 1;
        boundaryInfo_[eqIdx].isOutflow = 0;
187
188
        boundaryInfo_[eqIdx].isCouplingInflow = 0;
        boundaryInfo_[eqIdx].isCouplingOutflow = 0;
189
190
191
192
193
194
195
196

        Valgrind::SetDefined(boundaryInfo_[eqIdx]);
    }

    /*!
     * \brief Set a dirichlet boundary condition for a single primary
     *        variable
     *
197
198
199
200
     * \param pvIdx The index of the primary variable for which the
     *              Dirichlet condition should apply.
     * \param eqIdx The index of the equation which should used to set
     *              the Dirichlet condition
201
     */
202
    void setDirichlet(int pvIdx, int eqIdx)
203
204
205
    {
        boundaryInfo_[eqIdx].visited = 1;
        boundaryInfo_[eqIdx].isDirichlet = 1;
Andreas Lauser's avatar
Andreas Lauser committed
206
        boundaryInfo_[eqIdx].isNeumann = 0;
207
        boundaryInfo_[eqIdx].isOutflow = 0;
208
209
        boundaryInfo_[eqIdx].isCouplingInflow = 0;
        boundaryInfo_[eqIdx].isCouplingOutflow = 0;
210
211
212
213
214
215
216
217

        // update the equation <-> primary variable mapping
        eq2pvIdx_[eqIdx] = pvIdx;
        pv2eqIdx_[pvIdx] = eqIdx;

        Valgrind::SetDefined(boundaryInfo_[eqIdx]);
    }

218
219
220
    /*!
     * \brief Set a neumann boundary condition for a single a single
     *        equation.
221
222
223
     *
     * \param eqIdx The index of the equation on which the outflow
     *              condition applies.
224
225
226
227
228
229
230
     */
    void setOutflow(int eqIdx)
    {
        boundaryInfo_[eqIdx].visited = 1;
        boundaryInfo_[eqIdx].isDirichlet = 0;
        boundaryInfo_[eqIdx].isNeumann = 0;
        boundaryInfo_[eqIdx].isOutflow = 1;
231
232
        boundaryInfo_[eqIdx].isCouplingInflow = 0;
        boundaryInfo_[eqIdx].isCouplingOutflow = 0;
233
234
235
236

        Valgrind::SetDefined(boundaryInfo_[eqIdx]);
    }

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
    /*!
     * \brief Set a boundary condition for a single equation to coupling inflow.
     */
    void setCouplingInflow(int eqIdx)
    {
        boundaryInfo_[eqIdx].visited = 1;
        boundaryInfo_[eqIdx].isDirichlet = 0;
        boundaryInfo_[eqIdx].isNeumann = 0;
        boundaryInfo_[eqIdx].isOutflow = 0;
        boundaryInfo_[eqIdx].isCouplingInflow = 1;
        boundaryInfo_[eqIdx].isCouplingOutflow = 0;

        Valgrind::SetDefined(boundaryInfo_[eqIdx]);
    }

    /*!
     * \brief Set a boundary condition for a single equation to coupling outflow.
     */
    void setCouplingOutflow(int eqIdx)
    {
        boundaryInfo_[eqIdx].visited = 1;
        boundaryInfo_[eqIdx].isDirichlet = 0;
        boundaryInfo_[eqIdx].isNeumann = 0;
        boundaryInfo_[eqIdx].isOutflow = 0;
        boundaryInfo_[eqIdx].isCouplingInflow = 1;
        boundaryInfo_[eqIdx].isCouplingOutflow = 0;

        Valgrind::SetDefined(boundaryInfo_[eqIdx]);
    }

267
268
    /*!
     * \brief Set a dirichlet boundary condition for a single primary
269
270
271
272
273
     *        variable. 
     *
     * WARNING: This method assumes that the equation with the same
     * index as the primary variable to be set is used to specify the
     * Dirichlet condition. USE WITH _GREAT_ CARE!
274
275
276
277
278
279
280
281
282
283
284
285
     *
     * \param eqIdx The index of the equation which is assumed to be
     *              equal to the index of the primary variable
     */
    void setDirichlet(int eqIdx)
    {
        setDirichlet(eqIdx, eqIdx);
    }

    /*!
     * \brief Returns true if an equation is used to specify a
     *        dirichlet condition.
286
287
     *
     * \param eqIdx The index of the equation
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
     */
    bool isDirichlet(unsigned eqIdx) const
    { return boundaryInfo_[eqIdx].isDirichlet; };

    /*!
     * \brief Returns true if some equation is used to specify a
     *        dirichlet condition.
     */
    bool hasDirichlet() const
    {
        for (int i = 0; i < numEq; ++i)
            if (boundaryInfo_[i].isDirichlet)
                return true;
        return false;
    };

    /*!
     * \brief Returns true if an equation is used to specify a
     *        neumann condition.
307
308
     *
     * \param eqIdx The index of the equation
309
310
     */
    bool isNeumann(unsigned eqIdx) const
Andreas Lauser's avatar
Andreas Lauser committed
311
    { return boundaryInfo_[eqIdx].isNeumann; };
312
313
314
315
316
317
318
319

    /*!
     * \brief Returns true if some equation is used to specify a
     *        neumann condition.
     */
    bool hasNeumann() const
    {
        for (int i = 0; i < numEq; ++i)
Andreas Lauser's avatar
Andreas Lauser committed
320
            if (boundaryInfo_[i].isNeumann)
321
322
323
324
                return true;
        return false;
    };

325
326
327
    /*!
     * \brief Returns true if an equation is used to specify an
     *        outflow condition.
328
329
     *
     * \param eqIdx The index of the equation
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
     */
    bool isOutflow(unsigned eqIdx) const
    { return boundaryInfo_[eqIdx].isOutflow; };

    /*!
     * \brief Returns true if some equation is used to specify an
     *        outflow condition.
     */
    bool hasOutflow() const
    {
        for (int i = 0; i < numEq; ++i)
            if (boundaryInfo_[i].isOutflow)
                return true;
        return false;
    };

346
347
348
    /*!
     * \brief Returns true if an equation is used to specify an
     *        inflow coupling condition.
349
350
     *
     * \param eqIdx The index of the equation
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
     */
    bool isCouplingInflow(unsigned eqIdx) const
    { return boundaryInfo_[eqIdx].isCouplingInflow; };

    /*!
     * \brief Returns true if some equation is used to specify an
     *        inflow coupling condition.
     */
    bool hasCouplingInflow() const
    {
        for (int i = 0; i < numEq; ++i)
            if (boundaryInfo_[i].isCouplingInflow)
                return true;
        return false;
    };

    /*!
     * \brief Returns true if an equation is used to specify an
     *        outflow coupling condition.
370
371
     *
     * \param eqIdx The index of the equation
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
     */
    bool isCouplingOutflow(unsigned eqIdx) const
    { return boundaryInfo_[eqIdx].isCouplingOutflow; };

    /*!
     * \brief Returns true if some equation is used to specify an
     *        outflow coupling condition.
     */
    bool hasCouplingOutflow() const
    {
        for (int i = 0; i < numEq; ++i)
            if (boundaryInfo_[i].isCouplingOutflow)
                return true;
        return false;
    };

388
389
    /*!
     * \brief Returns the index of the equation which should be used
390
     *        for the Dirichlet condition of the pvIdx's primary
391
     *        variable.
392
393
394
     *
     * \param pvIdx The index of the primary variable which is be set
     *              by the Dirichlet condition.
395
396
397
398
399
400
     */
    unsigned dirichletToEqIndex(unsigned pvIdx) const
    { return pv2eqIdx_[pvIdx]; };

    /*!
     * \brief Returns the index of the primary variable which should
401
     *        be used for the Dirichlet condition given an equation
402
     *        index.
403
404
405
     *
     * \param eqIdx The index of the equation which is used to set
     *              the Dirichlet condition.
406
407
408
409
410
411
412
413
414
     */
    unsigned eqToDirichletIndex(unsigned eqIdx) const
    { return eq2pvIdx_[eqIdx]; };

private:
    // this is a bitfield structure!
    struct __packed__ {
        unsigned char visited : 1;
        unsigned char isDirichlet : 1;
Andreas Lauser's avatar
Andreas Lauser committed
415
        unsigned char isNeumann : 1;
416
        unsigned char isOutflow : 1;
417
418
        unsigned char isCouplingInflow : 1;
        unsigned char isCouplingOutflow : 1;
419
420
421
422
423
424
425
426
427
    } boundaryInfo_[numEq];

    unsigned char eq2pvIdx_[numEq];
    unsigned char pv2eqIdx_[numEq];
};

}

#endif