PROBLEMA.

1) Etapas de discretización
Iniciamos calculando las áreas de cada nodo, obteniendo para cada una A1= 0.45 m2, A2= 0.35 m2,
A3= 0.25 m2 y A4= 0.15 m2. Así mismo AA= 0.5 m2, AB= 0.4 m2, AC= 0.3 m2, AD= 0.2 m2 y AE= 0.1
m2
Asumiendo una variación de presión lineal entre los nodos A y E se obtiene, pA=p0 = 10,0 Pa, pB=
7,5 Pa, pC= 5,0 Pa, pD= 2,5 Pa y pE= 0,0 Pa.
Las ecuaciones que rigen el ejercicio son las siguientes:

• Conservación de Masa:
𝑑
(𝜌𝐴𝑢) = 0
𝑑𝑥
• Conservación de momento:
𝑑𝑢 𝑑𝑝
𝜌𝑢𝐴 = −𝐴
𝑑𝑥 𝑑𝑥
La forma discretizada de la ecuación de momento es:
Δ𝑝
(𝜌𝑢𝐴)𝑒 𝑢𝑒 − (𝜌𝑢𝐴)𝑤 𝑢𝑤 = Δ𝑉
Δ𝑥
Con Δ𝑝 = 𝑝𝑤 − 𝑝𝑒
La ecuación del momento discretizado para este problema unidimensional se puede escribir como:
∗
𝑎𝑝 𝑢𝑝∗ = 𝑎𝑊 𝑢𝑊 + 𝑎𝐸 𝑢𝐸∗ + 𝑆𝑢

Si utilizamos el esquema de upwind, los coeficientes se pueden obtener de:

𝑎𝑊 = 𝐷𝑤 + max⁡(𝐹𝑤 , 0)
𝑎𝐸 = 𝐷𝑒 + max⁡(0, −𝐹𝑒 )
 𝑎𝑃 = 𝑎𝑊 + 𝑎𝐸 + (𝐹𝑒 − 𝐹𝑤 )
Fw y Fe son caudales másicos a través de las caras oeste y este del volumen de control u.
Calculamos las velocidades de la cara necesarias para Fw y Fe a partir de promedios de valores de
velocidad en los nodos que se encuentran a ambos lados de la cara y utilizamos los valores
correctos del área de las caras oeste y este que se dan en la tabla anterior. Al comienzo de los
cálculos utilizamos el campo de velocidad inicial generado a partir del caudal másico estimado.
Para iteraciones posteriores utilizamos la velocidad corregida obtenida después de resolver la
ecuación de corrección de presión.
El término fuente Su contiene el gradiente de presión integrado sobre el volumen de control:
Δ𝑝 Δ𝑝 1
𝑆𝑢 = ∗ Δ𝑉 = ∗ 𝐴𝑎𝑣 Δ𝑥 = Δ𝑝 ∗ (𝐴𝑤 + 𝐴𝑒 )
Δ𝑥 Δ𝑥 2
Dado que la geometría tiene un área de sección transversal variable, utilizamos un área promediada
para calcular ∆V. A primera vista, esto parece una aproximación muy burda, pero es posible
demostrar que el orden de precisión de Su no es peor que la diferenciación ascendente utilizada para
los términos de flujo de momento.
En resumen, los coeficientes de las ecuaciones u discretizadas están dados por:
𝐹𝑤 = 𝜌𝐴𝑤 𝑢𝑤 ⁡𝑦⁡𝐹𝑒 = 𝜌𝐴𝑒 𝑢𝑒
𝑎𝑊 = 𝐹𝑤
𝑎𝐸 = 0
𝑎𝑃 = 𝑎𝑊 + 𝑎𝐸 + (𝐹𝑒 − 𝐹𝑤 )
1
𝑆𝑢 = Δ𝑝 ∗ (𝐴𝑤 + 𝐴𝑒 ) = Δ𝑝 ∗ 𝐴𝑃
2
Por su parte para la ecuación de presión, la forma discretizada de la ecuación de continuidad en esta
geometría unidimensional es:
(𝜌𝑢𝐴)𝑒 − (𝜌𝑢𝐴)𝑤 = 0

La ecuación de corrección de presión correspondiente es:

𝑎𝑃 𝑝𝑃′ = 𝑎𝑊 𝑝𝑊
′
+ 𝑎𝐸 𝑝𝐸′ + 𝑏′
Donde:
𝑎𝑊 = (𝜌𝑑𝐴)𝑤 ⁡𝑦⁡𝑎𝐸 = (𝜌𝑑𝐴)𝑒
𝑏 ′ = (𝐹𝑤∗ − 𝐹𝑒∗ )
En el algoritmo SIMPLE, las correcciones de presión p′ se utilizan para calcular las correcciones de
velocidad u′ y los campos de presión y velocidad corregidos usando:
𝑢′ = 𝑑(𝑝𝐼′ − 𝑝𝐼+1
′ )

