LBCouetteSolverInterface.h

// Copyright (C) 2016 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_SOLVERS_LBCOUETTESOLVERINTERFACE_H_
#define _MOLECULARDYNAMICS_COUPLING_SOLVERS_LBCOUETTESOLVERINTERFACE_H_
 
#include "coupling/interface/MacroscopicSolverInterface.h"
 
namespace coupling {
namespace solvers {
class LBCouetteSolverInterface;
}
} // namespace coupling

class coupling::solvers::LBCouetteSolverInterface : public coupling::interface::MacroscopicSolverInterface<3> {
public:
  LBCouetteSolverInterface(tarch::la::Vector<3, unsigned int> avgNumberLBCells, tarch::la::Vector<3, unsigned int> numberProcesses,
                           tarch::la::Vector<3, int> offsetMDDomain, tarch::la::Vector<3, unsigned int> globalNumberCouplingCells, unsigned int outerRegion = 1)
      : _avgNumberLBCells(avgNumberLBCells), _numberProcesses(numberProcesses), _offsetMDDomain(offsetMDDomain), _outerRegion(outerRegion),
        _globalNumberCouplingCells(globalNumberCouplingCells) {}
  ~LBCouetteSolverInterface() {}
 
  unsigned int getOuterRegion() override { return _outerRegion; }

  std::vector<unsigned int> getRanks(I01 idx) override {
    std::vector<unsigned int> ranks;
    // determine global index of cell in LB simulation
    tarch::la::Vector<3, int> globalLBCellIndex(idx.get() + _offsetMDDomain);
    // modify global LB cell index due to ghost layer
    for (int d = 0; d < 3; d++) {
#if (COUPLING_MD_DEBUG == COUPLING_MD_YES)
      std::cout << "LB cell index for global cell index " << idx << ": " << globalLBCellIndex << std::endl;
#endif
      globalLBCellIndex[d]--;
    }
    // loop over all neighbouring cells within a one-cell surrounding and detect
    // the respective ranks. IMPROVE: This currently only allows for simulations
    // with MD located inside the domain (no simulation across boundary)
    for (int z = -1; z < 2; z++) {
      for (int y = -1; y < 2; y++) {
        for (int x = -1; x < 2; x++) {
          // neighbour cell index
          const tarch::la::Vector<3, unsigned int> nbIndex(globalLBCellIndex[0] + x, globalLBCellIndex[1] + y, globalLBCellIndex[2] + z);
          // coordinates of process of neighbour celll
          const tarch::la::Vector<3, unsigned int> processCoordinates((nbIndex[0] / _avgNumberLBCells[0]), (nbIndex[1] / _avgNumberLBCells[1]),
                                                                      (nbIndex[2] / _avgNumberLBCells[2]));
          if (processCoordinates[0] < 0 || processCoordinates[0] > _numberProcesses[0] - 1)
            continue;
          if (processCoordinates[1] < 0 || processCoordinates[1] > _numberProcesses[1] - 1)
            continue;
          if (processCoordinates[2] < 0 || processCoordinates[2] > _numberProcesses[2] - 1)
            continue;
          // corresponding rank
          const unsigned int rank = processCoordinates[0] + _numberProcesses[0] * (processCoordinates[1] + processCoordinates[2] * _numberProcesses[1]);
 
          // if this rank is not part of the vector, push it back
          bool found = false;
          for (unsigned int i = 0; i < ranks.size(); i++) {
            found = found || (rank == ranks[i]);
          }
          if (!found) {
            ranks.push_back(rank);
          }
        }
      }
    }
#if (COUPLING_MD_DEBUG == COUPLING_MD_YES)
    std::cout << "Ranks for cell " << idx << ":";
    for (unsigned int i = 0; i < ranks.size(); i++) {
      std::cout << " " << ranks[i];
    }
    std::cout << std::endl;
#endif
    return ranks;
  }

  std::vector<unsigned int> getSourceRanks(I01 idx) override {
    // determine global index of cell in LB simulation
    tarch::la::Vector<3, int> globalLBCellIndex(idx.get() + _offsetMDDomain);
    // modify global LB cell index due to ghost layer
    for (int d = 0; d < 3; d++) {
#if (COUPLING_MD_DEBUG == COUPLING_MD_YES)
      std::cout << "LB cell index for global cell index " << idx << ": " << globalLBCellIndex << std::endl;
#endif
      globalLBCellIndex[d]--;
    }
    std::vector<unsigned int> ranks;
    // determine process coordinates and respective rank
    const tarch::la::Vector<3, unsigned int> processCoordinates(globalLBCellIndex[0] / _avgNumberLBCells[0], globalLBCellIndex[1] / _avgNumberLBCells[1],
                                                                globalLBCellIndex[2] / _avgNumberLBCells[2]);
    // if idx is outside of LB domain, return no ranks
    for (int d = 0; d < 3; d++)
      if (processCoordinates[d] < 0 || processCoordinates[d] > _numberProcesses[d] - 1)
        return ranks;
 
    const unsigned int rank = processCoordinates[0] + _numberProcesses[0] * (processCoordinates[1] + processCoordinates[2] * _numberProcesses[1]);
    ranks.push_back(rank);
#if (COUPLING_MD_DEBUG == COUPLING_MD_YES)
    std::cout << "Source rank for cell " << idx << ": " << ranks[0] << std::endl;
#endif
    return ranks;
  }
 
private:
  const tarch::la::Vector<3, unsigned int> _avgNumberLBCells;
  const tarch::la::Vector<3, unsigned int> _numberProcesses;
  const tarch::la::Vector<3, int> _offsetMDDomain;
  const unsigned int _outerRegion;
  const tarch::la::Vector<3, unsigned int> _globalNumberCouplingCells;
};
 
#endif // _MOLECULARDYNAMICS_COUPLING_SOLVERS_LBCOUETTESOLVERINTERFACE_H_