Tras el anterior
GOCE el geoide, para ver la relación entre el geoide y el valor de la gravedad trato un caso sencillito. Considero la Tierra perfectamente esférica, homogénea, salvo por una enorme cavidad subterránea también esférica, próxima a la superficie y sobre el plano ecuatorial, como se indica en la siguiente figura.
Para los cálculos reflejados en las figuras de abajo se ha usado \(r=50.000\ m\), \(D=R-2r\), y como valores de las constantes :
radio terrestre \(R=6.371.000\ m\)
constante de gravitación universal \(G=6\text{,}67 \ 10^{-11} m^3/(kg\ s^2)\)
masa de la Tierra \(M=5\text{,}97 \ 10^{24} kg\)
velocidad angular de la Tierra \(\omega=2 \pi/(23\text{,}9345\times 60\times 60) rad/s\).
La ventaja de semejante modelo es que el potencial gravitacional \(V\) es fácil de calcular: se pone el potencial de la esfera terrestre de radio \(R\) y masa \(M\), y se le
resta el que crearía la cavidad si estuviera llena, con una masa \(M (r/R)^3 \), con lo que en un punto fuera de la esfera terrestre tenemos, siendo \(m_r= {\left(\frac{ r }{R }\right)} ^3\) , es
\[V(P)=GM \left(\frac{ -1 }{d }+ \frac{ m_r }{e } \right)\] La superficie equipotencial que pasa a 200 metros sobre el polo norte, el punto en que \(\theta=90º\) y \(d=R+200\), es la mostrada en color azul en las gráficas de la siguiente figura, donde la ordenada es la altura sobre la superficie terrestre, \(d-R\).
En la gráfica de la derecha se ha considerado el potencial efectivo \(W = V + C\), que incluye el potencial centrífugo \[ C(P) = -\frac{ 1 }{2 }(\omega \cdot d \cdot cos\theta )^2 \] La curva equipotencial de W se muestra en color rojo.
Empezemos excluyendo de momento el efecto centrífugo. Si toda la esfera terrestre estuviese rellena homogéneamente la simetría provocaría que la superficie equipotencial pasando a 200m sobre el Polo Norte fuese una esfera perfecta. Debido a la existencia de ese nada despreciable hueco esférico de 50km de radio, centrado a una profundidad de 100km, la superficie equipotencial azul se hunde hasta casi 200m cerca de la oquedad, cuando \(\theta \approx 0º\), y sobresale muy poquito cuando \(\theta \approx 180º\). Lo que parece claro físicamente es que cuando \(\theta \approx 0º\) también deberá ser menor la gravedad \(g\) debido a la ausencia de masa en el hueco. En este caso por tanto un hundimiento de la superficie equipotencial conlleva disminución de \(g\).
Si ahora incluimos el potencial centrífugo \(C\), su influencia es muy marcada, y provoca una gran elevación de la superficie equipotencial de \(W=V+C\) en el ecuador ( \(\theta\approx 0º \ \text{ó} \ \theta \approx 180º\) ) respecto a los polos ( \(\theta\approx 90º \) ). La intuición física nos hace esperar ahora que esa elevación vaya aparejada con una disminución de \(g\), pues cuanto mayor sea la fuerza centrífuga más se contraresta la atracción gravitacional.
Bueno, en un caso es el hundimiento y en el otro la elevación de la superficie equipotencial lo que se corresponde con la disminución de \(g\). Lo que ilustra esta aparente contradicción es que \(g\) no va ligado a subidas o bajadas de la superficie equipotencial, sino a diferencia de distancias entre superficies equipotenciales contiguas. Para calcular aproximadamente \(g\) se buscan dos superficies equipotenciales próximas, por ejemplo las que pasan a 200m y a 300m sobre el Polo Norte. Los valores de \(g\) sin y con el efecto centrífugo se muestran en la siguiente figura. En adelante se expresarán siempre en \(m/s^2\).
Vemos que \(g\) disminuye tanto por el efecto de falta de masa de la oquedad, como por el efecto centrífugo, que es bastante más acusado. Este es el responsable en \(\theta = 180º\) de una disminución de 0,068 en \(g\) comparado con solamente 0,019 en \(\theta = 0º\) debido a la oquedad. Claro que la disminución centrífuga sale casi del doble de lo que uno esperaría considerando el efecto centrífugo sobre la superficie terrestre, \( \omega^2 R \approx 0,034 \). Menudo susto, un factor dos, las cifras que no cuadran, a ver dónde está el error. Bueno, en este caso no hay un 2 extra que sobre, es que \(g\) se calcula sobre la superficie equipotencial, que en el ecuador se eleva 11060 m, con lo que la atracción gravitacional allí se reduce respecto a la superficie terrestre en
\[
\frac{GM}{(R+200)^2} - \frac{GM}{(R+11060+200)^2} \approx 0\text{,}034
\]
¡Vaya! 0,034. De modo que la disminución de 0,068 desde el polo al ecuador es una contribución equitativa de efecto centrífugo "directo" y alejamiento del centro. Es curioso que ese factor 2 no es casual. Olvidando la oquedad, atendiendo al efecto centrífugo exclusivamente, el valor de \(g\) en \(\theta=90º\) viene dado por \(GM/d^2\), y para \(\theta=0º\) viene dado por \(GM/{d'}^2-\omega^2 d'\). Si se restan ambos valores teniendo en cuenta que estamos
sobre la misma superficie equipotencial, y que por tanto
\[
- \frac{GM}{d} = - \frac{GM}{d'} - \frac{1}{2} \omega^2 {d'}^2
\]
usamos la expresión de \(1/d\) en términos de \(1/d'\) y ...
\[
\frac{GM}{d^2} - \left(\frac{GM}{{d'}^2} - \omega^2 d' \right) = \frac{\omega^4 {d'}^4}{4GM} + \omega^2 d' + \omega^2 d' \approx 2 \omega^2 d'
\]
puesto que \( \omega^4 {d'}^4 / 4GM = 0\text{,}000029\).
La superficie terrestre tiene cierta plasticidad, y la forma abombada en el ecuador se adapta al efecto centrífugo. En el modelo sencillo utilizado la diferencia entre ecuador y polo norte para la superficie equipotencial es de 11.060m. En realidad la Tierra está muy bien descrita por un elipsoide de radio ecuatorial de 6.378.137 m y radio polar de 6.356.752 m, que supone una diferencia de
21.385 m. Ese elipsoide se diferencia del geoide en pocos cientos de metros. Una vez que hay redistribución de masa el propio potencial gravitacional cambia y no se calcula tan sencillamente.
En cuento a los efectos gravitacionales por exceso o defecto de masa, ya se ve que una oquedad tan inmensa e irreal como la supuesta conlleva un efecto relativamente "débil", de 0,019. No es de extrañar que en geofísica las anomalías de la gravedad a lo largo del terreno se midan en mGal, miligales, siendo un
Gal una peculiar designación de \(cm/s^2 = 10^{-2} m/s^2\). Es decir, \(1mGal = 10^{-5}m/s^2\), así de pequeña es la variación local de \(g\).
En el
sitio técnico de la misión GOCE hay una enorme cantidad de recursos, en inglés, como el de
la gravedad en detalle donde nos presentan entre otras unidades el Gal. El apartado de
cantidades del campo de gravedad (Gravity Field Quantities) contiene unas definiciones muy útiles, algunas de las cuales traduzco aquí:
Potencial gravitacional (V): Potencial generado por la atracción de las masas.
Potencial de gravedad (W): Suma del potencial gravitacional (V) y del potencial centrífugo (C) de la Tierra en rotación. Las diferencias entre dos puntos pueden observarse por nivelación.
Superficie equipotencial: Una superficie en la que W es constante. Sus puntos pueden determinarse localmente con mareógrafos, que definen el nivel medio del mar a escala regional.
Geoide: Superficie equipotencial que aproxima el nivel global medio del mar, es decir, una red global de mareógrafos y de señales de nivelación, tras sustracción de los componentes dinámicos. Puede considerarse como un hipotético océano en reposo.
Gravedad: La magnitud, g, del gradiente de W en la superficie terrestre y de V en el espacio. Puede observarse mediante una técnica absoluta (p.ej. experimento de caida libre) o relativa (como una diferencia) mediante un gravímetro de muelle.
He usado esta notación de V, W y C en las anteriores figuras.
En cuanto al código de Octave para sacar semejantes gráficas, primero va esto
#unidades mks
global G = 6.67e-11; #cte gravitacion unidades m^3/(kg s^2)
global MT=5.97e24; #masa de la Tierra en kg
#producto G por MT, unidad de longitud metros
global GMT = G*MT;
global RT = 6371000; #radio terrestre medio en m
global omegaT = (2*pi/(23.9345*60*60));#vel. angular Tierra tomando periodo rot. 23.9345h
function val = Vd(x,m,D,theta)
global GMT;
val = ((-1/x)+(m/sqrt(x^2+D^2-2*x*D*cos(theta))))*GMT;
endfunction
function val = fdz(xp,m,D,theta)
global GMT;
C=Vd(xp,m,D,pi/2);
fun = @(x) ((-1/x)+(m/sqrt(x^2+D^2-2*x*D*cos(theta))))*GMT - C;
[val, fval, info] = fzero(fun,xp+0.000001*rand);
endfunction
function val = fdzcf(xp,m,D,theta)
global GMT;
global omegaT;
C=Vd(xp,m,D,pi/2);
fun = @(x) (((-1/x)+(m/sqrt(x^2+D^2-2*x*D*cos(theta))))*GMT)+ (-0.5*(omegaT*x*cos(theta))^2) - C;
[val, fval, info] = fzero(fun,xp+0.000001*rand);
endfunction
R=RT;
r=50000;
D=RT-2*r;
m=(r/R)^3;
theta=(0:180) *pi/180;
res1=zeros(1,181);
x1=R+200;
C1PN=Vd(x1,m,D,pi/2);
for k=0:180
res1(k+1)=fdz(x1,m,D,theta(k+1));
endfor
res2=zeros(1,181);
x2=x1+100;
C2PN=Vd(x2,m,D,pi/2);
for k=0:180
res2(k+1)=fdz(x2,m,D,theta(k+1));
endfor
cfres1=zeros(1,181);
for k=0:180
cfres1(k+1)=fdzcf(x1,m,D,theta(k+1));
endfor
cfres2=zeros(1,181);
for k=0:180
cfres2(k+1)=fdzcf(x2,m,D,theta(k+1));
endfor
y después
subplot(1,2,1)
plot(theta*180/pi,res1-R,'b',90,200,'xk')
set(text(55,190,"equipotencial a 200m"),"color","black")
set(text(55,180,"sobre el Polo Norte"),"color","black")
set(text(55,170,"V (sin efecto centrífugo)"),"color","blue")
set (gca (), "xlim", [0, 180]);
set (gca (), "ylim", [0, 210]);
xlabel(strcat('\theta (grados)'));
ylabel(strcat('d-R (metros)'));
subplot(1,2,2)
plot(theta*180/pi,res1-R,'b',theta*180/pi,cfres1-R,'r',90,200,'xk')
set(text(50,11000,"equipotencial a 200m"),"color","black")
set(text(50,10500,"sobre el Polo Norte"),"color","black")
set(text(5,500,"V (sin efecto centrífugo)"),"color","blue")
set(text(25,10000,"W = V + C (con efecto centrífugo)"),"color","red")
set (gca (), "xlim", [0, 180]);
set (gca (), "ylim", [0, 11500]);
xlabel(strcat('\theta (grados)'));
y para el cálculo y gráfica de \(g\)
#calcula g para caso sin cf
gg=zeros(1,181);
alfa=zeros(1,181);
for k=2:180
alfa(k)=asin((res1(k+1)-res1(k-1))./(res1(k-1)*2));
gg(k)=(C2PN-C1PN)./((res2(k)-res1(k))*cos(alfa(k)));
endfor
gg(1)=(C2PN-C1PN)./(res2(1)-res1(1));
gg(181)=(C2PN-C1PN)./(res2(181)-res1(181));
#calcula g para caso con cf
ggcf=zeros(1,181);
alfa=zeros(1,181);
for k=2:180
alfa(k)=asin((cfres1(k+1)-cfres1(k-1))./(cfres1(k-1)*2));
ggcf(k)=(C2PN-C1PN)./((cfres2(k)-cfres1(k))*cos(alfa(k)));
endfor
ggcf(1)=(C2PN-C1PN)./(cfres2(1)-cfres1(1));
ggcf(181)=(C2PN-C1PN)./(cfres2(181)-cfres1(181));
subplot(1,1,1)
plot(theta*180/pi,gg,'b',theta*180/pi,ggcf,'r',90,gg(91),'xk')
set(text(120,9.813,"V (sin efecto centrífugo)"),"color","blue")
set(text(120,9.795,"W=V+C (con efecto centrífugo)"),"color","red")
title('gravedad en equipotencial a 200m sobre el Polo Norte');
set (gca (), "xlim", [0, 180]);
set (gca (), "ylim", [9.72, 9.82]);
xlabel(strcat('\theta (grados)'));
ylabel(strcat('g (m/s^2)'));
y ... se acabó. No siento los miligales.
