Transcript Document

Computing Waves in the Face
of Uncertainty
E. Bruce Pitman
Department of Mathematics
University at Buffalo
[email protected]
Nonlinearity and Randomness
in Complex Systems
1
Part of a large project investigating
geophysical mass flows

Interdisciplinary research project funded by NSF (ITR and EAR)

UB departments/people involved:

Mechanical engineering: A Patra, A Bauer, T Kesavadas, C
Bloebaum, A. Paliwal, K. Dalbey, N. Subramaniam, P. Nair, V.
Kalivarappu, A. Vaze, A. Chanda

Mathematics: E.B. Pitman, C Nichita, L. Le

Geology: M Sheridan, M Bursik, B.Yu, B. Rupp, A. Stinton, A. Webb,
B. Burkett

Geography (National Center for Geographic Information and Analysis): C
Renschler, L. Namikawa, A. Sorokine, G. Sinha

Center for Computational.Research M Jones, M. L. Green

Iowa State University E Winer
2
Mt. St. Helens, USA
3
Volcan Colima, Mexico
4
Atenquique, Mexico 1955
5
Atenquique, Mexico 1955
6
San Bernardino Mountain: Waterman Canyon
7
Guinsaugon. Phillipines, 02/16/06
Heavy rain sent a
torrent of earth, mud
and rocks down on the
village of Guinsaugon.
Phillipines, 02/16/06.
A relief official says
1,800 people are
feared dead.
8
Ruapehu, New Zealand
9
Pico de Orizaba, Mexico
Ballistic particle Simulations of pyroclastic flows and hazard map
at Pico de Orizaba -- hazard maps by Sheridan et. al.
10
“Hazard map” based on flow simulations
and input uncertainty characterizations
Regions for which
probability of flow > 1m for
initial volumes ranging from
5000 m3 to 108 m3 -- flow
volume distribution from
historical data
11
Introduction

Geophysical flows e.g. rock falls, debris flows, avalanches, volcanic
lava flows may have devastating consequences for the human
population

Need “what if …?” simulation tool to estimate hazards for formulating
public safety measures

We have developed TITAN2D

Simulate flows on natural terrain,

Be robust, numerically accurate and run efficiently on a large
variety of serial and parallel machines,

Quantify the effect of uncertain inputs

Have good visualization capabilities.
12
Goals of this talk
 Basic

mathematical modeling
Will not address extensions such as erosion, two
phase flows, that are important in the field
 Uncertainty
Quantification
Hyperbolic PDE system – poses special
difficulties for uncertainty computations
 Ultimate aim is Hazard Maps

13
Modeling
Savage , Hutter, Iverson, Denlinger, Gray, Pitman, …
Nonlinearity and Randomness
in Complex Systems
14
Modeling
Many models – complex physics is still not
perfectly represented !
Savage-Hutter Model
Iverson-Denlinger mixture theory Model
Pitman-Le Two-phase model
Debris Flows are hazardous mixture of soil,
rocks, clasts with interstitial fluid present
15
Micromechanics and Macromechanics

Characteristic length scales (from mm to Km)

e.g. for Mount St. Helens (mudflow –1985)
Runout distance  31,000 m
 Descent height  2,150 m
 Flow length(L)  100-2,000
 Flow thickness(H)  1-10 m
 Mean diameter of sediment material 0.001-10 m

(data from Iverson 1995, Iverson & Denlinger 2001)
16
Model Topography and Equations(2D)
z  s( x, y, t )
Upper free surface
Fs(x,t) = s(x,y,t) – z = 0,
flowing mass
ground
h  s b
z  b ( x, y )
Basal material surface
Fb(x,t) = b(x,y) – z = 0
Kinematic BC:
at Fs (x, t )  0 :  t Fs  v  Fs  0
at Fb (x, t )  0 :  t Fb  v  Fb  es
Iverson and Denlinger JGR, 2001; Pitman et. al. Phys. Fluids, 2003; Patra et. al, JVGR, 2005
17
Model System-Basic Equations
Solid Phase Only
The conservation laws for a continuum incompressible medium are:
 u  0
 t ρ 0u     ρ 0u  u     T  ρ 0 g
stress-strain rate relationship derived from Coulomb theory
[Aside: this system of equations is ill-posed (Schaeffer 1987)]
Boundary conditions for stress:
at F s (x, t )  0 :
Tsn s  0
at F b (x, t )  0 :
r
u
Tbn b  n b n b  Tbn b  r tan  n b  Tbn b
u




