Théorie de Darwin-de Sitter

Théorie des figures d'équilibre

La théorie des figures d'équilibre considérée ici résulte de nombreuses études traitant du problème de la forme d'équilibre de la Terre, en supposant que celle-ci soit causée par la seule force de pesanteur, à l'exclusion de forces de cohésion internes ou de forces électriques et magnétiques. La théorie présentée repose entièrement sur la notion de surface de niveau, et requiert donc l'existence d'un potentiel de pesanteur à l'intérieur du corps que l'on étudie. Cela exclut la présence d'une rotation différentielle sauf si celle-ci correspond à une rotation différentielle en couches cylindriques coaxiales.
Bien que ce problème soit vieux de plus de trois siècles (cf Modèle ellipsoïdal de la Terre), il n'en continue pas moins à être un problème actuel qui possède des applications intéressantes en géodésie et en géophysique. Il est fermement enraciné dans la théorie générale des figures de corps cosmiques en rotation, qui constitue elle-même un domaine de recherche active en astrophysique et en physique planétaire.^[1] Toutefois, pour des corps cosmiques en rotation lente — comme c'est le cas de la Terre — et pour autant que nous puissions admettre qu'il y a équilibre hydrostatique, la théorie des figures semble avoir atteint un niveau de finition acceptable.^[2] L'étude de la figure hydrostatique de la Terre n'a pas seulement exercée une profonde influence sur la naissance et le développement subséquent de la géodésie physique et de la géophysique, mais elle a aussi grandement contribué aux fondations de la physique de Newton, de l'hydrostatique, de la mécanique analytique, de l'analyse mathématique et de la physique mathématique pendant les XVIII^e et XIX^e siècles, et se trouve associée à la plupart des noms des mathématiciens et physiciens les plus distingués de cette époque.^[3]

Sommaire

1 Équations fondamentales
2 Relation constitutive
3 Rotateur permanent
4 Symétrie des figures d'équilibre
5 Théorie du champ interne de Clairaut, Laplace et Lyapunov
6 Fondement de la théorie des figures
7 Notes
8 Voir aussi ...

Équations fondamentales

Pour étudier la figure d'équilibre hydrostatique de la Terre, on se place dans le cadre général de la mécanique des milieux continus. Les équations de base d'un continuum sont, en formulation eulérienne :

l'équation de continuité, qui exprime la conservation de la masse :

$\frac{\partial \rho}{\partial t} + \sum_{k=1}^{3}\frac{\partial (\rho v_k)}{\partial x_k} = 0$ ,
l'équation de mouvement, qui exprime la conservation de la quantité de mouvement, ou impusion :

$\rho \frac{\partial v_i}{\partial t} + \rho \sum_{k=1}^{3} v_k \frac{\partial v_i}{\partial x_k} = f_i + \sum_{k=1}^{3} \frac{\partial T_{ik}}{\partial x_k}$ .

Dans ces équations, ρ(x₁, x₂, x₃) désigne la densité en un point spatial P repéré par ses coordonnées cartésiennes (x₁, x₂, x₃)), (v₁, v₂, v₃) le champ de vitesses au même point, (f₁, f₂, f₃) le champ de force volumique et T_ik, avec i ,k = 1, 2, 3, le tenseur des contraintes en P.

Dans les circonstances habituelles, que nous supposerons réalisées ici

l'équation de rotation qui exprime la conservation de la quantité de rotation, ou moment cinétique, implique que le tenseur des tensions est symétrique

$T i k = T k i$ ,

et se réduit donc à six composantes indépendantes au lieu de neuf : T₁₁, T₂₂, T₃₃, T₁₂ = T₂₁, T₂₃ = T₃₂, T₃₁ = T₁₃.

Relation constitutive

Afin de rendre le problème d'une figure d'équilibre traitable, il faut se donner une loi rhéologique (ou relation constitutive) permettant d'exprimer le champ des tensions en termes de champs dérivés du champ des vitesses ou des déplacements. En général, on suppose que le comportement rhéologique à long-terme du milieu continu étudié est celui d'un fluide parfait. Alors, le tenseur des tensions T_ij est isotrope et peut se décrire au moyen d'un champ de pression scalaire p, soit

T i j = - p δ i j

où δ_ij représente le symbole de substitution de Kronecker : δ_ij vaut 1 si i et j sont égaux, et 0 sinon.

