forked from ALaDyn/piccante
-
Notifications
You must be signed in to change notification settings - Fork 0
/
particle_species.h
178 lines (150 loc) · 6.88 KB
/
particle_species.h
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
/* Copyright 2014 - Andrea Sgattoni, Luca Fedeli, Stefano Sinigardi */
/*******************************************************************************
This file is part of piccante.
piccante is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
piccante is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with piccante. If not, see <http://www.gnu.org/licenses/>.
*******************************************************************************/
#ifndef __PARTICLE_SPECIES_H__
#define __PARTICLE_SPECIES_H__
#define _USE_MATH_DEFINES
#include <mpi.h>
#include "commons.h"
#include "structures.h"
#include "grid.h"
#include "current.h"
#include "em_field.h"
#define _VERY_SMALL_MOMENTUM 1.0e-5
class SPECIE{
public:
static const int allocsize = 1000;
particlesType type;
bool allocated;
int allocatedNp, Np, Ncomp;
double Z, A; //atomic number, mass number
double chargeSign, coupling; // charge/mass*(speed_of_light^2)
double Numerical2Physical_particles; // physical particles corresponding to one num. particle
double mass;
double minima[7], maxima[7]; //components minima and maxima
double totalMomentum[3], totalEnergy;
std::string name;
PLASMA plasma;
// double (*initial_profile)(double , double , double , PLASMA*);
SPECIEspectrum spectrum;
GRID *mygrid;
SPECIE();
SPECIE(GRID *grid);
~SPECIE();
void allocate_species();
void erase();
void reallocate_species();
SPECIE operator = (SPECIE &destro);
void creation();
void creationFromFile1D(std::string name);
void move_window();
void addMarker();
void setTestSpecies();
bool amIWithMarker();
void output(std::ofstream &ff);
static const int myWidth = 12;
static const int myNarrowWidth = 6;
std::string getName();
void init_output_diag(std::ofstream &ff);
void output_diag(int istep, std::ofstream &ff);
void init_output_extrems(std::ofstream &ff);
void output_extrems(int istep, std::ofstream &ff);
void init_output_stat(std::ofstream &fdiag, std::ofstream &fextrem);
void output_stat(int istep, std::ofstream &fdiag, std::ofstream &fextrem, std::ofstream &fspectrum);
void position_advance();
void position_pbc();
void position_parallel_pbc();
void position_obc();
void momenta_advance(EM_FIELD *ebfield);
void momentaStretchedAdvance(EM_FIELD *ebfield);
void momenta_advance_with_friction(EM_FIELD *ebfield, double lambda);
void current_deposition(CURRENT *current);
void add_momenta(double uxin, double uyin, double uzin);
void add_momenta(gsl_rng* ext_rng, double uxin, double uyin, double uzin, tempDistrib distribution);
void current_deposition_standard(CURRENT *current);
void currentStretchedDepositionStandard(CURRENT *current);
void density_deposition_standard(CURRENT *current);
void densityStretchedDepositionStandard(CURRENT *current);
void setParticlesPerCellXYZ(int numX, int numY, int numZ);
void setName(std::string iname);
double getKineticEnergy();
void computeKineticEnergyWExtrems();
void outputSpectrum(std::ofstream &fspectrum);
void dump(std::ofstream &f);
void dumpBigBuffer(std::ofstream &f);
void debugDump(std::ofstream &f);
void reloadDump(std::ifstream &f);
void reloadBigBufferDump(std::ifstream &f);
bool areEnergyExtremesAvailable();
uint64_t printParticleNumber();
//PUBLIC INLINE FUNCTIONS
#ifdef _ACC_SINGLE_POINTER
inline double &ru(int c, int np) { return pData[c + np*Ncomp]; }
inline double &r0(int np) { return pData[np*Ncomp + 0]; }
inline double &r1(int np) { return pData[np*Ncomp + 1]; }
inline double &r2(int np) { return pData[np*Ncomp + 2]; }
inline double &u0(int np) { return pData[np*Ncomp + 3]; }
inline double &u1(int np) { return pData[np*Ncomp + 4]; }
inline double &u2(int np) { return pData[np*Ncomp + 5]; }
inline double &w(int np) { return pData[np*Ncomp + 6]; }
inline uint64_t &marker(int np) { return *((uint64_t*)(pData + (np*Ncomp + 7))); }
//inline uint64_t &marker(int np) { return *((uint64_t*)(&dummy)); }
#else
inline double &ru(int c, int np) { return val[c][np]; }
inline double &r0(int np) { return val[0][np]; }
inline double &r1(int np) { return val[1][np]; }
inline double &r2(int np) { return val[2][np]; }
inline double &u0(int np) { return val[3][np]; }
inline double &u1(int np) { return val[4][np]; }
inline double &u2(int np) { return val[5][np]; }
inline double &w(int np) { return val[6][np]; }
inline uint64_t &marker(int np) { return *((uint64_t*)(&val[7][np])); }
//inline uint64_t &marker(int np) { return *((uint64_t*)(&dummy)); }
#endif
private:
#ifdef _ACC_SINGLE_POINTER
double *pData;
#else
double **val;
#endif
bool isTestSpecies;
double dummy;
int valSize;
int particlePerCell;
int particlePerCellXYZ[3];
uint64_t lastParticle;
double savedExtrema[14];
double savedEnergy;
bool energyExtremesFlag;
bool flagWithMarker;
void callWaterbag(gsl_rng* ext_rng, double p0_x, double p0_y, double p0_z, double uxin, double uyin, double uzin);
void callUnifSphere(gsl_rng* ext_rng, double p0, double uxin, double uyin, double uzin);
void callSupergaussian(gsl_rng* ext_rng, double p0, double alpha, double uxin, double uyin, double uzin);
void callMaxwell(gsl_rng* ext_rng, double temp, double uxin, double uyin, double uzin);
void callJuttner(gsl_rng* ext_rng, double a, double uxin, double uyin, double uzin);
void callSpecial(gsl_rng* ext_rng, double Ta);
void computeParticleMassChargeCoupling();
void setNumberOfParticlePerCell();
void setLocalPlasmaMinimaAndMaxima(double *plasmarmin, double *plasmarmax);
uint64_t getSumNewParticlesOfAllPreviousProcessors(int number);
int getNumberOfParticlesWithin(double plasmarmin[3], double plasmarmax[3]);
int getNumberOfParticlesWithinFromFile1D(double plasmarmin[3], double plasmarmax[3], std::string name);
void createParticlesWithinFrom(double plasmarmin[3], double plasmarmax[3], int oldNumberOfParticles, uint64_t disp);
void createStretchedParticlesWithinFrom(double plasmarmin[3], double plasmarmax[3], int oldNumberOfParticles, uint64_t disp);
void createParticlesWithinFromButFromFile1D(double plasmarmin[3], double plasmarmax[3], int oldNumberOfParticles, uint64_t disp, std::string name);
void createStretchedParticlesWithinFromButFromFile1D(double plasmarmin[3], double plasmarmax[3], int oldNumberOfParticles, uint64_t disp, std::string name);
void computeLorentzMatrix(double ux, double uy, double uz, double matr[16]);
void debug_warning_particle_outside_boundaries(double x, double y, double z, int nump);
};
#endif /* _PARTICLE_SPECIES_H_ */