PeTar
N-body code for collisional gravitational systems
soft_ptcl.hpp
Go to the documentation of this file.
1 #pragma once
2 #include"ptcl.hpp"
3 
4 class ForceSoft{
5 public:
6  PS::F64vec acc;
7  PS::F64 pot;
8 #ifdef KDKDK_4TH
9  PS::F64vec acorr;
10 #endif
11 #ifdef SAVE_NEIGHBOR_ID_IN_FORCE_KERNEL
12  PS::S64 id_ngb[4];
13 #endif
14  PS::S64 n_ngb;
15  static PS::F64 grav_const;
16  void clear(){
17  acc = 0.0;
18  pot = 0.0;
19  n_ngb = 0;
20 #ifdef SAVE_NEIGHBOR_ID_IN_FORCE_KERNEL
21  id_ngb[0] = id_ngb[1] = id_ngb[2] = id_ngb[3] = 0;
22 #endif
23  }
24 };
25 
26 class FPSoft: public Ptcl{
27 public:
28  PS::F64vec acc; // soft
29 #ifdef KDKDK_4TH
30  PS::F64vec acorr;
31 #endif
32  PS::F64 pot_tot; // soft + hard
33  PS::F64 pot_soft; // soft only
34 #ifdef EXTERNAL_POT_IN_PTCL
35  PS::F64 pot_ext; // external potential
36 #endif
37 #ifdef SAVE_NEIGHBOR_ID_IN_FORCE_KERNEL
38  PS::S64 id_ngb[4];
39 #endif
40  PS::S64 n_ngb;
41  PS::S32 rank_org;
42  PS::S32 adr;
43 // static PS::F64 r_out;
44 
45  FPSoft() {}
46 
48 
50  PS::F64vec getPos() const{
51  return pos;
52  }
53 
55 
57  PS::F64vec getVel() const{
58  return vel;
59  }
60 
61  //template<class Tptcl>
62  //FPSoft& operator = (const FPSoft& p) {
63  // Ptcl::DataCopy(p);
64  // acc = p.acc;
65  // pot_tot = p.pot_tot;
66  // n_ngb = p.n_ngb;
67  // return *this;
68  //}
69 
70  template<class Tptcl>
71  FPSoft(const Tptcl& p, const PS::S32 rank_, const PS::S32 adr_) {
72  Ptcl::DataCopy(p);
73  rank_org = rank_;
74  adr = adr_;
75  acc = 0;
76  pot_tot = 0;
77  pot_soft= 0;
78 #ifdef EXTERNAL_POT_IN_PTCL
79  pot_ext = 0;
80 #endif
81  n_ngb = 0;
82  }
83 
84  void copyFromForce(const ForceSoft & force){
85  acc = force.acc;
86  pot_tot = force.pot;
87  pot_soft= pot_tot;
88 #ifdef KDKDK_4TH
89  acorr = force.acorr;
90 #endif
91 #ifdef SAVE_NEIGHBOR_ID_IN_FORCE_KERNEL
92  for (int k=0; k<4; k++) id_ngb[k] = force.id_ngb[k];
93 #endif
94  n_ngb = force.n_ngb;
95  }
96 
97  PS::F64 getRSearch() const {
98 // return std::max(r_search * std::pow(2.0*mass/m_average,0.3333),r_search_min) * SAFTY_FACTOR_FOR_SEARCH + SAFTY_OFFSET_FOR_SEARCH;
99 // return this->r_search * SAFTY_FACTOR_FOR_SEARCH + SAFTY_OFFSET_FOR_SEARCH;
100  return this->r_search*SAFTY_FACTOR_FOR_SEARCH;
101  }
102 
103  void writeAscii(FILE* fp) const{
104  Ptcl::writeAscii(fp);
105  fprintf(fp, "%26.17e %26.17e %26.17e %26.17e %26.17e ",
106  this->acc.x, this->acc.y, this->acc.z, // 9-11
107  this->pot_tot, this->pot_soft);
108 #ifdef EXTERNAL_POT_IN_PTCL
109  fprintf(fp, "%26.17e ",this->pot_ext);
110 #endif
111  fprintf(fp, "%lld\n",this->n_ngb);
112  }
113 
114  void writeBinary(FILE* fp) const{
115  Ptcl::writeBinary(fp);
116 #ifdef EXTERNAL_POT_IN_PTCL
117  fwrite(&(this->acc), sizeof(PS::F64), 7, fp);
118 #else
119  fwrite(&(this->acc), sizeof(PS::F64), 6, fp);
120 #endif
121  }
122 
123  void readAscii(FILE* fp) {
124  Ptcl::readAscii(fp);
125  PS::S64 rcount=fscanf(fp, "%lf %lf %lf %lf %lf ",
126  &this->acc.x, &this->acc.y, &this->acc.z, // 9-11
127  &this->pot_tot, &this->pot_soft);
128  if (rcount<5) {
129  std::cerr<<"Error: FPSoft Data reading fails! requiring data number is 6, only obtain "<<rcount<<".\n";
130  std::cerr<<"Check your input data, whether the consistent features (interrupt mode and external mode) are used in configuring petar and the data generation\n";
131  abort();
132  }
133 #ifdef EXTERNAL_POT_IN_PTCL
134  rcount=fscanf(fp, "%lf ", &this->pot_ext);
135  if (rcount<1) {
136  std::cerr<<"Error: FPSoft Data reading fails! requiring data number is 1, only obtain "<<rcount<<".\n";
137  std::cerr<<"Check your input data, whether the consistent features (interrupt mode and external mode) are used in configuring petar and the data generation\n";
138  abort();
139  }
140 #endif
141  rcount=fscanf(fp, "%lld\n", &this->n_ngb);
142  if (rcount<1) {
143  std::cerr<<"Error: FPSoft Data reading fails! requiring data number is 1, only obtain "<<rcount<<".\n";
144  std::cerr<<"Check your input data, whether the consistent features (interrupt mode and external mode) are used in configuring petar and the data generation\n";
145  abort();
146  }
147  }
148 
149  void readBinary(FILE* fp) {
150  Ptcl::readBinary(fp);
151 #ifdef EXTERNAL_POT_IN_PTCL
152  size_t rcount = fread(&(this->acc), sizeof(PS::F64), 7, fp);
153  if (rcount<7) {
154  std::cerr<<"Error: Data reading fails! requiring data number is 7, only obtain "<<rcount<<".\n";
155  std::cerr<<"Check your input data, whether the consistent features (interrupt mode and external mode) are used in configuring petar and the data generation\n";
156  abort();
157  }
158 #else
159  size_t rcount = fread(&(this->acc), sizeof(PS::F64), 6, fp);
160  if (rcount<6) {
161  std::cerr<<"Error: Data reading fails! requiring data number is 6, only obtain "<<rcount<<".\n";
162  std::cerr<<"Check your input data, whether the consistent features (interrupt mode and external mode) are used in configuring petar and the data generation\n";
163  abort();
164  }
165 #endif
166  }
167 
168  void print(std::ostream & fout){
169  Ptcl::print(fout);
170  fout<<" acc= "<<acc
171  <<" pot_tot= "<<pot_tot
172  <<" pot_soft= "<<pot_soft
173 #ifdef EXTERNAL_POT_IN_PTCL
174  <<" pot_ext= "<<pot_ext
175 #endif
176  <<" N_b= "<<n_ngb;
177  }
178 
180 
184  static void printColumnTitle(std::ostream & _fout, const int _width=20) {
185  Ptcl::printColumnTitle(_fout, _width);
186  _fout<<std::setw(_width)<<"acc_soft.x"
187  <<std::setw(_width)<<"acc_soft.y"
188  <<std::setw(_width)<<"acc_soft.z"
189  <<std::setw(_width)<<"pot_tot"
190  <<std::setw(_width)<<"pot_soft"
191 #ifdef EXTERNAL_POT_IN_PTCL
192  <<std::setw(_width)<<"pot_ext"
193 #endif
194  <<std::setw(_width)<<"n_b";
195  }
196 
198 
203  static int printTitleWithMeaning(std::ostream & _fout, const int _counter=0, const int _offset=0) {
204  int counter = _counter;
205  counter = Ptcl::printTitleWithMeaning(_fout, counter, _offset);
206  counter++;
207  _fout<<std::setw(_offset)<<" "<<counter<<"-"<<counter+2<<". acc_soft.[x/y/z]: 3D soft (long-range) acceleration (0.0)\n";
208  counter+=3;
209  _fout<<std::setw(_offset)<<" "<<counter<<". pot_tot: total potential (0.0)\n";
210  counter++;
211  _fout<<std::setw(_offset)<<" "<<counter<<". pot_soft: soft potential (0.0)\n";
212  counter++;
213 #ifdef EXTERNAL_POT_IN_PTCL
214  _fout<<std::setw(_offset)<<" "<<counter<<". pot_ext: external potential (0.0)\n";
215  counter++;
216 #endif
217  _fout<<std::setw(_offset)<<" "<<counter<<". n_b: number of neighbors (0)\n";
218  return counter;
219  }
220 
222 
226  void printColumn(std::ostream & _fout, const int _width=20){
227  Ptcl::printColumn(_fout, _width);
228  _fout<<std::setw(_width)<<acc.x
229  <<std::setw(_width)<<acc.y
230  <<std::setw(_width)<<acc.z
231  <<std::setw(_width)<<pot_tot
232  <<std::setw(_width)<<pot_soft
233 #ifdef EXTERNAL_POT_IN_PTCL
234  <<std::setw(_width)<<pot_ext
235 #endif
236  <<std::setw(_width)<<n_ngb;
237  }
238 
240 
245  template <class Tpcm>
246  void printColumnWithOffset(Tpcm& _pcm, std::ostream & _fout, const int _width=20){
247  Ptcl::printColumnWithOffset(_pcm, _fout, _width);
248  _fout<<std::setw(_width)<<acc.x
249  <<std::setw(_width)<<acc.y
250  <<std::setw(_width)<<acc.z
251  <<std::setw(_width)<<pot_tot
252  <<std::setw(_width)<<pot_soft
253 #ifdef EXTERNAL_POT_IN_PTCL
254  <<std::setw(_width)<<pot_ext
255 #endif
256  <<std::setw(_width)<<n_ngb;
257  }
258 
260  void clearForce() {
261  acc = 0.0;
262  pot_tot = 0.0;
263  pot_soft = 0.0;
264 #ifdef EXTERNAL_POT_IN_PTCL
265  pot_ext = 0.0;
266 #endif
267  n_ngb = 0;
268  }
269 };
270 
271 class EPISoft{
272 public:
273  PS::S64 id;
274  PS::F64vec pos;
275  PS::F64 r_search;
276  PS::S32 rank_org;
277  PS::S32 type; // 0: orbital artificial particles; 1: others
278 #ifdef KDKDK_4TH
279  PS::F64vec acc;
280 #endif
281  static PS::F64 eps;
282  static PS::F64 r_out;
283  PS::F64vec getPos() const { return pos;}
284  void copyFromFP(const FPSoft & fp){
285  id = fp.id;
286  pos = fp.pos;
287  if (fp.group_data.artificial.isArtificial() && !fp.group_data.artificial.isCM() && fp.mass>0 ) type = 0;
288  else type = 1;
289 
290 #ifdef KDKDK_4TH
291  acc = fp.acc;
292 #endif
293  r_search = fp.r_search;
294  rank_org = fp.rank_org;
295  }
296  void print(std::ostream & fout=std::cout) const {
297  fout<<" id="<<id
298  <<" rank_org="<<rank_org
299  <<" pos="<<pos
300  <<" eps="<<eps;
301  }
302  PS::F64 getRSearch() const {
303 // return std::max(r_search * std::pow(2.0*mass/m_average,0.3333),r_search_min) * SAFTY_FACTOR_FOR_SEARCH + SAFTY_OFFSET_FOR_SEARCH;
304 // return r_search; * SAFTY_FACTOR_FOR_SEARCH + SAFTY_OFFSET_FOR_SEARCH;
305  return r_search*SAFTY_FACTOR_FOR_SEARCH;
306  }
307 
308 };
309 
310 
311 class EPJSoft{
312 public:
313  PS::S64 id;
314  PS::F64 mass;
315  PS::F64vec pos;
316  PS::F64vec vel;
317 #ifdef KDKDK_4TH
318  PS::F64vec acc;
319 #endif
320  PS::F64 r_in;
321  PS::F64 r_out;
322  PS::F64 r_search;
323  PS::F64 r_scale_next;
324  GroupDataDeliver group_data;
325  PS::S32 rank_org;
326  PS::S32 adr_org;
327 // static PS::F64 r_out;
328 // static PS::F64 m_average;
329 // static PS::F64 r_search_min;
330  void copyFromFP(const FPSoft & fp){
331  id = fp.id;
332  mass = fp.mass;
333  pos = fp.pos;
334  vel = fp.vel;
335 #ifdef KDKDK_4TH
336  acc = fp.acc;
337 #endif
338  r_in = fp.changeover.getRin();
339  r_out = fp.changeover.getRout();
340  r_scale_next = fp.changeover.r_scale_next;
341  r_search = fp.r_search;
342  group_data = fp.group_data;
343  rank_org = fp.rank_org;
344  adr_org = fp.adr;
345  }
346  PS::F64vec getPos() const { return pos; }
347  void setPos(const PS::F64vec & pos_new){ pos = pos_new;}
348  PS::F64 getCharge() const { return mass; }
349  PS::F64 getRSearch() const {
350 // return std::max(r_search * std::pow(2.0*mass/m_average,0.3333),r_search_min) * SAFTY_FACTOR_FOR_SEARCH + SAFTY_OFFSET_FOR_SEARCH;
351 // return r_search * SAFTY_FACTOR_FOR_SEARCH + SAFTY_OFFSET_FOR_SEARCH;
352  return r_search*SAFTY_FACTOR_FOR_SEARCH;
353  }
354 
355  PS::S64 getId() const {
356  return id;
357  }
358 
359  // FORDEBUG
360  void print(std::ostream & fout=std::cout) const {
361  fout<<" id="<<id
362  <<" rank_org="<<rank_org
363  <<" mass="<<mass
364  <<" pos="<<pos
365  <<" vel="<<vel
366  <<" r_search="<<r_search;
367  }
368  void clear(){
369  mass = 0.0;
370  pos = vel = 0.0;
371  r_in = r_out = 0.0;
372  r_search = 0.0;
373  r_scale_next = 1.0;
374  id = rank_org = adr_org = -1;
375  }
376 };
377 
EPJSoft::adr_org
PS::S32 adr_org
Definition: soft_ptcl.hpp:326
EPISoft
Definition: soft_ptcl.hpp:271
FPSoft::printColumn
void printColumn(std::ostream &_fout, const int _width=20)
print data of class members using column style
Definition: soft_ptcl.hpp:226
FPSoft
Definition: soft_ptcl.hpp:26
Ptcl::DataCopy
void DataCopy(const Tptcl &_p)
Definition: ptcl.hpp:58
FPSoft::copyFromForce
void copyFromForce(const ForceSoft &force)
Definition: soft_ptcl.hpp:84
FPSoft::writeAscii
void writeAscii(FILE *fp) const
Definition: soft_ptcl.hpp:103
ForceSoft::clear
void clear()
gravitational constant
Definition: soft_ptcl.hpp:16
EPJSoft::mass
PS::F64 mass
Definition: soft_ptcl.hpp:314
PIKG::S32
int32_t S32
Definition: pikg_vector.hpp:24
FPSoft::readAscii
void readAscii(FILE *fp)
Definition: soft_ptcl.hpp:123
EPJSoft::getPos
PS::F64vec getPos() const
Definition: soft_ptcl.hpp:346
galpy_pot_movie.fp
fp
Definition: galpy_pot_movie.py:131
FPSoft::writeBinary
void writeBinary(FILE *fp) const
Definition: soft_ptcl.hpp:114
ForceSoft
Definition: soft_ptcl.hpp:4
PIKG::F64vec
Vector3< F64 > F64vec
Definition: pikg_vector.hpp:167
EPISoft::rank_org
PS::S32 rank_org
Definition: soft_ptcl.hpp:276
EPJSoft::r_out
PS::F64 r_out
Definition: soft_ptcl.hpp:321
EPISoft::pos
PS::F64vec pos
Definition: soft_ptcl.hpp:274
Ptcl::printTitleWithMeaning
static int printTitleWithMeaning(std::ostream &_fout, const int _counter=0, const int _offset=0)
print column title with meaning (each line for one column)
Definition: ptcl.hpp:99
EPISoft::getPos
PS::F64vec getPos() const
Definition: soft_ptcl.hpp:283
FPSoft::printColumnTitle
static void printColumnTitle(std::ostream &_fout, const int _width=20)
print titles of class members using column style
Definition: soft_ptcl.hpp:184
GroupDataDeliver
group data delivery, used for two purpose
Definition: ptcl.hpp:17
FPSoft::n_ngb
PS::S64 n_ngb
Definition: soft_ptcl.hpp:40
EPISoft::getRSearch
PS::F64 getRSearch() const
Definition: soft_ptcl.hpp:302
FPSoft::readBinary
void readBinary(FILE *fp)
Definition: soft_ptcl.hpp:149
FPSoft::pot_tot
PS::F64 pot_tot
Definition: soft_ptcl.hpp:32
FPSoft::getRSearch
PS::F64 getRSearch() const
Definition: soft_ptcl.hpp:97
PIKG::F64
double F64
Definition: pikg_vector.hpp:17
ForceSoft::pot
PS::F64 pot
soft acceleration (c.m.: averaged force from orbital particles; tensor: c.m. is substracted)
Definition: soft_ptcl.hpp:7
EPJSoft::clear
void clear()
Definition: soft_ptcl.hpp:368
EPJSoft::getId
PS::S64 getId() const
Definition: soft_ptcl.hpp:355
EPISoft::id
PS::S64 id
Definition: soft_ptcl.hpp:273
EPISoft::copyFromFP
void copyFromFP(const FPSoft &fp)
Definition: soft_ptcl.hpp:284
Ptcl::writeBinary
void writeBinary(FILE *_fout) const
write class data with BINARY format
Definition: ptcl.hpp:158
EPISoft::type
PS::S32 type
Definition: soft_ptcl.hpp:277
Ptcl::printColumnTitle
static void printColumnTitle(std::ostream &_fout, const int _width=20)
print titles of class members using column style
Definition: ptcl.hpp:85
ParticleBase::pos
PS::F64vec pos
Definition: particle_base.hpp:24
Ptcl
Particle class.
Definition: ptcl.hpp:36
PIKG::S64
int64_t S64
Definition: pikg_vector.hpp:23
ParticleBase::vel
PS::F64vec vel
Definition: particle_base.hpp:25
SAFTY_FACTOR_FOR_SEARCH
const PS::F64 SAFTY_FACTOR_FOR_SEARCH
Definition: ptcl.hpp:8
ptcl.hpp
EPISoft::r_out
static PS::F64 r_out
Definition: soft_ptcl.hpp:282
FPSoft::acc
PS::F64vec acc
Definition: soft_ptcl.hpp:28
EPJSoft::getCharge
PS::F64 getCharge() const
Definition: soft_ptcl.hpp:348
FPSoft::getVel
PS::F64vec getVel() const
Get velocity (required for ARC::chain)
Definition: soft_ptcl.hpp:57
EPISoft::r_search
PS::F64 r_search
Definition: soft_ptcl.hpp:275
FPSoft::getPos
PS::F64vec getPos() const
Get position (required for ARC::chain)
Definition: soft_ptcl.hpp:50
EPJSoft::r_in
PS::F64 r_in
Definition: soft_ptcl.hpp:320
EPJSoft::getRSearch
PS::F64 getRSearch() const
Definition: soft_ptcl.hpp:349
EPJSoft::rank_org
PS::S32 rank_org
Definition: soft_ptcl.hpp:325
Ptcl::r_search
PS::F64 r_search
Definition: ptcl.hpp:38
Ptcl::writeAscii
void writeAscii(FILE *_fout) const
write class data with ASCII format
Definition: ptcl.hpp:142
EPJSoft
Definition: soft_ptcl.hpp:311
FPSoft::FPSoft
FPSoft(const Tptcl &p, const PS::S32 rank_, const PS::S32 adr_)
Definition: soft_ptcl.hpp:71
EPJSoft::print
void print(std::ostream &fout=std::cout) const
Definition: soft_ptcl.hpp:360
EPJSoft::r_scale_next
PS::F64 r_scale_next
Definition: soft_ptcl.hpp:323
FPSoft::clearForce
void clearForce()
clear force
Definition: soft_ptcl.hpp:260
Ptcl::printColumn
void printColumn(std::ostream &_fout, const int _width=20)
print data of class members using column style
Definition: ptcl.hpp:116
FPSoft::adr
PS::S32 adr
Definition: soft_ptcl.hpp:42
FPSoft::pot_soft
PS::F64 pot_soft
Definition: soft_ptcl.hpp:33
Ptcl::printColumnWithOffset
void printColumnWithOffset(Tpcm &_pcm, std::ostream &_fout, const int _width=20)
print data of class members with pos and vel offset using column style
Definition: ptcl.hpp:131
FPSoft::printTitleWithMeaning
static int printTitleWithMeaning(std::ostream &_fout, const int _counter=0, const int _offset=0)
print column title with meaning (each line for one column)
Definition: soft_ptcl.hpp:203
EPJSoft::pos
PS::F64vec pos
Definition: soft_ptcl.hpp:315
EPJSoft::id
PS::S64 id
Definition: soft_ptcl.hpp:313
EPISoft::eps
static PS::F64 eps
Definition: soft_ptcl.hpp:281
FPSoft::FPSoft
FPSoft()
Definition: soft_ptcl.hpp:45
EPJSoft::vel
PS::F64vec vel
Definition: soft_ptcl.hpp:316
Ptcl::readAscii
void readAscii(FILE *_fin)
read class data with ASCII format
Definition: ptcl.hpp:168
FPSoft::printColumnWithOffset
void printColumnWithOffset(Tpcm &_pcm, std::ostream &_fout, const int _width=20)
print data of class members with pos and vel offset using column style
Definition: soft_ptcl.hpp:246
FPSoft::print
void print(std::ostream &fout)
Definition: soft_ptcl.hpp:168
EPJSoft::group_data
GroupDataDeliver group_data
Definition: soft_ptcl.hpp:324
FPSoft::rank_org
PS::S32 rank_org
Definition: soft_ptcl.hpp:41
ForceSoft::acc
PS::F64vec acc
Definition: soft_ptcl.hpp:6
EPJSoft::setPos
void setPos(const PS::F64vec &pos_new)
Definition: soft_ptcl.hpp:347
ForceSoft::n_ngb
PS::S64 n_ngb
full potential
Definition: soft_ptcl.hpp:14
Ptcl::readBinary
void readBinary(FILE *_fin)
read class data with BINARY format
Definition: ptcl.hpp:194
ForceSoft::grav_const
static PS::F64 grav_const
neighbor number+1
Definition: soft_ptcl.hpp:15
EPJSoft::r_search
PS::F64 r_search
Definition: soft_ptcl.hpp:322
Ptcl::print
void print(std::ostream &_fout) const
Definition: ptcl.hpp:72
EPJSoft::copyFromFP
void copyFromFP(const FPSoft &fp)
Definition: soft_ptcl.hpp:330
EPISoft::print
void print(std::ostream &fout=std::cout) const
Definition: soft_ptcl.hpp:296