Una varilla con temperaturas fijas en sus extremos

Supongamos una varilla metálica de longitud L, conectada por sus extremos a dos focos de calor a temperaturas Ta y Tb respectivamente. Sea T0 la temperatura inicial de la varilla cuando se conectan los focos a los extremos de la varilla.

Para calcular la temperatura T(x,t) en un punto x de la varilla y en un isntante t, resolvemos la ecuación de la conducción del calor

1 α T t = 2 T x 2 0<x<Lt>0

con las siguientes condiciones de contorno: T(0,t)=Ta, T(L,t)=Tb

y la temperatura inicial T(x,0)=T0 constante en todos los puntos de la varilla

El estado estacionario

Al cabo de cierto tiempo, teóricamente infinito, que en la práctica depende del tipo de material que empleamos, se establece un estado estacionario en el que la temperatura de cada punto de la varilla no varía con el tiempo. Dicho estado está caracterizado por un flujo J constante de energía. La ley de Fourier establece que la temperatura variará linealmente con la distancia x al origen de la varilla.

T( x, )= T a +( T b T a ) x L

Solución de la ecuación de la conducción del calor

Definimos la función u(x,t)=T(x,t)-T(x,∞), en términos de esta nueva función, la ecuación de la conducción del calor, las condición inicial en el instante t=0 y las condiciones de contorno en x=0 y en x=L se escriben

1 α u t = 2 u x 2 0<x<Lt>0 contornou(0,t)=0u(L,t)=0 inicialu(x,0)= T 0 T a ( T b T a ) x L