𝑝 = 𝑝∗ + 𝑝′
𝑢 = 𝑢∗ + 𝑢′
 𝑚̇
Para generar un campo de velocidad inicial para este problema usamos 𝑢 = 𝜌𝐴 junto con las áreas
de la sección transversal en los nodos de velocidad para calcular la velocidad inicial. Así se obtiene
velocidades iniciales u1= 2.222 m/s, u2= 2.857 m/s, u3= 4.0 m/s y u4= 6.666 m/s

Ecuación para nodo 2:

m m
𝑢1 + 𝑢2 2.222 s + 2.857 s m
𝐹𝑤 = (𝜌𝑢𝐴𝐵 )𝑤 = 1.0 ∗ ∗ 0.4 = 1.0 ∗ ∗ 0.4 = 1.016⁡
2 2 s
m m
𝑢2 + 𝑢3 2.857⁡ s + 4.0⁡ s m
𝐹𝑒 = (𝜌𝑢𝐴𝐶 )𝑒 = 1.0 ∗ ∗ 0.3 = 1.0 ∗ ∗ 0.3 = 1.029⁡
2 2 s
𝑎𝑊 = 𝐹𝑤 = 1.016⁡
𝑎𝐸 = 0
𝑎𝑃 = 𝑎𝑊 + 𝑎𝐸 + (𝐹𝑒 − 𝐹𝑤 ) = 1.016 + 0 + (1.029 − 1.016) = 1.029
De donde podemos calcular d2:
𝐴2 0.35
𝑑2 = = = 0.340
𝑎𝑝 1.029

Obtenemos para Su:

𝑆𝑢 = Δ𝑃 ∗ 𝐴2 = (𝑝𝐵 − 𝑝𝐶 ) ∗ 𝐴2 = (7.5 − 5) ∗ 0.35 = 0.875
Así, la ecuación discretizada para el nodo 2 es:
𝑎𝑝 ∗ 𝑢2 = 𝐹𝑊 ∗ 𝑢1 + 𝑆𝑢

1.029 ∗ 𝑢2 = 1.016 ∗ 𝑢1 + 0.875

Ecuación para nodo 3:

m m
𝑢2 + 𝑢3 2.857 s + 4.0 s m
𝐹𝑤 = (𝜌𝑢𝐴𝐶 )𝑤 = 1.0 ∗ ∗ 0.3 = 1.0 ∗ ∗ 0.3 = 1.029⁡
2 2 s
m m
𝑢3 + 𝑢4 4.0⁡ s + 6.666 s m
𝐹𝑒 = (𝜌𝑢𝐴𝐷 )𝑒 = 1.0 ∗ ∗ 0.2 = 1.0 ∗ ∗ 0.2 = 1.067⁡
2 2 s
𝑎𝑊 = 𝐹𝑤 = 1.029⁡
𝑎𝐸 = 0
𝑎𝑃 = 𝑎𝑊 + 𝑎𝐸 + (𝐹𝑒 − 𝐹𝑤 ) = 1.029 + 0 + (1.067 − 1.029) = 1.067
De donde podemos calcular d3:
𝐴3 0.25
𝑑3 = = = 0.234
𝑎𝑝 1.067
Obtenemos para Su:
𝑆𝑢 = Δ𝑃 ∗ 𝐴3 = (𝑝𝐶 − 𝑝𝐷 ) ∗ 𝐴3 = (5 − 2,5) ∗ 0.25 = 0.625
Así, la ecuación discretizada para el nodo 3 es:
𝑎𝑝 ∗ 𝑢3 = 𝐹𝑊 ∗ 𝑢2 + 𝑆𝑢

1.067 ∗ 𝑢3 = 1.029 ∗ 𝑢2 + 0.625

Para los nodos 1 y 2 se sigue el siguiente procedimiento tal que:

Ecuación nodo 1:
Utilizando la ecuación de Bernoulli para calcular la presión estática en A:
1
𝑝𝐴 = 𝑝0 − (𝜌𝑢𝐴2 )
2
Y escribiendo uA de la siguiente manera:
𝐴1
𝑢𝐴 = 𝑢1
𝐴𝐴
Con lo que obtenemos para Bernoulli:

1 𝐴1 2
𝑝𝐴 = 𝑝0 − (𝜌⁡ (𝑢1 ) ⁡)
2 𝐴𝐴

Ahora escribimos la ecuación de momento discretizado para el volumen 1 de con ayuda del
esquema upwind:
𝐹𝑒 𝑢1 − 𝐹𝑤 𝑢𝐴 = (𝑝𝐴 − 𝑝𝐵 ) ∗ 𝐴1
En donde Fw:
𝐹𝑤 = 𝜌𝐴𝐴 𝑢𝐴 = 𝜌𝐴1 𝑢1
Obteniendo:

