motion_sweeper.h

// -*-c++-*-
#ifndef STORMM_MOTION_SWEEPER_H
#define STORMM_MOTION_SWEEPER_H
 
#include "copyright.h"
#include "Accelerator/hybrid.h"
#include "Constants/fixed_precision.h"
#include "DataTypes/stormm_vector_types.h"
#include "Synthesis/phasespace_synthesis.h"
 
namespace stormm {
namespace trajectory {
 
using card::HybridTargetLevel;
using numerics::default_momentum_scale_bits;
using numerics::default_com_scale_bits;
using numerics::default_inertia_scale_bits;
using synthesis::PhaseSpaceSynthesis;

struct MotionSweepWriter {
public:

  MotionSweepWriter(int nwu_in, const int4* work_units_in, const double* total_mass_in,
                    double com_scale_in, double mv_scale_in, double inrt_scale_in, llint* xcom_in,
                    llint* ycom_in, llint* zcom_in, llint* xcom_nxt_in, llint* ycom_nxt_in,
                    llint* zcom_nxt_in, int* xcom_ovrf_in, int* ycom_ovrf_in, int* zcom_ovrf_in,
                    int* xcom_nxt_ovrf_in, int* ycom_nxt_ovrf_in, int* zcom_nxt_ovrf_in,
                    llint* xmv_in, llint* ymv_in, llint* zmv_in, llint* xmv_nxt_in,
                    llint* ymv_nxt_in, llint* zmv_nxt_in, int* xmv_ovrf_in, int* ymv_ovrf_in,
                    int* zmv_ovrf_in, int* xmv_nxt_ovrf_in, int* ymv_nxt_ovrf_in,
                    int* zmv_nxt_ovrf_in, llint* rxmv_in, llint* rymv_in, llint* rzmv_in,
                    llint* rxmv_nxt_in, llint* rymv_nxt_in, llint* rzmv_nxt_in, int* rxmv_ovrf_in,
                    int* rymv_ovrf_in, int* rzmv_ovrf_in, int* rxmv_nxt_ovrf_in,
                    int* rymv_nxt_ovrf_in, int* rzmv_nxt_ovrf_in, llint* inrt_in,
                    llint* inrt_nxt_in, int* inrt_ovrf_in, int* inrt_nxt_ovrf_in);

  MotionSweepWriter(const MotionSweepWriter &original) = default;
  MotionSweepWriter(MotionSweepWriter &&original) = default;
  
  const int nwu;             
  const int4* work_units;    
  const double* total_mass;  
  const double com_scale;    
  const double mv_scale;     
  const double inrt_scale;   
 
  // Each of the following arrays has a dimension equal to the number of systems in the associated
  // synthesis.  The number of systems is not stored in this abstract as the correct index of each
  // array can be found in the work unit itself.
  llint* xcom;         
  llint* ycom;         
  llint* zcom;         
  llint* xcom_nxt;     
  llint* ycom_nxt;     
  llint* zcom_nxt;     
  int* xcom_ovrf;      
  int* ycom_ovrf;      
  int* zcom_ovrf;      
  int* xcom_nxt_ovrf;  
  int* ycom_nxt_ovrf;  
  int* zcom_nxt_ovrf;  
  llint* xmv;          
  llint* ymv;          
  llint* zmv;          
  llint* xmv_nxt;      
  llint* ymv_nxt;      
  llint* zmv_nxt;      
  int* xmv_ovrf;       
  int* ymv_ovrf;       
  int* zmv_ovrf;       
  int* xmv_nxt_ovrf;   
  int* ymv_nxt_ovrf;   
  int* zmv_nxt_ovrf;   
  llint* rxmv;         
  llint* rymv;         
  llint* rzmv;         
  llint* rxmv_nxt;     
  llint* rymv_nxt;     
  llint* rzmv_nxt;     
  int* rxmv_ovrf;      
  int* rymv_ovrf;      
  int* rzmv_ovrf;      
  int* rxmv_nxt_ovrf;  
  int* rymv_nxt_ovrf;  
  int* rzmv_nxt_ovrf;  
  llint* inrt;         
  llint* inrt_nxt;     
  int* inrt_ovrf;      
  int* inrt_nxt_ovrf;  
};

struct MotionSweepReader {
public:

  MotionSweepReader(int nwu_in, const int4* work_units_in, const double* total_mass_in,
                    double com_scale_in, double mv_scale_in, double inrt_scale,
                    const llint* xcom_in, const llint* ycom_in, const llint* zcom_in,
                    const llint* xcom_nxt_in, const llint* ycom_nxt_in, const llint* zcom_nxt_in,
                    const int* xcom_ovrf_in, const int* ycom_ovrf_in, const int* zcom_ovrf_in,
                    const int* xcom_nxt_ovrf_in, const int* ycom_nxt_ovrf_in,
                    const int* zcom_nxt_ovrf_in, const llint* xmv_in, const llint* ymv_in,
                    const llint* zmv_in, const llint* xmv_nxt_in, const llint* ymv_nxt_in,
                    const llint* zmv_nxt_in, const int* xmv_ovrf_in, const int* ymv_ovrf_in,
                    const int* zmv_ovrf_in, const int* xmv_nxt_ovrf_in, const int* ymv_nxt_ovrf_in,
                    const int* zmv_nxt_ovrf_in, const llint* rxmv_in, const llint* rymv_in,
                    const llint* rzmv_in, const llint* rxmv_nxt_in, const llint* rymv_nxt_in,
                    const llint* rzmv_nxt_in, const int* rxmv_ovrf_in, const int* rymv_ovrf_in,
                    const int* rzmv_ovrf_in, const int* rxmv_nxt_ovrf_in,
                    const int* rymv_nxt_ovrf_in, const int* rzmv_nxt_ovrf_in, const llint* inrt_in,
                    const llint* inrt_nxt_in, const int* inrt_ovrf_in,
                    const int* inrt_nxt_ovrf_in);
 
  MotionSweepReader(const MotionSweepWriter *w);
  
  MotionSweepReader(const MotionSweepWriter &w);

  MotionSweepReader(const MotionSweepReader &original) = default;
  MotionSweepReader(MotionSweepReader &&original) = default;
  
  const int nwu;             
  const int4* work_units;    
  const double* total_mass;  
  const double com_scale;    
  const double mv_scale;     
  const double inrt_scale;   
 
  // Each of the following arrays has a dimension equal to the number of systems in the associated
  // synthesis.  The number of systems is not stored in this abstract as the correct index of each
  // array can be found in the work unit itself.
  const llint* xcom;         
  const llint* ycom;         
  const llint* zcom;         
  const llint* xcom_nxt;     
  const llint* ycom_nxt;     
  const llint* zcom_nxt;     
  const int* xcom_ovrf;      
  const int* ycom_ovrf;      
  const int* zcom_ovrf;      
  const int* xcom_nxt_ovrf;  
  const int* ycom_nxt_ovrf;  
  const int* zcom_nxt_ovrf;  
  const llint* xmv;          
  const llint* ymv;          
  const llint* zmv;          
  const llint* xmv_nxt;      
  const llint* ymv_nxt;      
  const llint* zmv_nxt;      
  const int* xmv_ovrf;       
  const int* ymv_ovrf;       
  const int* zmv_ovrf;       
  const int* xmv_nxt_ovrf;   
  const int* ymv_nxt_ovrf;   
  const int* zmv_nxt_ovrf;   
  const llint* rxmv;         
  const llint* rymv;         
  const llint* rzmv;         
  const llint* rxmv_nxt;     
  const llint* rymv_nxt;     
  const llint* rzmv_nxt;     
  const int* rxmv_ovrf;      
  const int* rymv_ovrf;      
  const int* rzmv_ovrf;      
  const int* rxmv_nxt_ovrf;  
  const int* rymv_nxt_ovrf;  
  const int* rzmv_nxt_ovrf;  
  const llint* inrt;         
  const llint* inrt_nxt;     
  const int* inrt_ovrf;      
  const int* inrt_nxt_ovrf;  
};

class MotionSweeper {
public:

  MotionSweeper(const PhaseSpaceSynthesis *poly_ps,
                int momentum_bit_count_in = default_momentum_scale_bits,
                int center_of_mass_bit_count_in = default_com_scale_bits,
                int inertia_bit_count_in = default_inertia_scale_bits);
 
  MotionSweeper(const PhaseSpaceSynthesis &poly_ps,
                int momentum_bit_count_in = default_momentum_scale_bits,
                int center_of_mass_bit_count_in = default_com_scale_bits,
                int inertia_bit_count_in = default_inertia_scale_bits);

