Sistemas bi y tridimensionales (IV)

Una vez más, se encuentra de la ecuación anterior que si el coeficiente de Ti,j^n es negativo, existe la posibilidad de que se violen las leyes de la termodinámica, dando por resultado inestabilidad de la solución. Por lo tanto, se requiere que para que la solución numérica sea estable, el coeficiente Ti,j^n debe ser positivo. Entonces, una condición sencilla para la estabilidad de una solución numérica para un problema transitorio bidimensional con Δx = Δy es que

Sistemas bi y tridimensionales (III)

La expresión (5-50) y el miembro izquierdo de la ecuación (5-23) deben ser iguales para satisfacer la ley de conservación de la energía en un nodo interno Por lo tanto.

Sistemas bi y tridimensionales (II)

Una vez más, la base par ala formulación de las ecuaciones nodales es el principio de la energía. En la sección 5-4, la razón de energía neta que llega a un nodo interno, o frontera, se hizo igual a cero. Ahora,  a dicha razón de energía neta la hacemos igual a la razón de cambio de energía interna del elemento, esto es

Sistemas bi y tridimensionales (I)

Tanto la alerta rectangular como la placa de amplitud finita, L, fueron problemas transitorios unidimensionales. Las distribuciones de temperatura transitoria en sistemas bidimensionales tales como los cuerpos simétricos con respecto a un eje, o placas rectangulares resulta ser importantes por dos razones principales. Primera, dichos sistemas determinan el periodo de sobrecalentamiento - por ejemplo, el tiempo de residencia de una varilla métalica dentro de un horno antes de colocarla en un molino. Segunda, durante el estado no estacionario, dentro del cuerpo con frecuencia se crean grandes esfuerzos térmicos y la distribución de esfuerzo transitorio resulta ser de una particular importancia desde un punto de vista estructural. Algunas veces los caracteres transitorios surgen debido a la naturaleza de dependencia del tiempo de las condiciones en la frontera, como en el caso del retorno a la tierra de un vehículo espacial.

Formulación Implícita (II)

El miembro izquierdo representa el cambio de energía interna provocado por el flujo de calor asociado a los gradientes de temperatura "presentes" en el tiempo τ, según aparecen escritos en el miembro derecho. Es igualmente posible que el cambio en la energia interna se deba al flujo de calor asociado con los gradientes de temperatura en el tiempo (τ + Δτ), dando lugar a la ecuación siguiente:

La ecuación anterior contiene solamente una temperatura en el tiempo τ, esto es, Ti^n, mientras que todos los demás términos contienen temperaturas en el tiempo (τ + Δτ). Por lo tanto, las temperaturas en (τ + Δτ) no pueden resolverse explicitamente; en lugar de esto, se deben resolver las ecuaciones para todos los nodos simultáneamente, después de establecer las ecuaciones para los nodos internos y frontera. Dicha formulación implícita es estable independientemente del valor del incremento del tiempo, Δτ, que se elija. Un valor excesivamente grande de Δτ provocará grandes errores asociados inherentemente con el método de diferencias finitas. Sugerimos al lector que consulte las referencias 4, 5, 6 y 7 para una discusión acerca de la estabilidad, convergencia y errores en las formulaciones implícita y explícita.

Formulación Implícita (I)

El requisito de que se debe restringir el tamaño de Δτ para asegurar la estabilidad nos trae como consecuencia que algunas veces el incremento de tiempo debe ser muy pequeño, del orden de una fracción de segundo, especialmente cuando se involucra conducción transitoria en capas múltiples. La formulación implícita elimina esta restricción, pero involucra la solución de ecuaciones simúltaneas en cada incremento de tiempo de tiempo. Considere la ecuación (5-40), que reproducimos en seguida.

Programación de las ecuaciones - Problema 2 (IV)

