1 

Solving Stochastic PDEs with Approximate Gaussian Markov Random Fields using Different Programming Environments
This thesis is a study on the implementation of the Gaussian Markov Random Field (GMRF) for random sample generation and also the Multilevel Monte Carlo (MLMC) method to reduce the computational costs involved with doing uncertainty quantification studies. The GMRF method is implemented in different programming environments in order to evaluate the potential performance enhancements given varying levels of language abstraction. It is seen that the GMRF method can be used to generate Gaussian Fields with a Mat{\'e}rn type covariance function and reduces the computational requirements for large scale problems. Speedups of as much as 1000 can be observed when compared to the standard Cholesky Decomposition sample generation method, even for a relatively small problem size. The MLMC method was shown to be at least 6 times faster than the standard Monte Carlo method and the speedup increases with grid size. It is also seen that in any Monte Carlo type methods, a Krylov subspace type solver is almost always recommended together with a suitable preconditioner for robust sampling.
This thesis also studies the ease of implementation of these methods in varying levels of programming abstraction. The methods are implemented in different languages ranging from the most common language used by mathematicians (MATLAB), to the more performance oriented language (C++PETSc/MPI), and ends with one of the newest programming concept (ExaStencils). The GMRF method featured in this thesis also is one of the earliest application to be implemented in ExaStencils.

[PDF]
[Abstract]

2 

Unit Root Testing for AR(1) Processes
The purpose of this study is to investigate the asymptotics of a first order auto regressive unit root process, AR(1). The goal is to determine which tests could be used to test for the presence of a unit root in a first order auto regressive process. A unit root is present when the root of the characteristic equation of this process equals unity. In order to test for the presence of a unit root, we developed an understanding of the characteristics of the AR(1) process, such that the difference between a trend stationary process and a unit root process is clear.
The first test that will be examined is the DickeyFuller test. The estimator of this test is based on Ordinary Least Square Regression and a ttest statistic, which is why we have computed an ordinary least square estimator and the test statistic to test for the presence of unit root in the first order auto regressive process. Furthermore we examined the consistency of this estimator and its asymptotic properties. The limiting distribution of the test statistic is known as the DickeyFuller distribution. With a Monte Carlo approach, we implemented the DickeyFuller test statistic in Matlab and computed the (asymptotic) power of this test. Under the assumption of Gaussian innovations (or shocks) the limiting distribution of the unit root process is the same as without the normality assumption been made. When there is a reason to assume Gaussianity of the innovations, the Likelihood Ratio test can be used to test for a unit root.
The asymptotic power envelope is obtained with help of the Likelihood Ratio test, since the NeymanPearson lemma states that the Likelihood Ratio test is the point optimal test for simple hypotheses. By calculating the likelihood functions the test statistic was obtained, such that an explicit formula for the power envelope was found. Since each fixed alternative results in a different critical value and thus in a different unit root test, there is no uniform most powerful test available. Instead we are interested in asymptotically point optimal tests and we will analyze which of these point optimal tests is the overall best performing test. By comparing the asymptotic powercurve to the asymptotic power envelope for each fixed alternative we could draw a conclusion on which fixed alternative results in the overall best performing test.
On the basis of the results of this research, it can be concluded that there does not exist a uniform most powerful test, nonetheless we can define an overall best performing test.

[PDF]
[Abstract]

3 

Optimal boundary point control for linear elliptic equations
This master thesis is concerned with optimal boundary point and function control problems for linear elliptic equations subject to control constraints. The elliptic partial differential equation with Robin boundary condition is considered. The control is chosen as a linear combination of the Dirac delta functions in the point control problem. The weak formulation and optimality conditions are obtained for the function control problem. The main goal is to examine the existence of the weak solution and to derive optimality conditions of boundary point control. Introducing sufficient discretization methods such as a finite volume and a finite element methods, we obtain finite dimensional problem.
We apply efficient numerical methods including primaldual active set strategy, projected gradient and conjugate gradient methods. The test examples are presented clarifying the performance of numerical methods mentioned earlier. The conjugate gradient method is compared to the unconstrained matlab function QUADPROG and primaldual active set strategy is compared to the constrained QUADPROG. The projected conjugate gradient method is applied to improve the projected gradient method. Due to its importance, the sparse point and function control problems are studied. We apply primaldual active set strategy where the sparse parameter f is chosen differently. All of the results are presented in both point and function control problems.

[PDF]
[Abstract]

4 

Estimating the extreme value index for imprecise data
In extreme value theory the focus is on the tails of the distribution. The main focus is to estimate the tail distribution for a rounded data set. To estimate this tail distribution the extreme value index should be estimated, but due to the rounded data this extreme value index oscillates heavily. Therefore a correct estimate can not be obtained. By adding a small uniform stochast the rounded data can be smoothend out and in this way the oscillation can be cancelled.

[PDF]
[Abstract]

5 

Investigation of Different Solvers for Radiotherapy Treatment Planning Problems
Radiotherapy treatment planning involves solving inequality constrained minimization problems. The currently used interior point solver performs well, but is considered relatively slow. In this thesis we investigate two different solvers based on the logarithmic barrier method and Sequential Quadratic Programming (SQP) respectively. We argue that the behaviour of the logarithmic barrier solver is uncertain, thereby making it generally unreliable in this context. In addition we substantiate that the performance of the SQP solver is solid, but lacks efficiency in computing the minimizers of its related quadratic subproblems.
We conclude that without serious improvements, none of the solvers investigated are faster than the currently used interior point optimizer.

[PDF]
[Abstract]

6 

Anomaly detection for internet banking using supervised learning on high dimensional data
Nowadays, a high number of transactions are performed via internet banking. Rabobank processes more than 10 million transactions per day. Most of these transactions are (part of) normal behaviour. On the other hand, some transactions are considered to be out of the ordinary. These anomalous events occur relatively infrequently (less than 10 per day). Employees, that try to find these anomalous events, combine the transactions data, historical knowledge of the anomalous events and their expertise to detect and quantify them. Several types of anomalies are considered to be interesting and so they are labelled. These anomalies need to be detected, so they can be prevented in the future. The employees try to find events similar to known anomalies. Characteristics of anomalies change over time and employees also need to detect this slightly changed, but similar, behaviour. It is not our goal to detect completely new types of anomalies. In this thesis, the focus lies on finding events similar to the known anomalies. In order to assist these employees, a model that uses the transaction data and incorporates known anomalous events is built. Our model is able to score new incoming transactions and use these to update the model parameters. The scores can be returned to the employees to assist them in finding transactions that are similar to a particular type of anomaly. The AdaGrad algorithm with diagonal matrices is used. Also, l1regularization is used on the parameter to create a more sparse solution.

[PDF]
[Abstract]

7 

Simulating sprouting angiogenesis: using a new 3D substrate dependent cellbased model
Angiogenesis1 is the biological mechanism by which new blood vessels sprout from existing ones. It differs from vasculogenesis, which is the de novo growth of the primary vascular network from initially dispersed endothelial cells (ECs). Vasculogenesis is predominant in embryonic tissue whilst new vasculature in the adult body arises mostly from angiogenesis. ECs, lining the inside of blood vessels, react to different angiogenic stimuli and inhibitors. Among the stimuli is the vascular endothelial growth factor (VEGF) which is upregulated in tissue where the vascular structure is damaged or insufficiently developed to meet oxygen demand.
The identification of the processes involved in angiogenesis is quite recent and has stirred increased interest in therapeutic and clinical applications according to Carmeliet et al. [1]. One can think of tissue repair in wound beds, inhibition of growth of tumorous tissue or vascular reform during the female reproductive cycle. Rossiter et al. [2] showed that VEGF induced angiogenesis is crucial for wound healing in an experiment where wounds were inflicted upon normal and VEGFdeficient mice. New vasculature ensures supply of oxygen and lymphocytes and disposal of carbon dioxide and lactates, accelerating wound healing and tissue reconstruction. The increased creation of new vasculature around tumorous tissue is believed to follow the same process and inhibiting angiogenesis is therefore an important topic in clinical studies on cancer treatment.
Biochemical laboratory experiments can be hard, time consuming, expensive or unethical. Computational models can be used to provide an easy, quick and cheap way to get insights that would otherwise require laboratory experiments. The understanding of biological processes needs quantification and in this sense mathematical formulation of the relations involved becomes useful. Their mathematical interpretation and experimental verification is an iterative process resulting in better understanding of the process itself. Computer simulation will never make laboratory experiments obsolete, but it can provide guidance in targeting viable hypotheses before conducting in vitro or in vivo experiments.
Mathematical modeling of biological cellular processes dates back to the simulation by Glazier and Graner in 1992. They describe natural sorting behavior of different cell types [3] and different rearrangement patterns driven by the differential adhesion hypothesis [4]. This hypothesis states that cells of different types have specific potential energies upon adhesion, driving sorting behavior. In these simulations, the cellular Potts model2 (CPM) is used. A CPM for vasculogenesis based on this work was made byMerks et al. [5, 6] in which a layer of partial differential equations (PDEs) models the chemoattractants. Later, Merks added Vascular Endothelial cadherin (VEcadherin) caused contactinhibited chemotaxis to simulate angiogeniclike sprout formation [7]. From an initial clump of ECs in the model sprouting behavior appears. Merks postulates that both vasculogenesis and angiogenesis must be driven by the same principles. To produce these results, a generic library called the Tissue Simulation Toolkit (TST) was written in C++ starting from 2004 modeling the CPM described by Glazier et al. [4] in a generic way. Merks [7] extensively describes the advantages of a cell based approach over a continuum approach that is widely used in mathematical biology. Although his CPM is a nice method that increases insight in the angiogenic process, it is computationally heavy, limiting the scalability of the tractable problem domain.
Vermolen and Gefen [8] described tissue behavior using a semistochastic cellbased formalism to model the migration of cells in colonies in the context of wound healing, tumor growth, bone ingrowth and contraction formation. Movement of cells is assumed to be the result of a strain energy density working as a mechanical stimulus. Like the CPM, the model tracks displacement and viability of individual cells.
The aim of this study is to adapt this semistochastic cellbased formalism to describe angiogenesis, hence connecting this modeling approach to the subject ofMerks’ work. The need for such a model is clearly stated in the discussion of Vermolen’s work [9]. Thanks to the computational less heavy character in comparison with the CPM, we hope to be able to simulate larger areas to get a better glance at large scale behavior whilst still being able to benefit from the cellbased character of the model. We also improve the biochemical model for the degrading of the substrate by the cells and formulate all relevant parameters based on local properties.
The challenge is to translate the advantages of Merks’ CPM, like cell shape specific behavior, tracking of elongation patterns and cellcell contact behavior, to this new formalism without compromising the computational simplicity.
To verify our simulation results with biochemical experiments, this study is performed in collaboration with the Dermatology Department of the VU Medical Center. This department does in vitro laboratory research on many processes that occur in the skin, for example the role of endothelial cells during skin wound healing.
The first aim of this research is tomimic their in vitro angiogenesis sprouting assay using our computational model, simulating the response to different chemical stimuli like VEGF. Formulating a way to visually and numerically compare the laboratory work to the simulated results is key to making the model applicable in practice.

[PDF]
[Abstract]

8 

Investigation of Different Solvers for Radiotherapy Treatment Planning Problems
Radiotherapy treatment planning involves solving inequality constrained minimization problems. The currently used interior point solver performs well, but is considered relatively slow. In this thesis we investigate two different solvers based on the logarithmic barrier method and Sequential Quadratic Programming (SQP) respectively. We argue that the behaviour of the logarithmic barrier solver is uncertain, thereby making it generally unreliable in this context. In addition we substantiate that the performance of the SQP solver is solid, but lacks efficiency in computing the minimizers of its related quadratic subproblems. We conclude that without serious improvements, none of the solvers investigated are faster than the currently used interior point optimizer.

[PDF]
[Abstract]

9 

Optimal distributed point control of linear elliptic equations
In this master thesis, we consider the optimal function and point control problems governed by linear elliptic partial differential equations together with bilateral control constraints. The aim of this work is to choose a control function by linear combination of the Dirac delta and solve the optimal distributed point control problem. In this work, optimality systems of point and function control problems are derived by Lagrangian principle and reduced functional respectively. The optimality system is discretized by the finite element method (FEM) and finite volume method (FVM). We apply a semismooth Newton (SSN) method which is equivalent to a Primal dual active set strategy (PDASS) to solve the discretized optimality system. As a second solution method, we propose Projected gradient (PG) method for the same problem. For each method, we compare the result of FEM and FVM and give a preference. In order to have the best solution method, we compare the results of PDASS and PG methods to the results of matlab command QUADPROG. We have to solve two partial differential equations in every iteration, namely the state and the adjoint equations. Therefore, we develop a Multigrid Preconditioned Conjugate Gradient (MGPCG) method for solving the discretized optimality system as fast as possible. Finally, we consider the linear elliptic optimal control problems with L1 norms in the cost functional which results sparse control problem. Due to L1 norm, objective functional becomes nondfferentiable and the optimal controls are identically zero on large parts of the control domain. Using an appropriate smoothing of the nondifferentiable terms for the cost functional, we solve the optimal sparse control problems theoretically and numerically.

[PDF]
[Abstract]

10 

Efficiency Improvement of Panel Codes
Panel codes are used by the Maritime Reseach Institute Natherlands (MARIN) to compute flows around ships and propellers. These codes are based on Boundary Element Methods (BEM). A known drawback of BEM is that it forms dense linear system of equations that have to be solved. By improving the efficiency of the dense linear solver, the computational time required by panel codes can be significantly reduced. Since applications of panel codes in MARIN include automatic optimization, where a large number of hull forms or propeller geometries have to be evaluated, the reduction of computational time is important.
Four strategies were explored to improve the performance of the dense linear solver. First, to replace the current GMRES solver with IDR(s). Second, the updating of a fixed size block Jacobi preconditioner into a variable size block Jacobi preconditioner. Third, to use a hierarchical matrixvector multiplication in the solver instead of dense matrixvector multiplication. Lastly, to replace the block Jacobi preconditioner with a hierarchicalLU preconditioner. Out of the four strategies, the use of hierarchicalLU preconditioner was found to speed up the dense linear solver substantially, especially for large systems. The use of IDR(s) instead of GMRES is also recommended as it removes the problems introduced by the need to restart.
This report discusses the theory, implementation and test results obtained from the four strategies aforementioned. As a result of this project, the use of IDR(s) combined with hierarchicalLU preconditioner is recommended to be implemented in the panel codes.

[PDF]
[Abstract]

11 

Stochastic Modeling of Order Book Dynamics
In this project the order book model proposed by Cont et al. [10] is used as a starting point to model order book dynamics. This model nicely combines three desirable properties from earlier studies: it is easy to calibrate, it reproduces statistical properties of the order book and it allows to make analytical computations in the order book. The model is studied, calibrated and tested on realtime data from the London Stock Exchange. Possible improvements to the model are discussed and tested. A method to compute probabilities in the model will be presented: recovering densities by inverting continued fraction representations of Laplace transforms. This is also implemented and evaluated.

[PDF]
[Abstract]

12 

History Matching via Ensemble Kalman Filter of Norne Field
The very popular method, namely, Ensemble Kalman Filter is used to history match the Esegment of the Norne Field reservoir. A complete history matching workflow for a real reservoir requires to focus on a variety of issues involved. According to EnKF method there is a need to create initial ensemble and since the part of the field is considered the issue related to boundary conditions appears.

[PDF]
[Abstract]

13 

De geïmpliceerde boom en de scheefheid van BlackScholes

[PDF]

14 

Heuristics of heavytailed distributions and the Obesity index
In this thesis a new metric for the heavytailedness of dataset is proposed.

[PDF]
[Abstract]

15 

A rainmask for the ROBIN Lite radar system
In this thesis we describe the development of a new rainmask for the ROBIN Lite radar system.

[PDF]
[Abstract]

16 

Simulatie van wachtrijen met flexibele bedieningssnelheid; Simulation of queues with adaptable service speed
Three different models of singleserver queues with adaptable service speed are discussed and Monte Carlo simulations of these models are presented. Moreover, some results are given and some stochastic orderings between steadystate workload distributions are proven.

[PDF]
[Abstract]

17 

Dantzig Selector: Het schatten van een sparse parametervector met weinig metingen  Dantzig Selector: Statistical estimation of a sparse vector with undersampled data
Bij lineaire regressie geldt in het algemeen, dat het aantal metingen (n) groter is dan het aantal parameters (p). In praktische situaties kan het echter voorkomen, dat het aantal metingen veel kleiner is dan het aantal parameters. De kleinste kwadratenschatter kan in dat geval niet gebruikt worden om de parametervector te schatten, er is geen unieke oplossing voor de normaalvergelijkingen. Wanneer echter gegeven is dat een klein aantal parameters nietnul is, is het soms toch mogelijk om de parameters te schatten.
In 2005 introduceerden Emmanuel Candès en Terence Tao de Dantzig Selector (DS) als een schatter voor de parametervector in deze situatie. In dit verslag wordt uitgelegd hoe de Dantzig Selector werkt. Tevens wordt op theoretische en praktische wijze nagegaan hoe goed de DS werkt.

[PDF]
[Abstract]

18 

Modelling and simulation of bone implant healing
In this report a model for the ingrowth of a prosthesis in a bone will be formulated and simulated. First an overview is given of all the biological processes that are part of this healing process and what external factors can influence this process.
Then a first mathematical model is presented in which the mechanical stimuli (one of the external factors that can influence the healing process) are neglected. The model is solved by numerical means with both the finite volume method as the finite element method. For the finite element method a short introduction is first given to get familiar with this technique.
Results for this model will be presented, followed by a short discussion about the results and a conclusion.
Subsequently the previous model is extended to incorporate mechanical stimuli, this is done by combining it with an elasticity equation. A short introduction will also be given to the theory of linear elasticity. The model is then solved using finite element analysis and finally also the results for this extended model will be presented, together with once again a short discussion and conclusion.

[PDF]
[Abstract]

19 

Non parametric Bayesian belief nets (NPBBNs) versus ensemble Kalman filter (EnKF) in reservoir simulation
The thesis is concerned with applying a method that uses graphical models to solve a reservoir simulation problem. The focus is on the graphical models called Bayesian belief nets. A Bayesian belief net based approach is compared with an already popular technique used in reservoir simulation, the ensemble Kalman filter.

[PDF]
[Abstract]

20 

Rainwindinduced vibrations of cables

[PDF]