𝐴1 1 𝐴1 2
𝐹𝑒 𝑢1 − 𝐹𝑤 𝑢1 = ((𝑝0 − (𝜌⁡ (𝑢1 ) ⁡)) − 𝑝𝐵 ) ∗ 𝐴1
𝐴𝐴 2 𝐴𝐴

𝐴1 1 2 1
𝐴2
𝐹𝑒 𝑢1 − 𝐹𝑤 𝑢1 = ((𝑝0 − (𝜌𝑢1 ⁡ 2 ⁡)) − 𝑝𝐵 ) ∗ 𝐴1
𝐴𝐴 2 𝐴𝐴

1 𝐴 2
Utilizando 𝐹𝑤 = 𝜌𝐴𝐴 𝑢𝐴 = 𝜌𝐴1 𝑢1 y 𝑝𝐴 = 𝑝0 − 2 (𝜌⁡ (𝑢1 𝐴 1 ) ⁡):
𝐴

𝐴1 1 1 𝐴1
𝐹𝑒 𝑢1 − 𝐹𝑤 𝑢1 = ((𝑝0 − (𝜌𝑢1 𝐴1 𝑢1 ⁡ ⁡)) − 𝑝𝐵 ) ∗ 𝐴1
𝐴𝐴 2 𝐴𝐴 𝐴𝐴
 𝐴1 𝐹𝑤 𝐴12
𝐹𝑒 𝑢1 − 𝐹𝑤 𝑢1 = 𝑝0 𝐴1 − 𝑢1 ⁡ 2 − 𝑝𝐵 𝐴1
𝐴𝐴 2 𝐴𝐴

⁡𝐴1 𝐹𝑊 𝐴1 2
[𝐹𝑒 − 𝐹𝑤 + ∗ ( ) ] 𝑢1 = (𝑝0 − 𝑝𝐵 )𝐴1
𝐴𝐴 2 𝐴𝐴

Esto se puede transformar a:

𝐹𝑊 𝐴1 2 ⁡𝐴1 𝑎𝑛𝑡𝑖𝑔𝑢𝑜
[𝐹𝑒 + ∗ ( ) ] 𝑢1 = (𝑝0 − 𝑝𝐵 )𝐴1 + −𝐹𝑤 𝑢
2 𝐴𝐴 𝐴𝐴 1

Entonces:

𝐹𝑊 𝐴1 2
𝑎𝑃 = 𝐹𝑒 + ∗( )
2 𝐴𝐴

Calculando:
⁡𝐴1 0.45
𝑢𝐴 = 𝑢1 = 2.222 ∗ =2
𝐴𝐴 0.5
𝐹𝑤 = (𝜌𝑢𝐴)𝑤 = 𝑝𝑢𝐴 𝐴𝐴 = 1 ∗ 2 ∗ 0.5 = 1
𝑢1 + 𝑢2 2.222 + 2.857
𝐹𝑒 = (𝜌𝑢𝐴)𝑒 = 1 ∗ ∗ 0.4 = 1 ∗ ∗ 0.4 = 1.016
2 2
𝑎𝑊 = 0⁡
𝑎𝐸 = 0
Obteniendo para ap:

𝐹𝑊 𝐴1 2 1 0.45 2
𝑎𝑃 = 𝐹𝑒 + ∗ ( ) = 1.016 + ∗ ( ) = 1.421
2 𝐴𝐴 2 0.5
Con lo que:
𝐴1 0.45
𝑑1 = = = 0.317
𝑎𝑝 1.421

Por su parte aplicando p0=10 Pa y una velocidad inicial de u1=2.222 m/s obtenemos para Su:
𝐴1 0.45
𝑆𝑢 = (𝑝0 − 𝑝𝐵 )𝐴1 + 𝐹𝑤 ( ) ∗ 𝑢1 = (10 − 7.5) ∗ 0.45 + 1 ∗ ∗ 2.222 = 3.125
𝐴𝐴 0.5
Así, la ecuación discretizada para el nodo 1:
𝑎𝑝 ∗ 𝑢1 = 𝑆𝑢

1.421⁡𝑢1 = 3.125

Ecuación nodo 4:
 𝑢3 + 𝑢4
