PARALLEL LARGE-SCALE SIMULATION OF VISCOELASTIC FLUID FLOW - - PowerPoint PPT Presentation

PARALLEL LARGE-SCALE SIMULATION OF VISCOELASTIC FLUID FLOW INSTABILITIES

Mehmet SAHIN Astronautical Engineering Department, 17th International Workshop on Numerical Faculty of Aeronautics and Astronautics, Methods for non-Newtonian Flows Istanbul Technical University, 34469, 25 March 2012 — Blois, FRANCE Maslak/Isatanbul, TURKEY

Two-dimensional inertial instability

Flow Instabilities in Newtonian and non-Newtonian Fluids

Two-dimensional inertial instability Three-dimensional inertial instability Three-dimensional viscoelastic instability

The Linear Stability (Normal Mode) Analysis for the Three-Dimensional Oldroyd-B Fluid Past Periodic Array of Cylinders with L=2.5R

L

Sahin & Wilson JNNFM (2008).

The Linear Stability (Normal Mode) Analysis for the Three-Dimensional Oldroyd-B Fluid Past Periodic Array of Cylinders with L=4.0R

The linear stability analysis was not conclusive for L=4.0R due to the classical High Weissenberg Number problem.

Can we use direct numerical simulations with the log conformation in

• rder to investigate viscoelastic fluid flow instabilities?

Unstructured Finite Volume Formulation

The governing equations of an incompressible Oldroyd-B fluid can be written in dimensionless form as follows: Integrating the differential equations over an arbitrary irregular control volumes leads to

Numerical Discretization

(a) Two-dimensional dual volume (b) Three-dimensional dual volume The side centered finite volume method was initially used by Hwang (1995) and Rida et al. (1997) for the solution of the incompressible Navier-Stokes equations on unstructured triangular

• meshes. The present arrangement of the primitive variables leads to a stable numerical scheme

and it does not require any ad-hoc modifications in order to enhance pressure-velocity-stress

• coupling. The most appealing feature of this primitive variable arrangement is the availability of

very efficient multigrid solvers. Sahin JNNFM (2011).

The discrete contribution from the right cell is given for the momentum equation along the x-axis. The time derivation: The convective term

Numerical Discretization (Continued…)

The convective term The pressure term The viscous term

The extra stress term The discretization of the constitutive equation for the Oldroyd-B fluid

Numerical Discretization (Continued…)

The gradient terms are calculated by the use of Green’s theorem In order to evaluate the face values of the extra-stress tensor, an upwind least square interpolation is employed.

Any component of the extra stress tensor can be extrapolated to control volume boundaries using a Taylor series expansion The neighbouring cell center values may be expressed as

Numerical Discretization (Continued…)

This over determined system may be solved in a least square sense. The discretization of the above governing equations leads to

A time splitting scheme decouples the calculation of the extra stress tensor from the evaluation of the velocity and pressure fields by solving a generalized Stokes problem.

Iterative Methods

However, the classical iterative methods (Richardson, Gauss-Seidel, Jacobi, etc) and the multilevel methods can not be applied directly because of the zero block in the saddle point problem. Then we will apply the two-level non-nested geometric multigrid preconditioner (Sahin JNNFM, 2011) to solve the modified Stokes systems.

^ ^ ^ ^

The Basic Idea of the Multigrid Method

smoothing Fine Grid restriction prolongation (interpolation) Fine Grid Smaller Coarse Grid

The basic idea of the multigrid method is to carry out iterations on a fine grid and then progressively transfer these flow field variables and residuals to a series of coarser grids. On the coarser grids, the low frequency errors become high frequency ones and they can be easily annihilated by simple explicit methods.

High Weissenberg Number Problem

The relation between the conformation tensor and the extra stress tensor is given by The constitutive equation for the Oldroyd-B fluid in terms of the conformation tensor is given by The conformation tensor is a quantity that describes the internal microstructure of polymer molecules in a continuum level. The conformation tensor is symmetric and positive definite. Unless special care is taken, the conformation tensor may lose this property at high Weissenberg numbers and the numerical solution will soon diverge. Log conformation method is proposed by R. Fattal and R. Kupferman, J. Non- NewtonianFluid Mech. 123 (2004).

