Materials Process Design and Control Laboratory

STOCHASTIC MODELING OF MULTISCALE SYSTEMS NICHOLAS ZABARAS Materials Process Design and Control Laboratory Sibley School of Mechanical & Aerospace Engineering 188 Frank H T Rhodes Hall Cornell University Ithaca, NY 14853 Email: zabaras@cornell.edu URL: http://mpdc.mae.cornell.edu

OUTLINE • Motivation: coupling multiscaling and uncertainty analysis • Mathematical representation of uncertainty • Variational multiscale method (VMS) • Stochastic support method • Stochastic convection-diffusion equations • Computing PDFs of microstructures (Maximum Entropy) • Sparse grid collocation methods (Smolyak quadrature) • Natural convection on rough surfaces • Diffusion in stochastic heterogeneous media • Future research directions

NEED FOR UNCERTAINTY ANALYSIS • Uncertainty is everywhere From DOE From GE-AE website From Intel website From NIST Porous media Silicon wafer Aircraft engines Material process • Variation in properties, constitutive relations • Imprecise knowledge of governing physics, surroundings • Simulation based uncertainties (irreducible)

WHY UNCERTAINTY AND MULTISCALING ? • Uncertainties introduced across various length scales have a non-trivial interaction • Current sophistications – resolve macro uncertainties Micro Meso Macro • Imprecise boundary conditions • Initial perturbations • Use micro averaged models for resolving physical scales • Physical properties, structure follow a statistical description

UNCERTAINTY ANALYSIS TECHNIQUES • Monte-Carlo : Simple to implement, computationally expensive • Perturbation, Neumann expansions : Limited to small fluctuations, tedious for higher order statistics • Sensitivity analysis, method of moments : Probabilistic information is indirect, small fluctuations • Spectral stochastic uncertainty representation • Basis in probability and functional analysis • Can address second order stochastic processes • Can handle large fluctuations, derivations are general

RANDOM VARIABLES = FUNCTIONS ? • Math: Probability space (W, F, P) Sample space Probability measure Sigma-algebra • Random variable • : Random variable • A stochastic process is a random field with variations across space and time

SPECTRAL STOCHASTIC REPRESENTATION • A stochastic process = spatially, temporally varying random function CHOOSE APPROPRIATE BASIS FOR THE PROBABILITY SPACE GENERALIZED POLYNOMIAL CHAOS EXPANSION HYPERGEOMETRIC ASKEY POLYNOMIALS SUPPORT-SPACE REPRESENTATION PIECEWISE POLYNOMIALS (FE TYPE) KARHUNEN-LOÈVE EXPANSION SPECTRAL DECOMPOSITION SMOLYAK QUADRATURE, CUBATURE, LH COLLOCATION, MC (DELTA FUNCTIONS)

KARHUNEN-LOEVE EXPANSION ON random variables Mean function Stochastic process Deterministic functions • Deterministic functions ~ eigen-values , eigenvectors of the covariance function • Orthonormal random variables ~ type of stochastic process • In practice, we truncate (KL) to first N terms

GENERALIZED POLYNOMIAL CHAOS • Generalized polynomial chaos expansion is used to represent the stochastic output in terms of the input Stochastic input Askey polynomials in input Stochastic output Deterministic functions • Askey polynomials ~ type of input stochastic process • Usually, Hermite, Legendre, Jacobi etc.

SUPPORT-SPACE REPRESENTATION • Any function of the inputs, thus can be represented as a function defined over the support-space FINITE ELEMENT GRID REFINED IN HIGH-DENSITY REGIONS • SMOLYAK QUADRATURE • IMPORTANCE MONTE CARLO JOINT PDF OF A TWO RANDOM VARIABLE INPUT OUTPUT REPRESENTED ALONG SPECIAL COLLOCATION POINTS

VARIATIONAL MULTISCALE METHOD WITH ALGEBRAIC SUBGRID MODELLING • Application : deriving stabilized finite element formulations for advection dominant problems

VARIATIONAL MULTISCALE HYPOTHESIS EXACT SOLUTION COARSE SOLUTION INTRINSICALLY COUPLED SUBGRID SOLUTION H COARSE GRID RESOLUTION CANNOT CAPTURE FINE SCALE VARIATIONS THE FUNCTION SPACES FOR THE EXACT SOLUTION ALSO SHOW A SIMILAIR DECOMPOSITION • In the presence of uncertainty, the statistics of the solution are also coupled for the coarse and fine scales

