NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF - - PowerPoint PPT Presentation

NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF FLUIDS

Pedro Lima

CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBOA PORTUGAL

October 11, 2011

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 1 / 28

Joint work with: Luisa Morgado Department of Mathematics, Universidade de Tras-os-Montes e Alto Douro, Portugal

• G. Hastermann and E. Weinmuller

Institute for Analysis and Scientific Computing, Vienna University of Technology, Austria

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 2 / 28

Outline of the talk

1 Introduction 2 Existence of Solution 3 The singularities of the problem and the

associated one-parameter families of solutions

4 Shooting method based on asymptotic

expansions

5 Collocation method 6 Numerical results 7 Conclusions and future work Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 3 / 28

INTRODUCTION Physical interpretation:

The behavior of mixtures of fluids (for example: liquid-gas) is described by the Cahn-Hillard theory. Free volume energy: E(ρ, |∇ρ|2) = E0(ρ) + γ 2 |∇ρ|2, γ > 0, where ρ- density of the fluid. E0(ρ) - classical volume free energy γ - surface tension coefficient (independent from |∇ρ|).

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 4 / 28

Generalized Model

If we allow that the surface tension depends on ∇ρ, the free volume energy takes the form E(ρ, |∇ρ|) = E0(ρ) + c p |∇ρ|p, γ > 0, p > 1; in this case we obtain the following PDE: c div(|∇ρ|p−2∇ρ) = µ(ρ) − µ0; The operator in the left-hand side is the p-laplacian, where p > 1 (if p = 2 we obtain the classical laplacian). In the case of spherical bubbles, we obtain the radial ODE: r1−N rN−1|ρ′(r)|p−2ρ′(r) ′ = fp(ρ), (0 < r < ∞), where fp is a function with three real roots, whose specific form depends

• n p.

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 5 / 28

Right-Hand Side

In the classical laplacian case (p = 2), f2 is a third degree polynomial f2(ρ) = 4λ2(ρ − ξ)(ρ + 1)ρ, where ξ is a real parameter; In the degenerate laplacian case (p = 2), fp has the form fp(ρ) = 2pλ2(ρ − ξ)(ρ + 1)ρ|ρ − ξ|α|ρ + 1|α, where α = 0 in the case p ≤ 2; for p > 2 the value of α will be discussed later.

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 6 / 28

Boundary Conditions

• lim

r→0+ ρ(r)

• < ∞,

lim

r→0+ rρ′(r) = 0,

lim

r→∞ ρ(r) = ξ,

lim

r→∞ ρ′(r) = 0.

In the bubble case (if ξ > 0) , we search for a strictly increasing solution. In the droplet case (if ξ < −1), we search for a strictly decreasing solution.

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 7 / 28

References

F.dell’Isola, H.Gouin and P.Seppecher, ”Radius and Surface Tension of Microscopic Bubbles by Second Gradient Theory”, C.R.Acad. Sci. Paris, 320(Serie IIb), 211–216 (1995). F.dell’Isola, H.Gouin and G.Rotoli, ”Nucleation of Spherical Shell–Like Interfaces by Second Gradient Theory: Numerical Simulations”, Eur. J.

• Mech. B / Fluids 15, 545–568 (1996).

H.Gouin and G.Rotoli, ”An Analytical Approximation of Density Profile and Surface Tension of Microscopic Bubbles for Van der Waals Fluids”, Mechanics Research Communications 24, 255–260 (1997).

• N. Kim, L. Consiglieri and J.F.Rodrigues, On non-newtonian

incompressible fluids with phase transitions, Mathematical Methods in Applied Sciences, 29 1523–1541 .

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 8 / 28

Existence of Solution

Existence results for problems of this type can be found in:

• F. Gazzola, J. Serrin and M. Tang, Existence of ground states and free

boundary problems for quasilinear elliptic operators, Adv. Diff. Equ., 5 (2000) 1-30. From this work, it follows that, when p ≤ 2, for 0 < ξ < 1, the considered problem (choosing α = 0) has a bubble-type solution. For p > 2,existence of solution is guaranteed only if we choose α = p − 2 in the right hand side function.

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 9 / 28

The Singularity at r = 0

Initial conditions:

lim

r→0+ ρ(r) = ρ0

lim

r→0+ rρ′(r) = 0.

(1)

We assume that in the neighborhood of r = 0 the solution can be represented as ρ(r) = ρ0 + Crk(1 + o(1)), as r → 0+, (2)

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 10 / 28

Asymptotic approximation close to the origin

Proposition 3.1. Let N > 1 and p > 1. For each ρ0, the considered singular Cauhy problem has, in the neighborhood of r = 0, a unique holomorphic solution that can be represented by ρ (x, ρ0) = ρ0 + p − 1 p fp(ρ0) N

• 1

p−1

