▼ 

Ice flow, forced by gravity is governed by the Full Stokes (FS) equations, which are computationally intensive to solve due to their non-linearity. Therefore, it has been unavoidable to approximate the FS equations when modeling growth and collapse of an ice sheet-shelf system, which requires simulating many thousands of years. However, the most popular Shallow Ice Approximation (SIA) and Shallow Shelf Approximation (SSA) are only accurate in certain parts of an ice sheet, both excluding the grounding line where the ice starts floating. Using the Finite Element software Elmer/Ice, SIA and SSA are dynamically coupled to FS aiming to maintain high precision and reduce computation time. An existing coupling of SIA and FS, called ISCAL [Ahlkrona et al., 2013a], is investigated for robustness. It is shown that instabilities in the FS solution limit ISCAL’s robustness. A novel way of iteratively coupling SSA and FS has been implemented into the open source Finite Element software Elmer/Ice and applied to both 2D and 3D conceptual marine ice sheets. The SSA-FS coupling shows to be very accurate, for both diagnostic and prognostic runs (error in velocity respectively below 0.5% and 5%). Grounding line dynamics of the SSA-FS coupling are similar to the FS model under a periodical forcing in a simulation covering 3000 years. The current implementation does not yield speed up in 2D due to inefficient assembly of the matrices to be solved. In 3D, the cpu time is reduced to two thirds of the cpu time of the FS model despite inefficient assembly. The total number of FS iterations in the SSA-FS coupling is comparable to the FS model, showing a large potential of reducing computation time since computation time of the SSA model is up to 3% of the FS model’s computation time when applied to the same ice shelf ramp. In future research, the SSA-FS coupling can be combined with ISCAL, but this requires both efficient implementation of the SSA-FS coupling and improved stabilization methods for FS.

Applied Mathematics

Advisors/Committee Members: van Gijzen, Martin (mentor), Delft University of Technology (degree granting institution).

▼ 

In this thesis, the adaptation of the spectral clustering algorithm to the privacy preserving domain was investigated. The spectral clustering algorithm divides data points into clusters according to a measure of connectivity. A pivotal part of spectral clustering is the partial eigendecomposition of the graph Laplacian. Two numerical algorithms are used to approximate the eigenvectors of the Laplacian: the Lanczos algorithm and the QR algorithm. These numerical methods are adapted to work on Paillier encrypted data values. The main challenge is the fact that Paillier encryption can only be applied to a field of positive integers. Moreover, cryptographic protocols have to be invoked to perform non-linear operations. Also, the square root and division operations are computationally heavy in the privacy preserving domain. The numerical algorithms are adapted to overcome these challenges and be more suitable to work on encrypted values.

Applied Mathematics

Advisors/Committee Members: van Gijzen, Martin (mentor), Veugen, Thijs (mentor), Delft University of Technology (degree granting institution).

▼ 

Magnetic resonance imaging (MRI) is a medical diagnostic imaging technology characterized by its high cost and infrastructure demands. This makes it unsuitable for implementation in low- and medium-income countries. Its costs are directly related to the strength of the system’s magnetic field. Therefore, an MRI scanning technology dedicated for use in low- and medium-income countries has been developed by Delft University of Technology and Leiden University Medical Center that makes use of a low magnetic field. This low-field technology, that requires a smaller size magnet and yields lower resolution images, is already capable of diagnosing infant hydrocephalus, a frequently occurring condition in sub-Saharan Africa. For a successful implementation, it is crucial that this technology is integrated in an MRI system that adequately addresses challenges within the context of sub-Saharan Africa and is equipped with all necessary functionality to diagnose infants with hydrocephalus. From observations and interviews in the healthcare environment of Uganda, it was found that – in addition to high cost – poor infrastructure and challenges that occur during imaging procedures limit the implementation- and use-potential of current MRI systems and thus can affect the implementation of this new system. Poor pathways and unavailability of forklifts or cranes often cause breakage during installation and transport. Moreover, lack of space in hospitals hinders implementation of large systems that require dedicated, shielded rooms. During imaging, movement of the (infant) patient may affect the image, and often mothers are kept in the imaging room with the patient to calm her. Other challenges for the operation of these systems in this context include heat, humidity, dust, power outages and the presence of pests. A concept of a complete MRI system that integrates the developed technology and addresses these challenges is designed. This MRI system enables an ergonomic and safe operation, facilitates safe and comfortable insertion and removal of the patient from the machine, and allows mothers and operators to remain close and in visual contact with the patient. The entire concept takes up no more than 3 m2, can be installed in any room with a power socket and features a plug and play installation. It can be carried over rough pathways without the need for a forklift or crane. Furthermore, the concept features a bed with sliding mechanism that optimizes the space inside the system, allowing hydrocephalus patients with significant expansion of the head to be imaged. With this concept, the first steps are taken towards implementation of the technology in hospitals in low- and medium-income countries. Through further development and testing in close collaboration with stakeholders, the concept can be optimized to enable diagnosis and treatment of infants with hydrocephalus and eventually many others, and hereby contribute to increasing the (e)quality of healthcare around the world.

