Esfera iluminada por el Sol

Estado estacionario

Condición de contorno

Una esfera de radio a recibe un flujo de calor y pierde calor debido a la radiación, la condición de contorno en su superficie se escribe

T(r,θ) r | r=a =h(T(a,θ) T s )+ q K

Donde T(a,θ) es la temperatura de la superficie de la esfera, Ts es la temperatura ambiente que la rodea y cuyo valor supondremos que es cero. Explicaremos el significado del término q

Una esfera de radio a se encuentra a una distancia r del centro del Sol. Si r es muy grande, las ondas electromagnéticas procedentes del Sol serán aproximadamente planas, los rayos procedentes del Sol que llegan a media esfera serán paralelos tal como se muestra en la figura.

Denominamos I al vector que nos da la intensidad o energía emitida por el Sol en todas las frecuencias por unidad de tiempo y unidad de área. Lo pondremos, arbitrariamente, paralelo al eje Z. La dirección del vector superficie dS es radial.

El flujo de energía (energía por unidad de tiempo) que se recibe en la porción dS de la superficie de la esfera es el producto escalar I · dS =I·dScosθ .

La energía por unidad de área y unidad de tiempo que se recibe en la posición θ sobre la superficie de la esfera es q=Icosθ

La condición de contorno en la superficie de la esfera r=a, se escribe

T(r,θ) r | r=a =hT(a,θ)+f(θ),f(θ)={ pcosθ,0θ< π 2 0, π 2 θ<π

Las posiciones θ≥π/2 representan la cara no iluminada de la esfera

Solución de la ecuación de Laplace

Solución de la ecuación de Laplace en coordenadas esféricas es

T(r,θ)= n=0 ( A n r n + B n r n+1 ) P n (cosθ)

Vamos a calcular la temperatura en un punto interior de la esfera (r<a, θ). Dado que T(r, θ) no puede ser infinito cuando r→0, los coeficientes Bn=0 tienen anularse

Probamos la solución

T(r,θ)= n=0 A n ( r a ) n P n (cosθ)

Sustituimos en la condición de contorno

n=0 A n n a ( a a ) n1 P n (cosθ) =h( n=0 A n ( a a ) n P n (cosθ) )+f(θ)

Desarrollo en serie en términos de los polinomios de Legendre

Expresamos una función f(x) como combinación lineal de los polinomios de Legendre

f(x)= n=0 c n P n (x) c n = 2n+1 2 1 1 f(x) P n (x)dx ,n=0,1,2...

Como f(θ) tiene valores no nulos en el intervalo (0, π/2), f(x)=p·x, con x=cosθ, está definida en el intervalo (0,1)

c n = 2n+1 2 0 π/2 pcosθ· P n (cos θ)sinθ·dθ=p 2n+1 2 0 1 x P n (x)dx ,n=0,1,2...

Calculamos los seis primeros coeficientes

>> syms x;
>> for n=0:5
int(legendreP(n,x)*x, 0,1)*(2*n+1)/2
end
ans =1/4
ans =1/2
ans =5/16
ans =0
ans =-3/32
ans =0

c0=p/4, c1=p/2, c2=5p/6, c3=0, c4=-3p/32, c5=0, ...

Representamos la función f(θ)=cosθ, en el intervalo (0, π/2) y los cinco primeros términos del desarrollo en serie, f(θ)≈c0P0(cosθ)+c1P1(cosθ)+c2P2(cosθ)+c3P3(cosθ)+c4P4(cosθ)

N=4; %número de términos
%coeficienets c_n
c=zeros(1,N+1);
for n=0:N %inicial
    c(n+1)=integral(@(x) legendreP(n,x).*x, 0,1)*(2*n+1)/2;
end
th=linspace(0,pi/2,100);
x=cos(th);
y=zeros(1,length(x));
for n=1:N+1 
    y=y+c(n)*legendreP(n-1,x);
end
hold on
fplot(@(x) cos(x),[0,pi/2])
plot(th,y)
set(gca,'XTick',0:pi/12:pi/2)
set(gca,'XTickLabel',{'0','\pi/12','\pi/6','\pi/4','\pi/3','5\pi/12','\pi/2'})
hold off
legend('coseno','serie')
xlabel('\theta')
ylabel('f(\theta)')
grid on
title('Polinomios de Legendre')

Temperatura de los puntos de la esfera T(r,θ) en el estado estacionario

Conocidos los coeficientes cn, calculamos los coeficientes An

n=0 n a A n P n (cosθ) =h( n=0 A n P n (cosθ) )+ n=0 c n P n (cosθ) n=0 ( ( n a +h ) A n c n ) P n (cosθ) =0 A n = a n+ha c n

La temperatura de los puntos de la esfera T(r, θ) en el estado estacionario es

T(r,θ)= n=0 A n ( r a ) n P n (cosθ) = pa{ 1 4 + 1 2 1 1+ha ( r a ) P 1 (cosθ)+ 5 16 1 2+ha ( r a ) 2 P 2 (cosθ) 3 32 5 4+ha ( r a ) 4 P 4 (cosθ)+... }

Con los datos

Representamos la temperatura T(r, θ), en el intervalo 0≤θ<π. Para los siguientes valores de r: r=a en la superficie de la esfera, para r=3a/4, r=a/2, r=a/4. Los puntos θ≥π/2 corresponden a la zona no iluminada, la temperatura de estos puntos es más fría en la superficie que en el interior

N=20; %número de términos
%coeficienets c_n
c=zeros(1,N+1);
for n=0:N %inicial
    c(n+1)=integral(@(x) legendreP(n,x).*x, 0,1)*(2*n+1)/2;
end
th=linspace(0,pi,100);
x=cos(th);
R=1;
h=0.2; %radiación
hold on
for r=[1,3/4,1/2,1/4]*R
    y=zeros(1,length(x));
    for n=1:N+1 
        y=y+(R/(n-1+h*R))*c(n)*(r/R)^(n-1)*legendreP(n-1,x);
    end
    plot(th,y)
end
line([pi/2,pi/2],[0.9,1.8],'lineStyle','--','color','k')
set(gca,'XTick',0:pi/6:pi)
set(gca,'XTickLabel',{'0','\pi/6','\pi/3','\pi/2', '2\pi/3', '5\pi/3','\pi'})
hold off
xlabel('\theta')
ylabel('f(\theta)')
legend('R','3R/4','R/2','R/4')
grid on
title('Temperatura en la esfera')

Si cambiamos el parámetro h=0.05, de modo que la esfera emita menos radiación al exterior, las temperaturas de los puntos de la esfera se incrementan, el perfil de las temperaturas es similar

Referencias

N. S. Kosholyakov, M. M. Smirnov, E. B. Gliner. Differential equations of Mathematical Physics. North-Holland Publishing Company (1964) pp. 371-373