Research from home, a practical guide

Things have changed, to say the least. As a MSCA Fellow the first years of my research project consisted of travelling through Europe with seminars, training courses and conferences spread over multiple countries. But then COVID-19 struck and we’ve been con fined to our homes for more then a year now. The funny thing about doing a PhD in the ROMSOC project is that most of us are quite independent and don’t have an urgent necessity to work in an office.

With a charged laptop and noise cancelling headphones, I can work wherever, whenever. This turned out to not be completely true though. There seems to be a reason that we never abandoned offices. At fi rst this newfound freedom felt great. I could go for a run whenever I felt like it and be back at work with no lost time. If I wanted to make my favourite slow cooked dish, rendang, I would just prepare it in the morning and stir whilst waiting for a simulation run. Technically, as mathematicians we were already quite far from following a dress-code, but I feel that writing a paper in your pyjamas at the faculty might be a bridge too far. But when you work, eat and relax in the same room for a months on end, the lines do start to blur. Therefore, some caution is needed to keep your spirits up and those stress levels down.

Going further into the lockdown I found my schedule drifting off . Waking up a little bit later each day, consequently working a bit longer. Then not working that bit longer and feeling under-productive. Although I had almost complete freedom in choosing when I would do whatever I wanted, I am still a mathematician. As mathematicians, we are always thinking in structures and how to optimise them. More and more I started to experiment with structuring my day and kept the parts that I liked in place. I ended up with rising at 7, work at 7:30, an hour of exercise and lunch at 12, then work till 6 and make sure to relax in the evening to go right back at it again the next day. In the end this is not that diff erent from a regular office schedule, it seems that I reinvented the wheel.

I hope you enjoyed this instalment of our ROMSOC blog. Please find the references below to see the results we have produced since my last blog post [1, 2] to follow the progress of me and my fellow researchers please keep an eye on the ROMSOC website and follow us on Twitter and Facebook!

Marcus Bannenberg

References
[1] MWFM Bannenberg, A Ciccazzo, and M G ̈unther. Reduced order multirateschemes for coupled differential-algebraic systems. 2021.2
[2] MWFM Bannenberg, F Kasolis, M G ̈unther, and M Clemens. Maximumentropy snapshot sampling for reduced basis modelling.preprint, 2020.

MSCA Fellow and the Lockdown

In order to fully understand what it means to be an MSCA Fellow in the Lockdown, it is important to understand what it meant to be an MSCA Fellow before. Having a Ph.D. position in an MSCA funded project is considered as some kind of “special opportunity” among the young scientists, but almost no-one who hasn’t experienced it can tell exactly why, and what makes it so special.

Coming from a pure science background and from an economically challenging region, I didn’t have much of a travel experience or inter-field connections before my Ph.D. program in European Industrial Doctorate. During an interview, my supervisor warned me, that I would have to travel all around Europe to attend the different courses and interact with professionals from different applied or industrial backgrounds.

My first reaction was a polite chuckle, since, how in the world, could such a thing require a warning when this sounded like an enormous advantage of the project?! Little did I know, that this would be the beginning of eighteen months of my life, which were simultaneously one of the most amazing, but also one of the most stressful.

For those first eighteen months of the project, I felt like an anchor-less ship, always in motion. Almost every month, I had to travel somewhere. Sometimes it was a one-week course or training, sometimes a conference or a workshop, sometimes it was a visit back home to friends and family, sometimes it was a complete relocation from one country to another, but there was always somewhere to go, some travel to plan.

And yet, I was not alone. At all those courses I was meeting ten other fellows who had the same experience, the same amount of travel, the same feeling of constant motion. The same experience of being in so many towns and hearing so many languages, that even remembering which word was “hello” in a particular place was becoming a challenge. This is what was making us so “special”, this is what was simultaneously our curse and our blessing.

The world became smaller for us. Everything in Europe was within arm’s length, just one step to take. We had the feeling of having everything within the reach. There was no place we couldn’t have gone to, no knowledge or professional opinion that we could need and not be able to obtain. Everything was accessible!

And then suddenly, comes the pandemic and the lockdown. If our favorite food is not sold in our closest supermarkets, it is not accessible anymore, if we are not in the same town as our family, we can’t see them for months. Even if we live next to our offices and universities, they are out of reach. Even if we live next to our favorite leisure places, their doors are shut tight. We might not be disconnected, but we are isolated!

So how does an MSCA fellow feel under the lockdown? I think like an explorer, whose ship is damaged in the storm and who is stuck in a remote village. Who knows that other explorers are still out there, fighting the storm, and who is waiting for them to come to aid. Just like we are waiting for our fellow scientists, doctors, and pharmacists, who had been on the front lines of the battle against pandemics, to defeat it, and bring the world within the reach once again!

About the Author