Se calcula que el incremento Δτ, consistente con el criterio e estabilidad, tiene un valor de 74.8 segundos. La computadora entrega las temperaturas en la pared durante un intervalo de 29.94 minutos, cuyo valor se redondea a 30 minutos en la figura 5-17. El resultado nos muestra tambien el número de pasos de computación que se requieren para cada línea en el resultado. La figura 5-17 nos muestra la historia tiempo-temperatura de la placa. Resulta claro de la figura que, aun cuando toma tres horas y 42 minutos alcanzar el estado estacionario según se definio con anterioridad, ocurre un pequeño cambio en la temperatura después de los primeros 90 minutos.

Programación de las ecuaciones - Problema 2 (III)

Los pasos siguientes constituyen un bosquejo de las principales operaciones en el programa de la computadora.

1. Leer las propiedades físicas, y las condiciones iniciales y en la frontera, e imprimirlas.
2. Calcular Δx, α, y Δτ
3. Con Δτ  como incremento de tiempo, calcular el número de incrementos de tiempo, N, en 30 minutos.
4. Calcular Ti^n para N incrementos y guardar los Ti^ns.
5. Calcular el calor que se conduce a través de la pared y el calor que se transfiere por convección de la superficie derecha. Si las dos se encuentran dentro de un rango de uno por ciento, ir al paso (6). En caso contrario ir al paso
6. Entregar los valores de la temperatura y de la pérdida de calor y parar.

En seguida se dan el programa de computadora y las salidas.

Formulación de las ecuaciones - Problema 2 (II)

En el problema se incluye conducción de calor unidimensional no estacionario. El método para resolver el problema es el mismo que el que tratamos en la sección 5-7.1. Inicialmente, la pared de ladrillo tiene una temperatura de 20°C. Correspondiendo a esta condición en τ = 0, el valor del indice superior n en Ti^n se designa como uno. Así que

SE calcularán las temperaturas para intervalos sucesivos hasta que se alcanza el estado estacionario. En estado estacionario, el calor que se conduce a través de la pared, Qcond, es igual al calor que se transfiere por convección de la cara derecha, Qconv. Si estas dos cantidades se encuentran dentro de un rango de uno por ciento una de la otra, para fines prácticos, se consigue el estado estacionario.

Criterio de Estabilidad - Problema 2 (I)

Una pared larga de ladrillo, con espesor de 10 cm, se encuentra a una temperatura inicial de 20°C. Tiene una conductividad térmica de 0.5 W/m.K, una densidad de 1800kg/m³, y un calor específico de 840J/kg.K. En el tiempo τ = 0, la temperatura de la cara izquierda se eleva hasta 200°C mientras que la cara derecha se expone al aire cuya temperatura es de 40°C. El coeficiente convectivo de transferencia de calor en la cara derecha es de 50W/m²K. Escriba un programa de computadora para calcular la distribución de temperatura en la pared como función del tiempo. Estime el tiempo requerido para conseguir el estado estacionario.

Solución:

Datos: Se da una pared larga de ladrillo

Criterio de Estabilidad - Problema (III)

Para asegurar estabilidad de la solución númerica, requerimos una vez más que los coeficientes de Ti que aparecen en la ecuación (5-45), y de TM en la ecuación (5-47) no sean negativos. Esto nos lleva a las condiciones siguientes:

De estos dos requisitos la condición (5-48a) es más restrictiva que la condición (5-48) y, por lo tanto, debemos elegir valor de Δτ para satisfacer la condición (5-48a).

Criterio de Estabilidad - Problema (II)

Se puede encontrar la temperatura en el modo M, en el tiempo (τ + τΔ) igualando la energía neta que se transporta por utilidad de tiempo hacia el elemento asociado con el nodo M-ésimo [el miembro izquierdo de la ecuación (5-5)] con al razón de cambio de energia del elemento.

Cuando se resuelven las ecuaciones (5-45) a (5-47), nos dan los valores de las temperaturas en (τ + τΔ) para los siguientes nodos. Observe que dichas ecuaciones son algebraicas, pero no son ecuaciones simultáneas.

