← Back to branch summary

~ubuntu-branches/ubuntu/saucy/deal.ii/saucy

« back to all changes in this revision

Viewing changes to deal.II/include/grid/grid_generator.h

  • Committer: Bazaar Package Importer
  • Author(s): Adam C. Powell, IV
  • Date: 2009-05-08 23:13:50 UTC
  • Revision ID: james.westby@ubuntu.com-20090508231350-rrh1ltgi0tifabwc
Tags: upstream-6.2.0
Import upstream version 6.2.0

Show diffs side-by-side

added added

removed removed

Lines of Context:
deal.II/include/grid/grid_generator.h
 
1
//---------------------------------------------------------------------------
 
2
//    $Id: grid_generator.h 17866 2008-12-05 22:27:44Z bangerth $
 
3
//    Version: $Name$
 
4
//
 
5
//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008 by the deal.II authors
 
6
//
 
7
//    This file is subject to QPL and may not be  distributed
 
8
//    without copyright and license information. Please refer
 
9
//    to the file deal.II/doc/license.html for the  text  and
 
10
//    further information on this license.
 
11
//
 
12
//---------------------------------------------------------------------------
 
13
#ifndef __deal2__grid_generator_h
 
14
#define __deal2__grid_generator_h
 
15
 
 
16
 
 
17
#include <base/config.h>
 
18
#include <base/exceptions.h>
 
19
#include <base/point.h>
 
20
#include <base/table.h>
 
21
#include <grid/tria.h>
 
22
#include <map>
 
23
 
 
24
DEAL_II_NAMESPACE_OPEN
 
25
 
 
26
template <int dim, int spacedim> class Triangulation;
 
27
template <typename number> class Vector;
 
28
template <typename number> class SparseMatrix;
 
29
 
 
30
 
 
31
/**
 
32
 * This class provides a collection of functions for generating basic
 
33
 * triangulations. Below, we try to provide some pictures in order to
 
34
 * illustrate at least the more complex ones.
 
35
 *
 
36
 * Some of these functions receive a flag @p colorize. If this is
 
37
 * set, parts of the boundary receive different boundary numbers,
 
38
 * allowing them to be distinguished by application programs. See the
 
39
 * documentation of the functions for details.
 
40
 * 
 
41
 * Additionally this class provides a function
 
42
 * (@p laplace_transformation) that smoothly transforms a grid
 
43
 * according to given new boundary points. This can be used to
 
44
 * transform (simple-shaped) grids to a more complicated ones, like a
 
45
 * shell onto a grid of an airfoil, for example.
 
46
 *
 
47
 * No meshes for the codimension one case are provided at the moment.
 
48
 *
 
49
 *
 
50
 * @ingroup grid
 
51
 * @author Wolfgang Bangerth, Ralf Hartmann, Guido Kanschat, Stefan
 
52
 * Nauber, Joerg Weimar, Yaqi Wang, Luca Heltai, 1998, 1999, 2000, 2001, 2002,
 
53
 * 2003, 2006, 2007, 2008.
 
54
 */
 
55
class GridGenerator
 
