pml_time_harmonic_elasticity_tensor.h

Go to the documentation of this file.

// LIC// ====================================================================
// LIC// This file forms part of oomph-lib, the object-oriented,
// LIC// multi-physics finite-element library, available
// LIC// at http://www.oomph-lib.org.
// LIC//
// LIC// Copyright (C) 2006-2022 Matthias Heil and Andrew Hazel
// LIC//
// LIC// This library is free software; you can redistribute it and/or
// LIC// modify it under the terms of the GNU Lesser General Public
// LIC// License as published by the Free Software Foundation; either
// LIC// version 2.1 of the License, or (at your option) any later version.
// LIC//
// LIC// This library is distributed in the hope that it will be useful,
// LIC// but WITHOUT ANY WARRANTY; without even the implied warranty of
// LIC// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
// LIC// Lesser General Public License for more details.
// LIC//
// LIC// You should have received a copy of the GNU Lesser General Public
// LIC// License along with this library; if not, write to the Free Software
// LIC// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
// LIC// 02110-1301  USA.
// LIC//
// LIC// The authors may be contacted at oomph-lib@maths.man.ac.uk.
// LIC//
// LIC//====================================================================
// Header file for objects representing the fourth-rank elasticity tensor
// for linear elasticity problems
 
// Include guards to prevent multiple inclusion of the header
#ifndef OOMPH_PML_TIME_HARMONIC_ELASTICITY_TENSOR_HEADER
#define OOMPH_PML_TIME_HARMONIC_ELASTICITY_TENSOR_HEADER
 
// Config header generated by autoconfig
#ifdef HAVE_CONFIG_H
#include <oomph-lib-config.h>
#endif
 
#include "../generic/oomph_utilities.h"
#include "math.h"
#include <complex>
 
namespace oomph
{
  //=====================================================================
  //=====================================================================
  class PMLTimeHarmonicElasticityTensor
  {
    static const unsigned Index[3][3][3][3];
 
  protected:
    virtual inline std::complex<double> independent_component(
      const unsigned& i) const
    {
      return std::complex<double>(0.0, 0.0);
    }
 
 
    void range_check(const unsigned& i,
                     const unsigned& j,
                     const unsigned& k,
                     const unsigned& l) const
    {
      if ((i > 2) || (j > 2) || (k > 2) || (l > 2))
      {
        std::ostringstream error_message;
        if (i > 2)
        {
          error_message << "Range Error : Index 1 " << i
                        << " is not in the range (0,2)";
        }
        if (j > 2)
        {
          error_message << "Range Error : Index 2 " << j
                        << " is not in the range (0,2)";
        }
 
        if (k > 2)
        {
          error_message << "Range Error : Index 2 " << k
                        << " is not in the range (0,2)";
        }
 
        if (l > 2)
        {
          error_message << "Range Error : Index 4 " << l
                        << " is not in the range (0,2)";
        }
 
        // Throw the error
        throw OomphLibError(error_message.str(),
                            OOMPH_CURRENT_FUNCTION,
                            OOMPH_EXCEPTION_LOCATION);
      }
    }
 
 
    PMLTimeHarmonicElasticityTensor() {}
 
  public:
    virtual ~PMLTimeHarmonicElasticityTensor() {}
 
  public:
    std::complex<double> operator()(const unsigned& i,
                                    const unsigned& j,
                                    const unsigned& k,
                                    const unsigned& l) const
    {
      // Range check
#ifdef PARANOID
      range_check(i, j, k, l);
#endif
      return independent_component(Index[i][j][k][l]);
    }
  };
 
 
  //===================================================================
  //===================================================================
  class PMLTimeHarmonicIsotropicElasticityTensor
    : public PMLTimeHarmonicElasticityTensor
  {
    std::complex<double> C[4];
 
    static const unsigned StaticIndex[21];
 
    double Nu;
 
  public:
    PMLTimeHarmonicIsotropicElasticityTensor(const double& nu, const double& E)
      : PMLTimeHarmonicElasticityTensor(), Nu(nu)
    {
      // Set the three indepdent components
      C[0] = std::complex<double>(0.0, 0.0);
      double lambda = E * nu / ((1.0 + nu) * (1.0 - 2.0 * nu));
      double mu = E / (2.0 * (1.0 + nu));
      this->set_lame_coefficients(lambda, mu);
    }
 
    PMLTimeHarmonicIsotropicElasticityTensor(const double& nu)
      : PMLTimeHarmonicElasticityTensor(), Nu(nu)
    {
      // Set the three indepdent components
      C[0] = std::complex<double>(0.0, 0.0);
 
      // reference value
      double E = 1.0;
      double lambda = E * nu / ((1.0 + nu) * (1.0 - 2.0 * nu));
      double mu = E / (2.0 * (1.0 + nu));
      this->set_lame_coefficients(lambda, mu);
    }
 
    double nu() const
    {
      return Nu;
    }
 
    void update_constitutive_parameters(const double& nu, const double& E = 1.0)
    {
      // Set the three indepdent components
      C[0] = std::complex<double>(0.0, 0.0);
      double lambda = E * nu / ((1.0 + nu) * (1.0 - 2.0 * nu));
      double mu = E / (2.0 * (1.0 + nu));
      this->set_lame_coefficients(lambda, mu);
    }
 
 
    inline std::complex<double> independent_component(const unsigned& i) const
    {
      return C[StaticIndex[i]];
    }
 
 
  private:
    // Set the values of the lame coefficients
    void set_lame_coefficients(const double& lambda, const double& mu)
    {
      C[1] = lambda + 2.0 * mu;
      C[2] = lambda;
      C[3] = mu;
    }
  };
 
} // namespace oomph
#endif