Rotateur permanent

Définition d'un rotateur permanent

On peut simplifier le problème encore davantage en supposant que la figure de'équilibre constitue un rotateur permanent. On définit un corps en rotation permanente comme suit :^[4]

il est isolé dans l'espace et tourne autour d'un axe fixe ;
son mouvement est stationnaire dans un système d'inertie, et la densité de chaque point massique individuel reste constante au cours du mouvement ;
on peut négliger les effets dissipatifs, notamment les frictions ;
aucune force électromagnétique n'agit sur le corps.

Champ de vitesse d'un rotateur permanent

Constatons que l'on n'a pas encore spécifié la vitesse angulaire ω. Admettons maintenant que le corps tourne autour de l'axe (O x₃) et prenons l'origine O du système d'axes au centre de masse du corps. Considérons des coordonnées cylindriques (h, λ, z) définies en termes des coordonnées cartésiennes (x₁, x₂, x₃) par la loi de transformation :

$\begin{cases} h = (x_1^2 + x_2^2)^{\frac{1}{2}}\\ \lambda = \mathrm{arctan}(\frac{x_2}{x_1})\\ z = x_3\end{cases}$ .

Ainsi, dans un système d'inertie, les composantes de la vitesse prennent la forme

$\begin{cases} v_h = 0\\ v_\lambda = \omega h\\ v_z = 0\end{cases}$ .

Équation de continuité pour un rotateur permanent

Comme la densité d'un point massique appartenant à un rotateur permanent ne varie pas le long de sa trajectoire, nous avons

$\frac{\mathrm{d}\rho}{\mathrm{d}t} = \frac{\partial \rho}{\partial t} + v_j \frac{\partial \rho}{\partial x_j} = 0$ .

L'équation de continuité implique donc que le mouvement est indivergentiel^[5], c'est-à-dire

$\frac{\partial v_j}{\partial x_j} = 0$ .

Symétrie du champ de vitesse d'un rotateur permanent

En coordonnées cylindriques, cette équation devient^[6]

$\frac{\partial \omega}{\partial \lambda} = 0$ ,

ce qui implique qu'une rotation permanente est caractérisée par un champ de vitesse qui ne dépend pas de la variable λ^[7], soit

ω = ω(h, z)

Un cas particulier d'un rotateur permanent est un rotateur rigide, pour lequel la vitesse angulaire est constante :

ω = Constant

Le corps concerné tourne alors en bloc sur lui-même : on dit que la figure d'équilibre correspondante est en équilibre relatif.

Or, l'hypothèse d'une rotation permanente ne requiert pas seulement que dρ/dt = 0, mais aussi que ∂ρ/∂t = 0, d'où l'on tire au moyen de l'équation de continuité que

$v_j \frac{\partial \rho}{\partial x_j} = 0$ .

Cette dernière relation montre que pour une rotation permanente les lignes de courant sont orthogonales au gradient de densité local. Dans le système de coordonnées cylindriques que nous avons adopté, cette relation s'écrit

$\omega \frac{\partial \rho}{\partial \lambda} = 0$ .

Symétrie de la distribution de masse d'un rotateur permanent

En supposant que le corps tourne effectivement (ω ≠ 0) comme un rotateur permanent par rapport à un référentiel d'inertie, nous trouvons donc que

ρ = ρ(h, z)

On en déduit que la distribution de densité dans un rotateur permanent possède une symétrie axiale, et la même conclusion s'impose (voir ci-dessus) pour le champ de vitesse angulaire. Ainsi, ce champ ω = ω(h, z) ne possède pas de tourbillon^[8] dans le plan (O x₁ x₂), donc^[9]

$x_2 \frac{\partial \omega}{\partial x_1} - x_1 \frac{\partial \omega}{\partial x_2} = 0$ .

Équation de mouvement d'un rotateur permanent

En termes de coordonnées cartésiennes, le champ de vitesse peut donc s'écrire

$\begin{cases} v_1 = - \omega x_2\\ v_2 = \omega x_2\\ v_3 = 0\end{cases}$ .

En vertu des relations ∂v_j/∂x_j = 0 et v_j ∂ρ/∂x_j = 0,, la quantité ∂(ρ v_i v_j) / ∂_j devient simplement ρ v_j ∂v_i / ∂x_j. Dès lors, l'équation de mouvement prend la forme particulière