VARIATIONAL MULTISCALE BASICS DERIVE THE WEAK FORMULATION FOR THE GOVERNING EQUATIONS PROJECT THE WEAK FORMULATION ON COARSE AND FINE SCALES SOLUTION FUNCTION SPACES ARE NOW STOCHASTIC FUNCTION SPACES COARSE WEAK FORM FINE (SUBGRID) WEAK FORM ALGEBRAIC SUBGRID MODELS COMPUTATIONAL SUBGRID MODELS REMOVE SUBGRID EFFECTS IN THE COARSE WEAK FORM USING STATIC CONDENSATION APPROXIMATE SUBGRID SOLUTION NEED TECHNIQUES TO SOLVE STOCHASTIC PDEs MODIFIED MULTISCALE COARSE WEAK FORM INCLUDING SUBGRID EFFECTS

VMS – ILLUSTRATION [NATURAL CONVECTION] Mass conservation Momentum conservation Energy conservation Constitutive laws FINAL COARSE FORMULATION VMS HYPOTHESIS OBTAIN ASGS DERIVE WEAK FORM DERIVE SUBGRID DEFINE PROBLEM

FINAL COARSE FORMULATION VMS HYPOTHESIS OBTAIN ASGS DERIVE WEAK FORM DERIVE SUBGRID DEFINE PROBLEM WEAK FORM OF EQUATIONS • Energy function space • Test • Trial • Energy equation – Find such that, for all , the following holds • VMS hypothesis: Exact solution = coarse scale solution + fine scale (subgrid) solution

FINAL COARSE FORMULATION VMS HYPOTHESIS OBTAIN ASGS DERIVE WEAK FORM DERIVE SUBGRID DEFINE PROBLEM ENERGY EQUATION – SCALE DECOMPOSITION • Energy equation – Find and such that, for all and , the following holds • Coarse scale variational formulation • Subgrid scale variational formulation • These equations can be re-written in the strong form with assumption on regularity as follows

FINAL COARSE FORMULATION VMS HYPOTHESIS OBTAIN ASGS DERIVE WEAK FORM DERIVE SUBGRID DEFINE PROBLEM ELEMENT FOURIER TRANSFORM • Element Fourier transform RANDOM FIELD DEFINED OVER THE DOMAIN RANDOM FIELD DEFINED IN WAVENUMBER SPACE SPATIAL MESH • Addressing spatial derivatives NEGLIGIBLE FOR LARGE WAVENUMBERS  SUBGRID APPROXIMATION OF DERIVATIVE

ASGS [ALGEBRAIC SUBGRID SCALE] MODEL STRONG FORM OF EQUATIONS FOR SUBGRID CHOOSE AND APPROPRIATE TIME INTEGRATION ALGORITHM TIME DISCRETIZED SUBGRID EQUATION TAKE ELEMENT FOURIER TRANSFORM

FINAL COARSE FORMULATION VMS HYPOTHESIS OBTAIN ASGS DERIVE WEAK FORM DERIVE SUBGRID DEFINE PROBLEM MODIFIED COARSE FORMULATION • Assume the solution obeys the following regularity conditions • By substituting ASGS model in the coarse scale weak form • A similar derivation ensues for stochastic Navier-Stokes

FLOW PAST A CIRCULAR CYLINDER NO-SLIP TRACTION FREE RANDOM UINLET NO-SLIP INLET VELOCITY ASSUMED TO BE A UNIFORM RANDOM VARIABLE KARHUNEN-LOEVE EXPANSION YIELD A SINGLE RANDOM VARAIBLE THUS, GENERALIZED POLYNOMIAL CHAOS  LEGENDRE POLYNOMIALS (ORDER 3 USED) • Investigations: Vortex shedding, wake characteristics

FULLY DEVELOPED VORTEX SHEDDING • Mean pressure • First LCE coefficient • Second LCE coefficient • Wake region in the mean pressure is diffusive in nature • Also, the vortices do not occur at regular intervals [Karniadakis J. Fluids. Engrg]

