annular_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// Version 1.0; svn revision $LastChangedRevision: 1097 $
7 //LIC//
8 //LIC// $LastChangedDate: 2015-12-17 11:53:17 +0000 (Thu, 17 Dec 2015) $
9 //LIC//
10 //LIC// Copyright (C) 2006-2016 Matthias Heil and Andrew Hazel
11 //LIC//
12 //LIC// This library is free software; you can redistribute it and/or
13 //LIC// modify it under the terms of the GNU Lesser General Public
14 //LIC// License as published by the Free Software Foundation; either
15 //LIC// version 2.1 of the License, or (at your option) any later version.
16 //LIC//
17 //LIC// This library is distributed in the hope that it will be useful,
18 //LIC// but WITHOUT ANY WARRANTY; without even the implied warranty of
19 //LIC// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
20 //LIC// Lesser General Public License for more details.
21 //LIC//
22 //LIC// You should have received a copy of the GNU Lesser General Public
23 //LIC// License along with this library; if not, write to the Free Software
24 //LIC// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
25 //LIC// 02110-1301 USA.
26 //LIC//
27 //LIC// The authors may be contacted at oomph-lib@maths.man.ac.uk.
28 //LIC//
29 //LIC//====================================================================
30 #ifndef OOMPH_ANNULAR_DOMAIN_HEADER
31 #define OOMPH_ANNULAR_DOMAIN_HEADER
32 
33 // Generic oomph-lib includes
34 #include "../generic/quadtree.h"
35 #include "../generic/domain.h"
36 #include "../generic/geom_objects.h"
37 
38 namespace oomph
39 {
40 
41 //=================================================================
42 /// \short Annular domain
43 //=================================================================
44 class AnnularDomain : public Domain
45 {
46 
47 public:
48 
49  /// \short Constructor: Specify azimuthal fraction (1.0 is 360 degrees)
50  /// number of macro elements in azimuthal and radial direction,
51  /// inner radius and thickness. Rotate mesh by angle phi.
52  AnnularDomain(const double& azimuthal_fraction,
53  const unsigned &ntheta, const unsigned &nr,
54  const double& a, const double &h, const double& phi) :
55  Azimuthal_fraction(azimuthal_fraction),Inner_radius(a),
56  Thickness(h), Ntheta(ntheta), Nr(nr), Phi(phi)
57  {
58  const unsigned n_macro=ntheta*nr;
59  Macro_element_pt.resize(n_macro);
60 
61  // Create the macro elements
62  for (unsigned i=0;i<n_macro;i++)
63  {
65  }
66  }
67 
68  /// Broken copy constructor
70  {
71  BrokenCopy::broken_copy("AnnularDomain");
72  }
73 
74  /// Broken assignment operator
75  void operator=(const AnnularDomain&)
76  {
77  BrokenCopy::broken_assign("AnnularDomain");
78  }
79 
80 
81  /// Destructor: Kill all macro elements
83  {
84  unsigned n=nmacro_element();
85  for(unsigned i=0;i<n;i++)
86  {
87  delete Macro_element_pt[i];
88  }
89  }
90 
91 
92  /// \short Vector representation of the i_macro-th macro element
93  /// boundary i_direct (N/S/W/E) at time level t
94  /// (t=0: present; t>0: previous):
95  /// f(s).
96  void macro_element_boundary(const unsigned& t,
97  const unsigned& i_macro,
98  const unsigned& i_direct,
99  const Vector<double>& s,
100  Vector<double>& f);
101 
102 private:
103 
104  /// Azimuthal fraction
106 
107  /// Inner radius
108  double Inner_radius;
109 
110  /// Thickness
111  double Thickness;
112 
113  /// Number of macro elements in azimuthal direction
114  unsigned Ntheta;
115 
116  /// Number of macro elements in radial direction
117  unsigned Nr;
118 
119  /// Rotation angle
120  double Phi;
121 
122 };
123 
124 
125 /////////////////////////////////////////////////////////////////////////
126 /////////////////////////////////////////////////////////////////////////
127 /////////////////////////////////////////////////////////////////////////
128 
129 
130 
131 //=================================================================
132 /// \short Vector representation of the imacro-th macro element
133 /// boundary idirect (N/S/W/E) at time level t
134 /// (t=0: present; t>0: previous): f(s)
135 //=================================================================
137  const unsigned& t,
138  const unsigned& imacro,
139  const unsigned& idirect,
140  const Vector<double>& s,
141  Vector<double>& f)
142  {
143 
144  using namespace QuadTreeNames;
145 
146  // Get coordinates of macro element
147  unsigned i_theta=imacro%Ntheta;
148  unsigned i_r=(imacro-i_theta)/Ntheta;
149 
150  // Angle and radius limits
151  double theta_lo=Azimuthal_fraction*2.0*MathematicalConstants::Pi*
152  double(i_theta)/double(Ntheta);
153 
154  double theta_hi=Azimuthal_fraction*2.0*MathematicalConstants::Pi*
155  double(i_theta+1)/double(Ntheta);
156 
157  // Revert direction (convoluted -- don't ask. It mirrors what happens
158  // in the mesh...
159  theta_lo=-MathematicalConstants::Pi+
161  theta_hi=-MathematicalConstants::Pi+
163 
164  double r_lo=Inner_radius+Thickness*double(i_r)/double(Nr);
165  double r_hi=Inner_radius+Thickness*double(i_r+1)/double(Nr);
166 
167  // Actual radius and angle
168  double r=0.0;
169  double theta=0.0;
170 
171  // Which direction?
172  switch(idirect)
173  {
174  case N:
175 
176  theta=theta_lo+0.5*(s[0]+1.0)*(theta_hi-theta_lo);
177  r=r_hi;
178 
179  break;
180 
181  case S:
182 
183  theta=theta_lo+0.5*(s[0]+1.0)*(theta_hi-theta_lo);
184  r=r_lo;
185 
186  break;
187 
188  case W:
189 
190  theta=theta_lo;
191  r=r_lo+0.5*(s[0]+1.0)*(r_hi-r_lo);
192 
193  break;
194 
195  case E:
196 
197  theta=theta_hi;
198  r=r_lo+0.5*(s[0]+1.0)*(r_hi-r_lo);
199 
200  break;
201 
202  default:
203  std::ostringstream error_stream;
204  error_stream << "idirect is " << idirect
205  << " not one of N, S, W, E" << std::endl;
206 
207  throw OomphLibError(
208  error_stream.str(),
209  OOMPH_CURRENT_FUNCTION,
210  OOMPH_EXCEPTION_LOCATION);
211  }
212 
213  f[0]=r*cos(theta+Phi);
214  f[1]=r*sin(theta+Phi);
215 }
216 
217 }
218 
219 #endif
220 
unsigned nmacro_element()
Number of macro elements in domain.
Definition: domain.h:107
void broken_copy(const std::string &class_name)
Issue error message and terminate execution.
cstr elem_len * i
Definition: cfortran.h:607
const double Pi
50 digits from maple
Vector< MacroElement * > Macro_element_pt
Vector of pointers to macro elements.
Definition: domain.h:237
double Phi
Rotation angle.
double Inner_radius
Inner radius.
char t
Definition: cfortran.h:572
double Thickness
Thickness.
unsigned Nr
Number of macro elements in radial direction.
double Azimuthal_fraction
Azimuthal fraction.
~AnnularDomain()
Destructor: Kill all macro elements.
static char t char * s
Definition: cfortran.h:572
void operator=(const AnnularDomain &)
Broken assignment operator.
AnnularDomain(const double &azimuthal_fraction, const unsigned &ntheta, const unsigned &nr, const double &a, const double &h, const double &phi)
Constructor: Specify azimuthal fraction (1.0 is 360 degrees) number of macro elements in azimuthal an...
void broken_assign(const std::string &class_name)
Issue error message and terminate execution.
AnnularDomain(const AnnularDomain &)
Broken copy constructor.
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=...
Annular domain.
Base class for Domains with curvilinear and/or time-dependent boundaries. Domain boundaries are typic...
Definition: domain.h:71
unsigned Ntheta
Number of macro elements in azimuthal direction.