Stabilized finite elements for Bingham and Herschel-Bulkley confined flows. Part I : Formulation
Under a Creative Commons license
En este trabajo se presenta una metodología para la resolución de las ecuaciones de Navier-Stokes para los fluidos viscoplásticos de Bingham y de Herschel-Bulkley mediante el método de los elementos finitos mixtos estabilizados velocidad/presión. Se desarrolla una formulación teórica y se realiza la implementación computacional.
Los fluidos viscoplásticos se caracterizan por presentar una tensión de corte mínima, denominada tensión de fluencia. Por encima de esta tensión de corte mínima el fluido comienza a moverse. En caso de no superarse esta tensión de fluencia, el fluido se comporta como un cuerpo rígido o cuasirrígido, con velocidad de deformación nula.
Se presentan inicialmente las ecuaciones de Navier-Stokes para un fluido incompresible. Se incluye una revisión de los modelos reológicos viscoplásticos. Se hace una descripción detallada de los mismos. Se describen los modelos viscoplásticos regularizados de Papanastasiou. Se proponen modelos regularizados de doble viscosidad como alternativa a los comúnmente usados.
Se formula el modelo discreto, así como la formulación estabilizada con los métodos de subescalas algebraica (Algebraic SubGrid Scale [ASGS]), de subescalas ortogonales (Orthogonal Subgrid Scale [OSS]) y de subescalas ortogonales desacopladas (split -OSS).
La metodología descrita en este trabajo proporciona la base para el desarrollo de una herramienta computacional para estudiar flujos viscoplásticos confinados, muy comunes en la industria.
This work presents a methodology for the solution of the Navier-Stokes equations for Bingham and Herschel-Bulkley viscoplastic fluids using stabilized mixed velocity/pressure finite elements. The theoretical formulation is developed and implemented in a computer code.
Viscoplastic fluids are characterized by a minimum shear stress called yield stress. Above this yield stress, the fluid is able to flow. Below this yield stress, the fluid behaves as a quasi-rigid body, with zero strain-rate.
First, the Navier-Stokes equations for incompressible fluid are presented. A review of the viscoplastic rheological models is included, with a detailed description of these models. The regularized viscoplastic models due to Papanastasiou are described. Double viscosity regularized models are proposed as an alternative to the models commonly used.
The discrete model is developed, and the Algebraic SubGrid Scale (ASGS) stabilization method, the Orthogonal Subgrid Scale (OSS) method and the split orthogonal subscales method are introduced.
The methodology proposed in this work provides a computational tool to study confined viscoplastic flows, common in industry.
Método de los elementos finitos estabilizados ; Subescalas ortogonales OSS ; Incompresibilidad ; Fluidos viscoplásticos ; Modelo de Bingham ; Modelo de Herschel-Bulkley
Stabilized finite elements ; Orthogonal sub-grid scales OSS ; Incompressibility ; Viscoplastic fluid ; Bingham model ; Herschel-Bulkley model
En el presente trabajo se presenta la formulación continua y su correspondiente versión discreta para modelos de elementos finitos mixtos velocidad/presión de flujos confinados viscoplásticos de Bingham y de Herschel-Bulkley.
Los fluidos viscoplásticos de Bingham y de Herschel-Bulkley son fluidos no-newtonianos que se caracterizan por presentar una tensión de corte mínima, denominada «tensión de fluencia». Por encima de esta tensión de corte mínima el fluido comienza a moverse. En caso de no superar esta tensión de fluencia, el fluido se comporta como un cuerpo rígido o cuasirrígido, con velocidad de deformación nula.
En la industria, los fluidos de Bingham pueden modelar el comportamiento de las pinturas, de los plásticos, de productos alimenticios como la mayonesa y el kétchup, entre otros. Los fluidos de Herschel-Bulkley incluyen, por ejemplo, el comportamiento de las pastas, algunos geles y los fluidos de perforación. En el medio ambiente, estos fluidos pueden modelar flujos de detritos, entre otros.
El movimiento de los fluidos isotérmicos se describe mediante las ecuaciones de conservación de masa y momentum , representado por las ecuaciones de Navier-Stokes. Numerosos ensayos experimentales han demostrado que las ecuaciones de Navier-Stokes bajo condiciones isotérmicas describen exactamente el flujo incompresible de los fluidos. Las ecuaciones de Navier-Stokes requieren de una ecuación constitutiva para caracterizar el tipo de fluido. Esta ecuación define el valor de las tensiones en función de la dinámica del flujo y está asociada con la viscosidad del fluido (modelo reológico).
El modelo que dio inicio al estudio de los materiales viscoplásticos fue el modelo plástico de Bingham [1] , formulado por Eugene C. Bingham para describir el comportamiento de las pinturas. El modelo de Herschel-Bulkley [2] se considera como un modelo generalizado de Bingham, aunque ha sido menos estudiado que el modelo de Bingham.
Ambos fluidos exhiben una fuerte discontinuidad en su comportamiento reológico debido a la existencia de la tensión de fluencia que es difícil de tratar numéricamente. Para solventar este problema, autores como Bercovier y Engelman [3] , Tanner y Milthorpe [4] y Beris et al. [5] , entre otros, han propuesto diferentes formulaciones regularizadas. Tanner y Milthorpe fueron los primeros que simularon el problema utilizando un modelo de doble viscosidad aplicable a ambos fluidos. Beris y sus colegas centraron sus estudios en el fluido de Bingham, utilizando el criterio de Von Mises [6] en las zonas de no fluencia y el modelo ideal de Bingham en la zona de fluencia. En 1987, Papanastasiou [7] propuso un modelo regularizado aplicable tanto en las zonas de no fluencia como en las zonas de fluencia para estos 2 fluidos. Souza Mendes y Dutra (SMD) [8] han propuesto recientemente una modificación del modelo de Papanastasiou.
En el presente trabajo se proponen nuevos modelos regularizados para el fluido de Bingham y el fluido de Herschel-Bulkley como alternativa a los modelos regularizados comúnmente usados.
En el caso de los materiales viscoplásticos, el método numérico más utilizado es el método de los elementos finitos (MEF) [7] , [9] , [10] and [11] . Para abordar el problema de flujo incompresible mediante el MEF, se emplea la formulación mixta de velocidad/presión (u /p). La formulación estándar de Galerkin presenta 2 fuentes de inestabilidades.
La primera es la presencia del término convectivo en las ecuaciones de gobierno que puede resultar en oscilaciones numéricas en el campo de la velocidad. La segunda fuente de inestabilidad es la combinación inapropiada de espacios de interpolación para los campos de velocidad y presión. Esta falta de estabilidad produce oscilaciones numéricas en el campo de las presiones. Para que el problema discreto sea estable, los espacios de interpolación usados para la velocidad y la presión deben satisfacer la condición inf-sup de compatibilidad o condición de Babuška-Brezzi [12] . La formulación de igual interpolación lineal usada en este trabajo no cumple con la condición Babuška-Brezzi [13] .
En ambos casos el problema necesita estabilizarse para poder probar convergencia a la solución del problema. Los métodos de estabilización más usados en la actualidad están basados en los métodos de subescalas. Hughes fue el pionero en estos métodos de subescalas (SubGrid Scale [SGS]), proponiendo el método de estabilización de subescalas algebraicas (Algebraic SubGrid Scale stabilization method [ASGS]) [14] para una ecuación escalar de difusión-reacción. Codina [15] amplió esta aproximación algebraica aplicándola a sistemas escalares multidimensionales.
Posteriormente, Codina [15] propuso adoptar un espacio de subescalas ortogonales al espacio de los elementos finitos, fundamentando así el método de estabilización de subescalas ortogonales (Orthogonal Subscale Stabilization method [OSS]). El método OSS se ha aplicado al problema de Stokes, al problema de convección-difusión-reacción y a las ecuaciones de Navier-Stokes, entre otros [15] and [16] . La estabilización OSS ha sido reformulada en una nueva versión del método llamada estabilización split -OSS [17] , computacionalmente más ventajosa. Actualmente se usan en problemas muy variados, tanto de mecánica de fluidos [15] , [17] , [18] , [19] , [20] and [21] como de mecánica de sólidos [19] , [22] , [23] , [24] , [25] , [26] , [27] , [28] , [29] , [30] and [31] .
En la parte i de este trabajo se presenta el problema del flujo confinado con un amplio desarrollo de los modelos constitutivos para flujos viscoplásticos de Bingham y de Herschel-Bulkley, así como los modelos viscoplásticos regularizados propuestos en este trabajo. Finalmente, se presenta el modelo discreto incorporando los modelos de Bingham y de Herschel-Bulkley. Como métodos de estabilización se discuten los métodos ASGS, OSS y split -OSS.
El problema continuo de dinámica de fluidos incompresibles e isotérmicos puede resolverse completamente considerando las ecuaciones de Navier-Stokes.
Las ecuaciones de Navier-Stokes, planteadas inicialmente para fluidos newtonianos, pueden usarse conjuntamente con los modelos reológicos viscoplásticos de Bingham y Herschel-Bulkley en la ecuación constitutiva.
Considérese Ω, un dominio abierto y acotado dimensional de Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \mathbb{R}^d}
, donde d = 2 o 3 es el número de dimensiones del espacio, Γ = ∂Ω es su contorno, que puede ser dividido en el contorno con condiciones de Dirichlet (velocidad impuesta) Γd = ∂Ωd y el contorno con condiciones de Neumann (tracciones impuestas) Γn = ∂Ωn de forma que Γ = Γd ∪ Γn , [0, T ] es el intervalo de tiempo de análisis.
El problema de Navier-Stokes consiste en encontrar una velocidad u y una presión p tal que:
|
( 1) |
con
|
( 2) |
y
|
( 3) |
en Ω, t ∈ [0,t ], donde σ es el tensor de tensiones, τ es el tensor de tensiones desviadoras, f es el vector de fuerzas de volumen y ρ es la densidad del fluido.
A esas ecuaciones deben añadirse condiciones iniciales de la forma u = u0 en Ω, t0 = 0 y condiciones de contorno:
Condiciones de Dirichlet:
|
( 4) |
Condiciones de Neumann:
|
( 5) |
donde n es el vector unitario normal al contorno ∂ Ω. Por simplicidad se tomará en el contorno Γd la velocidad Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \bar{\mbox{u}}=0,t\in [0,T]}
. El vector t es el vector de tracción sobre el contorno con condiciones de Neumann.
Las líneas de corriente son curvas tangentes en cada punto al campo de velocidades. En un flujo estacionario, las líneas de corriente no varían con el tiempo, mientras que en flujo transitorio sí lo hacen.
Las líneas de corriente para un flujo bidimensional con un campo de velocidades u = (ux,uy ) coinciden con las líneas de nivel de la función ϕ , solución de la ecuación laplaciana:
|
( 6) |
con la condición de contorno ϕ = 0.
La ecuación constitutiva relaciona las tensiones con la presión y la velocidad de deformación. En el caso de los fluidos, esta relación se denomina también modelo reológico.
El tensor de tensiones σ se descompone en su parte volumétrica y desviadora como:
|
( 7) |
donde p es la presión, I es el tensor de identidad de segundo orden y τ es el tensor de las tensiones desviadoras.
Para un fluido newtoniano, usando la hipótesis de Stokes y usando la ecuación de incompresibilidad, el tensor desviador de tensiones se expresa como:
|
( 8) |
donde u es el vector de velocidades, μ es la viscosidad dinámica (constante en caso de fluido newtoniano) y ɛ (·) es el gradiente simétrico de la velocidad:
|
( 9) |
donde ▿u es el gradiente de la velocidad y (▿ut ) es la transpuesta del mismo.
El valor de la magnitud del tensor de la velocidad de deformación, γ, se toma como la raíz del segundo invariante del tensor simétrico:
|
( 10) |
La magnitud del tensor desviador o tensión efectiva, τ , se toma como la raíz del segundo invariante del tensor de las tensiones desviadoras:
|
( 11) |
De acuerdo con lo anterior, la ecuación (7) se puede escribir como:
|
( 12) |
Dependiendo de los valores de la viscosidad en función de la velocidad de deformación, Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \mu =\mu \left(\dot{\gamma }\right)}
, pueden distinguirse diferentes ecuaciones constitutivas que representan diferentes modelos reológicos no-newtonianos. Un amplio número de estos materiales pueden verse en Bird et al. [32] . Se estudian en este trabajo los modelos para fluidos viscoplásticos, en particular para el modelo de Bingham y de Herschel-Bulkley.
Eugene C. Bingham describió las pinturas con este modelo en 1919, publicado en su libro Fluidity and Plasticity[1] . El modelo fue analizado por Oldroyd [33] , Reiner [34] y Prager [35] . Los plásticos de Bingham requieren de una tensión de corte mínima, τy , a partir de la cual comienzan a moverse (fig. 1 ).
|
Figura 1. Curvas reológicas. Modelo bilineal. |
En el modelo de Bingham la viscosidad está dada por:
|
( 13) |
donde μ0 es la viscosidad plástica, y la viscosidad aparente Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \mu (\dot{\gamma })}
disminuye con el incremento en la magnitud de la velocidad de deformación Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \dot{\gamma }} ; τ es la magnitud del tensor de tensiones desviadoras.
En consecuencia, el tensor de tensiones desviadoras es:
|
( 14) |
Para definir si una partícula del fluido se mueve o no, es decir, si está en fluencia o no, se comprueba si la magnitud del tensor de tensiones desviadoras, τ , excede o no el valor de la tensión de fluencia, τy . Cuando la magnitud del tensor de tensiones del fluido, τ , supera la tensión de fluencia, el comportamiento es similar al de un fluido newtoniano; en caso contrario, el fluido no presenta deformaciones por corte.
4.1.1.1. Número de Bingham
Para estudiar flujos de Bingham se define un número adimensional denominado número de Bingham, Bn . El número de Bingham, sugerido por Bird et al. [32] , se define como:
|
( 15) |
donde V es una velocidad característica del flujo viscoplástico, H es una longitud característica y μ es la viscosidad del fluido de Bingham.
El número de Bingham relaciona la tensión de fluencia, τy , con la tensión ocasionada por una velocidad de deformación característica Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \left(\dot{\gamma }_0=\frac{V}{H}\right)}
.
En el caso de un fluido newtoniano el valor Bn es nulo, Bn = 0; en el límite opuesto, para fluidos en no fluencia (sólido) el número de Bingham puede tener valores muy altos, Bn → ∝ [10] .
4.1.1.2. Tensión adimensional de fluencia
En flujos viscoplásticos es conveniente mostrar los resultados en función de una tensión de fluencia adimensional, Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \tau _y^\ast }
, definida por Papanastasiou [7] como:
|
( 16) |
donde VN es una velocidad característica tomada como la velocidad promedio del líquido newtoniano de la misma viscosidad del modelo viscoplástico.
4.1.1.3. Número de Reynolds
El número de Reynolds, Rn , define si el régimen del flujo es laminar o turbulento. El empleado en este trabajo es el que Scott et al. [36] y Mitsoulis y Huilgol [37] , entre otros, usan en trabajos previos para flujo newtoniano. Para el flujo laminar de Bingham:
|
( 17) |
donde ρ es la densidad del fluido, VB es la velocidad promedio del flujo, H es una longitud característica y μ es la viscosidad para el fluido de Bingham.
En el modelo plástico de Herschel-Bulkley [2] se combinan la tensión de fluencia y la ley potencial. En el modelo de Herschel-Bulkley la viscosidad aparente está dada por:
|
( 18) |
Al igual que para el modelo de Bingham, los materiales de Herschel-Bulkley requieren de una tensión de corte mínima, τy , para que el material fluya. Para niveles de tensión por encima de la tensión de fluencia, el material fluye con una relación no lineal tensión-velocidad de deformación como un fluido pseudoplástico (n < 1) o dilatante (n > 1) determinado por el exponente de la ley de potencia (n) .
El tensor desviador resulta:
|
( 19) |
Si n = 1, se tiene el fluido de Bingham como caso particular [35] y el índice de consistencia es igual a la viscosidad plástica del material k = μ0 . Si la tensión de fluencia es nula, τy = 0, se recupera la ley potencial.
4.2.1.1. Número generalizado de Bingham
Para materiales que obedecen el modelo de Herschel-Bulkley se define el número generalizado de Bingham, Bn * o número de Oldroyd [33] , Od , como:
|
( 20) |
donde τy es la tensión de fluencia, H es una longitud característica, k es el índice de consistencia y n es el índice potencial. La velocidad V es una velocidad característica.
4.2.1.2. Número de Reynolds
El número de Reynolds usado en flujos de Herschel-Bulkley viene dado por la ley potencial [4] and [38] :
|
( 21) |
Los modelos reológicos ideales con tensión de fluencia presentan 2 problemas:
.
Estos problemas no constituyen una limitación en soluciones analíticas para problemas simples, pero sí constituyen un serio inconveniente de cara a la solución numérica [7] , [39] and [40] .
Para evitar estas dificultades y lograr una conveniente formulación computacional, se han propuestos diferentes modelos regularizados.
Un modelo muy utilizado es el modelo de doble viscosidad, inicialmente propuesto por Tanner y Milthorpe [4] . Otro de los modelos más usados en nuestros días es el modelo de Papanastasiou [7] . Souza Mendes y Dutra (SMD) [8] han propuesto recientemente una modificación del modelo de Papanastasiou.
El modelo introducido originalmente por Tanner y Milthorpe [4] es un modelo con doble viscosidad lineal que regulariza el fluido de Bingham.
Este modelo de doble viscosidad sustituye el comportamiento rígido del modelo ideal para valores de tensiones por debajo de la tensión de fluencia por una dependencia lineal entre la tensión y la velocidad de deformación. El modelo para el fluido de Bingham es:
|
( 22) |
donde μr es la viscosidad crítica y Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \dot{\gamma }_c}
es la velocidad de deformación crítica. Para esta velocidad de deformación crítica la tensión es:
|
( 23) |
Por tanto, la velocidad de deformación crítica es:
|
( 24) |
El valor de μr debe ser grande para aproximar el modelo ideal. Una buena aproximación recomendada por Beverly y Tanner es tomar Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle 300\leq \frac{\mu _r}{\mu _0}\leq 1000}
.
El inconveniente de este modelo es que para la viscosidad crítica se tiene una tensión crítica τc algo mayor que la tensión de fluencia, como se muestra en la figura 1 . Esta tensión crítica es:
|
( 25) |
Uno de los intentos por solventar la limitación debida a la singularidad de la viscosidad para Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \dot{\gamma }\rightarrow 0}
se debe a Papanastasiou [7] , que propuso una regularización exponencial para el término de la tensión de fluencia del modelo de Bingham. La misma idea se ha usado posteriormente con el modelo de Herschel-Bulkley.
La ventaja que presenta el modelo es que describe con una sola ecuación tanto las zonas de fluencia como las de no fluencia, mediante una función suavizada de la viscosidad que depende de la velocidad de deformación y de un parámetro de regularización (m) , modificando la viscosidad aparente Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \mu (\dot{\gamma })}
del modelo ideal de la manera:
|
( 26) |
En la figura 2 se puede apreciar la influencia del parámetro de regularización m .
|
Figura 2. Modelo regularizado de Papanastasiou para el fluido de Bingham con diferentes valores del parámetro de regularización, m. |
La viscosidad en la ecuación (26) está acotada cuando el gradiente de la velocidad de deformación tiende a cero. Desarrollando en serie de Taylor y despreciando los términos de más de segundo orden, se tiene que:
|
( 27) |
Para valores muy altos del parámetro de regularización, m , este valor límite de la viscosidad puede causar problemas numéricos. En este caso es aconsejable definir un valor de truncamiento μt < μmax para velocidades de deformación muy bajas Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \dot{\gamma }<\dot{\gamma }_c}
.
El modelo regularizado SMD de Souza Mendes y Dutra [8] es similar al modelo de Papanastasiou, pero la regularización exponencial afecta a todos los términos de la viscosidad:
|
( 28) |
Además, el parámetro de regularización m se sustituye por un parámetro reológico que depende de la viscosidad del fluido a cero cizallamiento ηo y la tensión de fluencia τy , m = η0 /τy .
El límite para la viscosidad cuando la velocidad de deformación es cero es:
|
( 29) |
Para el fluido de Herschel-Bulkley se proponen modelos regularizados análogos a los del fluido de Bingham.
El modelo de Tanner y Milthorpe [4] para el fluido de Herschel-Bulkley es:
|
( 30) |
donde μr es la viscosidad crítica y Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \dot{\gamma }_c}
es la velocidad de deformación crítica. La velocidad de deformación crítica se obtiene resolviendo la siguiente ecuación no lineal implícita:
|
( 31) |
La regularización propuesta por Papanastasiou (fig. 3 ) es también aplicable al modelo de Herschel-Bulkley. La viscosidad aparente queda definida como:
|
( 32) |
|
Figura 3. Modelo regularizado de Papanastasiou para el fluido de Herschel-Bulkley con diferentes valores del parámetro de regularización m. |
La influencia del parámetro m en el fluido de Herschel-Bulkley-Papanastasiou puede verse en la figura 3 .
Al igual que para el modelo de Bingham, el modelo de Herschel-Bulkley requiere en la implementación numérica un valor de truncamiento μt .
El valor límite de la viscosidad cuando la velocidad de deformación tiende a cero varía de acuerdo con el valor n . El límite para cada término de la ecuación (32) y el límite resultante para la viscosidad se muestran en la tabla 1 . Se observa que para fluidos pseudoplásticos (n < 1) la viscosidad no está acotada. En estos casos es imprescindible la aplicación del procedimiento de truncamiento.
Exponente | Términos de la viscosidad | Viscosidad | |
---|---|---|---|
n | Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \underset{\dot{\gamma }\rightarrow 0}{lim}k\dot{\gamma }^{n-1}} | Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \underset{\dot{\gamma }\rightarrow 0}{lim}\frac{\tau _y}{\dot{\gamma }}[1-exp(-m\dot{\gamma })]} | Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \underset{\dot{y}\rightarrow 0}{lim}\mu (\dot{\gamma })} |
Si n > 1 | 0 | mτy | mτy |
Si n = 1 | k = μ | mτy | μ + mτy |
Si n < 1 | ∞ | mτy | ∞ |
El modelo SMD (Souza-Mendez-Dutra) regulariza el modelo ideal de Herschel-Bulkley en la forma:
|
( 33) |
con m = η0 /τy
El valor de la viscosidad cuando la velocidad de deformación tiende a cero se muestra en la tabla 2 . La viscosidad está acotada para cualquier valor de n , lo cual es una ventaja en la implementación numérica.
Exponente | Términos de la viscosidad | Viscosidad | |
---|---|---|---|
n | Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \underset{\dot{\gamma }\rightarrow 0}{lim}k\dot{\gamma }^{n-1}[1-exp(-m\dot{\gamma })]} | Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \underset{\dot{\gamma }\rightarrow 0}{lim}\frac{\tau _y}{\dot{\gamma }}[1-exp(-m\dot{\gamma })]} | Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \underset{\dot{\gamma }\rightarrow 0}{lim}\mu (\dot{\gamma })]} |
Si n > 1 | 0 | mτy | mτy |
Si n = 1 | 0 | mτy | mτy |
Si n < 1 | 0 | mτy | mτy |
Estos modelos presentados han formado parte de los estudios de diferentes problemas resueltos por Papanastasiou [7] , Kelessidis et al. [41] , Westerberg et al. [42] y Dall’Onder dos Santos et al. [8] , entre otros.
Se proponen a continuación sendos modelos viscoplásticos regularizados para los fluidos de Bingham y de Herschel-Bulkley.
Ambos son modelos de doble viscosidad, basados en los modelos descritos anteriormente.
Presentan, por tanto, viscosidad constante, pero de diferente magnitud en las zonas de fluencia (por encima de la velocidad de deformación crítica) y en las de no fluencia (por debajo de la velocidad de deformación crítica).
Este modelo es idéntico al modelo bilineal, pero la viscosidad crítica, μr , se toma igual al valor límite regularizado de Papanastasiou y SMD correspondiente a Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \dot{y}=0}
, esto es, μr = mτy , en función del parámetro de regularización m ( fig. 4 ). Por tanto:
|
( 34) |
|
Figura 4. Modelo ideal de Bingham y modelo regularizado Bingham-DV. |
En este caso, el valor de la velocidad de deformación crítica es:
|
( 35) |
En la figura 5 se compara el modelo bilineal propuesto con los modelos de Papanastasiou y SMD para m = 10 s , 100 s con tensión de fluencia τy = 10 Pa y μ0 = 0.2 Pas .
|
Figura 5. Comparación entre el modelo bilineal propuesto para el fluido de Bingham (DV) con el modelo regularizado de Papanastasiou y el modelo SMD para m = 10, 100 s. |
De forma análoga, se propone un modelo regularizado de doble viscosidad para el fluido de Herschel-Bulkley, en el que la viscosidad crítica se toma como el valor límite de la viscosidad del modelo de Papanastasiou y SMD, μr = mτy , en función del parámetro de regularización m ( fig. 6 ). Esto es:
|
( 36) |
|
Figura 6. Modelo ideal de Herschel-Bulkley y modelos regularizados de Herschel-Bulkley-Papanastasiou y Hershel-Bulckley-DV, n >1. |
El valor crítico de la velocidad de deformación cuando n ≠ 1 ha de determinarse de forma iterativa imponiendo continuidad en las tensiones de corte (para Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \dot{\gamma }=\dot{\gamma }_c}
):
|
( 37) |
En esta sección se describe el modelo discreto de elementos finitos para las ecuaciones de Navier-Stokes, correspondiente al modelo continuo que se describe en la sección anterior.
La estrategia de discretización adoptada en este trabajo consiste en 2 pasos. El primero es discretizar las ecuaciones en el tiempo usando un esquema de integración de diferencias finitas, y luego, la aproximación de elementos finitos se desarrolla en el espacio. Este procedimiento desacopla errores que vienen de la discretización temporal y espacial. Es de destacar que la estrategia más común es al revés: primero se discretiza en el espacio y luego en el tiempo.
El segundo y tercer término de la ecuación (1) es el término convectivo y el de la ecuación constitutiva para fluidos no-newtonianos, en particular para los modelos viscoplástico de Bingham y de Herschel-Bulkley, con Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \mu =\mu (\dot{\gamma })}
. Ambos términos son no lineales y, por tanto, es necesaria la linealización del mismo. En este trabajo se linealizan el término convectivo y el término viscoso con el método de Picard por su robustez.
La presencia de una derivada temporal en la ecuación (1) precisa de un algoritmo de integración en el tiempo. La discretización temporal puede hacerse por diferencias finitas usando la regla trapezoidal generalizada. Este es el método de diferencias finitas más simple y de un solo paso. Incluye como caso particular el método de diferenciación hacia atrás de Euler Backward Differentiation Formula (BDF1), entre otras posibilidades.
Para la discretización espacial mediante el método de elementos finitos es necesario construir subespacios discretos Vh ⊂ V , y Qh ⊂ Q que aproximen los espacios continuos.
Sean Vh y Qh los espacios de elementos finitos para interpolar las funciones vectoriales (velocidad) y escalares (presión), respectivamente, y sea Ω una partición de elementos finitos Ω = ∪ Ωe , e = 1, …, nele , donde nele es el número de elementos.
En la formulación estándar de Galerkin se toman las funciones de test iguales a las funciones de forma, así que vh ∈ Vh y qh ∈ Qh .
El cálculo con elementos finitos de flujos incompresibles con la formulación estándar de Galerkin deben estabilizarse debido a que la formulación de igual interpolación lineal tanto para la velocidad como para la presión usada en este trabajo no cumple con la condición Babuška-Brezzi. La referencia [13] ofrece una descripción comparativa de varios de los métodos de estabilización propuestos en las últimas décadas.
Los métodos de estabilización más usados en la actualidad están basados en los métodos de subescalas [43] and [44]. Estos métodos consisten en descomponer la solución, por ejemplo, la velocidad u en 2 componentes Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \mbox{u}=\mbox{u}_h+\bar{\mbox{u}}}
; una componente uh , resuelta en la escala de la malla de elementos finitos considerada y una subescala Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \bar{\mbox{u}}} , que no puede ser capturada por la partición de elementos finitos y que se resuelve analíticamente. La aproximación particular usada para la escala submalla define el modelo numérico.
La solución Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \bar{\mbox{u}}}
en la escala fina se obtiene a partir del residuo de la solución en la escala gruesa. La solución de esta escala fina se localiza en el interior de cada elemento finito y se supone Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \bar{\mbox{u}}=0} en el contorno de los elementos.
El residuo de la ecuación de momento en la escala grande, Rh , resulta en:
|
( 38) |
La subescala de velocidad, Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \bar{\mbox{u}}}
, se aproxima de distinta forma en cada método de estabilización.
En este trabajo se utilizarán el método ASGS y el método OSS en los problemas para flujo confinado. El método split -OSS se utiliza para problemas grandes, como en flujos con superficie libre.
En el método de las subescalas algebraicas (ASGS), Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \bar{\mbox{u}}}
se toma proporcional al residuo Rh
|
( 39) |
donde τ1 es un parámetro numérico.
En el método de las subescalas ortogonales (OSS) la subescala Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \bar{\mbox{u}}}
se toma proporcional a la proyección ortogonal de dicho residuo Rh :
|
( 40) |
donde Ph es la proyección sobre el espacio de los elementos finitos y Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle P_h^{\perp }=\mbox{I}-P_n}
es la proyección ortogonal.
Comparando ambos métodos se observa que la diferencia radica en sustituir el término Rh de la versión ASGS por Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle P_h^{\perp }\mbox{R}_h}
en la versión OSS.
De acuerdo con la formulación anterior, el problema discreto con linealización de Picard, integración temporal BDF1, la estabilización ASGS consiste en:
Hallar Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \mbox{u}_h^{n+1}\mbox{ }y\mbox{ }p_h^{n+1}}
tales que
|
( 41) |
|
( 42) |
con los términos con segundas derivadas de las funciones de elementos finitos nulos para elementos lineales.
En el método OSS, en el caso de considerar distintas densidades, como por ejemplo en problemas con interfase entre fluidos inmiscibles, el residuo en puntos de integración de lados opuestos a la interfase varía fuertemente y en proporción a la densidad. En tales casos puede adoptarse una proyección modificada [45] :
|
( 43) |
Sustituyendo Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \mbox{R}_h^{n+1}}
y la proyección modificada para el método OSS queda que:
El problema discretizado con linealización de Picard, integración temporal BDF1 y estabilización OSS consiste en hallar Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \mbox{u}_h^{n+1}\mbox{ }y\mbox{ }p_h^{n+1}}
tales que:
|
( 44) |
|
( 45) |
El método split -OSS [15] and [45] separa la proyección del término convectivo y de la presión en 2 proyecciones, lo que permite una convergencia a la solución más rápida que en los métodos ASGS y OSS. Esto lo hace ventajoso para la resolución de problemas en 3 D. Este método resulta en:
|
( 46) |
|
( 47) |
para las iteraciones i = 1,2… hasta la convergencia, es decir, hasta que un +1,i −1 ≈ un +1,i y Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle p_h^{n+1,i}\approx p_h^{n+1,i-1}}
en la norma elegida.
El parámetro τ1 de las ecuaciones (39) y (40) se elige con el fin de obtener esquemas numéricos estables y velocidades de convergencia óptimas (ver [46] ). Este parámetro se calcula para cada elemento Ωe . Para τ1 se toma:
|
( 48) |
donde h es la longitud del elemento Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle e\mbox{ }\mbox{y}\mbox{ }\left|\mbox{u}^e\right|}
es la norma de la velocidad en el elemento e , c1 y c2 son coeficientes a elegir, μ y ρ son la viscosidad dinámica y densidad del fluido, respectivamente.
Se recomiendan los valores de c1 = 1 y c2 = 2 [46] .
Se presenta a continuación la versión matricial para los métodos ASGS, OSS y split -OSS. Se utiliza linealización de Picard y discretización en tiempo BDF1.
En la versión matricial para el método ASGS la proyección Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle P_h\left(\mbox{u}_h^{n+1}\cdot \nabla \mbox{u}_h^{n+1}+\frac{1}{\rho }\nabla p_h^{n+1}\right)}
se trata con un ciclo iterativo al igual que para la linealización del término convectivo. Se definen:
|
( 49) |
En lo que sigue se usa la notación compacta Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle (a,b)=\int_{\Omega }a\cdot b\mbox{ }\mbox{d}\Omega }
. En esta notación, la proyección de la ecuación (49) es la solución de:
|
( 50) |
para todo Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \mbox{v}_h^\ast \in \mbox{V}_h^\ast }
, donde Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \mbox{V}_h^\ast } es el espacio Vh ampliado con los vectores de funciones continuas asociados a los nodos del contorno.
El sistema algebraico resultante es:
|
( 51) |
|
( 52) |
donde U y P son los vectores de las incógnitas nodales para la velocidad u y la presión p , respectivamente. F es el vector de fuerzas nodales.
La ecuación (51) corresponde a la ecuación (41). La ecuación (52) corresponde a la ecuación (42).
Si se denotan los índices nodales a , b , los índices espaciales con i , j , la función de forma de los nodos a por Na y la función de forma de los nodos b por Nb , entonces las matrices de las ecuaciones anteriores, las cuales son válidas para los métodos restantes, son:
|
|
|
|
|
|
|
( 53) |
donde δij es la delta de Kronecker.
En la versión matricial para el método OSS la proyección Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle P_k\left(\mbox{u}_k^{n+1}\cdot \nabla \mbox{u}_k^{n+1}+\frac{1}{\rho }\nabla p_h^{n+1}\right)}
también se trata con un ciclo iterativo, al igual que para la linealización del término convectivo. Se definen:
|
( 54) |
En notación compacta, la proyección de la ecuación (54) es la solución de:
|
( 55) |
para todo Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \mbox{v}_h^\ast \in \mbox{V}_h^\ast }
, donde Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \mbox{V}_h^\ast } es el espacio Vh ampliado con los vectores de funciones continuas asociados a los nodos del contorno.
El sistema algebraico resultante es:
|
( 56) |
|
( 57) |
|
( 58) |
donde U , P y Y son los vectores de las incógnitas nodales para la velocidad u , la presión p y la proyección y , respectivamente. F es el vector de fuerzas nodales. Solo se consideran las fuerzas nodales del residuo en la proyección de los términos respectivos.
La ecuación (56) corresponde a la ecuación (44). La ecuación (57) corresponde a la ecuación (45). La ecuación (58) es la proyección del residuo en la ecuación de conservación de momentum .
Las matrices de las ecuaciones (56) a (58) no definidas anteriormente son:
|
|
|
|
( 59) |
Se presenta a continuación la versión matricial para el método split -OSS, algo más compleja que los métodos ASGS y OSS. Las proyecciones Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle P_h\left(\mbox{u}_h^{n+1}\cdot \nabla \mbox{u}_h^{n+1}\right)\mbox{ }y\mbox{ }P_h\left(\frac{1}{\rho }\nabla p_h^{n+1}\right)}
también se tratan con un ciclo iterativo al igual que para la linealización del término convectivo. Se definen:
|
( 60) |
y
|
( 61) |
En la notación compacta, las proyecciones de las ecuaciones (60) y (61) son la solución de:
|
( 62) |
|
( 63) |
para todo Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \mbox{v}_h^\ast \in \mbox{V}_h^\ast }
, donde Failed to parse (MathML with SVG or PNG fallback (recommended for modern browsers and accessibility tools): Invalid response ("Math extension cannot connect to Restbase.") from server "https://test.scipedia.com:8081/localhost/v1/":): {\textstyle \mbox{V}_h^\ast } es el espacio Vh ampliado con los vectores de funciones continuas asociados a los nodos del contorno.
El sistema algebraico resultante es:
|
( 64) |
|
( 65) |
|
( 66) |
|
( 67) |
donde U , P , Y y Z son los vectores de las incógnitas nodales para la velocidad u , la presión p y las proyecciones y y z , respectivamente. F es el vector de fuerzas nodales.
La ecuación (64) corresponde a la ecuación (46). La ecuación (65) corresponde a la ecuación (47). Las ecuaciones (66) y (67) son las proyecciones de los residuos de la ecuación de conservación de momentum y la ecuación de incompresibilidad, respectivamente.
La matriz de las ecuaciones (64) a (67) no definida anteriormente es:
|
( 68) |
En este trabajo se presentan el modelo continuo y la correspondiente formulación discreta para la resolución de las ecuaciones de Navier-Stokes para los flujos viscoplásticos de Bingham y de Herschel-Bulkley usando elementos finitos mixtos estabilizados con interpolación lineal tanto para la velocidad como para la presión.
Se tratan en detalle los fluidos viscoplásticos ideales y regularizados de Bingham y de Herschel-Bulkley. Se proponen, asimismo, nuevos modelos viscoplásticos regularizados para el fluido de Bingham y de Herschel-Bulkley.
Posteriormente, se presenta el modelo discreto de elementos finitos estabilizados para flujos confinados. Los métodos de estabilización usados son el método de estabilización de subescalas algebraicas (ASGS), el método de subescalas ortogonales (OSS) y el método de subescalas ortogonales desacopladas (split -OSS). El modelo discreto se ha extendido a los modelos viscoplásticos de Bingham y de Herschel-Bulkley.
Los autores agradecen el apoyo del Prof. R. Codina, mostrado a través de múltiples sugerencias y discusiones fructíferas, y a la Universidad de los Andes, en Venezuela, por su apoyo económico en el desarrollo de esta investigación, dentro de su programa de becarios en el exterior.
Published on 01/01/2015
Volume 18, Issue 1, 2016
Licence: CC BY-NC-SA license
Are you one of the authors of this document?