aggregatedTypeVectorAccessWrapper.hpp

/*
 * CoDiPack, a Code Differentiation Package
 *
 * Copyright (C) 2015-2024 Chair for Scientific Computing (SciComp), University of Kaiserslautern-Landau
 * Homepage: http://www.scicomp.uni-kl.de
 * Contact:  Prof. Nicolas R. Gauger (codi@scicomp.uni-kl.de)
 *
 * Lead developers: Max Sagebaum, Johannes Blühdorn (SciComp, University of Kaiserslautern-Landau)
 *
 * This file is part of CoDiPack (http://www.scicomp.uni-kl.de/software/codi).
 *
 * CoDiPack is free software: you can redistribute it and/or
 * modify it under the terms of the GNU General Public License
 * as published by the Free Software Foundation, either version 3 of the
 * License, or (at your option) any later version.
 *
 * CoDiPack is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty
 * of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
 *
 * See the GNU General Public License for more details.
 * You should have received a copy of the GNU
 * General Public License along with CoDiPack.
 * If not, see <http://www.gnu.org/licenses/>.
 *
 * For other licensing options please contact us.
 *
 * Authors:
 *  - SciComp, University of Kaiserslautern-Landau:
 *    - Max Sagebaum
 *    - Johannes Blühdorn
 *    - Former members:
 *      - Tim Albring
 */
#pragma once
 
#include <complex>
#include <vector>
 
#include "../../config.h"
#include "../../expressions/lhsExpressionInterface.hpp"
#include "../../misc/macros.hpp"
#include "../../tapes/misc/vectorAccessInterface.hpp"
#include "../../traits/computationTraits.hpp"
#include "../../traits/expressionTraits.hpp"
#include "../../traits/realTraits.hpp"
 
namespace codi {
 
  template<typename T_Type, typename = void>
  struct AggregatedTypeVectorAccessWrapper : public VectorAccessInterface<CODI_ANY, CODI_ANY> {
      CODI_STATIC_ASSERT(false && std::is_void<T_Type>::value,
                         "Instantiation of unspecialized AggregatedTypeVectorAccessWrapper.");
 
      using Type = CODI_DD(T_Type, CODI_ANY);  
  };
 
  template<typename T_Real, typename T_Identifier, typename T_InnerInterface>
  struct AggregatedTypeVectorAccessWrapperBase : public VectorAccessInterface<T_Real, T_Identifier> {
    public:
 
      using Real = CODI_DD(T_Real, double);           
      using Identifier = CODI_DD(T_Identifier, int);  
 
      using InnerInterface =
          CODI_DD(T_InnerInterface,
                  CODI_T(VectorAccessInterface<double, int>));  
 
    protected:
 
      InnerInterface& innerInterface;  
 
      std::vector<Real> lhs;     
      std::vector<Real> buffer;  
 
    public:
 
      AggregatedTypeVectorAccessWrapperBase(InnerInterface* innerInterface)
          : innerInterface(*innerInterface),
            lhs(innerInterface->getVectorSize()),
            buffer(innerInterface->getVectorSize()) {}
 
      /*******************************************************************************/
 
      size_t getVectorSize() const {
        return innerInterface.getVectorSize();
      }
 
      bool isLhsZero() {
        bool isZero = true;
 
        for (size_t curDim = 0; isZero && curDim < lhs.size(); curDim += 1) {
          isZero &= RealTraits::isTotalZero(lhs[curDim]);
        }
 
        return isZero;
      }
 
      /*******************************************************************************/
 
      void setLhsAdjoint(Identifier const& index) {
        getAdjointVec(index, lhs.data());
        this->resetAdjointVec(index);
      }
 
      void updateAdjointWithLhs(Identifier const& index, Real const& jacobian) {
        for (size_t curDim = 0; curDim < lhs.size(); curDim += 1) {
          Real update = jacobian * lhs[curDim];
          this->updateAdjoint(index, curDim, update);
        }
      }
 
      /*******************************************************************************/
 
      void setLhsTangent(Identifier const& index) {
        updateAdjointVec(index, lhs.data());
 
        for (size_t curDim = 0; curDim < lhs.size(); curDim += 1) {
          lhs[curDim] = Real();
        }
      }
 
      void updateTangentWithLhs(Identifier const& index, Real const& jacobian) {
        for (size_t curDim = 0; curDim < lhs.size(); curDim += 1) {
          lhs[curDim] += jacobian * this->getAdjoint(index, curDim);
        }
      }
 
      /*******************************************************************************/
 
