Newton Class Reference

#include <orsa_interaction.h>

Inheritance diagram for Newton:

Inheritance graph
[legend]
Collaboration diagram for Newton:

Collaboration graph
[legend]
List of all members.

Public Member Functions

 Newton ()
 Newton (const Newton &)
void Acceleration (const Frame &, std::vector< Vector > &)
double PotentialEnergy (const Frame &)
bool depends_on_velocity () const
Interactionclone () const
InteractionType GetType () const
void IncludeMultipoleMoments (const bool b)
bool IsIncludingMultipoleMoments () const
void IncludeRelativisticEffects (const bool b)
bool IsIncludingRelativisticEffects () const
void IncludeFastRelativisticEffects (const bool b)
bool IsIncludingFastRelativisticEffects () const
void SkipJPLPlanets (const bool b)
bool IsSkippingJPLPlanets () const

Protected Attributes

bool skip_JPL_planets

Constructor & Destructor Documentation

Newton  ) 
 

Definition at line 153 of file orsa_interaction.cc.

References orsa::GetC().

Referenced by Newton::clone().

00153                  : Interaction(), include_multipole_moments(false), include_relativistic_effects(false), include_fast_relativistic_effects(false), one_over_c2(1.0/(GetC()*GetC())) {
00154     skip_JPL_planets = false;
00155   }

Here is the call graph for this function:

Newton const Newton  ) 
 

Definition at line 157 of file orsa_interaction.cc.

References orsa::GetC().

00157                                  : Interaction(), include_multipole_moments(n.include_multipole_moments), include_relativistic_effects(n.include_relativistic_effects), include_fast_relativistic_effects(n.include_fast_relativistic_effects), one_over_c2(1.0/(GetC()*GetC())) {
00158     skip_JPL_planets = n.skip_JPL_planets;
00159   }

Here is the call graph for this function:


Member Function Documentation

void Acceleration const Frame ,
std::vector< Vector > & 
[virtual]
 

Implements Interaction.

Definition at line 199 of file orsa_interaction.cc.

References orsa::alpha_delta_meridian(), MappedTable::Distance(), MappedTable::DistanceVector(), MappedTable::DistanceVectorOverDistanceCube(), Legendre::dP2, Legendre::dP22, Legendre::dP3, Legendre::dP31, Legendre::dP32, Legendre::dP33, Legendre::dP4, Legendre::dP41, Legendre::dP42, Legendre::dP43, Legendre::dP44, orsa::ECLIPTIC, orsa::ExternalProduct(), Angle::GetRad(), Universe::GetReferenceSystem(), orsa::halfpi, Vector::LengthSquared(), MappedTable::load(), orsa::NONE, Vector::Normalized(), orsa::obleq_J2000(), MappedTable::OneOverDistance(), MappedTable::OneOverDistanceCube(), MappedTable::OneOverDistanceSquare(), ORSA_WARNING, Legendre::P2, Legendre::P22, Legendre::P3, Legendre::P31, Legendre::P32, Legendre::P33, Legendre::P4, Legendre::P41, Legendre::P42, Legendre::P43, Legendre::P44, Vector::rotate(), Vector::Set(), orsa::sincos(), and orsa::universe.