Integrated Product Design

Advisors/Committee Members: Diehl, Jan-Carel (mentor), Bakker, Martien (graduation committee), van Gijzen, Martin (graduation committee), Delft University of Technology (degree granting institution).

▼ For a proper three-dimensional reconstruction of histology serial sections, adjustment of the slices is necessary for combining serial sections. Deformations occur due to the sectioning and acquisition pro-cess of the microscopic analysis of histology. By reconstructing the deformations with a transformation of the sections, a mathematical correction on the images can be applied. By using image registration, a transformation function is searched for to minimize differences between the histology slices. Due to the non-linearity of the distortions, prior knowledge is required in order to have a solvable problem. Additional information in the form of elasticity regularisation is considered. Implementing the elastic regularisation with a finite element method, provides a continuous transformation function With the continuous function, in a natural way, alignment can be monitored for folding transformations. In this work, the (bi-)linear and (bi-)quadratic elements for the finite element method are implemented and compared with the finite difference method. It is observed that for the different kinds of elements, the (bi-)linear elements yield best results with the validity of the transformation. Moreover, the computa- tional costs for the bi-linear elements are the cheapest. Compared with the finite difference method, the differences in accuracy are not noteworthy but the computational time of the finite element method is longer. Furthermore, to steer the matching in an accurate direction, improvements are proposed by applying local stiffness of the elements or adding soft constraints on the volume of the elements. This results in significant improvements in the transformation. For these two approaches it is observed that local stiffness is more restrictive than volume-preserving. Solving the optimization problem, a Gauss-Newton method to search for descent directions is applied. A matrix-based and a matrix-free approach of elasticity regularisation is considered in the linear system of finding a descent direction. While the matrix-free approach decreases the memory usage, the computational costs are significantly increased. Advisors/Committee Members: van Gijzen, Martin (mentor), Modersitzki, Jan (mentor), Dubbeldam, Johan (graduation committee), Remis, Rob (graduation committee), Delft University of Technology (degree granting institution).

▼ The TU Delft and the LUMC are creating a low-field MRI scanner to use in third world countries. This type of MRI scanner has many advantages, but a downside is the amount of noise in the obtained images. A way to reduce noise is diffusion filtering. This thesis discusses the theory of some linear and nonlinear diffusion filtering methods and tests them on several test problems. The methods range from the basic linear method, the heat equation, to the more advanced Perona-Malik method. The results indeed show that the Perona-Malik method generates better results, sometimes when combined with a Gaussian kernel to decrease its ill-posedness. Two numerical methods have been compared for these anisotropic diffusion filtering methods: the Forward Time, Central Space method and the Additive Operator Splitting method. The advantage of the AOS method is the unconditional stability for al positive time step sizes, while FTCS implementation is only stable for time steps smaller than 0.25. The results for the AOS method were also slightly better than for the FTCS method and almost never lead to instability. The FTCS method showed instability more often. The parameter choice is of significant importance for the outcome. Several options for determining the gradient threshold parameter K, time step size and stopping time S have been investigated. A new estimation for the time step size and stopping time S have been proposed and compared. The results are images with good visual quality for both introduced methods. Also, an adaptive time step has been tested to overcome the instability for unstable methods and this seems to be succesful. Lastly, region growing segmentation has been applied to images obtained with the proposed stopping time S. Segmentation is an additional way to investigate the quality of these outcomes. For most test problems the partitioning resembles the partitioning of the original image. However, choosing certain parameters stays a challenge and can be of interest for further examination. It also remains a goal to investigate the discussed methods on data obtained with the actual TUD/LUMC MRI scanner. Advisors/Committee Members: van Gijzen, Martin (mentor), Vuik, Kees (graduation committee), Remis, Rob (graduation committee), Delft University of Technology (degree granting institution).

