General-Purpose Software for Large-Scale Bifurcation Analysis

1
General-Purpose Software for Large-Scale
Bifurcation Analysis
Andy Salinger, Eric Phipps Computer Science
Research Institute, Sandia National
Laboratories Albuquerque, NM, USA Tipping
Points in Complex Flows - Numerical Methods for
Bifurcation Analysis of Large-Scale Systems from
31 Oct 2011 through 4 Nov 2011
Sandia is a multiprogram laboratory operated by
Sandia Corporation, a Lockheed Martin
Company,for the United States Department of
Energys National Nuclear Security
Administration under contract DE-AC04-94AL85000.
2
Trilinos Algorithms and Enabling Technologies
for Large-Scale Applications

• Two-level design
• Self-contained packages (50)
• Leveraged common tools.
• Version Control
• Build System
• Test Harness

Nonlinear, Transient Optimization Solvers
Linear Eigen Solvers
Discretizations
Scalable Linear Algebra
Geometry, Meshing Load Balancing
Framework, Tools Interfaces
http//trilinos.sandia.gov
3
Trilinos Package Summary
Objective Package(s)
Linear algebra objects Epetra, Tpetra
Krylov solvers AztecOO, Belos, Komplex
ILU-type preconditioners AztecOO, IFPACK, ShyLU
Multilevel preconditioners ML, CLAPS
Eigenvalue problems Anasazi
Block preconditioners Teko
Direct sparse linear solvers Amesos (MUMPS, SuperLU, UMFPack, )
Load Balancing, Graph Algs Zoltan, Isorropia
Continuation/Bifurcation LOCA
Nonlinear system solver NOX
Time Integrators/DAEs Rythmos
Optimization Moocho, Aristos, Dakota (via TriKota)
Automatic Diffferentiation Sacado
Parameter List, Timers, Memory Teuchos
Uncertainty Quantification Stokhos, Dakota (via Trikota)
Abstract interfaces Thyra, EpetraExt
Node Kernels on New Architectures Kokkos (Cuda, Threads, OpenMP)
4
LOCA Provides Analysis Capabilities to
Large-Scale Applications
Pseudo-Arclength Continuation
Linear eigen-analysis through Anasazi (Thornquist
Lehoucq)
Multi-parameter continuation through Multifario
(Henderson)
Periodic orbit tracking (experimental)
Pitchfork
Bifurcation location and continuation (turning
point, pitchfork, and Hopf)
Hopf
5
Why Do We Need Stability Analysis Tools for
Large-Scale Applications

• Several powerful continuation and bifurcation
  analysis tools are available
• AUTO (Doedel et al)
• CONTENT (Kuznetsov et al)
• MATCONT (Govaerts et al)
• PyDSTool (Guckenheimer et al)
• Large-scale applications have specific
  requirements
• Massively parallel distributed memory
  architectures
• Complicated parallel data structures and sparse
  matrices
• Application-tuned linear algebra
• Limited derivative capabilities
• Tools and algorithms are needed that
• Do not change matrix sparsity or increase memory
  requirements
• Agnostic to linear algebra and architecture
• Can be incorporated into existing simulation
  codes (i.e., libraries)

6
Basic Defining Equations in LOCA
7
Bifurcations Discovered Through Eigenvalue
Analysis and Spectral Transformations
Generalized Cayley Transformation

• Eigenvalues/vectors approximated via Block
  Krylov-Schur Iterations (Anasazi Thornquist
  Lehoucq)

• Analogies to time integration can be used to pick
  transformation parameters (Lehoucq and Salinger,
  2001, Burroughs et al, 2004).

8
NOX Object-Oriented Nonlinear Solverin
Trilinos Pawlowski et al

• Methods
• Line Search
• Trust Region

Solver Layer

• Searches
• Directions

Linear Algebra
Abstract Layer
Application Interface
Linear Algebra
Epetra LAPACK
User Defined Thyra
Concrete Implementation Layer

• User Interface
• Residual
• Jacobian

9
LOCA Builton and around NOX

1. Step
2. Solve
3. Analyze
4. Predict
5. Stop?

Stepper Layer

• Eigensolvers
• Spectral transformations

• Methods
• Line Search
• Trust Region

Solver Layer

• Searches
• Directions

Continuation
Augmented Equations Layers
Bifurcation

Linear Algebra
Abstract Layer
Application Interface

• Mix-and-match between
• Continuation methods
• Predictor modules
• Step-size control modules
• Bifurcation modules
• Nonlinear solvers
• Linear solvers/preconditioners