$- \rho \omega^2 (x_i - \delta_{i3} x_3) = f_i - \frac{\partial p}{\partial x_i}$ .

Potentiel gravifique d'un rotateur permanent

Selon la définition d'un rotateur permanent énoncée plus haut, il est évident que le champ de force volumique f_i est simplement le champ de force engendré par l'autoattraction gravitationnelle, c'est-à-dire

f i = ρ g i

Il est bien connu que la force gravifique spécifique g_i en un point P est fournie par l'expression

$g_i (P) = G \int_{B} \frac{\rho(Q)}{d^2(P, Q)} e_i(P, Q) \;\mathrm{d}\tau(Q)$ ,

où d(P, Q) représente la distance entre le point attiré P et un quelconque point attirant Q du corps B. La quantité e_i désigne la i^e composante d'un vecteur unitaire porté par le segment de droite qui joint P à Q, et orienté de P à Q. L'élément de volume en Q est dénoté par dτ(Q). On sait de même que la force volumique g_i dérive en chaque point de l'espace d'une fonction scalaire V, appelé potentiel gravitationnel. Ainsi,

$g_i = - \frac{\partial V}{\partial x_i}$ .

Nous définissons ici ce potentiel V, en accord avec la convention généralement adoptée par les physiciens, comme l'énergie potentielle par unité de masse d'un élément matériel situé au point P. Cette énergie potentielle résulte de l'interaction gravitationnelle de l'élément matériel en P avec tous les autres éléments matériels composant le corps B, de sorte que

$V(P) = - G \int_{B} \frac{\rho(Q)}{d(P, Q)} \;\mathrm{d}\tau(Q)$ .

Selon cette définition, le potentiel gravifique est toujours négatif, et correspond physiquement à une énergie de liaison. Il faut remarquer que les astronomes et les géodésiens, et aussi de nombreux géophysiciens, adoptent une définition différente pour le potentiel gravifique, à savoir ils l'assimilent au travail qu'un élément matériel de masse unitaire accomplit dans le champ de force gravifique g_i en passant d'un point situé à l'infini au point P. Ce travail représente l'énergie potentielle par unité de masse changée de signe, et est donc positif. La convention des physiciens adoptée ici possède l'avantage, en retenant la correspondance directe entre potentiel et énergie potentielle, de faire en sorte qu'un point massique tend toujours à mouvoir vers des régions où le potentiel est plus bas. La constante G ci-dessus est la constante de gravitation.

Symétrie des figures d'équilibre

Lorsqu'on étudie la théorie des figures pour des corps cosmiques comme la Terre ou les planètes, c'est-à-dire pour des corps hétérogènes en rotation autour d'un axe polaire, la principale difficulté provient du fait qu'on ne peut pas admettre a priori une forme déterminée pour une surface de niveau, comme dans le cas de la théorie des figures d'équilibre des masses homogènes (par exemple, le sphéroïde de Maclaurin ou l'ellipsoïde de Jacobi) considérés en astronomie ou des ellipsoïdes de référence discutés en géodésie physique. Dès lors, l'examen des éléments de symétrie d'une figure d'équilibre permet de simplifier appréciablement la théorie. Dans la suite de cet article, nous allons considérer le cas de la Terre, mais la théorie générale s'applique aussi aux autres corps planétaires ainsi qu'aux étoiles simples.

