Theme
Organizers: Ray Spiteri (Saskatchewan) and Buks van Rensburg (York)
For
the first time at the CAIMS annual meeting, the two distinct communities
of computational PDE and discrete mathematics will be connected
in this theme. Discretizations of PDE lead to (large) discrete dynamical
systems, and a popular strategy of studying discrete dynamical systems
is to go in the opposite direction (ie, limits of large numbers,
continuum limits). As a consequence, we hope that linking the two
communities will lead to a crossfertilization of ideas. The session
will focus on scientific computing and numerical modeling of complex
systems such as percolation, polymers, Monte Carlo simulations of
protein folding and numerical simulation of biological and biomedical
systems including knotting and linking in biopolymers.
Speakers:
Dhavide Aruliah (University of Ontario
Institute of Technology
Robert Bridson (University of British
Columbia)
Matthew Emmett (University of North
Carolina)
Peter Grassberger (University of
Calgary)
Chen Greif (University of British
Columbia)
Ray Kapral (University of Toronto)
Vince Lyzinski (John Hopkins University)
Felicia Magpantay (McGill University)
Colin Macdonald (University of
Oxford)
Peter Minev (University of Alberta)
Enzo Orlandini (University of Padova)
George Patrick (University of Saskatchewan)
Thomas Prellberg (Queen Mary University,
London)
Oluwaseun Sharomi (Saskatchewan)
Erkan Tuzel (Worcester
Polytechnic Institute)
Buks van Rensburg (York)
Ivona
Bezakova (Rochester Institute of Technology)
(talk cancelled)
Daniel Stefankovic (University of Rochester) (talk cancelled)
Plenary
Speaker
Ted Kim (Santa
Barbara)
Subspace Simulation of Nonlinear Materials
Everyday
materials such as biological tissue exhibit geometric and constitutive
nonlinearities which are crucial in applications such as surgical
simulation and realistic human animation. However, these same
nonlinearities also make efficient finite element simulation
difficult, as computing nonlinear forces and their Jacobians
can be computationally intensive. Reduced order methods, which
limit simulations to a lowdimensional subspace, have the potential
to accelerate these simulations by orders of magnitude. In this
talk, I will describe a subspace method we have developed that
efficiently computes all the necessary nonlinear quantities by
evaluating them at a discrete set of "cubature" points,
an online integration method that can accelerate simulations even
when the subspace is not known a priori, and a domain decomposition
method that efficiently adds deformations to discretely partitioned,
articulated simulations. Using our methods, we have observed speedups
of one to four orders of magnitude.
Bio:
Theodore Kim joined the faculty of Media Arts and Technology Program
at the University of California, Santa Barbara, in the Fall of
2011 as an Assistant Professor. He conducts research in physicallybased
simulation for computer graphics, particularly in fluid dynamics,
solid mechanics, and multicore techniques. Techniques from his
research have been used in over a dozen movies. Previously, he
has been an Assistant Professor in Computer Science at the University
of Saskatchewan, a PostDoctoral Associate at Cornell University,
and a PostDoctoral Fellow at IBM TJ Watson Research Center. He
received his PhD in Computer Science from the University of North
Carolina at Chapel Hill in 2006.
Speakers
Dhavide
Aruliah (University of Ontario Institute of Technology)
Toward Accurate Methods for Forensic Bloodstain Pattern Analysis
At
present, bloodstain pattern analysis (BPA) is largely an empirical,
qualitative subdiscipline of forensic science. Forensic BPA specialists
analyze bloodstains found at crime scenes and attempt to infer
the location of the bloodletting impact (or impacts) that caused
to the bloodstains. Traditional BPA methods for reconstructing
crime scenes (notably stringing) ignore the effects of gravity
and aerodynamic drag on trajectories of blood droplets in flight.
Such assumptions may produce acceptable qualitative predictions
under certain conditions (e.g., when analyzing bloodstains caused
by droplets moving at high speeds due to discharge of a firearm).
However, in many circumstances, backtracking blood droplets along
straight lines from bloodstains to a static impact location are
misleading (e.g., when bloodstains are caused by assault with
a blunt instrument or an edged weapon). Our ultimate aim is to
develop software tools for quantitative analysis to support forensic
BPA analysts.
We
describe our framework for recording simulated bloodletting events
and extracting droplet flight trajectories. The simulations consist
of fake blood encased in ballistic gel being splattered by projectiles.
The resulting blood flight trajectories are recorded by a stereo
pair of highspeed cameras and the bloodstains are photographed
to provide a collection of video and static image data sets to
validate our inversion framework. We employ a sophisticated algorithm
incorporating background removal, segmentation, 2D motion tracking
and 3D reconstruction to extract meaningful flight trajectories
from the video data.
Robert
Bridson (University of British Columbia)
Generic implicit constrained dynamics and algebraic preconditioners
for graphics
The
emerging approach of "virtual practical effects" in
feature film production leverages artists' intuition for the real
world by simulating physics in a virtual world  instead of struggling
to control complex software models, artists set up effects as
they might in reality and let physics do the rest. This demands
robust numerical solvers which can couple diverse physics together
in unanticipated ways. We present a framework which unifies many
previous elements (from viscous incompressible fluids to rigid
body contact to inextensible cloth) and reduces to a sequence
of leastsquareslike problems. We further explore a new approach
to algebraic almostblackbox domain decomposition which shows
promise for tackling linear systems of this category.
Matthew
Emmett (University of North Carolina)
The Parallel Full Approximation Scheme in Space and Time (PFASST)
algorithm
The
Parallel Full Approximation Scheme in Space and Time (PFASST)
algorithm for parallelizing PDEs in time is presented. To efficiently
parallelize PDEs in time, the PFASST algorithm decomposes the
time domain into several time slices. After a provisional solution
is obtained using a relatively inexpensive time integration scheme,
the solution is iteratively improved using a deferred correction
scheme. To further improve parallel efficiency, the PFASST algorithm
uses a hierarchy of discretizations at different spatial and temporal
resolutions and employs an analog of the multigrid full approximation
scheme to transfer information between the discretizations.
Peter
Grassberger (University of Calgary)
"Who said that he understands percolation?"
Although
percolation theory was considered a mature subject several years
ago, recent progress has changed this radically. While "classical"
or "ordinary" percolation (OP) is a second order phase
transition between long range connectivity and disconnectedness
on diluted regular lattices or random graphs, examples have now
been found where this transition can range from infinite order
to first order. The latter is of particular importance in social
sciences, where first order percolation transitions show up as
a consequence of synergistic effects, and I will point out analogies
with the relationship between percolation and rough interfaces
in physics. Another case where first order percolation transitions
show up is interdependent networks, although first claims about
this have to be substantially modified  in some cases of interdependent
networks the transition is second order but in a new universality
class. A similar but even more unexpected result holds for socalled
"Achleoptas processes" that were originally claimed
to show first order transitions, but which actually show second
order transitions with a completely new type of finite size scaling.
Finally, I will present "agglomerative percolation"
(AP), a model originally introduced to understand the claim that
network renormalization can demonstrate the fractality of some
small world networks. Due to a spontaneously broken symmetry on
bipartite graphs, AP leads e.g. to different scaling behaviors
on square and triangular 2d lattices, in flagrant violations
of universality.
Chen
Greif (University of British Columbia)
Numerical Solution of SaddlePoint Linear Systems
Constrained
partial differential equations and optimization problems typically
require the need to solve special linear systems known as saddlepoint
systems. When the matrices are very large and sparse, iterative
methods must be used. A challenge here is to derive and apply
solution methods that exploit the properties and the structure
of the given discrete operators, and yield fast convergence while
imposing reasonable computer storage requirements. In this talk
I will provide an overview of solution techniques. In particular,
we will discuss effective preconditioners and their spectral properties,
bounds on convergence rates, and computational challenges.
Ray
Kapral (Toronto)
Mesoscopic Dynamics of Biopolymers and Protein Machines
The
dynamics of polymers and biopolymers in solution and in crowded
molecular environments which mimic some features of the interior
of a biochemical cell will be discussed. In particular, the dynamics
of protein machines that utilize chemical energy to effect cyclic
conformational changes will be described. The investigation of
the dynamics of such complex systems requires knowledge of the
time evolution on physically relevant long distance and time scales.
This often necessitates a coarse grained or mesoscopic treatment
of the dynamics. A particlebased mesoscopic dynamical method,
multiparticle collision dynamics, which conserves mass, momentum
and energy, has been constructed and utilized to study the dynamics
of these systems. The dynamics can be described by a Markov chain
model in the full phase space of the system and, using projection
operator or other methods, can be shown to lead to the continuum
NavierStokes or reactiondiffusion descriptions on long distance
and time scales.
Vince
Lyzinski (Johns Hopkins)
Strong Stationary Duality for Diffusions
Strong
stationary duality has had a wideranging impact on Markov chain
theory since its conception by Diaconis and Fill in 1990. Its
diverse applications range from perfect sampling extensions of
Markov Chain Monte Carlo to the establishment of cutoff phenomena
for wide classes of Markov chains. We extend the idea of strong
stationary duality to onedimensional diffusion processes and
in doing so recover some classical Markov chain results in the
diffusion setting. This is joint work with my PhD advisor James
Allen Fill.
Colin
Macdonald (University of Oxford)
Mathematics and Algorithms for Simple Computing on Surfaces
The
Closest Point Method is a simple numerical technique for solving
partial differential equations (PDEs) on curved surfaces or other
general domains. The method works by embedding the surface in
a higherdimensional space and solving the PDE in that space.
Numerically, we can use simple finite difference and interpolation
schemes on uniform Cartesian grids. This presentation provides
a brief overview of the Closest Point Method and outlines some
current results. In the spirit of a minisymposium that is examining
links between continuous and discrete computational mathematics,
I will discuss the issue of finding a narrow band surrounding
a complicated surface (this is useful for an efficient implementation)
and how we approach this (discrete) problem.
Joint
work with Tom Maerz, Ingrid von Glehn, and Yujia Chen (Oxford)
and Steve Ruuth (SFU).
Felicia
Magpantay (York
University)
Stability of backward Euler applied to a model state dependent
delay differential equation
We
consider the stability properties of the backward Euler method
for delay differential equations (DDEs) with respect to a test
equation. We consider two cases: (i) constant delay (ii) state
dependent delay. RungeKutta methods applied to DDEs have been
studied by many authors who have mainly considered stability regions
independent of the delay and/or require the stepsize to be a
submultiple of the delay (for the constant delay case). These
assumptions put too much restriction on the method and cannot
be used in the state dependent case. We direct attention to the
dependence of the stability regions to the delay function and
present results that use general step sizes. The techniques used
are derived from the method of Lyapunov functionals to directly
prove the stability of zero solution. We also prove the local
stability of backward Euler when the stepsize used in an important
case.
Joint
work with A.R. Humphries and N. Guglielmi.
Peter
Minev (University
of Alberta)
A Direction Splitting Algorithm for Flow Problems in Complex/Moving
Geometries
An
extension of the direction splitting method for the incompressible
NavierStokes equations proposed in [1], to flow problems in complex,
possibly time dependent geometries will be presented. The idea
stems from the idea of the fictitious domain/penalty methods for
flows in complex geometry. In our case, the velocity boundary
conditions on the domain boundary are approximated with a secondorder
of accuracy while the pressure subproblem is harmonically extended
in a fictitious domain such that the overall domain of the problem
is of a simple rectangular/parallelepiped shape.
The
new technique is still unconditionally stable for the Stokes problem
and retains the same convergence rate in both, time and space,
as the CrankNicolson scheme. A key advantage of this approach
is that the algorithm has a very impressive parallel performance
since it requires the solution of onedimensional problems only,
which can be performed very efficiently in parallel by a domaindecomposition
Schur complement approach. Numerical results illustrating the
convergence of the scheme in space and time will be presented.
Finally, the implementation of the scheme for particulate flows
will be discussed and some validation results for such flows will
be presented.
[1]
J.L. Guermond, P.D. Minev, A new class of massively parallel direction
splitting for the incompressible NavierStokes equations. Computer
Methods in Applied Mechanics and Engineering, 200 (2011), 2083â€“2093.
Enzo
Orlandini (Padova)
Dynamics of knotted polymers
Knots
are frequent in long polymer rings at equilibrium and it is now
well established that their presence can affect to some extent
the static properties of the chain. On the other hand, the presence
of topological constraints (knots) in circular and linear polymers
may influence also their dynamical properties. This has been indeed
shown in recent experiments where the motion of a single knotted
DNA has been followed within a viscous solution and in the presence
of a stretching force. These experiments raise interesting challenges
to the theoretical comprehension of the problem, an issue which
is still in its infancy. As a first step towards the understanding
of the mechanism underlying the mobility of a knot along the ring
backbone and its effect on the overall polymer dynamic we investigate,
by Monte Carlo and molecular dynamics simulations, the dynamics
of knotted rings under good solvent conditions. By using an algorithm
that detects position and size of the knot we are able to monitor
the motion of the topologically entangled sub region both in space
and along the ring backbone.
This
allows identifying in knotted rings a novel, slow topological
timescale, whose origin can be related to a selfreptation motion
of the knotted region.
For
open chains, knotted configurations do not represent an equilibrium
state any more. However, under suitable conditions, (for example
very tight knots or quite rigid chains) knotted metastable states
persist for a very long time and a statistical description of
their dynamical properties is then possible. By performing off
lattice molecular dynamic simulations of a semiflexible polymer
we estimate the average living time and the survival probability
as a function of the initial conditions (size of the initial knot)
and knot type. This analysis has been extended to the case in
which the linear polymer is subject to an external stretching
force.
George
Patrick (University of Saskatchewan)
Automatically
generated variational integrators
Many
fundamental physical systems have variational formulations, such
as mechanical systems in their Lagrangian formulation. Discretization
of the variational principles leads to (implicit) symplectic and
momentum preserving one step integration methods. However, such
methods can be very complicated.
I
will describe some advances in the basic theory of variational
integrators, and a software system called AUDELSI, which converts
any ordinary one step method into a variational integrator of
the same order.
Thomas
Prellberg (Queen Mary University London)
Rare event sampling with stochastic growth algorithms
We
discuss uniform sampling algorithms that are based on stochastic
growth methods, using sampling of extreme configurations of polymers
in simple lattice models as a motivation. We shall show how a
series of clever enhancements to a fiftyodd year old algorithm,
the Rosenbluth method, leads to a suite of cuttingedge algorithms
capable of uniform sampling of equilibrium statistical mechanical
systems of polymers in situations where competing algorithms failed
to perform well. Examples range from collapsed homopolymers near
sticky surfaces to models of protein folding.
Oluwaseun
Sharomi (University of Saskatchewan)
The Significance of Order in the Numerical Simulation of the
Bidomain Equation in Parallel
The
propagation of electrical activity in the heart can be modelled
by the bidomain equations. However to obtain clinically useful
data from the bidomain equations, they must be solved with many
millions of unknowns. Naturally, to obtain such data in real time
is the ultimate goal, but at present we are still an order of
magnitude or two away from being able to do so. The spectral/hp
element method can be considered to be a highorder extension
of the traditional “finite” or “spectral”
element methods, where convergence is not only possible through
reducing the mesh size h but also through increasing the local
polynomial order p of the basis functions used to expand the solution.
We are interested in evaluating the effectiveness of a highorder
method in serial against that of a loworder method in parallel.
We find that highorder methods in serial can outperform loworder
methods in parallel. These findings suggest software developers
for the bidomain equations should not forego the implementation
of highorder methods in their efforts to parallelize their solvers.
Erkan
Tuzel (Physics,
Worcester Polytechnic Institute)
Constrained Polymer Dynamics in a Mesoscale Solvent
While
modeling polymer solutions, the presence of multiple time scales,
such as the intermolecular bond potentials, makes quantitative
analysis of results difficult, and simulations costly. Here, we
show how these degrees of freedom can be replaced by rigid bond
constraints—as commonly done in Brownian dynamics—for
polymers embedded in a recently introduced hydrodynamic solvent
known as Stochastic Rotation Dynamics (SRD) (or MultiParticle
Collision Dynamics  MPCD). We then discuss applications of this
approach to various systems of biological interest.
EJ
Janse van Rensburg (York University)
Statistics of knotted lattice polygons
In
this talk I discuss the implementation of the GAS algorithm using
BFACF elementary moves which we implement to sample knotted polygons
in the cubic lattice. The GAS algorithm is an approximate enumeration
algorithm, and I show how it can be implemented to estimate the
number of distinct polygons of given length and fixed knot types.
The numerical results we obtain make it possible to examine the
scaling of knotted lattice polygons. For example, our data indicate
that unknotted polygons dominate cubic lattice polygon statistics
up to lengths of about 170,000 steps, thereafter polygons of knot
type the trefoil becomes more numerous. In addition, the relative
frequencies of various knot types can be determined  we found
that trefoil knots are about 28 times more likely to occur than
figure eight knots in long polygons.
Contributed
Talks
*Ken Roberts, Western University
Coauthors: S. R. Valluri, Muralikrishna Molli, M. ChandraSekhar,
K. Venkataramaniah, P. C. Deshmukh
A Study of Polylogarithmic Equations for Thermoelectric Nanomaterials
In
the design of thermoelectric nanomaterials, various expressions
arise which involve the polylogarithms. A computational study
of some of those equations has led to curious results, suggesting
additional properties to be explored, either the mathematics of
polylogarithms, or the statistical mechanics underlying the material
models which lead to polylogarithms.
We will present a progress report on our efforts to explore and
understand these relationships. There are possibly some insights
to be gained into the statistical mechanics of thermoelectric
materials, via utilizing polylogarithms of complex order, or via
generalizing polylogarithms to a form related to the Dirichlet
Lseries in analytic number theory.
Thomas
Humphries
Department of Mathematics and Statistics, Memorial University of
Newfoundland
Coauthors: Ronald Haynes, Lesley James
Simultaneous optimization of oil well placement and control using
a hybrid globallocal strategy
Two
important decisions in maximizing production from an oil field
are where to place injection and production wells, and how to
control the flow rates at these wells. In this presentation we
address the reservoir optimization problem using an automatic
approach that hybridizes two blackbox optimization methods: particle
swarm optimization (PSO) and generalized pattern search (GPS).
PSO provides a semirandom, global exploration of search space,
while GPS systematically explores local regions of space. We present
simulation results showing that this hybridized approach outperforms
the independent application of PSO to this problem.
Erin
Moulding, Department of Mathematics, UBC
Coauthors: Chen Greif, UBC, and Dominique Orban, École Polytechnique
Bounds on the Eigenvalues of 3x3 Block Matrices arising from
InteriorPoint Methods
Interiorpoint
methods feature prominently in the solution of constrained optimization
problems, and involve solving linear systems with a sequence of
3x3 block matrices which become increasingly illconditioned.
A review of the literature suggests that most practitioners reduce
these to either 2x2 block saddlepoint matrices or 1x1 block normal
equations. In this talk we explore whether it pays off to perform
such reductions. We use energy estimates to obtain bounds on the
eigenvalues of the unreduced matrix, which indicate that in terms
of spectral structure, it may be better to keep the matrix in
its original form.

Wayne
Enright ,
University of Toronto
Coauthors: Bo
Wang
Parameter Estimation for ODEs using a CrossEntropy Approach
Parameter
Estimation for ODEs and DDEs is an important topic in numerical
analysis. In this paper, we present a novel approach to address
this inverse problem. Crossentropy algorithms are general algorithm
which can be applied to solve global optimization problems. The
main steps of crossentropy methods are first to generate a set
of trial samples from a certain distribution, then to update the
distribution based on these generated sample trials. To overcome
the prohibitive computation of standard crossentropy algorithms,
we develop a modification combining local search techniques. The
modified crossentropy algorithm can speed the convergence rate
and improve the accuracy simultaneously.
Two different coding schemes (continuous coding and discrete coding)
are also introduced. Continuous coding uses a truncated multivariate
Gaussian to generate trial samples, while discrete coding reduces
the search space to a finite (but dense) subset of the feasible
parameter values and uses a Bernoulli distribution to generate
the trial samples (which are fixed point approximation of the
actual parameters) . Extensive numerical and real experiments
are conducted to illustrate the power and advantages of the proposed
methods. Compared to other existing stateoftheart approaches
on some benchmark problems for parameter estimation, our methods
have three main advantages: 1) They are robust to noise in the
data to be fitted; 2) They are not sensitive to the number of
observation points (in contrast to most existing approaches) ;
3) The modified versions exhibit faster convergence than the original
versions, thus they are more efficient, without sacrificing accuracy.