Log Conformation

The constitutive equation for the Oldroyd-B fluid in terms of the conformation tensor It is possible to decompose the gradient of divergence free velocity field into nonsymmetric and tensors, and symmetric tensor. Use eigen decomposition theorem Then let The evolution equation for becomes Then the conformation tensor is positive definite.

Mesh Generation (GAMBIT, CUBIT, ...) Mesh Partition (METIS library) Parallel Unstructured Finite Volume Code

Parallelization and Efficiency

Linear Solver with Two-Level Multigrid Preconditioner Post Processing (Tecplot) Parallel Unstructured Finite Volume Code Kroylov subspace methods (PETSc library) Preconditoners (PETSc library)

Computing Resources for Parallel Calculations

• SGI Altix 3000 (1300MHz, Itanium 2) with 32 nodes
• National Center for High Performance Computing of Turkey,
• TUBITAK ULAKBIM, High Performance and Grid Computing Center

TEST CASE I: Two-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder

The computational coarse mesh (M1) with 35,313 quadrilateral element (DOF=283,508). The ratio of the channel height to the cylinder diameter is 2.

TEST CASE I: Two-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder

The mesh convergence of Txx with mesh renement on the cylinder surface and in the cylinder wake at We = 0.7 with β=0.59 for an Oldroyd-B fluid.

TEST CASE I: Two-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder

The comparison of Txx with the numerical results of Yurun et al. (1999), Hulsen et

• al. (2005) and Afonso et al. (2009) on the cylinder surface and in the cylinder wake

at We = 0.7 with β=0.59 for an Oldroyd-B uid.

TEST CASE I: Two-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder

The mesh convergence of Txx with mesh renement on the cylinder surface and in the cylinder wake at We = 0.8 with β=0.59 for an Oldroyd-B fluid.

At this point, we are not sure whether the extra stress along the center line in the wake should exhibit exponential unbounded growth with time to infinity or leads to a time-dependent solution for the present two- dimensional calculations.

TEST CASE I: Two-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder

The RMS convergence for the extra stress tensor at We = 0.9 with β=0.59 for an Oldroyd-B fluid past a confined cylinder (∆t = 0.01).

TEST CASE II: Three-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder in a Rectangular Channel

The computational mesh with 582,400 hexahedral element (DOF=9,397,972). The ratio of the channel height to the cylinder diameter is 2. The ratio of channel width to the channel height is 5.

TEST CASE II: Three-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder in a Rectangular Channel

The computed Txx contours at We=0.7 for the Oldroyd-B fluid past a confined circular cylinder in a rectangular channel. Contour levels are 0, 0.1, 2, 4 and 8.

TEST CASE II: Three-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder in a Rectangular Channel

The computed isobaric surfaces at We=0.7 for the Oldroyd-B fluid past a confined circular cylinder in a rectangular channel.

TEST CASE II: Three-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder in a Rectangular Channel

The computed Txx contours at We=2.0 for the Oldroyd-B fluid past a confined circular cylinder in a rectangular channel. Contour levels are 0, 0.1, 4 and 8.

TEST CASE II: Three-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder in a Rectangular Channel

The computed isobaric surfaces at We=2.0 for the Oldroyd-B fluid past a confined circular cylinder in a rectangular channel.

TEST CASE II: Three-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder in a Rectangular Channel

We=0.0 We=2.0 The computed streamtraces for the Oldroyd-B fluid past a confined circular cylinder in a rectangular channel.

TEST CASE II: Three-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder in a Rectangular Channel

Three-dimensional experimental viscoelastic fluid flow instability taken from Gareth McKinley's Non-Newtonian Fluid Dynamics Research Group. The spanwise wave number of the three-dimensional instability is approximately equal to R which is in accord with that of the computational results.

TEST CASE II: Three-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder in a Rectangular Channel

Corner vortex (Shiang et al. 2000) The computed streamtraces indicating a hourseshoe vortice at z=4.99 plane at We=2.0 for the Oldroyd-B fluid past a confined circular cylinder in a rectangular channel at We=2.0.