Tout d'abord, nous supposons que la Terre est en rotation permanente, excluant ainsi de notre étude les phénomènes de marée ou des phénomènes électromagnétiques, par exemple. On a vu plus haut que sous ces conditions, le champ de vitesse angulaire et le champ de densité possèdent tous les deux une symétrie de révolution. Il s'ensuit que le potentiel gravifique d'un rotateur permanent possède lui-aussi la même symétrie axiale. Cette propriété de symétrie de $V(P) = - G \int_B \rho(Q) / d(P, Q)\; {\rm d} \tau(Q)$ peut être établie mathématiquement comme suit : soient (h, λ, z) et (h', λ', z') les coordonnées cylindriques du point potentié P et d'un point potentiant quelconque Q, respectivement. Les coordonnées cartésiennes correspondantes de P et Q, c'est-à-dire (x₁, x₂, x₃) et (x'₁, x'₂, x'₃), s'expriment au moyen des relations

$x_1 = h \cos \lambda, \qquad x_2 = h \sin \lambda, \qquad x_3 = z$ ,

$x'_1 = h' \cos \lambda', \qquad x'_2 = h' \sin \lambda', \qquad x'_3 = z'$ .

L'élément de volume en Q est ${\rm d} \tau(Q) = {\rm d} x'_1\;{\rm d} x'_2\;{\rm d} x'_3 = h'\;{\rm d} h'\; {\rm d}z'\;{\rm d} \lambda'$ , et la distance entre P et Q est

$d(P, Q) = \sqrt{h^2 + h'^2 - 2 h h' \cos(\lambda - \lambda') + (z - z')^2}$ .

Comme ρ(Q) ne dépend pas de l'angle λ', l'intégration sur λ' dans l'expression de V(P) peut être remplacée par une intégration sur la différence λ' – λ avec les mêmes limites d'intégration 0 et 2π que pour λ'. De la sorte, la variable λ disparaît de l'intégrale, et

V = V (h, z)

ensemble avec

ρ = ρ(h, z)

ω = ω(h, z)

Principe de symétrie de Pierre Curie

En fait, la symétrie axiale de V aurait pu être déduite de la symétrie axiale de ρ en se servant du principe de symétrie établi par le physicien français Pierre Curie (1859–1906). Ce principe peut s'énoncer de façon concise comme suit : Tout effet produit par une cause déterminée doit posséder au moins les mêmes éléments de symétrie mais peut être plus symétrique que la cause. Toute dissymétrie rencontrée dans un effet donné doit aussi exister dans la cause qui a produit l'effet. En corollaire, s'il existe une relation de réciprocité entre une cause et un effet, dans le sens que la cause peut tout aussi bien être considérée comme l'effet, et l'effet comme la cause, les deux — cause et effet — doivent posséder exactement les mêmes éléments de symétrie.^[10] Un exemple de ce cas est fourni par une relation entre les tensions et les déformations, à savoir ce que nous avons appelé plus haut une « relation constitutive ». Il semble que l'on n'ait pas encore trouvé jusqu'à présent d'expérience ni d'observation se rapportant à des phénomènes macrophysiques qui contredise le principe de Curie, et ce dernier devrait être considéré comme une loi universelle de la Nature au même titre que les lois de la thermodynamique, par exemple.^[11] Nous allons adopter ce point de vue ici pour déduire certaines propriétés des figures d'équilibre.

Tout d'abord, nous notons que la force d'attraction gravifique est centrale, et possède donc le plus haut degré de symétrie, celui de la sphère. Le principe de Curie requiert alors que l'arrangement des masses produit par la seule action de l'attraction gravifique doit correspondre à la symétrie sphérique : une quelconque configuration moins symétrique est exclue si l'on veut avoir équilibre. Une preuve mathématique du fait que la sphère est l'unique figure d'équilibre possible d'un corps fluide isolé au repos fut donnée par Lyapunov^[12] en 1884.

Ensuite, considérons l'action combinée d'une force d'attraction gravifique et d'une force axifuge. Cette dernière introduit une certaine asymétrie par rapport à la sphère par le fait que la force axifuge possède la symétrie d'un cylindre de révolution. Le principe de Curie requiert donc que la distribution d'équilibre des masses possède la même symétrie axiale, puisqu'il existe une relation de réciprocité entre la distribution de densité et le potentiel. C'est le cas étudié plus haut en termes mathématiques.

Équation d'état : barotropes et baroclines

Un point particulier doit retenir toute notre attention à ce stade. En effet, afin de compléter la formulation du problème, nous devons inclure une équation d'état appropriée. Sous forme générale, celle-ci peut s'écrire

$p = p(\rho, T, \varphi_1, ... \varphi_N)$ .

La lettre p désigne la pression, ρ est la densité, T est la température, et l'ensemble des variables φ₁ à φ_N décrit la composition chimique et minéralogique. Par définition, une structure baroclinique (ou, simplement, un barocline) est un système pour lequel une relation physique de cette forme existe. Toutefois, dans certaines circonstances spéciales qui sont approximativement remplies dans le cas de l'intérieur de la Terre,^[13] la pression peut simplement être reliée à la densité, à savoir

p = p (ρ)

Si tel est le cas, on dira qu'il s'agit d'une équation d'état barotropique et la structure correspondante est un barotrope. En première approximation, nous pouvons considérer la Terre comme un barotrope.

La distinction fondamentale entre modèles barotropiques et modèles barocliniques se trouve dans leurs stratifications respectives. Il est clair que les surfaces d'égale densité (surfaces isopycniques) et les surfaces d'égale pression (surfaces isobares) coïncident dans un barotrope. Par contre, dans les baroclines les surfaces d'égale densité coupent les surfaces d'égale pression, à moins que des conditions très spécifiques ne soient remplies.

Tenant compte de la symétrie axiale de ρ et V, l'équation de mouvement d'un rotateur permanent devient en coordonnées cylindriques

$\frac{1}{\rho}\frac{\partial \rho}{\partial h} = - \frac{\partial V}{\partial h} + \omega^2 h$ ,

$\frac{1}{\rho} \frac{\partial p}{\partial z} = - \frac{\partial V}{\partial z}$ ,

$\frac{1}{\rho h} \frac{\partial p}{\partial \lambda} = 0$ .

Le fait que ∂p/∂λ = 0 établit la symétrie de révolution du champ de pression, c'est-à-dire

p = p (h, z)

Équilibre hydrostatique relatif

On définit les composantes de la pesanteur en coordonnées cylindriques par les relations

$g_h = - \frac{\partial V}{\partial h} + h \omega^2$ ,

$g_z = - \frac{\partial V}{\partial z}$ ,

g λ = 0

Dès lors, dans un corps en rotation permanente, la pesanteur est partout orthogonale aux isobares, comme on le constate en se souvenant de la définition du gradient. Constatons que cette propriété est très générale : elle est valable tout aussi bien pour les baroclines que pour les barotropes, et qu'elle est valable pour une rotation uniforme comme pour une rotation différentielle, pour autant que le champ de vitesse angulaire conserve une symétrie axiale. En général, on définit le potentiel axifuge comme

$Z = -\frac{1}{2} \omega^2 h^2$ ,

mais il est évident que cette définition n'est valable que pour une rotation uniforme. Si tel est le cas, nous pouvons employer un système d'axes fixé dans le corps et tournant avec lui. Si chaque particule du corps est au repos par rapport à ce référentiel qui n'est pas d'inertie, on dit que le corps est en équilibre hydrostatique relatif ou, plus simplement, en équilibre hydrostatique. Le plus souvent, dans les théories géodésiques ou géophysiques, on admet que tel est le cas pour la Terre.

Au contraire, si la vitesse angulaire n'est pas constante d'un point à l'autre, il ne peut exister qu'un équilibre hydrodynamique. Une rotation différentielle de ce type joue probablement un rôle important dans les étoiles et les planètes géantes, mais semble exclue dans les planètes telluriques à cause d'un trop fort couplage visqueux ou électromagnétique entre les différentes strates. Les études des figures d'équilibre générales donnent souvent le plus d'importance à la rotation zonale, que l'on nomme circulation méridienne ou courants de von Zeipel, ou courants d'Eddington, ou encore courants de Jardetzky.^[14] Toutefois, l'étude des courants de convection ne fait généralement pas partie de la théorie des figures d'équilibre.

Théorème de Wavre-Poincaré

Les relations ci-dessus définissant les composantes de la pesanteur montrent qu'un potentiel axifuge Z, et donc une surface de niveau

U = V + Z = Constante

peut exister si et seulement si ω ne dépend pas de z, c'est-à-dire si la vitesse angulaire est constante sur des surfaces cylindriques centrées sur l'axe de rotation, autrement dit si ω = ω(h). Dans ce cas plus général que la rotation uniforme, on définit le potentiel axifuge par

$Z = - \int_0^h \omega^2(h') h'\;{\rm d} h'$ .

Nous constatons, en particulier, que pour une rotation uniforme, cette relation se réduit bien à la relation précédente comme il se doit.

Une conclusion importante qu'on peut tirer de l'existence d'un potentiel de rotation, et donc de l'existence d'un potentiel de pesanteur, est

$U(h,z) = V(h,z) - \int_0^h\omega^2(h') h'\;{\rm d} h'$ ,

de sorte que

$g_h = - \frac{\partial U}{\partial h}$ ,

$g_z = - \frac{\partial U}{\partial z}$ ,

g λ = 0

On donne à cette conclusion souvent le nom de théorème de Wavre-Poincaré^[15], et on peut l'énoncer sous la forme suivante : Pour un rotateur permanent, l'une quelconque des propriétés qui suivent entraînent automatiquement les trois autres, à savoir

la vitesse angulaire est constante sur des cylindres centrés sur l'axe de rotation ;
la pesanteur dérive d'un potentiel ;
la pesanteur est normal aux surfaces d'égale densité ;
les surfaces d'égale densité coïncident avec les surfaces d'égale pression.^[16]. En effet, selon les équations pour le gradient de la pression et pour les composantes de la pesanteur fournies plus haut, on peut écrire pour la différentielle totale

${\rm d} p = \frac{\partial p}{\partial h} {\rm d} h + \frac{\partial p}{\partial z} {\rm d} z + \frac{\partial p}{\partial \lambda} {\rm d} \lambda$

l'équation suivante :

${\rm d} p = \rho \left(g_h {\rm d} h + g_z {\rm d} z \right)$ ,

donc

${\rm d} U = - \frac{{\rm d} p}{\rho}$ .

Un déplacement arbitraire sur une surface de niveau fournit dU = 0, impliquant donc que dp = 0. Les surfaces isobares coïncident donc avec les surfaces de niveau, et p = p(U) ou U = U(p). Nous trouvons ainsi à cause de la relation sur dU que la densité est aussi constante sur surface de niveau. En conclusion, les isobares, les isopycniques et les surfaces de niveau coïncident toutes s'il existe un potentiel de pesanteur, et la force par unité de masse due à la pesanteur est normale à ces surfaces. En particulier, elle est normale aux isopycniques. Réciproquement, considérons un corps pour lequel les isobares et les isopycniques coïncident. On peut alors définir une fonction U(p), en termes de laquelle la relation donnant la différentielle de la pression devient

d U = - g h d h - g z d z

montrant que dU est une différentielle totale exacte. Par conséquent, les relations g_h(h,z) = – ∂U/∂h, g_z(h,z) = – ∂U/∂z, g_λ(h,z) = 0 doivent s'appliquer, et le champ de pesanteur dérive d'un potentiel. Finalement, supposons que la pesanteur est partout normale aux surfaces d'égale densité, soit gx∇ρ = 0. Comme g = – ∇U = ρ^–1 ∇p, la coïncidence des surfaces isobares et isopycniques est établie, et tous les énoncés du théorème de Wavre-Poincaré sont démontrés.

Baroclines et pseudo-barotropes

Il est important de se rendre compte que la condition ω = ω(h), ou la condition équivalente ∂ω/∂z = 0, ne requiert aucune connaissance préalable de l'équation d'état. Selon le théorème de Wavre-Poincaré, les isopycniques d'un barocline coïncideront donc avec ses isobares si ∂ω/∂z = 0. Toutefois, il est évident que de tels baroclines particuliers devraient être appelés de façon plus appropriée des pseudo-barotropes, parce qu'ils partagent la plupart de leurs propriétés avec les barotropes véritables, et se distinguent clairement des baroclines véritables.

Cette constatation nous permet de revenir aux considérations de symétrie. Tout d'abord, le principe de Curie nous apprend que les barotropes sont au moins aussi symétriques que les baroclines véritables, mais peuvent être plus symétriques.

Théorème de Lichtenstein

Le principe de Curie (sans symétrie brisée) permet d'inférer une autre propriété de symétrie très importante des figures d'équilibre : un quelconque rotateur permanent barotropique ou pseudo-barotropique doit posséder un élément de symétrie supplémentaire, à savoir un plan équatorial, c'est-à-dire un plan normal à l'axe de rotation qui passe par le centre de masse. Cet énoncé constitue le théorème de Lichtenstein.

Mais cette preuve physique du théorème de Lichtenstein dans toute sa généralité exigerait l'unicité de solution, et la démonstration mathématique rigoureuse n'est pas triviale du tout.^[17]

Le cas d'un barocline véritable est quelque peu plus compliqué, mais Dive a démontré en 1930 qu'un barocline en rotation différentielle est symétrique par rapport à un plan équatorial si la vitesse angulaire est partout une fonction biunivoque de la densité et de la distance à l'axe de rotation.^[18]

Théorie du champ interne de Clairaut, Laplace et Lyapunov

Considérons la Terre comme un barotrope tournant à une vitesse angulaire uniforme ω autour de l'axe Ox₃. Nous désirons trouver sa forme d'équilibre pour une distribution de densité ρ spécifiée le long d'une direction inclinée d'un certain angle par rapport à l'axe de rotation. Sous le vocable « figure d'équilibre » nous ne comprenons pas seulement la forme de la surface externe et du champ de pesanteur externe, mais aussi la forme d'une strate interne quelconque, c'est-à-dire la forme de toutes les surfaces isopycniques, isobares et équipotentielles qui coïncident toutes dans le cas étudié ici. Nous prendrons l'origine O du système de référence au centre de masse de la Terre, et nous dénoterons par r le rayon vecteur qui joint l'origine à un point arbitraire situé sur une surface de niveau donnée

U (r,θ) = V (r,θ) + Z (r,θ) = Constant

Comme il est coutumier, la lettre grecque θ (thêta) désigne la colatitude géocentrique, et r = |r| est la distance au centre de masse. En termes de coordonnées sphériques r, θ, λ le potentiel axifuge s'écrit

$Z(r, \theta) = - \frac{1}{2} r^2 \omega^2 \sin^2\theta$ ,

ou encore^[19],

$Z(r, \theta) - \frac{1}{3} r^2 \omega^2 \left[1 - P_2 (\cos \theta) \right]$ ,

P₂ désignant le polynôme de Legendre d'indice 2.

Les fonctions de figure

Pour un barotrope quelconque en rotation lente, et en particulier pour les modèles de Terre, les surfaces de niveau ne diffèrent que légèrement de la forme sphérique. Il est donc naturel, compte tenu des éléments de symétrie que doivent posséder ces surfaces, de rechercher l'équation décrivant la surface d'une couche quelconque sous la forme

$r(s, \theta) = s \left[1 + \sum_{n = 0}^\infty s_{2n} (s) P_{2n} (cos\theta) \right]$ .

La variable s désigne ici un paramètre judicieusement choisi devant permettre de reconnaître de quelle couche il s'agit — pour simplifier, on dira que s représente le nom de la couche considérée. Ainsi, on peut choisir s comme rayon équatorial a, ou comme rayon polaire c, de la couche en question. Dans la suite, nous choisissons s comme rayon équivolumétrique moyen de la couche, autrement dit, comme rayon de la sphère qui possède le même volume que celui contenu sous la couche.

De la sorte, le problème de trouver la forme d'équilibre est ramené au problème de déterminer un ensemble infini de fonctions s_2n qu'on appelle les « fonctions de figure ». Afin de pouvoir utiliser la formule précédente en pratique, il est nécessaire que l'ordre de grandeur des fonctions de figure décroisse suffisamment vite avec l'ordre n pour assurer une convergence raisonnablement rapide de la série figurant dans le membre de droite. Pour les modèles de Terre en particulier, tel est bien le cas comme il sera montré plus loin.

Si nous tronquons la série à n = 1, nous sommes en présence de l'approximation du premier ordre qui correspond au cas classique de la théorie de Clairaut. Si nous tronquons la série à n = 2, nous obtenons l'approximation du second ordre qui correspond à la théorie de Darwin–de Sitter. Pour la Terre, l'approximation d'ordre 2 correspond à négliger des effets d'environ 10^-5 ; comme le rayon moyen de la Terre est de 6 371 km, on commet donc des erreurs qui peuvent se chiffrer typiquement à quelques mètres. De telles erreurs sont incompatibles avec la précision des mesures géodésiques atteinte à l'heure actuelle, qui est fournissent le rayon de la Terre à quelques centimètres près. Depuis la fin des années 1970, divers auteurs^[20] ont considéré pour le cas de la figure hydrostatique de la Terre des approximations d'ordre 3 et en ont fourni les formules appropriées. En réalité, pour autant qu'on ne désire pas faire des calculs effectifs, il n'est guère plus difficile de développer la théorie des figures d'équilibre à un ordre général n = N que de se limiter à l'ordre n = 2 ou n = 3, par exemple.

Théorie de Clairaut-Radau

Si nous négligeons les termes en $\mathcal{O}(m^2)$ ou plus petits, les équations de figure se réduisent à la seule équation

$F_2 = - s_2 S_0 + S_2 + T_2 - \frac{1}{3} m = 0$ ,

et l'aplatissement devient

$f = - \frac{3}{2} s_2$ .

Combinant les deux, on obtient l'équation intégro-différentielle de Clairaut

$f S_0 + \frac{3}{2} (S_2 + T_2) - \frac{1}{2} m = 0$ .

Les quantités S₂ et T₂ ne doivent être évaluées qu'à l'ordre $\mathcal{O}(m)$ , et S₀ à l'ordre $\mathcal{O}(m)$ . Dès lors, avec $x = s / s 1$ et $\delta = \rho / \bar \rho$ , on trouve

$S_2(x) = - \frac{2}{5} \frac{1}{x^5} \int_{y=0}^x \delta(y) {\rm d}(y^5 f)$ ,

$T_2 (x) = - \frac{2}{5} \int_{y=x}^1 \delta(y) {\rm d} f$ .

Équation différentielle de Clairaut

Sous cette forme, l'équation de Clairaut est une équation intégro-différentielle linéaire pour l'aplatissement des strates internes de la Terre. Elle peut être résolue itérativement par la méthode générale indiquée plus haut pour les fonctions de figure. Toutefois, l'approche classique consiste à dériver l'équation de de Clairaut sous la forme intégro-différentielle par rapport à x et éliminer les quantités m et T₂ en faisant usage des deux dernières relations, obtenant ainsi

$\frac{5}{2} S_2 = \left(f - \frac{1}{3} x \frac{{\rm d} f}{{\rm d} x} \right) S_0$ ,

puis à multiplier ce résultat intermédiaire par x⁵ et dériver encore une fois par rapport à x. On obtient alors l'équation différentielle de Clairaut

$\frac{{\rm d}^2 f}{{\rm d} x^2} + 6 \frac{\delta}{x \; S_0} \frac{{\rm d} f}{{\rm d} x} + \frac{6}{x^2} \left( \frac{\delta}{S_0} - 1\right) f = 0$ .

On peut intégrer numériquement cette dernière équation différentielle ordinaire en imposant que f reste fini en x = 0 (au centre) et que f satisfasse en x = 1 (à la surface) à la condition

$\left. \frac{{\rm d} f}{{\rm d} x} \right|_1 + 2 f(1) = \frac{5}{2} m$ .

On aboutit à cette condition en multipliant l'équation intégro-différentielle de Clairaut ci-dessus par x⁵, en posant ensuite x = 1 et en notant que S₀(1) = 1.

Équation différentielle de Clairaut-Radau

Tout au long du XIX^e siècle, en l'absence des données fournies plus tard par la sismologie, les astronomes essayaient de déterminer la densité des couches internes de la Terre en imposant deux contraintes à diverses fonctions d'essai ρ = ρ(x) : (1) la fonction d'essai, par intégration sur le volume, devait fournir la valeur observée M de la masse de la Terre ; (2) en intégrant l'équation de Clairaut, la fonction d'essai devait fournir la valeur observée f de l'aplatissement géométrique en surface. Or, avant l'avènement des ordinateurs et, surtout, des microordinateurs personnels dans la deuxième moitié du XX^e siècle, la solution de l'équation différentielle de Clairaut s'obtenait par une intégration numérique qui devait se faire à la main et qui était fort fastidieuse.

Pour cette raison, l'élégante méthode de solution proposée en 1885^[21] par l'astronome français d'origine prussienne Rodolphe Radau (1835–1911), qui simplifie énormément l'intégration de l'équation de Clairaut, devint la méthode de résolution standard et reste jusqu'à présent d'une grande valeur dans les travaux théoriques ayant recours à l'aplatissement calculé au premier ordre.^[22]

La méthode de Radau consiste à transformer l'équation de Clairaut en introduisant la fonction sans dimension dite « paramètre (ou variable) de Radau »

$\eta(x) = \frac{x}{f} \frac{{\rm d} f}{{\rm d} x}$ ,

au moyen de laquelle on obtient l'équation de Clairaut-Radau

$\frac{{\rm d} x^5 \sqrt{1 + \eta} \;S_0}{{\rm d} x} = 5 x^4 S_0 \;F(\eta)$ ,

avec

$F(\eta) = \frac{1 + \frac{\eta}{2} - \frac{\eta^2}{10}}{\sqrt{ 1 + \eta}}$ .

Équation de Radau