: basal friction angle
18
Model System-Scaling
Scaling variables are chosen to reflect the shallowness
of the geophysical mass
x, y, z   Lx* , Ly * , Hz * , h  Hh *

(v x , v y )  gL v*x , v*y

H / L    1
 g *
 t , T  ρgH T*
t  

L


L – characteristic length in the downstream and cross-stream
directions (Ox,Oy)
H – characteristic length in normal direction to the flow (Oz)
Drop (most) terms of O()
19
Model System-Depth Average Theory
Depth average
s
s
s
1
1
1


u
dz
,
T
dz
,
ρu dz



hb
hb
hb
where
is the avalanche thickness
h( x, y, t )  s( x, y, t )  b( x, y)
z – dimension is removed from the problem - e.g. for
h (hvx ) (hv y )
the continuity equation:


e
t
x
y
s
where vx , and v y are the averaged lateral velocities defined as:
s
s
hvx   vx dz,
hv y   v y dz,
b
b
20
Modeling of Granular Stresses
Earth pressure coefficient is employed to relate
normal stresses
Tsxx  kapTszz
Shear stresses assumed proportional to normal
stresses
Tsxy
 vx 
 sin int kapTszz
  sgn 
 y 
Hydraulic assumption in normal direction
Tszz  g z h
21
Model System – 2D
Depth averaging and scaling: Hyperbolic System of
balance laws
h hv x hv y


 es
t
x
y
continuity
x momentum
2
2
hv x  (hv x  .5  k ap g z h ) hv y v x



t
x
y
 g x h  v x es 
1


 v 
1
hg z
  g z   v x2   h  tan bed   sgn  x   hk ap
sin int 
2
2


y

y
vx  v y 


x

vx
2
3
1. Gravitational driving force
2. Resisting force due to Coulomb friction at the base
3. Intergranular Coulomb force due to velocity gradients normal to the
direction of flow
22
Uncertainty
Dalbey, Patra
Nonlinearity and Randomness
in Complex Systems
23
Modeling and Uncertainty
“Why prediction of grain behavior is difficult in geophysical granular
systems””
“…there is no universal constitutive description of this
phenomenon as there is for hydraulics”
the variability of granular agglomerations is so large that
fundamental physics is not capable of accurately describing
the system and its variations
P. Haff (Powders and Grains ’97)
24
Uncertainty in Outputs of Simulations of
Geophysical Mass Flows
Model Uncertainty
Model Formulation: Assumptions and Simplifications
Model Evaluation: Numerical Approximation,
Solution strategies – error estimation
Data Uncertainty
propagation of input data uncertainty
25
Modeling Uncertainty
• Sources of Input Data Uncertainty
Initial conditions – flow volume and position
Bed and internal friction parameters
Terrain errors
Erosion and two phase model parameters
26
INPUT UNCERTAINTY PROPAGATION
Model inputs – material, loading and boundary data are always
uncertain
range of data and its distributions may be estimated
propagate input range and distribution to an output range and
distributions
e.g. maximum strain, maximum excursion
How does uncertain input produce a solution distribution?
27
Effect of different initial volumes
Left – block and Ash flow on Colima, V =1.5 x 105 m3
Right – same flow -- V = 8 x105 m3
28
Effect of initial position, friction angles
Figure shows output of
simulation from TITAN2D –
A) initial pile location, C) and D)
used different friction angles,
and, F) used a perturbed starting
location
29
Comparison of Models San Bernardino
Single phase model – low basal friction 4 deg!
Single phase model – water with
frictional dissipation term!
30
Comparison of Models
50% solid fraction
70% solid fraction
31
Quantifying Uncertainty -- Approach
Methods
• Monte Carlo (MC)
• Latin Hypercube Sampling (LHS)
}
• Polynomial Chaos (PC)
• Non Intrusive Spectral Projection (NISP)
―Polynomial Chaos Quadrature (PCQ)
Random sampling based
}
Functional
Approximation
• Stochastic Collocation
32
Quantifying Uncertainty -- MC Approach
• Monte Carlo (MC): random sampling of input pdf
• Moments can be computed by running averages e.g. mean
and standard deviation is given by:
1
 U ( ) 