𝐹𝑤 = (𝜌𝑢𝐴)𝑤 = 1 ∗ ∗ 0.2 = 1.067
2
𝐹𝑒 = (𝜌𝑢𝐴)4 = 1 ∗ 𝑢4 ∗ 𝐴4 = 1⁡(𝑝𝑎𝑟𝑎⁡𝑝𝑟𝑖𝑚𝑒𝑟𝑎⁡𝑖𝑡𝑒𝑟𝑎𝑐𝑖ó𝑛⁡𝑝𝑜𝑟⁡𝑙𝑎⁡𝑠𝑢𝑝𝑜𝑠𝑖𝑐𝑖ó𝑛)
𝑎𝑊 = 𝐹𝑤 = 1.067⁡
𝑎𝐸 = 0
𝑎𝑃 = 𝑎𝑊 + 𝑎𝐸 + (𝐹𝑒 − 𝐹𝑤 ) = 1.067 + 0 + (1.0 − 1.067) = 1.0
De donde podemos calcular d4:
𝐴4 0.15
𝑑4 = = = 0.15
𝑎𝑝 1.0

Obtenemos para Su:

𝑆𝑢 = Δ𝑃 ∗ 𝐴𝑎𝑣 = (𝑝𝐷 − 𝑝𝐸 ) ∗ 𝐴4 = (2.5 − 0) ∗ 0.15 = 0.375
Así, la ecuación discretizada para el nodo 4 es:
𝑎𝑝 ∗ 𝑢4 = 𝐹𝑊 ∗ 𝑢3 + 𝑆𝑢

1.0 ∗ 𝑢4 = 1.067 ∗ 𝑢3 + 0.375

Así, con todas las ecuaciones se obtiene para la primera iteración:

1.421⁡𝑢1 = 3.125
1.029 ∗ 𝑢2 = 1.016 ∗ 𝑢1 + 0.875
1.067 ∗ 𝑢3 = 1.029 ∗ 𝑢2 + 0.625
1.0 ∗ 𝑢4 = 1.067 ∗ 𝑢3 + 0.375
U1 (m/s) = 2.199 U3 (m/s) = 3.500
U2 (m/s) = 3.023 U4 (m/s) = 4.109

Corrección de la presión:
Presión en nodo B:
Iniciamos calculando:
𝑎𝑊 = (𝜌𝑑𝐴)1 = 1 ∗ 0.317 ∗ 0.45 = 0.143
𝑎𝐸 = (𝜌𝑑𝐴)2 = 1 ∗ 0.340 ∗ 0.35 = 0.119
𝐹𝑤∗ = (𝜌𝑢∗ 𝐴)1 = 1 ∗ 2.199 ∗ 0.45 = 0.990
𝐹𝑒∗ = (𝜌𝑢∗ 𝐴)2 = 1 ∗ 3.023 ∗ 0.35 = 1.058
Con lo que:
 𝑎𝑃 = 𝑎𝑊 + 𝑎𝐸 = 0.143 + 0.119 = 0.262
𝑏 ′ = 𝐹𝑤∗ − 𝐹𝑒∗ = 0.990 − 1.058 = −0.068
Así la ecuación de corrección de la presión en el nodo B es:
𝑎𝑝 ⁡𝑝𝐵′ = 𝑎𝑊 ⁡𝑝𝐴′ + 𝑎𝐸 𝑝𝐶′ + 𝑏′

0.262⁡𝑝𝐵′ = 0.143⁡𝑝𝐴′ + 0.119⁡𝑝𝐶′ − 0.068

Presión en nodo C:
Iniciamos calculando:
𝑎𝑊 = (𝜌𝑑𝐴)2 = 1 ∗ 0.340 ∗ 0.35 = 0.119
𝑎𝐸 = (𝜌𝑑𝐴)3 = 1 ∗ 0.234 ∗ 0.25 = 0.059
𝐹𝑤∗ = (𝜌𝑢∗ 𝐴)2 = 1 ∗ 3.023 ∗ 0.35 = 1.058
𝐹𝑒∗ = (𝜌𝑢∗ 𝐴)3 = 1 ∗ 3.500 ∗ 0.25 = 0.875
Con lo que:
𝑎𝑃 = 𝑎𝑊 + 𝑎𝐸 = 0.119 + 0.059 = 0.178
𝑏 ′ = 𝐹𝑤∗ − 𝐹𝑒∗ = 1.058 − 0.875 = 0.183
Así la ecuación de corrección de la presión en el nodo B es:
𝑎𝑝 ⁡𝑝𝐶′ = 𝑎𝑊 ⁡𝑝𝐵′ + 𝑎𝐸 𝑝𝐷′ + 𝑏′

0.178⁡𝑝𝐶′ = 0.119⁡𝑝𝐵′ + 0.059𝑝𝐷′ + 0.183

