StrainHelper.hpp

Go to the documentation of this file.

/*
 * ------------------------------------------------------------------------------------------------------------
 * SPDX-License-Identifier: LGPL-2.1-only
 *
 * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC
 * Copyright (c) 2018-2024 Total, S.A
 * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University
 * Copyright (c) 2023-2024 Chevron
 * Copyright (c) 2019-     GEOS/GEOSX Contributors
 * All rights reserved
 *
 * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details.
 * ------------------------------------------------------------------------------------------------------------
 */
 
#ifndef GEOS_PHYSICSSOLVERS_SOLIDMECHANICS_KERNELS_STRAINHELPER_HPP_
#define GEOS_PHYSICSSOLVERS_SOLIDMECHANICS_KERNELS_STRAINHELPER_HPP_
 
#include "common/DataTypes.hpp"
#include "common/GEOS_RAJA_Interface.hpp"
#include "finiteElement/FiniteElementDispatch.hpp"
#include "mesh/CellElementSubRegion.hpp"
#include "mesh/utilities/AverageOverQuadraturePointsKernel.hpp"
#include "physicsSolvers/solidMechanics/SolidMechanicsFields.hpp"
 
namespace geos
{
template< typename SUBREGION_TYPE,
          typename FE_TYPE >
class AverageStrainOverQuadraturePoints :
  public AverageOverQuadraturePointsBase< SUBREGION_TYPE,
                                          FE_TYPE >
{
public:
 
  using Base = AverageOverQuadraturePointsBase< SUBREGION_TYPE,
                                                FE_TYPE >;
 
  using Base::m_elementVolume;
  using Base::m_elemsToNodes;
  using Base::m_finiteElementSpace;
 
  AverageStrainOverQuadraturePoints( NodeManager & nodeManager,
                                     EdgeManager const & edgeManager,
                                     FaceManager const & faceManager,
                                     SUBREGION_TYPE const & elementSubRegion,
                                     FE_TYPE const & finiteElementSpace,
                                     fields::solidMechanics::arrayViewConst2dLayoutTotalDisplacement const displacement,
                                     fields::solidMechanics::arrayView2dLayoutStrain const avgStrain ):
    Base( nodeManager,
          edgeManager,
          faceManager,
          elementSubRegion,
          finiteElementSpace ),
    m_displacement( displacement ),
    m_avgStrain( avgStrain )
  {}
 
  struct StackVariables : Base::StackVariables
  {real64 uLocal[FE_TYPE::maxSupportPoints][3]; };
 
  GEOS_HOST_DEVICE
  void setup( localIndex const k,
              StackVariables & stack ) const
  {
    Base::setup( k, stack );
 
    for( localIndex a = 0; a < FE_TYPE::maxSupportPoints; ++a )
    {
      localIndex const localNodeIndex = m_elemsToNodes( k, a );
      for( int i = 0; i < 3; ++i )
      {
        stack.uLocal[a][i] = m_displacement[localNodeIndex][i];
      }
    }
 
    for( int icomp = 0; icomp < 6; ++icomp )
    {
      m_avgStrain[k][icomp] = 0.0;
    }
  }
 
  GEOS_HOST_DEVICE
  void quadraturePointKernel( localIndex const k,
                              localIndex const q,
                              StackVariables & stack ) const
  {
    //real64 const weight = FE_TYPE::transformedQuadratureWeight( q, stack.xLocal, stack.feStack ) / m_elementVolume[k];
 
    real64 dNdX[ FE_TYPE::maxSupportPoints ][3];
    real64 const detJxW = m_finiteElementSpace.template getGradN< FE_TYPE >( k, q, stack.xLocal, stack.feStack, dNdX );
    real64 strain[6] = {0.0};
    FE_TYPE::symmetricGradient( dNdX, stack.uLocal, strain );
 
    for( int icomp = 0; icomp < 6; ++icomp )
    {
      m_avgStrain[k][icomp] += detJxW*strain[icomp]/m_elementVolume[k];
    }
  }
 
  template< typename POLICY,
            typename KERNEL_TYPE >
  static void
  kernelLaunch( localIndex const numElems,
                KERNEL_TYPE const & kernelComponent )
  {
    forAll< POLICY >( numElems,
                      [=] GEOS_HOST_DEVICE ( localIndex const k )
    {
      typename KERNEL_TYPE::StackVariables stack;
 
      kernelComponent.setup( k, stack );
      for( integer q = 0; q < FE_TYPE::numQuadraturePoints; ++q )
      {
        kernelComponent.quadraturePointKernel( k, q, stack );
      }
    } );
  }
 
protected:
 
  fields::solidMechanics::arrayViewConst2dLayoutTotalDisplacement const m_displacement;
 
  fields::solidMechanics::arrayView2dLayoutStrain const m_avgStrain;
 
};
 
 
 
class AverageStrainOverQuadraturePointsKernelFactory
{
public:
 
  template< typename SUBREGION_TYPE,
            typename FE_TYPE,
            typename POLICY >
  static void
  createAndLaunch( NodeManager & nodeManager,
                   EdgeManager const & edgeManager,
                   FaceManager const & faceManager,
                   SUBREGION_TYPE const & elementSubRegion,
                   FE_TYPE const & finiteElementSpace,
                   fields::solidMechanics::arrayViewConst2dLayoutTotalDisplacement const displacement,
                   fields::solidMechanics::arrayView2dLayoutStrain const avgStrain )
  {
    AverageStrainOverQuadraturePoints< SUBREGION_TYPE, FE_TYPE >
    kernel( nodeManager, edgeManager, faceManager, elementSubRegion, finiteElementSpace,
            displacement, avgStrain );
 
    AverageStrainOverQuadraturePoints< SUBREGION_TYPE, FE_TYPE >::template
    kernelLaunch< POLICY >( elementSubRegion.size(), kernel );
  }
};
 
 
 
}
 
 
#endif /* GEOS_PHYSICSSOLVERS_SOLIDMECHANICS_KERNELS_STRAINHELPER_HPP_ */