Krylov subspace methods for approximating functions of symmetric positive definite matrices with applications to applied statistics and anomalous diffusion

Matrix function approximation is a current focus of worldwide interest and finds application in a variety of areas of applied mathematics and statistics. In this thesis we focus on the approximation of A^(-α/2)b, where A ∈ ℝ^(n×n) is a large, sparse symmetric positive definite matrix and b ∈ ℝ^n is a vector. In particular, we will focus on matrix function techniques for sampling from Gaussian Markov random fields in applied statistics and the solution of fractional-in-space partial differential equations. Gaussian Markov random fields (GMRFs) are multivariate normal random variables characterised by a sparse precision (inverse covariance) matrix. GMRFs are popular models in computational spatial statistics as the sparse structure can be exploited, typically through the use of the sparse Cholesky decomposition, to construct fast sampling methods. It is well known, however, that for sufficiently large problems, iterative methods for solving linear systems outperform direct methods. Fractional-in-space partial differential equations arise in models of processes undergoing anomalous diffusion. Unfortunately, as the fractional Laplacian is a non-local operator, numerical methods based on the direct discretisation of these equations typically requires the solution of dense linear systems, which is impractical for fine discretisations. In this thesis, novel applications of Krylov subspace approximations to matrix functions for both of these problems are investigated. Matrix functions arise when sampling from a GMRF by noting that the Cholesky decomposition A = LL^T is, essentially, a `square root' of the precision matrix A. Therefore, we can replace the usual sampling method, which forms x = L^(-T)z, with x = A^(-1/2)z, where z is a vector of independent and identically distributed standard normal random variables. Similarly, the matrix transfer technique can be used to build solutions to the fractional Poisson equation of the form ϕn = A^(-α/2)b, where A is the finite difference approximation to the Laplacian. Hence both applications require the approximation of f(A)b, where f(t) = t^(-α/2) and A is sparse. In this thesis we will compare the Lanczos approximation, the shift-and-invert Lanczos approximation, the extended Krylov subspace method, rational approximations and the restarted Lanczos approximation for approximating matrix functions of this form. A number of new and novel results are presented in this thesis. Firstly, we prove the convergence of the matrix transfer technique for the solution of the fractional Poisson equation and we give conditions by which the finite difference discretisation can be replaced by other methods for discretising the Laplacian. We then investigate a number of methods for approximating matrix functions of the form A^(-α/2)b and investigate stopping criteria for these methods. In particular, we derive a new method for restarting the Lanczos approximation to f(A)b. We then apply these techniques to the problem of sampling from a GMRF and construct a full suite of methods for sampling conditioned on linear constraints and approximating the likelihood. Finally, we consider the problem of sampling from a generalised Matern random field, which combines our techniques for solving fractional-in-space partial differential equations with our method for sampling from GMRFs.

Differential and numerical models of hysteretic systems with stochastic and deterministic inputs

Many deterministic models with hysteresis have been developed in the areas of economics, finance, terrestrial hydrology and biology. These models lack any stochastic element which can often have a strong effect in these areas. In this work stochastically driven closed loop systems with hysteresis type memory are studied. This type of system is presented as a possible stochastic counterpart to deterministic models in the areas of economics, finance, terrestrial hydrology and biology. Some price dynamics models are presented as a motivation for the development of this type of model. Numerical schemes for solving this class of stochastic differential equation are developed in order to examine the prototype models presented. As a means of further testing the developed numerical schemes, numerical examination is made of the behaviour near equilibrium of coupled ordinary differential equations where the time derivative of the Preisach operator is included in one of the equations. A model of two phenotype bacteria is also presented. This model is examined to explore memory effects and related hysteresis effects in the area of biology. The memory effects found in this model are similar to that found in the non-ideal relay. This non-ideal relay type behaviour is used to model a colony of bacteria with multiple switching thresholds. This model contains a Preisach type memory with a variable Preisach weight function. Shown numerically for this multi-threshold model is a pattern formation for the distribution of the phenotypes among the available thresholds.

On the limit cycles of a class of piecewise linear differential systems in R-4 with two zones

Fundação de Amparo à Pesquisa do Estado de São Paulo (FAPESP)

Optimal mean-variance control for discrete-time linear systems with Markovian jumps and multiplicative noises

In this paper, we consider the stochastic optimal control problem of discrete-time linear systems subject to Markov jumps and multiplicative noises under two criteria. The first one is an unconstrained mean-variance trade-off performance criterion along the time, and the second one is a minimum variance criterion along the time with constraints on the expected output. We present explicit conditions for the existence of an optimal control strategy for the problems, generalizing previous results in the literature. We conclude the paper by presenting a numerical example of a multi-period portfolio selection problem with regime switching in which it is desired to minimize the sum of the variances of the portfolio along the time under the restriction of keeping the expected value of the portfolio greater than some minimum values specified by the investor. (C) 2011 Elsevier Ltd. All rights reserved.

Results of analyses of fluid inclusions from ODP Site 735B, Leg 118 and 176

Microthermometric and isotopic analyses of fluid inclusions in primitive olivine gabbros, oxide gabbros, and evolved granitic material recovered from Ocean Drilling Program Hole 735B at the Southwest Indian Ridge provide new insights into the evolution of C-O-H-NaCl fluids in the plutonic foundation of the oceanic crust. The variably altered and deformed plutonic rocks span a crustal section of over 1500 m and record a remarkably complex magma-hydrothermal history. Magmatic fluids within this suite followed two chemically distinct paths during cooling through the subsolidus regime: the first path included formation of CO2+CH4+H2O+C fluids with up to 43 mole% CH4; the second path produced hypersaline brines that contain up to 50% NaCl equivalent salinities. Subsequent to devolatilization, respeciation of magmatic CO2, attendant graphite precipitation, and cooling from 800°C to 500°C promoted formation of CH4-enriched fluids. These fluids are characterized by average d13C(CH4) values of -27.1+/-4.3 per mil (N=45) with associated d13C(CO2) compositions ranging from -24.9 per mil to -1.9 per mil (N=39), and average dD values of exsolved vapor of -41+/-12 per mil (N=23). In pods, veins, and lenses of highly fractionated residual material, hypersaline brines formed during condensation and by direct exsolution in the absence of a conjugate vapor phase. Entrapped CO2+CH4+H2O-rich fluids within many oxide-bearing rocks and felsic zones are significantly depleted in 13C (with d13C(CO2) values down to about -25 per mil) and contain CO2 concentrations higher than those predicted by equilibrium devolatilization models. We hypothesize that lower effective pressures in high-temperature shear zones promoted infiltration of highly fractionated melts and compositionally evolved volatiles into focused zones of deformation, significantly weakening the rock strength. In felsic-rich zones, volatile build-up may have driven hydraulic fracturing of gabbroic wall rocks resulting in the formation of magmatic breccias. Comparison of isotopic compositions of fluids in plutonic rocks from 735B, the MARK area of the Mid-Atlantic Ridge, and the Mid-Cayman Rise indicate (1) that the carbon isotope composition of the lower oceanic crust may be far more heterogeneous than previously believed and (2) that carbon-bearing species in the oceanic crust and their distribution at depth are highly variable.

Chemical composition of glass inclusions of volcanic rocks from the Bonin Arc

The chemical composition of glass inclusions in phenocrystic plagioclase and pyroxene from Sites 792 and 793, drilled during Ocean Drilling Program Leg 126 in the Bonin Arc, is examined. Immiscible liquid, which is preserved as glass inclusions with unmixed textures in plagioclase, is observed in a high-magnesian andesite, which suggests an important role of liquid immiscibility in the fractionation of high-magnesian andesite. In other andesitic rocks (SiO2 = 57-60 wt%), such unmixed textures of glass inclusions in calcic plagioclase with a similar percentage of An (around 80%) is not found. The degree of fractionation and mixing of liquid are inferred from the glass composition in pyroxene.

Chemical composition of glass inclusions from ODP Holes 157-953C and 157-956B

We report S concentrations and relative proportions of (SO4)2- and S2- in OL- and CPX-hosted glass inclusions and in host glassy lapilli from Miocene basaltic hyaloclastites drilled north and south of Gran Canaria during ODP Leg 157. Compositions of glass inclusions and lapilli resemble those of subaerial Miocene shield basalts on Gran Canaria and comprise mafic to more evolved tholeiitic to alkali basalt and basanite (10.3-3.7 wt.% MgO, 44.5-56.9 wt.% SiO2). Glass inclusions fall into three groups based on their S concentrations: a high-sulfur group (1050 to 5810 ppm S), an intermediate-sulfur group (510 to 1740 ppm S), and a low-sulfur group (<500 ppm S). The most S-rich inclusions have the highest and nearly constant proportion of sulfur dissolved as sulfate determined by electron microprobe measurements of SKa peak shift. Their average S6+/S_total value is 0.75+/-0.09, unusually high for ocean island basalt magmas. The low-sulfur group inclusions have low S6+/S_total ratios (0.08+/-0.05), whereas intermediate sulfur group inclusions show a wide range of S6+/S_total (0.05-0.83). Glassy lapilli and their crystal-hosted glass inclusions with S concentrations of 50 to 1140 ppm S have very similar S6+/S_total ratios of 0.36+/-0.06 implying that sulfur degassing does not affect the proportion of (SO4)2- and S2- in the magma. The oxygen fugacities estimated from S6+/S_total ratios and from Fe3+/Fe2+ ratios in spinel inclusions range from NNO-1.1 to NNO+1.8. The origin of S-rich magmas is unclear. We discuss (1) partial melting of a mantle source at relatively oxidized fO2 conditions, and (2) magma contamination by seawater either directly or through magma interaction with seawater-altered Jurassic oceanic crust. The intermediate sulfur group inclusions represent undegassed or slightly degassed magmas similar to submarine OIB glasses, whereas the low-sulfur group inclusions are likely to have formed from magmas significantly degassed in near-surface reservoirs. Mixing of these degassed magmas with stored volatile-rich ones or volatile-rich magma replenishing the chamber filled by partially degassed magmas may produce hybrid melts with strongly varying S concentrations and S6+/S_total ratios.

Solving Differential Equations by Parallel Laplace Method with Assured Accuracy

The paper has been presented at the 12th International Conference on Applications of Computer Algebra, Varna, Bulgaria, June, 2006

Impulsive Partial Hyperbolic Functional Differential Equations of Fractional Order with State-Dependent Delay

MSC 2010: 26A33, 34A37, 34K37, 34K40, 35R11

Reconstruction of inclusions in photothermal imaging

Photothermal imaging allows to inspect the structure of composite materials by means of nondestructive tests. The surface of a medium is heated at a number of locations. The resulting temperature field is recorded on the same surface. Thermal waves are strongly damped. Robust schemes are needed to reconstruct the structure of the medium from the decaying time dependent temperature field. The inverse problem is formulated as a weighted optimization problem with a time dependent constraint. The inclusions buried in the medium and their material constants are the design variables. We propose an approximation scheme in two steps. First, Laplace transforms are used to generate an approximate optimization problem with a small number of stationary constraints. Then, we implement a descent strategy alternating topological derivative techniques to reconstruct the geometry of inclusions with gradient methods to identify their material parameters. Numerical simulations assess the effectivity of the technique.

The resistance of PRESENT-80 against related-key differential attacks

We examine the security of the 64-bit lightweight block cipher PRESENT-80 against related-key differential attacks. With a computer search we are able to prove that for any related-key differential characteristic on full-round PRESENT-80, the probability of the characteristic only in the 64-bit state is not higher than 2−64. To overcome the exponential (in the state and key sizes) computational complexity of the search we use truncated differences, however as the key schedule is not nibble oriented, we switch to actual differences and apply early abort techniques to prune the tree-based search. With a new method called extended split approach we are able to make the whole search feasible and we implement and run it in real time. Our approach targets the PRESENT-80 cipher however,with small modifications can be reused for other lightweight ciphers as well.