coordinate_series.h

// -*-c++-*-
#ifndef STORMM_COORDINATE_SERIES_H
#define STORMM_COORDINATE_SERIES_H
 
#include "copyright.h"
#include "Accelerator/gpu_details.h"
#include "Accelerator/hybrid.h"
#include "Accelerator/hybrid_util.h"
#include "Constants/behavior.h"
#include "Constants/fixed_precision.h"
#include "Constants/hpc_bounds.h"
#include "DataTypes/common_types.h"
#include "FileManagement/file_enumerators.h"
#include "FileManagement/file_listing.h"
#include "FileManagement/file_util.h"
#include "Math/rounding.h"
#include "amber_ascii.h"
#include "coordinateframe.h"
#include "phasespace.h"
#include "trajectory_enumerators.h"
 
namespace stormm {
namespace trajectory {
 
using card::deepCopy;
using card::deepRecast;
using card::default_hpc_format;
using card::GpuDetails;
using card::Hybrid;
using card::HybridTargetLevel;
using constants::CartesianDimension;
using data_types::isFloatingPointScalarType;
using data_types::isSignedIntegralScalarType;
using diskutil::detectCoordinateFileKind;
using diskutil::DataFormat;
using diskutil::DrivePathType;
using diskutil::getDrivePathType;
using diskutil::getTrajectoryFormat;
using numerics::default_trajpos_scale_bits;
using stmath::roundUp;
  
struct ValidCoordinateTypes {
 
  ValidCoordinateTypes(size_t double_id_in, size_t float_id_in, size_t short_id_in,
                       size_t int_id_in, size_t llint_id_in);

  ValidCoordinateTypes(const ValidCoordinateTypes &original) = default;
  ValidCoordinateTypes(ValidCoordinateTypes &&original) = default;
  ValidCoordinateTypes& operator=(const ValidCoordinateTypes &original) = default;
  ValidCoordinateTypes& operator=(ValidCoordinateTypes &&original) = default;
  
