Selected publications by Christina Christara
-
A high-order deferred correction method for the solution of
free boundary problems using penalty iteration, with an
application to American option pricing,
Dawei Wang, Kirill Serkh, and Christina Christara
-
Bilateral XVA pricing under stochastic default intensity:
PDE modelling and computation,
Yuwei Chen and Christina Christara
-
Penalty and Penalty-Like Methods for Nonlinear HJB PDEs,
Christina Christara and Ruining Wu
-
Penalty Methods for Bilateral XVA Pricing,
in European and American Contingent Claims by a PDE Model,
Yuwei Chen and Christina Christara
-
Analysis of quantization error in financial pricing
via finite difference methods,
Christina C. Christara and Nat Chun-Ho Leung
-
Partial Differential Equation pricing of contingent claims
under stochastic correlation,
Nat Chun-Ho Leung, Christina C. Christara and Duy-Minh Dang
-
Spread option pricing using ADI methods,
Vida Heidarpour-Dehkordi and Christina C. Christara
-
PDE option pricing with variable correlations,
Christina C. Christara and Nat Chun-Ho Leung
-
Option pricing in jump diffusion models with
quadratic spline collocation,
Christina C. Christara and Nat Chun-Ho Leung
-
An efficient numerical PDE approach for pricing foreign exchange
interest rate hybrid derivatives,
Duy Minh Dang, Christina C. Christara, Kenneth R. Jackson and Asif Lakhany
-
Graphics processing unit pricing
of exotic cross-currency interest rate derivatives with a
foreign exchange volatility skew model,
Duy Minh Dang, Christina C. Christara and Kenneth R. Jackson
-
A Highly Efficient Implementation on GPU Clusters of
PDE-Based Pricing Methods for
Path-Dependent Foreign Exchange Interest Rate Derivatives,
Duy Minh Dang, Christina C. Christara and Kenneth R. Jackson
-
An Efficient GPU-Based Parallel Algorithm for
Pricing Multi-Asset American Options,
Duy Minh Dang, Christina C. Christara and Kenneth R. Jackson
-
Adaptive and high-order methods for valuing American options,
Christina C. Christara and Duy Minh Dang
-
A Parallel Implementation on GPUs of ADI Finite Difference Methods
for Parabolic PDEs with Applications in Finance,
Duy Minh Dang, Christina C. Christara and Kenneth R. Jackson
-
A PDE Pricing Framework for Cross-Currency Interest Rate Derivatives,
Duy Minh Dang, Christina C. Christara, Kenneth R. Jackson and Asif Lakhany
-
Pricing Multi-Asset American Options on Graphics
Processing Units Using a PDE Approach,
Duy Minh Dang, Christina C. Christara and Kenneth R. Jackson
-
Quartic spline collocation for second-order boundary value problems,
Christina C. Christara and Guohong Liu
-
Quadratic spline collocation for one-dimensional
linear parabolic partial differential equations,
Christina C. Christara, Tong Chen and Duy Minh Dang
-
Quartic spline collocation for fourth-order boundary value problems,
Christina C. Christara, Ying Zhu and Jingrui Zhang
-
A High-Performance Method for the Biharmonic Dirichlet Problem on Rectangles,
Christina C. Christara and Jingrui Zhang
-
Optimal Quadratic and Cubic Spline Collocation on Nonuniform Partitions,
Christina C. Christara, and Kit Sun Ng
-
Adaptive Techniques for Spline Collocation,
Christina C. Christara, and Kit Sun Ng
-
GIA-induced secular variations in the Earth's long wavelength
gravity field: Influence of 3-D viscosity variations
(short communication),
Konstantin Latychev, Jerry X. Mitrovica, Mark E. Tamisiea, Jeroen Tromp,
Christina C. Christara and Robert Moucha
-
Glacial isostatic adjustment on 3-D Earth models: a finite-volume formulation,
Konstantin Latychev, Jerry X. Mitrovica, Jeroen Tromp, Mark E. Tamisiea,
Dimitri Komatitsch, and Christina C. Christara
-
Optimal Quadratic Spline Collocation on Non-Uniform Partitions,
Christina C. Christara, and Kit Sun Ng
-
Optimal Cubic Spline Collocation on Non-Uniform Partitions,
Christina C. Christara, and Kit Sun Ng
-
Quadratic Spline Galerkin Method for the Shallow Water Equations,
Anita T. Layton, Christina C. Christara, and Kenneth R. Jackson
-
Optimal Quadratic Spline Collocation Methods for the
Shallow Water Equations,
Anita T. Layton, Christina C. Christara, and Kenneth R. Jackson
-
Quadratic Spline Collocation Revisited:
Extension to Systems of Elliptic PDEs,
Christina C. Christara and Kit Sun Ng
-
Fast Fourier Transform Solvers and Preconditioners
for Quadratic Spline Collocation,
Christina C. Christara and Kit Sun Ng
-
An Efficient Transposition Algorithm for Distributed Memory Computers,
C. Christara, X. Ding and K. R. Jackson
-
Multigrid and Multilevel Methods for Quadratic Spline Collocation,
Christina C. Christara and Barry Smith
-
Parallel solvers for spline collocation equations, C. C. Christara
-
Quadratic spline collocation methods for
elliptic partial differential equations, C. C. Christara
-
High performance computing of elliptic partial differential equations
with spline collocation, C. C. Christara
-
Parallel computation of partial differential equations
on distributed memory machines, C. C. Christara
-
Multicolour orderings and iterative methods for elliptic equations,
C. C. Christara
-
Domain decomposition and incomplete factorisation methods for
partial differential equations, C. C. Christara
-
Schur complement preconditioned conjugate gradient methods for
spline collocation equations, C. C. Christara
-
Conjugate gradient methods for spline collocation equations, C. C. Christara
-
A Domain Decomposition Spline Collocation Method for Elliptic Partial
Differential Equations,
C. C. Christara and E. N. Houstis
-
A Parallel Spline Collocation - Capacitance Method for Elliptic Partial
Differential Equations,
C. C. Christara, E. N. Houstis and J. R. Rice
-
Geometry Decomposition based Methods for solving Elliptic Partial
Differential Equations,
C. C. Christara, A. Hadjidimos, E. N. Houstis, J. R. Rice and E. A. Vavalis
-
Quadratic Spline Collocation Methods for two point boundary value problems,
E. N. Houstis, C. C. Christara and J. R. Rice
-
Performance of Scientific Software,
E. N. Houstis, J. R. Rice, C. C. Christara and E. A. Vavalis
-
Scientific Computing by Numerical Methods,
C. C. Christara and K. R. Jackson
Supervised students' theses
-
Numerical PDE Methods for Pricing Bilateral XVA,
Yuwei Chen, Ph.D. thesis, University of Toronto,
February 2023, 123 pages.
-
Penalty Methods for Nonlinear Problems in Financial Option Pricing,
Ruining (Ray) Wu, MSc Research Paper,
February 2021, 121 pages.
-
Numerical Methods for Pricing Multi-Asset Options,
Yuwei Chen, M.Sc. thesis, University of Toronto,
January 2018, 75 pages.
-
PDE option pricing: analysis and application
to stochastic correlation,
Nat Chun-Ho Leung, Ph.D. thesis, University of Toronto,
June 2017, 119 pages.
-
Efficient and accurate numerical PDE methods
for pricing financial derivatives,
Mufan (Bill) Li, B.A.Sc. thesis, University of Toronto,
April 2015, 35 pages.
-
Modeling multi-factor financial derivatives
by a Partial Differential Equation approach
with efficient implementation on Graphics Processing Units,
Duy Minh Dang, Ph.D. thesis, University of Toronto,
August 2011, 228 pages.
-
A study of a discrete element method based granular dynamics solver,
Tina Yee, M.Sc. thesis, University of Toronto, 2009, 43 pages.
-
Multigrid and Spline Collocation Methods on Non-uniform Grids,
Guoyu Wang, M.Sc. thesis, University of Toronto, September 2008, 98 pages.
-
Bi-quartic Spline Collocation Methods for Fourth-order
Boundary Value Problems with an Application to
the Biharmonic Dirichlet Problem,
Jingrui Zhang, Ph.D. thesis,
University of Toronto, January 2008, 160 pages.
-
Adaptive Finite Difference Methods for Valuing American Options,
Duy Minh Dang, M.Sc. Thesis
-
Quartic Spline Collocation Methods For Second-Order Two-Point
Boundary Value ODE Problems,
Guohong Liu, M.Sc. thesis
-
Pricing Convertible Bonds with Dividend Protection
subject to Credit Risk Using a Numerical PDE Approach,
Qingkai Mo, M.Sc. Thesis
-
An Efficient Algorithm Based on Quadratic Spline Collocation and
Finite Difference Methods for Parabolic Partial Differential Equations,
Tong Chen, M.Sc. Thesis
-
Pricing Convertible Bonds using Partial Differential Equations,
Lucy Xingwen Li, M.Sc. Thesis
-
Spline Collocation on Adaptive Grids and Non-Rectangular Domains,
Kit Sun Ng, Ph.D. Thesis
-
Multigrid and Cubic Spline Collocation Methods for Advection Equations,
Zheng Zeng, M.Sc. Thesis
-
Towards a high-performance method for the biharmonic Dirichlet problem,
Jingrui Zhang, Research Paper (equivalent to M.Sc. Thesis)
-
Quartic-Spline Collocation Methods for Fourth-Order Two-Point
Boundary Value Problems,
Ying Zhu, M.Sc. Thesis
-
High-Order Spatial Discretization Methods for the Shallow Water
Equations,
Anita W. Tam, Ph.D. Thesis, co-supervised by Ken Jackson
-
Quadratic Spline Collocation Methods for Systems of Elliptic PDEs,
Kit Sun Ng, M.Sc. Thesis
-
Numerical Solution of the Shallow-Water Equations on Distributed
Memory Systems,
Xiaoliang (Lloyd) Ding, M.Sc. Thesis, co-supervised by Ken Jackson
-
Fast Fourier transform solvers for quadratic spline collocation,
Athanasia Constas, M.Sc. Thesis
Title: |
A high-order deferred correction method for the solution of
free boundary problems using penalty iteration, with an
application to American option pricing
|
Authors: | Dawei Wang, Kirill Serkh, and Christina Christara
|
Date: | May 2022, updated November 2022, March 2023 |
Pages: | 41 |
Note: | Submitted |
Keywords: | Green's functions; deferred correction;
free boundary problems; immersed interface methods;
finite differences; American options;
Linear Complementarity Problem; penalty iteration; derivative jumps
|
File: |
HOFB.pdf
|
Abstract:
This paper presents a high-order deferred correction algorithm combined with
penalty iteration for solving free and moving boundary problems, using a
fourth-order finite difference method. Typically, when free boundary
problems are solved on a fixed computational grid, the order of the solution
is low due to the discontinuity in the solution at the free boundary, even
if a high order method is used. Using a detailed error analysis, we observe
that the order of convergence of the solution can be increased to
fourth-order by solving successively corrected finite difference systems,
where the corrections are derived from the previously computed lower order
solutions. The corrections are applied solely to the right-hand side,
and leave the finite difference matrix unchanged.
The penalty iterations converge quickly given a good initial guess.
We demonstrate the accuracy and efficiency of our algorithm using
several examples. Numerical results show that our algorithm gives up to
fourth-order convergence for both the solution and the free boundary
location. We also test our algorithm on the challenging American put option
pricing problem. Our algorithm gives the expected high-order convergence.
Title: |
Bilateral XVA pricing under stochastic default intensity:
PDE modelling and computation
|
Authors: | Yuwei Chen and Christina Christara |
Date: | June 2022 |
Pages: | 32 |
Note: | A later version is to be published in
Applied Numerical Mathematics, vol 185, pages 236-259, March 2023.
link to journal page
|
Keywords: | Partial Differential Equations, Black-Scholes,
Crank-Nicolson finite differences, mean reversion CIR process,
stochastic default intensity, asymptotic approximation,
nonlinear iteration, penalty iteration method.
|
File: |
sto_inten_xva_v3.pdf
|
Abstract:
Adjusting derivative prices to take into account default risk
has attracted the attention of several researchers and practitioners,
especially after the 2007-2008 financial crisis.
We derive a novel partial differential equation (PDE) model
for derivative pricing including the adjustment for default risk,
assuming that the default risk of one of
the counterparties (the buyer) follows a
Cox-Ingersoll-Ross (CIR) process, while the other party
has constant default risk.
The time-dependent PDE derived is of Black-Scholes type and
involves two ``space'' variables, namely the asset price
and the buyer default intensity,
as well as a nonlinear source term.
We formulate boundary conditions appropriate
for the default intensity variable.
The numerical solution of the PDE is based
on standard finite differences, and a penalty-like iteration
for handling the nonlinearity.
We also develop and analyze a novel asymptotic approximation formula
for the adjusted price of derivatives,
resulting in a very efficient, accurate, and convenient
for practitioners formula.
We present numerical results that indicate
stable second order convergence for the 2D PDE solution
in terms of the discretization size.
We compare the effectiveness of the 2D PDE and asymptotic approximations.
We study the effect of various numerical and market parameters
to the values of the adjusted prices
and to the accuracy of the computed solutions.
Title: |
Penalty and Penalty-Like Methods for Nonlinear HJB PDEs
|
Authors: | Christina Christara and Ruining Wu |
Date: | June 2022 |
Pages: | 27 |
Note: | A later version is published in
Applied Mathematics and Computation, vol 425, pages 1-19, 2022.
link to journal page
|
Keywords: | Partial Differential Equations, Black-Scholes,
nonlinear iteration, finite differences, Crank-Nicolson,
control problem, Hamilton-Jacobi-Bellman (HJB) equation,
transaction costs, stock borrowing fees, penalty methods, policy iteration.
|
File: |
hjb-penalty.pdf
|
Abstract:
There are numerous financial problems that can be posed as
optimal control problems, leading to Hamilton-Jacobi-Bellman
or Hamilton-Jacobi-Bellman-Issacs equations.
We reformulate these problems as nonlinear PDEs,
involving max and/or min terms of the unknown function,
and/or its first and second spatial derivatives.
We suggest efficient numerical methods for handling the nonlinearity
in the PDE through an adaptation of the discrete penalty method
[Forsyth, Vetzal 2002] that gives rise to tridiagonal penalty matrices.
We formulate a penalty-like method for the use with European exercise rights,
and extend this to American exercise rights resulting
in a double-penalty method.
We also use our findings to improve the policy iteration algorithms
described in [Forsyth, Labahn 2007].
Numerical results are provided showing clear second-order convergence,
and where applicable, we prove the convergence of our algorithms.
Title: | Penalty Methods for Bilateral XVA Pricing
in European and American Contingent Claims by a PDE Model
|
Authors: | Yuwei Chen and Christina C. Christara |
Date: | Nov 2019, Rev. Feb 2020
|
Pages: | 22
|
Keywords: |
Partial Differential Equations, Black-Scholes,
credit risk, default risk, credit valuation adjustment,
call option, put option, forward contract,
nonlinear iteration, finite differences, Crank-Nicolson,
control problem, Hamilton-Jacobi-Bellman (HJB) equation.
|
Note: | A version of this paper will appear in the
Journal of Computational Finance
|
File: | ccc/xva_full.pdf
|
Abstract:
Taking into account default risk
in the valuation of financial derivatives
has become increasingly important,
especially after the 2007-2008 financial crisis.
Under some assumptions, the valuation of financial derivatives,
including a value adjustment to account for default risk
(the so-called XVA),
gives rise to a nonlinear PDE \cite{Christoph Burgard and Mats Kjaer, 2011}.
We propose numerical methods for handling the nonlinearity in this PDE,
the most efficient of which are discrete penalty iteration methods.
We first formulate a penalty iteration method for the case of European
contingent claims, and study its convergence.
We then extend the method to the case of American contingent claims,
resulting in a double-penalty iteration.
We also propose boundary conditions
and their discretization for the XVA PDE problem
in the cases of call and put options, as well as
the case of a forward contract.
Numerical results demonstrate the effectiveness of our methods.
Title: | Analysis of quantization error in financial pricing
via finite difference methods
|
Authors: | Christina C. Christara and Nat Chun-Ho Leung |
Date: | July 2017, rev. March 2018
|
Pages: | 27
|
Keywords: |
non-smooth initial conditions, option pricing,
numerical solution, partial differential equation,
convection-diffusion equations, Fourier analysis, finite difference methods
|
Note: | A version of this paper is published in SIAM J. Numerical Analysis, 56:3, 2018, pp. 1731-1757
|
File: | ccc/paper-error.pdf
|
Abstract:
In this paper, we study the error of a second order finite
difference scheme for the one-dimensional convection-diffusion
equation. We consider non-smooth initial conditions commonly
encountered in financial pricing applications. For these initial
conditions, we establish the explicit expression of the
quantization error, which is loosely defined as the error
of the numerical solution due to the placement of the point of
non-smoothness on the numerical grid. Based on our analysis,
we study the issue of optimal placement of such
non-smoothness points on the grid, and the effect of smoothing
operators on quantization errors.
Title: | Partial Differential Equation pricing of contingent claims
under stochastic correlation
|
Authors: | Nat Chun-Ho Leung, Christina C. Christara and Duy-Minh Dang |
Date: | October 2017
|
Pages: | 30
|
Keywords: |
stochastic correlation, option pricing,
numerical solution, asymptotic solution, partial differential equation
|
Note: | A version of this paper appeared in
SIAM J. Scientific Computing, Vol 40, Issue 1, pages B1-B31 (2018).
|
File: | ccc/paper_sto_corrL.pdf
|
Abstract:
In this paper, we study a partial differential equation (PDE) framework
for option pricing
where the underlying factors exhibit stochastic correlation,
with an emphasis on computation.
We derive a multi-dimensional time-dependent PDE
for the corresponding pricing problem, and present
a numerical PDE solution. We prove a stability result,
and study numerical issues regarding the boundary conditions used.
Moreover, we develop and analyze an asymptotic analytical approximation
to the solution, leading to a novel computational asymptotic approach
based on quadrature with a perturbed transition density.
Numerical results are presented to verify second order convergence
of the numerical PDE solution and to demonstrate its agreement
The effect of certain problem parameters to the PDE solution,
as well as to the asymptotic approximation solution, is also studied.
Title: | Spread option pricing using ADI methods
|
Authors: | Vida Heidarpour-Dehkordi and Christina C. Christara |
Date: | July 2017
|
Pages: | 17
|
Keywords: |
Modified Craig-Sneyd, Alternating Direction Implicit method,
two-dimensional Black-Scholes, American option, spread option, exchange option,
analytical approximation, numerical PDE solution, penalty iteration
|
Note: | Published in
International J. Numerical Analysis and Modelling,
15:3, 2018. pp. 353-369.
|
File: | ccc/spreadADI.pdf
|
Abstract:
Spread option contracts are becoming increasingly important,
as they frequently arise in the energy derivative markets,
e.g. exchange electricity for oil.
In this paper, we study the pricing of European and American
spread options.
We consider the two-dimensional Black-Scholes
Partial Differential Equation (PDE),
use finite difference discretization in space
and consider Crank-Nicolson (CN)
and Modified Craig-Sneyd (MCS) Alternating Direction Implicit (ADI)
methods for timestepping.
In order to handle the early exercise feature arising in American options,
we employ the discrete penalty iteration method,
introduced and studied
in Forsyth and Vetzal (2002), for one-dimensional PDEs
discretized in time by the CN method.
The main novelty of our work is the incorporation
of the ADI method into the discrete penalty iteration method,
in a highly efficient way,
so that it can be used for two or higher-dimensional problems.
The results from spread option pricing are compared with those
obtained from the closed-form approximation formulae
of Kirk (1995), Venkatramanan and Alexander (2011),
Monte Carlo simulations, and
the Brennan-Schwartz ADI Douglas-Rachford method,
as implemented in MATLAB.
In all spread option test cases we considered, including American ones,
our ADI-MCS method, implemented on appropriate non-uniform grids,
gives more accurate prices and Greeks than the MATLAB ADI method.
Title: | PDE option pricing with variable correlations
|
Authors: | Christina C. Christara and Nat Chun-Ho Leung |
Date: | May 2017
|
Pages: | 17
|
Keywords: |
stochastic correlation, option pricing, regime switching,
numerical solution, partial differential equation, system of PDEs,
exchange option, basket option, Greeks.
|
Note: | Proceedings of the
International Conference on
Computational and Informational Sciences and Engineering,
honoring Professor Elias N. Houstis, University of Thessaly publications,
September 2018, pp. 43-57.
|
File: | ccc/paperRS.pdf
|
Abstract:
Modelling correlation between financial quantities
is important in the accurate pricing of financial derivatives.
In this paper,
we introduce some stochasticity in correlation,
by considering a regime-switching correlation model,
in which the transition rates between regimes are given.
We present a derivation of the associated
Partial Differential Equation (PDE) problem.
The problem involves a system of $\lsc$ PDEs,
where $\lsc$ is the number of regimes.
We formulate a finite difference method for the solution
of the PDE system, and numerically demonstrate that
it converges with second order.
We study the effect of certain model parameters on the computed prices.
We compare prices from this model, associated PDE and method with those from
a stochastic correlation model, associated PDE and method
in \cite{Emmerich06, Leung2017, LeungChristara2017}
and discuss advantages and disadvantages.
Title: | Option pricing in jump diffusion models with
quadratic spline collocation
|
Authors: | Christina C. Christara and Nat Chun-Ho Leung |
Date: | May 2014; revised May 2015
|
Pages: | 15
|
Keywords: |
quadratic spline, collocation, finite difference,
European option, American option, partial integro-differential equation,
Merton model, Kou model, calculation of Greeks
|
Note: | Published in
Applied Mathematics and Computation, Vol 279, 10 April 2016, pp 28-42.
|
File: | ccc/paper-jump.pdf,
|
Abstract:
In this paper, we develop a robust numerical method in pricing options,
when the underlying asset follows a jump diffusion model.
We demonstrate that, with the quadratic spline collocation method,
the integral approximation in the pricing PIDE is intuitively simple,
and comes down to the evaluation of the probabilistic moments
of the jump density. The calculation of Greeks is also simple.
When combined with a Picard iteration scheme, the pricing problem
can be solved efficiently.
We present the method and the numerical results from pricing
European and American options with Merton's and Kou's models.
Title: | An efficient numerical PDE approach for pricing foreign exchange
interest rate hybrid derivatives
|
Authors: | Duy Minh Dang, Christina C. Christara, Kenneth R. Jackson and
Asif Lakhany |
Date: | March 2012; revised May 2013
|
Pages: | 54
|
Keywords: |
Power-Reverse Dual-Currency (PRDC) swaps, Target Redemption (TARN),
knockout, Partial Differential Equation (PDE), finite differences,
non-uniform grids.
|
Note: | Published in
Journal of Computational Finance, Vol 18, Issue 4, 2015, pp 39-93.
|
File: | prdcTARN.pdf,
Link to the paper on the SSRN
|
Abstract:
We discuss efficient pricing methods via a Partial Differential Equation (PDE)
approach for long-dated foreign exchange (FX) interest rate hybrids
under a three-factor multi-currency pricing model with FX volatility skew.
The emphasis of the paper is on
Power-Reverse Dual-Currency (PRDC) swaps with popular exotic features,
namely knockout and FX Target Redemption (FX-TARN).
Challenges in pricing these derivatives via a PDE approach
arise from the high-dimensionality of the model PDE,
as well as from the complexities in handling the exotic features,
especially in the case of the FX-TARN provision,
due to its path-dependency.
Our proposed PDE pricing framework
for FX-TARN PRDC swaps is based on partitioning
the pricing problem into several independent pricing
sub-problems over each time period of the swap's tenor structure,
with possible communication at the end of the time period.
Each of these pricing sub-problems can be viewed as
equivalent to a knockout PRDC swap with a known time-dependent barrier,
and requires a solution of the model PDE, which, in our case,
is a time-dependent parabolic PDE in three space dimensions.
Finite difference schemes on non-uniform grids are used
for the spatial discretization of the model PDE,
and the Alternating Direction Implicit (ADI)
timestepping methods are employed for its time discretization.
Numerical examples illustrating the convergence properties
and efficiency of the numerical methods are provided.
Title: | Graphics processing unit pricing
of exotic cross-currency interest rate derivatives with a
foreign exchange volatility skew model,
|
Authors: | Duy Minh Dang, Christina C. Christara and Kenneth R. Jackson
|
Date: | February 2011
|
Pages: | 17
|
Keywords: |
Power Reverse Dual Currency (PRDC) Swaps, Bermudan Cancelable,
Partial Differential Equation (PDE), Alternating Direction Implicit (ADI),
Finite Differences, Graphics Processing Units (GPUs), GPU Clusters, MPI,
Parallel Computing
|
Note: | Published in
Journal of Concurrency and Computation: Practice and Experience (JCCPE)
Volume 26, Issue 9, 2014, pp.~1609--1625.
|
File: | prdc.gpu.bc.mod.pdf,
Link to the paper on the SSRN
|
Abstract:
We present a Graphics Processing Unit (GPU)
parallelization of the computation of the price of exotic
cross-currency interest rate derivatives via a Partial Differential
Equation (PDE) approach. In particular, we focus on the GPU-based
parallel pricing of long-dated foreign exchange (FX) interest rate hybrids,
namely Power Reverse Dual Currency (PRDC) swaps
with Bermudan cancelable features.
We consider a three-factor pricing model with
FX volatility skew which results in a
time-dependent parabolic PDE in three spatial dimensions.
for the spatial discretization of the PDE,
and the Alternating Direction Implicit (ADI) technique is
employed for the time discretization.
We then exploit the parallel architectural features
of GPUs together with the Compute Unified Device
Architecture (CUDA) framework to design and implement an
efficient parallel algorithm for pricing PRDC swaps.
Over each period of the tenor structure,
we divide the pricing of a Bermudan cancelable PRDC swap
into two independent pricing subproblems,
each of which can efficiently be solved on a GPU via
a parallelization of the ADI timestepping technique.
Numerical results indicate that GPUs can provide
significant increase in performance over CPUs when pricing PRDC swaps.
An analysis of the impact of the FX skew on such derivatives is provided.
Title: | A Highly Efficient Implementation on GPU Clusters of
PDE-Based Pricing Methods for
Path-Dependent Foreign Exchange Interest Rate Derivatives
|
Authors: | Duy Minh Dang, Christina C. Christara and Kenneth R. Jackson
|
Date: | Submitted March 2013, revised April 2013.
|
Pages: | 17
|
Keywords: |
Power Reverse Dual Currency (PRDC) Swaps, Bermudan Cancelable,
Partial Differential Equation (PDE), Alternating Direction Implicit (ADI),
Finite Differences, Graphics Processing Units (GPUs), GPU Clusters, MPI,
Parallel Computing
|
Note: | Published in
Springer Lecture Notes in Computer Science,
Volume 7975, 2013, pp.~107--126,
Proceedings of
The 13th International Conference on Computational Science
and Its Applications (ICCSA 2013),
Ho Chi Minh City, Vietnam, 24--27 June 2013.
|
File: | dang-ccc-krj-iccsa13.pdf,
Link to the paper on the SSRN
|
Abstract:
We present a highly efficient parallelization of the computation of
the price of exotic cross-currency interest rate derivatives with
path-dependent features via a Partial Differential Equation (PDE)
approach. In particular, we focus on the parallel pricing on Graphics
Processing Unit (GPU) clusters of long-dated foreign exchange (FX)
interest rate derivatives, namely Power-Reverse Dual-Currency (PRDC)
swaps with FX Target Redemption (FX-TARN) features under a three-factor
model. Challenges in pricing these derivatives via a PDE approach
arise from the high-dimensionality of the model PDE, as well as from
path-dependency of the FX-TARN feature. The PDE pricing framework for
FX-TARN PRDC swaps is based on partitioning the pricing problem into
several independent pricing sub-problems over each time period of the
swap's tenor structure, with possible communication at the end of the
time period. Finite difference methods on non-uniform grids are used for
the spatial discretization of the PDE, and the Alternating Direction
Implicit (ADI) technique is employed for the time discretization. Our
implementation of the pricing procedure on a GPU cluster involves
(i) efficiently solving each independent sub-problems on a GPU via a
parallelization of the ADI timestepping technique, and (ii) utilizing
MPI for the communication between pricing processes at the end of the
time period of the swap's tenor structure. Numerical results showing
the efficiency of the parallel methods are provided.
Title: | An Efficient GPU-Based Parallel Algorithm for
Pricing Multi-Asset American Options
|
Authors: | Duy Minh Dang, Christina C. Christara and Kenneth R. Jackson
|
Date: | March 2011
|
Pages: | 14
|
Keywords: | American Option, Multi-Asset, Penalty Method,
Alternating Direction Implicit Approximate Factorization
(ADI-AF), Graphics Processing Units, GPUs,
Parallel Computing, Finite Difference
|
Note: | A revised version of this paper (with a slightly different title)
appeared in the
Journal of Concurrency and Computation: Practice and Experience (JCCPE)
special issue on high-performance computational finance, Vol. 24, No. 8, June 2012,
pp. 849-866.
|
File: |
Link to the paper on the SSRN
|
Abstract:
We develop highly-efficient parallel
Partial Differential Equation (PDE) based pricing methods
on Graphics Processing Units (GPUs) for multi-asset American options.
Our pricing approach is built upon a combination
of a discrete penalty approach for the
linear complementarity problem arising due to the free boundary
and a GPU-based parallel Alternating Direction Implicit
Approximate Factorization technique with finite differences
on uniform grids for the solution of the linear algebraic system
arising from each penalty iteration.
A timestep size selector implemented efficiently on GPUs
is used to further increase the efficiency of the methods.
We demonstrate the efficiency and accuracy of the parallel numerical
methods by pricing American options written on three assets.
Title: | Adaptive and high-order methods for valuing American options
|
Authors: | Christina C. Christara and Duy Minh Dang |
Date: | December 2009
|
Pages: | 24
|
Keywords: | adaptive mesh selection, error equidistribution,
quadratic splines, collocation, finite differences,
European option, American option, penalty method
|
Note: | A later version of this paper is published in the
Journal of Computational Finance, Vol 14, Issue 4, 2011, pages 73-113.
|
File: | ccc/americanOption.pdf
|
Abstract:
We develop space-time adaptive and high-order methods for
valuing American options using a partial differential equation
approach.
The linear complementarity problem arising due to
the free boundary is handled by a penalty method.
Both finite difference and finite element methods are considered
for the space discretization of the PDE, while classical finite differences,
such as Crank-Nicolson, are used for the time discretization.
The high-order discretization in space is based on
an optimal finite element collocation method,
the main computational requirements of which are the solution of one
tridiagonal linear system at each time step, while the resulting
errors at the grid points and midpoints of the space partition are
fourth-order. To control the space error, we use adaptive
gridpoint distribution based on an error equidistribution principle.
A time stepsize selector is used to further increase the efficiency
of the methods. Numerical examples show that our methods converge fast
and provide highly accurate options prices, Greeks, and early exercise
boundaries.
Title: | A Parallel Implementation on GPUs of ADI Finite Difference
Methods for Parabolic PDEs with Applications in Finance
|
Author: | Duy Minh Dang, Christina C. Christara and Kenneth R. Jackson
|
Date: | March 2010; revised December 2010; published 2011
|
Pages: | 21
|
Keywords: | Alternating Direction Implicit, ADI,
Partial Differential Equation, PDE, Graphics Processing
Units, GPUs, parallel computing, finite difference,
multi-asset options, multi-currency swaps
|
Note: | Appeared in the
Canadian Applied Mathematics Quarterly (CAMQ),
Vol. 17, No. 4, 2009, pp. 627-660,
which was published in 2011, for a belated 30-year anniversary issue of CAMQ.
|
File: | Link to the paper on the SSRN
|
Abstract:
We study a parallel implementation on a Graphics Processing Unit (GPU)
of Alternating Direction Implicit (ADI) time-discretization methods for
solving time-dependent parabolic Partial Differential Equations (PDEs)
in three spatial dimensions with mixed spatial derivatives in a variety
of applications in computational finance. Finite differences on uniform
grids are used for the spatial discretization of the PDEs. As examples,
we apply the GPU-based parallel methods to price European rainbow and
European basket options, each written on three assets. Numerical results
showing the efficiency of the parallel methods are provided.
Title: | A PDE Pricing Framework for Cross-Currency Interest Rate
Derivatives
|
Author: | Duy Minh Dang, Christina C. Christara, Kenneth R. Jackson
and Asif Lakhany |
Date: | November 2009
|
Pages: | 10
|
Keywords: | Power Reverse Dual Currency Swaps, Bermudan Cancelable,
Partial Differential Equation, PDE,
Alternating Direction Implicit, ADI,
Generalized Minimal Residual, GMRES,
Fast Fourier Transform, FFT
|
Note: | Submitted to the
"2010 Workshop on Computational Finance and Business Intelligence",
International Conference In Computational Science (ICCS),
Amsterdam, May 31 to June 2, 2010.
|
File: | Link to the paper on the SSRN
|
Abstract:
We propose a general framework for efficient pricing via a Partial
Differential Equation (PDE) approach of cross-currency interest rate
derivatives under the Hull-White model. In particular, we focus on pricing
long-dated foreign exchange (FX) interest rate hybrids, namely Power
Reverse Dual Currency (PRDC) swaps with Bermudan cancelable features. We
formulate the problem in terms of three correlated processes that
incorporate FX skew via a local volatility function. This formulation
results in a parabolic PDE in three spatial dimensions. Finite
difference methods are used for the spatial discretization of the
PDE. The Crank-Nicolson method and the Alternating Direction Implicit
(ADI) method are considered for the time discretization. In the former
case, the preconditioned Generalized Minimal Residual (GMRES) method is
employed for the solution of the resulting block banded linear system at
each time step, with the preconditioner solved by Fast Fourier Transform
(FFT) techniques. Numerical examples illustrating the convergence
properties of the methods are provided.
Title: | Pricing Multi-Asset American Options on Graphics
Processing Units Using a PDE Approach
|
Authors: | Duy Minh Dang, Christina C. Christara and Kenneth R. Jackson
|
Date: | September 2010
|
Pages: | 8
|
Keywords: | American Option, Multi-Asset, Penalty Method,
Alternating Direction Implicit Approximate Factorization
(ADI-AF), Graphics Processing Units, GPUs,
Parallel Computing, Finite Difference
|
Note: | Proceedings of the
Third Workshop on High Performance Computational Finance
at SC10, The International Conference for High Performance Computing,
Networking, Storage, and Analysis, sponsored by ACM/IEEE,
New Orleans, USA, 13-19 November 2010.
|
File: |
DOI:10.1109/WHPCF.2010.5671831
ccc/american.gpu.ieee.pdf,
link to pdf file in ieeexplore
|
Abstract:
We develop highly efficient parallel pricing methods on Graphics
Processing Units (GPUs) for multi-asset American options via a Partial
Differential Equation (PDE) approach. The linear complementarity problem
arising due to the free boundary is handled by a penalty method. Finite
difference methods on uniform grids are considered for the space
discretization of the PDE, while classical finite differences, such
as Crank-Nicolson, are used for the time discretization. The discrete
nonlinear penalized equations at each timestep are solved using a
penalty iteration. A GPU-based parallel Alternating Direction Implicit
Approximate Factorization technique is employed for the solution of
the linear algebraic system arising from each penalty iteration. We
demonstrate the efficiency and accuracy of the parallel numerical methods
by pricing American options written on three assets.
Title: |
Quartic spline collocation for second-order boundary value problems
|
Authors: | Christina C. Christara and Guohong Liu |
Date: | September 2009
|
Pages: | 8
|
Keywords: | sixth order convergence, quartic splines,
second-order BVP, collocation
|
Note: | Proceedings of the
9th HERCMA Conference on Computer Mathematics and Applications,
Athens University of Economics and Business,
September 23-26, 2009, Athens, Greece, pp. 1-8.
|
File: | ccc/paper-qrtbvp2.pdf
|
Abstract:
Collocation methods based on quartic splines are presented for
second-order two-point boundary value problems.
In order to obtain a uniquely solvable linear system
for the degrees of freedom of the quartic spline collocation
approximation, in addition to the boundary conditions
specified by the problem, extra boundary or near-boundary
conditions are introduced.
Non-optimal (fourth-order) and optimal (sixth-order)
quartic-spline methods are considered.
The theoretical behavior of the collocation methods
is verified by numerical experiments.
The extension of the methods to two-dimensional problems
is briefly considered.
Abstract:
New methods for solving
general linear parabolic partial differential equations (PDEs)
in one space dimension are developed.
The methods combine
quadratic-spline collocation for the space discretization and
classical finite differences, such as Crank-Nicolson,
for the time discretization.
The main computational requirements of the most efficient method
are the solution of one tridiagonal linear system
at each time step, while
the resulting errors at the gridpoints and midpoints
of the space partition are fourth order.
The stability and convergence properties of some of the new methods
are analyzed for a model problem.
Numerical results demonstrate the stability and accuracy of the methods.
Adaptive mesh techniques are introduced in the space dimension,
and the resulting method is applied to the American put option
pricing problem, giving very competitive results.
Title: | Quartic spline collocation for fourth-order boundary value problems
|
Authors: | Christina C. Christara, Ying Zhu and Jingrui Zhang
|
Date: | September 2008
|
Pages: | 6
|
Note: | Proceedings of the 2008 Numerical Analysis conference,
September 1-5, 2008, Kalamata, Greece, pp. 62-67.
|
File: | ccc/paper-qrtbvp4.pdf
|
Abstract:
We consider the numerical solution
of linear fourth-order boundary value problems (BVPs)
in one and two dimensions,
by methods based on quartic splines and the collocation
methodology. The discretization error is sixth order
on the gridpoints and midpoints of a uniform partition
and fifth order globally in the uniform norm.
For the linear systems arising from the
discretization of certain biharmonic problems
by quartic spline collocation, we develop fast solvers
based on Fourier transforms, leading to asymptotically
almost optimal solution techniques.
Title: | A High-Performance Method for the Biharmonic Dirichlet
Problem on Rectangles
|
Authors: | Christina C. Christaramand Jingrui Zhang |
Date: | July 2007
|
Pages: | 3
|
Note: | Proceedings of the 2007 International Conference
On Preconditioning Techniques
for Large Sparse Matrix Problems in
Scientific and Industrial Applications,
July 9-12, 2007 Météopole, Toulouse, France, pp. 35-37.
|
File: | Please contact ccc SYMBOL_AT cs DOT toronto DOT edu
|
Abstract:
We propose a fast solver for the linear system resulting from the
application of a sixth-order Bi-Quartic Spline Collocation method
to the biharmonic Dirichlet problem on a rectangle. The fast
solver is based on Fast Fourier Transforms (FFTs) and
preconditioned GMRES (PGMRES) and has complexity $O(n^2\log(n))$
on an $n \times n$ uniform partition. The FFTs are applied to two
auxiliary problems with different boundary conditions on the two
vertical sides of the domain, while the PGMRES is applied to a
problem related to the two vertical boundaries. We show that the
number of PGMRES iterations required to reduce the relative
residual to a certain tolerance is independent of the gridsize $n$.
Numerical experiments verify the effectiveness of the solver.
Title: | Optimal Quadratic and Cubic Spline Collocation on Nonuniform Partitions
|
Authors: | Christina C. Christara, and Kit Sun Ng
|
Date: | 2006
|
Pages: | 31
|
Note: | Appeared in Computing, 76:3-4, 2006, pp. 227-257.
|
File: | SCnonuni.pdf,
DOI
|
Keywords: |
spline collocation; second-order two-point boundary value problem;
error bounds; optimal order of convergence; spline interpolation.
|
Abstract:
We develop optimal quadratic and cubic spline collocation methods
for solving linear second-order two-point boundary value problems
on non-uniform partitions.
To develop optimal non-uniform partition methods,
we use a mapping function from uniform to non-uniform partitions
and develop expansions of the error at the non-uniform collocation
points of some appropriately defined spline interpolants.
The existence and uniqueness of the spline collocation approximations
are shown, under some conditions.
Optimal global and local orders of convergence of the
spline collocation approximations
and derivatives are derived,
similar to those of the respective methods for uniform partitions.
Numerical results on a variety of problems,
including a boundary layer problem, and a non-linear problem,
verify the optimal convergence of the methods,
even under more relaxed conditions than those assumed by theory.
Title: | Adaptive Techniques for Spline Collocation
|
Authors: | Christina C. Christara, and Kit Sun Ng
|
Date: | 2006
|
Pages: | 19
|
Note: | Appeared in Computing, 76:3-4, 2006, pp. 259-277.
|
File: | adapt.pdf,
DOI
|
Keywords: |
spline collocation; second-order two-point boundary value problem;
error bounds; optimal order of convergence; adaptive grid;
grid size estimator; error estimator; spline interpolation.
|
Abstract:
We integrate optimal quadratic and cubic spline collocation methods
for second-order two-point boundary value problems
with adaptive grid techniques, and grid size and error estimators.
Some adaptive grid techniques are based on the construction
of a mapping function that maps uniform to non-uniform points,
placed appropriately to minimize a certain norm of the error.
One adaptive grid technique for cubic spline collocation
is mapping-free and resembles the technique used in COLSYS (COLNEW).
Numerical results on a variety of problems,
including problems with boundary or interior layers,
and singular perturbation problems
indicate that, for most problems,
the cubic spline collocation method requires less computational
effort for the same error tolerance, and has equally
reliable error estimators,
when compared to Hermite piecewise cubic collocation.
Comparison results with quadratic spline collocation
are also presented.
Title: | GIA-induced secular variations in the Earth's long wavelength
gravity field: Influence of 3-D viscosity variations
(short communication),
|
Authors: | Konstantin Latychev, Jerry X. Mitrovica, Mark E. Tamisiea,
Jeroen Tromp, Christina C. Christara and Robert Moucha
|
Date: | November 2005
|
Pages: | 6
|
Note: | Published in Earth and Planetary Science Letters,
Vol. 240, Issue 2, Pages 322-327
|
File: | DOI link ,
journal page
|
Keywords: |
geopotential harmonics, glacial isostatic adjustment, 3-D structure
|
Abstract:
Predictions of present day secular variations in the Earth's
long wavelength geopotential driven by glacial isostatic adjustment (GIA)
have previously been analyzed to infer the radial profile of
mantle viscosity and to constrain ongoing cryospheric mass balance.
These predictions have been based on spherically symmetric Earth models.
We explore the impact of lateral variations in mantle viscosity
using a new finite-volume formulation for computing the response
of 3-D Maxwell viscoelastic Earth models.
The geometry of the viscosity field is constrained from seismic-to-mographic
images of mantle structure, while the amplitude of the lateral viscosity
variations is tuned by a free parameter in the modeling.
We focus on the zonal J_l harmonics for degrees l = 2, ..., 8
and demonstrate that large-scale lateral viscosity variations
of two to three orders of magnitude have a modest, 5-10%,
impact on predictions of J_2. In contrast, predictions of
higher degree harmonics show a much greater sensitivity
to lateral variation in viscosity structure. We conclude that
future analyses of secular trends (for degree l > 2) estimated
from ongoing (GRACE, CHAMP) satellite missions must incorporate
GIA predictions based on 3-D viscoelastic Earth models.
Title: | Glacial isostatic adjustment on 3-D Earth models:
a finite-volume formulation
|
Authors: | Konstantin Latychev, Jerry X. Mitrovica, Jeroen Tromp,
Mark E. Tamisiea, Dimitri Komatitsch, Christina C. Christara
|
Date: | May 2005
|
Pages: | 24
|
Note: | Published in Geophysical Journal International,
Vol. 161, Issue 2, Pages 421-444
|
File: | DOI link ,
journal page
|
Keywords: |
crystal motions, finite volumes, glacial isostatic adjustment,
3-D viscoelastic Earth models, parallel computation,
domain decomposition
|
Abstract:
We describe and present results from a finite-volume (FV)
parallel computer code for forward modelling the Maxwell viscoelastic response
of a 3-D, self-gravitating, elastically compressible Earth
to an arbitrary surface load. We implement a conservative,
control volume discretization of the governing equations
using a tetrahedral grid in Cartesian geometry and a low-order,
linear interpolation. The basic starting grid honours all major
radial discontinuities in the Preliminary Reference Earth Model (PREM),
and the models are permitted arbitrary spatial variations in viscosity
and elastic parameters. These variations may be either continuous
or discontinuous at a set of grid nodes forming a 3-D surface
within the (regional or global) modelling domain.
In the second part of the paper, we adopt the FV methodology
and a spherically symmetric Earth model to generate a suite
of predictions sampling a broad class of glacial isostatic adjustment (GIA)
data types (3-D crustal motions, long-wavelength gravity anomalies).
These calculations, based on either a simple disc load history
or a global Late Pleistocene ice load reconstruction (ICE-3G),
are benchmarked against predictions generated using the traditional
normal-mode approach to GIA. The detailed comparison provides a guide
for future analyses (e.g. what grid resolution is required to obtain
a specific accuracy?) and it indicates that discrepancies in predictions
of 3-D crustal velocities less than 0.1 mm yr^{-1} are generally obtainable
for global grids with ~3 x 10^6 nodes; however, grids of higher resolution
are required to predict large-amplitude (>1 cm yr^{-1}) radial velocities
in zones of peak post-glacial uplift (e.g. James bay) to the same level
of absolute accuracy. We conclude the paper with a first application
of the new formulation to a 3-D problem. Specifically, we consider
the impact of mantle viscosity heterogeneity on predictions of
present-day 3-D crustal motions in North America. In these tests,
the lateral viscosity variation is constructed, with suitable scaling,
from tomographic images of seismic S-wave heterogeneity, and
it is characterized by approximately 2 orders of magnitude (peak-to-peak)
lateral variations within the lower mantle and 1 order of magnitude
variations in the bulk of the upper mantle (below the asthenosphere).
We find that the introduction of 3-D viscosity structure has a profound
impact on horizontal velocities; indeed, the magnitude of the perturbation
(of order 1 mm yr^{-1}) is as large as the prediction generated
from the underlying (1-D) radial reference model and it far exceeds
observational uncertainties currently being obtained from space-geodetic
surveying. The relative impact of lateral viscosity variations
on predicted radial motions is significantly smaller.
Title: | Optimal Quadratic Spline Collocation on Non-Uniform Partitions
|
Authors: | Christina C. Christara, and Kit Sun Ng
|
Date: | June 2003, revised 2004
|
Pages: | 28
|
Note: | A revised version appeared in Computing.
|
File: |
ccc/QSCnonuni.ps.gz
|
Keywords: |
spline collocation; second-order two-point boundary value problem;
error bounds; optimal order of convergence; adaptive grid;
grid size estimator; error estimator; spline interpolation.
|
Abstract:
Quadratic and Cubic Spline Collocation (QSC and CSC) methods
of optimal orders of convergence have been developed for the solution
of linear second-order two-point Boundary Value Problems (BVPs)
discretized on uniform partitions.
In this paper, we extend the optimal QSC methods to non-uniform
partitions, and in a companion paper
we do the same for CSC.
To develop optimal non-uniform partition QSC methods,
we use a mapping function from uniform to non-uniform
partitions and develop expansions of the error at the non-uniform collocation
points of some appropriately defined spline interpolants.
The existence and uniqueness of the QSC approximations
are shown, under some conditions.
The $j$th derivative of the QSC approximation
is $O(h^{3-j})$ globally, and
$O(h^{4-j})$ locally on certain points, for $j \ge 0$.
Title: | Optimal Cubic Spline Collocation on Non-Uniform Partitions
|
Authors: | Christina C. Christara, and Kit Sun Ng
|
Date: | June 2003, revised 2004
|
Pages: | 17
|
Note: | A revised version appeared in Computing.
|
File: |
ccc/CSCnonuni.ps.gz
|
Keywords: |
spline collocation; second-order two-point boundary value problem;
error bounds; optimal order of convergence; adaptive grid;
grid size estimator; error estimator; spline interpolation.
|
Abstract:
In a recent paper,
non-uniform partition
Quadratic Spline Collocation (QSC) methods
of optimal order of convergence were developed for
second-order two-point Boundary Value Problems (BVPs).
In this paper, we develop optimal non-uniform
Cubic Spline Collocation (CSC) methods.
The formulation and analysis of the methods
are based, as in the QSC case, on
a mapping function from uniform to non-uniform partitions.
However, the CSC method is turned into a mapping-free method.
The $j$th derivative of the CSC approximation
is $O(h^{4-j})$ globally, for $j \ge 0$, and
$O(h^{5-j})$ locally on certain points, for $j > 0$.
The non-uniform partition CSC method is integrated with
adaptive grid techniques, and grid size and error estimators.
One adaptive grid technique is based on the construction
of a mapping function. Another is mapping-free
and resembels the technique used in COLSYS (COLNEW) \cite{AMR95}.
Numerical results on a variety of problems,
including problems with boundary or interior layers,
and singular perturbation problems
indicate that, for most problems,
the CSC method require less computational
effort for the same error tolerance, and has equally
reliable error estimators,
when compared to Hermite piecewise cubic collocation (HPCC).
Title: | Quadratic Spline Galerkin Method for the Shallow Water Equations
|
Authors: | Anita T. Layton, Christina C. Christara, and Kenneth R. Jackson
|
Date: | November 2002
|
Pages: | 30
|
Note: | Mathematics and Computers in Simulation, 71:3, 2006, pp. 175-186
|
File: |
sweqsg02.ps.gz
|
Keywords: | numerical weather prediction; finite element;
semi-Lagrangian semi-implicit; Rossby stability; staggered grids
|
Abstract:
Currently in most global meteorological applications,
the spectral transform method or low-order
finite difference/finite element methods are used.
The spectral transform method, which yields high-order approximations,
requires Legendre transforms.
The Legendre transforms have a computational complexity of
$\mathcal O(N^3)$, where $N$ is the number of subintervals in one dimension,
and thus render the spectral transform method unscalable.
In this study, we present an alternative numerical method for solving
the shallow water equations (SWEs) on a sphere in spherical coordinates.
In this implementation, the SWEs are discretized in time using the
two-level semi-Lagrangian semi-implicit method, and in space on
staggered grids using the quadratic spline Galerkin method.
We show that, when applied to a simplified version of the
SWEs, the method yields a neutrally stable solution for the
meteorologically significant Rossby waves.
Moreover, we demonstrate that the Helmholtz equation arising from
the discretization and solution of the
SWEs should be derived algebraically rather than analytically, in
order for the method to be stable with respect to the Rossby waves.
These results are verified numerically using Boyd's equatorial wave equations
\cite{boyd1} with initial conditions chosen to generate a soliton.
Title: | Optimal Quadratic Spline Collocation Methods for the
Shallow Water Equations
|
Authors: | Anita T. Layton, Christina C. Christara, and Kenneth R. Jackson
|
Date: | November 2002
|
Pages: | 38
|
Note: | Mathematics and Computers in Simulation, 71:3, 2006, pp. 187-205
|
File: |
sweqsc02.ps.gz
|
Keywords: | numerical weather prediction; finite element;
semi-Lagrangian semi-implicit; Rossby stability; staggered grids
|
Abstract:
In this study, we present numerical methods, based on the optimal
quadratic spline collocation (OQSC) methods,
for solving the shallow water equations (SWEs) in spherical coordinates.
A quadratic spline collocation method approximates the solution of a
differential problem by a quadratic spline.
In the standard formulation, the quadratic spline is computed by making the
residual of the differential equations zero at a set of collocation points;
the resulting error is second order, while the error associated with
quadratic spline interpolation is fourth order locally at certain points
and third order globally.
The OQSC methods generate approximations of the same order as
quadratic spline interpolation.
In the one-step OQSC method, the discrete differential operators are
perturbed to eliminate low-order error terms, and a high-order
approximation is computed using the perturbed operators.
In the two-step OQSC method, a second-order approximation is generated
first, using the standard formulation, and then a high-order
approximation is computed in a second phase by perturbing the right
sides of the equations appropriately.
In this implementation, the SWEs are discretized in
time using the semi-Lagrangian semi-implicit scheme, which allows large
timesteps while maintaining numerical stability, and in space using the
OQSC methods.
The resulting methods are efficient and yield stable and accurate
representation of the meteorologically important Rossby waves.
Moreover, by adopting the Arakawa C-type grid, the methods also
faithfully capture the group velocity of inertia-gravity waves.
Title: | Quadratic Spline Collocation Revisited:
Extension to Systems of Elliptic PDEs,
|
Authors: | Christina C. Christara and Kit Sun Ng
|
Date: | July 2002
|
Pages: | 32
|
Note: | Depart. of Computer Science, University of Toronto,
Technical Report 318/2001, rev. July 2002.
|
File: |
ccc/sys.ps.gz
|
Abstract:
We revisit the optimal Quadratic Spline Collocation (QSC) methods
for single elliptic PDEs developed in \cite{Chri94}
and prove some new results, including a discrete maximum principle.
We extend QSC to
systems of two coupled linear second-order PDEs in two dimensions.
The system of PDEs is treated as one entity; no decoupling is applied.
We study the matrix properties of the linear system arising
from the discretization of systems of two PDEs by QSC.
We give sufficient conditions under which the QSC linear system
is uniquely solvable and
develop analytic formulae for the eigenvalues and eigenvectors
of the QSC matrix from $2 \times 2$ Helmholtz operators with
constant coefficients.
Numerical results demonstrate that the QSC methods are optimal,
that is, fourth order locally at certain points and third order globally,
even for some problems that do not satisfy the conditions of the analysis.
The QSC methods are compared to conventional second-order discretization
methods and are shown to produce smaller approximation errors
in the same computation time, while they achieve the same accuracy
in less time.
The formulation of the optimal two-step QSC method for
a $2 \times 2$ system of PDEs can be extended in a straightforward way to
a $n \times n$ system of PDEs.
Title: | Fast Fourier Transform Solvers and Preconditioners
for Quadratic Spline Collocation
|
Authors: | Christina C. Christara and Kit Sun Ng
|
Date: | November 2000
|
Pages: | 33
|
Note: | Depart. of Computer Science, University of Toronto,
Technical Report 317/2000.
A revised version to appear in BIT, 42:4 (December 2002).
|
File: | ccc/precond.pdf,
ccc/precond.ps.gz
|
DOI: | 10.1023/A:1021944218806
|
Journal: | PDF file (fulltext),
page
|
Abstract:
Quadratic Spline Collocation (QSC) methods of optimal order of convergence
have been recently developed for the solution of elliptic Partial
Differential Equations (PDEs). In this paper, linear solvers based on
Fast Fourier Transforms (FFT) are developed for the solution of the QSC
equations. The complexity of the FFT solvers is O(N^2 log N),
where N is the gridsize in one dimension.
These direct solvers can handle PDEs with coefficients
in one variable or constant, and Dirichlet, Neumann and
periodic boundary conditions. General variable coefficient PDEs
are handled by preconditioned iterative solvers. The preconditioner
is the QSC matrix arising from a constant coefficient PDE.
The convergence analysis of the preconditioner is presented.
It is shown that, under certain conditions, the convergence rate
is independent of the gridsize.
The preconditioner is solved by FFT techniques,
and integrated with one-step or acceleration methods,
giving rise to asymptotically almost optimal linear solvers,
with complexity O(N^2 log N).
Numerical experiments verify the effectiveness of the solvers
and preconditioners, even on problems more general than the
analysis assumes.
The development and analysis of FFT solvers and preconditioners is extended
to QSC equations corresponding to systems of elliptic PDEs.
Title: | An Efficient Transposition Algorithm for Distributed Memory Computers
|
Authors: | C. Christara, X. Ding and K. R. Jackson
|
Date: | October 1999
|
Pages: | 19
|
Note: | Proceedings of the 13th Annual International
Symposium on High Performance Computing Systems and Applications
(HPCS'99), Queen's University, Kingston, Ontario, Canada, 13-16 June 1999.
|
File: | transpose.99b.ps.gz
|
Abstract:
Data transposition is required in many numerical applications. When
implemented on a distributed-memory computer, data transposition
requires all-to-all communication, a time consuming operation. The Direct
Exchange algorithm, commonly used for this task, is inefficient if the
number of processors is large. We investigate a series of more
sophisticated techniques: the Ring Exchange, Mesh Exchange and Cube
Exchange algorithms. These data transposition schemes were
incorporated into a parallel solver for the shallow-water equations. We
compare the performance of these schemes with that of the Direct
Exchange Algorithm and the MPI all-to-all communication routine,
MPI_AllToAll. The numerical experiments were performed on a Cray T3E
computer with 512 processors and on an ethernet-connected cluster of
36 Sun workstations. Both the analysis and the numerical results indicate
that the more sophisticated Mesh and Cube Exchange algorithms perform
better than either the simpler well-known Direct and Ring Exchange
schemes or the MPI_AllToAll routine. We also generalize the Mesh and
Cube Exchange algorithms to a d-dimensional mesh algorithm, which can
be viewed as a generalization of the standard hypercube data
transposition algorithm.
Title: | Multigrid and Multilevel Methods for Quadratic Spline Collocation
|
Authors: | Christina C. Christara and Barry Smith
|
Date: | February 1997
|
Pages: | 21
|
Note: | Appeared in BIT 37:4 (1997), pp. 781-803.
|
File: | multign.pdf,
multign.ps.gz
|
DOI: | 10.1007/BF02510352
|
Journal: | PDF file (fulltext),
page
|
Abstract:
Multigrid methods are developed and analyzed
for quadratic spline collocation equations arising
from the discretization of one-dimensional second-order
differential equations.
The rate of convergence of the two-grid method
integrated with a damped Richardson relaxation scheme
as smoother is shown to be faster than 1/2,
independently of the step-size.
The additive multilevel versions of the
algorithms are also analyzed.
The development of quadratic spline collocation
multigrid methods is extended to two-dimensional
elliptic partial differential equations.
Multigrid methods for quadratic spline collocation methods
are not straightforward: because the basis functions used
with quadratic spline collocation are not nodal basis functions,
thus the design of efficient restriction and extension
operators is nontrivial.
Experimental results, with V-cycle and full multigrid,
indicate that suitably chosen multigrid
iteration is a very efficient solver for the
quadratic spline collocation equations.
Title: | Parallel solvers for spline collocation equations
|
Author: | C. C. Christara
|
Date: | January 1995
|
Pages: | 35
|
Note: | Appeared in Advances in Engineering Software, 27:1/2, 1996, pp. 71-89
|
File: | ccc/civil.ps.gz
|
Abstract:
A variety of solvers for the spline collocation equations
arising from the discretisation of elliptic Partial Differential
Equations (PDEs) are considered.
The convergence properties of semi-iterative
and Krylov subspace acceleration methods
applied to the system of spline collocation equations are studied.
The preconditioners tested include incomplete factorisation
and SSOR for both the natural and multicolour orderings,
domain decomposition based on
Schur complement methods with nonoverlapping subdomains,
or Schwarz methods with overlapping subdomains,
and multigrid methods.
The parallelisation of some of the above
iterative methods is studied and their advantages
and disadvantages discussed.
The communication requirements of the methods are discussed
when the methods are implemented on distributed memory machines.
Results which show that spline collocation methods
are very competitive for the solution of PDEs
are presented.
Title: | Quadratic spline collocation methods for
elliptic partial differential equations
|
Author: | C. C. Christara
|
Date: | March 1994 (replaces Tech. Rep. DCS-135, 1990)
|
Pages: | 34
|
Note: | A slightly revised version of this paper appeared in BIT, 34:1, 1994, pp. 33-61
|
File: | ccc/cs-90-135.ps.gz
|
DOI: | 10.1007/BF01935015
|
Journal: | PDF file (fulltext),
page
|
Abstract:
We consider Quadratic Spline Collocation (QSC) methods for linear
second order elliptic Partial Differential Equations (PDEs).
The standard formulation of these methods leads to non-optimal approximations.
In order to derive optimal QSC approximations, high order perturbations
of the PDE problem are generated. These perturbations can be applied
either to the PDE problem operators or to the right sides,
thus leading to two different formulations of optimal QSC methods.
The convergence properties of the QSC methods are studied.
Optimal $O(h sup 3-j )$ global error estimates for the $j$-th partial
derivative are obtained for a certain class of problems. Moreover,
$O(h sup 4-j )$ error bounds for the $j$-th partial derivative are
obtained at certain sets of points.
Results from numerical experiments verify the theoretical behaviour
of the QSC methods. Performance results also show that the QSC methods
are very effective from the computational point of view.
The QSC methods have been implemented efficiently
on parallel machines.
Title: | High performance computing of elliptic partial differential equations
with spline collocation
|
Author: | C. C. Christara
|
Date: | June 1994
|
Pages: | 13
|
Note: | Appeared in Proceedings of the
1994 Computational Structures Technology conference (CST94),
August 30 - September 1, 1994, Athens, Greece, sponsored by Civil-Comp,
Vol. A, pp. 23-34
|
File: | ccc/highperf.ps.gz
|
Abstract:
We consider a variety of solvers for the spline collocation equations
arising from the discretisation of elliptic Partial Differential Equations.
We study the convergence properties of semi-iterative
and conjugate gradient acceleration methods
applied to the system of spline collocation equations.
The preconditioners tested include incomplete factorisation
and SSOR for both the natural and multicolour orderings,
domain decomposition based on
Schur complement methods with nonoverlapping subdomains,
or Schwarz methods with overlapping subdomains,
and multigrid methods.
We study the parallelisation of some of the above
iterative methods and discuss their advantages
and disadvantages.
The communication requirements of the methods are discussed
when the methods are implemented on distributed memory machines.
Results which show that spline collocation methods
are very competitive for the solution of BVPs are presented.
Title: | Parallel computation of partial differential equations
on distributed memory machines
|
Author: | C. C. Christara
|
Date: | June 1992
|
Pages: | 13
|
Note: | Appeared in Advances in Computer Methods for
Partial Differential Equations - VII, pp. 142-149
R. Vichnevetsky, D. Knight, G. Richter eds.
|
File: | ccc/paral.ps.gz
|
Abstract:
We review a number of parallel methods for the solution
of elliptic Boundary Value Problems (BVPs).
The methods considered are of coarse or medium grain parallelism
and fixed grid discretisation.
The cases where the degree of parallelism grows
linearly with the size of the problem in all its
dimensions or in fewer dimensions are considered.
We also consider the parallelism in both the discretisation of the BVP phase
and the solution of the resulting linear system phase.
We discuss the assignment of data to processors
and the communication requirements of the methods,
when implemented on distributed memory machines.
Based on a general model for a distributed memory machine,
we develop a framework in order to classify
the communication requirements of the methods considered.
We study the communication complexity of the methods
on a general distributed machine model,
and present performance results on the iPSC/2 hypercube.
Title: | Multicolour orderings and iterative methods for elliptic equations
|
Author: | C. C. Christara
|
Date: | June 1991
|
Pages: | 11
|
Note: | Appeared in Proceedings of the Supercomputing Symposium '91 (SS '91),
June 3-5, 1991, Fredericton, NB, Canada, pp. 221-230
|
File: | ccc/colour.ps.gz
|
Abstract:
In this paper we study the performance and parallel implementation
of iterative methods applied to linear systems arising from the
discretisation of linear elliptic Boundary Value Problems (BVPs)
which are ordered according to multicolour orderings.
We view the discretisation of a BVP as a process of
generating a set of equations between stencils.
These equations, when ordered according to some indexing of the stencil,
result in a sparse linear system.
We discuss the effect of indexing of the stencils
on the sparsity pattern of the linear system,
on the convergence of iterative methods applied to the system,
as well as on the parallel implementation of these methods.
We discuss which stencil indexings are preferable in order to minimise
communication in message-passing parallel architectures.
Title: | Domain decomposition and incomplete factorisation methods
for partial differential equations
|
Author: | C. C. Christara
|
Date: | May 1991
|
Pages: | 11
|
Note: | Appeared in Proceedings of the sixth Distributed Memory Computing
Conference (DMCC6), April 28 - May 1, 1991, Portland, OR, U.S.A., pp. 369-377
|
File: | ccc/iccg.ps.gz
|
Abstract:
In this paper we develop and study a method which tries to combine
the merits of Domain Decomposition (DD)
and Incomplete Cholesky preconditioned Conjugate Gradient method (ICCG)
for the parallel solution of linear elliptic
Partial Differential Equations (PDEs) on rectangular domains.
We first discretise the PDE problem, using Spline Collocation,
a method of Finite Element type based on smooth splines.
This gives rise to a sparse linear system of equations.
The ICCG method provides us with a very efficient,
but not straightforward parallelisable
linear solver for such systems.
On the other hand, DD methods
are very effective for elliptic PDEs.
A combination of DD and ICCG methods, in which
the subdomain solves are carried out with ICCG,
leads to efficient and highly parallelisable solvers.
We implement this hybrid DD\-ICCG method on a hypercube,
discuss its parallel efficiency, and show results from
experiments on configurations with up to 32 processors.
We apply a totally local communication scheme
and discuss its performance on the iPSC/2 hypercube.
A similar approach can be used with other PDE discretisation methods.
Title: | Schur complement preconditioned conjugate gradient methods for
spline collocation equations
|
Author: | C. C. Christara
|
Date: | May 1991
|
Pages: | 14
|
Note: | An earlier version of this paper appeared in
Computer architecture news, 18:3, 1990, pp. 108-120, and
Proceedings of the 1990 International Conference on Supercomputing,
June 1990, Amsterdam, the Netherlands, sponsored by ACM
|
File: | ccc/scmpcg.ps.gz
|
Abstract:
We are interested in the efficient solution of linear second order
Partial Differential Equation (PDE) problems on rectangular domains.
The PDE discretisation scheme used is of Finite Element type
and is based on quadratic splines and the collocation methodology.
We integrate the Quadratic Spline Collocation (QSC) discretisation scheme
with a Domain Decomposition (DD) technique.
We develop DD motivated orderings of the QSC equations and unknowns
and apply the Preconditioned Conjugate Gradient (PCG) method
for the solution of the Schur Complement (SC) system.
Our experiments show that the SC-PCG-QSC method
in its sequential mode is very efficient compared to standard
direct band solvers for the QSC equations.
We have implemented the SC-PCG-QSC method on the iPSC/2 hypercube
and present performance evaluation results for up to 32 processors
configurations.
We discuss a type of nearest neighbour communication scheme,
in which the amount of data transfer per processor does not
grow with the number of processors.
The estimated efficiencies are at the order of 90%.
Title: | Conjugate gradient methods for spline collocation equations
|
Author: | C. C. Christara
|
Date: | May 1991
|
Pages: | 9
|
Note: | An earlier version of this paper appeared in Proceedings of the
fifth Distributed Memory Computing Conference (DMCC5),
April 8-12, 1990, Charleston, SC, U.S.A., pp. 550-558
|
File: | ccc/pcgcoef.ps.gz
|
Abstract:
We study the parallel computation of linear second order elliptic
Partial Differential Equation (PDE) problems in rectangular domains.
We discuss the application of Conjugate Gradient (CG)
and Preconditioned Conjugate Gradient (PCG) methods to the linear
system arising from the discretisation of such problems using
quadratic splines and the collocation discretisation methodology.
Our experiments show that the number of iterations
required for convergence of CG-QSC
(Conjugate Gradient applied to Quadratic Spline Collocation equations)
grows linearly with the square root of the number of equations.
We implemented the CG and PCG methods for the solution of the
Quadratic Spline Collocation (QSC) equations on the iPSC/2 hypercube
and present performance evaluation results
for up to 32 processors configurations.
Our experiments show efficiencies of the order of 90%,
for both the fixed and scaled speedups.
Title: |
A Domain Decomposition Spline Collocation Method for Elliptic Partial
Differential Equations
|
Author(s): | C. C. Christara and E. N. Houstis
|
Date: | March 1989
|
Pages: | 7
|
Note: | Proceedings of the fourth Conference on Hypercubes, Concurrent Computers
and Applications (HCCA4),
March 6-8, 1989, Monterey, CA, U.S.A., pp. 1267-1273.
|
Abstract:
We consider the parallel solution of a discretization
system obtained by approximating the solution of the second order
elliptic boundary value problem
Lu == a uxx + b uxy + c uyy + d ux + e uy + fu = g in OMEGA == [ax,bx] X [ay,by]
(1.1)
Bu == alpha u + beta un = g_0 on partial OMEGA == boundary of OMEGA
(1.2)
using the quadratic spline collocation method.
The input functions in the partial differential equation (PDE) problem (1.1),
(1.2) are assumed to be functions of $x,y$ in
$C^1 [OMEGA]$.
In [Hous88], we have studied the computational performance of various
direct and iterative methods for solving such discrete systems on sequential
machines.
An overview of spline collocation methods is presented in Section 2.
The objective of this paper is to develop and analyze parallel solvers on
shared and non-shared memory MIMD architectures for solving elliptic
quadratic spline collocation equations.
The class of domain decomposition or substructuring methods has recently been
applied by many researchers to solve efficiently PDE discretization
systems.
Their build in natural parallelism can be easily mapped to various MIMD
architectures.
Our focus in this study is limited to a parallel formulation of the
capacitance method for spline collocation equations.
According to the basic idea of this method, the computational effort of
solving the original system is reduced to the solution of the so called
capacitance matrix problem.
For positive definite systems, the capacitance system is solved with
preconditioning conjugate gradient (PCG).
In the case of quadratic spline collocation equations,
in addition to the PCG method,
we have developed a scaled Jacobi scheme for its solution,
which is formulated and analyzed in Section 3.
The parallelization of the capacitance - Jacobi solver, is presented in
Section 4 together with its complexity and performance analysis on a
hypercube and bus architectures.
Finally, the performance on the MIMD machines (NCUBE/7 and SEQUENT) are
reported in Section 5.
Title: |
A Parallel Spline Collocation - Capacitance Method for Elliptic Partial
Differential Equations
|
Author(s): | C. C. Christara, E. N. Houstis and J. R. Rice
|
Date: | July 1988
|
Pages: | 10
|
Note: |
Proceedings of the 1988 International Conference on Supercomputing (ICS88),
July 1988, Saint-Malo, France, pp. 375-384.
|
Abstract:
We consider the integration of a domain decomposition technique with a
new quadratic spline collocation discretization scheme for solving
second order elliptic boundary value problems on rectangles.
The domain decomposition method is based on the capacitance matrix
technique. Due to the limitations of existing methods for solving the
corresponding capacitance
problem, we develop and analyze iterative methods for its solution.
The optimum partitioning and mapping of the underlying computation is
studied on hypercube architectures.
A numerical realization of this method is presented on NCUBE/7 (128
processors) and its comparative efficiency is measured.
The resulting parallel quadratic spline collocation-capacitance
method is seen to be efficient in achieving accurate solutions and
in using parallel architectures.
Title: |
Geometry Decomposition based Methods for solving Elliptic Partial
Differential Equations
|
Author(s): |
C. C. Christara, A. Hadjidimos, E. N. Houstis, J. R. Rice and E. A. Vavalis
|
Date: | May 1988
|
Pages: | 8
|
Note: |
Comp. Methods in Flow Analysis, (H. Niki and M. Kawahara, eds),
1988, Univ. of Okayama, Japan, pp. 175-182.
|
Abstract:
Title: |
Quadratic Spline Collocation Methods for two point boundary value problems
|
Author(s): | E. N. Houstis, C. C. Christara and J. R. Rice
|
Date: | April 1988
|
Pages: | 18
|
Note: |
International Journal for Numerical Methods in Engineering, 26:4,
April 1988, pp. 935-952.
|
Abstract:
A new collocation method based on quadratic splines is presented for second
order two point boundary value problems.
First, $O(h^4)$ approximations to the first and second derivative
of a function are derived using a quadratic spline interpolant of $u$.
Then these approximations are used to define an $O(h^4)$
perturbation of the given boundary value problem.
Second, the perturbed problem is used to define a collocation
approximation at interval midpoints
for which an optimal $O(h^{3-j})$ global estimate for the $j$th
derivative of the error is derived.
Further, $O(h^{4-j})$ error bounds for the $j$th derivative are
obtained for certain superconvergence points.
It should be observed that standard collocation at midpoints gives
$O(h^{2-j})$ bounds.
Results from numerical experiments are reported that verify the
theoretical behavior of the method.
Title: |
Performance of Scientific Software
|
Author(s): |
E. N. Houstis, J. R. Rice, C. C. Christara and E. A. Vavalis
|
Date: | 1988
|
Pages: | 23
|
Note: |
Mathematical Aspects of Scientific Software,
The IMA Volumes in Mathematics and its applications, Vol 14,
Springer Verlag, 1988, pp. 123-155.
|
Abstract:
We review performance methodologies used for the evaluation of scientific
software in von Neumann architectures. A prototype evaluation facility for
second
order elliptic partial differential equation (PDE) solvers is described which
realizes the main objectives of these methodologies. Finally, the results
of an evaluation study for a new class of spline collocation solvers for
elliptic PDEs are presented.
Title: |
Scientific Computing by Numerical Methods
|
Author(s): | C. C. Christara and K. R. Jackson
|
Date: | 1996
|
Pages: | 121
|
File: | survey.pdf
|
Note: |
A different version of this article appeared in
Encyclopaedia of Applied Physics, Vol 17, 1996, pp. 1-79,
and in Mathematical Tools for Physicists.
|
Abstract:
Numerical methods are an indispensable tool in solving many problems
that arise in applied physics.
In this article, we briefly survey a few of the most common mathematical
problems and review some numerical methods to solve them.
As can be seen from the table of contents,
the topics covered in this survey are those that appear
in most introductory numerical methods books.
However, our discussion of each topic is more brief
than is normally the case in such texts
and we frequently provide references to more advanced topics that
are not usually considered in introductory books.
We originally intended to include a section on the computation of
eigenvalues and eigenvectors of matrices,
but, because of time and space constraints,
we omitted it.
The interested reader should consult an
introductory numerical methods book,
such as those listed in the last section,
or an advanced text, such as
\cite{Golub-Van-Loan, Parlett, Saad, Stewart, Wilkinson}.
Supervised students' theses
Title: |
Bilateral XVA pricing under stochastic default intensity:
PDE modelling and computation,
|
Author: | Yuwei Chen
|
Date: | February 2023
|
Pages: | 123
|
Note: | Ph.D. thesis, Computer Science Department,
University of Toronto
|
File: |
ywchen-23-phd.pdf
|
Abstract:
Adjusting derivative prices to take into account default risk
has attracted the attention of several researchers and practitioners,
especially after the 2007-2008 financial crisis.
The thesis is a study, via numerical Partial Differential Equation (PDE)
approaches, of the modeling and computation
of valuation adjustment problems in financial derivative pricing
if we consider the bilateral counterparty default intensities.
Under some assumptions, the valuation of financial derivatives,
including a value adjustment to account for bilateral default risk
(XVA), assuming constant default intensities
gives rise to a nonlinear PDE \cite{burgard2011partial}.
We formulate a penalty iteration method to handle the nonlinearity
in this PDE, study its convergence, and
extend it to American derivatives,
resulting in a double-penalty iteration.
We also propose boundary conditions
and their discretization for the XVA PDE problem
in the cases of call, put options, and forward contracts.
Numerical results demonstrate the effectiveness of our methods.
Then, we derive a novel PDE
for derivative pricing including XVA,
assuming that the default risk of one of
the counterparties follows a
Cox-Ingersoll-Ross (CIR) process, while the other party
has constant default risk.
The time-dependent PDE derived is of Black-Scholes type and
involves two ``space'' variables, namely the asset price
and the buyer default intensity,
as well as a nonlinear source term.
We formulate boundary conditions appropriate
for the default intensity variable.
The numerical solution of the PDE uses a penalty-like iteration
for handling the nonlinearity.
We also develop a novel asymptotic approximation formula
for the adjusted price of derivatives,
resulting in a very efficient, accurate, and convenient
for practitioners formula.
Numerical results indicate
stable second order convergence for the 2D PDE solution
in terms of the discretization size and convergence of order at least
1.5 for the asymptotic approximation in terms of inverse of
the mean-reversion rate.
We extend the PDE model we developed to price European derivatives
including XVA, considering stochastic counterparty default intensity,
to American derivatives.
We also extend the asymptotic approximation
to the American put option problem.
A key step is to derive the asymptotic approximation to
the free boundary.
We present numerical experiments in order to study
the accuracy and effectiveness of the 2D PDE and asymptotic approximations.
Title: | Penalty Methods for Nonlinear Problems
in Financial Option Pricing
|
Authors: | Ruining (Ray) Wu
|
Date: | February 2021 |
Pages: | 121 |
Note: | MSc Research Paper, Computer Science Department,
University of Toronto |
File: | rwu-21-msc.pdf
|
Abstract:
We study the pricing of nonlinear problems in computational finance by
numerical Partial Differential Equation (PDE) methods. We examine the numerical
solution obtained via the use of Policy Iteration methods from
Hamilton-Jacobi-Bellman equations and introduce new penalty-like iterative
methods to solve the nonlinear equations. We consider the following problems:
Unequal Borrowing and Lending rates, Stock Borrowing Fees, Stock Borrowing Fees
with American exercise rights, Uncertain Volatility, and Transaction cost
models. We demonstrate second-order convergence of the solution and of the
Greeks, and for the nonlinear problems we study the number of nonlinear
iterations taken as a proxy to the total computational cost. Furthermore, the
effects of various details such as different boundary conditions and nonuniform
discretization grids are studied. Finally, where applicable, we prove
monotonicity of the new penalty-like methods that we have introduced and
demonstrate that our numerical implementations uphold this property.
Title: |
Numerical Methods for Pricing Multi-Asset Options
|
Author: | Yuwei Chen
|
Date: | January 2018
|
Pages: | 75
|
Note: | M.Sc. thesis, Computer Science Department,
University of Toronto
|
File: |
ywchen-18-msc.pdf
|
Abstract:
We consider the pricing of two-asset European and American options
by numerical Partial Differential Equation (PDE) methods,
and compare the results with certain analytical formulae.
Two cases of options are tested: exchange option and spread option.
For exchange options, the analytical formula considered
is the (exact) Margrabe formula. For spread options,
we consider Kirk's formula and the formula by Li, Deng and Zhou.
In pricing European two-asset options, the basic numerical PDE model
is the two-dimensional Black-Scholes PDE.
Different boundary conditions are considered, and the effect
of them on the solution at various points of the grid is studied.
Furthermore, various types of non-uniform grids are considered,
aiming at reducing the error at certain areas of the grid.
Moreover, the effect of the truncated domain on the
PDE approximation is studied.
We also discuss the effect of certain problem parameters,
such as the length of maturity time, and the values of volatility
and correlation, to the accuracy and convergence of the PDE approximations.
The experiments indicate that the numerical PDE computed price
and Greeks are second-order, for appropriately chosen grid
discretizations.
In the American pricing problem, the discrete penalty method
is applied to the linear complementarity problem involving
the two-dimensional Black-Scholes PDE and additional constraints.
The convergence of the American Spread put option approximation
computed with penalty iteration remains second-order, with
the number of penalty iterations per timestep remaining small (2-3).
We also consider an iterative method with preconditioning techniques
for solving the arising large sparse linear system at each timestep,
and show that this solution technique is asymptotically optimal.
Title: |
PDE option pricing: analysis and application
to stochastic correlation,
|
Author: | Nat Chun-Ho Leung
|
Date: | June 2017
|
Pages: | 119
|
Note: | Ph.D. Thesis, Depart. of Computer Science,
University of Toronto, 2017.
|
File: |
nat-17-phd.pdf
|
Keywords: |
quantization error, Fourier analysis, digital options, call and put options,
stochastic correlation, spread options, quanto options, max options,
partial differential equation, localization, boundary conditions,
convergence, stability.
|
Abstract:
This thesis is a study of numerical Partial Differential Equation (PDE)
methods in financial derivatives pricing. The first part of the thesis
is concerned with the behaviour of a numerical PDE solution when the initial
condition is not smooth. The second part of the thesis develops
computational PDE methods for option pricing problems with stochastic
correlation.
In the first part of this thesis, we provide an analysis of the error
arising from a non-smooth initial condition when solving a pricing
problem with a finite difference method.
We build our framework on the sharp error estimate in Giles and Carter (2006),
and study three types of non-smoothness that are of financial interest.
We show that the error of the numerical solution under
Crank-Nicolson-Rannacher timestepping with central spatial differences
can be decomposed into two components, respectively a second order
error resulting from the approximation to the heat kernel by a
discrete operator, and a quantization error that depends on the
positioning of non-smoothness on the grid. We discuss how this positioning
affects the quality of the numerical solution, and the possibility of
an optimal placement of the non-smoothness in the mesh.
We also study explicitly the effect of
smoothing on the approximation error.
The second part of the thesis focuses on the pricing of European
options using a stochastic correlation model. We derive a
time-dependent PDE for the pricing problem under stochastic
correlation, and develop computational approaches for its solution.
The first approach is a finite difference scheme. We study
such issues as localization of domain, boundary conditions and stability
of the numerical scheme. The second approach
is an asymptotic solution of the PDE, appropriate for cases
when the correlation process exhibits fast mean reversion and when
a numerical PDE solution is considered costly.
Numerical experiments demonstrate the effectiveness of our methods,
and the agreement among the two solutions and Monte Carlo simulations.
We also experimentally demonstrate the effect of smoothing on the
numerical solution, and study the effect of certain problem
parameters on the approximate solution.
Title: |
Efficient and Accurate Numerical PDE Methods for
Pricing Financial Derivatives
|
Author: | Mufan (Bill) Li
|
Date: | April 2015
|
Pages: | 35
|
Note: | B.A.Sc. Thesis, Depart. of Electrical and Computer Eng.,
Univ. of Toronto, 2015.
|
File: |
mufanLi-15-basc.pdf
|
Keywords: |
American options, (implicit) penalty iteration, operator splitting
(explicity penalty) method, convergence, accuracy, complexity, efficiency,
Greek calculation.
|
Abstract:
The main difficulty in pricing American options comes from
the early exercise right,
creating a non-linear constraint on the Black-Scholes PDE.
Under a finite difference discretization of the PDE,
the price of an American can be approximated,
with several techniques to properly handle the American Constraint.
While both an iterative penalty method and
a direct operator splitting method are convergent,
the efficiency and quality require a comprehensive study.
Using Crank-Nicolson time stepping and non-uniform grids,
the methods are compared in numerical experiments.
The criteria include order of convergence, efficiency, and complexity.
Title: |
Modeling Multi-Factor Financial Derivatives by a
Partial Differential Equation Approach with
Efficient Implementation on Graphics Processing Units
|
Author: | Duy Minh Dang
|
Date: | August 2011
|
Pages: | 229
|
Note: | Ph.D. Thesis, Computer Science Dept., Univ. of Toronto, 2011.
|
File: |
Duy_Minh_Dang_PhD_Thesis.pdf
|
Keywords: |
cross-currency interest rate derivatives,
Power Reverse Dual Currency (PRDC) swaps,
Bermudan cancelability, knockout, and FX Target Redemption,
Alternating Direction Implicit (ADI) timestepping methods,
multi-GPU platforms, multi-asset American options.
|
Abstract:
This thesis develops efficient modeling frameworks via a Partial
Differential Equation (PDE) approach for multi-factor financial
derivatives, with emphasis on three-factor models, and studies
highly efficient implementations of the numerical methods on novel
high-performance computer architectures, with particular focus on Graphics
Processing Units (GPUs) and multi-GPU platforms/clusters of GPUs. Two
important classes of multi-factor financial instruments are considered:
cross-currency/foreign exchange (FX) interest rate derivatives and
multi-asset options.
For cross-currency interest rate derivatives, the focus of the thesis
is on Power Reverse Dual Currency (PRDC) swaps with three of the most
popular exotic features, namely Bermudan cancelability, knockout, and FX
Target Redemption. The modeling of PRDC swaps using one-factor Gaussian
models for the domestic and foreign interest short rates, and a one-factor
skew model for the spot FX rate results in a time-dependent parabolic PDE
in three space dimensions. Our proposed PDE pricing framework is based
on partitioning the pricing problem into several independent pricing
subproblems over each time period of the swap's tenor structure, with
possible communication at the end of the time period. Each of these
subproblems requires a solution of the model PDE. We then develop a
highly efficient GPU-based parallelization of the Alternating Direction
Implicit (ADI) timestepping methods for solving the model PDE. To
further handle the substantially increased computational requirements
due to the exotic features, we extend the pricing procedures to multi-GPU
platforms/clusters of GPUs to solve each of these independent subproblems
on a separate GPU. Numerical results indicate that the proposed GPU-based
parallel numerical methods are highly efficient and provide significant
increase in performance over CPU-based methods when pricing PRDC swaps.
An analysis of the impact of the FX volatility skew on the price of PRDC
swaps is provided.
In the second part of the thesis, we develop efficient pricing algorithms
for multi-asset options under the Black-Scholes-Merton framework, with
strong emphasis on multi-asset American options. Our proposed pricing
approach is built upon a combination of (i) a discrete penalty approach
for the linear complementarity problem arising due to the free boundary
and (ii) a GPU-based parallel ADI Approximate Factorization technique for
the solution of the linear algebraic system arising from each penalty
iteration. A timestep size selector implemented efficiently on GPUs is
used to further increase the efficiency of the methods. We demonstrate
the efficiency and accuracy of the proposed GPU-based parallel numerical
methods by pricing American options written on three assets.
Title: |
A study of a discrete element method based granular dynamics solver
|
Author: | Tina Siu Yee
|
Date: | October 2009
|
Pages: | 43
|
Note: | M.Sc. Thesis, Computer Science Dept., Univ. of Toronto, 2009.
|
File: | ccc/tyee-09-msc.pdf
|
Keywords: |
molecular dynamics, sand, powder, particle simulation, ODE solver, Runge-Kutta,
friction, dampening, piling, stiction
|
Abstract:
Granular dynamics is the dynamics of a large set of small particles (grains).
Convincing simulation of natural granular phenomena
(e.g. dust, sand or powders) is a challenging mathematical
and computational problem.
Our observation is that the more realistically the collection
of grains approaches its static state, the more natural
the simulation appears. This study focuses on the simulation
of sets of grains as the set approaches its static state.
The method begins with a discrete element (also referred to
as molecular dynamics) model of the inter-particle contacts
within the granular collection. Inertia terms (friction/dampening)
are added to the basic contact model to facilitate static piling.
An examination of the different contact models on the formation
of the final static state and a discussion of the numerical
consequences of each model is presented.
The discrete element approach demonstrates to be a straightforward
and natural way to model many granular behaviors.
Its versatility makes it possible to use it to build a
general-purpose granular solver.
Title: |
Multigrid and Spline Collocation Methods on Non-uniform Grids
|
Author: | Guoyu Wang
|
Date: | September 2008
|
Pages: | 98
|
Note: | M.Sc. Thesis, Computer Science Dept., Univ. of Toronto, 2008.
|
File: | Please contact ccc SYMBOL_AT cs DOT toronto DOT edu
|
Keywords: |
quadratic splines, cubic splines, fourth order of convergence,
collocation, multigrid, extension operator, restriction operator,
preconditioning
|
Abstract:
This thesis describes numerical methods for the solution of
one-dimensional second-order differential equations on
non-uniform grids. The methods combine the multigrid and spline
collocation methodologies. Multigrid methods are presented.
Appropriate restriction and extension operators are developed for
quadratic and cubic splines on non-uniform grids. Multigrid
methods are applied to quadratic and cubic spline collocation equations
arising from the discretization of one-dimensional second-order
differential equations. The convergence analysis of a two-grid
method integrated with a damped Richardson relaxation scheme as
smoother shows that the rate of convergence is faster than
$\frac{1}{2}$ for a model second-order equation with
homogeneous Dirichlet boundary conditions. An adaptive
technique for estimating the minimal acceptable coarse grid
size is proposed. Numerical experiments using full multigrid
methods are performed on various one-dimensional second-order
differential equations discretized on non-uniform grids, and
the conditions under which they are convergent and efficient
are studied. The numerical results confirm our analysis.
Title: | Bi-quartic Spline Collocation Methods for Fourth-order
Boundary Value Problems with an Application to
the Biharmonic Dirichlet Problem
|
Author: | Jingrui Zhang
|
Date: | January 2008
|
Pages: | 160
|
Note: | Ph.D. Thesis, Computer Science Dept., Univ. of Toronto, 2008.
|
File: | ccc/jingrui-08-phd.ps.gz
|
Keywords: |
sixth order of convergence, quartic splines, collocation,
fast Fourier Transform (FFT), GMRES, preconditioning
|
Abstract:
We present optimal bi-quartic spline collocation methods for
the numerical solution of fourth-order boundary value problems
on rectangular domains, and apply a particular instance
of these methods to the biharmonic Dirichlet problem.
The methods are based on quartic splines and the collocation
discretization methodology with midpoints of a uniform partition
as the collocation points. While the standard quartic spline
method usually provides second-order approximations, the two
bi-quartic spline collocation methods developed in this thesis are
observed to produce approximations which are of sixth order
at grid points and midpoints, and of fifth order
at other points. Both are based on high order perturbations of the
differential and boundary conditions operators. The one-step
(extrapolated) method forces the approximations to satisfy a
perturbed problem, while the three-step (deferred-correction)
method proceeds in three steps and perturbs the right sides of the
linear equations appropriately.
The three-step bi-quartic spline collocation method is then
applied to the biharmonic Dirichlet problem and a fast solver for
the resulting linear systems is developed. The fast solver
consists of solving an auxiliary biharmonic problem with Dirichlet
and second derivative boundary conditions along the two opposite
boundaries, and a one-dimensional problem related to the two
opposite boundaries. The linear system arising from the bi-quartic
spline collocation discretization of the auxiliary problem is
solvable by fast Fourier transforms. The problem related to the
two opposite boundaries is solved by preconditioned GMRES
(PGMRES). We present a detailed analysis of the performance of
PGMRES by studying bounds for the eigenvalues of the
preconditioned matrix. We show that the number of PGMRES
iterations is independent of the grid size $n$. As a result, the
complexity of the fast solver is $O(n^2\log(n))$. Numerical
experiments from a variety of problems, including practical
applications and problems more general than the analysis assumes,
verify the accuracy of the discretization scheme and the
effectiveness of the fast solver.
Title: | Adaptive Finite Difference Methods for Valuing American Options
|
Author: | Duy Minh Dang
|
Date: | September 2007
|
Pages: | 128
|
Note: | M.Sc. Thesis, Computer Science Dept., Univ. of Toronto, 2007.
|
File: | dmdang-07-msc.pdf
|
Keywords: |
adaptive mesh selection, monitor function, finite differences, Crank-Nicolson,
European option, American option, penalty method, projected SOR
|
Abstract:
We develop space-time adaptive methods for valuing
American options with strong emphasis on American put options.
We examine the application of adaptive techniques
to the Black-Scholes partial differential equation
problem associated with an American put option
in the context of non-uniform second-order finite differences.
At certain timesteps, we obtain a redistribution of
the spatial points based on a monitor function that
attempts to equidistribute the error.
The proposed finite difference discretization on non-uniform grids
and redistribution of the spatial points lead to
linear complementarity problems with $\Mcal$-matrices.
The Projected Successive Over-relaxation
and a penalty method are considered to handle the free boundaries.
We study and compare the accuracy and efficiency
of the considered methods.
A complete proof of convergence and uniqueness
of the projected SOR method under certain conditions is
also presented.
Title: | Quartic Spline Collocation Methods For Second-Order
Two-Point Boundary Value ODE Problems
|
Author: | Guohong Liu
|
Date: | August 2007
|
Pages: | 130
|
Note: | M.Sc. Thesis, Computer Science Dept., Univ. of Toronto, 2007.
|
File: | ghliu-07-msc.ps.gz
|
Keywords: |
sixth order of convergence, quartic splines, collocation,
boundary conditions
|
Abstract:
Collocation methods based on quartic splines are presented for
second-order two-point boundary value problems. In addition to the
boundary conditions specified by the problem, extra boundary
conditions are introduced in order to uniquely determine the
quartic spline approximation. The standard quartic spline
collocation method gives fourth order bounds. Two optimal methods,
namely the extrapolated (one-step) and the deferred-correction
(two-step) methods, are formulated based on appropriate extra
boundary conditions and an appropriate perturbation of the
operators of the differential equation, boundary conditions and
extra boundary conditions. The convergence analysis and error
bounds are developed using a Green's functions approach. The
analysis shows that the maximum discrete error is of sixth order,
and the maximum global error is of fifth order for the optimal
methods.
The properties of the matrices arising from the optimal methods
for a certain class of problems are studied. Non-optimal
collocation methods based on different extra boundary conditions
are also investigated. Problems with layers are also considered,
and a grid refinement technique is presented. The theoretical
behavior of the collocation methods is verified by numerical
results.
Title: | Pricing Convertible Bonds with Dividend Protection
subject to Credit Risk Using a Numerical PDE Approach
|
Author: | Qingkai Mo
|
Date: | April 2006
|
Pages: | 92
|
Note: | M.Sc. Thesis, Computer Science Dept., Univ. of Toronto, 2006.
|
File: | ccc/moqk-06-msc.pdf
|
Keywords: |
convertible bond, dividend protection, credit risk,
Conversion Ratio Adjustment, Dividend Pass-Thru, finite differences,
Crank-Nicolson, Projected Successive Over-Relaxation (PSOR), penalty method.
|
Abstract:
We develop a pricing model for convertible bonds with dividend protection
subject to credit risk by extending the models developed by
Tsiveriotis and Fernandes (TF), and by Ayache, Forsyth and Vetzal (AFV).
We consider two techniques to incorporate the dividend protection feature:
Conversion Ratio Adjustment and Dividend Pass-Thru.
We apply finite difference methods to discretize the PDEs associated
with our models, and study the Projected Successive Over-Relaxation (PSOR)
and penalty methods to handle the free boundaries.
We compare these two methods in terms of convergence rate,
number of iterations per time step and computation time
for pricing convertible bonds without dividends.
Finally, we apply the penalty method, the better of the two methods,
to solve the systems arising from our models for convertible bonds
with dividend protection.
We examine the convergence rates and discuss the difference
between the results from the extended TF and AFV models,
with both dividend protection techniques.
Title: | An Efficient Algorithm Based on Quadratic Spline Collocation and
Finite Difference Methods for Parabolic Partial Differential Equations
|
Author: | Tong Chen
|
Date: | September 2005
|
Pages: | 92
|
Note: | M.Sc. Thesis, Computer Science Dept., Univ. of Toronto, 2005.
|
File: | ccc/tchen-05-msc.pdf
|
Keywords: |
semi-explicit semi-implicit method, optimal order of convergence,
Crank-Nicolson, relaxed conditions for stability
|
Abstract:
An efficient algorithm which combines
quadratic spline collocation methods (QSC) for the space discretization and
classical finite difference methods (FDMs), such as Crank-Nicolson,
for the time discretization
to solve general linear parabolic partial differential equations
has been studied.
By combining QSC and finite differences, a form of the
approximate solution of the problem at each time step can be obtained; thus
the value of the approximate solution and its derivatives
can be easily evaluated at any point of the space domain for each time step.
There are two typical ways for solving this problem:
(a) using QSC in its standard formulation, which has low accuracy
$\mathcal{O}(h^2)$ and low computational work.
More precisely, it requires the solution of a tridiagonal
linear system at each time step;
(b) using optimal QSC, which has high accuracy
$\mathcal{O}(h^4)$ and requires the solution of
either two tridiagonal linear systems or
an almost pentadiagonal linear system at each time step.
A new technique is introduced here which has
the advantages of the above two techniques;
more precisely,
it has high accuracy $\mathcal{O}(h^4)$
and almost the same low computational work as the standard QSC.
Title: | Pricing Convertible Bonds using Partial Differential Equations
|
Author: | Lucy Xingwen Li
|
Date: | September 2005
|
Pages: | 81
|
Note: | M.Sc. Thesis, Computer Science Dept., Univ. of Toronto, 2005.
|
File: | lucy-05-msc.pdf
|
Keywords: |
Coupon payments, call, put, convert, discontinuities,
implicit methods, successive over-relaxation, penalty method
|
Abstract:
A Convertible Bond (CB) is a corporate debt security that gives the holder
the right to exchange future coupon payments and principal repayment for a
prescribed number of shares of equity. Thus, it has both an equity part
and a fixed-income part, and may contain some additional features,
such as callability and puttability.
In this paper, we study the model for valuing Convertible Bonds
with credit risk originally developed by Kostas Tsiveriotis and
Chris Fernandes (TF). The Convertible Bond is a derivative of the
stock price, and the pricing model developed by TF is based on a
free boundary value problem associated with a pair of parabolic
Partial Differential Equations (PDEs) with discontinuities at
the time points when there is a coupon payment, or when the bond
is converted, or when it is called back (purchased) by the issuer,
or when it is put (sold) to the issuer.
We explore the possible derivation of the TF model and study the
convergence of several numerical methods for solving the free
boundary value problem associated with the TF model.
In particular, we consider the Successive Over-Relaxation (SOR) iteration
and a penalty method for solving the linear complementarity problem
used to handle the free boundary.
Special emphasis is given to the effectiveness of the numerical scheme.
Title: | Spline Collocation on Adaptive Grids and Non-Rectangular Domains
|
Author: | Kit Sun Ng
|
Date: | June 2005
|
Pages: | 187
|
Note: | Ph.D. Thesis, Computer Science Dept., Univ. of Toronto, 2005.
|
File: | ngkit-05-phd.ps.gz
|
Keywords: |
quadratic splines, cubic splines, optimal convergence, non-uniform grid,
adaptive methods, gridsize estimator, error estimator,
L-shaped domains, T-shaped domains, skipped grid, moving mesh PDE.
|
Abstract:
Optimal spline collocation is a computationally cost effective
finite element method that has been relatively recently developed.
Although there are many important results
for optimal spline collocation, we believe the method can further
be extended to solve a series of more relevant and difficult problems.
We first consider using optimal spline collocation on non-uniform partitions
to solve one-dimensional linear second-order Boundary Value Problems (BVPs)
that have solutions with regions of high variations, such as boundary layers.
Optimal quadratic spline collocation (QSC) is extended to
non-uniform partitions by using a mapping function from uniform
to non-uniform partitions.
Optimal cubic spline collocation (CSC) on non-uniform partitions
is also examined. The formulation and analysis of the CSC methods
are based, as in the QSC case, on a mapping function from uniform
to non-uniform partitions. However, this method is converted into
a mapping-free method and is integrated with adaptive techniques
that include error and grid size estimators.
Spline collocation is then considered for solving two-dimensional
linear second-order BVPs. We extend optimal spline collocation
to solving BVPs on L- and T-shaped domains, which occur often in practice.
The CSC method can handle general boundary conditions, which
had previously never been done, but our approach requires a
third step to achieve optimal convergence.
We also consider two adaptive approaches for solving problems on rectangular
domains whose solutions have rapid variations. As in the one-dimensional case,
such problems present a significant challenge to obtain accurate approximations
at minimal computational cost.
We first implement a static method, which involves adding collocation
points to regions where the solution has sharp variations, with CSC that uses
skipped grids, which are grids that have several levels of refinement,
in certain regions of the domain.
An optimal skipped grid method is developed by converting the
second-order BVP into a system of three first-order BVPs.
We then implement a dynamic method, which involves moving the collocation
points to regions where the solution has sharp variations, with CSC
that incorporates the {\it Moving Mesh Partial Differential Equation} (MMPDE).
The MMPDE method is further integrated into an adaptive method
that can solve BVPs within a specific error tolerance.
Title: | Multigrid and Cubic Spline Collocation Methods for Advection Equations
|
Author: | Zheng Zeng
|
Date: | April 2005
|
Pages: | 111
|
Note: | M.Sc. Thesis, Computer Science Dept., Univ. of Toronto, 2005.
|
File: |
ccc/zengz-05-msc.ps.gz
|
Keywords: |
Shallow water equations, semi-Lagrangian methods, algebraic elimination,
analytic elimination, systems of PDEs
|
Abstract:
This thesis describes numerical methods for the solution of
advection equations in one-dimensional space. The methods are
based on combining the multigrid and cubic spline collocation (CSC)
methodologies.
Multigrid methods for cubic splines are presented.
Appropriate restriction and extension operators are developed for
cubic splines that satisfy various boundary conditions (BCs),
i.e., Dirichlet, Neumann, periodic and general BCs.
The multigrid methods are applied to CSC equations
arising from the discretization of one-dimensional
second-order differential equations.
The convergence analysis of a two-grid method integrated
with a damped Richardson relaxation scheme as smoother shows the
rate of convergence is equal to or faster than $1/2$
for a model second order equation with periodic BCs.
Multigrid methods for cubic splines are then extended to the solution of
one-dimensional shallow water equations (SWEs).
The SWEs are discretized in time with a semi-Lagrangian
semi-implicit scheme and in space with CSC methods.
At each time step, there are three different space discretization approaches:
Analytic Elimination - Discretization (AED),
Discretization -Algebraic Elimination (DAE),
and Discretization - No Elimination (DNE).
We discuss advantages and disadvantages of each, and
develop new numerical methods based on multigrid and CSC
for solving the SWEs discretized with
either the AED or the DNE approaches.
We compare our methods with Jacobi's iterative method applied
to the same problems.
Through comparison and convergence discussion, we show that our
multigrid methods for CSC are convergent and efficient.
The numerical results confirm our analysis.
Title: | Towards a high-performance method for the biharmonic Dirichlet problem
|
Authors: | Jingrui Zhang
|
Date: | April 2004
|
Pages: | 77
|
Note: | Research paper (equivalent to M.Sc. Thesis),
Computer Science Dept., Univ. of Toronto, 2004.
|
File: | ccc/jingrui-04-msc.ps.gz
|
Keywords: |
spline collocation; fourth-order two-point boundary value problem;
optimal order of convergence; FFT; GMRES; preconditioning;
eigenvalues; eigenvectors
|
Abstract:
The biharmonic Dirichlet problem appears in many applications.
Therefore, there is a lot of interest in developing
high-performance methods for its solution. Two important
components of performance are accuracy and efficiency. In this
research, we combine a sixth order discretization method and a
Fast Fourier Transform (FFT)-based solver for the solution of the
biharmonic Dirichlet problem.
The discretization of the problem uses the collocation methodology
and quartic splines. Quartic Spline Collocation (Quartic SC)
methods that produce sixth order approximations have been recently
developed for the solution of fourth order Boundary Value Problems
(BVPs). In this paper, we apply an adjustment to the quartic
spline basis functions so that the Quartic SC matrix has more
favorable properties. Particular attention is given to the
one-dimensional harmonic problem. We study the properties of two
Quartic SC matrices, more specifically, one that corresponds to
Dirichlet and Neumann conditions (the standard harmonic problem)
and another that corresponds to Dirichlet and second derivative
conditions. The latter is solvable by FFTs and used as
preconditioner for the first. We develop formulae for the
eigenvalues and eigenvectors of the preconditioned matrix. We show
that three iterations suffice to obtain convergence for the GMRES
preconditioned method. Numerical experiments verify the
effectiveness of the solver and preconditioner. We also discuss
the application of the preconditioned GMRES method to the
two-dimensional biharmonic Dirichlet problem.
Title: | Quartic-Spline Collocation Methods for Fourth-Order Two-Point
Boundary Value Problems
|
Author: | Ying Zhu
|
Date: | April 2001
|
Pages: | 73
|
Note: | M.Sc. Thesis, Computer Science Dept., Univ. of Toronto, 2001.
|
File: |
ccc/yz-01-msc.ps.gz
|
Abstract:
This thesis presents numerical methods for the solution of
general linear fourth-order boundary value problems in one dimension.
The methods are based on quartic splines, that is, piecewise quartic
polynomials with C^3 continuity, and the collocation discretization
methodology with the midpoints of a uniform partition being the
collocation points.
The standard quartic-spline collocation method is shown
to be second order. Two sixth-order quartic-spline collocation methods
are developed and analyzed. They are both based on a high order
perturbation of the differential equation and boundary conditions
operators. The one-step method forces the approximation to satisfy
a perturbed problem, while the three-step method proceeds in three steps
and perturbs the right sides of the equations appropriately.
The error analysis follows the Green's function approach and shows
that both methods exhibit optimal order of convergence, that is,
they are locally sixth order on the gridpoints and midpoints
of the uniform partition, and fifth order globally.
The properties of the matrices arising from a restricted class
of problems are studied. Analytic formulae for the eigenvalues
and eigenvectors are developed, and related to those arising
from quadratic-spline collocation matrices.
Numerical results verify the orders of convergence predicted by the
analysis.
Title: | High-Order Spatial Discretization Methods for the Shallow Water
Equations
|
Author: | Anita W. Tam
|
Date: | February 2001
|
Pages: | 160
|
Note: | Ph.D. Thesis, Computer Science Dept., Univ. of Toronto, 2001.
|
File: |
tam-01-phd.ps.gz
|
Abstract:
We present new numerical methods for the shallow water equations on a
sphere in spherical coordinates. In our implementation, the equations
are discretized in time with the two-level semi-Lagrangian
semi-implicit (SLSI) method, and in space on a staggered grid with the
quadratic spline Galerkin (QSG) and the optimal quadratic spline
collocation (OQSC) methods. When discretized on a uniform spatial
grid, the solutions are shown through numerical experiments to be
fourth-order in space locally at the nodes and midpoints of the spatial
grids, and third-order globally.
We also show that, when applied to a simplified version of the shallow
water equations, each of our algorithms yields a neutrally stable
solution for the meteorologically significant Rossby waves. Moreover,
we demonstrate that the Helmholtz equation associated with the shallow
water equations should be derived algebraically rather than
analytically in order for the algorithms to be stable with respect to
the Rossby waves. These results are verified numerically using Boyd's
equatorial wave equations with initial conditions to generate a
soliton.
We then analyze the performance of our methods on various staggered
grids --- the A-, B-, and C-grids. From an eigenvalue analysis of our
simplified version of the shallow water equations, we conclude that,
when discretized on the C-grid, our algorithms faithfully capture the
group velocity of inertia-gravity waves. Our analysis suggests that
neither the A- nor B-grids will produce such good results. Our
theoretical results are supported by numerical experiments, in which we
discretize Boyd's equatorial wave equations using different staggered
grids and set the initial conditions for the problem to generate
gravitation modes instead of a soliton. With both the A- and B-grids,
some waves are observed to travel in the wrong direction, whereas, with
the C-grid, gravity waves of all wavelengths propagate in the correct
direction.
Title: | Quadratic Spline Collocation Methods for Systems of Elliptic PDEs
|
Authors: | Kit Sun Ng
|
Date: | April 2000
|
Pages: | 67
|
Note: | M.Sc. Thesis, Computer Science Departmnent,
University of Toronto, 2000.
|
File: | ccc/ngkit-00-msc.ps.gz
|
Abstract:
We consider Quadratic Spline Collocation (QSC) methods
for solving systems of two linear second-order PDEs in two dimensions.
Optimal order approximation to the solution is obtained,
in the sense that the convergence order
of the QSC approximation is the same as the
order of the quadratic spline interpolant.
We study the matrix properties of the linear system arising
from the discretization of systems of two PDEs by QSC.
We give sufficient conditions under which the QSC linear system
is uniquely solvable
and the optimal order of convergence for the QSC approximation
is guaranteed.
We develop fast direct solvers based on Fast Fourier Transforms (FFTs)
and iterative methods using multigrid or FFT preconditioners
for solving the above linear system.
Numerical results demonstrate that the QSC methods are fourth
order locally on certain points and third order globally,
and that the computational complexity of the linear solvers developed
is almost asymptotically optimal.
The QSC methods are compared to conventional second order discretization
methods and are shown to produce smaller approximation errors
in the same computation time, while they achieve the same
accuracy in less time.
Title: | Numerical Solution of the Shallow-Water Equations
on Distributed Memory Systems
|
Author: | Xiaoliang (Lloyd) Ding
|
Date: | October 1998
|
Pages: | 73
|
Note: | M.Sc. Thesis, Computer Science Departmnent,
University of Toronto, 1998.
|
File: | ding.98.msc.ps.gz
|
Abstract:
The shallow-water equations are often used as a mathematical model when
numerical methods for solving weather or climate prediction problems
are tested. This thesis studies the performance and scalability of
numerical methods for the shallow-water equations on distributed memory
systems. Time integration of the numerical methods is based on a
two-time-level semi-implicit semi-Lagrangian scheme. To solve the
Helmholtz problem that arises at each time-step, a fast direct solver
based on FFTs is used. This technique requires a data transposition,
which is a time consuming operation on a distributed memory system.
The data transposition requires an all-to-all type communication. The
Direct Exchange algorithm, commonly used for this task, is inefficient
if the number of processors is large. We investigate a series of more
sophisticated techniques: the Ring Exchange, Mesh Exchange and Cube
Exchange algorithms. We compare the performance of these techniques
with that of the Message Passing Interface collective communication
routine. We perform an isoefficiency analysis to measure the
scalability of the parallel system. The numerical experiments are
carried out on a Cray T3E with 512 processors and on an
ethernet-connected cluster of 36 workstations.
Both the analysis and our numerical results indicate that the Cube
Exchange and Mesh Exchange algorithms perform better than the other
methods. However, even with the more sophisticated techniques, the
parallel algorithm for the shallow-water equations remains essentially
unscalable.
Title: | Fast Fourier transform solvers for quadratic spline collocation
|
Author: | A. Constas
|
Date: | July 1996
|
Pages: | 64
|
Note: | M.Sc. Thesis, Computer Science Dept., Univ. of Toronto, 1996.
|
File: | ccc/constas-96-msc.ps.gz
|
Abstract:
The solution of linear systems arising from the discretisation
of second order elliptic Boundary Value Problems (BVPs)
defined in rectangular domains
by solvers employing Fast Fourier Transforms (FFTs) is studied.
The types of boundary conditions considered are
homogeneous Dirichlet, periodic, homogeneous Neumann or
combinations of these.
Two-dimensional and three-dimensional problems are considered.
Quadratic spline collocation methods are used
for the discretisation of the BVPs.
The FFT solvers are first developed for differential operators
with constant coefficients and even order derivatives.
Diagonal scalings of these solvers are then used
as preconditioners for the solution of linear systems
arising from the discretisation of general second order
linear elliptic BVPs.
Several acceleration methods are tested.
Results from numerical experiments
that demonstrate the asymptotically optimal convergence
and computational performance of the FFT solvers are presented.