adJointImpl.hpp

Go to the documentation of this file.

/** \file dHODataOperators.hpp
  * \brief Operators managing HO geometry derivatives
 
*/
 
#ifndef __ADJOINT_IMPL_HPP__
#define __ADJOINT_IMPL_HPP__
 
namespace MoFEM {
 
//! [declarations of templates]
 
template <int SPACE_DIM, IntegrationType I, typename OpBase>
struct OpGetCoFactorImpl;
 
//! [declarations of templates]
 
//! [definitions of templates]
 
template <int SPACE_DIM, typename OpBase>
struct OpGetCoFactorImpl<SPACE_DIM, IntegrationType::GAUSS, OpBase>
    : public OpBase {
 
  OpGetCoFactorImpl(boost::shared_ptr<MatrixDouble> jac_ptr,
                    boost::shared_ptr<MatrixDouble> diff_jac_ptr,
                    boost::shared_ptr<VectorDouble> diff_out_ptr
 
                    )
      : OpBase(NOSPACE, OpBase::OPSPACE), jacPtr(jac_ptr),
        diffJacPtr(diff_jac_ptr), diffOutPtr(diff_out_ptr) {}
 
  MoFEMErrorCode doWork(int side, EntityType type,
                        EntitiesFieldData::EntData &data);
 
protected:
  boost::shared_ptr<MatrixDouble> jacPtr;
  boost::shared_ptr<MatrixDouble> diffJacPtr;
  boost::shared_ptr<VectorDouble> diffOutPtr;
};
 
//! [definitions of templates]
 
//! [implementations of templates functions]
 
template <int SPACE_DIM, typename OpBase>
MoFEMErrorCode
OpGetCoFactorImpl<SPACE_DIM, IntegrationType::GAUSS, OpBase>::doWork(
    int side, EntityType type, EntitiesFieldData::EntData &data) {
  MoFEMFunctionBegin;
  auto nb_integration_pts = OpBase::getGaussPts().size2();
  diffOutPtr->resize(nb_integration_pts, false);
  FTENSOR_INDEX(SPACE_DIM, i);
  FTENSOR_INDEX(SPACE_DIM, j);
  auto t_jac = getFTensor2FromMat<SPACE_DIM, SPACE_DIM>(*jacPtr);
  auto t_diff_jac = getFTensor2FromMat<SPACE_DIM, SPACE_DIM>(*diffJacPtr);
  auto t_diff_out = getFTensor0FromVec(*diffOutPtr);
  for (int gg = 0; gg != nb_integration_pts; gg++) {
    auto t_det = determinantTensor(t_jac);
    FTensor::Tensor2<double, SPACE_DIM, SPACE_DIM> t_inv_jac;
    CHKERR invertTensor(t_jac, t_det, t_inv_jac);
    // note that integration weight is 
    t_diff_out = (t_inv_jac(j, i) * t_diff_jac(i, j)) * t_det; 
    ++t_jac;
    ++t_diff_jac;
    ++t_diff_out;
  }
  MoFEMFunctionReturn(0);
}
 
//! [implementations of templates functions]
 
}; // namespace MoFEM
 
#endif