N MC
N MC
U ( )
i 1
i
  U   U 
2
Central Limit Theorem :
2


N MC
Computationally expensive.
Estimated computational time for 10-3
error in sample TITAN calculation on
64 processors ~ 217 days
33
Latin Hypercube Sampling -- MMC
McKay 1979, Stein 1987, …
1.
2.
For each random direction (random variable or input), divide that direction
into Nbin bins of equal probability;
Select one random value in each bin;
3.
Divide each bin into 2 bins of equal probability; the random value chosen
above lies in one of these sub-bins;
4.
Select a random value in each sub-bin without one;
5.
Repeat steps 3 and 4 until desired level of accuracy is obtained.
34
Functional Approximations
In these approaches we attempt to compute an
approximation of the output pdf based on
functional approximations of the input pdf
Prototypical method of this is the Karhunen
Loeve expansion
35
Quantifying Uncertainty -- Approach
Wiener ’34, Xiu and Karniadakis’02
• Polynomial Chaos (PC):
approximate pdf as the
truncated sum of infinite
number of orthogonal
polynomials yi
y
 g ( y (t );  ( ))
t
y (t )   yi (t )yi ( )
   jyj ( )
• Multiply by ym and integrate to use orthogonality
37
Chaos solver
38
PC for Burger’s equation
Let = kyk U= Ui yI
@
Ui
@
t
i
= à Uk
@
Ul
k @
x
i=1..n k=1..n
l
+ ÷k
@2Ul
k @x 2
l
Coupled across all
Equations m=1..n
Multiply by ψm and integrate
@Um R 2
m dø =
@t
@Ul
Uk @x
R
k l
+
m dø
@2Ul
÷k @x 2
R
k l
m dø
39
Polynomial Chaos Quadrature
Instead of Galerkin projection, integrate by
quadrature weights
Analogy with


Non-Intrusive Spectral Projection
Stochastic Collocation
Leads to a method that has the simplicity of MC
sampling and cost of PC
Can directly compute all moment integrals
Efficiency degrades for large number of random
variables
40
NISP
Replace integration with quadrature and interchange
order of integration of time and stochastic dimension
41
Polynomial Chaos Quadrature (PCQ)
42
Quantifying Uncertainty -- Approach
PCQ: a simple deterministic sampling method with sample points
chosen based on an understanding of PC and quadrature rules
 makes PC computationally
feasible for non-
linear non-polynomial forms
 easy to implement and parallelize; statistics
obtained directly
can use random variables from multiple
distributions simultaneously
 difficult
to find sample points for
very high orders
“curse of dimensionality” -samples required grows
exponentially as a function of
number of RV
44
Quantification of Uncertainty
Test Problem
Application to flow at Volcan Colima
Starting location, and,
Initial volume
are assumed to be random variables distributed according to
assumption, or available data
45
Test Problem
Burgers equation
Figure shows statistics
of time required to reach
steady state for randomly
positioned shock in
initial condition; PCQ
converges much faster
than Monte Carlo
46
Quantifying Uncertainty
Starting location
Gaussian with std.
deviation of 150m
Mean Flow
Flow from starting
locations 3 std. dev
away
Mean Flow
Flow from starting
locations 3 std. dev
away
48
Application to Volcan Colima
Initial volume 
uniformly
distributed from
1.57x106
to 1.57x107
Mean and standard
deviation
of flow spread
computed with
MC and PCQ
Monte Carlo
PCQ
49
Mean Flow for Volcan Colima for initial volume uncertainty
50
Mean+3std. dev for Volcan Colima -- initial volume uncertainty
51
“Hazard Map” for Volcan Colima
Probability of flow
Exceeding 1m for
Initial volume ranging
From 5000 to 108 m3
And basal friction from
28 to 35 deg
52
Conclusions
PCQ is an attractive methodology for
determining the solution distribution as a
consequence of uncertainty
Find full pdf
Curse of dimensionality still strikes
MC, LH, NISP, Point Estimate methods, PCQ –
which to use depends on the problem at hand
54
Conclusions
How to handle uncertainty in terrain? In the models?
More work to integrate PCQ into output functionals that
prove valuable
All developed software is available free and open
source from www.gmfg.buffalo.edu
Software can be accessed on the Computational Grid
(DOE Open Science Grid) at http://grid.ccr.buffalo.edu
55