Swarm-NG  1.1
stop_on_ejection.hpp
Go to the documentation of this file.
1 /*************************************************************************
2  * Copyright (C) 2011 by Eric Ford and the Swarm-NG Development Team *
3  * *
4  * This program is free software; you can redistribute it and/or modify *
5  * it under the terms of the GNU General Public License as published by *
6  * the Free Software Foundation; either version 3 of the License. *
7  * *
8  * This program is distributed in the hope that it will be useful, *
9  * but WITHOUT ANY WARRANTY; without even the implied warranty of *
10  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
11  * GNU General Public License for more details. *
12  * *
13  * You should have received a copy of the GNU General Public License *
14  * along with this program; if not, write to the *
15  * Free Software Foundation, Inc., *
16  * 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. *
17  ************************************************************************/
18 
25 #pragma once
26 
27 #include <limits>
28 #include "../swarm/types/config.hpp"
29 
30 namespace swarm {
31  namespace monitors {
32 
42 struct stop_on_ejection_params {
43  double rmax;
44  bool deactivate_on, log_on, verbose_on;
47  stop_on_ejection_params(const config &cfg)
48  {
49  rmax = cfg.optional("rmax",std::numeric_limits<float>::max());
50  deactivate_on = cfg.optional("deactivate_on_ejection",false);
51  log_on = cfg.optional("log_on_ejection",false);
52  verbose_on = cfg.optional("verbose_on_ejection",false);
53  }
54 };
55 
64 template<class log_t>
65 class stop_on_ejection {
66  public:
67  typedef stop_on_ejection_params params;
68 
69  private:
70  params _params;
71  bool condition_met;
72 
73  ensemble::SystemRef& _sys;
74  log_t& _log;
75 
76  public:
77 
78  GPUAPI bool is_deactivate_on() { return _params.deactivate_on; };
79  GPUAPI bool is_log_on() { return _params.log_on; };
80  GPUAPI bool is_verbose_on() { return _params.verbose_on; };
81  GPUAPI bool is_any_on() { return is_deactivate_on() || is_log_on() || is_verbose_on() ; }
82  GPUAPI bool is_condition_met () { return ( condition_met ); }
83  GPUAPI bool need_to_log_system ()
84  { return (is_log_on() && is_condition_met() ); }
85  GPUAPI bool need_to_deactivate ()
86  { return ( is_deactivate_on() && is_condition_met() ); }
87 
88  GPUAPI void log_system() { log::system(_log, _sys); }
89 
90  template<class T>
91  static GENERIC int thread_per_system(T compile_time_param){
92  return 1;
93  }
94 
95  template<class T>
96  static GENERIC int shmem_per_system(T compile_time_param) {
97  return 0;
98  }
99 
100  public:
101 
108  GPUAPI bool test_body(const int& b) {
109 
110  double x,y,z,vx,vy,vz; _sys[b].get(x,y,z,vx,vy,vz);
111  // double r = sqrt(_sys[b].distance_to_origin_squared()); // WARNING: Deceiving function name
112  double r = _sys[b].distance_to_origin(); // WARNING: Deceiving function name
113  if( r < _params.rmax ) return false;
114  double rdotv = x*vx+y*vy+z*vz;
115  if( rdotv <= 0. ) return false;
116 
117  bool stopit = false;
118  double speed_sq = _sys[b].speed_squared(); // WARNING: Deceiving function name
119  double epp = 0.5*speed_sq*r/_sys[b].mass()-1.;
120  if( fabs(epp) < 1e-4 ) {
121  double energy = 0.5*speed_sq-_sys[b].mass()/r;
122  if(is_verbose_on() )
123  lprintf(_log, "Orbit is nearly parabolic: _sys=%d, bod=%d, T=%lg r=%lg energy=%lg energy*r/GM=%lg.\n"
124  , _sys.number(), b, _sys.time() , r, energy, epp );
125  stopit = true;
126  }else if ( epp > 0 ){
127  double energy = 0.5*speed_sq-_sys[b].mass()/r;
128  if(is_verbose_on() )
129  lprintf(_log, "Orbit is hyperbolic: _sys=%d, bod=%d, T=%lg r=%lg energy=%lg energy*r/GM=%lg.\n"
130  , _sys.number(), b, _sys.time() , r, energy, epp );
131  // TODO: Make sure that planet is not near another body
132  // This is very unlikely to be an issue, provided that rmax
133  // is set to be well beyond the initial semi-major axes
134  stopit = true;
135  }
136 
137  return stopit;
138  }
139 
140  GPUAPI void operator () (const int thread_in_system)
141  {
142  pass_one(thread_in_system);
143  pass_two(thread_in_system);
144  if(need_to_log_system() && (thread_in_system==0) )
145  log_system();
146  }
147 
148 #if 0
149  // GPUAPI void operator () () {
150  GPUAPI void operator () (int thread_in_system)
151  {
152  if(!is_any_on()) return;
153  bool stopit = false;
154 
155  if(thread_in_system==0)
156  {
157  // Check each body other than central star
158  for(int b = 1; b < _sys.nbod(); b++)
159  stopit = stopit || test_body(b);
160 
161  if(stopit)
162  {
163  if(is_log_on())
164  log::system(_log, _sys);
165  if(is_deactivate_on())
166  _sys.set_disabled();
167  }
168  }
169  }
170 #endif
171 
173  GPUAPI stop_on_ejection(const params& p,ensemble::SystemRef& s,log_t& l)
174  :_params(p),_sys(s),_log(l){}
175 
177  GPUAPI bool pass_one (int thread_in_system)
178  {
179  bool need_full_test = false;
180  condition_met = false;
181  if(is_any_on()&&(thread_in_system==0))
182  {
183  // Check each body other than central star
184  for(int b = 1; b < _sys.nbod(); b++)
185  condition_met = condition_met || test_body(b);
186 
187  if( condition_met && is_log_on() )
188  { need_full_test = true; }
189 
190 
191  }
192  return need_full_test;
193  }
194 
196  GPUAPI int pass_two (int thread_in_system)
197  {
198  if(is_condition_met() && is_deactivate_on() )
199  { _sys.set_disabled(); }
200  return _sys.state();
201  }
202 
203 
204 };
205 
206 }
207 }