c_loops.hh

Go to the documentation of this file.

// Copyright Epic Games, Inc. All Rights Reserved.
// Modified version of Voro++'s source file
 
// Voro++, a 3D cell-based Voronoi library
//
// Author   : Chris H. Rycroft (LBL / UC Berkeley)
// Email    : chr@alum.mit.edu
// Date     : August 30th 2011
 
#ifndef VOROPP_C_LOOPS_HH
#define VOROPP_C_LOOPS_HH
 
#include "Misc/AssertionMacros.h"
 
#include <cstddef>
 
#include "config.hh"
 
namespace voro {
 
enum c_loop_subset_mode {
    sphere,
    box,
    no_check
};
 
class particle_order {
    public:
        int *o;
        int *op;
        int size;
        particle_order(int init_size=init_ordering_size)
            : o(new int[size_t(init_size)<<1]),op(o),size(init_size) {}
        ~particle_order() {
            delete [] o;
        }
        inline void add(int ijk,int q) {
            if(op==o+size) add_ordering_memory();
            *(op++)=ijk;*(op++)=q;
        }
    private:
        void add_ordering_memory();
};
 
class c_loop_base {
    public:
        const int nx;
        const int ny;
        const int nz;
        const int nxy;
        const int nxyz;
        const int ps;
        double **p;
        int **id;
        int *co;
        int i;
        int j;
        int k;
        int ijk;
        int q;
        template<class c_class>
        c_loop_base(c_class &con) : nx(con.nx), ny(con.ny), nz(con.nz),
                        nxy(con.nxy), nxyz(con.nxyz), ps(con.ps),
                        p(con.p), id(con.id), co(con.co) {}
        inline void pos(double &x,double &y,double &z) {
            double *pp=p[ijk]+ps*q;
            x=*(pp++);y=*(pp++);z=*pp;
        }
        inline void pos(int &pid,double &x,double &y,double &z,double &r) {
            pid=id[ijk][q];
            double *pp=p[ijk]+ps*q;
            x=*(pp++);y=*(pp++);z=*pp;
            r=ps==3?default_radius:*(++pp);
        }
        inline double x() {return p[ijk][ps*q];}
        inline double y() {return p[ijk][ps*q+1];}
        inline double z() {return p[ijk][ps*q+2];}
        inline int pid() {return id[ijk][q];}
};
 
class c_loop_all : public c_loop_base {
    public:
        template<class c_class>
        c_loop_all(c_class &con) : c_loop_base(con) {}
        inline bool start() {
            i=j=k=ijk=q=0;
            while(co[ijk]==0) if(!next_block()) return false;
            return true;
        }
        inline bool inc() {
            q++;
            if(q>=co[ijk]) {
                q=0;
                do {
                    if(!next_block()) return false;
                } while(co[ijk]==0);
            }
            return true;
        }
    private:
        inline bool next_block() {
            ijk++;
            i++;
            if(i==nx) {
                i=0;j++;
                if(j==ny) {
                    j=0;k++;
                    if(ijk==nxyz) return false;
                }
            }
            return true;
        }
};
 
class c_loop_block_range : public c_loop_base {
    public:
        template<class c_class>
        c_loop_block_range(c_class &con) : c_loop_base(con) {}
        inline bool start() {
            ijk = ijk_start;
            checkSlow(ijk >= 0 && ijk < nxyz);
            k = ijk / nxy;
            int ijkt = ijk - nxy * k;
            j = ijkt / nx;
            i = ijkt - j * nx;
            q=0;
            while(co[ijk]==0) if(!next_block()) return false;
            return true;
        }
        void setup_range(int ijk_start_in, int ijk_stop_in)
        {
            ijk_start = ijk_start_in;
            ijk_stop = ijk_stop_in;
 
            check(ijk_start >= 0 && ijk_start < nxyz && ijk_stop >= 0 && ijk_stop <= nxyz);
            checkSlow(ijk_start < ijk_stop);
        }
        inline bool inc() {
            q++;
            if(q>=co[ijk]) {
                q=0;
                do {
                    if(!next_block()) return false;
                } while(co[ijk]==0);
            }
            return true;
        }
    private:
        int ijk_start = 0, ijk_stop = 0;
 