      void getAdjointVec(Identifier const& index, Real* const vec) {
        for (size_t curDim = 0; curDim < lhs.size(); curDim += 1) {
          vec[curDim] = this->getAdjoint(index, curDim);
        }
      }
 
      Real const* getAdjointVec(Identifier const& index) {
        getAdjointVec(index, buffer.data());
        return buffer.data();
      }
 
      void updateAdjointVec(Identifier const& index, Real const* const vec) {
        for (size_t curDim = 0; curDim < lhs.size(); curDim += 1) {
          this->updateAdjoint(index, curDim, vec[curDim]);
        }
      }
 
      /*******************************************************************************/
 
      bool hasPrimals() {
        return innerInterface.hasPrimals();
      }
  };
 
  template<typename T_Type, typename = void>
  struct AggregatedTypeVectorAccessWrapperFactory {
    public:
      using Type = CODI_DD(T_Type, CODI_ANY);  
 
      using RType = AggregatedTypeVectorAccessWrapper<Type>;  
 
      template<typename Real, typename Identifier>
      static RType* create(VectorAccessInterface<Real, Identifier>* access) {
        return new RType(access);
      }
 
      static void destroy(RType* access) {
        delete access;
      }
  };
 
#ifndef DOXYGEN_DISABLE
  template<typename T_InnerType>
  struct AggregatedTypeVectorAccessWrapper<std::complex<T_InnerType>>
      : public AggregatedTypeVectorAccessWrapperBase<
            std::complex<typename T_InnerType::Real>, std::complex<typename T_InnerType::Identifier>,
            VectorAccessInterface<typename T_InnerType::Real, typename T_InnerType::Identifier>> {
    public:
 
      using InnerType = CODI_DD(T_InnerType, CODI_DEFAULT_LHS_EXPRESSION);  
      using Type = std::complex<InnerType>;                                 
 
      using InnerInterface = VectorAccessInterface<
          typename InnerType::Real,
          typename InnerType::Identifier>;  
 
      using Real = std::complex<typename InnerType::Real>;              
      using Identifier = std::complex<typename InnerType::Identifier>;  
 
      using Base =
          AggregatedTypeVectorAccessWrapperBase<Real, Identifier, InnerInterface>;  
 
      AggregatedTypeVectorAccessWrapper(InnerInterface* innerInterface) : Base(innerInterface) {}
 
      /*******************************************************************************/
 
      void resetAdjoint(Identifier const& index, size_t dim) {
        Base::innerInterface.resetAdjoint(std::real(index), dim);
        Base::innerInterface.resetAdjoint(std::imag(index), dim);
      }
 
      void resetAdjointVec(Identifier const& index) {
        Base::innerInterface.resetAdjointVec(std::real(index));
        Base::innerInterface.resetAdjointVec(std::imag(index));
      }
 
      Real getAdjoint(Identifier const& index, size_t dim) {
        return Real(Base::innerInterface.getAdjoint(std::real(index), dim),
                    Base::innerInterface.getAdjoint(std::imag(index), dim));
      }
 
      void updateAdjoint(Identifier const& index, size_t dim, Real const& adjoint) {
        Base::innerInterface.updateAdjoint(std::real(index), dim, std::real(adjoint));
        Base::innerInterface.updateAdjoint(std::imag(index), dim, std::imag(adjoint));
      }
 
      /*******************************************************************************/
 
      void setPrimal(Identifier const& index, Real const& primal) {
        Base::innerInterface.setPrimal(std::real(index), std::real(primal));
        Base::innerInterface.setPrimal(std::imag(index), std::imag(primal));
      }
 
      Real getPrimal(Identifier const& index) {
        return Real(Base::innerInterface.getPrimal(std::real(index)), Base::innerInterface.getPrimal(std::imag(index)));
      }
 
      VectorAccessInterface<Real, Identifier>* clone() const {
        return new AggregatedTypeVectorAccessWrapper(Base::innerInterface.clone());
      }
  };
 
  template<typename T_Type>
  struct AggregatedTypeVectorAccessWrapperFactory<T_Type, ExpressionTraits::EnableIfLhsExpression<T_Type>> {
    public:
      using Type = CODI_DD(T_Type, CODI_DEFAULT_LHS_EXPRESSION);  
 
      using RType = VectorAccessInterface<typename Type::Real, typename Type::Identifier>;
 
      static RType* create(RType* access) {
        return access;
      }
 
      static void destroy(RType* access) {
        CODI_UNUSED(access);
 
        // Do nothing
      }
  };
#endif
}