00199                                                                {
00200     
00201     const unsigned int size = f.size();
00202     
00203     if (size < 2) return;
00204     
00205     a.resize(size);
00206     /* 
00207        for (unsigned int i=0;i<size;++i)
00208        a[i].Set(0.0,0.0,0.0);
00209     */
00210     
00211     zero_mass.resize(size);
00212     for (unsigned int i=0;i<size;++i) {
00213       zero_mass[i] = f[i].has_zero_mass();
00214     }
00215     
00216     skip.resize(size);
00217     if (skip_JPL_planets) {
00218       for (unsigned int k=0;k<size;++k) {
00219         skip[k] = (f[k].JPLPlanet() != NONE);
00220       }
00221     } else {
00222       for (unsigned int k=0;k<size;++k) {
00223         skip[k] = false;
00224       }
00225     }
00226     
00227     const bool only_newton = !(include_multipole_moments || include_relativistic_effects || include_fast_relativistic_effects);
00228     
00229     if (only_newton) {
00230       for (unsigned int i=0;i<size;++i) {
00231         a[i].Set(0.0,0.0,0.0);
00232       }
00233       //
00234       fast_newton_acc(f,a);
00235       //
00236       return;
00237     }
00238     
00239     if (include_relativistic_effects && include_fast_relativistic_effects) {
00240       ORSA_WARNING("Both the accurate and the fast version of the Relativistic corrections are activated!");
00241     }
00242     
00243     if (a_newton.size() != size) {
00244       a_newton.resize(size);
00245       a_multipoles.resize(size);
00246       a_relativity.resize(size);
00247     }
00248     //
00249     for (unsigned int i=0;i<size;++i) {
00250       a_newton[i].Set(0.0,0.0,0.0);
00251       a_multipoles[i].Set(0.0,0.0,0.0);
00252       a_relativity[i].Set(0.0,0.0,0.0);
00253     }
00254     
00255     // mod_f initialization needed!
00256     Frame mod_f = f;
00257     //
00258     if (include_relativistic_effects) {
00259       // a_newton.resize(size);
00260       //
00261       // The RelativisticBarycenter is _very_ close to the normal Barycenter,
00262       // and the accuracy of a "double" is often barely enough to distinguish them.
00263       // I'll use the normal Barycenter...
00264       //
00265       const Vector barycenter          = f.Barycenter();
00266       const Vector barycenter_velocity = f.BarycenterVelocity();
00267       // 
00268       // mod_f already equalt to f
00269       // mod_f = f;
00270       //
00271       for (unsigned int k=0; k<f.size(); ++k) {
00272         mod_f[k].AddToPosition(-barycenter);
00273         mod_f[k].AddToVelocity(-barycenter_velocity);
00274       }
00275       //
00276       mapped_table.load(mod_f,skip_JPL_planets);
00277     } else {
00278       mapped_table.load(f,skip_JPL_planets);
00279     }
00280     
00281     // note: mapped(i,j) = something(i) - something(j)
00282     
00283     /* 
00284        {
00285        // test        
00286        for (unsigned int i=0;i<size-1;++i) {
00287        for (unsigned int j=i+1;j<size;++j) {
00288        if (zero_mass[i] && zero_mass[j]) continue;
00289        if (skip[i] && skip[j]) continue;
00290        if (fabs((mapped_table.Distance(i,j)-f[i].DistanceVector(f[j]).Length())/(mapped_table.Distance(i,j))) > 1.0e-10) {
00291        fprintf(stderr,
00292        "PROBLEM-1!! d(%10s,%10s): mapped=%20.12f num=%20.12f\n",
00293        f[i].name().c_str(),
00294        f[j].name().c_str(),
00295        mapped_table.Distance(i,j),
00296        f[i].DistanceVector(f[j]).Length());
00297        }
00298        if (fabs((mapped_table.DistanceVector(i,j).x+f[i].DistanceVector(f[j]).x)/mapped_table.DistanceVector(i,j).x) > 1.0e-10) {
00299        fprintf(stderr,
00300        "PROBLEM-2!! d(%10s,%10s): mapped=%20.12f num=%20.12f\n",
00301        f[i].name().c_str(),
00302        f[j].name().c_str(),
00303        mapped_table.DistanceVector(i,j).x,
00304        f[i].DistanceVector(f[j]).x);
00305        }
00306        }
00307        }
00308        }
00309     */
00310     
00311     if (include_multipole_moments) {
00312       axis.resize(size);
00313       x_axis.resize(size);
00314       R1.resize(size);
00315       R2.resize(size);
00316       R3.resize(size);
00317       R4.resize(size);
00318       //
00319       {
00320         Angle alpha,delta,meridian;
00321         Vector planet_axis;
00322         Vector planet_x_axis;
00323         for (unsigned int i=0;i<size;++i) {
00324           if (f[i].JPLPlanet() == NONE) continue;
00325           alpha_delta_meridian(f[i].JPLPlanet(),f,alpha,delta,meridian);
00326           //
00327           planet_axis.Set(0,0,1);
00328           planet_x_axis.Set(1,0,0);
00329           // glRotated(90.0+(180.0/pi)*W.GetRad(),0,0,1);
00330           // no halfpi needed!
00331           planet_axis.rotate(meridian.GetRad(),0,0);
00332           planet_x_axis.rotate(meridian.GetRad(),0,0);
00333           // glRotated(90.0-(180.0/pi)*delta.GetRad(),1,0,0);
00334           planet_axis.rotate(0,halfpi-delta.GetRad(),0);
00335           planet_x_axis.rotate(0,halfpi-delta.GetRad(),0);
00336           // glRotated(90.0+(180.0/pi)*alpha.GetRad(),0,0,1);
00337           planet_axis.rotate(halfpi+alpha.GetRad(),0,0);
00338           planet_x_axis.rotate(halfpi+alpha.GetRad(),0,0);
00339           //
00340           if (universe->GetReferenceSystem() == ECLIPTIC) {
00341             // glRotated(-(180.0/pi)*obleq_J2000().GetRad(),1,0,0);
00342             planet_axis.rotate(0,-obleq_J2000().GetRad(),0);
00343             planet_x_axis.rotate(0,-obleq_J2000().GetRad(),0);
00344           }
00345           //
00346           axis[i] = planet_axis;
00347           x_axis[i] = planet_x_axis;
00348         }
00349       }
00350       //
00351       {
00352         double R;
00353         for (unsigned int i=0;i<size;++i) {
00354           R = f[i].radius();
00355           if (R1[i] != R) {
00356             R1[i] = R;
00357             R2[i] = R*R;
00358             R3[i] = R2[i]*R;
00359             R4[i] = R3[i]*R;
00360           }
00361         }
00362       }
00363     }
00364     
00365     for (unsigned int i=0;i<size-1;++i) {
00366       for (unsigned int j=i+1;j<size;++j) {
00367         
00368         if (zero_mass[i] && zero_mass[j]) continue;
00369         
00370         if (skip[i] && skip[j]) continue;
00371         
00372         if (mapped_table.Distance(i,j) < (std::numeric_limits<double>::min() * 1.0e3)) {
00373           ORSA_WARNING("two objects in the same position! (%s and %s)",f[i].name().c_str(),f[j].name().c_str());
00374           continue;
00375         }
00376         
00377         a_newton[i] += mapped_table.DistanceVectorOverDistanceCube(j,i) * f[j].mu();
00378         a_newton[j] += mapped_table.DistanceVectorOverDistanceCube(i,j) * f[i].mu();
00379       }
00380     }
00381     
00382     if (include_relativistic_effects) {
00383       
00384       for (unsigned int i=0;i<size-1;++i) {     
00385         for (unsigned int j=i+1;j<size;++j) {
00386           
00387           /* 
00388              for (unsigned int i=0;i<size;++i) {        
00389              for (unsigned int j=0;j<size;++j) {
00390           */
00391           
00392           // if (i == j) continue;
00393           
00394           if (zero_mass[i] && zero_mass[j]) continue;
00395           
00396           if (skip[i] && skip[j]) continue;
00397           
00398           // Now mod_f should be the Frame in the relativistic barycenter reference frame
00399           
00400           // Relativistic interaction, from:
00401           // "Explanatory Supplement to the Astronomical Almanac",
00402           // Edited by P. Kenneth Seidelmann, U.S. Naval Observatory,
00403           // 2nd edition, 1992, sect. 5.211, pag. 281.
00404           
00405           double sum_ik = 0.0;
00406           double sum_jk = 0.0;
00407           for (unsigned int k=0; k<size; ++k) {
00408             if (i != k) sum_ik += f[k].mu()*mapped_table.OneOverDistance(i,k);
00409             if (j != k) sum_jk += f[k].mu()*mapped_table.OneOverDistance(j,k);
00410           }
00411           
00412           a_relativity[i] += one_over_c2*f[j].mu()*mapped_table.DistanceVectorOverDistanceCube(j,i)*( // missing Newton term...
00413                                                                                                      - 4.0*sum_ik
00414                                                                                                      - 1.0*sum_jk
00415                                                                                                      + mod_f[i].velocity().LengthSquared()
00416                                                                                                      + 2.0*mod_f[j].velocity().LengthSquared()
00417                                                                                                      - 4.0*mod_f[i].velocity()*mod_f[j].velocity()
00418                                                                                                      - 1.5*pow(mapped_table.DistanceVector(i,j)*mod_f[j].velocity()*mapped_table.OneOverDistance(i,j),2)
00419                                                                                                      + 0.5*mapped_table.DistanceVector(j,i)*a_newton[j]
00420                                                                                                      )
00421             + one_over_c2*f[j].mu()*mapped_table.OneOverDistanceCube(i,j)*((mapped_table.DistanceVector(i,j)*(4.0*mod_f[i].velocity()-3.0*mod_f[j].velocity()))*(mod_f[i].velocity()-mod_f[j].velocity()))
00422             + 3.5*one_over_c2*f[j].mu()*mapped_table.OneOverDistance(i,j)*a_newton[j];
00423           
00424           a_relativity[j] += one_over_c2*f[i].mu()*mapped_table.DistanceVectorOverDistanceCube(i,j)*( // missing Newton term...
00425                                                                                                      - 4.0*sum_jk
00426                                                                                                      - 1.0*sum_ik
00427                                                                                                      + mod_f[j].velocity().LengthSquared()
00428                                                                                                      + 2.0*mod_f[i].velocity().LengthSquared()
00429                                                                                                      - 4.0*mod_f[j].velocity()*mod_f[i].velocity()
00430                                                                                                      - 1.5*pow(mapped_table.DistanceVector(j,i)*mod_f[i].velocity()*mapped_table.OneOverDistance(j,i),2)
00431                                                                                                      + 0.5*mapped_table.DistanceVector(i,j)*a_newton[i]
00432                                                                                                      )
00433             + one_over_c2*f[i].mu()*mapped_table.OneOverDistanceCube(j,i)*((mapped_table.DistanceVector(j,i)*(4.0*mod_f[j].velocity()-3.0*mod_f[i].velocity()))*(mod_f[j].velocity()-mod_f[i].velocity()))
00434             + 3.5*one_over_c2*f[i].mu()*mapped_table.OneOverDistance(j,i)*a_newton[i];
00435         }
00436       }
00437     }
00438     
00439     if (include_fast_relativistic_effects) {
00440       for (unsigned int i=0;i<size-1;++i) {     
00441         for (unsigned int j=i+1;j<size;++j) {
00442           
00443           // if (i == j) continue;
00444           
00445           if (zero_mass[i] && zero_mass[j]) continue;
00446           
00447           if (skip[i] && skip[j]) continue;
00448           
00449           // tests...
00450           // good old code...
00451           
00452           // const Vector r = mod_f[i].position() - mod_f[j].position();
00453           const Vector r = mapped_table.DistanceVector(i,j);
00454           //
00455           const Vector v  = mod_f[i].velocity() - mod_f[j].velocity();
00456           const double v2 = v.LengthSquared();
00457           
00458           // see http://www.projectpluto.com/relativi.htm
00459           //
00460           // acc = acc_{newton} - \frac{G M}{r^3 c^2} ( (\frac{4 G M}{r} - v^2) \vec(r) + 4 (v \cdot r) \vec(v) )
00461           //
00462           
00463           a_relativity[i] += one_over_c2*f[j].mu()*mapped_table.OneOverDistanceCube(i,j)*((4.0*f[j].mu()*mapped_table.OneOverDistance(i,j)-v2)*r + 
00464                                                                                           4.0*(r*v)*v);
00465           a_relativity[j] -= one_over_c2*f[i].mu()*mapped_table.OneOverDistanceCube(i,j)*((4.0*f[i].mu()*mapped_table.OneOverDistance(i,j)-v2)*r +
00466                                                                                           4.0*(r*v)*v);
00467         }
00468       }
00469     }
00470     
00471     if (include_multipole_moments) {
00472       
00473       /* for (unsigned int i=0;i<size-1;++i) {
00474          for (unsigned int j=i+1;j<size;++j) {
00475       */
00476       // scan all i's and j's, becaues i is always the extended body,
00477       // while j is always the point mass body
00478       for (unsigned int i=0;i<size;++i) {
00479         for (unsigned int j=0;j<size;++j) {
00480           
00481           if (i == j) continue;
00482           
00483           if (zero_mass[i] && zero_mass[j]) continue;
00484           
00485           if (skip[i] && skip[j]) continue;
00486           
00487           /* 
00488              const Vector      d = mapped_table.DistanceVector(i,j);
00489              const Vector unit_d = d.Normalized();
00490              
00491              const double l2 = mapped_table.Distance2(i,j);
00492              const double l3 = mapped_table.Distance3(i,j);
00493              const double l4 = mapped_table.Distance4(i,j);
00494           */
00495           
00496           if ((R1[i] > 0.0) && (fabs(f[i].J2()) > 0.0) && (f[i].JPLPlanet() != NONE)) {
00497             
00498             // cerr << "figure: extended body: " << f[i].name() << "  point-mass body: " << f[j].name() << endl;
00499             
00500             const Vector unit_x_local = mapped_table.DistanceVector(j,i).Normalized();
00501             const Vector unit_y_local = ExternalProduct(axis[i],unit_x_local);
00502             const Vector unit_z_local = ExternalProduct(unit_x_local,unit_y_local);
00503             
00504             const Vector tmp_y_axis = ExternalProduct(axis[i],x_axis[i]);
00505             const double tmp_uy_y   = unit_y_local*tmp_y_axis;
00506             const double tmp_uy_mx  = unit_y_local*x_axis[i];
00507             const double lambda = atan2(tmp_uy_mx,tmp_uy_y);
00508             //
00509             double s1l,c1l; orsa::sincos(lambda,s1l,c1l);
00510             double s2l,c2l; orsa::sincos(lambda,s2l,c2l);
00511             double s3l,c3l; orsa::sincos(lambda,s3l,c3l);
00512             double s4l,c4l; orsa::sincos(lambda,s4l,c4l);
00513             
00514             const double theta = acos(unit_x_local*axis[i]);
00515             //
00516             double st,ct; orsa::sincos(theta,st,ct);
00517             const double cosec_theta = 1.0/st;      
00518             
00519             const double ar  = R1[i]*mapped_table.OneOverDistance(i,j);
00520             const double ar2 = ar*ar;
00521             const double ar3 = ar*ar2;
00522             const double ar4 = ar*ar3;
00523             
00524             const Legendre l(axis[i]*unit_x_local);
00525             
00526             // cerr << "x: " << axis[i]*unit_x_local << endl;
00527             
00528             // const double base_coefficient = - f[j].mu() * mapped_table.OneOverDistanceSquare(i,j);
00529             const double base_coefficient = - mapped_table.OneOverDistanceSquare(i,j);
00530             
00531             const double acc_x_local = 
00532               base_coefficient *
00533               ( // zonal harmonics
00534                f[i].J2()*ar2*3.0*l.P2 +
00535                f[i].J3()*ar3*4.0*l.P3 + 
00536                f[i].J4()*ar4*5.0*l.P4 +
00537                // tesseral harmonics
00538                ar2*(-3.0*l.P22*(f[i].C22()*c2l) ) +
00539                ar3*(-4.0*l.P31*(f[i].C31()*c1l+f[i].S31()*s1l) +
00540                     -4.0*l.P32*(f[i].C32()*c2l+f[i].S32()*s2l) +
00541                     -4.0*l.P33*(f[i].C33()*c3l+f[i].S33()*s3l) ) +
00542                ar4*(-5.0*l.P41*(f[i].C41()*c1l+f[i].S41()*s1l) +
00543                     -5.0*l.P42*(f[i].C42()*c2l+f[i].S42()*s2l) +
00544                     -5.0*l.P43*(f[i].C43()*c3l+f[i].S43()*s3l) + 
00545                     -5.0*l.P44*(f[i].C44()*c4l+f[i].S44()*s4l) ) );
00546             
00547             const double acc_y_local = 
00548               base_coefficient *
00549               ( // zonal harmonics
00550                // no zonal harmonics along the y direction
00551                // tesseral harmonics
00552                ar2*(2.0*cosec_theta*l.P22*(-f[i].C22()*s2l) ) +
00553                ar3*(1.0*cosec_theta*l.P31*(-f[i].C31()*s1l+f[i].S31()*c1l) +
00554                     2.0*cosec_theta*l.P32*(-f[i].C32()*s2l+f[i].S32()*c2l) +
00555                     3.0*cosec_theta*l.P33*(-f[i].C33()*s3l+f[i].S33()*c3l) ) +
00556                ar4*(1.0*cosec_theta*l.P41*(-f[i].C41()*s1l+f[i].S41()*c1l) +
00557                     2.0*cosec_theta*l.P42*(-f[i].C42()*s2l+f[i].S42()*c2l) +
00558                     3.0*cosec_theta*l.P43*(-f[i].C43()*s3l+f[i].S43()*c3l) + 
00559                     4.0*cosec_theta*l.P44*(-f[i].C44()*s4l+f[i].S44()*c4l) ) );
00560             
00561             const double acc_z_local = 
00562               base_coefficient *
00563               ( // zonal harmonics
00564                f[i].J2()*ar2*(-st)*l.dP2 +
00565                f[i].J3()*ar3*(-st)*l.dP3 + 
00566                f[i].J4()*ar4*(-st)*l.dP4 +
00567                // tesseral harmonics
00568                ar2*(st*l.dP22*(f[i].C22()*c2l) ) +
00569                ar3*(st*l.dP31*(f[i].C31()*c1l+f[i].S31()*s1l) +
00570                     st*l.dP32*(f[i].C32()*c2l+f[i].S32()*s2l) +
00571                     st*l.dP33*(f[i].C33()*c3l+f[i].S33()*s3l) ) +
00572                ar4*(st*l.dP41*(f[i].C41()*c1l+f[i].S41()*s1l) +
00573                     st*l.dP42*(f[i].C42()*c2l+f[i].S42()*s2l) +
00574                     st*l.dP43*(f[i].C43()*c3l+f[i].S43()*s3l) + 
00575                     st*l.dP44*(f[i].C44()*c4l+f[i].S44()*s4l) ) );
00576             
00577             const Vector local_acc =
00578               unit_x_local * acc_x_local +
00579               unit_y_local * acc_y_local +
00580               unit_z_local * acc_z_local;
00581             
00582             // cerr << "F: i: " << f[i].name() << " j: " << f[j].name() << " ax: " << acc_x_local << " tot: " << local_acc.Length() << endl;
00583             
00584             a_multipoles[i] += f[j].mu() * local_acc;
00585             
00586             // reaction
00587             a_multipoles[j] -= f[i].mu() * local_acc;
00588             
00589           }
00590         }
00591       }
00592     }
00593     
00594     for (unsigned int i=0;i<size;++i) {
00595       a[i] = a_newton[i] + (a_multipoles[i] + a_relativity[i]);
00596     }
00597     
00598     // done
00599   }