Criterio de Estabilidad - Problema (I)

Plantee las ecuaciones nodales para la aleta rectangular recta, suponiendo una temperatura inicial uniforme, To, a lo largo de la aleta, con la temperatura en la base de la aleta elevada súbitamente a T* en el tiempo τ = 0. La superficie de la aleta está cediendo calor por convección hacia el aire ambiente, cuya temperatura es T∞.

Solución


Formulación de las ecuaciones:

Criterio de Estabilidad (IV)

Cuando comparamos las condiciones (5-42) y (5-43), encontramos que esta última expresión nos da como resultado un incremento del tiempo menor. Por lo tanto, usamos la condición (5-43) para seleccionar el incremento de tiempo Δτ.

Si β = 2, i.e., αΔτ = (Δx)², entonces se obtiene de la ecuación (5-40a)

O bien, que la temperatura (τ + Δτ) en un nodo es igual a la media aritmética de la temperatura τ en los nodos vecinos. La ecuación (5-44) nos lleva de inmediato a una solución gráfica conocida como gráfica de Schmidt. Con ayuda de las computadoras, uno puede generar soluciones mucho más exactas que las que se obtienen por medio de las gráficas de Schmidt.

Propiedades de radiación básica (IV)

Puesto que la potencia emisiva de un cuerpo real es igual a su emisividad por la potencia emisiva de un cuerpo negro a la misma temperatura, tenemos que ε1eb1 = energía radiante emitida por la superficie 1 por unidad de tiempo y de área.

Para el equilibrio térmico del cuerpo 1, debemos tener

α1G1 = ε1eb1

Finalmente, combinando las ecuaciones (6-5) y (6-5a), concluimos que

α1 = ε1

Por lo tanto, se puede establecer que la absortividad de un cuerpo es igual a su emisividad si las temperaturas de la fuente de la radiación incidente, y del cuerpo, son iguales. Este hecho se conoce como la ley de Kirchhoff 

Criterio de Estabilidad (III)

En las referencias 4 y 5 se discuten criterios más avanzados para la estabilidad de soluciones numéricas de ecuaciones diferenciales parciales. El criterio sencillo, dado anteriormente, restringe de inmediato el valor del incremento de tiempo, Δτ, que puede usarse una vez que se fija la red (y, por tanto, Δx). Un examen semejante de la ecuación (5-41a) indica que el coeficiente de TM^n debe ser positivo. Entonces, el criterio de estabilidad se transforma en

Criterio de Estabilidad (II)

De la fisica del problema, sabemos que la temperatura en un punto cualquiera de la placa se debe encontrar en el rango de 100°F a 400°F. Aún así, la tabla anterior nos señala que de τ = 3Δτ en adelante, obtenemos valores de temperaturas que están en violación directa de las leyes de la termodinámica.

Consecuentemente, no se debe permitir un coeficiente negativo de Ti^n si se desea obtener una solución con significado físico. Dicho requisito implica automáticamente que

Criterio de Estabilidad (I)

Un examen minucioso de la ecuación (5-40a) revela que si β es menor que 2, el coeficiente de Ti^n será negativo. Un coeficiente negativo de Ti^n puede llevar a valores fluctuantes de temperaturas al repetir los cálculos para Δτs sucesivos. Pueden tener tan violentas las fluctuaciones que se violan las leyes de la termodinámica, y el número de tolerancia para la computadora para los números máximo o mínimo queda excedido.

Por ejemplo, considere una placa cuya temperatura original es de 100°F. A partir de τ = 0, haga que la cara izquierda de la placa se mantenga a 400°F y que la cara de la derecha se aisle. Al seleccionar un valor de β, observamos que se incluyen la difusividad térmica α, el espaciamiento nodal Δx, y el intervalo Δτ. El valor de α está fijo para la placa dada. Si seleccionamos β = 1, entonces sólo uno de los dos parámetros, Δx y Δτ, se puede seleccionar independientemente. Para β = 1, la ecuación (5-40a) se transforma en:

Sistema unidimensional (IV)

La temperatura en el nodo 1 será fija, mientras que en el nodo M cambiará continuamente hasta alcanzar un estado estacionario. Esta forma de las ecuaciones de diferencia se conoce como la forma explícita, ya que las temperaturas Ti^(n+1), correspondientes al tiempo ( τ + Δτ) se pueden resolver explicitamente; solamente las temperaturas Ti^n que corresponden al tiempo τ, aparecen en los miembros derecho de las ecuaciones (5-40a) y (5-41a). Esto se debe al uso de las Ti^n en el balance de energía para todos los nodos.

Sistema unidimensional (III)

La ecuación (5-4d) expresa la temperatura en el nodo i en el tiempo (τ + Δτ) en términos de la temperatura en el tiempo τ. En nuestro caso, ya que todas las temperaturas para τ = 0 son iguales a To y conocidas, es posible calcular las temperaturas en el tiempo Δτ en todos los nodos internos con un valor prescrito para β. Una vez que se han calculado las temperaturas en todos los nodos para el tiempo Δτ, se usan sus valores como entradas cuando se calculan las temperaturas en el tiempo 2Δτ, etc. Ya que la solución se obtiene sucesivamente en los tiempos Δτ, 2 Δτ, ..... se llama al procedimiento computacional un procedimiento en marcha.

La ecuación para el nodo frontera 1 es

Sistema unidimensional (II)

La energía que se conduce hacia el nodo i a través de un área transversal, A, durante un tiempo unitario está dado por

en cuya expresión el índice superior, n, representa las temperaturas en el tiempo τ = nΔτ. Utilizando un índice superior para representar el tiempo con el fin de hacer énfasis en la naturaleza no estable del problema y para mantener una distribución entre las coordenadas de espacio y la coordenada de tiempo. Estrictamente hablando, también se debe utilizar un índice superior para Q. Sin embargo, esto haría solamente engorrosa la notación y no aportaría efecto alguno al análisis. Por lo tanto, no usaremos un índice superior para Q. En vista de la naturaleza no estable del problema, existe un cambio de energía interna, ΔUi en el tiempo Δτ que está dado por

Sistema unidimensional (I)

Considere una placa de amplitud, L, en la dirección x y sin gradientes de temperatura en las direcciones y y z, suponga que la placa se encuentra inicialmente a una temperatura uniforme, To. A partir del tiempo x = 0 en adelante, se mantiene la cara izquierda a una temperatura estable T*, y se expone la cara derecha a una pérdida de calor por convección. Nos interesa obtener una solución numérica a este problema. Se divide la amplitud de una placa, L, en (M-1) partes iguales, dando por resultado M planos nodales, como se muestra en la figura 5-16. El espacio de separación de nodos en la red es exactamente el mismo que en la situación de estado estacionario. Cada uno de los nodos internos tiene material de amplitud (Δx/2) asociado a cada uno de ellos a cada lado de su plano central, mientras que los nodos frontera tienen material de amplitud (Δx/2) sólo en un lado.

Sistema en estado no estacionario

Cuando se trata con problemas en estado no estacionario, utilizando métodos numéricos, es necesario usar una red de diferencia finita para el tiempo así como para las variables de espacio. De hecho, se consigue la solución de "presenta" a "futuro" utilizando pequeños incrementos de tiempo. En cuanto más pequeño sea el incremento de tiempo, mayor es la exactitud que se logra. Si se usa un incremento de tiempo excesivamente grande, se puede tener inestabilidad numérica; o, en otras palabras, se puede obtener una solución sin sentido. Más adelante se discutirán las razones para dichos posibles comportamientos.

