fsi.cc
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 #include "fsi.h"
31 #include "integral.h"
32 #include "algebraic_elements.h"
33 
34 namespace oomph
35 {
36 
37 //======================================================================
38 /// Namespace for "global" FSI functions
39 //======================================================================
40 namespace FSI_functions
41 {
42 
43  /// Strouhal number St = a/(UT) for application of no slip condition.
44  /// Initialised to 1.0.
46 
47  /// Apply no-slip condition for N.St. on a moving wall node
48  /// u = St dR/dt, where the Strouhal number St = a/(UT) is defined by
49  /// FSI_functions::Strouhal_for_no_slip and is initialised to 1.0.
50  /// Note: This requires the x,y,[z] velocity components to be stored
51  /// in nodal values 0,1,[2]. This is the default for all currently
52  /// existing Navier-Stokes elements. If you use any others,
53  /// use this function at your own risk.
55  {
56  // Find out spatial dimension of node
57  unsigned ndim=node_pt->ndim();
58  // Assign velocity
59  for (unsigned i=0;i<ndim;i++)
60  {
61  node_pt->set_value(i,Strouhal_for_no_slip*(node_pt->dposition_dt(i)));
62  }
63  }
64 
65 }
66 
67 
68 //====================================================================
69 /// Flag that allows the suppression of warning messages
70 //====================================================================
72 
73 //====================================================================
74 /// \short Function to describe the local dofs of the element. The ostream
75 /// specifies the output stream to which the description
76 /// is written; the string stores the currently
77 /// assembled output that is ultimately written to the
78 /// output stream by Data::describe_dofs(...); it is typically
79 /// built up incrementally as we descend through the
80 /// call hierarchy of this function when called from
81 /// Problem::describe_dofs(...)
82 //====================================================================
84 describe_local_dofs(std::ostream& out,const std::string& current_string) const
85 {
86  // Call the standard finite element classification.
87  FiniteElement::describe_local_dofs(out,current_string);
88  describe_solid_local_dofs(out,current_string);
89  //Find the number of external field data
90  const unsigned n_external_field_data = nexternal_interaction_field_data();
91  //Now loop over the field data again to assign local equation numbers
92  for(unsigned i=0;i<n_external_field_data;i++)
93  {
94  std::stringstream conversion;
95  conversion<<" of External Interaction Field Data "<<i<<current_string;
96  std::string in(conversion.str());
98  }
99 
100  //Find the number of external geometric data
101  unsigned n_external_geom_data = nexternal_interaction_geometric_data();
102 
103  //Now loop over the field data again assign local equation numbers
104  for(unsigned i=0;i<n_external_geom_data;i++)
105  {
106  std::stringstream conversion;
107  conversion<<" of External Interaction Geometric Data "<<i<<current_string;
108  std::string in(conversion.str());
110  }
111 }
112 
113  //================================================================
114  /// Compete build of element: Assign storage -- pass the Eulerian
115  /// dimension of the "adjacent" fluid elements and the
116  /// number of local coordinates required to parametrise
117  /// the wall element. E.g. for a FSIKirchhoffLoveBeam,
118  /// bounding a 2D fluid domain ndim_fluid=2 and nlagr_solid=1.
119  /// Note that, by default, we're only providing storage for fluid
120  /// loading on the "front" of the element. Call
121  /// FSIWallElement::enable_fluid_loading_on_both_sides()
122  /// to enable fluid loading on the back, too.
123  //====================================================================
124  void FSIWallElement::setup_fsi_wall_element(const unsigned& nlagr_solid,
125  const unsigned& ndim_fluid)
126  {
127 
128  // Set storage for underlying GeomObject
129  set_nlagrangian_and_ndim(nlagr_solid, ndim_fluid);
130 
131  // Set source element storage - one interaction
132  this->set_ninteraction(1);
133  }
134 
135 
136  //==================================================================
137  /// \short Allow element to be loaded by fluid on both
138  /// sides. (Resizes containers for lookup schemes and initialises
139  /// data associated with elements at the "back" of the FSIWallElement
140  /// to NULL.
141  //==================================================================
143  {
144  // Both faces are loaded
146 
147  // Set source element storage - two interactions
148  this->set_ninteraction(2);
149  }
150 
151 
152  //==================================================================
153  /// Get the contribution to the load vector provided by
154  /// the adjacent fluid element: Pass number of integration point
155  /// in solid element,
156  /// and the unit normal vector (pointing into the fluid on the "front"
157  /// of the FSIWallElement) and return the load vector.
158  /// Note that the load is non-dimensionalised on the wall-stress scale,
159  /// i.e. it is obtained by computing the traction (on the fluid stress-scale)
160  /// from the adjacent fluid element and then multiplying it by
161  /// the stress-scale-ratio \f$ Q. \f$.
162  /// If element is loaded on both sides, the fluid load computed
163  /// from the reverse element is \b subtracted, because it's
164  /// computed with the same normal vector.
165  //==================================================================
166  void FSIWallElement::fluid_load_vector(const unsigned& intpt,
167  const Vector<double>& N,
168  Vector<double>& load)
169  {
170 
171  // Size of load vector
172  unsigned n_load=load.size();
173 
174  // Initialise
175  for (unsigned i=0;i<n_load;i++) load[i]=0.0;
176 
177  // Create vector for fluid load
178  Vector<double> fluid_load(n_load);
179 
180  // Loop over front and back if required: Get number of fluid-loaded faces
181  unsigned n_loaded_face=2;
182  if (Only_front_is_loaded_by_fluid) n_loaded_face=1;
183 
184  for (unsigned face=0;face<n_loaded_face;face++)
185  {
186  //Get the local coordinate in the fluid element (copy
187  // operation for Vector)
188  Vector<double> s_adjacent(external_element_local_coord(face,intpt));
189 
190  //Call the load function for adjacent element
191  FSIFluidElement* el_f_pt=dynamic_cast<FSIFluidElement*>
192  (external_element_pt(face,intpt));
193  if (el_f_pt!=0)
194  {
195  el_f_pt->get_load(s_adjacent,N,fluid_load);
196  }
197  else
198  {
199 #ifdef PARANOID
201  {
202  std::ostringstream warning_stream;
203  warning_stream
204  << "Info: No adjacent element set in FSIWallElement.\n\n"
205  << "Note: you can disable this message by setting \n "
206  << "FSIWallElement::Dont_warn_about_missing_adjacent_fluid_elements"
207  << "\n to true or recompiling without PARANOID.\n";
208  OomphLibWarning(warning_stream.str(),
209  "FSIWallElement::fluid_load_vector()",
210  OOMPH_EXCEPTION_LOCATION);
211  }
212 #endif
213  }
214 
215  // Adjust sign of fluid traction. If normal on the front
216  // points into the fluid, the normal at the back points away
217  // from it.
218  double sign=1.0;
219  if (face==1) sign=-1.0;
220 
221  //The load is scaled by the stress-scale ratio Q
222  for(unsigned i=0;i<n_load;i++)
223  {
224  load[i] += fluid_load[i]*sign*q();
225  }
226 
227  } // end of loop over faces
228 
229  }
230 
231 
232  //==================================================================
233  /// Update the nodal positions in all fluid elements that affect
234  /// the traction on this FSIWallElement
235  //==================================================================
237  {
238  // Don't update elements repeatedly
239  std::map<FSIFluidElement*,bool> done;
240 
241  // Number of integration points
242  unsigned n_intpt=integral_pt()->nweight();
243 
244  // Loop over front and back if required: Get number of fluid-loaded faces
245  unsigned n_loaded_face=2;
246  if (Only_front_is_loaded_by_fluid) n_loaded_face=1;
247  for (unsigned face=0;face<n_loaded_face;face++)
248  {
249 
250  // Loop over all integration points in wall element
251  for (unsigned iint=0;iint<n_intpt;iint++)
252  {
253  // Get fluid element that affects this integration point
254  FSIFluidElement* el_f_pt=dynamic_cast<FSIFluidElement*>
255  (external_element_pt(face,iint));
256 
257  // Is there an adjacent fluid element?
258  if (el_f_pt!=0)
259  {
260  // Have we updated its positions yet?
261  if (!done[el_f_pt])
262  {
263  // Update nodal positions
264  el_f_pt->node_update();
265  done[el_f_pt]=true;
266  }
267  }
268  }
269  }
270  }
271 
272 
273 //=================================================================
274 /// Static default value for the ratio of stress scales
275 /// used in the fluid and solid equations (default is 1.0)
276 //=================================================================
278 
279 
280 //=======================================================================
281 /// Overload the function that must return all field data involved
282 /// in the interactions from the external (fluid) element. It allows
283 /// the velocity degrees of freedom to be ignored if we want to
284 /// ignore the shear stresses when computing the Jacobian.
285 //=======================================================================
287  Vector<std::set<FiniteElement*> > const &external_elements_pt,
288  std::set<std::pair<Data*,unsigned> > &paired_interaction_data)
289  {
290  //Loop over each inteaction
291  const unsigned n_interaction = this->ninteraction();
292  for(unsigned i=0;i<n_interaction;i++)
293  {
294  //Loop over each element in the set
295  for(std::set<FiniteElement*>::const_iterator it=
296  external_elements_pt[i].begin();
297  it != external_elements_pt[i].end(); it++)
298  {
299  //Cast the element the specific fluid element
300  FSIFluidElement *external_fluid_el_pt =
301  dynamic_cast<FSIFluidElement*>(*it);
302 
303  //Only add pressure load if so desired
305  {
306  // Add the "pressure" load data
307  external_fluid_el_pt->identify_pressure_data(paired_interaction_data);
308  }
309  else
310  {
311  // Add the "direct" load data (usually velocities and pressures)
312  // to the set
313  external_fluid_el_pt->identify_load_data(paired_interaction_data);
314  }
315  }
316  } //End of loop over interactions
317  }
318 
319  //===================================================================
320  /// Function that must return all geometric data involved
321  /// in the desired interactions from the external element
322  //===================================================================
324  Vector<std::set<FiniteElement*> > const &external_elements_pt,
325  std::set<Data*> &external_geometric_data_pt)
326  {
327  //If we are ignoring the shear stress, then we can ignore all
328  //external geometric data
330  {
331  return;
332  }
333  else
334  {
335  //Loop over each inteaction
336  const unsigned n_interaction = this->ninteraction();
337  for(unsigned i=0;i<n_interaction;i++)
338  {
339  //Loop over each element in the set
340  for(std::set<FiniteElement*>::const_iterator it=
341  external_elements_pt[i].begin();
342  it != external_elements_pt[i].end(); it++)
343  {
344  (*it)->identify_geometric_data(external_geometric_data_pt);
345  }
346  }
347  }
348  }
349 
350 
351 }
void node_update_adjacent_fluid_elements()
Update the nodal positions in all fluid elements that affect the traction on this FSIWallElement...
Definition: fsi.cc:236
Vector< double > & external_element_local_coord(const unsigned &interaction_index, const unsigned &ipt)
Access function to get source element's local coords for specified interaction index at specified int...
void describe_local_dofs(std::ostream &out, const std::string &current_string) const
Function to describe the local dofs of the element. The ostream specifies the output stream to which ...
Definition: fsi.cc:84
void identify_all_geometric_data_for_external_interaction(Vector< std::set< FiniteElement * > > const &external_elements_pt, std::set< Data * > &external_geometric_data_pt)
Function that must return all geometric data involved in the desired interactions from the external e...
Definition: fsi.cc:323
virtual void get_load(const Vector< double > &s, const Vector< double > &N, Vector< double > &load)=0
Compute the load vector that is applied by current element (at its local coordinate s) onto the adjac...
cstr elem_len * i
Definition: cfortran.h:607
void setup_fsi_wall_element(const unsigned &nlagr_solid, const unsigned &ndim_fluid)
Setup: Assign storage – pass the Eulerian dimension of the "adjacent" fluid elements and the number o...
Definition: fsi.cc:124
Nodes are derived from Data, but, in addition, have a definite (Eulerian) position in a space of a gi...
Definition: nodes.h:852
void set_ninteraction(const unsigned &n_interaction)
Set the number of interactions in the element This function is usually called in the specific element...
void apply_no_slip_on_moving_wall(Node *node_pt)
Apply no-slip condition for N.St. on a moving wall node, u = St dR/dt, where the Strouhal number St =...
Definition: fsi.cc:54
virtual void identify_pressure_data(std::set< std::pair< Data *, unsigned > > &paired_pressure_data)=0
Add to the set paired_pressure_data pairs containing.
virtual void identify_load_data(std::set< std::pair< Data *, unsigned > > &paired_load_data)=0
Add to the set paired_load_data pairs containing.
bool Only_front_is_loaded_by_fluid
Is the element exposed to (and hence loaded by) fluid only on its "front"? True by default...
Definition: fsi.h:541
void describe_solid_local_dofs(std::ostream &out, const std::string &current_string) const
Classifies dofs locally for solid specific aspects.
Definition: elements.cc:6548
virtual void describe_dofs(std::ostream &out, const std::string &current_string) const
Function to describe the dofs of the Node. The ostream specifies the output stream to which the descr...
Definition: nodes.cc:905
void set_value(const unsigned &i, const double &value_)
Set the i-th stored data value to specified value. The only reason that we require an explicit set fu...
Definition: nodes.h:267
Integral *const & integral_pt() const
Return the pointer to the integration scheme (const version)
Definition: elements.h:1900
void identify_all_field_data_for_external_interaction(Vector< std::set< FiniteElement * > > const &external_elements_pt, std::set< std::pair< Data *, unsigned > > &paired_iteraction_data)
Overload the function that must return all field data involved in the interactions from the external ...
Definition: fsi.cc:286
bool Ignore_shear_stress_in_jacobian
Set this flag to true to ignore shear stress component of load when calculating the Jacobian...
Definition: fsi.h:553
void fluid_load_vector(const unsigned &intpt, const Vector< double > &N, Vector< double > &load)
Get FE Jacobian by systematic finite differencing w.r.t. nodal positition Data, internal and external...
Definition: fsi.cc:166
unsigned nexternal_interaction_geometric_data() const
Return the number of geometric Data items that affect the external interactions in this element: i...
unsigned nexternal_interaction_field_data() const
Return the number of Data items that affect the external interactions in this element. This includes e.g. fluid velocities and pressures in adjacent fluid elements in an FSI problem.
void set_nlagrangian_and_ndim(const unsigned &n_lagrangian, const unsigned &n_dim)
Set # of Lagrangian and Eulerian coordinates.
Definition: geom_objects.h:188
virtual unsigned nweight() const =0
Return the number of integration points of the scheme.
FiniteElement *& external_element_pt(const unsigned &interaction_index, const unsigned &ipt)
Access function to source element for specified interaction index at specified integration point...
void enable_fluid_loading_on_both_sides()
Allow element to be loaded by fluid on both sides. (Resizes containers for lookup schemes and initial...
Definition: fsi.cc:142
std::string string(const unsigned &i)
Return the i-th string or "" if the relevant string hasn't been defined.
unsigned ninteraction() const
Return the number of interactions in the element.
unsigned ndim() const
Return (Eulerian) spatial dimension of the node.
Definition: nodes.h:992
static bool Dont_warn_about_missing_adjacent_fluid_elements
Static flag that allows the suppression of warning messages.
Definition: fsi.h:228
double Strouhal_for_no_slip
Strouhal number St = a/(UT) for application of no slip condition. Initialised to 1.0.
Definition: fsi.cc:45
double dposition_dt(const unsigned &i) const
Return the i-th component of nodal velocity: dx/dt, either directly or via hanging node representatio...
Definition: nodes.cc:2542
virtual void describe_local_dofs(std::ostream &out, const std::string &current_string) const
Function to describe the local dofs of the element[s]. The ostream specifies the output stream to whi...
Definition: elements.cc:1680
const double & q() const
Return the ratio of the stress scales used to non-dimensionalise the fluid and solid equations...
Definition: fsi.h:263
static double Default_Q_Value
Static default value for the ratio of stress scales used in the fluid and solid equations (default is...
Definition: fsi.h:534
The FSIFluidElement class is a base class for all fluid finite elements that apply a load (traction) ...
Definition: fsi.h:67
virtual void node_update()
Update the positions of all nodes in the element using each node update function. The default impleme...
Definition: elements.cc:4865