Here is the call graph for this function:

Interaction * clone  )  const [virtual]
 

Implements Interaction.

Definition at line 161 of file orsa_interaction.cc.

References Newton::Newton().

00161                                     {
00162     return new Newton(*this);
00163   }

Here is the call graph for this function:

bool depends_on_velocity  )  const [inline, virtual]
 

Reimplemented from Interaction.

Definition at line 250 of file orsa_interaction.h.

00250                                             { 
00251       return (include_relativistic_effects || include_fast_relativistic_effects); 
00252     }

InteractionType GetType  )  const [inline, virtual]
 

Implements Interaction.

Definition at line 257 of file orsa_interaction.h.

References orsa::NEWTON.

00257                                     {
00258       return NEWTON;
00259     }

void IncludeFastRelativisticEffects const bool  b  )  [inline]
 

Definition at line 300 of file orsa_interaction.h.

Referenced by OrsaFile::Read().

00300                                                       {
00301       include_fast_relativistic_effects = b;
00302     }

void IncludeMultipoleMoments const bool  b  )  [inline]
 

Definition at line 271 of file orsa_interaction.h.

Referenced by OrsaFile::Read().

00271                                                {
00272       include_multipole_moments = b;
00273     }

void IncludeRelativisticEffects const bool  b  )  [inline]
 