Al trabajar con sistemas en estado estacionario, sólo una vez se requiere una solución de ecuaciones simultáneas. En el caso de problemas en estado no estacionario, se deben tener nuevas soluciones a ecuaciones algebraicas, o bien a ecuaciones simultáneas, para cada incremento de tiempo sucesivo, hasta que se alcanza el tiempo requerido. Si se obtiene una solución en estado estacionario siguiendo un transitorio, la misma se obtiene numéricamente cuando el proceso computacional se lleva a cabo par aun número suficientemente grande de incrementos de tiempo adicionales.

Al plantear las ecuaciones para todos los puntos nodales, se usa una vez más el principio de conservación de la energía. El flujo neto de calor hacia dentro de un elemento durante un pequeño intervalo, Δτ, se hace igual al cambio de energía interna de dicho elemento, incorporando así tanto la temperatura presente como la temperatura futura en el nodo. Ahora discutiremos un caso unidimensional de estado no estacionario.

Ecuaciones de diferencias a partir de la ecuación diferencial de energía (V)

La ecuación anterior nos muestra que el uso del método de diferencia finita, basado en la formulación de la diferencia central, incluye un error que es del orden de [(Δx)^2/12]. La expansión en series de Taylor para T(i+1), j nos muestra también que la aproximación por diferencia finita de la primera derivada, ( ∂T/ ∂x)p, incluye un error del orden de (Δx/2). Así que, cuando la razón de flujo de calor, -kA(dT/dx),se evalúa utilizando la técnica de diferencia finita, ésta incluirá un error del orden de (Δx/2). Por otra parte, cuando se obtiene una distribución de temperatura utilizando la misma técnica, se incluye un error del orden de [(Δx)^2/12]

Ecuaciones de diferencias a partir de la ecuación diferencial de energía (IV)

No resulta sorprendente que esta expresión sea la misma que la ecuación (5-24a). Es recomendable efectuar un balance de energía física en los nodos de frontera para obtener ecuaciones de diferencias adecuadas. En las secciones 5-4 y 5-4.1 se consideraron muchas de las condiciones en la frontera que se encuentran con más frecuencia en la práctica, utilizando técnicas de balance de energía.

Ahora examinaremos el orden de magnitud del error que se introduce en la aproximación por diferencia finita. Para este fin, expandamos

Sumando y rearreglando estas dos ecuaciones, tenemos

en cuya expresión se substituye T(i,j) por T(xi,yj)

Ecuaciones de diferencias a partir de la ecuación diferencial de energía (III)

Entonces, la ecuación diferencial que gobierna el proceso, V²T =0, toma la siguiente forma (la forma de diferencia finita)

Ecuaciones de diferencias a partir de la ecuación diferencial de energía (II)

Observe que la primera cantidad que aparece dentro del paréntesis es la aproximación por diferencias centrales para (δT/δx) en A, mientras que el segundo término es la aproximación por diferencias centrales para (δT/δx) en B, como se muestra en la figura 5-15. Por lo tanto.

Ecuaciones de diferencias a partir de la ecuación diferencial de energía (I)

REsulta interesante deducir las ecuaciones de diferencias de una manera más formal, aunque menos física. La ecuación que gobierna el proceso de conducción de calor bidimensional en estado estacionario, con propiedades constantes y sin fuentes de calor es la siguiente.

Sistemas Tridimensionales en estado estacionario

El procedimiento para trabajar numéricamente los problemas de conducción de calor tridimensional es exactamente el mismo que el que seguimos para casos bidimensionales; la única complicación adicional consiste en que se debe hacer una red tridimensional. Para el caso de un problema tridimensional, habrá seis nodos vecinos para cada nodo interno, según se ilustra en la figura 5-14.

o bien
T1 + T2 + T3 +T4 + T5+ T6 -6TP =0

Se puede escribir las ecuaciones para los nodos frontera efectuando balances de energía adecuados.