Presión en nodo D:
Iniciamos calculando:
𝑎𝑊 = (𝜌𝑑𝐴)3 = 1 ∗ 0.234 ∗ 0.25 = 0.059
𝑎𝐸 = (𝜌𝑑𝐴)4 = 1 ∗ 0.15 ∗ 0.15 = 0.023
𝐹𝑤∗ = (𝜌𝑢∗ 𝐴)3 = 1 ∗ 3.500 ∗ 0.25 = 0.875
𝐹𝑒∗ = (𝜌𝑢∗ 𝐴)4 = 1 ∗ 4.109 ∗ 0.15 = 0.616
Con lo que:
𝑎𝑃 = 𝑎𝑊 + 𝑎𝐸 = 0.059 + 0.023 = 0.082
𝑏 ′ = 𝐹𝑤∗ − 𝐹𝑒∗ = 0.875 − 0.616 = 0.259
Así la ecuación de corrección de la presión en el nodo D es:
𝑎𝑝 ⁡𝑝𝐷′ = 𝑎𝑊 ⁡𝑝𝐶′ + 𝑎𝐸 𝑝𝐸′ + 𝑏′

0.082⁡𝑝𝐷′ = 0.059⁡𝑝𝐶′ + 0.023⁡𝑝𝐸′ + 0.259

Podemos considerar que en cada nivel de iteración la presión estática pA se fija temporalmente por
p0 y el valor actual de u*1, justificando así el uso de pA′ = 0,0. Así, las correcciones de presión se
establecen en cero para ambos nodos A y B, por lo que PA’=0 y PB’=0.
Así para la primera iteración se obtiene:
0.262⁡𝑝𝐵′ = 0.143⁡𝑝𝐴′ + 0.119⁡𝑝𝐶′ − 0.068
0.178⁡𝑝𝐶′ = 0.119⁡𝑝𝐵′ + 0.059𝑝𝐷′ + 0.183

0.082⁡𝑝𝐷′ = 0.059⁡𝑝𝐶′ + 0.259

p’A (Pa)= 0 p’D (Pa)= 6.208
p’B (Pa)= 1.639 p’E (Pa)= 0
p’C (Pa)= 4.175
Las presiones nodales corregidas ahora se calculan utilizando estas correcciones de presión:
𝑝𝐵 = 𝑝𝐵∗ + 𝑝𝐵′ = 7.5 + 1.639 = 9.139
𝑝𝐶 = 𝑝𝐶∗ + 𝑝𝐶′ = 5 + 4.175 = 9.175

𝑝𝐷 = 𝑝𝐷∗ + 𝑝𝐷′ = 2.5 + 6.208 = 8.708

Las velocidades corregidas al final de la primera iteración son:
𝑚
𝑢1 = 𝑢1∗ + 𝑑1 (𝑝𝐴′ − 𝑝𝐵′ ) = 2.199 + 0.317 ∗ (0 − 1.639) = 1.68
𝑠
𝑚
𝑢2 = 𝑢2∗ + 𝑑2 (𝑝𝐵′ − 𝑝𝐶′ ) = 3.023 + 0.34 ∗ (1.639 − 4.17461) = 2.16
𝑠
𝑚
𝑢2 = 𝑢2∗ + 𝑑2 (𝑝𝐵′ − 𝑝𝐶′ ) = 3.023 + 0.34 ∗ (1.639 − 4.17461) = 2.16
𝑠
𝑚
𝑢3 = 𝑢3∗ + 𝑑3 (𝑝𝐶′ − 𝑝𝐷′ ) = 3.501 + 0.234 ∗ (4.17461 − 6.208) = 3.024
𝑠
𝑚
𝑢4 = 𝑢4∗ + 𝑑4 (𝑝𝐷′ − 𝑝𝐸′ ) = 4.109 + 0.15 ∗ (6.208 − 0) = 5.04
𝑠
Y para A se puede calcular la presión mediante:

1 𝐴1 2 1 0.45 2
𝑝𝐴 = 𝑝0 − 𝜌𝑢12 ( ) = 10 − ∗ 1 ∗ 1.682 ∗ ( ) = 8.857
2 𝐴𝐴 2 0.5

2) Esquema de ecuaciones obtenidos previamente.

El algoritmo para este punto se encuentra desarrollado en el archivo de Excel que se adjunto en la
entrega.

• También puedo ingresar al archivo dando click AQUÍ.

 3) Distribución de velocidad y presión en dirección x a lo largo de la tobera.

Velocidad
4,5000

4,0000

3,5000

3,0000

2,5000
Algoritmo
2,0000
Bernoulli
1,5000

1,0000

0,5000

0,0000
0 0,5 1 1,5 2 2,5

Presión
12

10

6 Bernoulli
Algoritmo
4

0
0 0,1 0,2 0,3 0,4 0,5 0,6

4) Comparación de los resultados de la malla de 4 volúmenes con la forma teórica de

Bernoulli.
Distribución de velocidad y presión:
Después de 19 iteraciones se tiene un error menor a 10-5, arrojando los siguientes valores:
U1 (m/s) = 1.3321 U3 (m/s) = 2.3978
U2 (m/s) = 1.7127 U4 (m/s) = 3.9963
 p’A (Pa) = 9.2813 p’D (Pa) = 5.7488