r

p p−1

• 1 + y1r

p p−1 + o

• x

p p−1

• ,

(3) where y1 can be determined analytically.

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 11 / 28

Singularity at Infinity

As r → ∞ we introduce the variable substitution ρ(r) = ξ + r

1−N (p−1)2 z(r).

(4) In the new variable z we obtain an asymptotically autonomous equation. In order to analyse the asymptotic behavior of the solutions, we can consider the autonomous equation: (p − 1)z′′

∞(r) = 2pλ2 z∞(r)p−1ξp−1(ξ + 1)

z′

∞(r)p−2

. (5) We search for a solution of (5) in the form z∞(r) = c exp(τr), (6) where c and τ are constants.

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 12 / 28

Asymptotic Expansion at Infinity

Subsituting in the equation, we obtain: z∞(r) = c exp  −

p

• 2pλ2 (1 + ξ)ξp−1

p − 1 r   . (7) Then, the solution of the non-autonomous equation can be expressed in the form of a Lyapunov series: z(r) =

∞

∑

k=1

bkCk(r)e−τkr, (8) where the functions Ck can be determined by solving a set of linear ODEs.

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 13 / 28

Asymptotic Expansion at Infinity

We have obtained an asypmptotic expression of the solutions which satisfy the condition lim

r→∞ z(r) = lim r→∞ z′(r) = 0.

(9) In the old variable ρ, we obtain the asymptotic expression ρ(r) = ξ − bC1(r)rae−τr(1 + o(1)), r → ∞. (10) We must compute the value of b for which the solution satisfies the prescribed boundary condition close to 0.

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 14 / 28

Numerical Approximation - Shooting Method

R - bubble radius (ρ(R) = 0). r0 - initial approximation of R. First Auxiliary Problem ρ−(r) - monotone solution on [δ, r0], which satisfies the boundary conditions ρ (δ) = ρ0 + p − 1 p fp(ρ0) N

• 1

p−1

δ

p p−1

• 1 + y1δ

p p−1

• ,

(11) ρ(r0) = 0. (12) Second Auxiliary Problem ρ+(r) - monotone solution on [r0, r∞], which satisfies the boundary conditions (12) and ρ(r∞) = ξ − bra

∞C1(r∞)e−τr∞.

(13)

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 15 / 28

Numerical Approximation - Shooting Method

For a given value of r0, each of the auxiliary problems can be solved

• numerically. Then we construct the global solution:

ρ(r) = ρ−(r), if δ ≤ r ≤ r0; ρ+(r), if r0 ≤ r ≤ r∞. (14) Let ∆(r0) = ρ′

+(r0) − ρ′ −(r0).

The true value of r0 is computed from the condition that ∆(r0) = 0.

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 16 / 28

Numerical Approximation - Collocation Method

The problem was also solved by a package of programms developped at the Vienna University of Technology. By means of a variable substitution the equation over the unbounded domain is reduced to a system of equations over the interval [0, 1]. The solution is approximated by a piecewise polynomial function which satisfies the differential equation at the collocation points. In order to deal with very high derivatives at the origin and close to R, in the case p = 2 a non-uniform mesh is used.

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 17 / 28

Variational Formulation- Energy Integral

The considered BVP can be considered as a variational problem. In this case, we consider the minimization of the energy integral: J := J(ρ) :=

∞

• p−1ρ′(r)p + fp(ρ(r))
• rN−1dr

(15) The convergence of this integral is a necessary condition for the solvability

• f the considered BVP.

J(ρ) was computed numericallly, along with the solution, and the obtained approximations are displayed in this section.

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 18 / 28

Numerical Results - Tables

Table: Values for p = 4

ξ J errJ ρ0 R 0.1 0.0000110 1.2430e-06

• 0.105224707

1.4526123 0.2 0.0002406 3.3812e-05

• 0.222387362

1.3167518 0.3 0.0016604 2.5849e-04

• 0.352716078

1.3022994 0.4 0.0074125 1.1950e-03

• 0.496734089

1.3559345 0.5 0.0272362 4.3152e-03

• 0.651936824

1.4827317 0.6 0.0936560 1.3972e-02

• 0.806770828

1.7227242 0.7 0.3324319 4.5475e-02

• 0.932061496

2.1824727 0.8 1.3849031 1.7231e-01

• 0.992219569

3.1738130

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 19 / 28

Numerical Results - Tables

Table: Values for p = 3.5

ξ J errJ ρ0 R 0.1 0.0000663 7.1726e-06

• 0.1253685

1.731116 0.2 0.0009179 1.2943e-04

• 0.2622131

1.534249 0.3 0.0048775 7.6815e-04

• 0.4104203

1.501633 0.4 0.0182754 2.9618e-03

• 0.5677667

1.557478 0.5 0.0592836 9.3047e-03