Linear Algebra
Concrete Implementation Layer
Epetra LAPACK
User Defined Thyra

• User Interface
• Residual
• Jacobian

• Update parameters
• Mass matrix

10
Block Elimination Algorithm for Turning Point
(fold) Tracking Uses 4 Solves

• Full Newton Algorithm

• Turning Point Bifurcation

• Block Elimination Algorithm

11
Modified Turning Point Bordering Algorithm
Solve 5 bordered systems of equations using QR
approach Then
12
Minimally Augmented Turning Point Formulation
Given and , let then There are
constants such that Standard
formulation Note for Newtons method 3
linear solves per Newton iteration (5 for
modified bordering)! For symmetric problems
reduces to 2 solves.
13
Flow Calculations Performed with Sandia CFD
codes and Trilinos solvers

• MPSalsa, Charon, Albany codes
• Incompressible Navier-Stokes
• Heat and Mass Transfer, Reactions, variable
  properties
• Unstructured Finite Element
• Galerkin/Least-Squares Q1Q1
• Analytic Jacobian matrix in distributed sparse
  storage
• Compute with AD
• Fully Coupled Newton Method
• GMRES with ILUT or MultI-Level Preconditioners

14
Frank-Kamenetskii Explosion Model(230K hex
elements, 1.1M unknowns, 128 cores)

• Scenario
• Continuous Stirred Tank Reactor
• Exothermic Cehmical Reaction
• Cooling at Walls
• The stirring breaks!
• Will natural convection prevent explosion?

Arc-length Continuation
15
Frank-Kamenetskii Explosion ModelTurning Point
Location

• ILU(k) fill factor 1
• ILU(k) overlap 2
• Max Krylov space 1000
• MS Bordering
• Minimally augmented

16
Frank-Kamenetskii Explosion ModelTurning Point
Continuation
Method Continuation Steps Failed Steps Nonlinear Iterations Linear Solves Linear Iterations Total Time (hrs)
Moore-Spence Mod. Bordering 49 5 290 2012 472968 11.0
Min. Augmented 38 4 214 810 154027 6.9
17
Analysis Beyond Simulation
LOCA
18
Natural Convection Instability in 8x1 and 8x1x1
Cavities
19
Stability Analysis of Impinging Jets Pawlowski,
Salinger, Shadid, Mountziaris (2005)
Region I
Region II
Region III
H
H
P
20
Rayleigh-Benard in 5x5x1 cavity with Bifurcation
Tracking
Pitchfork
Hopf
Pitchfork
Codimension 2 Bifurcation near Pr0.0434, Ra2106
Eigenvectors at Hopf
21
Hydromagnetic Rayleigh-Bernard ProblemPawlowski,
Shadid
B0
g

• Parameters
• Q B02 (Chandresekhar number)
• Ra (Rayleigh number)

• Buoyancy driven instability initiates flow at
  high Ra numbers.
• Increased values of Q delay the onset of flow.
• Domain 1x20

22
Resistive, Extended MHD Equations
Extended MHD Model in Residual Form
Involution
23
Hydro-Magnetic Rayleigh-Bernard Stability Direct
Determination of Linear Stability and Nonlinear
Equilibrium Solutions (Steady State Solves)
Leading Eigenvector at Bifurcation Point, Ra
1945.78, Q10

• 2 Direct-to-steady-state solves at a given Q
• Arnoldi method using Cayley transform to
  determine approximation to 2 eigenvalues with
  largest real part
• Simple linear interpolation to estimate Critical
  Ra

24
Design (Two-Parameter) Diagram
Vx
Ra

• No flow does not equal no-structure
  pressure and magnetic fields must adjust/balance
  to maintain equilibrium.
• LOCA can perform continuation of bifurcation

25
Critical Mode is different for various Q values
Leading mode is 26 cells

• Analytic solution is on an infinite domain with
  two bounding surfaces (top and bottom)
• Multiple modes exist, mostly differentiated by
  number of cells/wavelength.
• Therefore tracking the same eigenmode does not
  give the stability curve!!!
• Periodic BCs will not fix this issue.