p’B (Pa) = 8.3571 p’E (Pa) = 0
p’C (Pa) = 7.6575

Con la ecuación de Bernoulli se hubiese obtenido:

1 2 1 𝑚2
𝑝0 = 𝑝𝑁 + 𝜌𝑢𝑁 = 𝑝𝑁 + 𝜌 ⁡
2 2 (𝜌𝐴𝑁 )2
De los datos del problema tenemos p0 = 10 Pa, ρ = 1 kg/m3 y N = E, por lo que AN = AE = 0,1 m2,
lo que produce K = 0,44721 kg/s. Con esto las presiones serían:
U1 (m/s) = 0.9938 pA (Pa)= 9.6
U2 (m/s) = 1.2778 pB (Pa)= 9.375
U3 (m/s) = 1.7889 pC (Pa)= 8.889
U4 (m/s) = 2.9814 pD (Pa)= 7.5
pE (Pa)= 0
Si tomamos como reales los valores calculados con Bernoulli se obtiene que los porcentajes de error
de los resultados del algoritmo son:
U1 (m/s) = 34,04% p’A (Pa)= 3,32%
U2 (m/s) = 34,04% p’B (Pa)= 10,86%
U3 (m/s) = 34,04% p’C (Pa)= 13,85%
U4 (m/s) = 34,04% p’D (Pa)= 23,35%

Se observa que los resultados obtenidos con el algoritmo son muy alejados de los reales obtenidos
mediante la ecuación de Bernoulli.

5) Modificación del número de volúmenes en el algoritmo.

El algoritmo para este punto se encuentra desarrollado en el archivo de Excel que se adjuntó en la
entrega; para el cual se tomó 8 volúmenes en su desarrollo.

• También puedo ingresar al archivo dando click AQUÍ.

6) Análisis y conclusiones.
El análisis revela que, inicialmente, el caudal másico calculado para la malla de cinco nodos
muestra una sobreestimación significativa con respecto al valor exacto. Este exceso disminuye
gradualmente conforme se refinan las mallas con más nodos, lo que sugiere una mejora progresiva
hacia la solución exacta. Por ejemplo, al emplear mallas de 15, 30 y 100 nodos, los caudales
másicos obtenidos son muchos más precisos con respecto a los valores obtenidos mediante la
ecuación de bernoulli, demostrando una reducción del error conforme aumenta la densidad de la
malla. Este patrón indica que el refinamiento sistemático de la malla puede reducir eficazmente las
discrepancias en la solución. Además, al aumentar aún más la cantidad de nodos, se logra una
convergencia del caudal másico calculado muy cercana al valor exacto.

PREGUNTAS - MÉTODOS DE SOLUCIÓN

a) ¿Cuáles son las diferencias entre resolver un problema de flujo de fluidos

analíticamente en comparación con resolverlo numéricamente? ¿Cuáles son las
ventajas y desventajas de cada método?
Sabemos que para resolver un problema de flujo de fluidos tenemos varios métodos para llegar a
ello. Primero nos enfocamos en la solución analítica, esta se basa en el uso de ecuaciones y
soluciones exactas en forma de expresiones matemáticas explícitas. Ahora para el caso de resolverlo
numéricamente, nos basamos en utilizar métodos de discretización para las ecuaciones de flujo,
dando así soluciones aproximadas a través de simulaciones numéricas, que se consigues haciendo
conjeturas y probando si el problema se resuelve lo suficientemente bien como para detenerse. En
muchos casos, se utilizan ambos métodos en conjunto, utilizando soluciones analíticas para validar
y comprender los resultados numéricos.

• Ventajas y Desventajas de la resolución analítica.

VENTAJAS DESVENTAJAS
- Soluciones exactas - Difícilmente aplicable a problemas
complejos con geometrías irregulares o
- Comprensión teórica
condiciones de frontera complicadas.
- Bajos costos computacionales - Simplificaciones Necesarias (Asumir
flujo laminar, incompresible, etc…)
que pueden no representar la realidad.

• Ventajas y Desventajas de la resolución numérica.

VENTAJAS DESVENTAJAS
- Variedad de métodos de discretización - Errores numéricos; es decir, soluciones
aproximadas y no exactas
- Aplicable a casos y geometrías
complejas - Requiere altos costos computacionales
- Útil en aplicaciones prácticas y
problemas reales donde las condiciones
 ideales no se cumplen.
- Permite visualizar el comportamiento
del flujo a través de gráficos y
simulaciones, facilitando la
interpretación de los resultados.

b) ¿Cuáles son las principales ventajas y desventajas de la discretización de las