• 0.7265184

1.705310 0.6 0.1854613 2.7005e-02

• 0.8679368

1.991638 0.7 0.6077452 8.0602e-02

• 0.9629840

2.537943 0.8 2.3553203 2.8465e-01

• 0.9973361

3.720690

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 20 / 28

Numerical Results - Tables

Table: Values for p = 3

ξ J errJ ρ0 R 0.1 0.0004458 4.6772e-05

• 0.1553041

2.1128421 0.2 0.0038033 5.4274e-04

• 0.3196498

1.8217047 0.3 0.0153870 2.4759e-03

• 0.4899689

1.7632890 0.4 0.0483425 7.9119e-03

• 0.6589481

1.8249067 0.5 0.1387958 2.1982e-02

• 0.8117346

2.0065574 0.6 0.3950271 5.6167e-02

• 0.9260635

2.3617865 0.7 1.1925252 1.5376e-01

• 0.9852478

3.0364177 0.8 4.3028403 5.0486e-01

• 0.9994262

4.4585079

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 21 / 28

Numerical Results - Tables

Table: Values for p = 2.5

ξ J errJ ρ(0) R errsol 0.2 0.02033902 6.657e-03

• 0.5679509

2.6858283 1.404e-04 0.3 0.09596951 2.202e-02

• 0.7709138

2.582680 1.925e-04 0.4 0.27328747 5.659e-02

• 0.9033395

2.721913 4.524e-04 0.5 0.74289829 5.523e-02

• 0.9000092

2.433356 5.304e-04

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 22 / 28

Errors and Condition Numbers

Table: Convergence and condition numbers for (p, ξ) = (3, 0.5)

h errz1

• rd

ρ0 R cond

• rd

2−3 3.80449e−03 8.60

• 0.98958

2.203714 3.16e+06

• 2.03

2−4 9.82483e−06 3.03

• 0.98879

2.197152 1.29e+07

• 1.55

2−5 1.20483e−06 1.34

• 0.9887820

2.197105 3.78e+07

• 1.41

2−6 4.76089e−07 1.49

• 0.9887808

2.197106 1.00e+08

• 1.71

2−7 1.69122e−07 1.51

• 0.98878037

2.19710606 3.28e+08

• 2.60

2−8 5.93217e−08 1.51

• 0.98878020

2.19710605 1.99e+09

• 3.66

2−9 2.08948e−08

• 0.98878014

2.19710605 2.51e+10

• 3.22

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 23 / 28

Solution Graphics

2 4 6 8 10 12 −1 −0.5 0.5 p = 2 2 4 6 8 −1 −0.5 0.5 p = 3.5 2 4 6 8 10 12 −1 −0.5 0.5 p = 3 2 4 6 8 −1 −0.5 0.5 p = 4 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8

ξ

Figure: Solutions ρ(r) for r ∈ (0, ∞)

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 24 / 28

Bubble Radius

0.2 0.4 0.6 0.8 xi 1 2 3 4 5 6 7 R

Figure: The dependence of the bubble radius on ξ, for p = 2 (blue), p = 3 (red) and p = 4 (green)

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 25 / 28

Density at the center

0.2 0.4 0.6 0.8 xi 1.0 0.8 0.6 0.4 0.2 rho

Figure: Dependence of ρ0 on ξ, for p = 2 (blue), p = 3 (green) and p = 4 (red)

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 26 / 28

Conclusions and Future Work

By analysing the associate singular Cauchy problems we were able to describe the behavior of the solutions near the singularities. It was proved that a solution of the considered boundary value problem exists if and only if 0 < ξ < 1, for p > 1. Numerical algorithms based on the shooting method are simple and work efficiently for a large number of values of ξ . Numerical algorithms based on collocation methods provide higher accuracy, specially when ξ is close to 0 or to 1.

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 27 / 28

Conclusions and Future Work

The results obtained for the p-laplacian confirm that many of the properties of the original model can be extended to the general one. In particular, for each value of p there is a minimal bubble radius, which is attained for a certain value of ξ. The density of the gas at the bubble centre (ρ0) tends to −1 as ξ tends to 1 . Comparing with the case p = 2, the BVP with degenerate laplacian is harder to compute, due to the fact that the second and higher derivatives of the solution become unbounded at the

• rigin.

A significant drop of the convergence order of the numerical methods is observed, in comparison with the case of the classical laplacian. In the future we intend to complete the numerical investigation

• f the computational methods.

Pedro Lima (CENTRO DE MATEM´ ATICA E APLICAC ¸ ˜ OES INSTITUTO SUPERIOR T´ ECNICO UNIVERSIDADE T´ ECNICA DE LISBO NUMERICAL SOLUTION OF THE DENSITY PROFILE EQUATION FOR MIXTURES OF October 11, 2011 28 / 28