TEST CASE III: Three-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder in a Channel with a Spanwise Periodicity of W=2R

The computational mesh with 1,279,200 hexahedral elements (DOF=20,605,860). The ratio of the channel height to the cylinder diameter is 2. The ratio of channel width to the channel height is 0.5. The periodic boundary condition is applied in the spanwise direction. The calculation at We=2.0 requires approximately 2-3 days on the Karadeniz machine with 128 processors (Intel Xeon 5550) at UYBHM.

TEST CASE III: Three-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder in a Channel with a Spanwise Periodicity of W=2R

. The computed streamtraces at We=2.0. Color contours show the magnitude of the w-velocity component. IMPORTANT: The flow is not symmetric according to the y=0 plane! The particles left at y=0 plane at the upstream does not necessarily end up at y=0 plane at the downstream.

TEST CASE III: Three-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder in a Channel with a Spanwise Periodicity of W=2R

. The computed streamtraces at We=2.0. Color contours show the magnitude of the w-velocity component. IMPORTANT: The flow is not symmetric according to the y=0 plane! The particles left at y=0 plane at the upstream does not necessarily end up at y=0 plane at the downstream.

TEST CASE III: Three-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder in a Channel with a Spanwise Periodicity of W=2R

. The computed streamtraces on the cylinder surface (r=1.01R) at We=2.0. Color contours show the magnitude of the v-velocity component.

Separation lines

TEST CASE III: Three-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder in a Channel with a Spanwise Periodicity of W=2R

The computed Txx contours at We=2.0 for the Oldroyd-B fluid past a confined circular cylinder in a rectangular channel. Contour levels are 10, 20, 40, 80 and 160.

TEST CASE III: Three-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder in a Channel with a Spanwise Periodicity of W=2R

The computed isobaric surfaces at We=2.0 for the Oldroyd-B fluid past a confined circular cylinder in a rectangular channel. Extremely low values of pressure are

• bserved behind the cylinder (in particular, along the separation lines). The total

drag value is 133.79 (slightly higher than that of the Stokes flow) and the 2D value reported by Hulsen (JNNFM, 2005) is 135.84.

TEST CASE III: Three-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder in a Channel with a Spanwise Periodicity of W=2R

. [a] [b] [c] The computed streamtraces at x=1.3R plane with u-velocity contours [a], v-velocity contours [b] and w-velocity contours [c] at We=2.0.

TEST CASE III: Three-Dimensional Oldroyd-B Fluid Past a Confined Circular Cylinder in a Channel with a Spanwise Periodicity of W=2R

. The u-velocity component variation with x/R on the y=0 symmetry plane at We=2.0.

• A parallel unstructured finite volume technique has been developed

for parallel large-scale viscoelastic fluid flow computations.

• The present arrangement of the primitive variables leads to a stable numerical

scheme and it does not require any ad-hoc modifications in order to enhance the pressure-velocity-stress coupling.

• The time stepping algorithm used decouples the calculation of the extra stresses

from the evaluation of the velocity and pressure fields by solving a generalized

Conclusions

from the evaluation of the velocity and pressure fields by solving a generalized Stokes problem.

• The most appealing feature of present primitive variable arrangement is the

availability of very efficient multigrid solvers for the Stokes system.

• The log-conformation representation has been implemented in order improve the

limiting Weissenberg numbers.

• These are the first numerical results indicating the three-dimensional viscoelastic

instabilities seen in the experimental works of McKinley et al. (1993), Liu (1997) and Shiang et al. (2000) for an Oldroyd-B fluid past a confined cylinder in a channel at We=2.0.

• The present numerical calculations correctly predict the spanwise wave number.

The author gratefully acknowledge the use of the Chimera machine at the Faculty of Aeronautics and Astronautics at ITU, the computing resources provided by the National Center for High Performance Computing of Turkey (UYBHM) under grant number 10752009 and the computing facilities at TUBITAK ULAKBIM, High Performance and Grid Computing Center. The author would also like to thank PETSc team for their helpful suggestions during the implementation of the present two-level preconditoner.

Acknowledgements