Definition at line 291 of file orsa_interaction.h.

Referenced by OrsaFile::Read().

00291                                                   {
00292       include_relativistic_effects = b;
00293     }

bool IsIncludingFastRelativisticEffects  )  const [inline]
 

Definition at line 303 of file orsa_interaction.h.

Referenced by OrsaFile::Write().

00303                                                     {
00304       return include_fast_relativistic_effects;
00305     }

bool IsIncludingMultipoleMoments  )  const [inline]
 

Definition at line 274 of file orsa_interaction.h.

Referenced by OrsaFile::Write().

00274                                              {
00275       return include_multipole_moments;
00276     }

bool IsIncludingRelativisticEffects  )  const [inline]
 

Definition at line 294 of file orsa_interaction.h.

Referenced by OrsaFile::Write().

00294                                                 {
00295       return include_relativistic_effects;
00296     }

bool IsSkippingJPLPlanets  )  const [inline, inherited]
 

Definition at line 98 of file orsa_interaction.h.

References Interaction::skip_JPL_planets.

Referenced by Evolution::Integrate().

00098                                       {
00099       return skip_JPL_planets;
00100     }

double PotentialEnergy const Frame  )  [virtual]
 

Implements Interaction.

Definition at line 601 of file orsa_interaction.cc.

References Vector::IsZero(), Vector::Length(), and ORSA_WARNING.