Buscamos una solución de la forma u(x,t)=F(xG(t), variables separadas

1 α 1 G(t) dG(t) dt = 1 F(x) d 2 F(x) d x 2 = ω 2

Integramos la primera ecuación diferencial

dG(t) dt +α ω 2 G(t)=0 G(t)=G(0)·exp( α ω 2 t )

Integramos la segunda ecuación diferencial

d 2 F(x) d x 2 + ω 2 F(x)=0

Es una ecuación diferencial similar a la de un Movimiento Armónico Simple, cuya solución es F(x)=Asin(ωx)+Bcos(ωx)

Condiciones de contorno

La solución completa u(x,t) es la superposición de los productos de las funciones Fn(xGn(x) para cada ωn=nπ/L

u(x,t)= n=1 A n sin( ω n x )exp( α ω n 2 t )

Condición inicial

Solamente, queda por determinar los coeficientes An, identificando la solución para t=0 con la condición inicial u(x,0)

u(x,0) = n = 1 A n sin ( n π L x )   u(x,0) =   T 0 T a T b T a L x

Multiplicamos ambos miembros por sin(mπx/L) e integramos entre 0 y L

0 L u(x,0) sin ( m π L x )   d x = n = 1 A n 0 L sin ( n π L x ) sin ( m π L x ) d x   

Calculamos la integral de la derecha, haciendo el cambio de variable zx/L, dz= πdx/L

0 L sin( nπ x L )sin( mπ x L )dx = L π 0 π sin(nz)sin(mz)dz

Integramos dos veces por partes

( 1 n 2 m 2 ) 0 π sin(nz)sin(mz)dz = 1 m sin(mz)cos(nz) n m 2 sin(nz)cos(mz) | 0 π =0

Cuando m=n

0 π sin 2 (nz) dz= 1 2 0 π ( 1cos(2nz) ) dz= 1 2 ( z 1 2n sin(2nz) ) | 0 π = π 2

Calculamos el coeficiente An conociendo la función u(x,0)

A n = 2 L 0 L u(x,0) sin( nπ x L )dx

Sabiendo que u(x,0)=a+bx, donde a=T0-Ta y b=-(Tb-Ta)/L

A n = 2 L 0 L ( a+bx )sin( nπ x L ) dx= 2 L { a L nπ cos( nπ x L )+b( xL nπ cos( nπ x L )+ L 2 n 2 π 2 cos( nπ x L ) ) } | 0 L = { 2 nπ ( 2a+bL )nimpar 2bL nπ npar

Finalmente,

  A n = { 2 n π ( 2 T 0 T a T b ) n  impar 2 n π ( T b T a ) n  par

Solución completa

La temperatura en cualquier punto de la varilla x, en un instante t, se compone de la suma de un término proporcional a x y de una serie rápidamente convergente que describe el estado transitorio.

u(x,t)= n=1 A n sin( nπ x L )exp( α n 2 π 2 L 2 t ) T(x,t)=u(x,t)+T(x,)= T a +( T b T a ) x L + n=1 A n sin( nπ x L )exp( α n 2 π 2 L 2 t )

El valor de α=K/(ρc) nos da una medida de la rapidez con la que el sistema alcanza el estado estacionario. Cuanto mayor sea α antes se alcanza el estado estacionario

Ejemplo

Metal Densidad, ρ Calor específico,c Conductividad térmica, K 1/α=ρc/K
Aluminio 2700 880 209.3 11352
Acero 7800 460 45 79733
Cobre 8900 390 389.6 8909
Latón 8500 380 85.5 37778
Plata 10500 230 418.7 5768
Plomo 11300 130 34.6 42457

Fuente: Koshkin N. I., Shirkévich M. G.. Manual de Física Elemental. Editorial Mir 1975. págs 36, 74-75, 85-86

Creamos la función temperatura, que admite como parámetros: las temperaturas de los extremos de la varilla, Ta y Tb, la temperatura inicial T0, la longitud L de la varilla, el parámetro 1/α, y el tiempo t en el que queremos calcular la distribución de temperaturas a lo largo de la varilla. Devuelve un vector x de puntos a lo largo de la varilla y sus correspondientes temperturas en un vector T.

function [x,T]=temperatura(Ta,Tb,T0,L,a2,t)
   x=linspace(0,L,100);
   T=zeros(length(x),1);
   if(t==0)
       T=T0;
       return;
   end
   k=zeros(2,1);  
   k(1)=2*(2*T0-(Ta+Tb))/pi;
   k(2)=2*(Tb-Ta)/pi;
   cte=pi*pi/(a2*L^2);
   
   T=Ta+(Tb-Ta)*x/L;
   n=0;
   v=1;
   while(v>0.01)
        for j=1:2
            n=n+1;
            v=exp(-cte*n^2*t);
            T=T+v*k(j)*sin(pi*n*x/L)/n;
        end
   end
end

Una serie infinita tenemos que interrumpirla en algún término de acuerdo con algún criterio. Un posible criterio sería que la exponencial

exp( α n 2 π 2 L 2 t )<0.01

sea menor que un valor prefijado, por ejemplo 0.01.

Creamos un script que establece

Defina el vector de los instantes t, en segundos, en el que queremos representar la distribución de temperaturas a lo largo de la varilla.

Llame a la función temperatura y represente mediante el comando plot cada una de las graficas T(x,t) en la misma ventana. Como ejercicio se sugiere, además:

Ta=100; %temperatura en el extremo izquierdo
Tb=30; %temperatura en el extremo derecho
T0=0; %temperatura inicial de la varilla
alfa=11352; %Coeficiente, 1/alfa del aluminio
L=0.5; % longitud de la varilla

hold on
axis([0 0.5 -5 100]);
for t=[2 30 60 90 200]
    [x,T]=temperatura(Ta,Tb,T0,L,alfa,t);
    plot(x,T,'displayName',num2str(t));
end
title('Evolución de la temperatura de una varilla')
xlabel('x')
ylabel('T')
legend('-DynamicLegend','location','northeast')
grid on
hold off  

Referencias

Puig Adam P., Curso teórico-práctico de ecuaciones diferencias aplicado a la Física y Técnica. Biblioteca Matemática (1950), págs. 300-303