VELOCITIES AND FFT FFT YIELDS A MEAN SHEDDING FREQUENCY OF 0.162 FFT SHOWS A DIFFUSE BEHAVIOR IMPLYING CHANGING SHEDDING FREQUENCIES MEAN VELOCITY AT NEAR WAKE REGION EXHIBITS SUPERIMPOSED FREQUENCIES

VARIATIONAL MULTISCALE METHOD WITH EXPLICIT SUBGRID MODELLING FOR MULTISCALE DIFFUSION IN HETEROGENEOUS RANDOM MEDIA

FINAL COARSE FORMULATION VMS HYPOTHESIS AFFINE CORRECTION DERIVE WEAK FORM COARSE-TO-SUBGRID MAP DEFINE PROBLEM MODEL MULTISCALE HEAT EQUATION in on in THE DIFFUSION COEFFICIENT K IS HETEROGENEOUS AND POSSESSES RAPID RANDOM VARIATIONS IN SPACE • OTHER APPLICATIONS • DIFFUSION IN COMPOSITES • FUNCTIONALLY GRADED MATERIALS FLOW IN HETEROGENEOUS POROUS MEDIA  INHERENTLY STATISTICAL DIFFUSION IN MICROSTRUCTURES

FINAL COARSE FORMULATION VMS HYPOTHESIS AFFINE CORRECTION DERIVE WEAK FORM COARSE-TO-SUBGRID MAP DEFINE PROBLEM STOCHASTIC WEAK FORM such that, for all : Find • Weak formulation • VMS hypothesis Exact solution Subgrid solution Coarse solution

EXPLICIT SUBGRID MODELLING: IDEA DERIVE THE WEAK FORMULATION FOR THE GOVERNING EQUATIONS PROJECT THE WEAK FORMULATION ON COARSE AND FINE SCALES COARSE WEAK FORM FINE (SUBGRID) WEAK FORM COARSE-TO-SUBGRID MAP  EFFECT OF COARSE SOLUTION ON SUBGRID SOLUTION AFFINE CORRECTION  SUBGRID DYNAMICS THAT ARE INDEPENDENT OF THE COARSE SCALE LOCALIZATION, SOLUTION OF SUBGRID EQUATIONS NUMERICALLY FINAL COARSE WEAK FORMULATION THAT ACCOUNTS FOR THE SUBGRID SCALE EFFECTS

FINAL COARSE FORMULATION VMS HYPOTHESIS AFFINE CORRECTION DERIVE WEAK FORM COARSE-TO-SUBGRID MAP DEFINE PROBLEM SCALE PROJECTION OF WEAK FORM such that, for all Find and and • Projection of weak form on coarse scale • Projection of weak form on subgrid scale EXACT SUBGRID SOLUTION COARSE-TO-SUBGRID MAP SUBGRID AFFINE CORRECTION

SPLITTING THE SUBGRID SCALE WEAK FORM • Subgrid scale weak form • Coarse-to-subgrid map • Subgrid affine correction

NATURE OF MULTISCALE DYNAMICS ASSUMPTIONS: NUMERICAL ALGORITHM FOR SOLUTION OF THE MULTISCALE PDE COARSE TIME STEP SUBGRID TIME STEP 1 1 ũC ūC Coarse solution field at end of time step Coarse solution field at start of time step ûF

REPRESENTING COARSE SOLUTION ELEMENT COARSE MESH RANDOM FIELD DEFINED OVER THE ELEMENT FINITE ELEMENT PIECEWISE POLYNOMIAL REPRESENTATION USE GPCE TO REPRESENT THE RANDOM COEFFICIENTS • Given the coefficients , the coarse scale solution is completely defined

FINAL COARSE FORMULATION VMS HYPOTHESIS AFFINE CORRECTION DERIVE WEAK FORM COARSE-TO-SUBGRID MAP DEFINE PROBLEM COARSE-TO-SUBGRID MAP ELEMENT COARSE MESH ANY INFORMATION FROM COARSE TO SUBGRID SOLUTION CAN BE PASSED ONLY THROUGH BASIS FUNCTIONS THAT ACCOUNT FOR FINE SCALE EFFECTS INFORMATION FROM COARSE SCALE COARSE-TO-SUBGRID MAP