▼ It has become clear that the global energy system needs to shift away from fossil fuels towards clean, renewable energy sources. Experience in the past decades has shown that at-plate solar heat panels can play a role in the energy system of the future, in particular in solar collector fields for district heat generation. This thesis concerns the dynamical modelling of such solar district heat plants (SDHPs). As the upswing of SDHPs occurred rather recently, research on their optimal control is ongoing. The main challenge for control is to adapt to fluctuations in solar radiation and other inputs in an optimal way, ensuring a high and stable outlet water temperature to the grid while minimizing flow fluctuations. Many modelling efforts of single collectors as well as full solar heat fields have been reported, although mostly in the form of detailed physical models. In this thesis, a new approach for describing the dynamics of a large at-plate solar field is proposed. We develop a continuous-discrete stochastic state space model suitable for prediction and control. This model form combines knowledge from physics and information from data, thereby allowing for relatively simple formulations while modelling complex dynamics. Retaining a physically meaningful model formulation has additional advantages for model development, as analysis of residuals and outputs can provide information on suitable model extensions. A basic model was formulated and systematically extended in a forward selection procedure, using a.o. likelihood ratio tests. The model development was based on May 2017 measurements from an Aalborg CSP solar heat plant in the municipality of Solrød , Denmark. The models are implemented using the R-package CTSM-R, which includes parameter estimation methods based on maximum likelihood estimation and the extended Kalman filter. It is found that the developed model is suitable for very short term (minutes to hours ahead) to short term (day ahead) prediction, as needed for control and heat output forecasting of a SDHP. It includes several new aspects compared to existing models, the most notable being non-parametric modelling of shading effects and a split of total radiation into diffuse and direct components. Including these elements improves model predictions considerably, and allows for asymmetric panel effciency over the day. Detailed analysis of the model's predictive performance is provided, including a comparison to current ISO standard model and the current Solrød control scheme, as well as a cross-validation on data from different seasons. In future work, the model's performance when using input predictions from weather forecasts should be tested. The model should further be used in a model predictive control scheme in order to improve current SDHP control strategies. This would lead to a smoother pump operation, thereby reducing electricity consumption, costs, and greenhouse gas emissions. Advisors/Committee Members: Heemink, Arnold (mentor), Bacher, Peder (mentor), van Gijzen, Martin (graduation committee), Delft University of Technology (degree granting institution).

▼ Combating air pollution has proven to be a difficult task for countries with rapidly developing economies. Poor air quality can be hazardous to people doing any outdoor activities. So being able to make accurate, short term air quality predictions can be very useful. However, making these predictions has proven to be quite difficult, since there are a lot of different physical and chemical processes involved in the emission and transport of the various aerosols that contribute to air pollution. So instead of the more traditional Chemical Transport Models (CTMs) we will be using neural networks in order to make predictions of one of these aerosols, PM2.5. In particular, we will be using a Long Short Term Memory (LSTM) network. In addition, we will include the simulations results from a CTM, LOTOS-EUROS, as input data to the LSTM network to improve the performance of the neural network. One of the main drawbacks of the LSTM approach is that whenever the PM2.5 concentration changes a lot, the predictions made by the LSTM network take some time to change as well, causing a visible time delay when looking at the measurements and predictions in the same time series plot. We will also try a simpler type of neural network, a Feedforward Neural Network (FNN) and compare its performance to that of LSTM. We found that using the simulation data does indeed improve the LSTM network. Not only in terms of the loss function used by the neural network and, but in particular in the amount gross overestimations by the network, which we use to quantify the LSTM time delay problem. We also found that FNN outperforms the LSTM approach, in particular on samples of high PM2.5 concentrations, which we argue is primarily caused by a low amount of samples in our dataset. Advisors/Committee Members: Lin, Hai Xiang (mentor), Heemink, Arnold (mentor), van Gijzen, Martin (graduation committee), Delft University of Technology (degree granting institution).