Giorgi Rukhaia is an early-stage researcher (ESR3) in the ROMSOC project. He is affiliated with the Institut national de recherche en informatique et en automatique (INRIA) in Paris, France and working in collaboration with Signify (former Philips Lighting), the Netherlands. In his research project he is working  on an Optimal Transportation computational approach and the development of numerical algorithms for the inverse free-form optical surfaces design using extended sources.

The curious case of modeling acoustics in fluid flow

The significant contributions and problems posed by the peers over time have all contributed to the development of scientific knowledge. Being curious about the history of scientific evolution can generally help develop a new organic perspective. It provides an invaluable lesson by leading us through the reasoning and choices along the path of scientific progress. The last couple blog posts ( viz., Origin and need for coupled models in acoustics – ROMSOC and Mathematical modelling acoustic measurements – ROMSOC ) have shed light into the motivation and approach followed in coupled systems modelling in acoustic industry, while this post intends to delve deep into the world of modeling acoustics in fluid flow. Describing fluid flow models while also coupled with acoustics is definitely harder in ‘simpler’ terms, I assume the reader is somewhat familiar to fluid flow models.

Almost always, fluid flow in mathematical models begins with the Navier-Stokes(NS) equations and involves some  reformulation/restructuring, simplification along with/without change in frame of reference. Lighthill and Curle (1950-60s) pioneered the development in this regard especially in the field of aeroacoustics providing models for jet noise prediction with a rearranged NS equation in the form of a linear wave equation for a fluid medium at rest. The model primarily assumes a equivalent-source analogy- suggesting that the source term may always be modeled approximately. Subsequent works by Ffowcs Williams, Hawkings, Phillips and Liley extended this approach by restructuring the NS equations into a more general form of convective wave equation. These models could help predict sound generation in simple geometries and flow profiles, but complex cases required by the industry were still out of reach. A multitude of empirical and semi-empirical approaches for specific applications followed, although failing to address a holistic analogy for flow acoustics.

Another favored modelling approach was to make some simplifications on the much harsher NS equation and trying to tame a simpler set of equations hoping that they cater to the complex cases. Euler equations for one, are simplified NS equations with the effects of viscosity and heat transfer omitted. Convenient linear approximations on density, pressure and velocity fields around the mean flow give the Linearized Euler Equations (ironically since the equations themselves are non-linear). These equations received interest especially in fluid-structure interaction problems since they could capture refractional effects and reflections at solid boundaries (especially relevant for aero-industry); however, they are often used for 2D geometries and behave well only for low amplitude acoustic perturbations. The model is still sensitive to small changes and there are further simplified LEE models e.g. Bogey-Bailly-Juve model which ignores some more terms to make the model simpler.

I would like to shed light into one of the less-known models – the Galbrun’s equation. This less-used model utilizes a different form of the NS equations within a mixed frame of reference. It assumes a known stationary underlying Eulerian flow and models the Lagrangian fluctuations of acoustics. The displacement formulation of it is given as,

Here,  and  are the density and speed of sound of the underlying flow with a pressure field  while,  is the material derivative. The model offers advantages in coupling problems, especially the compatibility of boundary conditions using Lagrangian terms, and offers an option to modularize the problems of fluid flow and the acoustic fluctuations. It also comes with a few setbacks – the equation is notoriously unstable with traditional finite elements and needs ‘regularization’ for effective use as was demonstrated by Bonnet-Ben Dhia and others [REFs]. Addressing these challenges could provide for a powerful tool to integrate acoustics in flowing medium, helping develop better noise predictions in fluid flow and the ROMSOC project hopes to address this for the real-world applications by providing the open-access knowledge and tools. Keep following this blog for more information into a diverse set of mathematical technologies that ROMSOC hopes to develop for industry.

 

References:

Howe, M. (1998). Acoustics of Fluid-Structure Interactions (Cambridge Monographs on Mechanics). Cambridge: Cambridge University Press. doi:10.1017/CBO9780511662898

Goldstein, M. E. (2003). A generalized acoustic analogy. Journal of Fluid Mechanics, 488, 315–333. https://doi.org/10.1017/S0022112003004890

Billson, M., Eriksson, L.-E., & Davidson, L. (2005). Acoustic Source Terms for the Linearized Euler Equations in Conservative Form. AIAA Journal, 43(4), 752–759. https://doi.org/10.2514/1.12858

Maeder, M. (2020). 90 Years of Galbrun’s Equation: An Unusual Formulation for Aeroacoustics and Hydroacoustics in Terms of the Lagrangian Displacement. 39.

 

 

About the author:

Ashwin Nayak is an early-stage researcher (ESR2) in the Reduced Order Modeling Simulation Optimization and Coupled methods (ROMSOC) project. He is affiliated with the Technological Institute of Industrial Mathematics (ITMATI), Spain and working in collaboration with Microflown Technologies B.V., Netherlands. His research aims to develop coupled mathematical models to enable acoustical measurements in fluid-porous media interactions.