ecuaciones gobernantes mediante el método de diferencias finitas, método de
elementos finitos y método de volúmenes finitos?
La discretización de ecuaciones rectoras utilizando diferentes métodos numéricos, como el método
de diferencias finitas, el método de elementos finitos y el método de volúmenes finitos, tiene cada
uno su propio conjunto de ventajas y desventajas:
1. Método de diferencias finitas:
Ventajas:
• Fácil de implementar y comprender, lo que lo convierte en una buena opción para
principiantes.
• Puede incorporar aproximaciones de diferenciación de orden superior en cuadrículas
regulares para aumentar la precisión.
• Muy adecuado para cuadrículas estructuradas y problemas con geometrías regulares.
• Eficiente computacionalmente para problemas con geometrías regulares.
• Es menos susceptible a errores de convergencia que otros métodos.
Desventajas:
• Las propiedades de conservación no se aplican de manera inherente y requieren un cuidado
especial para su mantenimiento.
• Menos flexible que otros métodos para refinar localmente la malla.
• Puede requerir una distribución de cuadrícula uniforme para obtener resultados precisos, lo
que dificulta el trabajo con cuadrículas no uniformes.
• La precisión puede ser limitada, especialmente para problemas con geometrías complejas o
cuando se requieren soluciones de alta resolución.
• Dificultad para manejar geometrías complejas: Adaptar la malla a geometrías complejas
puede ser difícil y laborioso.

2. Método de elementos finitos:

Ventajas:
• Muy adecuado para geometrías complejas y dominios irregulares debido a su flexibilidad en
la generación de malla.
• Puede manejar varios tipos de condiciones de contorno más fácilmente.
• Proporciona una buena aproximación para soluciones con precisión de orden superior.
• Buena precisión, especialmente para problemas estructurales.
• Es muy flexible y se puede aplicar a una amplia gama de problemas de CFD.
Desventajas:
• Más complejo de implementar en comparación con el método de diferencias finitas.
• Requiere más recursos computacionales debido a la necesidad de resolver grandes sistemas
de ecuaciones.
• Puede que no conserve propiedades como la masa o la energía a menos que se tengan en
cuenta específicamente.
• Menos intuitivo para problemas de flujo que el método de volúmenes finitos.
• La convergencia de la solución puede ser más difícil de garantizar que con las diferencias
finitas.

3. Método de volumen finito:

Ventajas:
• Conserva naturalmente propiedades como masa, momento y energía, lo que lo hace
adecuado para problemas de flujo de fluidos.
• Puede manejar rejillas no estructuradas, lo que permite una mayor flexibilidad en la
generación de mallas.
• Puede manejar geometrías complejas mejor que FDM.
• Buena precisión y estabilidad para una amplia gama de problemas de flujo.
• Es menos susceptible a errores de convergencia que el MEF.
Desventajas:
• Es más desafiante desarrollar aproximaciones de diferenciación de orden superior en tres
dimensiones debido a la necesidad de interpolación e integración.
• Requiere un tratamiento cuidadoso de las condiciones de contorno para garantizar las
propiedades de conservación.
• Puede ser computacionalmente más costoso en comparación con el método de diferencias
finitas para ciertos problemas.
• Puede ser más difícil de implementar que FDM, especialmente en 3D y menos preciso que
FEM para problemas estructurales o de elasticidad.
• Puede ser difícil manejar problemas con difusión dominante.

c) Describa la manera en que se acoplan las ecuaciones de conservación de masa y

cantidad de movimiento en el Método de Volúmenes Finitos mediante el esquema
SIMPLEC.

Usando la cuadrícula escalonada que se muestra en la Fig. b y c, las ecuaciones de volumen finito
para u y u, respectivamente. Reescribiendo la ecuación del momento u para el volumen de control
centrado en e para mostrar explícitamente el término de presión se obtiene:

𝑎𝑒 𝑢𝑒 = ∑ 𝑎𝑛𝑏 𝑢𝑛𝑏 + 𝐴𝑒 (𝑝𝑃 − 𝑝𝐸 ) + 𝑏𝑒 (1)

donde p es la presión, A, es el área de la cara del volumen de control P en e, y:

𝑎𝑒 = ∑ 𝑎𝑛𝑏 ⁡(2)

Para cualquier distribución de presión estimada p*, las velocidades u* obtenidas al resolver las
ecuaciones de momento satisface:

∗
𝑎𝑒 𝑢𝑒∗ = ∑ 𝑎𝑛𝑏 𝑢𝑛𝑏 + 𝑏𝑒 + 𝐴𝑒 (𝑝𝑃∗ − 𝑝𝐸∗ )(3)