56
{
 
57
  public:
 
58
                                     /**
 
59
                                      * Initialize the given
 
60
                                      * triangulation with a hypercube
 
61
                                      * (line in 1D, square in 2D,
 
62
                                      * etc) consisting of exactly one
 
63
                                      * cell. The hypercube volume is
 
64
                                      * the tensor product interval
 
65
                                      * <i>[left,right]<sup>dim</sup></i>
 
66
                                      * in the present number of
 
67
                                      * dimensions, where the limits
 
68
                                      * are given as arguments. They
 
69
                                      * default to zero and unity,
 
70
                                      * then producing the unit
 
71
                                      * hypercube.
 
72
                                      *
 
73
                                      * @image html hyper_cubes.png
 
74
                                      *
 
75
                                      * See also
 
76
                                      * subdivided_hyper_cube() for a
 
77
                                      * coarse mesh consisting of
 
78
                                      * several cells. See
 
79
                                      * hyper_rectangle(), if
 
80
                                      * different lengths in different
 
81
                                      * ordinate directions are
 
82
                                      * required.
 
83
                                      *
 
84
                                      * @note The triangulation needs to be
 
85
                                      * void upon calling this
 
86
                                      * function.
 
87
                                      */
 
88
    template <int dim, int spacedim>
 
89
    static void hyper_cube (Triangulation<dim,spacedim>  &tria,
 
90
                            const double        left = 0.,
 
91
                            const double        right= 1.);
 
92
 
 
93
                                     /**
 
94
                                      * Same as hyper_cube(), but
 
95
                                      * with the difference that not
 
96
                                      * only one cell is created but
 
97
                                      * each coordinate direction is
 
98
                                      * subdivided into
 
99
                                      * @p repetitions cells. Thus,
 
100
                                      * the number of cells filling
 
101
                                      * the given volume is
 
102
                                      * <tt>repetitions<sup>dim</sup></tt>.
 
103
                                      *
 
104
                                      * If spacedim=dim+1 the same
 
105
                                      * mesh as in the case
 
106
                                      * spacedim=dim is created, but
 
107
                                      * the vertices have an
 
108
                                      * additional coordinate =0. So,
 
109
                                      * if dim=1 one obtains line
 
110
                                      * along the x axis in the xy
 
111
                                      * plane, and if dim=3 one
 
112
                                      * obtains a square in lying in
 
113
                                      * the xy plane in 3d space.
 
114
                                      *
 
115
                                      * @note The triangulation needs
 
116
                                      * to be void upon calling this
 
117
                                      * function.
 
118
                                      */
 
119
    template <int dim>
 
120
    static void subdivided_hyper_cube (Triangulation<dim>  &tria,
 
121
                                       const unsigned int  repetitions,
 
122
                                       const double        left = 0.,
 
123
                                       const double        right= 1.);
 
124
 
 
125
                                     /**
 
126
                                      * Create a coordinate-parallel
 
127
                                      * brick from the two
 
128
                                      * diagonally opposite corner
 
129
                                      * points @p p1 and @p p2.
 
130
                                      *
 
131
                                      * If the @p colorize flag is
 
132
                                      * set, the
 
133
                                      * @p boundary_indicators of the
 
134
                                      * surfaces are assigned, such
 
135
                                      * that the lower one in
 
136
                                      * @p x-direction is 0, the
 
137
                                      * upper one is 1. The indicators
 
138
                                      * for the surfaces in
 
139
                                      * @p y-direction are 2 and 3,
 
140
                                      * the ones for @p z are 4 and
 
141
                                      * 5. Additionally, material ids
 
142
                                      * are assigned to the cells
 
143
                                      * according to the octant their
 
144
                                      * center is in: being in the right half
 
145
                                      * plane for any coordinate
 
146
                                      * direction <i>x<sub>i</sub></i>
 
147
                                      * adds 2<sup>i</sup>. For
 
148
                                      * instance, the center point
 
149
                                      * (1,-1,1) yields a material id 5.
 
150
                                      * 
 
151
                                      * @note The triangulation needs to be
 
152
                                      * void upon calling this
 
153
                                      * function.
 
154
                                      */
 
155
    template <int dim, int spacedim>
 
156
    static void hyper_rectangle (Triangulation<dim,spacedim> &tria,
 
157
                                 const Point<spacedim>       &p1,
 
158
                                 const Point<spacedim>       &p2,
 
159
                                 const bool                  colorize = false);
 
160
 
 
161
                                     /**
 
162
                                      * Create a coordinate-parallel
 
163
                                      * parallelepiped from the two
 
164
                                      * diagonally opposite corner
 
165
                                      * points @p p1 and @p p2. In
 
166
                                      * dimension @p i,
 
167
                                      * <tt>repetitions[i]</tt> cells are
 
168
                                      * generated.
 
169
                                      * 
 
170
                                      * To get cells with an aspect
 
171
                                      * ratio different from that of
 
172
                                      * the domain, use different
 
173
                                      * numbers of subdivisions in
 
174
                                      * different coordinate
 
175
                                      * directions. The minimum number
 
176
                                      * of subdivisions in each
 
177
                                      * direction is
 
178
                                      * 1. @p repetitions is a list
 
179
                                      * of integers denoting the
 
180
                                      * number of subdivisions in each
 
181
                                      * coordinate direction.
 
182
                                      * 
 
183
                                      * If the @p colorize flag is
 
184
                                      * set, the
 
185
                                      * @p boundary_indicators of the
 
186
                                      * surfaces are assigned, such
 
187
                                      * that the lower one in
 
188
                                      * @p x-direction is 0, the
 
189
                                      * upper one is 1. The indicators
 
190
                                      * for the surfaces in
 
191
                                      * @p y-direction are 2 and 3,
 
192
                                      * the ones for @p z are 4 and
 
193
                                      * 5.  Additionally, material ids
 
194
                                      * are assigned to the cells
 
195
                                      * according to the octant their
 
196
                                      * center is in: being in the right half
 
197
                                      * plane for any coordinate
 
198
                                      * direction <i>x<sub>i</sub></i>
 
199
                                      * adds 2<sup>i</sup>. For
 
200
                                      * instance, the center point
 
201
                                      * (1,-1,1) yields a material id 5.
 
202
                                      *
 
203
                                      * @note The triangulation needs to be
 
204
                                      * void upon calling this
 
205
                                      * function.
 
206
                                      *
 
207
                                      * @note For an example of the
 
208
                                      * use of this function see the
 
209
                                      * @ref step_28 "step-28"
 
210
                                      * tutorial program.
 
211
                                      */
 
212
    template <int dim>
 
213
    static
 
214
    void
 
215
    subdivided_hyper_rectangle (Triangulation<dim>              &tria,
 
216
                                const std::vector<unsigned int> &repetitions,
 
217
                                const Point<dim>                &p1,
 
218
                                const Point<dim>                &p2,
 
219
                                const bool                      colorize=false);
 
220
 
 
221
                                     /**
 
222
                                      * Like the previous
 
223
                                      * function. However, here the
 
224
                                      * second argument does not
 
225
                                      * denote the number of
 
226
                                      * subdivisions in each
 
227
                                      * coordinate direction, but a
 
228
                                      * sequence of step sizes for
 
229
                                      * each coordinate direction. The
 
230
                                      * domain will therefore be
 
231
                                      * subdivided into
 
232
                                      * <code>step_sizes[i].size()</code>
 
233
                                      * cells in coordinate direction
 
234
                                      * <code>i</code>, with widths
 
235
                                      * <code>step_sizes[i][j]</code>
 
236
                                      * for the <code>j</code>th cell.
 
237
                                      *
 
238
                                      * This function is therefore the
 
239
                                      * right one to generate graded
 
240
                                      * meshes where cells are
 
241
                                      * concentrated in certain areas,
 
242
                                      * rather than a uniformly
 
243
                                      * subdivided mesh as the
 
244
                                      * previous function generates.
 
245
                                      *
 
246
                                      * The step sizes have to add up
 
247
                                      * to the dimensions of the hyper
 
248
                                      * rectangle specified by the
 
249
                                      * points @p p1 and @p p2.
 
250
                                      */
 
251
    template <int dim>
 
252
    static
 
253
    void
 
254
    subdivided_hyper_rectangle(Triangulation<dim>                      &tria,
 
255
                               const std::vector<std::vector<double> > &step_sizes,
 
256
                               const Point<dim>                        &p_1,
 
257
                               const Point<dim>                        &p_2,
 
258
                               const bool                              colorize);
 
259
 
 
260
                                     /**
 
261
                                      * Like the previous function, but with
 
262
                                      * the following twist: the @p
 
263
                                      * material_id argument is a
 
264
                                      * dim-dimensional array that, for each
 
265
                                      * cell, indicates which material_id
 
266
                                      * should be set. In addition, and this
 
267
                                      * is the major new functionality, if the
 
268
                                      * material_id of a cell is <tt>(unsigned
 
269
                                      * char)(-1)</tt>, then that cell is
 
270
                                      * deleted from the triangulation,
 
271
                                      * i.e. the domain will have a void
 
272
                                      * there.
 
273
                                      */
 
274
    template <int dim>
 
275
    static
 
276
    void
 
277
    subdivided_hyper_rectangle (Triangulation<dim>                       &tria,
 
278
                                const std::vector< std::vector<double> > &spacing,
 
279
                                const Point<dim>                         &p,
 
280
                                const Table<dim,unsigned char>           &material_id,
 
281
                                const bool                               colorize=false);
 
282
    
 
283
                                     /**
 
284
                                      * A parallelogram. The first
 
285
                                      * corner point is the
 
286
                                      * origin. The <tt>dim</tt>
 
287
                                      * adjacent points are the
 
288
                                      * one-dimensional subtensors of
 
289
                                      * the tensor provided and
 
290
                                      * additional points will be sums
 
291
                                      * of these two vectors.
 
292
                                      * Colorizing is done according
 
293
                                      * to hyper_rectangle().
 
294
                                      *
 
295
                                      * @note This function is
 
296
                                      * implemented in 2d only.
 
297
                                      *
 
298
                                      * @note The triangulation needs to be
 
299
                                      * void upon calling this
 
300
                                      * function.
 
301
                                      */
 
302
    template <int dim>
 
303
    static void
 
304
    parallelogram(Triangulation<dim>&  tria,
 
305
                  const Tensor<2,dim>& corners,
 
306
                  const bool           colorize=false);
 
307
                  
 
308
                              
 
309
                                     /**
 
310
                                      * Hypercube with a layer of
 
311
                                      * hypercubes around it. The
 
312
                                      * first two parameters give the
 
313
                                      * lower and upper bound of the
 
314
                                      * inner hypercube in all
 
315
                                      * coordinate directions.
 
316
                                      * @p thickness marks the size of
 
317
                                      * the layer cells.
 
318
                                      *
 
319
                                      * If the flag colorize is set,
 
320
                                      * the outer cells get material
 
321
                                      * id's according to the
 
322
                                      * following scheme: extending
 
323
                                      * over the inner cube in
 
324
                                      * (+/-) x-direction: 1/2. In y-direction
 
325
                                      * 4/8, in z-direction 16/32. The cells
 
326
                                      * at corners and edges (3d) get
 
327
                                      * these values bitwise or'd.
 
328
                                      *
 
329
                                      * Presently only available in 2d
 
330
                                      * and 3d.
 
331
                                      *
 
332
                                      * @note The triangulation needs to be
 
333
                                      * void upon calling this
 
334
                                      * function.
 
335
                                      */
 
336
    template <int dim>
 
337
    static void enclosed_hyper_cube (Triangulation<dim> &tria,
 
338
                                     const double      left = 0.,
 
339
                                     const double      right= 1.,
 
340
                                     const double      thickness = 1.,
 
341
                                     const bool        colorize = false);
 
342
    
 
343
                                     /**
 
344
                                      * Initialize the given
 
345
                                      * triangulation with a
 
346
                                      * hyperball, i.e. a circle or a
 
347
                                      * ball around <tt>center</tt>
 
348
                                      * with given <tt>radius</tt>.
 
349
                                      *
 
350
                                      * In order to avoid degenerate
 
351
                                      * cells at the boundaries, the
 
352
                                      * circle is triangulated by five
 
353
                                      * cells, the ball by seven
 
354
                                      * cells. The diameter of the
 
355
                                      * center cell is chosen so that
 
356
                                      * the aspect ratio of the
 
357
                                      * boundary cells after one
 
358
                                      * refinement is optimized.
 
359
                                      *
 
360
                                      * This function is declared to
 
361
                                      * exist for triangulations of
 
362
                                      * all space dimensions, but
 
363
                                      * throws an error if called in
 
364
                                      * 1d.
 
365
                                      *
 
366
                                      * @note The triangulation needs to be
 
367
                                      * void upon calling this
 
368
                                      * function.
 
369
                                      */    
 
370
    template <int dim>
 
371
    static void hyper_ball (Triangulation<dim> &tria,
 
372
                            const Point<dim>   &center = Point<dim>(),
 
373
                            const double      radius = 1.);
 
374
 
 
375
                                     /**
 
376
                                      * This class produces a half
 
377
                                      * hyper-ball around
 
378
                                      * <tt>center</tt>, which
 
379
                                      * contains four elements in 2d
 
380
                                      * and 6 in 3d. The cut plane is
 
381
                                      * perpendicular to the
 
382
                                      * <i>x</i>-axis.
 
383
                                      *
 
384
                                      * The boundary indicators for the final 
 
385
                                      * triangulation are 0 for the curved boundary and
 
386
                                      * 1 for the cut plane. 
 
387
                                      *
 
388
                                      * The appropriate
 
389
                                      * boundary class is 
 
390
                                      * HalfHyperBallBoundary, or HyperBallBoundary.
 
391
                                      * 
 
392
                                      * @note The triangulation needs to be
 
393
                                      * void upon calling this
 
394
                                      * function.
 
395
                                      */
 
396
    template <int dim>
 
397
    static void half_hyper_ball (Triangulation<dim> &tria,
 
398
                                 const Point<dim>   &center = Point<dim>(),
 
399
                                 const double      radius = 1.);
 
400
 
 
401
                                     /**
 
402
                                      * Create a cylinder around the
 
403
                                      * x-axis.  The cylinder extends
 
404
                                      * from <tt>x=-half_length</tt> to
 
405
                                      * <tt>x=+half_length</tt> and its
 
406
                                      * projection into the
 
407
                                      * @p yz-plane is a circle of
 
408
                                      * radius @p radius.
 
409
                                      *
 
410
                                      * In two dimensions, the
 
411
                                      * cylinder is a rectangle from
 
412
                                      * <tt>x=-half_length</tt> to
 
413
                                      * <tt>x=+half_length</tt> and
 
414
                                      * from <tt>y=-radius</tt> to
 
415
                                      * <tt>y=radius</tt>.
 
416
                                      *
 
417
                                      * The boundaries are colored
 
418
                                      * according to the following
 
419
                                      * scheme: 0 for the hull of the
 
420
                                      * cylinder, 1 for the left hand
 
421
                                      * face and 2 for the right hand
 
422
                                      * face.
 
423
                                      *
 
424
                                      * @note The triangulation needs to be
 
425
                                      * void upon calling this
 
426
                                      * function.
 
427
                                      */
 
428
    template <int dim>
 
429
    static void cylinder (Triangulation<dim> &tria,
 
430
                          const double      radius = 1.,
 
431
                          const double      half_length = 1.);
 
432
    
 
433
                                     /**
 
434
                                      * Initialize the given
 
435
                                      * triangulation with a hyper-L
 
436
                                      * consisting of exactly
 
437
                                      * <tt>2^dim-1</tt> cells. It
 
438
                                      * produces the hypercube with
 
439
                                      * the interval [<i>left,right</i>] without
 
440
                                      * the hypercube made out of the
 
441
                                      * interval [<i>(a+b)/2,b</i>].
 
442
                                      *
 
443
                                      * @image html hyper_l.png
 
444
                                      *
 
445
                                      * The triangulation needs to be
 
446
                                      * void upon calling this
 
447
                                      * function.
 
448
                                      *
 
449
                                      * This function is declared to
 
450
                                      * exist for triangulations of
 
451
                                      * all space dimensions, but
 
452
                                      * throws an error if called in
 
453
                                      * 1d.
 
454
                                      *
 
455
                                      * @note The triangulation needs to be
 
456
                                      * void upon calling this
 
457
                                      * function.
 
458
                                      */
 
459
    template <int dim>
 
460
    static void hyper_L (Triangulation<dim> &tria,
 
461
                         const double      left = -1.,
 
462
                         const double      right= 1.);
 
463
    
 
464
                                     /**
 
465
                                      * Initialize the given
 
466
                                      * Triangulation with a hypercube
 
467
                                      * with a slit. In each
 
468
                                      * coordinate direction, the
 
469
                                      * hypercube extends from @p left
 
470
                                      * to @p right.
 
471
                                      *
 
472
                                      * In 2d, the split goes in
 
473
                                      * vertical direction from
 
474
                                      * <tt>x=(left+right)/2,
 
475
                                      * y=left</tt> to the center of
 
476
                                      * the square at
 
477
                                      * <tt>x=y=(left+right)/2</tt>.
 
478
                                      *
 
479
                                      * In 3d, the 2d domain is just
 
480
                                      * extended in the
 
481
                                      * <i>z</i>-direction, such that
 
482
                                      * a plane cuts the lower half of
 
483
                                      * a rectangle in two.
 
484
                                      
 
485
                                      * This function is declared to
 
486
                                      * exist for triangulations of
 
487
                                      * all space dimensions, but
 
488
                                      * throws an error if called in
 
489
                                      * 1d.
 
490
                                      *
 
491
                                      * @note The triangulation needs to be
 
492
                                      * void upon calling this
 
493
                                      * function.
 
494
                                      */
 
495
    template <int dim>
 
496
    static void hyper_cube_slit (Triangulation<dim> &tria,
 
497
                                 const double      left = 0.,
 
498
                                 const double      right= 1.,
 
499
                                 const bool colorize = false);
 
500
    
 
501
                                     /**
 
502
                                      * Produce a hyper-shell,
 
503
                                      * the region between two
 
504
                                      * spheres around <tt>center</tt>,
 
505
                                      * with given
 
506
                                      * <tt>inner_radius</tt> and
 
507
                                      * <tt>outer_radius</tt>.
 
508
                                      *
 
509
                                      * If the flag @p colorize is @p
 
510
                                      * true, then the outer boundary
 
511
                                      * will have the id 1, while the
 
512
                                      * inner boundary has id zero. If
 
513
                                      * the flag is @p false, both
 
514
                                      * have id zero.
 
515
                                      *
 
516
                                      * In 2D, the number
 
517
                                      * <tt>n_cells</tt> of elements
 
518
                                      * for this initial triangulation
 
519
                                      * can be chosen arbitrarily. If
 
520
                                      * the number of initial cells is
 
521
                                      * zero (as is the default), then
 
522
                                      * it is computed adaptively such
 
523
                                      * that the resulting elements
 
524
                                      * have the least aspect ratio.
 
525
                                      *
 
526
                                      * In 3D, only two different
 
527
                                      * numbers are meaningful, 6 for
 
528
                                      * a surface based on a
 
529
                                      * hexahedron and 12 for the
 
530
                                      * rhombic dodecahedron.
 
531
                                      *
 
532
                                      * @image html hypershell3d-6.png
 
533
                                      * @image html hypershell3d-12.png
 
534
                                      *
 
535
                                      * This function is declared to
 
536
                                      * exist for triangulations of
 
537
                                      * all space dimensions, but
 
538
                                      * throws an error if called in
 
539
                                      * 1d. It is also currently not
 
540
                                      * implemented in 3d.
 
541
                                      *
 
542
                                      * @note The triangulation needs to be
 
543
                                      * void upon calling this
 
544
                                      * function.
 
545
                                      */
 
546
    template <int dim>
 
547
    static void hyper_shell (Triangulation<dim>   &tria,
 
548
                             const Point<dim>     &center,
 
549
                             const double        inner_radius,
 
550
                             const double        outer_radius,
 
551
                             const unsigned int  n_cells = 0,
 
552
                             bool colorize = false);
 
553
    
 
554
                                     /**
 
555
                                      * Produce a half hyper-shell,
 
556
                                      * i.e. the space between two
 
557
                                      * circles in two space
 
558
                                      * dimensions and the region
 
559
                                      * between two spheres in 3d,
 
560
                                      * with given inner and outer
 
561
                                      * radius and a given number of
 
562
                                      * elements for this initial
 
563
                                      * triangulation.  However,
 
564
                                      * opposed to the previous
 
565
                                      * function, it does not produce
 
566
                                      * a whole shell, but only one
 
567
                                      * half of it, namely that part
 
568
                                      * for which the first component
 
569
                                      * is restricted to non-negative
 
570
                                      * values. The purpose of this
 
571
                                      * class is to enable
 
572
                                      * computations for solutions
 
573
                                      * which have rotational
 
574
                                      * symmetry, in which case the
 
575
                                      * half shell in 2d represents a
 
576
                                      * shell in 3d.
 
577
                                      *
 
578
                                      * If the number of
 
579
                                      * initial cells is zero (as is
 
580
                                      * the default), then it is
 
581
                                      * computed adaptively such that
 
582
                                      * the resulting elements have
 
583
                                      * the least aspect ratio.
 
584
                                      *
 
585
                                      * At present, this function only
 
586
                                      * exists in 2d.
 
587
                                      *
 
588
                                      * @note The triangulation needs to be
 
589
                                      * void upon calling this
 
590
                                      * function.
 
591
                                      */
 
592
    template <int dim>
 
593
    static void half_hyper_shell (Triangulation<dim>   &tria,
 
594
                                  const Point<dim>     &center,
 
595
                                  const double        inner_radius,
 
596
                                  const double        outer_radius,
 
597
                                  const unsigned int  n_cells = 0);
 
598
    
 
599
                                     /**
 
600
                                      * Produce a domain that is the space
 
601
                                      * between two cylinders in 3d, with
 
602
                                      * given length, inner and outer radius
 
603
                                      * and a given number of elements for
 
604
                                      * this initial triangulation. If @p
 
605
                                      * n_radial_cells is zero (as is the
 
606
                                      * default), then it is computed
 
607
                                      * adaptively such that the resulting
 
608
                                      * elements have the least aspect
 
609
                                      * ratio. The same holds for @p
 
610
                                      * n_axial_cells.
 
611
                                      *
 
612
                                      * @note Although this function
 
613
                                      * is declared as a template, it
 
614
                                      * does not make sense in 1D and
 
615
                                      * 2D.
 
616
                                      *
 
617
                                      * @note The triangulation needs
 
618
                                      * to be void upon calling this
 
619
                                      * function.
 
620
                                      */
 
621
    template <int dim>
 
622
    static void cylinder_shell (Triangulation<dim>   &tria,
 
623
                                const double        length,
 
624
                                const double        inner_radius,
 
625
                                const double        outer_radius,
 
626
                                const unsigned int  n_radial_cells = 0,
 
627
                                const unsigned int  n_axial_cells = 0);
 
628
 
 
629
                                     /** 
 
630
                                      * This class produces a square
 
631
                                      * on the <i>xy</i>-plane with a
 
632
                                      * circular hole in the middle,
 
633
                                      * times the interval [0.L]
 
634
                                      * (only in 3d).
 
635
                                      *
 
636
                                      *  @image html cubes_hole.png
 
637
                                      *  
 
638
                                      * It is implemented in 2d and
 
639
                                      * 3d, and takes the following
 
640
                                      * arguments:
 
641
                                      * 
 
642
                                      * @arg @p inner_radius: size of the
 
643
                                      *    internal hole 
 
644
                                      * @arg @p  outer_radius: size of the
 
645
                                      *    biggest enclosed cylinder
 
646
                                      * @arg @p L: extension on the @p z-direction
 
647
                                      * @arg @p repetitions: number of subdivisions
 
648
                                      *      along the @p z-direction
 
649
                                      * @arg @p colorize: wether to assign different 
 
650
                                      *     boundary indicators to different faces.
 
651
                                      *    The colors are given in lexicographic 
 
652
                                      *    ordering for the flat faces (0 to 3 in 2d, 
 
653
                                      *    0 to 5 in 3d) plus the curved hole 
 
654
                                      *    (4 in 2d, and 6 in 3d).
 
655
                                      *    If @p colorize is set to false, then flat faces 
 
656
                                      *    get the number 0 and the hole gets number 1.
 
657
                                      */
 
658
    template<int dim>
 
659
    static void hyper_cube_with_cylindrical_hole (Triangulation<dim> &triangulation, 
 
660
                                                const double inner_radius = .25,
 
661
                                                const double outer_radius = .5,
 
662
                                                const double L = .5,
 
663
                                                const unsigned int repetition = 1,
 
664
                                                const bool colorize = false);
 
665
 
 
666
                                     /**
 
667
                                      * Produce a ring of cells in 3D that is
 
668
                                      * cut open, twisted and glued together
 
669
                                      * again. This results in a kind of
 
670
                                      * moebius-loop.
 
671
                                      *
 
672
                                      * @param tria        The triangulation to be worked on.
 
673
                                      * @param n_cells     The number of cells in the loop. Must be greater than 4.
 
674
                                      * @param n_rotations The number of rotations (Pi/2 each) to be performed before glueing the loop together.
 
675
                                      * @param R           The radius of the circle, which forms the middle line of the torus containing the loop of cells. Must be greater than @p r.
 
676
                                      * @param r           The radius of the cylinder bend together as loop.
 
677
                                      */
 
678
    static void moebius (Triangulation<3,3>&  tria,
 
679
                         const unsigned int   n_cells,
 
680
                         const unsigned int   n_rotations,
 
681
                         const double         R,
 
682
                         const double         r);
 
683
    
 
684
                                     /**
 
685
                                      * This function transformes the
 
686
                                      * @p Triangulation @p tria
 
687
                                      * smoothly to a domain that is
 
688
                                      * described by the boundary
 
689
                                      * points in the map
 
690
                                      * @p new_points. This map maps
 
691
                                      * the point indices to the
 
692
                                      * boundary points in the
 
693
                                      * transformed domain.
 
694
                                      *
 
695
                                      * Note, that the
 
696
                                      * @p Triangulation is changed
 
697
                                      * in-place, therefore you don't
 
698
                                      * need to keep two
 
699
                                      * triangulations, but the given
 
700
                                      * triangulation is changed
 
701
                                      * (overwritten).
 
702
                                      *
 
703
                                      * In 1d, this function is not
 
704
                                      * currently implemented.
 
705
                                      */
 
706
    template <int dim>
 
707
    static void laplace_transformation (Triangulation<dim> &tria,
 
708
                                        const std::map<unsigned int,Point<dim> > &new_points);
 
709
 
 
710
                                     /**
 
711
                                      * Exception
 
712
                                      */
 
713
    DeclException0 (ExcInvalidRadii);
 
714
                                     /**
 
715
                                      * Exception
 
716
                                      */
 
717
    DeclException1 (ExcInvalidRepetitions,
 
718
                    int,
 
719
                    << "The number of repetitions " << arg1
 
720
                    << " must be >=1.");
 
721
                                     /**
 
722
                                      * Exception
 
723
                                      */
 
724
    DeclException1 (ExcInvalidRepetitionsDimension,
 
725
                    int,
 
726
                    << "The vector of repetitions  must have " 
 
727
                    << arg1 <<" elements.");
 
728
 
 
729
  private:
 
730
                                     /**
 
731
                                      * Perform the action specified
 
732
                                      * by the @p colorize flag of
 
733
                                      * the hyper_rectangle()
 
734
                                      * function of this class.
 
735
                                      */
 
736
    template <int dim, int spacedim>
 
737
    static void colorize_hyper_rectangle (Triangulation<dim,spacedim> &tria);
 
738
 
 
739
                                     /**
 
740
                                      * Perform the action specified
 
741
                                      * by the @p colorize flag of
 
742
                                      * the
 
743
                                      * subdivided_hyper_rectangle()
 
744
                                      * function of this class. This
 
745
                                      * function is singled out
 
746
                                      * because it is dimension
 
747
                                      * specific.
 
748
                                      */
 
749
    template <int dim>
 
750
    static
 
751
    void
 
752
    colorize_subdivided_hyper_rectangle (Triangulation<dim> &tria,
 
753
                                         const Point<dim>   &p1,
 
754
                                         const Point<dim>   &p2,
 
755
                                         const double        epsilon);
 
756
 
 
757
                                     /**
 
758
                                      * Assign boundary number zero to
 
759
                                      * the inner shell boundary and 1
 
760
                                      * to the outer.
 
761
                                      */
 
762
    template<int dim>
 
763
    static
 
764
    void
 
765
    colorize_hyper_shell (Triangulation<dim>& tria,
 
766
                          const Point<dim>& center,
 
767
                          const double inner_radius,
 
768
                          const double outer_radius);
 
769
    
 
770
                                     /**
 
771
                                      * Solve the Laplace equation for
 
772
                                      * @p laplace_transformation
 
773
                                      * function for one of the
 
774
                                      * @p dim space
 
775
                                      * dimensions. Externalized into
 
776
                                      * a function of its own in order
 
777
                                      * to allow parallel execution.
 
778
                                      */
 
779
    static
 
780
    void
 
781
    laplace_solve (const SparseMatrix<double>          &S,
 
782
                   const std::map<unsigned int,double> &m,
 
783
                   Vector<double>                      &u);
 
784
};
 
785
 
 
786
 
 
787
 
 
788
 
 
789
DEAL_II_NAMESPACE_CLOSE
 
790
 
 
791
#endif