What makes a solver a good one?

Solving large linear systems in real-time applications

In many practical examples solving large linear systems in real-time is the key issue. This is the case, e.g., for our ROMSOC project “Real Time Computing Methods for Adaptive Optics”, where we are dealing with a rapidly changing atmosphere of the earth. Real-time usually refers to a time frame within milliseconds. For our specific test configuration, it is 2 ms.

Before talking about efficient solvers for linear systems we need to define what is a linear system and what is large. A linear system is a system of linear equations Ax = b, where A is a matrix of size n x n, b is the right hand side vector of size n and x is the desired solution. In this context, large refers to an n up to several hundred thousand.  The larger the size of the matrix, the harder the system is to solve in real-time. Note, that in our setting we are dealing with a dense square matrix A.

Now that we have stated the mathematical problem formulation, we need to fix some performance criteria to compare the different ways to solve a linear system. What does it mean to have a good or efficient solver? There are plenty of ways to solve such systems, either directly or iteratively. However, there are some quantities that can be checked in advance no matter what specific solver you want to analyze in order to be able to choose a suitable one. In the end, efficient always refers somehow to the run-time of the algorithm. But the run-time heavily depends on the hardware and the implementation, which is usually not done in a preliminary state. Hence, it would be nice to have some criteria that can be checked in advance before implementing the whole method.

In mathematics, solvers are often compared using the computational complexity, which simply counts the number of operations needed to solve the problem. Usually, this quantity is given asymptotically in the order of n. This becomes clearer if we look at some examples: Adding two vectors of length n is of linear complexity O(n), because we need to perform n additions to obtain the result. A matrix vector multiplication is of quadratic complexity O(n2). However, this criterion does not take into account if a matrix is sparse or dense. Note that we say a matrix is sparse, if it has a considerable number of zero entries. A more precise measure in this case is the number of floating point operations (FLOPs). When performing a dense matrix vector multiplication, we require 2n2-n FLOPs to get the solution, however, a sparse matrix vector multiplication needs only 2nnz2-nnz FLOPS. Here, nnz refers to the non-zero entries of A. The above two quantities are suitable for a first, theoretical performance analysis. A next step is to look at parallelization possibilities of the method. We call a method parallelizable, if the whole algorithm or at least the main parts can be executed in parallel, i.e., without a dependency on intermediate results. The big advantage here is that the run-time is considerably reduced while the number of FLOPs stays the same. Most of the real-time applications require massively parallel algorithms in order to be able to solve the problem in the desired time frame. The efficiency of parallelization depends on the hardware in use. GPUs have an enormous computational throughput compared to CPUs and perform extremely well for highly parallelizable algorithms. Some applications can also benefit a lot from a so called matrix-free representation, i.e., the matrix A is formulated as a linear function rather than a sparse or dense matrix. Let us illustrate this on the example of the matrix A=vvT. Saving the matrix requires n2 units of memory and performing a matrix vector multiplication with a vector x needs 2n2-n FLOPs. However, if we save only the vector v instead of the whole matrix this requires just n units of storage and multiplying x first by vT and then by v needs only 3n FLOPs. In this example we used the term memory usage, which denotes the units of memory required to store a matrix or vector. This is another important criterion, since the memory of hardware is limited and storing matrices can become quite memory intensive.

Altogether, we listed here 5 criterions based on whom you can decide for a solver how well it is suited for your application. Nevertheless, in the end you have to implement and optimize the algorithm for your specific hardware, which can become a very crucial and time demanding task.

About the author

My name is Bernadett Stadler and I studied Industrial Mathematics at the Johannes Kepler University in Linz. Since May 2018 I am part of the ROMSOC program as a PhD student. I am involved in a cooperation of the university in Linz and the company Microgate in Bolzano. My research addresses the improvement of the image quality of extremely large telescopes. The goal is to enable a sharper view of more distant celestial objects.

Optimal shape design of air ducts

Shape optimization has proved its indispensability from both theoretical and application points of view in many real world applications [1, 2, 4], such as drag reduction of aircrafts, designs of pipelines, bridges, engine components, implants for blood vessels, etc. For many optimal design problems in engineering or biomedical sciences, an optimal shape of a component or a device is determined to minimize some suitable objectives subject to fluid flows. In the automotive industry particularly, the optimal shape design of several components of combustion engines, such as air ducts, piston, crankshaft, valves, etc., is crucial for their performance optimization. Additionally, the possible design is restricted by some geometric constraints. To describe this problem mathematically, we use the stationary Navier-Stokes equation (NSEs)

to model the air flow inside a tube illustrated by the following figure:

Here f and fin are the given source density and in flow profi le, u and p are the velocity and the kinematic pressure, n is the unit outward normal vector. The uniformity of the flow leaving the outlet is an important design criterion of automotive air ducts to enhance the eciency of distributing the air flow [2]. To achieve this criterion, we minimize a cost functional capturing the distance between the velocity (or only its normal component) and a given desired velocity in the outlet. Another criterion engineers also want to minimize is the power dissipated by air ducts (and any fluid dynamic devices in general) [2]. This dissipated power can be computed as the net inward flux of energy through the boundary of the considered tube. Taking into account both objectives, we consider a mixed cost functional as a convex combination of the cost functionals above. A typical PDE-constrained shape optimization problem [1, 4] can be established by finding an admissible shape to minimize the mixed cost functional under the given Navier-Stokes equation. After setting up our optimization problem, the optimal shape of the tube can be analyzed with some suitable methods based on shape derivatives. Furthermore, large Reynolds number flows are known to be unstable and computationally challenging. We will explore the world of turbulence models [3] to enhance our con guration for the shape optimization problem in the next blog.


References

  1. Mohammadi, Bijan; Pironneau, Olivier. Applied shape optimization for fluids. Second edition. Numerical Mathematics and Scienti c Computation. Oxford University Press, Oxford, 2010. xiv+277 pp. ISBN: 978-0-19-954690-9.
  2. Othmer, C. A continuous adjoint formulation for the computation of topological and surface sensitivities of ducted flows. Internat. J. Numer. Methods Fluids 58 (2008), no. 8, 861{877.
  3. Pope, Stephen B. Turbulent flows. Cambridge University Press, Cambridge, 2000. xxxiv+771 pp. ISBN: 0-521-59886-9.
  4. Sokolowski, Jan; Zolesio, Jean-Paul. Introduction to shape optimization. Shape sensitivity analysis. Springer Series in Computational Mathematics, 16. Springer-Verlag, Berlin, 1992. ii+250 pp. ISBN: 3-540-54177-2.

About the author

Hong Nguyen is an Early-Stage Researcher within the ROMSOC project and a Ph.D. student at the Weierstrass Institute for Applied Analysis and Stochastics (WIAS) in Berlin. He is working in collaboration with Math.Tec GmbH in Vienna (Austria) on optimal shape design of air ducts in combustion engines.

Artificial Neural Network and Data-driven techniques: Scientific Computing in the era of emerging technologies

by Nirav Vasant Shah

Rapidly changing technological and scientific environments pose the challenges, which require “adaptability” of our knowledge with the skills of the future. One of the most important emerging technology of the recent times is the Artificial Intelligence (AI). The Artificial Neural Network (ANN), foundation of AI, has entered into the field of scientific computing and is set to become new paradigm in modeling and computation. The concept of ANN is inspired from the massive, hierarchical neural network in the human brain. The salient feature of the ANN is its ability to mimic the analytical relationship between input and output. This is possible by adjusting the (hyper)parameters of the ANN during training procedure against known input and output pairs. The use of ANN is giving rise to new field called data-driven techniques. The aim of data-driven techniques is to allow the system to learn and to train from the data rather than to create and solve a system of equations.

Neural network schematic

The data-driven techniques have brought alternative approaches in the field of numerical analysis. Usually, classical numerical analysis methods require to solve a system of equations, based on an alternate form of governing physical equations. In the data-driven techniques, the system aims to compute the solution field as minimizer of a loss function. As an example, one of the interesting variant of data-driven techniques is Physics Informed Deep Learning, described in [3],[4]. In this approach, the ANN is used to minimize the residual based on the governing conservation equation. It was demonstrated that such ANN can be used either to compute the solution field or to learn the parameters of the system. The successful results proved that Physics Informed Deep Learning can become an alternative to classical methods.

In Model Order Reduction, another area of scientific computing, various data-driven approaches have been proposed as an alternative to the projection based approaches [2]. The biggest difference between projection based approaches and data-driven approaches is the non-intrusive nature of the latter. Non-intrusive methods are preferrable in the context of industrial applications. Besides, assembly of system of equations was eliminated adding to the significant speedup of the reduced order model.

Despite the promising results and important advantages, data-driven approaches, in their current form, cannot become fully reliable alternative to classical methods such as Finite Element Method. The classical methods have rigorous mathematical background and have been successfully applied to the problems of varying degree of complexity. Besides, as a recent SIAM news article [1] explained, the deep learning methods face some of the major challenges in the context of scientific computing. For example, generation of data and training of ANN can be very expensive or the integration of deep learning libraries with data from other open-source computing libraries may remain a challenging task. The trend that needs to be observed is whether the developments in the areas of mathematical theories, algorithmic efficiency and increasing hardware capabilities, will make the data-driven techniques robust and more acceptable within the scientific computing community. As an Early Stage Researcher, it is a matter of clear professional choice : whether to invest into a quickly emerging field with many known unknowns or to opt for classical approaches due to their reliability and availability of references. A choice that lays foundation for the skills of the future or a choice that matches with seemingly sustainable and reliable career path.

1] Aparna Chandramowlishwaran. Artificial intelligence and highperformance computing: The drivers of tomorrow’s science. SIAM news. https://sinews.siam.org/Details-Page/artificial-intelligence-andhigh-performance-computing-the-drivers-of-tomorrows-science, October 2020.