  size_t double_id;  
  size_t float_id;   
  size_t short_id;   
  size_t int_id;     
  size_t llint_id;   
};
  
template <typename T> struct CoordinateSeriesWriter {

  CoordinateSeriesWriter(int natom_in, int nframe_in, UnitCellType unit_cell_in, int gpos_bits_in,
                         double gpos_scale_in, double inv_gpos_scale_in, T* xcrd_in, T* ycrd_in,
                         T* zcrd_in, double* umat_in, double* invu_in, double* boxdim_in);

  CoordinateSeriesWriter(const CoordinateSeriesWriter &original) = default;
  CoordinateSeriesWriter(CoordinateSeriesWriter &&original) = default;
  CoordinateSeriesWriter& operator=(const CoordinateSeriesWriter &other) = default;
 
  const int natom;               
  const int nframe;              
  const UnitCellType unit_cell;  
  const int gpos_bits;           
  const double gpos_scale;       
  const double inv_gpos_scale;   
  T* xcrd;                       
  T* ycrd;                       
  T* zcrd;                       
  double* umat;                  
  double* invu;                  
  double* boxdim;                
};
  
template <typename T> struct CoordinateSeriesReader {

  CoordinateSeriesReader(int natom_in, int nframe_in, UnitCellType unit_cell_in, int gpos_bits_in,
                         double gpos_scale_in, double inv_gpos_scale_in, const T* xcrd_in,
                         const T* ycrd_in, const T* zcrd_in, const double* umat_in,
                         const double* invu_in, const double* boxdim_in);
  CoordinateSeriesReader(const CoordinateSeriesWriter<T> &csw);

  CoordinateSeriesReader(const CoordinateSeriesReader &original) = default;
  CoordinateSeriesReader(CoordinateSeriesReader &&original) = default;
  CoordinateSeriesReader& operator=(const CoordinateSeriesReader &other) = default;
 
  const int natom;               
  const int nframe;              
  const UnitCellType unit_cell;  
  const int gpos_bits;           
  const double gpos_scale;       
  const double inv_gpos_scale;   
  const T* xcrd;                 
  const T* ycrd;                 
  const T* zcrd;                 
  const double* umat;            
  const double* invu;            
  const double* boxdim;          
};

template <typename T>  class CoordinateSeries {
public:

  explicit CoordinateSeries(int natom_in = 0, int nframe_in = 0,
                            UnitCellType unit_cell_in = UnitCellType::NONE,
                            int globalpos_scale_bits_in = 0,
                            HybridFormat format_in = default_hpc_format);
 
  explicit CoordinateSeries(const std::string &file_name, int atom_count_in = 0,
                            CoordinateFileKind file_kind = CoordinateFileKind::UNKNOWN,
                            const std::vector<int> &frame_numbers = {},
                            int replica_count = 1, UnitCellType unit_cell_in = UnitCellType::NONE,
                            int globalpos_scale_bits_in = 0,
                            HybridFormat format_in = default_hpc_format);
 
  explicit CoordinateSeries(PhaseSpace *ps, int nframe_in, int globalpos_scale_bits_in = 0,
                            HybridFormat format_in = default_hpc_format);
 
  explicit CoordinateSeries(const PhaseSpace &ps, int nframe_in, int globalpos_scale_bits_in = 0,
                            HybridFormat format_in = default_hpc_format);
 
  explicit CoordinateSeries(CoordinateFrame *cf, int nframe_in, int globalpos_scale_bits_in = 0,
                            HybridFormat format_in = default_hpc_format);
 
  explicit CoordinateSeries(const CoordinateFrame &cf, int nframe_in,
                            int globalpos_scale_bits_in = 0,
                            HybridFormat format_in = default_hpc_format);

  CoordinateSeries(const CoordinateSeries &original) = default;
  CoordinateSeries(CoordinateSeries &&original) = default;
  CoordinateSeries& operator=(const CoordinateSeries &other) = default;
  CoordinateSeries& operator=(CoordinateSeries &&other) = default;

  explicit CoordinateSeries(const CoordinateSeries &original, HybridFormat format_in);
  
  template <typename Toriginal>
  CoordinateSeries(const CoordinateSeries<Toriginal> &original, int globalpos_scale_bits_in = 0,
                   HybridFormat format_in = default_hpc_format);

  HybridFormat getFormat() const;
  
  int getAtomCount() const;

  int getFrameCount() const;

  int getFrameCapacity() const;

  UnitCellType getUnitCellType() const;

  int getFixedPrecisionBits() const;
  
  template <typename Treport> std::vector<Treport>
  getInterlacedCoordinates(int frame_index, int globalpos_bits_out = -1,
                           HybridTargetLevel tier = HybridTargetLevel::HOST) const;
 
  template <typename Treport> std::vector<Treport>
  getInterlacedCoordinates(int frame_index, int low_index, int high_index,
                           int globalpos_bits_out = -1,
                           HybridTargetLevel tier = HybridTargetLevel::HOST) const;

  std::vector<double> getBoxSpaceTransform(int frame_index,
                                           HybridTargetLevel tier = HybridTargetLevel::HOST) const;

  std::vector<double> getInverseTransform(int frame_index,
                                          HybridTargetLevel tier = HybridTargetLevel::HOST) const;

  std::vector<double> getBoxDimensions(int frame_index,
                                       HybridTargetLevel tier = HybridTargetLevel::HOST) const;
 
  const Hybrid<double>& getBoxDimensions() const;

  void extractFrame(CoordinateFrame *cf, size_t frame_index, const GpuDetails &gpu) const;
 
  void extractFrame(PhaseSpace *ps, size_t frame_index, TrajectoryKind kind,
                    CoordinateCycle time_point, const GpuDetails &gpu) const;
 
  void extractFrame(PhaseSpace *ps, size_t frame_index, TrajectoryKind kind,
                    const GpuDetails &gpu) const;
 
  void extractFrame(PhaseSpace *ps, size_t frame_index, const GpuDetails &gpu) const;
 
  void extractFrame(CoordinateFrame *cf, size_t frame_index, HybridTargetLevel tier) const;
 
  void extractFrame(PhaseSpace *ps, size_t frame_index, TrajectoryKind kind,
                    CoordinateCycle time_point, HybridTargetLevel tier) const;
 
  void extractFrame(PhaseSpace *ps, size_t frame_index, TrajectoryKind kind,
                    HybridTargetLevel tier) const;
 
  void extractFrame(PhaseSpace *ps, size_t frame_index, HybridTargetLevel tier) const;

  CoordinateFrame exportFrame(size_t frame_index,
                              HybridTargetLevel tier = HybridTargetLevel::HOST) const;

  PhaseSpace exportPhaseSpace(size_t frame_index,
                              HybridTargetLevel tier = HybridTargetLevel::HOST) const;
  
  void exportToFile(const std::string &file_name,
                    CoordinateFileKind output_kind = CoordinateFileKind::AMBER_CRD,
                    PrintSituation expectation = PrintSituation::UNKNOWN, int low_index = 0,
                    int high_index = 0, HybridTargetLevel tier = HybridTargetLevel::HOST) const;

  const Hybrid<T>& getCoordinateReference(CartesianDimension dim) const;

  const Hybrid<T>* getCoordinatePointer(CartesianDimension dim) const;

  const Hybrid<double>& getBoxTransforms() const;

  const Hybrid<double>* getBoxTransformPointer() const;

  const Hybrid<double>& getInverseTransforms() const;

  const Hybrid<double>* getInverseTransformPointer() const;

  const Hybrid<double>* getBoxDimensionPointer() const;

  const Hybrid<T>* getFramesHandle(CartesianDimension dim) const;
  Hybrid<T>* getFramesHandle(CartesianDimension dim);

  const Hybrid<double>* getBoxTransformsHandle() const;
  Hybrid<double>* getBoxTransformsHandle();

  const Hybrid<double>* getInverseTransformsHandle() const;
  Hybrid<double>* getInverseTransformsHandle();

  const Hybrid<double>* getBoxDimensionsHandle() const;
  Hybrid<double>* getBoxDimensionsHandle();
  
  const CoordinateSeries<T>* getSelfPointer() const;

  const CoordinateSeriesReader<T> data(HybridTargetLevel tier = HybridTargetLevel::HOST) const;
  CoordinateSeriesWriter<T> data(HybridTargetLevel tier = HybridTargetLevel::HOST);

  const CoordinateSeriesReader<void>
  templateFreeData(HybridTargetLevel tier = HybridTargetLevel::HOST) const;
 
  CoordinateSeriesWriter<void> templateFreeData(HybridTargetLevel tier = HybridTargetLevel::HOST);
  
#ifdef STORMM_USE_HPC
  const CoordinateSeriesReader<T> deviceViewToHostData() const;
  CoordinateSeriesWriter<T> deviceViewToHostData();

  const CoordinateSeriesReader<void> deviceViewToTemplateFreeHostData() const;
  CoordinateSeriesWriter<void> deviceViewToTemplateFreeHostData();
  
  void upload();

  void download();
#endif
  
  void importFrame(const CoordinateFrameReader &cfr, int atom_start, int atom_end,
                   int frame_index = -1);
  void importFrame(const CoordinateFrameReader &cfr, int frame_index = -1);
  void importFrame(const CoordinateFrameWriter &cfw, int atom_start, int atom_end,
                   int frame_index = -1);
  void importFrame(const CoordinateFrameWriter &cfw, int frame_index = -1);
  void importFrame(const CoordinateFrame &cf, int atom_start, int atom_end, int frame_index = -1);
  void importFrame(const CoordinateFrame &cf, int frame_index = -1);
  void importFrame(const CoordinateFrame *cf, int atom_start, int atom_end, int frame_index = -1);
  void importFrame(const CoordinateFrame *cf, int frame_index = -1);
  void importFrame(const PhaseSpace &ps, int atom_start, int atom_end, int frame_index = -1,
                   TrajectoryKind kind = TrajectoryKind::POSITIONS,
                   CoordinateCycle orientation = CoordinateCycle::WHITE);
  void importFrame(const PhaseSpace &ps, int frame_index = -1,
                   TrajectoryKind kind = TrajectoryKind::POSITIONS,
                   CoordinateCycle orientation = CoordinateCycle::WHITE);
  void importFrame(const PhaseSpace *ps, int atom_start, int atom_end, int frame_index = -1,
                   TrajectoryKind kind = TrajectoryKind::POSITIONS,
                   CoordinateCycle orientation = CoordinateCycle::WHITE);
  void importFrame(const PhaseSpace *ps, int frame_index = -1,
                   TrajectoryKind kind = TrajectoryKind::POSITIONS,
                   CoordinateCycle orientation = CoordinateCycle::WHITE);

  void importFromFile(const std::string &file_name,                       
                      CoordinateFileKind file_kind = CoordinateFileKind::UNKNOWN,
                      const std::vector<int> &frame_numbers = {}, int replica_count = 1,
                      int frame_index_start = -1);

  void reserve(const int new_frame_capacity);

  void shrinkToFit();

  void resize(int new_frame_count);
  void resize(int new_frame_count, const CoordinateFrameReader &cfr, int atom_start = 0,
              int atom_end = 0);
  void resize(int new_frame_count, const CoordinateFrameWriter &cfw, int atom_start = 0,
              int atom_end = 0);
  void resize(int new_frame_count, const CoordinateFrame &cf, int atom_start = 0,
              int atom_end = 0);
  void resize(int new_frame_count, CoordinateFrame *cf, int atom_start = 0, int atom_end = 0);
  void resize(int new_frame_count, const PhaseSpace &ps, int atom_start = 0, int atom_end = 0);
  void resize(int new_frame_count, PhaseSpace *ps, int atom_start = 0, int atom_end = 0);
  
  void pushBack(const CoordinateFrameReader &cfr, int atom_start = 0, int atom_end = 0);
  void pushBack(const CoordinateFrameWriter &cfw, int atom_start = 0, int atom_end = 0);
  void pushBack(const CoordinateFrame &cf, int atom_start = 0, int atom_end = 0);
  void pushBack(CoordinateFrame *cf, int atom_start = 0, int atom_end = 0);
  void pushBack(const PhaseSpace &ps, int atom_start = 0, int atom_end = 0);
  void pushBack(PhaseSpace *ps, int atom_start = 0, int atom_end = 0);
  
private:
  HybridFormat format;                  
  int atom_count;                       
  int frame_count;                      
  int frame_capacity;                   
  int globalpos_scale_bits;             
  UnitCellType unit_cell;               
  double globalpos_scale;               
  double inverse_globalpos_scale;       
  Hybrid<T> x_coordinates;              
  Hybrid<T> y_coordinates;              
  Hybrid<T> z_coordinates;              
  Hybrid<double> box_space_transforms;  
  Hybrid<double> inverse_transforms;    
  Hybrid<double> box_dimensions;        

  void allocate(int new_frame_capacity);
};

template <typename T> CoordinateSeriesWriter<T> restoreType(CoordinateSeriesWriter<void> *rasa);
 
template <typename T>
CoordinateSeriesWriter<T> restoreType(const CoordinateSeriesWriter<void> &rasa);
 
template <typename T>
const CoordinateSeriesReader<T> restoreType(const CoordinateSeriesReader<void> *rasa);
 
template <typename T>
const CoordinateSeriesReader<T> restoreType(const CoordinateSeriesReader<void> &rasa);

template <typename Torig, typename Tnew>
CoordinateSeries<Tnew> changeCoordinateSeriesType(const CoordinateSeries<Torig> &cs,
                                                  int globalpos_scale_bits_in);
 
} // namespace trajectory 
} // namespace stormm
 
#include "coordinate_series.tpp"
 
#endif