La solución teórica 1D para una deposición frontal de potencia constante q (q=2000 W/m2) en un medio semiinfinito es (Problem 2 in Table 5 of Energy deposition in unbounded media, in Heat conduction):

>

eqT:=T(x,t)-T0=(q*x/k)*Eta;Eta:=exp(-eta^2)/(sqrt(Pi)*eta)+erf(eta)-1;eqeta:=eta=x/sqrt(4*a*t);eqT0:=T(0,t)-T0=2*q*sqrt(t/(Pi*rho*c*k));plot(Eta,eta=0..1.2,'Eta'=0..5);eqT_:=subs(eqeta,eqT):eqTp:=dT/dx=diff(rhs(eqT_),x):eqTp0:=(dT/dx)[x=0]=convert(series(rhs(eqTp),x=0,2),polynom);T_T0:=subs(eqeta,rhs(eqT));T_T0_:=subs(eqq_,dat2,dat,SI0,T_T0);pl2:=plot([seq(subs(t=10*i,T_T0_),i=1..10)],x=0..subs(dat2,SI0,L),'T_T0'=0..25);plot(subs(eqq_,dat2,SI0,rhs(eqT0)),t=0..100,T0_C=0..25);

i.e. se trata de una solución de semejanza en la variable eta, con la función Eta(eta) representada en la figura de arriba (comprobamos que en x=0 se cumple q=-k*dT/dx).

En el gráfico T(x) se ha representado la evolución T(x,t)-T0 a lo largo del espesor para diversos tiempos (entre cada línea transcurren 10 s) para la placa 2 (la gruesa) hasta los tc=100 s antes calculados. No es comparable este último perfil con el perfil lineal del caso estacionario del apartado anterior (que iba de T-T0=16 ºC a 0, porque aquí no se fija la Text-T0=0 (aparte de que aquí se envían todos los 2000 W/m2 hacia este lado).

En el gráfico T(t), último, se ha representado la evolución T(0,t)-T0 en función del tiempo (i.e. la evolución de la temperatura del calefactor). Nótese que la misma resistencia eléctrica puede usarse para calefacción y como termómetro (a partir de la relación entre resistividad y temperatura), por lo que si se miden a la vez q y T(0,t), la ecuación eqT0 anterior permite determinar la conductividad térmica del material (en realidad solo la efusividad, e=sqrt(k·rho·c)).

c) Evolución de los flujos de calor y de los saltos de temperaturas con el tiempo, para cada placa, aplicado el modelo semi-infinito a cada una.

Se sabe que el reparto de flujos de calor en el contacto entre sólidos semi-infinitos es proporcional a las efusividades térmicas, e=sqrt(rho·c·k), como se desprende la la eqT0 anterior; i.e.: como la T(x=0) ha de ser la misma, q1/e1=q2/e2 y q1+q2=q, de dionde:

>	q1:=q*sqrt(rho1*c1*k1)/(sqrt(rho1*c1*k1)+sqrt(rho2*c2*k2));q2:=q*sqrt(rho2*c2*k2)/(sqrt(rho1*c1*k1)+sqrt(rho2*c2*k2));q_:=rhs(eqq_);q1_:=subs(q=q_,rho1=rho,c1=c,k1=k,dat1,rho2=rho,c2=c,k2=k,dat2,q1);q2_:=subs(q=q_,rho1=rho,c1=c,k1=k,dat1,rho2=rho,c2=c,k2=k,dat2,q2);

	(2)

i.e. de los 2000 W/m2, 775 W/m2 van a 1 (entran en la placa delgada) y 1225 W/m2 van a 2 (entran en la placa gruesa), y este reparto no varía con el tiempo mientras sea válido el modelo doblemente semi-infinito. En valores absolutos, de los 20 W, 12,3 W van hacia la derecha y 7,7 W hacia la izquierda. Como comprobación, aplicamos la solución semiinfinita a cada lado:

>	eqT0_1:=evalf(subs(q=q1_,dat1,SI0,rhs(eqT0)));eqT0_2:=evalf(subs(q=q2_,dat2,SI0,rhs(eqT0)));plot([eqT0_1,eqT0_2],t=0..1e2,Tmax_T0=0..15);

efectivamente, coinciden, y vemos que al cabo de unos 100 s ha subido unos 15 ºC la temperatura en el contacto. Si calculamos y dibujamos las temperatures exteriores en función del tiempo:

>	plot([subs(x=1e-6,q=q1_,dat1,dat,SI0,rhs(eqT_)),subs(x=L,q=q1_,dat1,dat,SI0,rhs(eqT_))],t=0..300,T1_T0=0..30);plot([subs(x=1e-6,q=q2_,dat2,dat,SI0,rhs(eqT_)),subs(x=L,q=q2_,dat2,dat,SI0,rhs(eqT_))],t=0..300,T2_T0=0..30);

donde se aprecia que el salto en la placa 1 va creciendo, pero en la 2 se estabiliza. Podemos calcular estos saltos con el modelo estacionario y los flujos conocidos:

>	eqDT:=q=k*DT/L;DT_:=q*L/k;DT1_:=subs(q=q1_,dat1,dat,DT_);DT2_:=subs(q=q2_,dat2,dat,DT_);

	(3)

i.e. entre las caras de la placa 1, DT tiende a unos 31 ºC (en la figura ya llega a unos 20 ºC al cabo de 300 s), mientras que entre las caras de la placa 2 (de mayor espesor, pero más difusiva), DT tiende a unos 12 ºC (como se aprecia en la figura)..

d) Simulación numérica de la evolución transversal de las temperaturas, T(x,t), al conectar el calentador desde un instante inicial en el que todo está a 0 ºC, suponiendo que las caras externas son adiabáticas, y comparación con el resultado del modelo semi-infinito.

El caso de bordes adiabáticos se ha resuelto con Matlab, y se ha obtenido:

y con el modelo doblemente semi-infinito sería:

>

T_T0:=subs(eqeta,rhs(eqT));T1_T0:=subs(q=q1_,dat1,dat,SI0,T_T0);pl1:=plot([seq(subs(x=-x,t=10*i,T1_T0),i=1..10)],x=-subs(dat1,SI0,L)..0,'T_C'=0..15):T2_T0:=subs(q=q2_,dat2,dat,SI0,T_T0);pl2:=plot([seq(subs(t=10*i,T2_T0),i=1..10)],x=0..subs(dat2,SI0,L),'T_C'=0..15):display({pl1,pl2});T1_T0_100_L:=evalf(subs(t=100,x=L,dat1,SI0,T1_T0));T1_T0_100_0:=evalf(subs(t=100,x=1e-9,T1_T0));T2_T0_100_0:=evalf(subs(t=100,x=1e-9,T2_T0));T2_T0_100_L:=evalf(subs(t=100,x=L,dat2,SI0,T2_T0));

	(4)

i.e., con el modelo doblemente semi-infinito, al cabo de 100 s se alcanza una Tmax-T0=13.8 ºC (en el calefactor), con calentamientos en los extremos de T2-T0=4,9 ºC a la izquierda de la placa 1, y T1-T0=0,5 ºC a la derecha de la placa 2, mientras que con el modelo de bordes aislados era Tmax-T0=15,2 ºC, T2-T0=9,7 ºC y T1-T0=1,3 ºC. Nótese que todo se calienta más en el caso de bordes aislados, pero en el contacto apenas se nota (Tmax-T0=15,2 ºC frente a 13,8 ºC), mientras que en los bordes la diferencia es mucho mayor (e.g. a la derecha de la placa 2 es T2-T0=9,7 ºC frente a los T1-T0=4,9 ºC del modelo semi-infinito).

e) Evolución T(x,t) si se apaga el calentador tras 100 s de funcionamiento.

Resuelto con Matlab para bordes aislados.

i.e. si los extremos están aislados, toda la energía queda almacenada en las placas, que adquieren una temperatura uniforme de Tmean-T0=Q*Dt/(mc)=8 ºC (nótese la rapidez conque se amortigua el pico en el contacto, x=0).

En realidad (y con el modelo semiinfinito) toda la energía se disipará y las placas quedarían a la T inicial.

f) Estimación de las pérdidas laterales en las placas, debido al menor tamaño del calentador.

En el caso más desfavorable (más caliente), la placa 1 tiene una temperatura media en sus 4 mm de espesor de unos 5 ºC sobre la Tb=0 ºC, y la placa 2 de unos 12 ºC sobre la Tb. Pero esto es solo un modelo 1D que será muy bueno en el centro del calentador de 0,1·0,1 m2, pero irá perdiendo validez a medida que nos alejemos a lo largo del eje y (o z). Una estimación de las pérdidas de calor en esa dirección es (tomando como longitud característica Lh=longitud del heater):

>

Qy:=k*Ay*(Tm-T0)/Lh;Ay=Lh*L;Lh_:=subs(dat,SI0,sqrt(A))*m_;Ay1_:=subs(dat1,Lh_*L);'Ay1_'=%*1e6*mm_^2/m_^2;Ay2_:=subs(dat2,Lh_*L);'Ay2_'=%*1e6*mm_^2/m_^2;;Tm1:=5*K_+T0;Qy1_:=subs(dat1,dat,k*Ay1_*(Tm1-T0)/Lh_);Tm2:=T0+12*K_;Qy2_:=subs(dat2,dat,k*Ay2_*(Tm2-T0)/Lh_);Qlost:=4*(Qy1+Qy2);Qlost_:=4*(Qy1_+Qy2_);

	(5)

i.e.  el modelo 1D era muy bueno ya que las pérdidas laterales (sumando ambas placas y contando con las 4 direcciones), apenas son de 0,5 W frente a las 20 W del calentador. Nótese que casi todas las pérdidas laterales son debidas a la placa 2 porque se calienta más.

>