time_harmonic_elasticity_tensor.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-2023 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 // Header file for objects representing the fourth-rank elasticity tensor
27 // for linear elasticity problems
28 
29 // Include guards to prevent multiple inclusion of the header
30 #ifndef OOMPH_TIME_HARMONIC_ELASTICITY_TENSOR_HEADER
31 #define OOMPH_TIME_HARMONIC_ELASTICITY_TENSOR_HEADER
32 
33 // Config header generated by autoconfig
34 #ifdef HAVE_CONFIG_H
35 #include <oomph-lib-config.h>
36 #endif
37 
38 #include "../generic/oomph_utilities.h"
39 
40 namespace oomph
41 {
42  //=====================================================================
43  /// A base class that represents the fourth-rank elasticity tensor
44  /// \f$E_{ijkl}\f$ defined such that
45  /// \f[\tau_{ij} = E_{ijkl} e_{kl},\f]
46  /// where \f$e_{ij}\f$ is the infinitessimal (Cauchy) strain tensor
47  /// and \f$\tau_{ij}\f$ is the stress tensor. The symmetries of the
48  /// tensor are such that
49  /// \f[E_{ijkl} = E_{jikl} = E_{ijlk} = E_{klij}\f]
50  /// and thus there are relatively few independent components. These
51  /// symmetries are included in the definition of the object so that
52  /// non-physical symmetries cannot be accidentally imposed.
53  //=====================================================================
55  {
56  /// Translation table from the four indices to the corresponding
57  /// independent component
58  static const unsigned Index[3][3][3][3];
59 
60  protected:
61  /// Member function that returns the i-th independent component of the
62  /// elasticity tensor
63  virtual inline double independent_component(const unsigned& i) const
64  {
65  return 0.0;
66  }
67 
68 
69  /// Helper range checking function
70  /// (Note that this only captures over-runs in 3D but
71  /// errors are likely to be caught in evaluation of the
72  /// stress and strain tensors anyway...)
73  void range_check(const unsigned& i,
74  const unsigned& j,
75  const unsigned& k,
76  const unsigned& l) const
77  {
78  if ((i > 2) || (j > 2) || (k > 2) || (l > 2))
79  {
80  std::ostringstream error_message;
81  if (i > 2)
82  {
83  error_message << "Range Error : Index 1 " << i
84  << " is not in the range (0,2)";
85  }
86  if (j > 2)
87  {
88  error_message << "Range Error : Index 2 " << j
89  << " is not in the range (0,2)";
90  }
91 
92  if (k > 2)
93  {
94  error_message << "Range Error : Index 2 " << k
95  << " is not in the range (0,2)";
96  }
97 
98  if (l > 2)
99  {
100  error_message << "Range Error : Index 4 " << l
101  << " is not in the range (0,2)";
102  }
103 
104  // Throw the error
105  throw OomphLibError(error_message.str(),
106  OOMPH_CURRENT_FUNCTION,
107  OOMPH_EXCEPTION_LOCATION);
108  }
109  }
110 
111 
112  /// Empty Constructor
114 
115  public:
116  /// Empty virtual Destructor
118 
119  public:
120  /// Return the appropriate independent component
121  /// via the index translation scheme (const version).
122  double operator()(const unsigned& i,
123  const unsigned& j,
124  const unsigned& k,
125  const unsigned& l) const
126  {
127  // Range check
128 #ifdef PARANOID
129  range_check(i, j, k, l);
130 #endif
131  return independent_component(Index[i][j][k][l]);
132  }
133  };
134 
135 
136  //===================================================================
137  /// An isotropic elasticity tensor defined in terms of Young's modulus
138  /// and Poisson's ratio. The elasticity tensor is assumed to be
139  /// non-dimensionalised on some reference value for Young's modulus
140  /// so the value provided to the constructor (if any) is to be
141  /// interpreted as the ratio of the actual Young's modulus to the
142  /// Young's modulus used to non-dimensionalise the stresses/tractions
143  /// in the governing equations.
144  //===================================================================
147  {
148  // Storage for the independent components of the elasticity tensor
149  double C[4];
150 
151  // Translation scheme between the 21 independent components of the general
152  // elasticity tensor and the isotropic case
153  static const unsigned StaticIndex[21];
154 
155  public:
156  /// Constructor. Passing in the values of the Poisson's ratio
157  /// and Young's modulus (interpreted as the ratio of the actual
158  /// Young's modulus to the Young's modulus (or other reference stiffness)
159  /// used to non-dimensionalise stresses and tractions in the governing
160  /// equations).
161  TimeHarmonicIsotropicElasticityTensor(const double& nu, const double& E)
163  {
164  // Set the three indepdent components
165  C[0] = 0.0;
166  double lambda = E * nu / ((1.0 + nu) * (1.0 - 2.0 * nu));
167  double mu = E / (2.0 * (1.0 + nu));
168  this->set_lame_coefficients(lambda, mu);
169  }
170 
171  /// Constructor. Passing in the value of the Poisson's ratio.
172  /// Stresses and tractions in the governing equations are assumed
173  /// to have been non-dimensionalised on Young's modulus.
176  {
177  // Set the three indepdent components
178  C[0] = 0.0;
179 
180  // reference value
181  double E = 1.0;
182  double lambda = E * nu / ((1.0 + nu) * (1.0 - 2.0 * nu));
183  double mu = E / (2.0 * (1.0 + nu));
184  this->set_lame_coefficients(lambda, mu);
185  }
186 
187  /// Update parameters: Specify values of the Poisson's ratio
188  /// and (optionally) Young's modulus (interpreted as the ratio of the actual
189  /// Young's modulus to the Young's modulus (or other reference stiffness)
190  /// used to non-dimensionalise stresses and tractions in the governing
191  /// equations).
192  void update_constitutive_parameters(const double& nu, const double& E = 1.0)
193  {
194  // Set the three indepdent components
195  C[0] = 0.0;
196  double lambda = E * nu / ((1.0 + nu) * (1.0 - 2.0 * nu));
197  double mu = E / (2.0 * (1.0 + nu));
198  this->set_lame_coefficients(lambda, mu);
199  }
200 
201 
202  /// Overload the independent coefficient function
203  inline double independent_component(const unsigned& i) const
204  {
205  return C[StaticIndex[i]];
206  }
207 
208 
209  private:
210  // Set the values of the lame coefficients
211  void set_lame_coefficients(const double& lambda, const double& mu)
212  {
213  C[1] = lambda + 2.0 * mu;
214  C[2] = lambda;
215  C[3] = mu;
216  }
217  };
218 
219 } // namespace oomph
220 #endif
cstr elem_len * i
Definition: cfortran.h:603
An OomphLibError object which should be thrown when an run-time error is encountered....
A base class that represents the fourth-rank elasticity tensor defined such that.
virtual double independent_component(const unsigned &i) const
Member function that returns the i-th independent component of the elasticity tensor.
static const unsigned Index[3][3][3][3]
Translation table from the four indices to the corresponding independent component.
void range_check(const unsigned &i, const unsigned &j, const unsigned &k, const unsigned &l) const
Helper range checking function (Note that this only captures over-runs in 3D but errors are likely to...
double operator()(const unsigned &i, const unsigned &j, const unsigned &k, const unsigned &l) const
Return the appropriate independent component via the index translation scheme (const version).
virtual ~TimeHarmonicElasticityTensor()
Empty virtual Destructor.
An isotropic elasticity tensor defined in terms of Young's modulus and Poisson's ratio....
double independent_component(const unsigned &i) const
Overload the independent coefficient function.
static const unsigned StaticIndex[21]
Translation scheme for the isotropic elasticity tensor.
TimeHarmonicIsotropicElasticityTensor(const double &nu)
Constructor. Passing in the value of the Poisson's ratio. Stresses and tractions in the governing equ...
void update_constitutive_parameters(const double &nu, const double &E=1.0)
Update parameters: Specify values of the Poisson's ratio and (optionally) Young's modulus (interprete...
TimeHarmonicIsotropicElasticityTensor(const double &nu, const double &E)
Constructor. Passing in the values of the Poisson's ratio and Young's modulus (interpreted as the rat...
void set_lame_coefficients(const double &lambda, const double &mu)
//////////////////////////////////////////////////////////////////// ////////////////////////////////...