[2] Jan S. Hesthaven and Stefano Ubbiali. Non-intrusive reduced order modeling of nonlinear problems using neural networks. Journal of Computational Physics. http://www.sciencedirect.com/science/article/pii/S0021999118301190,363:55 – 78, 2018.

[3] Maziar Raissi, Paris Perdikaris, and George Em Karniadakis. Physics informed deep learning (part i): Data-driven solutions of nonlinear partial differential equations. arXiv preprint arXiv:1711.10561, 2017.

[4] Maziar Raissi, Paris Perdikaris, and George Em Karniadakis. Physics informed deep learning (part ii): Data-driven discovery of nonlinear partial differential equations arXiv preprint arXiv:1711.10566, 2017.

About the author

Nirav Vasant Shah is an Early-Stage Researcher (ESR10) within the ROMSOC project. He is a PhD student at the Scuola Internazionale Superiore di Studi Avanzati di Trieste (SISSA) in Trieste (Italy). He is working in collaboration with ArcelorMittal, the world’s leading steel and mining company, in Asturias (Spain) and Technological Institute for Industrial Mathematics (ITMATI) in Santiago de Compostela (Spain) on the mathematical modelling of thermo-mechanical phenomena arising in blast furnace hearth with application of model reduction techniques.

What’s the Extension in the Extended Finite Element Method?

After reading about Umberto’s amazing cycling experience across Italy in the last blog entry, we sadly had to admit that summer was really over. Then it is really time to sit again at your school desk and deal with a more technical topic. Therefore, this time I would like to introduce you to the basics of the Extended Finite Element Method or X-FEM. What does the extension in X-FEM stand for? Or, if you rather, what is the X factor of X-FEM?

Originally proposed by A. Hansbo and P. Hansbo [1], X-FEM is an advanced numerical technique that extends the range of applications of the more classical Finite Element Method (FEM), thanks to its geometric flexibility and the preservation of its accuracy even in complex topological configurations.
Before explaining the reasons for this gain in X-FEM, I first need to recall the main concepts in FEM theory. Used in most real-engineering problems, FEM is a numerical method used to solve Partial Differential Equations (PDEs) which model the dynamics of a physical system in a properly discretized domain. Indeed, the domain is typically divided into triangles (in 2D) or into tetrahedra (in 3D), called finite elements, which form all together the computational mesh. In the simplest case of linear FEM, the vertices of the finite elements are the points (or nodes) of the mesh where the solution is computed by solving a corresponding linear system. Then, such nodal solutions are combined with basis functions to represent the solution over the whole domain. Therefore, each mesh node corresponds to a degree of freedom (or dof) of the system and, in general, the more points are considered in the spatial discretization, the more accutate the results should be.

So how is the extension of X-FEM implemented in this framework? The key for the answer to this question resides in the enriched nature of the finite elements of X-FEM, which allows for solutions with internal discontinuities thanks to a proper treatment of the associated dofs.

But let‘s take another step back to better understand the context of my application of this numerical technique. In my Ph.D. research, I decided to adopt X-FEM to simulate the complex dynamics occurring in the so-called Wave Membrane Blood Pump, developed at CorWave SA (Clichy, France). Indeed, this novel pumping technology is based on the mutual interaction between the blood flow and an undulating polymer membrane, which results in an effective propulsion working against the adverse pressure gradient between the left ventricle and the aorta. Therefore, we are technically dealing with a 3D Fluid-Structure Interaction (FSI) problem with a thin immersed membrane that undergoes to large displacements. The X-FEM formulation for this challenging class of FSI problems was proposed in [2].

In FSI problems we need at least two meshes, one for the fluid and one for the structure. Since XFEM is an unfitted mesh technique, the two meshes lay on two different levels: in the background, there is the fluid mesh, that is fixed in time; while in the foreground, the structure mesh is free to move and intersect the underlying fluid mesh (see Figure 1). This means that some fluid mesh elements are cut by the structure mesh and perhaps split into multiple parts, where we need to represent the fluid solution. These fluid elements are called split elements. Notice that the split elements change in time because the structure mesh moves and intersects different fluid elements in different ways.

Illustration 1: Figure 1