SOLVING FOR THE COARSE-TO-SUBGRID MAP START WITH THE WEAK FORM APPLY THE MODELS FOR COARSE SOLUTION AND THE C2S MAP AFTER SOME ASSUMPTIONS ON TIME STEPPING THIS IS DEFINED OVER EACH ELEMENT, IN EACH COARSE TIME STEP

BCs FOR THE COARSE-TO-SUBGRID MAP INTRODUCE A SUBSTITUTION CONSIDER AN ELEMENT

FINAL COARSE FORMULATION VMS HYPOTHESIS AFFINE CORRECTION DERIVE WEAK FORM COARSE-TO-SUBGRID MAP DEFINE PROBLEM SOLVING FOR SUBGRID AFFINE CORRECTION START WITH THE WEAK FORM • WHAT DOES AFFINE CORRECTION MODEL? • EFFECTS OF SOURCES ON SUBGRID SCALE • EFFECTS OF INITIAL CONDITIONS CONSIDER AN ELEMENT IN A DIFFUSIVE EQUATION, THE EFFECT OF INITIAL CONDITIONS DECAY WITH TIME. WE CHOOSE A CUT-OFF • To reduce cut-off effects and to increase efficiency, we can use the quasistatic subgrid assumption

FINAL COARSE FORMULATION VMS HYPOTHESIS AFFINE CORRECTION DERIVE WEAK FORM COARSE-TO-SUBGRID MAP DEFINE PROBLEM MODIFIED COARSE SCALE FORMULATION • We can substitute the subgrid results in the coarse scale variational formulation to obtain the following • We notice that the affine correction term appears as an anti-diffusive correction • Often, the last term involves computations at fine scale time steps and hence is ignored

DIFFUSION IN A RANDOM MICROSTRUCTURE • A MIXTURE MODEL IS USED AS AN EXAMPLE OF GENERATING A HETEROGENEOUS DISTRIBUTION OF CONDUCTIVITY • WE ASSUME THAT THE DIFFUSION COEFFICIENTS OF INDIVIDUAL CONSTITUENTS ARE NOT KNOWN EXACTLY THE INTENSITY OF THE GRAY-SCALE IMAGE IS MAPPED TO THE CONCENTRATIONS DARKEST DENOTES b PHASE LIGHTEST DENOTES a PHASE

RESULTS AT TIME = 0.05 FIRST ORDER GPCE COEFF MEAN SECOND ORDER GPCE COEFF RECONSTRUCTED FINE SCALE SOLUTION (VMS) FULLY RESOLVED GPCE SIMULATION

RESULTS AT TIME = 0.2 FIRST ORDER GPCE COEFF MEAN SECOND ORDER GPCE COEFF RECONSTRUCTED FINE SCALE SOLUTION (VMS) FULLY RESOLVED GPCE SIMULATION

HIGHER ORDER TERMS AT TIME = 0.2 FOURTH ORDER GPCE COEFF THIRD ORDER GPCE COEFF FIFTH ORDER GPCE COEFF RECONSTRUCTED FINE SCALE SOLUTION (VMS) FULLY RESOLVED GPCE SIMULATION

SUPPORT-SPACE – STOCHASTIC GALERKIN - Joint probability density function of the inputs - The input support-space denotes the regions where input joint PDF is strictly positive Triangulation of the support-space Any function can be represented as a piecewise polynomial on the triangulated support-space - Function to be approximated - Piecewise polynomial approximation over support-space L2 convergence – (mean-square) Error in approximation is penalized severely in high input joint PDF regions. We use importance based refinement of grid to avoid this h = mesh diameter for the support-space discretization q = Order of interpolation

IMPLEMENTATION OF SUPPORT-SPACE A stochastic process can be interpreted as a random variable at each spatial point Two-level grid approach • Support-space grid • Mesh dense in regions of high input joint PDF Spatial domain Spatial grid • There is finite element interpolation at both spatial and random levels • Each spatial location handles an underlying support-space grid • Highly OOP structure Element

