# Presentation: 122319

## Advanced Software for the Calculation of Thermochemistry, Kinetics and Dynamicshttp://www.mcs.anl.gov/scidac/beskinetics

**ANL BES-SciDAC program has two inter-related projects**
**Parallelization of Cumulative Reaction Probabilities (****CRP****)**
**Parallel Implementation of **
**	****Subspace Projection Approximate Matrix  (****SPAM****) method**
**CRP *****(Stephen Gray and Al Wagner)***
**computationally intensive core of reaction rate constants**
_**mathematical kernel**_** (all matrices are sparse with some structure):**
**	****- method 1:****	****- iterative eigensolve (imbedded iterative linearsolves)**
**	****  ****	****	****- clever preconditioning important**
**	****	****	****- portability based on ANL PETSc library of kernels**
**	****- method 2:****	****- Chebyschev propagation (=> matrix vector multiplies)**
**	****  ****	****	****- novel finite difference representation (helps parallelize) **
_**programming issues**_**:**
**	****- parallelization **
**	****- exploiting data structure (i.e., preconditioning)**

## Advanced Software for the Calculation of Thermochemistry, Kinetics and Dynamics

**SPAM *****(Ron Shepard and Mike Minkoff)***
**novel iterative method to solve general matrix equations**
**eigensolve****	****	****- linear solve****	****   - nonlinear solve**
**applications are widespread **
**	****- in chemistry: CRP, electronic structure (SCF, MRSDCI,...)**
_**mathematical kernel**_**:**
**	****- related to Davidson, multigrid, and conjugate gradient methods**
**	****- subspace reduction (requiring usual matrix vector multiplies) **
**	****- projection operator decomposition of matrix vector product **
**	****- substitution of user-supplied approximate matrix **
**	****  ****in computationally intensive part of decomposition**
**	****- sequence of approx. matrices => multilevel method**
_**programming issues**_**:**
**	****- generalization of approach (only done for eigensolve)**
**	****- incorporation into libraries (connected to TOPS project part at ANL)**
**	****- test of efficacy in realistic applications (e.g., CRP)**

## Advanced Software for the Calculation of Thermochemistry, Kinetics & DynamicsParallelization of Cumulative Reaction Probabilities (CRP) Stephen Gray and Al Wagner (Argonne National Laboratory)

**Cumulative Reaction Probabilities (CRP)**
**	****- computational core of reaction rate constants**
**	****- exact computation computational intensive**
**	****- approximate computation underlies **
**	****  ****all major reaction rate theories in use**
**	****=> efficient exact CRP code will**
**	****	****- give exact rates (if the computed forces are accurate)**
**	****	****- calibrate ubiquitous approximate rate methods**
**Two methods:**
_**Time Independent**_** (Miller and Manthe, 1994, and others)**
_**Time Dependent **_**(Zhang and Light, 1996, and others)**
**Highly parallel approaches to both methods being pursued**

## Parallelization of Cumulative Reaction Probabilities Time Independent Approach

**• ****N(E) = ********k ****p****k****(E,J)**
**  ****where p****k****(E,J) = eigenvalues of Probability Operator:**
**   ****	****P(E) = 4 ********r****1/2**** (H+i********-E)****-1********p**** (H-i********-E)****-1********r****1/2**
**	****where i********x**** = absorbing diagonal potentials**
**                      ****= imaginary potentials**
**	****	****   ****H = hamiltonian (differential operator)**
**	****	****      ****=> for realistic problems**
**	****	****           ****size ~10****5****x10****5**** or much larger**
**	****	****          ****number of eigenvalues < 100**
**• ****iterative approach**
**	****macrocycle of iteration for eigenvalue**
**	****microcycle of iteration for**
**	****    ****action of Green’s function (H-i********-E)****-1**** **
**	****    ****=> linear solve**
**	****	****	**** **

## Parallelization of Cumulative Reaction Probabilities Time Independent Approach

**• ****PETSc: data structure, GMRES linear solve, preconditioners**
**  ****USER: Lanczos method for eigensolve**
**• ****Future: user supplied preconditioners**

**Parallelization of ****C****umulative ****R****eaction ****P****robabilities **_**Time Independent Approach**_

## Parallelization of Cumulative Reaction Probabilities Time Independent Approach

**Performance:**
**model problem**
** ****Optional number of dimensions**
**  • ****Eckhart potential along rxn coord.**
**  • ****parabolic potential perpendicular**
**    ****to reaction coord.**
** ****DVR representation of H**
**Algorithm options**
** ****Diagonal preconditioner**
** ****other PETSc preconditioners slower**
**Computers**
** ****NERSC SP**
** ****others include SGI Power Chanllenge,**
**  ****Cray T3E**

**Performance vs. processors**
**for model with increasing dimensions**

## Parallelization of Cumulative Reaction Probabilities Time Independent Approach

