VelocityGradientRelaxation.h

// Copyright (C) 2015 Technische Universitaet Muenchen
// This file is part of the Mamico project. For conditions of distribution
// and use, please see the copyright notice in Mamico's main folder, or at
// www5.in.tum.de/mamico
#ifndef _MOLECULARDYNAMICS_COUPLING_VELOCITYGRADIENTRELAXATION_H_
#define _MOLECULARDYNAMICS_COUPLING_VELOCITYGRADIENTRELAXATION_H_
 
#include "coupling/MomentumInsertion.h"
#include "coupling/cell-mappings/ComputeMomentumMapping.h"
#include "coupling/cell-mappings/VelocityGradientRelaxationMapping.h"
#include "coupling/interface/MDSolverInterface.h"
 
namespace coupling {
template <class LinkedCell, unsigned int dim> class VelocityGradientRelaxation;
template <class LinkedCell, unsigned int dim> class VelocityGradientRelaxationTopOnly;
} // namespace coupling

template <class LinkedCell, unsigned int dim> class coupling::VelocityGradientRelaxation : public coupling::MomentumInsertion<LinkedCell, dim> {
public:
  VelocityGradientRelaxation(double relaxationParam, coupling::interface::MDSolverInterface<LinkedCell, dim>* const mdSolverInterface,
                             const coupling::datastructures::CouplingCellWithLinkedCells<LinkedCell, dim>* const couplingCells)
      : coupling::MomentumInsertion<LinkedCell, dim>(mdSolverInterface), _couplingCells(couplingCells), _relaxationParam(relaxationParam) {}

  virtual ~VelocityGradientRelaxation() {}

  unsigned int getTimeIntervalPerMomentumInsertion() const override { return 1; }

  void insertMomentum(coupling::datastructures::CouplingCellWithLinkedCells<LinkedCell, dim>& cell, I02 idx) const override {
    coupling::cellmappings::ComputeMomentumMapping<LinkedCell, dim> momentumMapping(coupling::MomentumInsertion<LinkedCell, dim>::_mdSolverInterface);
    tarch::la::Vector<dim, double> oldVelocity(0.0);
    cell.iterateConstCells(momentumMapping);
    // set current average velocity within this cell
    oldVelocity = momentumMapping.getMeanVelocity();
    // set new momentum (based on velocity stored in microscopic
    // momentum-buffer)
    coupling::cellmappings::VelocityGradientRelaxationMapping<LinkedCell, dim> velocityGradientRelaxation(
        _relaxationParam, oldVelocity, idx, coupling::MomentumInsertion<LinkedCell, dim>::_mdSolverInterface, _couplingCells);
    cell.iterateCells(velocityGradientRelaxation);
  }
 
  void setInnerImposition(bool enable) override {
    throw std::runtime_error(std::string("coupling::AdditiveMomentumInsertion::setInnerImposition not implemented"));
  }
 
protected:
  const coupling::datastructures::CouplingCellWithLinkedCells<LinkedCell, dim>* const _couplingCells;
  const double _relaxationParam;
};

template <class LinkedCell, unsigned int dim> class coupling::VelocityGradientRelaxationTopOnly : public coupling::VelocityGradientRelaxation<LinkedCell, dim> {
public:
  VelocityGradientRelaxationTopOnly(double relaxationParam, coupling::interface::MDSolverInterface<LinkedCell, dim>* const mdSolverInterface,
                                    const coupling::datastructures::CouplingCellWithLinkedCells<LinkedCell, dim>* const couplingCells)
      : coupling::VelocityGradientRelaxation<LinkedCell, dim>(relaxationParam, mdSolverInterface, couplingCells) {}

  virtual ~VelocityGradientRelaxationTopOnly() {}

  void insertMomentum(coupling::datastructures::CouplingCellWithLinkedCells<LinkedCell, dim>& cell, I02 idx) const override {
    coupling::cellmappings::ComputeMomentumMapping<LinkedCell, dim> momentumMapping(coupling::MomentumInsertion<LinkedCell, dim>::_mdSolverInterface);
    tarch::la::Vector<dim, double> oldVelocity(0.0);
    cell.iterateConstCells(momentumMapping);
    // set current average velocity within this cell
    oldVelocity = momentumMapping.getMeanVelocity();
 
    // set new momentum (based on velocity stored in microscopic
    // momentum-buffer)
    coupling::cellmappings::VelocityGradientRelaxationTopOnlyMapping<LinkedCell, dim> velocityGradientRelaxation(
        coupling::VelocityGradientRelaxation<LinkedCell, dim>::_relaxationParam, oldVelocity, idx,
        coupling::MomentumInsertion<LinkedCell, dim>::_mdSolverInterface, coupling::VelocityGradientRelaxation<LinkedCell, dim>::_couplingCells);
    cell.iterateCells(velocityGradientRelaxation);
  }
};
#endif // _MOLECULARDYNAMICS_COUPLING_VELOCITYGRADIENTRELAXATION_H_