pmd.hpp

#include <vector>
#include <array>
#include <cmath>
#include <algorithm>
#include "mpi.h"
 
/* Constants for the random number generator */
#define D2P31M 2147483647.0
#define DMUL 16807.0
 
class Atom {
public:
    double type;      
    bool isResident; 
 
    double x;        
    double y;        
    double z;        
 
    double ax;        
    double ay;        
    double az;        
 
    double vx;        
    double vy;        
    double vz;        
 
    std::array<int, 6> shiftCount; 
 
    Atom(); 
};
 
class SubSystem {
public:
    int pid; 
    int n; 
    int nglob; 
    double comt; 
    std::array<double, 3> al; 
    std::array<int, 3> vid; 
    std::array<int, 3> myparity; 
    std::array<int, 6> nn; 
    std::vector<std::vector<double> > sv; 
    std::array<double, 3> vSum, gvSum;
    std::vector<Atom> atoms; 
 
    std::array<int, 3> vproc{}, 
    InitUcell{}; 
    double Density, 
    InitTemp, 
    DeltaT; 
    int StepLimit, 
    StepAvg; 
 
    double kinEnergy, 
    potEnergy, 
    totEnergy, 
    temperature; 
 
    SubSystem();
 
    void InitNeighborNode(std::array<int, 3> vproc);
 
    void Update(double DeltaT);
 
    void Kick(double DeltaT);
 
 
    void AtomCopy();
 
    std::vector<int> AtomMove(bool);
 
    void ShiftAtoms();
 
    void WrapAtoms();
 
    int bbd(Atom atom, int ku);
 
    int bmv(Atom atom, int ku);
 
    void EvalProps(int stepCount);
 
    void WriteXYZ(int step);
 
    static double Dmod(double a, double b) {
        int n;
        n = (int) (a / b);
        return (a - b * n);
    }
 
    static double RandR(double *seed) {
        *seed = Dmod(*seed * DMUL, D2P31M);
        return (*seed / D2P31M);
    }
 
    static void RandVec3(double *p, double *seed) {
        double x = 0, y = 0, s = 2.0;
        while (s > 1.0) {
            x = 2.0 * RandR(seed) - 1.0;
            y = 2.0 * RandR(seed) - 1.0;
            s = x * x + y * y;
        }
        p[2] = 1.0 - 2.0 * s;
        s = 2.0 * sqrt(1.0 - s);
        p[0] = s * x;
        p[1] = s * y;
    }
 
};