**Future Preconditioners:**
**• ****SPAM (see next poster)****	****	****             ****great ****	****	****       ****	****         scales**
**	****     ****Exact H****	****	****	****            performance ****	****	****	****         poorly**
**	****	****	****	****	****	****	****	****	****        ****- storage**
**	****	****	****	****	****	****	****	****	****        ****- banded**
**	****	****	****	****	****	**** ****	****	****	****          ****matrix**
**	****	****	****	****	****	**** ****	****	****	****          ****inversion**
**• ****Sparse optimal similarity transforms (Poirier) **
**	****IF Q =                         THEN find optimal Q such that QHQT =**
**	****	****	****	****	****	*****optimal block diagonal***** ->**

**fat band**
**precond.**
**SPAM**
**precond.**

**can SPAM get scalable performance?**

## Parallelization of Cumulative Reaction Probabilities Time Dependent Approach

**N(E) can be found from time dependent transition state wavepackets (TSWP)**
**(see Zhang and Light, J. Chem. Phys. 104, 6184 (1996))**
**N(E) = ********M**** N****i****(E)  **
**           ****where N****i****(E) ******** Im<********i****(E)|F********i****(E)>**
**	****                    ****where F      = differential flux operator**
**	****	****	****  ********i****(E) ********  ********exp(iEt)********i****(x,t) dt**
**	****	****	****	****where TSWP ********i****(x,t) from i********/********t ********i****(x,t) = H ********i****(x,t)**
**	****	****	****	****            ****where H is Schroedinger Eq. Operator**
**Work**
**• ****M different TSWPs (each independent of other)**
**• ****each TSWP **
**  ****- propagated over time**
**  ****- N****t**** time steps**
**  ****- each time step propagation dominated by H ********i**** multiply**
**  ****=> CPU work = M N****t**** (work of H ********i**** multiply)**

## Parallelization of Cumulative Reaction Probabilities Time Dependent Approach

**Solution Strategy**
**• ****Trivial parallelization over M wavepackets**
**  ****Nontrivial parallelization for wavepacket propagation**
**• ****Real Wavepackets (TS-RWP)**
**  ****(K. M. Forsythe and S. K. Gray, J. Chem. Phys. 112, 2623 (2000))**
**	****- half the storage, twice as fast relative to complex wavepackets**
**	****- Chebyshev iteration for propagation**
**	****- H********i**** work -> action of second order differential operators**

## Parallelization of Cumulative Reaction Probabilities Time Dependent Approach

**Solution strategy (more)**
**• ****Dispersion Fitted Finite Difference (DFFD)**
**  ****(S. K. Gray and E.M. Goldfield, J. Chem. Phys. 115, 8331 (2001))**
**	****- finite difference evaluation of action of differential operators**
**	****- optimized constants to reproduce dispersion relation**
**	****  ****(dispersion related momentum to kinetic energy)**
**	****- different optimized constants for selected propagation error ******
**  ****	********x.error for 3D H+H2**
**   ****Reaction Probability vs.**
**   ****order of finite difference**
**	****- parallelize via decomposition of ********i**** in x**
**	****  ****DFFD => less edge effects inducing processor communications **

## Parallelization of Cumulative Reaction Probabilities Time Dependent Approach

**Preliminary Results**
**• ****OH+H2 reaction **
**  ****- 6 dimensions**
**  ****- zero total ang. mom.**
** **
**• ****Wavepacket propagated **
**  ****reaction probability**
**• ****RWP + DFFD propagation **
**  ****technique**
**• ****up to 4 WinterHawk processor**
**  ****(MPI with shared memory off)**

## Advanced Software for the Calculation of Thermochemistry, Kinetics & Dynamics Subspace Projected Approximate Matrix: SPAMRon Shepard and Mike Minkoff (Argonne National Laboratory)

**New iterative method to solving general matrix equations**
**Eigenvalue**
**Linear**
**Non linear**
_**Method**_
**Based on subspace reduction, projection operators, decomposition, sequence of one or more approximate matrices**
_**Extensions**_
**Demonstrated for symmetric real eigensolves**
**Demonstrations on linear and non-linear equations planned**
_**Applications**_
**Applicable to problems with convergent**
**	****	****sequences of physical or numerical approximations**
**Parallelizable, multi-level library implementations via TOPS**

## Subspace Projected Approximate Matrix: SPAM Method

**• ****Subspace iterative solution to eigenvalue, linear, and nonlinear problems**
**	****e.g. eigenvalue problem  **(**H **–** ******j**** **)**v****j**** **=** 0**
**•**** ****Subspace iterative solutions have form  v****j**** = X****n**** c****j**** **
**	****	****where:****	****X****n**** = {x****1****,x****2****,...,x****n****}**
**• ****c****j**** is solved in a subspace:  H****n ****c****j**** = ********n****j ****c****j**** **
**	****	****where:****	**** H****n**** = (W****n****)****T ****X****n**** **
**	****	****	****where****: W****n ****= H X****n****    <---****for N>>n, where all the work is**
**• ****New x****n+1**** vector from residual:  r****n+1**** = (W****n****-********n****j ****X****n****)c****j**
**• ****SPAM gives more accurate or faster converging way to get x****n+1 ****in 3 steps**

## Subspace Projected Approximate Matrix: SPAM Method

**	****P****n ****x****n+1**** (=X****n****[X****n****]****T ****x****n+1****)****	****	****Q****n ****x****n+1**** (= x****n+1**** - P****n ****x****n+1****)**
**     ****where**
**       ****n trial vectors     n+1 trial vector  exact         matrix-vector         subspace**
**   ****(already processed)****	****   x****n+1****	****      matrix            product                matrix**
**                ****(X****n****) ****	****	**** ****	**

**S****ubspace ****P****rojected ****A****pproximate ****M****atrix: ****SPAM**** **_**Method**_

## Subspace Projected Approximate Matrix: SPAM Method

**Step 2: Decomposition and Approximation of H r****n+1**
**H x****n+1**** = (P****n****+Q****n****) H (P****n****+Q****n****) x****n+1**
**	****          ****cheap subspace operations****	****	****           expensive full space**
**	****	****	****	****	****	****                   ****operation**

**        ****cheap approximation**
**     ****to full matrix**

**Step 3: Solve approximate subspace problem**
**	****- solution is x****n+1**

## Subspace Projected Approximate Matrix: SPAM Method

**SPAM properties:**
**projection operators => convergence from any approx. matrix**
**multi-level SPAM with dynamic tolerances**
**	****	****Q****n****HQ****n**** approximated by ****Q****(1)n****H****(1)****Q****(1)n**** **
**	****	****	****Q****(1)n****H****(1)****Q****(1)n**** approximated by ****Q****(2)n****H****(2)****Q****(2)n**
**	****	****	****	****Q****(2)n****H****(2)****Q****(2)n**** approximated by ****Q****(3)n****H****(3)****Q****(3)n****...**
**tens of lines of code added to existing iterative subspace eigensolvers**
**highly parallelizable**
**applicable to any subspace problem**

## Subspace Projected Approximate Matrix: SPAM Method

**SPAM properties (broad view):**
**Relation to other subspace methods (e.g., Davidson)**
**More flexible (sequence of approx. matrices - no sequence => SPAM= Davidson)**
**If iteration tolerances are correct, always no worse than Davidson **
**Relation to multigrid methods**
**SPAM sequences of approx. matrices ~ multigrid sequences of approx. grids**
**Subspace method => solution vector composed of multiple vectors **
**Multigrid methods => single solution vector that is updated**
**Relation to preconditioned conjugate gradient (PCG) methods**
**SPAM has multiple vectors and approximations always improved by projection**
**	****PCG has a fixed single preconditioner and single vector improved by projection**
**Deep injection of physical insight into numerics**
**Application experts can design physical approximation sequences**
**SPAM maps approximation sequences onto numerics sequences**
**Projection operators continually improve the approximations**
**	****=> coarse approximations can still be numerically useful**

## Subspace Projected Approximate Matrix: SPAMExtensions

**Mathematical Extensions**
**Eigensolves**
**Symmetric real**
**Done with many applications **
**(http://chemistry.anl.gov/chem-dyn/Section-C-RonShepard.htm)**
**Code available **
**	****    ****(ftp: ftp.tcg.anl.gov/pub/spam/{README,spam.tar.Z})**
**Generalized symmetric planned**
**Generalized complex hermitian planned**
**General complex non-hermitian planned**
**Linear solves planned**
**Nonlinear solves planned**
**Formal connection to multigrid and conjugate gradient methods in progress**

## Subspace Projected Approximate Matrix: SPAMApplications

**broad view:**
**Any iterative problem solved in a subspace**
**	****with a user-supplied cheap approximate matrix (TOPS connection)**
**What is a cheap approximate matrix?**
**	****Sparser**
**	****Smaller underlying basis**
**	****Lower-order expansion of matrix elements**
**	****Coarser underlying grid **
**	****Lower-order difference equation**
**	****Tensor-product approximation **
**	****Operator approximation**
**           ****More highly parallelized approximation**
**Terascale Optimal PDE Siumlations (TOPS) connection**
**Basic parallelized multi-level code with user-supplied approx. matrix**
**Template approximate matrices stored in library**

## Subspace Projected Approximate Matrix: SPAMApplications

_**Two Chemistry Applications**_**:**
**Cumulative Reaction Probability ****(other poster)**
**Electronic structure:**
_**Non Chemistry Tensor Product Model**_
**• ****4x4 tensor products + tridiagonal perturbation**
**• ****Approx matrix is 4x4 tensor**

**exact matrix**

**exact matrix**

**approx. matrix**

**approx. matrix**