CAPTURING UNSTABLE EQUILIBRIUM • Computational details – 1600 bilinear elements for spatial grid • Time of simulation – 1.5 non-dimensional units • Rayleigh number – uniformly distributed random variable between 1530 and 1870 (10% fluctuation about 1700) • Prandtl number – 6.95 • Time stepping – 0.002 non-dimensional units • Support-space grid – One-dimensional with ten linear elements Cold wall Insulated Insulated Hot wall • Simulation about the critical Rayleigh number – conduction below, convection above • Both GPCE and support-space methods are used separately for addressing the problem • Failure of Generalized polynomial chaos approach • Support-space method – evaluation and results against a deterministic simulation

FAILURE OF THE GPCE Y-vel X-vel Mean X- and Y- velocities determined by GPCE yields extremely low values !! (Gibbs effect) X-vel Y-vel X- and Y- velocities obtained from a deterministic simulation with Ra = 1870 (the upper limit)

PREDICTION BY SUPPORT-SPACE METHOD Y-vel X-vel Mean X- and Y- velocities determined by support-space method at a realization Ra=1870 X-vel Y-vel X- and Y- velocities obtained from a deterministic simulation with Ra = 1870 (the upper limit)

SPARSE GRID COLLOCATION • If the number of random inputs is large (dimension D ~ 10 or higher), the number of grid points to represent an output on the support-space mesh increases exponentially • GPCE for very high dimensions yields highly coupled equations and ill-conditioned systems (relative magnitude of coefficients can be drastically different) • Instead of relying on piecewise interpolation, series representations, can we choose collocation points that still ensure accurate interpolations of the output (solution)

SMOLYAK ALGORITHM: SPARSE GRIDS Full tensor product grid: 289 points Example of using sparse grids to build interpolating functions: Discontinuous functions Left to right: Improving interpolation depth For an error around 2x10-2: Required number of points using sparse grids 3300 Required number of points using full tensor products: 32769 Sparse Grid: 65 points Number of points required to construct interpolating functions reduces combinatorially. Reduction more significant as the number of dimensions increases

SMOLYAK ALGORITHM LET OUR BASIC 1D INTERPOLATION SCHEME BE SUMMARIZED AS IN MULTIPLE DIMENSIONS, THIS CAN BE WRITTEN AS TO REDUCE THE NUMBER OF SUPPORT NODES WHILE MAINTAINING ACCURACY WITHIN A LOGARITHMIC FACTOR, WE USE SMOLYAK METHOD IDEA IS TO CONSTRUCT AN EXPANDING SUBSPACE OF COLLOCATION POINTS THAT CAN REPRESENT PROGRESSIVELY HIGHER ORDER POLYNOMIALS IN MULTIPLE DIMENSIONS A FEW FAMOUS SPARSE QUADRATURE SCHEMES ARE AS FOLLOWS: CLENSHAW CURTIS SCHEME, MAXIMUM-NORM BASED SPARSE GRID AND CHEBYSHEV-GAUSS SCHEME

SMOLYAK ALGORITHM Extensively used in statistical mechanics Provides a way to construct interpolation functions based on minimal number of points Univariate interpolations to multivariate interpolations Uni-variate interpolation Multi-variate interpolation Smolyak interpolation Some degradation in accuracy Maximal reduction when the function is assumed to be smooth D = 10

SPARSE GRID COLLOCATION METHOD Solution Methodology PREPROCESSING Compute list of collocation points based on number of stochastic dimensions, N and level of interpolation, q Compute the weighted integrals of all the interpolations functions across the stochastic space (wi) Use any validated deterministic solution procedure. Completely non intrusive Solve the deterministic problem defined by each set of collocated points POSTPROCESSING Compute moments and other statistics with simple operations of the deterministic data at the collocated points and the preprocessed list of weights Std deviation of temperature: Natural convection

USING THE COLLOCATION METHOD FOR HIGHER DIMENSIONS • Flow through heterogeneous random media Alloy solidification, thermal insulation, petroleum prospecting Look at natural convection through a realistic sample of heterogeneous material Square cavity with free fluid in the middle part of the domain. The porosity of the material is taken from experimental data1 Left wall kept heated, right wall cooled Numerical solution procedure for the deterministic procedure is a fractional time stepping method 1. Reconstruction of random media using Monte Carlo methods, Manwat and Hilfer, Physical Review E. 59 (1999)

Materials Process Design and Control Laboratory