  MotionSweeper(const MotionSweeper &original);
  MotionSweeper(MotionSweeper &&original);
  MotionSweeper& operator=(const MotionSweeper &other);  
  MotionSweeper& operator=(MotionSweeper &&other);  
  
  int getSystemCount() const;

  int getWorkUnitCount() const;

  CoordinateCycle getCyclePosition() const;

  int getMomentumBitCount() const;

  int getCenterOfMassBitCount() const;

  int getInertialTensorBitCount() const;

  double3 getNetVelocity(int idx, HybridTargetLevel tier = HybridTargetLevel::HOST) const;

  double3 getCenterOfMass(int idx, HybridTargetLevel tier = HybridTargetLevel::HOST) const;

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

  double3 getAngularVelocity(int idx, HybridTargetLevel tier = HybridTargetLevel::HOST) const;

  void updateCyclePosition();
  void updateCyclePosition(CoordinateCycle time_point);

  MotionSweepWriter data(CoordinateCycle orientation,
                         HybridTargetLevel tier = HybridTargetLevel::HOST);
 
  MotionSweepWriter data(HybridTargetLevel tier = HybridTargetLevel::HOST);
 
  const MotionSweepReader data(CoordinateCycle orientation,
                               HybridTargetLevel tier = HybridTargetLevel::HOST) const;
 
  const MotionSweepReader data(HybridTargetLevel tier = HybridTargetLevel::HOST) const;
  
#ifdef STORMM_USE_HPC
  void uploadWorkUnits();

  void downloadWorkUnits();
  
  void uploadAll();
  
  void downloadAll();
#endif
  
private:
  
  CoordinateCycle cycle_position;   
  int center_of_mass_bit_count;     
  int momentum_bit_count;           
  int inertia_bit_count;            
  Hybrid<llint> xcom_white;         
  Hybrid<llint> ycom_white;         
  Hybrid<llint> zcom_white;         
  Hybrid<llint> xcom_black;         
  Hybrid<llint> ycom_black;         
  Hybrid<llint> zcom_black;         
  Hybrid<int> xcom_overflow_white;  
  Hybrid<int> ycom_overflow_white;  
  Hybrid<int> zcom_overflow_white;  
  Hybrid<int> xcom_overflow_black;  
  Hybrid<int> ycom_overflow_black;  
  Hybrid<int> zcom_overflow_black;  
  Hybrid<double> total_mass;        
  Hybrid<llint> xmv_white;          
  Hybrid<llint> ymv_white;          
  Hybrid<llint> zmv_white;          
  Hybrid<llint> xmv_black;          
  Hybrid<llint> ymv_black;          
  Hybrid<llint> zmv_black;          
  Hybrid<int> xmv_overflow_white;   
  Hybrid<int> ymv_overflow_white;   
  Hybrid<int> zmv_overflow_white;   
  Hybrid<int> xmv_overflow_black;   
  Hybrid<int> ymv_overflow_black;   
  Hybrid<int> zmv_overflow_black;   
  Hybrid<llint> rxmv_white;         
  Hybrid<llint> rymv_white;         
  Hybrid<llint> rzmv_white;         
  Hybrid<llint> rxmv_black;         
  Hybrid<llint> rymv_black;         
  Hybrid<llint> rzmv_black;         
  Hybrid<int> rxmv_overflow_white;  
  Hybrid<int> rymv_overflow_white;  
  Hybrid<int> rzmv_overflow_white;  
  Hybrid<int> rxmv_overflow_black;  
  Hybrid<int> rymv_overflow_black;  
  Hybrid<int> rzmv_overflow_black;  

  Hybrid<llint> inertial_tensor_white;
  Hybrid<llint> inertial_tensor_black;

  Hybrid<int> inertial_tensor_overflow_white;
  Hybrid<int> inertial_tensor_overflow_black;

  Hybrid<llint> llint_data;

  Hybrid<int> int_data;
  
  int work_unit_count;
  
  Hybrid<int4> work_units;

  PhaseSpaceSynthesis *poly_ps_ptr;
 
  // Allocate memory in ARRAY-kind Hybrid objects targeted by POINTER-kind Hybrid accumulators.
  // Memory for work units and total system masses is allocated separately in the constructor.
  void allocate();
};
  
} // namespace trajectory
} // namespace stormm
#endif