refineable_line_spectral_element.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: 1182 $
7 //LIC//
8 //LIC// $LastChangedDate: 2016-05-20 16:50:20 +0100 (Fri, 20 May 2016) $
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_REFINEABLE_LINE_SPECTRAL_ELEMENT_HEADER
31 #define OOMPH_REFINEABLE_LINE_SPECTRAL_ELEMENT_HEADER
32 
33 
34 // Config header generated by autoconfig
35 #ifdef HAVE_CONFIG_H
36  #include <oomph-lib-config.h>
37 #endif
38 
39 
40 // oomph-lib headers
41 #include "mesh.h"
42 #include "algebraic_elements.h"
44 #include "Qspectral_elements.h"
46 
47 namespace oomph
48 {
49 
50  //==========================================================================
51  /// Refineable version of LineElements that add functionality for spectral
52  /// Elements.
53  //==========================================================================
54  template<>
55  class RefineableQSpectralElement<1> : public virtual RefineableQElement<1>
56  {
57 
58  public:
59 
60  /// Constructor
62  {
63 #ifdef LEAK_CHECK
64  LeakCheckNames::RefineableQSpectralElement<1>_build+=1;
65 #endif
66  }
67 
68  /// Broken copy constructor
70  {
71  BrokenCopy::broken_copy("RefineableQSpectralElement<1>");
72  }
73 
74  /// Broken assignment operator
75 //Commented out broken assignment operator because this can lead to a conflict warning
76 //when used in the virtual inheritence hierarchy. Essentially the compiler doesn't
77 //realise that two separate implementations of the broken function are the same and so,
78 //quite rightly, it shouts.
79  /*void operator=(const RefineableQSpectralElement<1>&)
80  {
81  BrokenCopy::broken_assign("RefineableQSpecralElement<1>");
82  }*/
83 
84  /// Destructor
86  {
87 #ifdef LEAK_CHECK
88  LeakCheckNames::RefineableQSpectralElement<1>_build-=1;
89 #endif
90  }
91 
92  /// The only thing to add is rebuild from sons
93  void rebuild_from_sons(Mesh* &mesh_pt)
94  {
95  // The timestepper should be the same for all nodes and node 0 should
96  // never be deleted.
97  if(this->node_pt(0)==0)
98  {
99  throw OomphLibError("The vertex node (0) does not exist",
100  OOMPH_CURRENT_FUNCTION,
101  OOMPH_EXCEPTION_LOCATION);
102  }
103 
104  TimeStepper* time_stepper_pt = this->node_pt(0)->time_stepper_pt();
105 
106  // Determine number of history values stored
107  const unsigned ntstorage = time_stepper_pt->ntstorage();
108 
109  // Allocate storage for local coordinates
110  Vector<double> s_fraction(1), s(1);
111 
112  // Determine the number of nodes in the element
113  const unsigned n_node = this->nnode_1d();
114 
115  // Loop over the nodes in the element
116  for(unsigned n=0;n<n_node;n++)
117  {
118  // Get the fractional position of the node in the direction of s[0]
119  s_fraction[0] = this->local_one_d_fraction_of_node(n,0);
120 
121  // Determine the local coordinate in the father element
122  s[0] = -1.0 + 2.0*s_fraction[0];
123 
124  // If the node has not been built
125  if(this->node_pt(n)==0)
126  {
127  // Has the node been created by one of its neighbours?
128  bool is_periodic = false;
129  Node* created_node_pt =
130  this->node_created_by_neighbour(s_fraction,is_periodic);
131 
132  // If it has, set the pointer
133  if(created_node_pt!=0)
134  {
135  // If the node is periodic
136  if(is_periodic)
137  {
138  throw OomphLibError(
139  "Cannot handle periodic nodes in refineable spectral elements",
140  OOMPH_CURRENT_FUNCTION,
141  OOMPH_EXCEPTION_LOCATION);
142  }
143  // Non-periodic case, just set the pointer
144  else
145  {
146  this->node_pt(n) = created_node_pt;
147  }
148  }
149  // Otherwise, we need to build it
150  else
151  {
152  // First, find the son element in which the node should live
153 
154  // Find coordinate in the son
155  Vector<double> s_in_son(1);
156  using namespace BinaryTreeNames;
157  int son=-10;
158  // If s_fraction is between 0 and 0.5, we are in the left son
159  if(s_fraction[0] < 0.5)
160  {
161  son = L;
162  s_in_son[0] = -1.0 + 4.0*s_fraction[0];
163  }
164  // Otherwise we are in the right son
165  else
166  {
167  son = R;
168  s_in_son[0] = -1.0 + 4.0*(s_fraction[0]-0.5);
169  }
170 
171  // Get the pointer to the son element in which the new node
172  // would live
173  RefineableQSpectralElement<1>* son_el_pt =
174  dynamic_cast<RefineableQSpectralElement<1>*>(
175  this->tree_pt()->son_pt(son)->object_pt());
176 
177  // In 1D we should never be rebuilding an element's vertex nodes
178  // (since they will never be deleted), so throw an error if we
179  // appear to be doing so
180 #ifdef PARANOID
181  if(n==0 || n==n_node-1)
182  {
183  std::string error_message =
184  "I am trying to rebuild one of the (two) vertex nodes in\n";
185  error_message +=
186  "this 1D element. It should not have been possible to delete\n";
187  error_message +=
188  "either of these!\n";
189 
190  throw OomphLibError(
191  error_message,
192  OOMPH_CURRENT_FUNCTION,
193  OOMPH_EXCEPTION_LOCATION);
194  }
195 #endif
196 
197  // With this in mind we will always be creating normal "bulk" nodes
198  this->node_pt(n) = this->construct_node(n,time_stepper_pt);
199 
200  // Now we set the position and values at the newly created node
201 
202  // In the first instance use macro element or FE representation
203  // to create past and present nodal positions.
204  // (THIS STEP SHOULD NOT BE SKIPPED FOR ALGEBRAIC ELEMENTS AS NOT
205  // ALL OF THEM NECESSARILY IMPLEMENT NONTRIVIAL NODE UPDATE
206  // FUNCTIONS. CALLING THE NODE UPDATE FOR SUCH ELEMENTS/NODES WILL
207  // LEAVE THEIR NODAL POSITIONS WHERE THEY WERE (THIS IS APPROPRIATE
208  // ONCE THEY HAVE BEEN GIVEN POSITIONS) BUT WILL NOT ASSIGN SENSIBLE
209  // INITIAL POSITIONS!)
210 
211  // Loop over history values
212  for(unsigned t=0;t<ntstorage;t++)
213  {
214  // Allocate storage for the previous position of the node
215  Vector<double> x_prev(1);
216 
217  // Get position from son element -- this uses the macro element
218  // representation if appropriate
219  son_el_pt->get_x(t,s_in_son,x_prev);
220 
221  // Set the previous position of the new node
222  this->node_pt(n)->x(t,0) = x_prev[0];
223 
224  // Allocate storage for the previous values at the node
225  // NOTE: the size of this vector is equal to the number of values
226  // (e.g. 3 velocity components and 1 pressure, say)
227  // associated with each node and NOT the number of history values
228  // which the node stores!
229  Vector<double> prev_values;
230 
231  // Get values from son element
232  // Note: get_interpolated_values() sets Vector size itself.
233  son_el_pt->get_interpolated_values(t,s_in_son,prev_values);
234 
235  // Determine the number of values at the new node
236  const unsigned n_value = this->node_pt(n)->nvalue();
237 
238  // Loop over all values and set the previous values
239  for(unsigned v=0;v<n_value;v++)
240  {
241  this->node_pt(n)->set_value(t,v,prev_values[v]);
242  }
243  } // End of loop over history values
244 
245  // Add new node to mesh
246  mesh_pt->add_node_pt(this->node_pt(n));
247 
248  } // End of case where we build the node ourselves
249 
250  // Check if the element is an algebraic element
251  AlgebraicElementBase* alg_el_pt =
252  dynamic_cast<AlgebraicElementBase*>(this);
253 
254  // If so, throw error
255  if(alg_el_pt!=0)
256  {
257  std::string error_message =
258  "Have not implemented rebuilding from sons for";
259  error_message +=
260  "Algebraic Spectral elements yet\n";
261 
262  throw
263  OomphLibError(error_message,
264  "RefineableQSpectralElement<1>::rebuild_from_sons()",
265  OOMPH_EXCEPTION_LOCATION);
266  }
267 
268  } // End of if this node has not been built
269  } // End of loop over nodes in element
270  }
271 
272 
273  /// Overload the nodes built function
274  virtual bool nodes_built()
275  {
276  const unsigned n_node = this->nnode();
277  for(unsigned n=0;n<n_node;n++)
278  {
279  if(node_pt(n)==0) { return false; }
280  }
281  // If we get to here, OK
282  return true;
283  }
284 
285  };
286 
287 } // End of namespace
288 
289 #endif
RefineableQSpectralElement(const RefineableQSpectralElement< 1 > &dummy)
Broken copy constructor.
void broken_copy(const std::string &class_name)
Issue error message and terminate execution.
char t
Definition: cfortran.h:572
Nodes are derived from Data, but, in addition, have a definite (Eulerian) position in a space of a gi...
Definition: nodes.h:852
unsigned ntstorage() const
Return the number of doubles required to represent history (one for steady)
Definition: timesteppers.h:562
void add_node_pt(Node *const &node_pt)
Add a (pointer to a) node to the mesh.
Definition: mesh.h:602
virtual bool nodes_built()
Overload the nodes built function.
static char t char * s
Definition: cfortran.h:572
void rebuild_from_sons(Mesh *&mesh_pt)
The only thing to add is rebuild from sons.
std::string string(const unsigned &i)
Return the i-th string or "" if the relevant string hasn't been defined.
long RefineableQElement< 2 > _build
Base class for time-stepping schemes. Timestepper provides an approximation of the temporal derivativ...
Definition: timesteppers.h:219
A general mesh class.
Definition: mesh.h:74
virtual ~RefineableQSpectralElement()
Broken assignment operator.