La aproximación SIMPLE llega a 𝑎𝑒 𝑢′𝑒 = ∑ 𝑎𝑛𝑏 𝑢′𝑛𝑏 + 𝐴𝑒 (𝑝𝑃∗ − 𝑝𝐸∗ ), donde posteriormente asume
∑ 𝑎𝑛𝑏 𝑢′𝑛𝑏 puede ignorarse, mientras que se retiene en el lado izquierdo el término 𝑎𝑒 𝑢′𝑒 que es
igual ∑ 𝑎𝑛𝑏 𝑢′𝑒 , que de acuerdo a SIMPLEC, este es de similar magnitud al ignorado, por lo que
considera que ignorar ese término es incoherente.
Se usa la corrección de la presión p' =p -p* y u' = u - u*. La relación entre p' y u' se obtiene restando
la ecuación (3) de la ecuación (1):

(𝑎𝑒 − ∑ 𝑎𝑛𝑏 ) 𝑢′𝑒 = ∑ 𝑎𝑛𝑏 ⁡(𝑢′𝑛𝑏 − 𝑢𝑒′ ) + 𝐴𝑒 (𝑝𝑃∗ − 𝑝𝐸∗ )⁡(4)

Así, en SIMPLEC, el término ∑ 𝑎𝑛𝑏 ⁡(𝑢′𝑛𝑏 − 𝑢𝑒′ ) se desprecia ya que se iguala a 0, ya que
∑ 𝑎𝑛𝑏 𝑢′𝑛𝑏 y ∑ 𝑎𝑛𝑏 𝑢′𝑒 son de similar magnitud. Reemplazando u' por u-u* y adoptando la
aproximación SIMPLEC, la ecuación (4) se convierte en:
 (𝑎𝑒 − ∑ 𝑎𝑛𝑏 ) (𝑢𝑒 − 𝑢∗ 𝑒 ) = 𝐴𝑒 (𝑝𝑃∗ − 𝑝𝐸∗ )

𝐴𝑒
(𝑢𝑒 − 𝑢∗ 𝑒 ) = (𝑝∗ − 𝑝𝐸∗ )
𝑎𝑒 − ∑ 𝑎𝑛𝑏 𝑃
𝐴𝑒
𝑢𝑒 − 𝑢 ∗ 𝑒 = (𝑝∗ − 𝑝𝐸∗ )
𝑎𝑒 − ∑ 𝑎𝑛𝑏 𝑃
𝐴𝑒
𝑢𝑒 = 𝑢 ∗ 𝑒 + (𝑝∗ − 𝑝𝐸∗ )
𝑎𝑒 − ∑ 𝑎𝑛𝑏 𝑃
𝑢𝑒 = 𝑢𝑒∗ + 𝑑𝑒 (𝑝𝑃′ − 𝑝𝐸′ )
Donde:
𝐴𝑒
𝑑𝑒 =
𝑎𝑒 − ∑ 𝑎𝑛𝑏

d) ¿Por qué son importantes los esquemas de upwind para un flujo fuertemente
convectivo?
El término "upwind" se refiere a un tipo de métodos de discretización numérica para ecuaciones
diferenciales parciales hiperbólicas, donde se utilizan variables ascendentes para calcular las
derivadas en un campo de flujo. Específicamente, en la dinámica de fluidos computacional, el
esquema de diferenciación upwind se emplea para resolver problemas de convección-difusión.
Este método se presenta como una alternativa al esquema de diferenciación central, el cual, aunque
es más preciso, no permite modelos rápidos y estables y no converge en presencia de flujos
convectivos fuertes. Este problema se puede resolver con el esquema upwind, que está diseñado
para manejar flujos convectivos intensos con efectos de difusión reducidos. Por lo tanto, el esquema
de diferenciación central se usa solo para flujos de baja velocidad, mientras que para flujos de alta
velocidad se prefiere el upwind.
Las principales características del esquema upwind son las siguientes:
✓ Estabilidad: Es estable cuando se calculan flujos de convección con |Pe|>2
✓ Conservación: Utiliza expresiones coherentes para calcular los flujos a través de las caras
de las celdas (formulación conservadora)
✓ Acotamiento: Los coeficientes de la ecuación discretizada son siempre positivos y
satisfacen los requisitos de acotamiento
✓ Transportabilidad: Tiene en cuenta la dirección del flujo, por lo que la transportabilidad
está integrada en la formulación
✓ Exactitud: Se basa en la fórmula de diferenciación backward, por lo que la exactitud es sólo
de primer orden sobre la base del error de truncamiento de la serie de Taylor
Los esquemas upwind son estables, cumplen con el requisito de transportabilidad, pero tienen una
precisión de primer orden y, por lo tanto, son propensos a errores de difusión numérica. Estos
errores pueden minimizarse utilizando discretizaciones de orden superior. Otro inconveniente de
este esquema es que produce resultados incorrectos cuando el flujo no está alineado con las líneas
de la malla.