External field effect of modified Newtonian dynamics in the Solar system

0.34 MB·English

by  Luc Blanchet

Mon.Not.R.Astron.Soc.000,1–15(2010) Printed4January2011 (MNLATEXstylefilev2.2)
External field effect of modified Newtonian dynamics
in the Solar system
1 Luc Blanchet,1⋆ J´erome Novak2†
1
1 RεC , Institut d’Astrophysique de Paris — UMR 7095 du CNRS, Universit´e Pierre et Marie Curie,
0 9G8bis boOulevard Arago, 75014 Paris, France
2
2Laboratoire Universet Th´eories, Observatoire de Paris, CNRS, Universit´e Paris Diderot,

n 5 place Jules Janssen, 92190 Meudon, France
a
J

3
4January2011

]
O
ABSTRACT
C
The Modified Newtonian Dynamics (MOND) has been formulated as a modification
h. of the Poissonequation for the Newtonian gravitationalfield. This theory generically
p predicts a violation of the strong version of the equivalence principle, and as a result
- the gravitational dynamics of a system depends on the external gravitational field in
o
which the system is embedded. This so-called external field effect has been recently
r
t shown to imply the existence of an anomalous quadrupolar correction, along the di-
s
rection of the external galactic field, in the gravitational potential felt by planets in
a
[ the Solar System. In this paper we confirm the existence of this effect by a numerical
integration of the MOND equation in the presence of an external field, and compute

2 the secular precession of the perihelion of planets induced by this effect. We find
v
that the precession effect is rather large for outer gaseous planets, and in the case of
9
Saturn is comparable to published residuals of precession obtained by Saturn range
4
tracking data. The effect is much smaller for inner planets, but in the case of the
3
Earthit appears to be in conflictfor most ofthe MOND functions µ(y) with the very
1
. good constraint on the perihelion precession obtained from Jupiter VLBI data. The
0
MONDfunctionsthatarecompatiblewiththisconstraintappeartohaveaveryrapid
1
transition from the MONDian regime to the Newtonian one.
0
1 Key words: planets andsatellites:general,galaxies:kinematicsanddynamics,meth-
:
v ods:numerical, dark matter
i
X
r
a
1 INTRODUCTION field and g = g its ordinary Euclidean norm. The modi-
| |
fication of the Poisson equation is encoded in the MOND
1.1 Motivation
function µ(y) of the single argument y g/a , where
0
The Modified Newtonian Dynamics (MOND) –Milgrom a0=1.2 10−10m/s2 denotestheMONDc≡onstantacceler-
×
(1983a,b,c)– is a successful alternative framework for in- ation scale. The MOND function interpolates between the
terpreting the galactic rotation curves and the empirical MOND regime corresponding to weak gravitational fields
Tully-Fisher relation without relying on dark matter halos y=g/a0 1, for which it behavesas µ(y)=y+o(y), and
≪
(see Sanders & McGaugh 2002 for a review). At the non- the Newtonian strong-field regime y 1, where µ reduces
≫
relativistic level, the best modified-gravity formulation of to 1 so that we recover theusual Newtonian gravity.
MONDisthemodifiedPoissonequationoriginallyproposed RelativisticextensionsofMONDmodifyinggeneralrel-
by Bekenstein & Milgrom (1984), ativityhavebeenproposed(Bekenstein2004;Sanders2005)
and extensively studied (see Bekenstein 2005 for a review).
g
∇ µ ∇U = 4πGρ, (1) Alternatively the MOND equation (1) can be reinterpreted
·(cid:20) (cid:18)a0(cid:19) (cid:21) − as a modification of dark matter rather than gravity by
where ρ is the density of ordinary (baryonic) matter, U is invoking a mechanism of “gravitational polarisation” — a
the gravitational potential, g = ∇U is the gravitational gravitationalanalogue oftheelectricpolarisation (Blanchet
2007);thisyieldstotheconceptofdipolardarkmatterwhich
hasbeenformulatedasarelativisticmodelinstandardgen-
eral relativity by Blanchet & Le Tiec (2008, 2009).
⋆ E-mail:[email protected]
E-mail:[email protected] An important consequence of the non-linearity of
†
2 Luc Blanchet & J´erˆome Novak
Eq. (1) in the MOND regime, is that the gravitational dy- by
namics of a system is influenced (besides the well-known
1
tidal force) by the external gravitational environment in δu= 2xixjQij, (3)
which the system is embedded. This is known as the exter-
nalfieldeffect(EFE),whichhasnon-trivialimplicationsfor where xi is the distance from the Solar System barycentre,
non-isolated gravitatingsystems.TheEFEwasconjectured and Qij is a trace-free quadrupolemoment of thetype
to explain the dynamics of open star clusters in our galaxy
1
(dMariklgmroamtte1r9d8e3sap,bit,ec)t,hesiinncveoltvheedywdeoaknionttesrhnoawlacecveildeernactieonosf Qij =Q2(cid:16)eiej− 3δij(cid:17), (4)
(i.e. below a0). The EFE effect shows that the dynamics of with Q2 being the quadrupole coefficient and e = (ei) the
these systems should actually be Newtonian as a result of preferred direction of the galactic centre (i.e. e = ge/ge).
theirimmersion in thegravitational fieldof theMilky Way. The anomalous term (3) is a harmonic solution of the
TheEFEisarigorouspredictionoftheBekenstein-Milgrom Laplace equation, i.e. ∆δu=0.
equation (1), and is best exemplified bythe asymptotic be- Theradialdependenceofthisanomaly is r2 andcan
∝
haviour of the solution of (1) far from a localised matter thus be separated from a quadrupolar deformation due to
distribution(say,theSolarSystem),inthepresenceofacon- theSun’soblatenesswhichdecreaseslike r−3.Thisresult
∝
stant external gravitational field g (the field of the Milky is valid wheneverr is much less than theMONDtransition
e
Way).At large distances r= x we have (Bekenstein distance for the Solar System, defined by r0 = GM/a0
| |→∞
& Milgrom 1984) with M the mass of the Sun and a0 the MONDpaccelera-
tion scale. This radius corresponds to the transition region
GM/µ 1
U =g x+ e + , (2) where the Newtonian acceleration becomes of the order of
e· r 1+λesin2θ O(cid:18)r2(cid:19) theMONDaccelerationa0andtherefore,MONDeffectsbe-
where M is the mass opf the localised matter distribution, comedominant.Wehaver0 7100AUsotheeffect(3)–(4)
≃
where θ is the azimuthal angle from the direction of the is valid in a large volume around the Sun including all the
external field g (see also Fig. 1), and where we denote planets (recall that Neptune’sorbit is at 30AU).
e
µ µ(y ) and λ y µ′/µ , with y = g /a and µ′ = In addition to (3) we have also the usual MOND ef-
dµe(≡y )/dye .Intheepr≡esenecee ofetheexterenalfieeld,0theMOeND fect.Inthespecialcaseofsphericalsymmetry,theequation
e e
internalpotentialu U ge xshowsaNewtonian-likefall- (1) reduces to µ(g/a0)g=gN, with gN the ordinary Newto-
off r−1 at large d≡istan−ces·but with an effective gravita- niangravitationalfield.ForaMONDfunctionbehavinglike
tciaosne∼a,litcoexnhstibanittsGa/nµoen.-1spHhoewriecvaelrd,ecfoonrtmraartyiotnoatlhoengNtehwetdonirieacn- tµh(ey)N∼ew1to−niǫayn−fiγelwdhiesnmyod→ifie∞db(yi.eδ.gNin thǫea0r(erg/ior0n)2rγ−≪2.rF0o)r,
≃
tion oftheexternalfield.Thefact thattheexternalfieldg γ =1 this corresponds to a constant acceleration similar to
e
doesnotdisappearfromtheinternaldynamicscanbeinter- the Pioneer anomaly, and for γ >2 the effect is very small
pretedasaviolationofthestrongversionoftheequivalence and the motion of planets is almost unaffected. This effect
principle.2 Forthereader’sconveniencewederivetheresult is independent of the presence of the external field, and is
(2) in the AppendixA. spherical, socan bedistinguished from thequadrupolarde-
In a recent paper, Milgrom (2009) has shown that the formation (3) induced by the external field; we neglect this
imprintoftheexternalgalacticfieldg ontheSolarSystem effect here.
e
(duetoaviolationofthestrongequivalenceprinciple)shows
up not only asymptotically, but also in the inner regions of
the system, where it may have implications for the motion 1.2 Summary
of planets. This is somewhat unexpected because gravity is
Inthepresentpaperweshallconfirmtheexistenceoftheef-
strong there (we have g a ) and the dynamics should
≫ 0 fect(3)–(4)andshallderiveitbyanindependentnumerical
beNewtonian.However,becauseofthespecialpropertiesof
integration of the equation (1) using an elliptic solver orig-
theequation(1),thesolutionwillbegivenbysomemodified
inally built for numerical relativity purposes.3 The magni-
Poisson-type integral, and the dynamics in the strong-field
tudeofthequadrupolecoefficientQ dependsonthedimen-
2
region will be affected by the anomalous behaviour in the
sionless ratio between the external field g and the MOND
e
asymptoticweak-fieldregion.Thustheresultsapplyonlyto
acceleration a , say
0
the nonlinear Poisson equation (1), not to modified inertia
g
formulations of MOND. The anomaly expresses itself as an η= e , (5)
anomalousquadrupolarcontributiontotheMONDianinter- a0
nal field u, as compared to the Newtonian potential, given and on the particular MOND interpolating function µ in
use.FortheMilkyWayfieldattheleveloftheSunwehave
g 1.9 10−10m/s2 which happens to be slightly above
1 RecallthatintheabsenceoftheexternalfieldtheMONDpo- e ≃ ×
theMOND scale, i.e. η 1.6. Ourcalculation implemented
tential behaves like U ∼ −√GMa0lnr, showing that there is for various standard cho≃ices for theµ-function gives
noescapevelocityfromanisolatedsystem(Famaeyetal.2007).
Howeversincenoobjectistrulyisolatedtheasymptoticbehaviour 2.1 10−27 s−2 .Q .4.1 10−26 s−2. (6)
2
ofthepotentialisalwaysgivenby(2),intheapproximationwhere × ×
theexternal fieldisconstant.
2 Howevertheweakversionoftheequivalenceprinciple,thatall
testparticlesintheMONDgravitational fieldhaveuniversalac- 3 The code is based on the library lorene, available from the
celerationa=g,remainssatisfied(Bekenstein&Milgrom1984). websitehttp://www.lorene.obspm.fr.
MOND effects in the Solar system 3
Thequadrupolecoefficient Q isfoundtobepositivecorre- ofgravity(belowa ),maynotbeextrapolatedwithoutmod-
2 0
spondingtoaprolateelongationalongthequadrupolaraxis. ification to the strong field of the Solar System. For in-
Weshall verify numerically that Q is indeed constant over stanceit hasbeen argued byFamaey & Binney(2005) that
2
alargeportionoftheinnerSolarSystem,andstartsdecreas- a MOND interpolating function µ which performs well at
ing when the distance becomes comparable to the MOND fittingtherotation curvesof galaxies isgiven byµ defined
1
transition scale r , at a few thousands AU from the Sun. by (34a) below. However this function has a rather slow
0
We shall also check that Q scales approximately like the transition to the Newtonian regime, given by µ 1 y−1
2 1
∼ −
inverse square root of the central mass, in agreement with when y = g/a , which is already excluded by Solar
0
→ ∞
dimensional analysis which shows that we should approxi- System observations (Sereno & Jetzer 2006). Indeed such
mately have Q = q a /r , where q (η) is a dimensionless slow fall-off y−1 predicts a constant supplementary accel-
2 2 0 0 2
−
quadrupolecoefficientdependingontheratioηbetweenthe eration directed toward the Sun δg =a (i.e. a “Pioneer”
N 0
externalfieldandtheMONDacceleration. Wefindthatfor effect), which is ruled out because not seen from the mo-
η 1.6, q ranges over0.02.q .0.36, in good agreement tionofplanets.Thusitcouldbethatthetransitionbetween
2 2
≃
with thevalues found by Milgrom (2009). MONDandtheNewtonianregimeismorecomplicatedthan
We then study some consequences of the anomalous ismodelledinEq.(1).Thisisalsotrueforthemodifieddark
quadrupole moment for the motion of planets, notably the matter model (Blanchet & Le Tiec 2008, 2009) which may
secular advances of planetary perihelia. Applying standard onlygiveaneffectivedescriptionvalidintheweakfieldlimit
perturbation equationsof celestial mechanics weobtain the andcannotbeextrapolatedasitstandstotheSolarSystem.
anomalous precession rate Whilelooking atMOND-likeeffectsintheSolarSystemwe
shouldkeepthepreviousprovisoinmind.Thepotentialcon-
dω˜ Q √1 e2
∆2 =(cid:28)dt(cid:29)= 2 4n− h1+5cos(2ω˜)i. (7) flwiictthwtheefiEnadrthheroerbwititahlptrheeceSsosliaorn)Smysatyemnodtynneacmesiscasri(lnyoitnavbally-
Tosimplify,wehaveassumedthatthedirectionofthegalac- idate those theories if they are not “fundamental” theories
ticcentreliesintheplaneoftheecliptic(thisisonlyapprox- but rather“phenomenological” models.
imatelycorrect)andthatallplanetsmoveinthisplane.We The plan of this paper is the following. In Sec. 2, we
defined the origin of the precession angle ω˜ to be the di- derivean expression of thesolution of theMOND equation
rection of the galactic centre. The orbital frequency of the near the origin in terms of multipole moments. Section 3 is
planet is n = 2π/P (P is the orbital period), and a, e and devotedtonumericaltechniquesandresults,withadescrip-
ω˜ =ω+Ωdenotethestandard orbital elements. Theeffect tion oftheparticularformulation andassumptions usedfor
seems to be the most interesting for Saturn for which we the numerical integration in Sec. 3.1, of the various tests
find passedbythecode(Sec.3.2),andthenumericalresultsand
values for the first multipoles of the potential in Sec. 3.3.
0.3mas/cy.∆ .5.8mas/cy. (8)
2 Section 4 details the consequences of this modified gravita-
Those values are within the range of published residuals tional potential on the orbits of Solar-system planets, with
for the precession of Saturn as obtained from global fits of firstthederivationoftheperturbationequationsinSec.4.1
the Solar System dynamics in Pitjeva (2005) and Fienga andnumericalvaluesfortheplanetaryprecessions(Sec.4.2).
et al. (2009). However we shall find that in the case of the Finally, Section 5 summarises the results and gives some
Earth the predicted perihelion precession may be greatly concluding remarks.
constrained by the best estimates of postfit residuals ob-
tained thanks to the Jupiter VLBI data using the INPOP
planetary ephemerides(Fienga etal. 2009). Thisis turnre- 2 MULTIPOLAR EXPANSION OF THE
duces the possible choices for the MOND function µ(y) to MONDIAN POTENTIAL
thosethatexhibitasharptransitionfromtheNewtonianto
InthispaperweshallsolvethemodifiedPoissonequation(1)
MONDian regime.
withtheboundaryconditiongivenbytheconstantexternal
To conclude, we present very accurate numerical com-
gravitational field g (defining a preferred spatial direction
putationsperformedwithinthewell-definedtheory(1),and e
denotede=g /g ),i.e.thatthegravitationalfieldg=∇U
comparetheresultswithrecentobservationsofthemotionof e e
planetsin theSolar System.Wefindthat theobservational asymptotes to limr→∞g =ge. The external field ge should
consistentlyobeyaMONDequation,buthereweshallsim-
constraints are rather strong, and may even conflict with
plyneedtoassumethat g isconstant overtheentireSolar
someofthepredictions.Thusourresultsprovidemotivation e
System.TheMONDianphysicistmeasuresfromthemotion
forinvestigatingmoresystematicallypossibleanomalousef-
of planets relatively to the Sun the internal gravitational
fects in the Solar System predicted by alternative theories.
potential u definedby
In particular, the external field effect associated with the
preferred direction of the galactic centre could be seen as u=U g x, (9)
e
− ·
a typical one arising in generic attempts to modify gravity
with motivation coming from dark matter and/or dark en- whichissuchthatlimr→∞u=0.Contrarytowhathappens
in theNewtonian case, theexternalfield g does not disap-
ergy.Indeedweexpectthatgenericdeparturesfromgeneral e
pearfromthegravitationalfieldequation(1)andwewantto
relativity will violate the strong equivalence principle (Will
investigatenumericallyitseffect.Theanomalydetectedbya
1993).
NewtonianphysicistwithrespecttotheMONDianphysicist
OntheotherhandletuscautiouslyremarkthatMOND
is thedifference of internal potentials,
and more sophisticated theories such as TeVeS (Bekenstein
2004), which are intendedto describe theweak field regime δu=u u , (10)
N
−
4 Luc Blanchet & J´erˆome Novak
where u denotes the ordinary Newtonian potential gener- Because the integration in (15) is limited to the domain
N
ated by the same ordinary matter distribution ρ, and thus r > ε and ∂ (1/r) is symmetric-trace-free (STF) there [in-
L
solution of the Poisson equation ∆u = 4πGρ with the deed ∆(1/r) = 0], we deduce that the multipole moments
N
−
boundaryconditionthatlimr→∞uN =0.HenceuN isgiven QL themselves are STF. This can also be immediately in-
by thestandard Poisson integral ferred from the fact that ∆δu = 0 when r 6 ε, hence the
d3x′ multipole expansion (14) must be a homogeneous solution
uN(x,t)=GZ x x′ ρ(x′,t). (11) of the Laplace equation which is regular at the origin, and
| − | is therefore necessarily made solely of STF tensors of type
A short calculation shows that the anomaly obeys the or- xˆL. Hence we can replace xL in (14) by its STF projection
dinary Poisson equation ∆δu = −4πGρpdm, where ρpdm is xˆL. It is clear from the formula (15) that the MONDian
thedensity of “phantom dark matter” definedby
gravitational field (for r > r ) can influence the near-zone
0
ρ = 1 ∇ (χ∇U) , (12) expansion of thefield when r 0.
pdm 4πG · With the expression (12)→for thephantom dark matter
wherewe denoteχ µ 1.The phantomdark matterrep- and the MOND equation (1), we can further transform Q
≡ − L
resents the mass density that a Newtonian physicist would as
attribute to dark matter. In the model by Blanchet & Le
1 1
Tasietche(2d0e0n8s,it2y00o9f)ptohlearpishaatinotnomofdsaormkemdaipttoelrarisdianrtkermpraettteedr QL=−4π Zr>εd3xh4πGρ+∆Ui∂L(cid:18)r(cid:19) . (16)
medium and the coefficient χ represents the “gravitational Approximating the central matter distribution as being
susceptibility” of thisdark matter medium. spherically symmetric (i.e. ignoring the “back-reaction” of
ThePoissonequation∆δu= 4πGρ istobesolved the non-spherical anomalous field δu on the matter which
pdm
−
withtheboundaryconditionthatlimr→∞δu=0;hencethe generates the field), we see that the first term is non-zero
solution is given by thePoisson integral only in the monopolar case l = 0 where it reduces in the
δu(x,t)=GZ xd3xx′′ ρpdm(x′,t). (13) tlihmeitorεig→in.0Otno mthienuosthtehrehNaenwdtothneiasnecpoontdentteiarml evinalu(1a6t)edcaant
| − | betransformed as a surface integral. Wefind
Weemphasiseherethat,contrarytotheNewtonian(linear)
1 r=∞ 1 1
case, the knowledge of thematter density distribution does Q = u (0)δ + d2S U∂ ∂ U∂ .
not allow to obtain any analytic solution for the potential; L − N l,0 4π Zr=ε ih iL(cid:18)r(cid:19)− i L(cid:18)r(cid:19)i
(17)
however, we can still infer the structure of the multipolar
Ournotationmeansthatthesurfaceintegraliscomposedof
expansion near the origin, and the moments will be com-
two contributions, an inner one at r =ε (denoted Qε) and
puted numerically. We can check that the phantom dark L
matter behaves like r−3 when r , so the integral (13) an outer one at infinity r = ∞ (say Q∞L). In Eq. (17) we
→ ∞ implicitly assume that the limit ε 0 is to be taken after
is perfectly convergent.
→
evaluating theinner surface integral.
Wewanttodiscusstheinfluenceoftheexternalgalactic
Inserting U = u+g x into the surface integral at
fieldintheinnerpartoftheSolarSystemwherethegravita- infinity Q∞ we find that iet·contributes only to the dipolar
tionalfield isstrong(g≫a0);thusµtendstoonethere,so termandrLeducestotheexternalgalacticfield;thusQ∞ =gi
χ tends to zero. For thediscussion in this section we adopt i e
theextremecasewhereχisexactly zeroinaneighbourhood withallotherQ∞L’sbeingzero.Ontheotherhandtheinner
surfaceintegralQε isfoundtohaveawell-definedlimitwhen
of the origin, say for r 6 ε, so that there is no phantom L
ε 0 given by Qε = ( )l(∂ˆ U)(0), where we recall that
dark matter for r 6 ε; in later sections devoted to the full ∂ˆ→denotes the STLF par−t of ∂L = ∂ ∂ . We then find
numerical integration we shall still make this assumption L L i1··· il
that the galactic field in U = u+g x cancels the dipole
by posing χ = 0 inside the Sun (in particular we shall al- e
·
term in thesurface integral at infinity,so that theresult is
ways neglect the small MOND effect at the centre of the
Sun where gravity is vanishingly small). If ρpdm = 0 when QL= uN(0)δl,0+( )l(∂ˆLu)(0). (18)
r6εwecandirectlyobtainthemultipolarexpansionofthe − −
anomalous term (13) about the origin by Taylor expanding TheinternalMONDianpotentialuadmitsthereforethefol-
theintegrand when r= x 0. In thisway we obtain lowingSTFmultipolarexpansion(equivalenttoanear-zone
| |→ expansion when r 0)
+∞( )l →
δu=Xl=0 −l! xLQL, (14) u=uN−uN(0)++∞l1!xL(∂ˆLu)(0). (19)
where themultipole momentsnear the origin are given by4 Xl=0
InAppendixBwederivean alternativeproofofthisresult.
1
QL =GZ d3xρpdm∂L(cid:18)r(cid:19) , (15) Atthisstage we haveelucidated thestructureof themulti-
r>ε pole expansion of the anomaly δu near the origin, but still
weneedtoresorttoanumericalintegrationofthenon-linear
4 Our notation is as follows: L = i1···il denotes a multi-index MOND equation in order to obtain quantitative values for
composedoflmultipolarspatialindicesi1,···,il(rangingfrom1
t∂o/∂3x);i;∂LxL==∂i1xi·1··∂ilxiilsitshethperopdruocdtucotflofpalrstpiaaltidaelripvoastiitvieosns∂ix≡i; ∂ˆL.Thedecompositionof∂L intermsofSTFcomponents ∂ˆL is
···
similarlynL=ni1 nil =xL/rl istheproductoflunitvectors given by (B2)–(B3) below. In the case of summed-up (dummy)
···
ni=xi/r;thesymmetric-trace-free(STF)projectionisindicated multi-indices L, we do not write the l summations from 1 to 3
withahat,forinstancexˆL STF[xL],andsimilarlyfornˆL and overtheirindices.
≡
MOND effects in the Solar system 5
3 NUMERICAL INTEGRATION OF THE
z
MOND EQUATION
3.1 Formulation and implementation
M g
e
Using spherical spatial coordinates r,θ,ϕ , we consider a
{ }
θ star represented by a spherically symmetric distribution of
x
matter ρ(r) obtained from a hydrostatic equilibrium model
inNewtoniantheory(polytrope).Thesemodelsareobtained
by integrating the hydrostatic spherical equilibrium equa-
r tion(relatingthepressurepandthegravitationalfield)and
the Newtonian equation for the gravitational field (relat-
ing the gravitational field and the density distribution ρ),
together with a polytropic equation of state (EOS) of the
formp=Kργ,whereK andγ aretwoconstants.Asitshall
beshownlaterin Sec.3.2, theresultsdonotdependon the
y
specificEOSusedtoobtain thishydrostaticequilibrium,so
wedonotdiscusstheparticularEOSusedforobtainingthe
density distribution. It also means that we neglect the ef-
fect of MOND theory on the structure of the star. This is
ϕ not a severe restriction, since inside the star the gravita-
tional field is much higher in amplitude than the constant
a , making µ 1 and Newtonian theory is recovered up to
0
x ∼
veryhighaccuracy.WethensolvetheMONDequation(1),
withtheboundaryconditionsgivenbytheconstantgalactic
Figure 1. Setting of the spherically symmetric star and the
asymptoticgalacticfieldge,usingsphericalcoordinates r,θ,ϕ . gravitational field ge (see Fig. 1),i.e. limr→∞g=ge=gee.
{ } In order to be closer to a Poisson-like form of the partial
differential equation (PDE), we are solving in terms of the
internal potential u = U ge x such that limr→∞u = 0.
Definingh=∇u=g g −,theM· ONDequation(1)becomes
e
−
thefollowing PDE:
the multipole moments. Section 3 will be devoted to this 1 µ′(g/a )
∆u= 4πGρ 0 gigj∂ h , (22)
task. µ(g/a0)(cid:26)− − a0g i j(cid:27)
FinallyweshowthatthedipolemomentQ (withl=1)
i
isactually zero(Milgrom 2009).Thisfollows from theweak withµ′ dµ/dy beingthederivativeofthefunctionµwith
≡
equivalence principle satisfied by the Bekenstein-Milgrom respect to its argument y=g/a0.
theory. As a consequence, the total acceleration of the cen- This PDE is solved numerically, using the library
treofmassoftheSolarSysteminthegalacticexternalfield lorenewhichimplementsspectralmethods(forareviewin
g doesnotdependontheSolarSystem’sinternaldynamics the case of numerical relativity see Grandcl´ement & Novak
e
andissimplygivenbyg (foradetailedproofseeBekenstein 2009) in spherical coordinates. In our case (see Fig. 1), the
e
& Milgrom 1984). The centre of mass of the Solar System fieldsdonotdependontheazimuthalangleϕandtheprob-
does not differ much from that of the Sun, so we deduce lem is axisymmetric. Since Eq. (22) is a non-linear elliptic
that the “self-acceleration” of the Sun, i.e. the total accel- PDE,thealgorithmusedisthefixedpointiterationmethod
erationduetotheinternalpotentialuwhenintegratedover in which one starts from an initial guess u0(r,θ), here the
thewholevolumeoftheSun,mustbe(approximately)zero: solution of theNewtonian Poisson equation∆u0 = 4πGρ.
−
Knowing u (r,θ) at a given iteration step n, one com-
n
putes the non-linear source terms in the right-hand-side of
d3xρ∂ u=0. (20) Eq.(22), say σ(u,ρ),toobtain a newvalueof thepotential
i
Z
u , solving a linear Poisson equation
n+1
WehaveobviouslythesameresultfortheNewtonianpoten- ∆u˜n+1 =σ(un,ρ), (23)
tialu inNewtoniangravitysothesamewillbetrueforthe
N and using a relaxation technique
anomalous field δu = u u . In other words the phantom
N
−
dark matter which is the source for the anomaly exerts no u =λu˜ +(1 λ)u , (24)
n+1 n+1 n
−
netforceontheSunandweget,forasphericallysymmetric
star, with λ (0,1] being the relaxation parameter (often taken
∈
tobe0.5).Thisiterationisstoppedwhentherelativevaria-
tionofthepotentialubetweentwosuccessivestepsbecomes
Q = 1 d3xρ∂ δu=0. (21) lower than a given threshold (usually 10−12).
i −M Z i The linear equation (23) is solved decomposing each
field on a truncated base of spherical harmonics Ym(θ,ϕ)
l
In Sec. 3.3 we shall numerically verify that the dipole mo- (with m = 0 because the problem is axisymmetric) in
ment Q is indeed zero within our numerical error bars. the θ-direction, and Chebyshev polynomials T (r) in the
i n
6 Luc Blanchet & J´erˆome Novak
r-direction. For this last coordinate, a multi-domain tech- reads as5
nique is used, with linear mappings of the r coordinate to +∞ ( )l
theinterval[−1,1],exceptforthelastdomainwhereamap- δu(r,θ)= (2l−1)!!rlQl(r)Pl(cosθ), (28)
ping of the type s = 1/r allows to treat spatial infinity in Xl=0 −
our computational domain and to impose the right bound-
where P is the usual Legendre polynomial and θ is defined
l
ary conditions at r (see Grandcl´ement & Novak 2009
→∞ in Fig. 1. Although from the considerations of Sec. 2 the
fordetails).Moreover,thismulti-domaintechniqueimposes
multipole moments Q should be approximately constant
l
that thepotential u, as well as its first radial derivative, be
within the MOND transition radius r , here we compute
0
continuousacross the stellar surface.
them directly from the numerical solution of (1) and shall
Finally, the only modification made to the exact equa-
obtain their dependence on r; we shall check numerically
tion(22)inordertobenumericallyintegratedisthesetting
that Q (r) and Q (r) are indeed almost constant in a large
2 3
µ(g/a )=1 inside thestar. Indeed,thereare two problems
0 spheresurroundingtheSolarsystem.Withourdefinitionthe
here.First, theMOND gravitational field g(r)does not ad-
monopole, quadrupole and octupole pieces in the internal
mit a regular Taylor expansion near the origin, as does the
field are given by
density in Eq. (B6) and the associated Newtonian gravita-
tional field which follows as δu1 = rQ1 cosθ, (29a)
−
1 1
δu = r2Q cos2θ , (29b)
+∞ r2k+1 2 2 2(cid:18) − 3(cid:19)
g =4πG (∆kρ)(0). (25)
N Xk=0(2k+3)(2k+1)! δu3 =−61r3Q3 cosθ(cid:18)cos3θ− 35cosθ(cid:19) . (29c)
From dimensional arguments, it is possible to see that Q
Thiscanbeimmediatelyseenfromthefactthatinspherical shouldscaleasQ a /rl−1 [see forinstanceEqs.(33) andl
symmetry the MOND equation (1) reduces to µ(g/a0)g = (37)] and thereforle∼6 0 0
g . Using (25) we see that the expansion of g when r 0
N
→
starts with a term proportional to the square root of r and r l
δu Qrl a r . (30)
is therefore not regular; more precisely we have ∼Xl l ∼ 0 0Xl (cid:18)r0(cid:19)
This last expression shows that higher-order multipoles
g= 4πGρ0a0r 1+ (r2) . (26) should have lower influence on Solar-system planets, for
r 3 h O i which r r0. This shall be exemplified with the influence
≪
of the octupole, as compared to that of the quadrupole, in
With our polynomial decomposition of fields in the radial Sec. 4.2.
direction, where we decompose in Chebyshev polynomials
T (r), we therefore get an incompatibility when trying to
n 3.2 Tests of the code
solve the Poisson equation (23) in the vicinity of the origin
atthecentreoftheSun.Thesecondproblemisthat,atthe Afirst test of ourcodeis thestandard comparison between
very centre of the star g 0 and therefore µ 0, which the Laplace operator applied to the potential u, and the
→ →
makesthedivisionpresentintheright-hand-sideofEq.(22) right-hand-sideofEq.(22).Inallresultsdisplayedhere,the
difficult tohandle numerically. maximumofthiserrorhasalwaysremainedlowerthan10−11
Here we have chosen to avoid both problems by im- andshallnotbeconsideredasaninterestingerrorindicator.
posing µ = 1 in the star so it is entirely Newtonian. This Next, an indication of the error comes from the use of the
approximationmayproduceanoticeablechangeinthevalue Gauss theorem:
wofhµe(rye)ρo0nilsytwhiethciennatrsaplhdeernesoiftyraadniudsRr⊙.is4π3tGah0ρe0r∼ad1iu0s−1o5fRth⊙e, I µ(cid:18)ag0(cid:19)g·d2S=−4πGM, (31)
Sun. It is thus supposed to have a completely negligible ef-
whereM isthemassofthestar,computedfrominitialdata
fect on theresults.
ataveryhighprecision.Then,anothercheck(whichishow-
The result which is finally used to study the influence
evernot a priori independent) is the asymptotic behaviour
of the modification of the Newtonian gravity on the orbits
of thepotential u(r,θ) when r ,which is
of planets is the value of the trace-free multipole moments →∞
QL defined in Sec. 2. In the case where all the multipole GM 1
u= + , (32)
moments are induced by the presence of the external field rµ 1+λ sin2θ O(cid:18)r2(cid:19)
g inthepreferreddirectione,thesituationisaxisymmetric e e
e p
andallthemomentsQ willhavetheiraxispointinginthat
L
direction e, and we can define the multipole coefficients Ql 5 Here, the expansion is defined for any range in radius r, con-
by trarytoSec. 2andAppendixB,wherethemultipoleswereonly
numbersdefinedforr 0andnotfunctions ofr.
6 Theradialdependen→ceofthecontributionofthel-thmultipole
QL=QleˆL, (27) moment isthe sameas that of the usual MONDspherical effect
corresponding to a MOND function behaving like µ 1 ǫy−γ
∼ −
when y , with γ = (l+1)/2 (see the discussion at the end
where eˆL denotes the STF part of the product of l unit of Sec. 1→.1)∞. But of course the effect we are considering here is
vectors eL = ei1 eil. Then the multipole expansion (14) non-spherical.
···
MOND effects in the Solar system 7
µ(y) = y / (1+y2)1/2
0.01
Test using the Gauss theorem at infinity
Test on the asymptotic behavior of u(r,θ)
Equation of state 1
Equation of state 2
0.001 10000
or
Relative err0.0001 3nsity [kg / m] 1100000
e
d
10
1e-05
0.1 1 10
η = g / a
e 0
1
0 1e+08 2e+08 3e+08 4e+08 5e+08 6e+08
radius [m]
Figure2.Relativeerrorofthecodeonthetwotestsdescribedin
Sec.3.2,fora1M⊙ sphericalmassdistribution,forthestandard
Figure3.Twodensityprofilesfortestingtheresultdependence
MOND function µ2(y) [Eq. (34b)] and the standard value a0 =
1.2 10−10m.s−2. The abscissa gives the ratio (5) between the ontheequationofstate.Bothprofilesyieldaone-solarmassstar
× inhydrostaticequilibrium.Thetwoprofilesyieldverylittle(same
norm of the asymptotic gravitational field ge and the parameter levelastheoverallnumericalaccuracy 10−3)differenceonour
a0.Notethatbothcurvesdocoincide. ∼
results.
3.3 Results on the induced quadrupole and
where µ = µ(y ) and λ = y µ′/µ with y = g /a (see octupole
e e e e e e e e 0
the proof in Appendix A). The two tests are not really in-
Inwhatfollows,unlessspecified,themassofthestaristhat
dependent because the Gauss theorem is obtained by inte-
of the Sun. As a first result, we show in Fig. 4 the profile
grating the asymptotic field over a sphere at infinity. The
of the quadrupole induced by the MOND theory through
main interest of these tests is that while the field is com-
the function Q (r) defined in Eq. (29b). We find that this
puted asymptotically on a sphere at infinity, the mass M 2
quadrupole is decreasing from the star’s neighbourhood to
is obtained “locally” by integrating numerically thedensity
zero, on a typical scale of 10000 astronomical units (AU).
overthevolumeof thestar. Weemphasise that in ourcode
However, this quadrupole can be considered as almost con-
the correct asymptotic behaviour (32) comes out directly
stant closer to the Sun, as it has a relative variation lower
and needs not to be imposed by hands at the beginning of
than 10−4 within 30 AU (see the zoomed region in Fig. 4).
thecalculation.
Weshallthereforerefertothequadrupoleasasimplenum-
From the use of the mapping s = 1/r (see Sec. 3.1),
ber,notedQ (0)orsimplyQ ,whenevaluatingitsinfluence
it is numerically possible to check the angular dependence 2 2
ontheorbitsofSolar-systemplanets.Ournumericalresults
of the quantity r u when r and thus to estimate
× → ∞ for thequadrupole are given in Table 1.
independently the error on the potential u. Both tests —
Wenowbrieflystudyothermultipolemoments,namely
Gausstheorem(31) andasymptoticbehaviour(32)—have
for l = 1 and l = 3. The profiles for the dipole Q (r) and
been performed and have given accuracy levels of 10−3– 1
10−4. Some examples are given in Fig. 2, for a s∼equence octupoleQ3(r),definedby(29a)and(29c),aredisplayedin
Fig. 5. Wefind that near the Sun and the solar planets the
of different values for the asymptotic gravitational field g .
e dipole can be considered as zero, since it is three orders of
Both tests give the same numbers, up to five digits. In all
magnitudelowerthanthetypicalvalueforanaccelerationin
results shown hereafter, the accuracy level defined by these
two tests is better than 10−2. theproblem(i.e.a0 orge).Indeedweexpectondimensional
analysis that Q should scale with the MOND acceleration
Among the parameters entering the numerical model, 1
a [see (30)],
wehavecheckedthatthedetailsofthedensityprofileρ(r)do 0
notcontributetothevalueofthequadrupolemomentQ2(r). Q1 =a0q1(η), (33)
As an illustration, we give here the quadrupole obtained
with two differentdensity profiles, displayed in Fig. 3, both where the dimensionless coefficient q depends on the ratio
1
givingastarofexactlyonesolarmass.Withthesameother η betweeng anda asdefinedbyEq.(5).Ourvaluesgiven
e 0
parameters[choicesofµ(y),a andg ],therelativedifference in Table 1 show that q is actually very small. It means
0 e 1
in theinducedquadrupoleis4 10−3,while theerror indi- that Q is lower than the numerical error and confirms the
1
cators give a relative numerica×l uncertainty of 2 10−3. analytical proof given in Sec. 2 that the dipole moment is
∼ ×
Wethereforeconsiderthattheparticularformofthedensity approximately zero. On the other hand, we find that the
profile used to model the star (i.e., the choice of EOS) has octupoletendstoanon-zerovalue,becausethecorrespond-
little influenceon theresults. ing q , which has also an influence on the orbits of planets
3
8 Luc Blanchet & J´erˆome Novak
Figure4.Leftpanel:profileofQ2(r)inthesolarsystem,forastandardchoiceoffunctionµ1(y)[seeEq.(34a)],a0=1.2 10−10m.s−2
andge=1.9 10−10 m.s−2.TheMONDtransitionradiusr0=pGM/a0 isshownbyadash-dottedlineatr 7100AU×.Rightpanel:
× ≃
zoom ofthecentral region(r650AU),wherethequadrupoleisalmostconstant (thespikesforr 0areduetothehighvalueofthe
→
monopoleterm,fromwhichitisnumericallydifficulttoextractthequadrupole). Notethedifferenceinthey-axisscales.
Figure5.ProfilesofdipoleQ1(r)(left)andoctupoleQ3(r)(right)inthesolarsystem,withsamesettingsasthoseofFig.4.TheMOND
transition radius r0 =pGM/a0 isshown by a dash-dotted lineat r 7100 AU. Numericallywe get Q1(0) 9 10−14m.s−2 a0
≃ | |∼ × ≪
whichshowsthatthedipoleisinfactzero.
(see the next section), is found to be of the order of one. general notation for any positive integer n:
Numerical valuesof theoctupole arealso given in Table 1. y
µ (y)= . (35)
The dependence of the quadrupole moment Q2 upon n √n1+yn
the value of the galactic gravitational field is displayed in
We include also the function µ having an exponentially
exp
Fig. 6, for four different coupling functions µ(y). Here we
fast transition to the Newtonian regime. The fourth choice
consider four cases widely used in the literature:
µ ismotivatedbythetheoryTeVeS(Bekenstein2004).
TeVeS
y Oneshouldnotethatnoneofthesefunctions,exceptmaybe
µ (y)= , (34a)
1 1+y the fourth one, derives from a fundamental physical princi-
y ple.
µ (y)= , (34b)
2 One observes that, up to the standard value of g =
1+y2 e
1.9 10−10 m.s−2 thefourcurvesare quiteclose, giving an
µ (y)=1p e−y, (34c) ×
exp − intervalof valuesfor Q2:
µTeVeS(y)= √√11++44yy−+11. (34d) 2.2×10−26 s−2 .Q2.4.1×10−26 s−2. (36)
However, for other choices of the MOND function µ(y), Q
2
Thefunctionµ hasbeenshowntoyieldgoodfitsofgalactic cantakelowervaluesdownto2.1 10−27 s−2.Inallcases,we
1
×
rotation curves (Famaey & Binney 2005). However because findthatQ ispositive,whichmeansaprolatedeformation
2
of its slow transition to the Newtonian regime it is a priori of the field toward the galactic centre. Two profiles of µ(y)
incompatible with Solar System observations. The function seem to give a maximum for Q near the standard value
2
µ isgenerallycalledthe“standard”choiceandwasusedin of g , namely the standard choice µ and the exponential
2 e 2
fitssuchasthoseofBegeman etal.(1991). Wehereusethe choice µ . For the two other choices of µ(y), we have not
exp
MOND effects in the Solar system 9
a = 1.2 x 10-10 m/s2
1e-24 0 a = 1.2 x 10-10 m/s2, µ(x) = x/(1+x)
0
µ1(y) 1e-25
µ(y)
2 Numerical results
µexp(y) M(-1/2) law
-2Quadrupole Q(0) [s]211ee--2265 µTeVeS(y) -102=1.9x10 m/s) -2Quadrupole Q(0) [s]2468eee---222666
1e-270.1 η = 1ge/a0 standard value (ge 10 2e-26 2 4 6 8
Mass of the star [solar masses]
Figure 6. Quadrupole in the solar system Q2 Q2(0) as a
function of the external galactic gravitational fiel≡d ge, for four Figure 8. Quadrupole moment Q2 as a function of the stellar
different expressions of the couplingfunction µ(y) [seeEqs. (34) massM.TheM−1/2-lawisrecoveredasexpectedfromEq.(37a).
foradescription].Thestandardvalueofge=1.9 10−10 m.s−2
×
isshownbyathindash-double-dotted line.
givevalues of thedipole, quadrupoleand octupole near the
Sun (r 6 50 AU), where they are observed to be constant,
a = 1.2 x 10-10 m/s2 for the four different expressions (34) of the interpolating
0
function.
We have also explored different values of a at fixed
0
g , with the results that Q was increasing for small a as
e 2 0
a power-law, reaching a maximum value for a 10−10
0
-2Quadrupole Q(0) [s]21e-26 ts1vinhia0ner−cFi9eeqimdguo.ant·dh8srge−uew2pmnhoeaaelrenrsaesdlm,ogotifrnhmoteuehpnnneadtsrscltioiescwonurntlelyreacaroled,vxesetpctrhraeeeercda.t.sMsRiTnte−hghs1i.aus/tlF2tistsindhnaaeeorlplteqye≃u,nadadwisdseuepnrrulcpahpeyraoievso−ldeeef
moment (and similarly the octupole moment) should scale
like (Milgrom 2009)7
a
Q = 0 q (η), (37a)
1e-27 10 20 2 r0 2
index n a
Q = 0 q (η), (37b)
3 r2 3
0
Figure 7. Quadrupole moment Q2 as afunction of the index n where r0 = GM/a0 is the MOND transition radius and
oftheMONDfunctionµn(y)definedbyEq.(35). q2, q3 denotepsome dimensionless coefficients depending on
the ratio η = g /a and on the choice of the interpolating
e 0
function µ. Having obtained the expected dependence on
been able to increase the value of η = g /a further than
e 0 the mass in Fig. 8 is another check that our code behaves
10,becausetheerrorindicatorswouldbecometoolargeand
correctly, since we are able to recover the analytic results
the results could not be trusted any longer. Note that in
knownforthissystem.Thevaluesofthedimensionlessmul-
Fig. 6, a has been kept fixed, while we have varied g . We
0 e tipole coefficients q and q , which are expectedly close to
have further explored the dependence of the EFE induced 2 3
unity,are reported in Table 1.
quadrupoleQ on thetypeof MOND function.
2
We have used several different functions of typeµ , as
n
defined in Eq. (35), and the results for the variation of Q
2
depending on n are displayed in Fig. 7. From this figure,
one can notice that the value of Q decreases with n, that
2 7 To compare with the results of Milgrom (2009) for the
iswithafastertransitionfromtheweak-fieldregime(where
quadrupole, one should note the different conventions in use,
µ(y) y) to the strong field regime (µ(y) 1). We have namelyqMilgrom= QBN andqMilgrom(η)= 2qBN(η).Ourre-
been∼unable to study higher values of n and∼to determine sultsfortijhequadru−poleijareingoodagreemen−tw3it2hvaluesgiven
a possible limit for Q2 as n goes to infinity. In Table 1, we byMilgrom(2009).
10 Luc Blanchet & J´erˆome Novak
Table1.NumericalvaluesofthedipoleQ1,thequadrupoleQ2andtheoctupoleQ3togetherwiththeassociateddimensionlessquantities
q1,q2 and q3 defined by Eqs. (33) and (37). All values are given near the Sun. We use different choices of the function µ(y) defined in
Eqs.(34)and(35).Thesevalueshavebeenobtainedfora0=1.2 10−10m s−2 andge=1.9 10−10m s−2 (andM =1M⊙).
× · × ·
MONDfunction µ1(y) µ2(y) µ5(y) µ20(y) µexp(y) µTeVeS(y)
Q1 [ms−2] 8.9 10−14 9.8 10−14 1.1 10−13 1.2 10−13 9.7 10−14 3.5 10−14
q1 −7.4× 10−4 −8.2× 10−4 −9.2× 10−4 − 1×0−3 −8.1× 10−4 −2.9× 10−4
− × − × − × − − × − ×
Q2 [s−2] 3.8 10−26 2.2 10−26 7.4 10−27 2.1 10−27 3.0 10−26 4.1 10−26
q2 ×0.33 ×0.19 6.5× 10−2 1.8× 10−2 ×0.26 ×0.36
× ×
Q3 [m−1s−2] 1.2 10−40 9.3 10−41 4.9 10−42 2.3 10−42 1.2 10−40 1.1 10−40
q3 − ×1.1 − ×0.87 −4.6× 10−2 −2.1× 10−2 − ×1.1 − ×1
− − − × − × − −
4 EFFECT ON THE DYNAMICS OF SOLAR dΩ 1 ∂R
= . (38f)
SYSTEM PLANETS dt a2nsinI√1 e2 ∂I
−
Weshall beinterested only in secular effects,so we av-
4.1 Perturbation equations erage these equations over one orbital period P; denoting
thetime average by bracketswe have
In this section, we investigate the consequence for the dy-
namics of inner planets of theSolar System of thepresence ∂R 1 P ∂R 1 2π ∂R
= dt = dℓ . (39)
of an abnormal multipole moment QL oriented toward the (cid:28)∂cA(cid:29) P Z0 ∂cA 2π Z0 ∂cA
direction e of the galactic centre, in the sense of Eq. (27).
In particular this implies that we shall always have
Recall that the domain of validity of this anomaly is ex-
∂R/∂ℓ = 0 hence the perturbation equation (38a) gives
pected to enclose all the inner Solar System (for distances h i
da/dt =0 (at first order in perturbation theory). Indeed,
r . r0 ≈ 7100 AU), with the quadrupole coefficient being hapertuirbativeforceofthetypeF =∇Risconservative,so
constant up to say 50 AU (see Fig. 4). As we have seen,
the energy of the orbit is conserved in average and there is
theanomalyinducesaperturbationontheNewtoniangrav-
nosecular change in the semi-major axis.
itational potential, namely u=u +δu, where the pertur-
N
Letusnowapplytheperturbationequations(38)tothe
bation function δu is given for various multipole moments
specific case of the perturbation function corresponding to
by Eqs. (29). Following the standard practice of celestial
thequadrupoleanomaly, namely
mechanics we denote the perturbation by R δu. The
≡
perturbation function R is such that the perturbing force 1 1
R=δu = r2Q cos2θ . (40)
(or,rather,acceleration) actingontheNewtonianmotionis 2 2 2(cid:18) − 3(cid:19)
F =∇R.
Here cosθ = e n, with n = x/r being the unit direction
The unperturbed Keplerian orbit of a planet around ·
of the planet and r2 =x2+y2+z2, where (x,y,z) are the
the Sun is described by six orbital elements (say c with
A
coordinatesoftheplanetinanabsoluteGalilean framecen-
A = 1, ,6). For these we adopt the semi-major axis a,
··· tred on the focus of the unperturbedKeplerian ellipse, and
theeccentricity e,theinclination I of theorbital plane, the
with respect to which the orbital elements a,e,I,ℓ,ω,Ω
mean anomaly ℓ defined by ℓ = n(t T) where n = 2π/P { }
− are defined. In the following, to simplify the presentation,
(n is themean motion, P is theorbital period and T is the
we shall choose the x-direction of this Galilean frame to be
instant of passage at the perihelion), the argument of the
thedirectionofthegalacticcentree=g /g .Thatis,weas-
perihelion ω (or angular distance from ascending node to e e
sumethattheorigin of thelongitudeof theascendingnode
perihelion), and the longitude of theascending node Ω. We
Ωliesinthedirectionofthegalacticcentre.Thismeansthat
alsousethelongitudeoftheperiheliondefinedbyω˜ =ω+Ω.
cosθ=x/r so that
TheperturbationisafunctionR(c )oftheorbitalele-
A
mentsoftheunperturbedKeplerianellipse.Withourchoice R= Q2 2x2 y2 z2 . (41)
for the orbital elements {cA} = {a,e,I,ℓ,ω,Ω} the per- 6 (cid:16) − − (cid:17)
turbation equations of celestial mechanics read as (see e.g.
We express the planet’s absolute coordinates (x,y,z) in
Brouwer & Clemence 1961)
terms of the orbital elements a,e,I,ℓ,ω,Ω by perform-
{ }
ingas usualthreesuccessive frame rotations with angles Ω,
da 2 ∂R
= , (38a) I andω,toarriveat theframe(u,v,w)associated with the
dt an ∂ℓ
motion, where (u,v) is in the orbital plane, with u in the
de = √1−e2 1 e2∂R ∂R , (38b) direction of the perihelion and v oriented in the sense of
dt ea2n (cid:20)p − ∂ℓ − ∂ω(cid:21) motion at perihelion. The unperturbed coordinates of the
dI 1 ∂R ∂R planet in this frame are
= cosI , (38c)
dt a2n√1 e2sinI (cid:20) ∂ω − ∂Ω(cid:21)
u = a(cosU e) , (42a)
−
dℓ 1 ∂R 1 e2∂R −
= n 2a + − , (38d) v = a 1 e2sinU, (42b)
dt − a2n(cid:20) ∂a e ∂e(cid:21) −
p
w = 0, (42c)
dω 1 √1 e2∂R cosI ∂R
= − , (38e)
dt a2n(cid:20) e ∂e − sinI√1 e2 ∂I (cid:21) where U denotes the eccentric anomaly, related to ℓ by the
−