▼ 

Structure Optimization has been an important subject with many applications for centuries. In the last sixty years, numerical optimization has facilitated large advancements in this field. One of the areas in Structure Optimization is Topology Optimization, which is used for Additive Manufacturing purposes. In this thesis we explore Static and Dynamic Topology Optimization. In the optimization problems matrix equations of the type Ax = b, with A sparse and badly conditioned, are accelerated using deflation techniques in addition to preconditioning. We have applied several iterative methods, preconditioners, and deflation types to the topology optimization problems. The static problem concerns compliance optimization of a two-dimensional MBB-beam. The deflation type that reduced the number of iterations the most without introducing large overhead costs was rigid body modes deflation, divided over element squares. It was found that dividing the rigid body modes vectors over density based regions did not reduce the number of iterations. The dynamic problem concerns eigenvalue optimization of a three-dimensional moving wafer stage that is used for laser-printing computer chips. The optimization formulation contains a shifted eigenvalue problem that is solved using model order reduction. In the computations of bases for the reduction matrix equations appear, to which deflation techniques were applied. All deflation types reduced the number of iterations and needed time to solve the matrix equations. The best deflation type was using rigid body modes (RBM) divided over element cubes combined with eigenvectors from the previous iteration. The next best was the same combination without the division over cubes. Using eigenvectors or RBM separately were the least effective deflation types. There is an optimal amount of element cubes to use when dividing the RBM. The tests on a few grid sizes suggested a quadrupling of the amount when the grid size doubles, but more research is needed to really identify a relation between the grid size and optimal amount. When increasing the grid size to a level where parallel computing on a cluster was required, the deflation type using element cubes could not be used due to complications with parallel implementation. Using the deflation type RBM and eigenvectors, reductions by a factor of 1.60 and 1.75 in the total needed time for 150 optimization iterations were achieved for grid sizes 120x120x20 and 180x180x30, respectively. In the first case the objective function of the optimized design converged further in the same amount of iterations when using deflation. future research could include using the promising deflation type RBM in cubes + eigenvectors in parallel computing to obtain an even larger time reduction in the optimization of the wafer stage with large grid sizes.

Applied Mathematics

Advisors/Committee Members: van Gijzen, Martin (mentor), Astudillo Rengifo, Reinaldo (mentor), Langelaar, Matthijs (graduation committee), van Horssen, Wim (graduation committee), Delft University of Technology (degree granting institution).

▼ 

The numerical simulation of brittle failure with nonlinear finite element analysis (NLFEA) remains a challenge due to robustness issues. These problems are attributed to the softening material behaviour and the iterative nature of the Newton-Raphson type methods used in NLFEA. However, robust numerical simulations become increasingly important, for example due to recent developments in Groningen.  To address these issues, sequentially linear analysis (SLA) was developed which exploits the fact that a linear analysis is inherently stable. By assuming a stepwise material degradation the nonlinear response of a structure can be approximated with a sequence of linear analyses. Although this approach has been proven to be effective for several case studies, the numerical performance is still a problem that has to be solved. After every linear analysis, a single element is damaged resulting in incremental damage. As a result, the system of equations only changes locally between these linear analyses. Traditional solution techniques do not exploit this property and calculate a matrix factorisation every linear analysis, resulting in high computational times per analysis step. Since SLA typically requires many linear analyses to obtain the desired structural response, this leads to unacceptable analysis times. The aim of this thesis is to improve the computational performance of SLA by developing numerical solution techniques which exploit the incremental approach of SLA. To this extend, the following methods have been developed. A direct solution technique has been developed which is based on the Woodbury matrix identity. This identity allows for the numerically cheap computation of the inverse of a low-rank corrected matrix. In this approach, the expensive matrix factorisation does not have to be calculated every linear analysis step. Instead, the old factorisation can be reused along with some additional matrix- and vector multiplications and solving a significantly smaller linear system of equations. An optimal strategy is derived to determine at which point a new factorisation should be calculated. An improved preconditioner for the conjugate gradient (CG) method has been developed. Instead of an incomplete factorisation, the complete factorisation is used as a preconditioner which reduces the number of required CG iterations significantly. The point at which too many CG iterations are required and a new factorisation is necessary, is determined using the same strategy as the first method. From numerical experiments it follows that both methods perform significantly better than the direct solution method, especially for large 3-dimensional problems. The best performance is achieved using the Woodbury matrix identity resulting in the solver no longer being the dominant factor in SLA. Furthermore, significantly larger problems are not solvable in time frames in which previously only smaller problems were solved.

