quarter_circle_sector_domain.h
Go to the documentation of this file.
1// LIC// ====================================================================
2// LIC// This file forms part of oomph-lib, the object-oriented,
3// LIC// multi-physics finite-element library, available
4// LIC// at http://www.oomph-lib.org.
5// LIC//
6// LIC// Copyright (C) 2006-2022 Matthias Heil and Andrew Hazel
7// LIC//
8// LIC// This library is free software; you can redistribute it and/or
9// LIC// modify it under the terms of the GNU Lesser General Public
10// LIC// License as published by the Free Software Foundation; either
11// LIC// version 2.1 of the License, or (at your option) any later version.
12// LIC//
13// LIC// This library is distributed in the hope that it will be useful,
14// LIC// but WITHOUT ANY WARRANTY; without even the implied warranty of
15// LIC// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16// LIC// Lesser General Public License for more details.
17// LIC//
18// LIC// You should have received a copy of the GNU Lesser General Public
19// LIC// License along with this library; if not, write to the Free Software
20// LIC// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
21// LIC// 02110-1301 USA.
22// LIC//
23// LIC// The authors may be contacted at oomph-lib@maths.man.ac.uk.
24// LIC//
25// LIC//====================================================================
26// Include guards
27#ifndef OOMPH_QUARTER_CIRCLE_SECTOR_DOMAIN_HEADER
28#define OOMPH_QUARTER_CIRCLE_SECTOR_DOMAIN_HEADER
29
30// Generic oomph-lib includes
31#include "../generic/quadtree.h"
32#include "../generic/domain.h"
33#include "../generic/geom_objects.h"
34
35namespace oomph
36{
37 //=================================================================
38 /// Circular sector as domain. Domain is bounded by
39 /// curved boundary which is represented by a GeomObject. Domain is
40 /// parametrised by three macro elements.
41 //=================================================================
43 {
44 public:
45 /// Constructor: Pass boundary object and start and end coordinates
46 /// and fraction along boundary object where outer ring is divided.
47 QuarterCircleSectorDomain(GeomObject* boundary_geom_object_pt,
48 const double& xi_lo,
49 const double& fract_mid,
50 const double& xi_hi)
51 : Xi_lo(xi_lo),
52 Fract_mid(fract_mid),
53 Xi_hi(xi_hi),
54 Wall_pt(boundary_geom_object_pt),
56 {
57 // There are three macro elements
58 unsigned nmacro = 3;
59
60 // Resize
61 Macro_element_pt.resize(nmacro);
62
63 // Create macro elements
64 for (unsigned i = 0; i < nmacro; i++)
65 {
67 }
68 }
69
70
71 /// Broken copy constructor
73
74 /// Broken assignment operator
75 void operator=(const QuarterCircleSectorDomain&) = delete;
76
77 /// Destructor: empty; cleanup done in base class
79
80 /// Typedef for function pointer for function that squashes
81 /// the outer two macro elements towards
82 /// the wall by mapping the input value of the "radial" macro element
83 /// coordinate to the return value
84 typedef double (*BLSquashFctPt)(const double& s);
85
86
87 /// Function pointer for function that squashes
88 /// the outer two macro elements towards
89 /// the wall by mapping the input value of the "radial" macro element
90 /// coordinate to the return value
92 {
93 return BL_squash_fct_pt;
94 }
95
96
97 /// Function that squashes the outer two macro elements towards
98 /// the wall by mapping the input value of the "radial" macro element
99 /// coordinate to the return value
100 double s_squashed(const double& s)
101 {
102 return BL_squash_fct_pt(s);
103 }
104
105 /// Vector representation of the i_macro-th macro element
106 /// boundary i_direct (N/S/W/E) at time level t
107 /// (t=0: present; t>0: previous):
108 /// f(s). Note that the local coordinate \b s is a 1D
109 /// Vector rather than a scalar -- this is unavoidable because
110 /// this function implements the pure virtual function in the
111 /// Domain base class.
112 void macro_element_boundary(const unsigned& t,
113 const unsigned& i_macro,
114 const unsigned& i_direct,
115 const Vector<double>& s,
116 Vector<double>& f);
117
118
119 private:
120 /// Lower limit for the (1D) coordinates along the wall
121 double Xi_lo;
122
123 /// Fraction along wall where outer ring is to be divided
124 double Fract_mid;
125
126 /// Upper limit for the (1D) coordinates along the wall
127 double Xi_hi;
128
129 /// Pointer to geometric object that represents the curved wall
131
132 /// Function pointer for function that squashes
133 /// the outer two macro elements towards
134 /// the wall by mapping the input value of the "radial" macro element
135 /// coordinate to the return value
137
138 /// Default for function that squashes
139 /// the outer two macro elements towards
140 /// the wall by mapping the input value of the "radial" macro element
141 /// coordinate to the return value: Identity.
142 static double default_BL_squash_fct(const double& s)
143 {
144 return s;
145 }
146
147 /// Boundary of top left macro element zeta \f$ \in [-1,1] \f$
148 void r_top_left_N(const unsigned& t,
149 const Vector<double>& zeta,
150 Vector<double>& f);
151
152 /// Boundary of top left macro element zeta \f$ \in [-1,1] \f$
153 void r_top_left_W(const unsigned& t,
154 const Vector<double>& zeta,
155 Vector<double>& f);
156
157 /// Boundary of top left macro element zeta \f$ \in [-1,1] \f$
158 void r_top_left_S(const unsigned& t,
159 const Vector<double>& zeta,
160 Vector<double>& f);
161
162 /// Boundary of top left macro element zeta \f$ \in [-1,1] \f$
163 void r_top_left_E(const unsigned& t,
164 const Vector<double>& zeta,
165 Vector<double>& f);
166
167 /// Boundary of bottom right macro element zeta \f$ \in [-1,1] \f$
168 void r_bot_right_N(const unsigned& t,
169 const Vector<double>& zeta,
170 Vector<double>& f);
171
172 /// Boundary of bottom right macro element zeta \f$ \in [-1,1] \f$
173 void r_bot_right_W(const unsigned& t,
174 const Vector<double>& zeta,
175 Vector<double>& f);
176
177 /// Boundary of bottom right macro element zeta \f$ \in [-1,1] \f$
178 void r_bot_right_S(const unsigned& t,
179 const Vector<double>& zeta,
180 Vector<double>& f);
181
182 /// Boundary of bottom right macro element zeta \f$ \in [-1,1] \f$
183 void r_bot_right_E(const unsigned& t,
184 const Vector<double>& zeta,
185 Vector<double>& f);
186
187 /// Boundary of central box macro element zeta \f$ \in [-1,1] \f$
188 void r_centr_N(const unsigned& t,
189 const Vector<double>& zeta,
190 Vector<double>& f);
191
192 /// Boundary of central box macro element zeta \f$ \in [-1,1] \f$
193 void r_centr_E(const unsigned& t,
194 const Vector<double>& zeta,
195 Vector<double>& f);
196
197 /// Boundary of central box macro element zeta \f$ \in [-1,1] \f$
198 void r_centr_S(const unsigned& t,
199 const Vector<double>& zeta,
200 Vector<double>& f);
201
202 /// Boundary of central box macro element zeta \f$ \in [-1,1] \f$
203 void r_centr_W(const unsigned& t,
204 const Vector<double>& zeta,
205 Vector<double>& f);
206 };
207
208
209 /// //////////////////////////////////////////////////////////////////////
210 /// //////////////////////////////////////////////////////////////////////
211 /// //////////////////////////////////////////////////////////////////////
212
213
214 //=================================================================
215 /// Vector representation of the imacro-th macro element
216 /// boundary idirect (N/S/W/E) at time level t (t=0: present; t>0: previous):
217 /// f(s)
218 //=================================================================
220 const unsigned& t,
221 const unsigned& imacro,
222 const unsigned& idirect,
223 const Vector<double>& s,
225 {
226#ifdef WARN_ABOUT_SUBTLY_CHANGED_OOMPH_INTERFACES
227 // Warn about time argument being moved to the front
229 "Order of function arguments has changed between versions 0.8 and 0.85",
230 "QuarterCircleSectorDomain::macro_element_boundary(...)",
231 OOMPH_EXCEPTION_LOCATION);
232#endif
233
234 // Which macro element?
235 // --------------------
236 switch (imacro)
237 {
238 using namespace QuadTreeNames;
239
240 // Macro element 0: Central box
241 case 0:
242
243 // Which direction?
244 if (idirect == N)
245 {
246 r_centr_N(t, s, f);
247 }
248 else if (idirect == S)
249 {
250 r_centr_S(t, s, f);
251 }
252 else if (idirect == W)
253 {
254 r_centr_W(t, s, f);
255 }
256 else if (idirect == E)
257 {
258 r_centr_E(t, s, f);
259 }
260 else
261 {
262 std::ostringstream error_stream;
263 error_stream << "idirect is " << idirect << " not one of N, S, E, W"
264 << std::endl;
265
266 throw OomphLibError(error_stream.str(),
267 OOMPH_CURRENT_FUNCTION,
268 OOMPH_EXCEPTION_LOCATION);
269 }
270
271 break;
272
273 // Macro element 1: Bottom right
274 case 1:
275
276 // Which direction?
277 if (idirect == N)
278 {
279 r_bot_right_N(t, s, f);
280 }
281 else if (idirect == S)
282 {
283 r_bot_right_S(t, s, f);
284 }
285 else if (idirect == W)
286 {
287 r_bot_right_W(t, s, f);
288 }
289 else if (idirect == E)
290 {
291 r_bot_right_E(t, s, f);
292 }
293 else
294 {
295 std::ostringstream error_stream;
296 error_stream << "idirect is " << idirect << " not one of N, S, E, W"
297 << std::endl;
298
299 throw OomphLibError(error_stream.str(),
300 OOMPH_CURRENT_FUNCTION,
301 OOMPH_EXCEPTION_LOCATION);
302 }
303
304 break;
305
306 // Macro element 2:Top left
307 case 2:
308
309 // Which direction?
310 if (idirect == N)
311 {
312 r_top_left_N(t, s, f);
313 }
314 else if (idirect == S)
315 {
316 r_top_left_S(t, s, f);
317 }
318 else if (idirect == W)
319 {
320 r_top_left_W(t, s, f);
321 }
322 else if (idirect == E)
323 {
324 r_top_left_E(t, s, f);
325 }
326 else
327 {
328 std::ostringstream error_stream;
329 error_stream << "idirect is " << idirect << " not one of N, S, E, W"
330 << std::endl;
331
332 throw OomphLibError(error_stream.str(),
333 OOMPH_CURRENT_FUNCTION,
334 OOMPH_EXCEPTION_LOCATION);
335 }
336
337 break;
338
339 default:
340
341 // Error
342 std::ostringstream error_stream;
343 error_stream << "Wrong imacro " << imacro << std::endl;
344
345 throw OomphLibError(
346 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
347 }
348 }
349
350
351 //=================================================================
352 /// Northern edge of top left macro element \f$ s \in [-1,1] \f$
353 //=================================================================
355 const Vector<double>& s,
357 {
358 Vector<double> x(1);
359
360 // Coordinate along wall
361 x[0] = Xi_hi + (Fract_mid * (Xi_hi - Xi_lo) - Xi_hi) * 0.5 * (s[0] + 1.0);
362
363 Wall_pt->position(t, x, f);
364 }
365
366
367 //=================================================================
368 /// Western edge of top left macro element \f$s \in [-1,1] \f$
369 //=================================================================
371 const Vector<double>& s,
373 {
374 Vector<double> x(1);
375
376 // Top left corner
377 Vector<double> r_top(2);
378 x[0] = Xi_hi;
379
380 Wall_pt->position(t, x, r_top);
381
382 f[0] = 0.0;
383 f[1] = 0.5 * r_top[1] * (1.0 + s_squashed(0.5 * (s[0] + 1.0)));
384 }
385
386
387 //=================================================================
388 /// Southern edge of top left macro element \f$ s \in [-1,1] \f$
389 //=================================================================
391 const Vector<double>& s,
393 {
394 Vector<double> x(1);
395
396 // Top left corner
397 Vector<double> r_top(2);
398 x[0] = Xi_hi;
399
400 Wall_pt->position(t, x, r_top);
401
402
403 // Bottom right corner
404 Vector<double> r_bot(2);
405 x[0] = 0.0;
406
407 Wall_pt->position(t, x, r_bot);
408
409 f[0] = 0.5 * r_bot[0] * 0.5 * (s[0] + 1.0);
410 f[1] = 0.5 * r_top[1];
411 }
412
413
414 //=================================================================
415 /// Eastern edge of top left macro element \f$ s \in [-1,1] \f$
416 //=================================================================
418 const Vector<double>& s,
420 {
421 Vector<double> x(1);
422
423 // Top left corner
424 Vector<double> r_top(2);
425 x[0] = Xi_hi;
426
427 Wall_pt->position(t, x, r_top);
428
429 // Bottom right corner
430 Vector<double> r_bot(2);
431 x[0] = Xi_lo;
432
433 Wall_pt->position(t, x, r_bot);
434
435 // Halfway along wall
436 Vector<double> r_half(2);
437 x[0] = Xi_lo + Fract_mid * (Xi_hi - Xi_lo);
438
439 Wall_pt->position(t, x, r_half);
440
441 f[0] = 0.5 * (r_bot[0] + s_squashed(0.5 * (s[0] + 1.0)) *
442 (2.0 * r_half[0] - r_bot[0]));
443 f[1] = 0.5 * (r_top[1] + s_squashed(0.5 * (s[0] + 1.0)) *
444 (2.0 * r_half[1] - r_top[1]));
445 }
446
447
448 //=================================================================
449 /// Northern edge of bottom right macro element
450 //=================================================================
452 const Vector<double>& s,
454 {
455 r_top_left_E(t, s, f);
456 }
457
458 //=================================================================
459 /// Western edge of bottom right macro element
460 //=================================================================
462 const Vector<double>& s,
464 {
465 Vector<double> x(1);
466
467 // Top left corner
468 Vector<double> r_top(2);
469 x[0] = Xi_hi;
470
471 Wall_pt->position(t, x, r_top);
472
473 // Bottom right corner
474 Vector<double> r_bot(2);
475 x[0] = Xi_lo;
476
477 Wall_pt->position(t, x, r_bot);
478
479 f[0] = 0.5 * r_bot[0];
480 f[1] = 0.5 * r_top[1] * 0.5 * (s[0] + 1.0);
481 }
482
483 //=================================================================
484 /// Southern edge of bottom right macro element
485 //=================================================================
487 const Vector<double>& s,
489 {
490 Vector<double> x(1);
491
492 // Bottom right corner
493 Vector<double> r_bot(2);
494 x[0] = Xi_lo;
495 Wall_pt->position(t, x, r_bot);
496
497
498 f[0] = 0.5 * r_bot[0] * (1.0 + s_squashed(0.5 * (s[0] + 1.0)));
499 f[1] = 0.0;
500 }
501
502 //=================================================================
503 /// Eastern edge of bottom right macro element
504 //=================================================================
506 const Vector<double>& s,
508 {
509 Vector<double> x(1);
510
511 // Coordinate along wall
512 x[0] = Xi_lo + (Fract_mid * (Xi_hi - Xi_lo) - Xi_lo) * (s[0] + 1.0) * 0.5;
513
514 Wall_pt->position(t, x, f);
515 }
516
517
518 //=================================================================
519 /// Northern edge of central box
520 //=================================================================
522 const Vector<double>& s,
524 {
525 r_top_left_S(t, s, f);
526 }
527
528
529 //=================================================================
530 /// Eastern edge of central box
531 //=================================================================
533 const Vector<double>& s,
535 {
536 r_bot_right_W(t, s, f);
537 }
538
539
540 //=================================================================
541 /// Southern edge of central box
542 //=================================================================
544 const Vector<double>& s,
546 {
547 Vector<double> x(1);
548
549 // Bottom right corner
550 Vector<double> r_bot(2);
551 x[0] = Xi_lo;
552 Wall_pt->position(t, x, r_bot);
553
554 f[0] = 0.5 * r_bot[0] * 0.5 * (s[0] + 1.0);
555 f[1] = 0.0;
556 }
557
558
559 //=================================================================
560 /// Western edge of central box
561 //=================================================================
563 const Vector<double>& s,
565 {
566 Vector<double> x(1);
567
568 // Top left corner
569 Vector<double> r_top(2);
570 x[0] = Xi_hi;
571 Wall_pt->position(t, x, r_top);
572
573 f[0] = 0.0;
574 f[1] = 0.5 * r_top[1] * 0.5 * (s[0] + 1.0);
575 }
576
577
578} // namespace oomph
579
580#endif
static char t char * s
Definition: cfortran.h:568
cstr elem_len * i
Definition: cfortran.h:603
char t
Definition: cfortran.h:568
Base class for Domains with curvilinear and/or time-dependent boundaries. Domain boundaries are typic...
Definition: domain.h:67
Vector< MacroElement * > Macro_element_pt
Vector of pointers to macro elements.
Definition: domain.h:301
/////////////////////////////////////////////////////////////////////
Definition: geom_objects.h:101
virtual void position(const Vector< double > &zeta, Vector< double > &r) const =0
Parametrised position on object at current time: r(zeta).
An OomphLibError object which should be thrown when an run-time error is encountered....
An OomphLibWarning object which should be created as a temporary object to issue a warning....
QMacroElement specialised to 2 spatial dimensions.
Circular sector as domain. Domain is bounded by curved boundary which is represented by a GeomObject....
void r_top_left_S(const unsigned &t, const Vector< double > &zeta, Vector< double > &f)
Boundary of top left macro element zeta .
QuarterCircleSectorDomain(const QuarterCircleSectorDomain &)=delete
Broken copy constructor.
void r_centr_N(const unsigned &t, const Vector< double > &zeta, Vector< double > &f)
Boundary of central box macro element zeta .
void r_centr_W(const unsigned &t, const Vector< double > &zeta, Vector< double > &f)
Boundary of central box macro element zeta .
void r_bot_right_S(const unsigned &t, const Vector< double > &zeta, Vector< double > &f)
Boundary of bottom right macro element zeta .
GeomObject * Wall_pt
Pointer to geometric object that represents the curved wall.
double s_squashed(const double &s)
Function that squashes the outer two macro elements towards the wall by mapping the input value of th...
double(* BLSquashFctPt)(const double &s)
Typedef for function pointer for function that squashes the outer two macro elements towards the wall...
void operator=(const QuarterCircleSectorDomain &)=delete
Broken assignment operator.
void r_top_left_W(const unsigned &t, const Vector< double > &zeta, Vector< double > &f)
Boundary of top left macro element zeta .
void r_bot_right_W(const unsigned &t, const Vector< double > &zeta, Vector< double > &f)
Boundary of bottom right macro element zeta .
void r_top_left_E(const unsigned &t, const Vector< double > &zeta, Vector< double > &f)
Boundary of top left macro element zeta .
void macro_element_boundary(const unsigned &t, const unsigned &i_macro, const unsigned &i_direct, const Vector< double > &s, Vector< double > &f)
Vector representation of the i_macro-th macro element boundary i_direct (N/S/W/E) at time level t (t=...
BLSquashFctPt & bl_squash_fct_pt()
Function pointer for function that squashes the outer two macro elements towards the wall by mapping ...
void r_centr_S(const unsigned &t, const Vector< double > &zeta, Vector< double > &f)
Boundary of central box macro element zeta .
double Xi_hi
Upper limit for the (1D) coordinates along the wall.
void r_bot_right_E(const unsigned &t, const Vector< double > &zeta, Vector< double > &f)
Boundary of bottom right macro element zeta .
QuarterCircleSectorDomain(GeomObject *boundary_geom_object_pt, const double &xi_lo, const double &fract_mid, const double &xi_hi)
Constructor: Pass boundary object and start and end coordinates and fraction along boundary object wh...
void r_bot_right_N(const unsigned &t, const Vector< double > &zeta, Vector< double > &f)
Boundary of bottom right macro element zeta .
static double default_BL_squash_fct(const double &s)
Default for function that squashes the outer two macro elements towards the wall by mapping the input...
double Fract_mid
Fraction along wall where outer ring is to be divided.
void r_centr_E(const unsigned &t, const Vector< double > &zeta, Vector< double > &f)
Boundary of central box macro element zeta .
BLSquashFctPt BL_squash_fct_pt
Function pointer for function that squashes the outer two macro elements towards the wall by mapping ...
~QuarterCircleSectorDomain()
Destructor: empty; cleanup done in base class.
void r_top_left_N(const unsigned &t, const Vector< double > &zeta, Vector< double > &f)
Boundary of top left macro element zeta .
double Xi_lo
Lower limit for the (1D) coordinates along the wall.
//////////////////////////////////////////////////////////////////// ////////////////////////////////...