Now, here’s the trick of X-FEM. The split elements are enriched – or extended – by duplicating the corresponding dofs. Hence, we can use a set of dofs (i, j, k) to integrate the solution in one fluid visible portion of the split element, and the second set of dofs (i’, j’, k’) to integrate the solution in the other portion (see Figure 2). In this way, we end up with a simple way to embed a discontinuous solution in the same (extended) finite element using standard basis functions for the integration in space. Otherwise, in FEM framework, we would need to switch to discontinuous modal basis functions directly defined on the polyhedra generated by the intersetion, making the problem very stiff in three dimensions.

Illsutration 2: Figure 2

However, I must add that even the implementation of X-FEM is not simple, because we need to compute the intersections between the fluid and the structure meshes at each timestep and duplicate the dofs accordingly. Indeed, to my knowledge, this is the first time that X-FEM is applied on a real 3D industrial problem. Nevertheless, we validated the X-FEM-based numerical model against experimental data, proving that it is a reliable predictive tool for the hydraulic performance of the pump device [3]. In Figure 3, I report a snapshot from an X-FEM simulation of the FSI in the wave membrane blood pump, where we visualize the blood velocity field and the displcament field of the elastic wave membrane.

Illustration 3: Figure 3

I hope I have “extended” a bit your knowledge. Sorry for the bad joke.

[1] A. Hansbo and P. Hansbo. An unfitted finite element method, based on Nitsche’s method, for elliptic interface problems. Computer methods in applied mechanics and engineering, 2002

[2] S. Zonca, C. Vergara, and L. Formaggia. An unfitted formulation for the interaction of an incompressible fluid with a thick structure via an XFEM/DG approach. SIAM Journal on Scientific Computing, 2018

[3] M. Martinolli, J. Biasetti, S. Zonca, L. Polverelli, C. Vergara, Extended Finite Element Method for Fluid- Structure Interaction in Wave Membrane Blood Pumps, submitted, MOX-Report 39/2020, 2020

Marco Martinolli

Bardonecchia – Otranto, stories of a ROMSOC fellow in cycling across Italy

August just ended and we are coming back to our research. But summer is still all around us and the mind easily go back to the last months with vivid memories of sunny days. So, I will take a break from any technicality and will tell you the story of my summer adventure.

It is June and we are finally allowed to get out. Being locked at home for months makes any relationship stiff and rusty. So one of the first thing I want to do is to meet again my dearest friend Paolo. Together we already cycled from Savona (north of Italy) to Lisbon years ago and we know each other so well that in few minutes the familiar feeling is back. Summer is approaching so we go back to the adventures we had together. Soon we get excited and start thinking about possible projects. We have only one week, not too bad but still too little time for great projects. Apparently. Some nice Alps trail, France, the Balkans. The ideas are many. But looking at the maps Paolo spots two points: Bardonecchia – Otranto and a straight line connecting them. The westernmost and the easternmost towns of Italian Peninsula. Already excited, we check the mileage. It is more than 1500km. We have seven days. Here comes the math: it is more than 200km per day. Non of us have ever made 200km in a day. I spot a twinkle in his eyes, the brain is already cycling. We fix the date.

One month later, the 25th of July, we are looking at the sunrise in Bardonecchia with fully packed bikes ready to cycle. From now on, our main thought is “keep pedalling”. The more we cycle, the souther we get. Every day, the accent changes and so the food, the attitude, the landscape. We go so fast through Italy that it is an overdose of sun and beauty.
We organized the trip to visit friends we do not meet since a long time. Then, when we arrive every evening, a friend is waiting for us. They immediately understand how exhausted we are after 10 hours of cycling. Every time, our host is full of care understanding our needs and thought. The feeling is to arrive at home every day.
As the days pass, the fatigue on our body increases, the days get hotter. Now we need to rest at midday because the thermometer goes over 40 degrees. In a single day, we drink over 20 water bottle trying to keep some water in the body. But we do not stop. Now we have crossed the Apennines, the number of kilometers per day increases but there are no more climbs to do. Orientation is not a problem anymore, we follow the sea line. The sixth day, we cycle in the sunshine the trabocchi coast. Every kilometer there is an ancient wood building going over the sea. They were built and used by farmers to fish because they did not know how to sail. I have never seen a clearer water and a more beautiful coast.

Reaching Barletta for the last stay before Otranto, we now we are going to success. We know that soon will be time to celebrate. The last day, we go so fast and happy that we arrive several hours earlier than expected. The town council is waiting for us and we arrive shouting and crying for the joy. The adventure is over. It was the adventure of a ROMSOC fellow.

About the author

Umberto Emil Morelli is an early stage researcher within ROMSOC. He is doing his PhD in reduced order modeling for boundary conditions estimation in heat transfer problems at the Technological Institute for Industrial Mathematics (ITMATI) in collaboration with Danieli & C. Officine Meccaniche SpA and the Scuola Internazionale Superiore di Studi Avanzati (SISSA).

Applications of mathematical optimization in railway planning