Applied Mathematics

Advisors/Committee Members: van Gijzen, Martin (mentor), Schreppers, G.M.A (graduation committee), Rots, Jan (graduation committee), van Horssen, Wim (graduation committee), Delft University of Technology (degree granting institution).

▼ 

We have tried to create a bin-packing algorithm that assigns items from a customer-order to totes such that the amount of totes is minimized. Analyzing the bin-packing algorithm that was used before this thesis had been written, taught us that xxx.xx% of the customer-orders was packed non optimal. In this thesis four algorithms are applied to Picnic data. The order in which the algorithms assign items to a tote has major consequences for the solutions. Eight different ways to order the items are combined with each algorithm, resulting in 32 different tote-calculations. Out of those 32 tote-calculations, the Best Fit Algorithm with items ordered in decreasing normalized values generates the best results. Remarkably, ordering items randomly also gives good solutions. This brought us to introducing a new method, where each customer-order is calculated at most eight times, each time shuffling the items before rerunning the algorithm and remembering the better solution. This heuristic is optimal for xxx.xx% of the customer-orders.

Applied Mathematics

Advisors/Committee Members: de Oliveira Filho, Fernando (mentor), Aardal, Karen (graduation committee), van Gijzen, Martin (graduation committee), Gorte, Frank (graduation committee), Delft University of Technology (degree granting institution).

▼ 

In this thesis a deterministic wave model is used to reconstruct and predict the sea surface motion from FMCW (Frequency Modulated Continuous Wave) radar data, produced by Radac. The deterministic model that is used to do this is based on the linear wave theory. The radar is looking horizontally straight towards the waves in 5 separate beam directions of -40,-20,0, 20 and 40 degrees. Using the FMCW principle the backscatterd signal is converted into velocity and spatial range information. After some compensations (current for example) this velocity data can be treated as horizontal component of the orbital velocity of the wave. By using a least-squares solving approach (the trust-region reflective algorithm) on these orbital velocities and the expression that holds for them in the linear wave theory the model can be fitted to the measurements. The result of the least squares solver consists of a set of parameters for wave amplitude, phase and frequency. With these parameters the deterministic motion of the sea surface can be computed. This method is tested using artificial data and a generalized one directional case (using information from 1 beam under assumption of infinitely long-crested waves). For the experiments with artificial data consisting of waves with Hs = 2 meters (significant waveheight) the results are promising. A prediction time of 30 seconds over a range of 150 meters with an average error of 15 cm in the one directional model (fitted on 10 second data over 384 meters) can be achieved. For the multi directional model this lies between 20 and 30 seconds with an average error of 25 cm, depending on the spreading of the waves. Experiments with real data show less impressive results, an accurate reconstruction of the surface can be given, but the predictive capability is very limited.

Applied Mathematics | Applied Physics

Advisors/Committee Members: Heemink, Arnold (mentor), van der Vlugt, Rolf (mentor), Dubbeldam, Johan (graduation committee), van Gijzen, Martin (graduation committee), Delft University of Technology (degree granting institution).

▼ 