Referenced by JPLPlanetsNewton::PotentialEnergy(), and Relativistic::PotentialEnergy().

00601                                                {
00602     
00603     if (f.size() < 2) return(0.0);
00604     
00605     double energy = 0.0;
00606     
00607     unsigned int i,j;
00608     
00609     Vector d;
00610     
00611     double l;
00612     
00613     for (i=1;i<f.size();++i) {
00614       
00615       if (f[i].mu()==0) continue;
00616       
00617       for (j=0;j<i;++j) {
00618         
00619         if (f[j].mu()==0) continue;
00620         
00621         d = f[i].DistanceVector(f[j]);
00622         
00623         l = d.Length();
00624         
00625         if (d.IsZero()) {
00626           ORSA_WARNING("two objects in the same position! (%s and %s)",f[i].name().c_str(),f[j].name().c_str());
00627           continue;
00628         }
00629         
00630         energy -= f[i].mu()*f[j].mass()/l;
00631         
00632       } 
00633     }  
00634     
00635     return (energy);
00636   }

Here is the call graph for this function:

void SkipJPLPlanets const bool  b  )  [inline, inherited]
 

Definition at line 95 of file orsa_interaction.h.

References Interaction::skip_JPL_planets.

Referenced by OrsaFile::Read().

00095                                       {
00096       skip_JPL_planets = b;
00097     }


Member Data Documentation

bool skip_JPL_planets [protected, inherited]
 

Definition at line 102 of file orsa_interaction.h.

Referenced by Interaction::IsSkippingJPLPlanets(), and Interaction::SkipJPLPlanets().


The documentation for this class was generated from the following files:
Generated on Tue Jan 11 15:29:17 2005 for liborsa by  doxygen 1.4.0