Regardless of the mode, transportation – as an industry –  faces numerous planning challenges, beyond these seen in other industries. They frequently entail thousands of interdependent decisions, with an element of uncertainty. This is why for many years, optimization techniques have found extensive use in that branch of economy. In this post, I would like to present some of the most prominent applications of optimization techniques in transportation, with a particular focus on rail freight transportation. I will also mention some works relating to each of the applications discussed. This choice is due to my own research focusing on that mode of transport.

Let us start from the broad perspective of the railway network of a country. With a history starting in the 19th century, these were frequently built by numerous companies, serving their own interests. Over time, they grew bigger and bigger, and gradually found themselves under the control of one entity, usually a state-owned company which maintains them and sells the access slots (just like at an airport). Now, as the transportation needs change, some of the lines are shut down while others are refurbished. On some occasions, completely new lines are built.

To assist decision-making while considering an extension to a railway network, policy-makers can use tools like those described in Andreas Bärmann’s [1]. Andreas uses forecasts of the demand for access to the railway network at different points in Germany and comes up with a model which highlights the sections of the network which will experience the greatest increase in demand. To make the problem solve efficiently, he inter alia develops a decomposition technique, made to the measure of the studied multi-period network design problem.

Going further, the network manager needs to distribute the access rights to the network. This is required in order to – in the most general terms – maintain safety in the network and assure that as many trains as possible can access the railway line. In particular, a certain distance – called headway – needs to be kept between two consecutive trains. Further, especially on the electrified sections of the network, an appropriate train schedule can help save energy thanks to recuperation.

One of the most prominent papers dealing with the described problem was written by Alberto Caprara, Matteo Fischetti and Paolo Toth [2]. They develop a graph theoretic formulation for the problem of scheduling trains on a single, one-way track linking two major stations, with a number of intermediate stations in between using a directed multigraph. They then use Lagrangian relaxation to arrive at an efficient heuristic solution algorithm, and demonstrate the results on instances supplied by two Italian railway network managers. More recently, my colleagues from the chair, Prof. Alexander Martin, Dr. Andreas Bärmann (who are the supervisors of my thesis) and others [3] together with VAG (transportation authority of Nuremberg metropolitan area) developed an approach to optimal timetabling aiming at minimizing the energy consumption of the trains through more energy-efficient driving as well as by increasing the usability of recuperated energy from braking. This solution could soon be used by VAG in their timetabling process. This work is part of a project on energy-saving via real-time driver assistant systems for interconnected trains in the ADA-Lovelace-Center in Nuremberg, Germany [4]. It is a research centre for artificial itelligence (AI) founded jointly by the Fraunhofer Institute for Integrated Circuits (IIS), the Friedrich-Alexander University Erlangen-Nürnberg (FAU) and the Ludwig-Maximilian University München (LMU) [5]. Its mission is to push methodological advancements in AI and to act as a platform for AI collaboration between research and industry.

Another important application of mathematical optimization pertains to the decisions relating to the production resources of a railway carrier, most important of which are locomotives and drivers. On one hand, one needs to ensure that both an appropriate locomotive (e.g. suitable to the line, powerful enough) and an appropriate driver (e.g. licensed to drive the mentioned locomotive on the train’s route) is attributed to the train, on the other hand, each driver has a number of working time constraints, and locos require maintenance. Other assets considered frequently include wagons and passenger train personnel.

A prominent paper pertaining to locomotive scheduling is the one by Ahuja et al. [6]. They study a problem of assigning locomotive consists to pre-planned trains, aiming at supplying each train with appropriate power. Another one, by Jütte et al. [7], uses column-generation techniques to develop and implement a crew-scheduling system at DB Schenker, the largest European rail freight carrier. My own research focuses on combining both problems (posed differently) into one, and solving them “from one hand”.

All of the abovementioned applications are just the “tip of the iceberg”. In this post, I only covered three macro-stages of planning in the railway industry, and showed how can mathematical optimization be of help there. Other applications could include e.g. train control, locomotive maintenance scheduling, shunting, train dispatch and delay management. Vast majority of those, as well as others, were all described in [8], which should serve as an reference point to the interested reader.

[1] ] Bärmann, A., 2015. Solving Network Design Problems via Decomposition, Aggregation and Approximation. Springer Spektrum. https://doi.org/10.1007/978-3-658-13913-1

[2] Caprara, A., Fischetti, M., Toth, P., 2002. Modeling and Solving the Train Timetabling Problem. Operations Research 50, 851–861. https://doi.org/10.1287/opre.50.5.851.362

[3] Bärmann, A., Gemander, P., Martin, A., Merkert, M., Nöth, F., 2019. Energy-Efficient Timetabling in a German Underground System. Preprint. http://www.optimization-online.org/DB_FILE/2020/04/7728.pdf