Delft University of Technology (TU Delft), Leiden University Medical Center (LUMC), Pennsylvania State University (PSU) and Mbarara University of Science and Technology (MUST) have an ongoing collaboration to create an affordable, portable and simplified version of the magnetic resonance imaging (MRI) scan for the CURE children’s hospital to diagnose children with hydrocephalus (water on the brain). As opposed to the conventional MRI scan, the low-field MRI prototype uses permanent magnets to create a magnetic field in the order of Milliteslas (mT). A downside of the low-field MRI application is the difficulty with spatial encoding due to small variations in the strength of magnetic field. This is a major problem for image reconstruction. The purpose of this research was to implement a deep learning (DL) network to overcome two of the major bottlenecks in image reconstruction for low-field MRI. These are the lack of real measured data for DL purposes, and the signal model associated with the low-field MRI. For DL purposes we generated synthetic data and acquired measured data. Each dataset consists of samples and each sample consist of an image and the corresponding signal. Due to technical limitations the measured dataset is small, 53 samples. To partially circumvent the problem, the data set was augmented to a total of 1908 samples. In addition, we used Transfer learning, which is a powerful method that applies knowledge gained from one problem to a different but related problem. We present three image reconstruction techniques, Model I, II, and III, based on convolutional and feedforward neural networks, which take MR signal data as input and directly and quickly outputs an image. We demonstrated that DL generates high quality images using synthetic data. In addition, we showed that Model III needs less training to reconstructs good quality images compared to Models I and III, respectively. Finally, Models I and III were unsuccessfully applied to real measured data. However, this study shows that neural networks are able to find a mapping between signal and image, therefore this idea can be extended to work on real measured data.

Applied Mathematics

Advisors/Committee Members: van Gijzen, Martin (mentor), de Leeuw den Bouter, Merel (graduation committee), Remis, Rob (graduation committee), Lin, Hai Xiang (graduation committee), Delft University of Technology (degree granting institution).

▼ The Wanda software package developed by Deltares can used for simulating both steady state and transient fluid flow in pipeline systems. Steady state simulations are used for initial system design and transient flow simulations are used for doing water hammer analysis in pipeline systems. For both types of simulations a system consisting of both linear and non-linear equations needs to be solved for the main unknown quantities, flowrate and head. This system is solved by linearising the equations using theNewton-Raphson method and solving the resulting system of linear equations. Currently, this is done by using a matrix solver from the proprietary IMSL numerical library which requires a paid license. The problem is that this solver sometimes either crashes or gets stuck in an infinite loop when dealing with singular matrices, while the proprietary nature of the library only allows for limited troubleshooting. The solution method therefore requires a redesign which should improve its robustness and maintainability. On the other hand, no ground should be yielded in terms of solution accuracy and performance. The singularity of the matrices are caused by quantities being underdetermined either due to user error in network design or phase changes such as a closing valve. In this report, a graph-theoretic approach is taken to detect these structural singularities in the form of determining a maximum size matching in a graph representing the system of equations. This approach gives information about which quantity is undetermined where in the pipeline system. As an alternative, condition number estimation is implemented. Furthermore, the IMSL library is replaced by LAPACK, which is a lightweight and versatile alternative. The open source nature of LAPACK and its permissive license ensure its maintainability. Since thematrices are banded, the band version of the LAPACK algorithms can be used. The new solutionmethod is compared to IMSL and evaluated in terms of robustness, accuracy and performance. The graph-theoretic method resolves the robustness issues of the IMSL-based method, while showing great performance. The information it gives is used to either correct the matrix and continue the simulation or output an appropriate error message, ensuring a user-friendly experience. Condition number estimation is too slow while also not being useful for further matrix correction purposes and is therefore disregarded. To improve the accuracy of LAPACK, iterative refinement is used. The maximum relative error measured over all the test cases was about 2.5%, resulting from an ill-conditioned case. Overall, the accuracy compared to IMSL is good, so that users will not be able to notice large difference in solution quality between the two solution methods. To improve the performance of LAPACK, Reverse Cuthill-Mckee (RCM) is applied to reduce the bandwidth of the matrices. Using several test cases, the LAPACK performance is shown to be similar to that of IMSL when using RCM. In that respect, the transition from IMSL… Advisors/Committee Members: van Gijzen, Martin (mentor), van der Zwan, S (graduation committee), Vuik, Kees (graduation committee), Heemink, Arnold (graduation committee), Delft University of Technology (degree granting institution).