4000
3000
Ra
2000
Leading mode is 20 cells
Q
Mode 20 Cells Q100, Ra4017
Mode 26 Cells Q100, Ra3757
26
Scaling studies
Cores Fine Mesh Level 0 Unkns. Intermed. Level 1 Unkns. Intermed. Level 2 Unkns. Coarse Level 3 Unkns. Newton Iters. Avg. No. Linear Its. / Newton Total Sim. Time (min.)
24,000 1.05 billion 23.3M .5M 11.2K 18 86 33
27
LOCA has impacted several application areas
without flow
Pattern Formation in Swift-Hohenberg Eqs
Avitabile, Sanstede
TeraHz Resonance in Quantum Tunneling Diode
Kelley
Super-Conductivity Transitions in Ginzburg-Landau
Schlomer, Vanroose
Phase Transitions in a Confined Fluid Frink
28
How do I attribute the successes that LOCA has
had?

1. Algorithmic research in large-scale bifurcations
2. Science demonstrations
3. LOCA is hooked up to Trilinos linear solvers

10 15 75
What broader lesson is there?
Build software in independent-yet-interoperable
components
29
Build Software in Independent-yet-Interoperable
Components
Continuation, Bifurcation, Stability Analysis
30
Build Software in Independent-yet-Interoperable
Components
Transient sensitivity analysis
31
Build Software in Independent-yet-Interoperable
Components
CVD Reactor Optimization
32
Agile Components
Analysis Tools (black-box)
Composite Physics
Particle Code Tools
MultiPhysics Coupling
Data Structures
Optimization
Solution Control
Neighbor Search / Sort
Parameter Studies
Utilities
System Models
UQ (non-invasive)
Input File Parser
PostProcessing
System UQ
VV, Calibration
Parameter List
Visualization
OUU, Reliability
I/O Management
Embedded Verification
Mesh Tools
Computational Steering
Feature Extraction
Memory Management
Mesh I/O
Data Reduction
Communicators
Analysis Tools (embedded)
Inline Meshing
Model Reduction
Runtime Compiler
Partitioning
Nonlinear Solver
MultiCore Parallelization Tools
Load Balancing
Time Integration
Mesh Database
Adaptivity
Continuation
Remeshing
Mesh Database
Sensitivity Analysis
Geometry Database
Grid Transfers
Stability Analysis
Software Quality
Mesh Quality
Solution Database
Constrained Solves
Version Control
Optimization
Local Fill
Regression Testing
UQ Solver
Build System
Discretizations
Physics Fill
Backups
Discretization Library
Linear Algebra
Element Level Fill
Mailing Lists
Field Manager
Data Structures
Material Models
Unit Testing
Iterative Solvers
Derivative Tools
Bug Tracking
Objective Function
Direct Solvers
UQ / PCE Propagation
Performance Testing
Constraints
Eigen Solver
Code Coverage
Error Estimates
Sensitivities
Preconditioners
Porting
MMS Source Terms
Matrix Partitioning
Derivatives
Web Pages
Adjoints
Architecture- Dependent Kernels
Release Process
33
Albany Code DemonstratingComponent-based Code
Design
Analysis Tools
Main()
Optimization
UQ
Mesh Database
Application
Exodus File
Problem Discretization
Solvers w/ Sensitivities
Nonlinear Model
Inline Mesher
Nonlinear
Optimization
Quality Improvement
Transient
Continuation
Stochastic Galerkin
Bifurcation
Load Balancing
Application
Linear Solve
PDE Assembly is Templated for AD, PCE
Linear Solvers / Preconditioners
Field Manager
Iterative
Domain Decomp
PDE TERMS
ManyCore Node Kernels
Block Iterative
Direct
MultiLevel
Discretization
Mulit-Core
Eigensolve
SchurComp
Accelerators
34
Applications in Albany are born with
Transformational Analysis Capabilities

• QCAD Quantum dot design tool.
• Optimization of gate voltages

• LCM Platform for RD in mechanics
• Load Stepping, AD of material models

ThermoElectrostatics ?Shape Optimization with
Embedded UQ
Std Deviation
2-Param Optimum
Initial Mesh
35
Summary

• LOCA and Trilinos provide powerful simulation and
  analysis capabilities
• Continuation, bifurcation, and linear stability
  analysis
• Scalable linear algebra
• Optimization, time integration, automatic
  differentiation, uncertainty quantification,
  discretization,
• Missing Capabilities (formerly future work)
• More generic algorithms for bordered matrix
  solves
• Much is hardwired to our Epetra format
• Periodic Orbit tracking beyond initial attempt
• Automated initial guess generation for null
  vectors
• Better documentation, examples, error checking,
  etc.
• Current Passion
• Component-based code design with Embedded
  Analysis in Mind from the beginning