[4] https://www.scs.fraunhofer.de/en/focus-projects/ada-center/ai-transport-mobility.html

[5] https://www.scs.fraunhofer.de/en/focus-projects/ada-center.html

[6] Ahuja, R.K., Liu, J., Orlin, J.B., Sharma, D., Shughart, L.A., 2005. Solving Real-Life Locomotive-Scheduling Problems. Transportation Science 39, 503–517. https://doi.org/10.1287/trsc.1050.0115

[7] Jütte, S., Albers, M., Thonemann, U.W., Haase, K., 2011. Optimizing Railway Crew Scheduling at DB Schenker. Interfaces 41, 109–122. https://doi.org/10.1287/inte.1100.0549

[8] Borndörfer, R., Klug, T., Lamorgese, L., Mannino, C., Reuther, M., Schlechte, T. (Eds.), 2018. Handbook of Optimization in the Railway Industry, International Series in Operations Research & Management Science. Springer International Publishing, Cham. https://doi.org/10.1007/978-3-319-72153-8

About the author:

Jonasz Staszek is an early-stage researcher (ESR7) within ROMSOC. He is working with Friedrich-Alexander-Universität Erlangen-Nürnberg in collaboration with DB Cargo Polska on mixed-integer optimization models for an optimal resource allocation and investment for border-crossing transport services in rail traffic.

Model hierarchy for the reduced order modelling

It is essential to be aware of the financial risk associated with the invested product. In my last blog, I explained that one has to perform costly numerical simulations to understand the risks associated with a financial product.
The financial instruments are valuated via the dynamics of short-rate models, based on the convection-diffusion-reaction partial differential equations (PDEs). The choice of the short-rate model depends on the underlying financial instrument. Some of the prominent financial models are the one-factor Hull-White model, the shifted Black Karasinski model, and the two-factor Hull-White model. These models are calibrated based on several thousand simulated yield curves that generate a high dimensional parameter space. In short, to perform the risk analysis, the financial model is needed to be solved for such a high dimensional parameter space, which requires efficient algorithms. My Ph. D. aims to develop an approach that will perform such a computationally costly task as fast as possible but with a reliable outcome. Thus, I am developing a model order reduction approach.
Now, let us dive into the main topic. The model hierarchy simplifies the process of obtaining a reduced-order model (ROM) and is shown in the Fig. 1.

Figure 1: Model hierarchy to obtain a reduced order model.

It starts by discretizing the partial differential equation using either the finite difference method (FDM) or the finite element method (FEM). The discretized model is known as the full order model (FOM). To obtain a reduced-order model, one has to compute a reduced-order basis (ROB). To do so, the proper orthogonal decomposition (POD) approach relies on the method of snapshots. The snapshots are obtained by solving the FOM for some training parameters, which are further combined into a single matrix know as the snapshot matrix. Subsequently, the ROB is obtained by computing the truncated singular value decomposition (SVD) of the snapshot matrix. Finally, the FOM is projected onto the ROB to get the reduced-order model. One can easily infer that the quality of the ROB is based on the selection of training parameters. In this work, the training parameters are chosen based on either the greedy sampling or the adaptive greedy sampling approaches.

The greedy sampling technique selects the parameters at which the error between the ROM and the FOM is maximum. Further, the ROB is obtained using these selected parameters. The calculation of the relative error is expensive, so instead, the residual error associated with the ROM is used as an error estimator. However, it is not reasonable to compute an error estimator for the entire parameter space. This problem forces us to select a pre-defined parameter set as a subset of the high dimensional parameter space to train the greedy sampling algorithm. We usually select this pre-defined subset randomly. But, a random selection may neglect the crucial parameters within the parameter space. Thus, to surmount this problem, I implemented an adaptive greedy sampling approach. The algorithm chooses the most suitable parameters adaptively using an optimized search based surrogate modeling. This approach evades the cost of computing the error estimator for each parameter within the parameter space and instead uses a surrogate model to locate the training parameter set. I have built the surrogate model using two approaches: (i) the principal component regression model (ii) machine learning model. See our pre-print available on Arxiv for more details [1].

Figure 2: Greedy sampling (left) and adaptive greedy sampling (right) algorithms

References:
[1] A. Binder, O. Jadhav, and V. Mehrmann. Model order reduction for parametric high dimensional models in the analysis of financial risk. Technical report, 2020. https://arXiv:2002.11976.

About the author

Onkar Sandip Jadhav is an early-stage researcher (ESR6) in the ROMSOC project. He is a PhD student at the Technische Universität Berlin (Germany) and is working in collaboration with MathConsult GmbH (Austria). In his research he is working on a parametric model order reduction (MOR) approach aiming to develop a new MOR methodology for high-dimensional convection-diffusion reaction PDEs arising in computational finance with the goal to reduce the computational complexity in the analysis of financial instruments.