        inline bool next_block() {
            ijk++;
            if (ijk >= ijk_stop) return false;
            i++;
            if(i==nx) {
                i=0;j++;
                if(j==ny) {
                    j=0;k++;
                }
            }
            return true;
        }
};
 
class c_loop_subset : public c_loop_base {
    public:
        c_loop_subset_mode mode;
        template<class c_class>
        c_loop_subset(c_class &con) : c_loop_base(con), ax(con.ax), ay(con.ay), az(con.az),
            sx(con.bx-ax), sy(con.by-ay), sz(con.bz-az), xsp(con.xsp), ysp(con.ysp), zsp(con.zsp),
            xperiodic(con.xperiodic), yperiodic(con.yperiodic), zperiodic(con.zperiodic) {}
        void setup_sphere(double vx,double vy,double vz,double r,bool bounds_test=true);
        void setup_box(double xmin,double xmax,double ymin,double ymax,double zmin,double zmax,bool bounds_test=true);
        void setup_intbox(int ai_,int bi_,int aj_,int bj_,int ak_,int bk_);
        bool start();
        inline bool inc() {
            do {
                q++;
                while(q>=co[ijk]) {q=0;if(!next_block()) return false;}
            } while(mode!=no_check&&out_of_bounds());
            return true;
        }
    private:
        const double ax,ay,az,sx,sy,sz,xsp,ysp,zsp;
        double px,py,pz,apx,apy,apz;
        double v0,v1,v2,v3,v4,v5;
        int ai,bi,aj,bj,ak,bk;
        int ci,cj,ck,di,dj,dk,inc1,inc2;
        const bool xperiodic,yperiodic,zperiodic;
        inline int step_mod(int a,int b) const {return a>=0?a%b:b-1-(b-1-a)%b;}
        inline int step_div(int a,int b) const {return a>=0?a/b:-1+(a+1)/b;}
        inline int step_int(double a) const {return a<0?int(a)-1:int(a);}
        void setup_common();
        bool next_block();
        bool out_of_bounds();
};
 
class c_loop_order : public c_loop_base {
    public:
        particle_order &vo;
        int *cp;
        int *op;
        template<class c_class>
        c_loop_order(c_class &con,particle_order &vo_)
        : c_loop_base(con), vo(vo_), nx(con.nx), nxy(con.nxy) {}
        inline bool start() {
            cp=vo.o;op=vo.op;
            if(cp!=op) {
                ijk=*(cp++);decode();
                q=*(cp++);
                return true;
            } else return false;
        }
        inline bool inc() {
            if(cp==op) return false;
            ijk=*(cp++);decode();
            q=*(cp++);
            return true;
        }
    private:
        const int nx;
        const int nxy;
        inline void decode() {
            k=ijk/nxy;
            int ijkt=ijk-nxy*k;
            j=ijkt/nx;
            i=ijkt-j*nx;
        }
};
 
class c_loop_all_periodic : public c_loop_base {
    public:
        template<class c_class>
        c_loop_all_periodic(c_class &con) : c_loop_base(con), ey(con.ey), ez(con.ez), wy(con.wy), wz(con.wz),
            ijk0(nx*(ey+con.oy*ez)), inc2(2*nx*con.ey+1) {}
        inline bool start() {
            i=0;
            j=ey;
            k=ez;
            ijk=ijk0;
            q=0;
            while(co[ijk]==0) if(!next_block()) return false;
            return true;
        }
        inline bool inc() {
            q++;
            if(q>=co[ijk]) {
                q=0;
                do {
                    if(!next_block()) return false;
                } while(co[ijk]==0);
            }
            return true;
        }
    private:
        int ey;
        int ez;
        int wy;
        int wz;
        int ijk0;
        int inc2;
        inline bool next_block() {
            i++;
            if(i==nx) {
                i=0;j++;
                if(j==wy) {
                    j=ey;k++;
                    if(k==wz) return false;
                    ijk+=inc2;
                } else ijk++;
            } else ijk++;
            return true;
        }
};
 
class c_loop_order_periodic : public c_loop_base {
    public:
        particle_order &vo;
        int *cp;
        int *op;
        template<class c_class>
        c_loop_order_periodic(c_class &con,particle_order &vo_)
        : c_loop_base(con), vo(vo_), nx(con.nx), oxy(con.nx*con.oy) {}
        inline bool start() {
            cp=vo.o;op=vo.op;
            if(cp!=op) {
                ijk=*(cp++);decode();
                q=*(cp++);
                return true;
            } else return false;
        }
        inline bool inc() {
            if(cp==op) return false;
            ijk=*(cp++);decode();
            q=*(cp++);
            return true;
        }
    private:
        const int nx;
        const int oxy;
        inline void decode() {
            k=ijk/oxy;
            int ijkt=ijk-oxy*k;
            j=ijkt/nx;
            i=ijkt-j*nx;
        }
};
 
}
 
#endif