TALK
CANCELLED
Ivona
Bezakova (RIT)
Counting and sampling minimum
cuts in weighted planar graphs
We
will discuss two minimum cut problems in weighted planar graphs:
minimum sourcesink cuts and contiguous minimum singlesourcemultisink
cuts.
A
sourcesink cut is a set of vertices containing the source vertex
and not the sink vertex (or, in the case of multiple sinks, not
containing any of the sink vertices). A cut is minimum if the sum
of the weights of the cut edges, connecting a vertex in the cut
set with a vertex outside the cut set, is the smallest possible.
A cut is contiguous if the cut set can be separated from the remaining
vertices by a simply connected planar region whose boundary intersects
only the cut edges.
We
will present an O(n^2) algorithm counting all minimum sourcesink
cuts in weighted planar graphs, where n is the number of vertices.
We will also sketch an O(n^3) algorithm counting all contiguous
minimum singlesourcemultisink cuts. In both cases, having completed
the counting part, subsequent sampling is very fast: a uniformly
random cut can be produced in additional linear time.
The
counting algorithms share a common outline. First, we reduce the
problem to the problem of counting a different type of cuts in an
unweighted planar directed acyclic graph (these cuts can also be
thought of as maximal antichains in the corresponding partially
ordered set). These cuts correspond to certain cycles in the planar
dual graph and we employ dynamic programming to count them. We will
discuss the first algorithm in detail and briefly sketch the issues
encountered by the contiguity requirement.
Minimum
sourcesink and contiguous minimum singlesourcemultisinks cuts
have applications in computer vision and medical imaging where the
underlying graph often forms a 2D grid (lattice).
Based
on joint works with Adam Friedlander and Zach Langley
Daniel
Stefankovic (Rochester)
Connection
between counting and sampling
Counting
problems arise in a wide variety of areas, for example, computer
science, enumerative combinatorics, and statistical physics (estimating
the value of a partition function is a counting problem). I will
talk about the following aspect of approximate counting: given a
sampling algorithm, how can one efficiently translate it into a
counting algorithm?
back
to top
