00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025 #ifndef _ORSA_INTEGRATOR_H_
00026 #define _ORSA_INTEGRATOR_H_
00027
00028 #include "orsa_frame.h"
00029 #include "orsa_error.h"
00030 #include "orsa_interaction.h"
00031
00032 #include <string>
00033
00034 namespace orsa {
00035
00036 enum IntegratorType {
00037 STOER=1,
00038 BULIRSCHSTOER=2,
00039 RUNGEKUTTA=3,
00040 DISSIPATIVERUNGEKUTTA=4,
00041 RA15=5,
00042 LEAPFROG=6
00043 };
00044
00045 inline void convert(IntegratorType &it, const unsigned int i) {
00046 switch(i) {
00047 case 1: it = STOER; break;
00048 case 2: it = BULIRSCHSTOER; break;
00049 case 3: it = RUNGEKUTTA; break;
00050 case 4: it = DISSIPATIVERUNGEKUTTA; break;
00051 case 5: it = RA15; break;
00052 case 6: it = LEAPFROG; break;
00053
00054 default:
00055 ORSA_ERROR("conversion problem: i = %i",i);
00056 break;
00057 }
00058 }
00059
00060 inline std::string label(const IntegratorType it) {
00061 std::string s;
00062 switch (it) {
00063 case STOER: s="Stoer"; break;
00064 case BULIRSCHSTOER: s="Bulirsch-Stoer"; break;
00065 case RUNGEKUTTA: s="Runge-Kutta 4th order"; break;
00066 case DISSIPATIVERUNGEKUTTA: s="Dissipative Runge-Kutta 4th order"; break;
00067 case RA15: s="Everhart's RADAU 15th order"; break;
00068 case LEAPFROG: s="Leapfrog 2nd order"; break;
00069 }
00070 return s;
00071 }
00072
00073
00074 class Integrator {
00075 public:
00076 virtual void Step(const Frame &, Frame &, Interaction *) = 0;
00077 virtual ~Integrator() { };
00078
00079 public:
00080 virtual Integrator * clone() const = 0;
00081
00082 public:
00083
00084 UniverseTypeAwareTimeStep timestep;
00085
00086 protected:
00087 UniverseTypeAwareTimeStep timestep_done;
00088
00089 public:
00090
00091 double accuracy;
00092 unsigned int m;
00093
00094 public:
00095 virtual bool can_handle_velocity_dependant_interactions() const { return false; }
00096
00097 public:
00098 inline IntegratorType GetType() const { return type; }
00099
00100 protected:
00101 IntegratorType type;
00102 };
00103
00104 void make_new_integrator(Integrator**, const IntegratorType);
00105
00106
00107
00108
00109
00110 class FixedTimestepIntegrator : public Integrator {
00111
00112 };
00113
00114 class MultistepIntegrator : public FixedTimestepIntegrator {
00115
00116 };
00117
00118 class VariableTimestepIntegrator : public Integrator {
00119
00120 };
00121
00122
00123
00124
00125
00126
00127 class ModifiedMidpoint : public MultistepIntegrator {
00128
00129 };
00130
00131
00132 class Stoer : public MultistepIntegrator {
00133 public:
00134 Stoer();
00135 Stoer(int);
00136 Stoer(const Stoer &);
00137 ~Stoer();
00138
00139 public:
00140 void Step(const Frame&, Frame&, Interaction*);
00141
00142 public:
00143 Integrator * clone() const;
00144 };
00145
00146
00147
00148
00149
00150
00151
00152
00153
00154
00155
00156
00157
00158
00159
00160
00161
00162
00163
00164
00165
00166
00167
00168
00169
00170
00171
00172
00173
00174
00175
00176
00177
00178
00179 class RungeKutta : public FixedTimestepIntegrator {
00180 public:
00181 RungeKutta();
00182 RungeKutta(const RungeKutta &);
00183 ~RungeKutta();
00184
00185 public:
00186 void Step(const Frame&, Frame&, Interaction*);
00187
00188 public:
00189 Integrator * clone() const;
00190 };
00191
00192 class DissipativeRungeKutta : public FixedTimestepIntegrator {
00193 public:
00194 DissipativeRungeKutta();
00195 DissipativeRungeKutta(const DissipativeRungeKutta &);
00196 ~DissipativeRungeKutta();
00197
00198 public:
00199 bool can_handle_velocity_dependant_interactions() const { return true; }
00200
00201 public:
00202 void Step(const Frame&, Frame&, Interaction*);
00203
00204 public:
00205 Integrator * clone() const;
00206 };
00207
00208
00209
00210
00211
00212
00213
00214 class Radau15 : public VariableTimestepIntegrator {
00215 public:
00216 Radau15();
00217 Radau15(const Radau15 &);
00218 ~Radau15();
00219
00220 private:
00221 void init();
00222 void Bodies_Mass_or_N_Bodies_Changed(const Frame&);
00223
00224 public:
00225 bool can_handle_velocity_dependant_interactions() const { return true; }
00226
00227 public:
00228 void Step(const Frame&, Frame&, Interaction*);
00229
00230 public:
00231 Integrator * clone() const;
00232
00233 private:
00234
00235
00236
00237
00238 double h[8];
00239
00240
00241
00242
00243
00244 double xc[8],vc[7];
00245
00246
00247 double r[28],c[21],d[21],s[9];
00248
00249 std::vector< std::vector<double> > g,b,e;
00250
00251
00252 unsigned int nv,niter;
00253
00254 std::vector<double> x,v,a,x1,v1,a1;
00255
00256
00257 std::vector<double> mass;
00258
00259
00260 std::vector<Vector> acc;
00261
00262
00263 unsigned int size;
00264 };
00265
00266 class Leapfrog : public FixedTimestepIntegrator {
00267 public:
00268 Leapfrog();
00269 Leapfrog(const Leapfrog &);
00270 ~Leapfrog();
00271
00272 public:
00273 void Step(const Frame&, Frame&, Interaction*);
00274
00275 public:
00276 Integrator * clone() const;
00277 };
00278
00279 }
00280
00281 #endif // _ORSA_INTEGRATOR_H_
00282