publications-journals

2023

I. Kicic, P. Vlachas, G. Arampatzis, M. Chatzimanolakis, L. Guibas, and P. Koumoutsakos, “Adaptive learning of effective dynamics for online modeling of complex systems,” Computer Methods in Applied Mechanics and Engineering, vol. 415, pp. 116204, 2023. Publisher's VersionAbstract

Predictive simulations are essential for applications ranging from weather forecasting to material design. The veracity of these simulations hinges on their capacity to capture the effective system dynamics. Massively parallel simulations predict the systems dynamics by resolving all spatiotemporal scales, often at a cost that prevents experimentation. On the other hand, reduced order models are fast but often limited by the linearization of the system dynamics and the adopted heuristic closures. We propose a novel systematic framework that bridges large scale simulations and reduced order models to extract and forecast adaptively the effective dynamics (AdaLED) of multiscale systems. AdaLED employs an autoencoder to identify reduced-order representations of the system dynamics and an ensemble of probabilistic recurrent neural networks (RNNs) as the latent time-stepper. The framework alternates between the computational solver and the surrogate, accelerating learned dynamics while leaving yet-to-be-learned dynamics regimes to the original solver. AdaLED continuously adapts the surrogate to the new dynamics through online training. The transitions between the surrogate and the computational solver are determined by monitoring the prediction accuracy and uncertainty of the surrogate. The effectiveness of AdaLED is demonstrated on three different systems - a Van der Pol oscillator, a 2D reaction–diffusion equation, and a 2D Navier–Stokes flow past a cylinder for varying Reynolds numbers (400 up to 1200), showcasing its ability to learn effective dynamics online, detect unseen dynamics regimes, and provide net speed-ups. To the best of our knowledge, AdaLED is the first framework that couples a surrogate model with a computational solver to achieve online adaptive learning of effective dynamics. It constitutes a potent tool for applications requiring many computationally expensive simulations.

E. Papadopoulou, G. W. Kim, P. Koumoutsakos, and G. Kim, “Molecular dynamics analysis of water flow through a multiply connected carbon nanotube channel,” Current Applied Physics, vol. 45, pp. 64-71, 2023. Publisher's VersionAbstract
The filling process of nanoconduits is an active research topic. In this study, we use molecular dynamics simulations to identify the filling process of water molecules in a multiply connected carbon nanotube (MCCNT). For water permeation, a local change in the channel cross-section affects the water filling of MCCNTs because it may lead to irregularities in the permeation profile. A decrease in hydrogen bonds at the junctions of the structure characterizes the permeability of MCCNTs. In contrast to pristine CNTs, the complex nanochannel exhibits a different imbibition profile due to the energy changes at the junction. Next, we examine the local water density and velocity patterns in MCCNT channels to understand how junction regions affect steady-state water transport. We find that there is congestion and irregularities in steady water flow density and velocity profiles. Through this study, we expect to develop effective channels with more complex geometries for water purification and drug delivery
P. Karnakov, S. Litvinov, and P. Koumoutsakos, “Flow reconstruction by multiresolution optimization of a discrete loss with automatic differentiation,” The European Physical Journal E, vol. 46, no. 59, 2023. Publisher's VersionAbstract
We present a potent computational method for the solution of inverse problems in fluid mechanics. We consider inverse problems formulated in terms of a deterministic loss function that can accommodate data and regularization terms. We introduce a multigrid decomposition technique that accelerates the convergence of gradient-based methods for optimization problems with parameters on a grid. We incorpo- rate this multigrid technique to the Optimizing a DIscrete Loss (ODIL) framework. The multiresolution ODIL (mODIL) accelerates by an order of magnitude the original formalism and improves the avoidance of local minima. Moreover, mODIL accommodates the use of automatic differentiation for calculating the gradients of the loss function, thus facilitating the implementation of the framework. We demonstrate the capabilities of mODIL on a variety of inverse and flow reconstruction problems: solution reconstruction for the Burgers equation, inferring conductivity from temperature measurements, and inferring the body shape from wake velocity measurements in three dimensions. We also provide a comparative study with the related, popular Physics-Informed Neural Networks (PINNs) method. We demonstrate that mODIL has three to five orders of magnitude lower computational cost than PINNs in benchmark problems including simple PDEs and lid-driven cavity problems. Our results suggest that mODIL is a very potent, fast and consistent method for solving inverse problems in fluid mechanics.
L. Amoudruz, A. Economides, G. Arampatzis, and P. Koumoutsakos, “The stress-free state of human erythrocytes: Data- driven inference of a transferable RBC model,” Biophysical Journal, 2023. Publisher's VersionAbstract

The stress-free state (SFS) of red blood cells (RBCs) is a fundamental reference configuration for the calibration of computational models, yet it remains unknown. Current experimental methods cannot measure the SFS of cells without affecting their mechanical properties, whereas computational postulates are the subject of controversial discussions. Here, we introduce data-driven estimates of the SFS shape and the visco-elastic properties of RBCs. We employ data from single- cell experiments that include measurements of the equilibrium shape of stretched cells and relaxation times of initially stretched RBCs. A hierarchical Bayesian model accounts for these experimental and data heterogeneities. We quantify, for the first time, the SFS of RBCs and use it to introduce a transferable RBC (t-RBC) model. The effectiveness of the proposed model is shown on predictions of unseen experimental conditions during the inference, including the critical stress of transitions between tum- bling and tank-treading cells in shear flow. Our findings demonstrate that the proposed t-RBC model provides predictions of blood flows with unprecedented accuracy and quantified uncertainties.

2022

M. Chatzimanolakis, P. Weber, and P. Koumoutsakos, “Vortex separation cascades in simulations of the planar flow past an impulsively started cylinder up to Re = 100 000,” Journal of Fluid Mechanics, vol. 953, 2022. Publisher's VersionAbstract
Direct numerical simulations of the flow past an impulsively started cylinder at high Reynolds numbers (25k–100k) reveal an intriguing portrait of unsteady separation. Vorticity generation and vortex shedding entails a cascade of separation events on the cylinder surface that are reminiscent of Kelvin–Helmholtz instabilities. Primary vortices roll up along the cylinder surface as a result of instabilities of the initially attached vortex sheets, followed by vortex eruptions, creation of secondary vorticity and formation of dipole structures that are subsequently ejected from the surface of the cylinder. We analyse the vortical structures and their relationship to the forces experienced by the cylinder. This striking cascade of vortex instabilities may serve as reference for reduced-order models of flow separation and as guide for flow control of separated flows at high Reynolds numbers.
S. Cantero-Chinchilla, et al., “Robust optimal sensor configuration using the value of information,” Structural Control and Health Monitoring, 2022. Publisher's VersionAbstract

Sensing is the cornerstone of any functional structural health monitoring technology, with sensor number and placement being a key aspect for reliable monitoring. We introduce for the first time a robust methodology for optimal sensor configuration based on the value of information that accounts for (1) uncertainties from updatable and nonupdatable parameters, (2) variability of the objective function with respect to nonupdatable parameters, and (3) the spatial correlation between sensors. The optimal sensor configuration is obtained by maximizing the expected value of information, which leads to a cost-benefit analysis that entails model parameter uncertainties. The proposed methodology is demonstrated on an application of structural health monitoring in plate-like structures using ultrasonic guided waves. We show that accounting for uncertainties is critical for an accurate diagnosis of damage. Furthermore, we provide critical assessment of the role of both the effect of modeling and measurement uncertainties and the optimization algorithm on the resulting sensor placement. The results on the health monitoring of an alu- minum plate indicate the effectiveness and efficiency of the proposed methodology in discovering optimal sensor configurations.

P. R. Vlachas, G. Arampatzis, C. Uhler, and P. Koumoutsakos, “Multiscale simulations of complex systems by learning their effective dynamics,” Nat Mach Intell, 2022. Publisher's VersionAbstract
Predictive simulations of complex systems are essential for applications ranging from weather forecasting to drug design. The veracity of these predictions hinges on their capacity to capture effective system dynamics. Massively parallel simulations predict the system dynamics by resolving all spatiotemporal scales, often at a cost that prevents experimentation, while their findings may not allow for generalization. On the other hand, reduced-order models are fast but limited by the frequently adopted linearization of the system dynamics and the utilization of heuristic closures. Here we present a novel systematic framework that bridges large-scale simulations and reduced-order models to learn the effective dynamics of diverse, complex systems. The framework forms algorithmic alloys between nonlinear machine learning algorithms and the equation-free approach for modelling complex systems. Learning the effective dynamics deploys autoencoders to formulate a mapping between fine- and coarse-grained representations and evolves the latent space dynamics using recurrent neural networks. The algorithm is validated on benchmark problems, and we find that it outperforms state-of-the-art reduced-order models in terms of predictability, and large-scale simulations in terms of cost. Learning the effective dynamics is applicable to systems ranging from chemistry to fluid mechanics and reduces the computational effort by up to two orders of magnitude while maintaining the prediction accuracy of the full system dynamics. We argue that learning the effective dynamics provides a potent novel modality for accurately predicting complex systems.
J. Lipkova, B. Menze, B. Wiestler, P. Koumoutsakos, and J. Lowengrub, “Modelling glioma progression, mass effect and intracranial pressure in patient anatomy,” J. R. Soc. Interface, vol. 19, no. 188, 2022. Publisher's VersionAbstract
Increased intracranial pressure is the source of most critical symptoms in patients with glioma, and often the main cause of death. Clinical interventions could benefit from non-invasive estimates of the pressure distribution in the patient's parenchyma provided by computational models. However, existing glioma models do not simulate the pressure distribution and they rely on a large number of model parameters, which complicates their calibration from available patient data. Here we present a novel model for glioma growth, pressure distribution and corresponding brain deformation. The distinct feature of our approach is that the pressure is directly derived from tumour dynamics and patient-specific anatomy, providing non-invasive insights into the patient's state. The model predictions allow estimation of critical conditions such as intracranial hypertension, brain midline shift or neurological and cognitive impairments. A diffuse-domain formalism is employed to allow for efficient numerical implementation of the model in the patient-specific brain anatomy. The model is tested on synthetic and clinical cases. To facilitate clinical deployment, a high-performance computing implementation of the model has been publicly released.
J. Bae and P. KOUMOUTSAKOS, “Scientific multi-agent reinforcement learning for wall-models of turbulent flows,” Nature Communications, vol. 13, 2022. Publisher's VersionAbstract

The predictive capabilities of turbulent flow simulations, critical for aerodynamic design and weather prediction, hinge on the choice of turbulence models. The abundance of data from experiments and simulations and the advent of machine learning have provided a boost to turbulence modeling efforts. However, simulations of turbulent flows remain hindered by the inability of heuristics and supervised learning to model the near-wall dynamics. We address this challenge by introducing scientific multi-agent reinforcement learning (SciMARL) for the discovery of wall models for large-eddy simulations (LES). In SciMARL, discretization points act also as cooperating agents that learn to supply the LES closure model. The agents self- learn using limited data and generalize to extreme Reynolds numbers and previously unseen geometries. The present simulations reduce by several orders of magnitude the computa- tional cost over fully-resolved simulations while reproducing key flow quantities. We believe that SciMARL creates unprecedented capabilities for the simulation of turbulent flows.

P. R. Vlachas, J. Zavadlav, M. Praprotnik, and P. Koumoutsakos, “Accelerated Simulations of Molecular Systems through Learning of Effective Dynamics,” J. Chem. Theory Comput. vol. 18, no. 1, pp. 538-549, 2022. Publisher's VersionAbstract
Simulations are vital for understanding and predicting the evolution of complex molecular systems. However, despite advances in algorithms and special purpose hardware, accessing the time scales necessary to capture the structural evolution of biomolecules remains a daunting task. In this work, we present a novel framework to advance simulation time scales by up to 3 orders of magnitude by learning the effective dynamics (LED) of molecular systems. LED augments the equation-free methodology by employing a probabilistic mapping between coarse and fine scales using mixture density network (MDN) autoencoders and evolves the non-Markovian latent dynamics using long short-term memory MDNs. We demonstrate the effectiveness of LED in the Müller–Brown potential, the Trp cage protein, and the alanine dipeptide. LED identifies explainable reduced-order representations, i.e., collective variables, and can generate, at any instant, all-atom molecular trajectories consistent with the collective variables. We believe that the proposed framework provides a dramatic increase to simulation capabilities and opens new horizons for the effective modeling of complex molecular systems.

2021

K. Daniilidis, et al., “Robotics in the AI era: a vision for a hellenic robotics initiative,” Found. Trends Mach. Learn. vol. 9, no. 3, pp. 201-265, 2021. Publisher's VersionAbstract

In January 2021, the Hellenic Institute of Advanced Study (HIAS) assembled a panel including world leading roboticists from the Hellenic diaspora, who volunteered their scientific expertise to provide a vision for Robotics in Greece. This monograph, entitled “Robotics in the Artificial Intelligence (AI) era,” will hopefully trigger a dialogue towards the development of a national robotics strategy. Our vision is that Robotics in the AI era will be an essential technology of the future for the safety and security of the Hellenic nation, its environment and its citizens, for modernizing its economy towards Industry 4.0, and for inspiring and educating the next generation workforce for the challenges of the 21st century. To contribute towards making this vision a reality, after reviewing global trends in robotics and assessing the Greek robotics ecosystem, we arrived at the following key findings and recommendations: Firstly, we think that Greece should develop a national Hellenic Robotics Initiative that serves as the nation’s long-term vision and strategy across the entire Greek robotics ecosystem. Also, certain societal drivers should be key in the areas of focus. Safety and security is an area of national importance necessitating a national initiative, while agrifood, maritime and logistics provide opportunities for internationally leading innovation.

We recommend the establishment of a mission-driven, government-funded, organization advancing unmanned vehicles in societal drivers of national importance, and Greece should leverage its unique geography and become a living testbed of robotics innovation turning the country into a development site for exportable technologies. In our opinion, universities should create Centers of Excellence in robotics and AI as well as consider innovation-leading research institutes such as the Italian Institute of Technology. We recommend investing in robotics education using Maker Spaces in order to prepare the workforce with 21st century skills to become Industry 4.0 innovators. Furthermore, we believe that the government should collect, measure, and analyze data on the robotics industry, robotics uses, labor shifts, and brain gain and promote awareness via a Hellenic Robotics Day. And finally, the government should regulate robot safety without stifling innovation, provide safe experimentation areas and mechanisms for certifying safety of locally developed robots. This monograph has many additional suggestions that enhance the above main recommendations. As authors, we advocate bringing the robotics ecosystem together in order to sharpen and expand these findings to an ambitious, long-term, and detailed national strategy and roadmap for robotics in the AI era.

L. Amoudruz and P. Koumoutsakos, “Independent Control and Path Planning of Microswimmers with a Uniform Magnetic Field,” Advanced Intelligent Systems, vol. 2100183, 2021. Publisher's VersionAbstract
Artificial bacteria flagella (ABFs) are magnetic helical microswimmers that can be remotely controlled via a uniform, rotating magnetic field. Previous studies have used the heterogeneous response of microswimmers to external magnetic fields for achieving independent control. Herein, analytical and reinforcement learning control strategies for path planning to a target by multiple swimmers using a uniform magnetic field are introduced. The comparison of the two algorithms shows the superiority of reinforcement learning in achieving minimal travel time to a target. The results demonstrate, for the first time, the effective independent navigation of realistic microswimmers with a uniform magnetic field in a viscous flow field.
P. Gunnarson, I. Mandralis, G. Novati, P. Koumoutsakos, and J. Dabiri, “Learning efficient navigation in vortical flow fields,” Nature Communications, vol. 12, no. 1, 2021. Publisher's VersionAbstract

Efficient point-to-point navigation in the presence of a background flow field is important for robotic applications such as ocean surveying. In such applications, robots may only have knowledge of their immediate surroundings or be faced with time-varying currents, which limits the use of optimal control techniques. Here, we apply a recently introduced Reinforcement Learning algorithm to discover time-efficient navigation policies to steer a fixed- speed swimmer through unsteady two-dimensional flow fields. The algorithm entails input- ting environmental cues into a deep neural network that determines the swimmer’s actions, and deploying Remember and Forget Experience Replay. We find that the resulting swimmers successfully exploit the background flow to reach the target, but that this success depends on the sensed environmental cue. Surprisingly, a velocity sensing approach significantly out- performed a bio-mimetic vorticity sensing approach, and achieved a near 100% success rate in reaching the target locations while approaching the time-efficiency of optimal navigation trajectories. 

 

S. M. Martin, D. Wälchli, G. Arampatzis, A. E. Economides, P. Karnakov, and P. Koumoutsakos, “Efficient and scalable software framework for Bayesian uncertainty quantification and stochastic optimization,” Comput. Method. Appl. M. 2021. Publisher's VersionAbstract

We present Korali, an open-source framework for large-scale Bayesian uncertainty quantification and stochastic optimization. The framework relies on non-intrusive sampling of complex multiphysics models and enables their exploitation for optimization and decision-making. In addition, its distributed sampling engine makes efficient use of massively-parallel architectures while introducing novel fault tolerance and load balancing mechanisms. We demonstrate these features by interfacing Korali with existing high-performance software such as APHROS, LAMMPS (CPU-based), and MIRHEO (GPU-based) and show efficient scaling for up to 512 nodes of the CSCS Piz Daint supercomputer. Finally, we present benchmarks demonstrating that Korali outperforms related state-of-the-art software frameworks.

E. Papadopoulou, J. Zavadlav, R. Podgornik, M. Praprotnik, and P. Koumoutsakos, “Tuning the Dielectric Response of Water in Nanoconfinement through Surface Wettability,” ACS Nano, 2021. Publisher's VersionAbstract

The tunable polarity of water can be exploited in emerging technologies including catalysis, gas storage, and green chemistry. Recent experimental and theoretical studies have shown that water can be rendered into an effectively apolar solvent under nanoconfinement. We furthermore demonstrate, through molecular simulations, that the static dielectric constant of water can be modified by changing the wettability of the confining material. We find the out-of-plane dielectric response to be highly sensitive to the level of confinement and can be reduced up to 40× , in accordance with experimental data. By altering the surface wettability from superhydrophilic to super- hydrophobic, we observe a 36% increase for the out-of-plane and a 31% decrease for the in-plane dielectric constants. Our findings demonstrate the feasibility of tunable water polarity, a phenomenon with great potential for scientific and technological impact.

E. Papadopoulou, C. M. Megaridis, J. H. Walther, and P. Koumoutsakos, “Nanopumps without Pressure Gradients: Ultrafast Transport of Water in Patterned Nanotubes,” J. Phys, Chem. B, 2021. Publisher's VersionAbstract
The extreme liquid transport properties of carbon nanotubes present new opportunities for surpassing conventional technologies in water filtration and purification. We demonstrate that carbon nanotubes with wettability surface patterns act as nanopumps for the ultrafast transport of picoliter water droplets without requiring externally imposed pressure gradients. Large-scale molecular dynamics simulations evidence unprecedented speeds and accelerations on the order of 1010 g of droplet propulsion caused by interfacial energy gradients. This phenomenon is persistent for nanotubes of varying sizes, stepwise pattern configurations, and initial conditions. We present a scaling law for water transport as a function of wettability gradients through simple models for the droplet dynamic contact angle and friction coefficient. Our results show that patterned nanotubes are energy-efficient nanopumps offering a realistic path toward ultrafast water nanofiltration and precision drug delivery.
I. Mandralis, P. Weber, G. Novati, and P. Koumoutsakos, “Learning swimming escape patterns for larval fish under energy constraints,” Phys. Rev. Fluids, vol. 6, pp. 093101, 2021. Publisher's VersionAbstract
Swimming organisms can escape their predators by creating and harnessing unsteady flow fields through their body motions. Stochastic optimization and flow simulations have identified escape patterns that are consistent with those observed in natural larval swimmers. However, these patterns have been limited by the specification of a particular cost function and depend on a prescribed functional form of the body motion. Here, we deploy reinforcement learning to discover swimmer escape patterns for larval fish under energy constraints. The identified patterns include the C-start mechanism, in addition to more energetically efficient escapes. We find that maximizing distance with limited energy requires swimming via short bursts of accelerating motion interlinked with phases of gliding. The present, data efficient, reinforcement learning algorithm results in an array of patterns that reveal practical flow optimization principles for efficient swimming and the methodology can be transferred to the control of aquatic robotic devices operating under energy constraints.
K. Larson, et al., “Data-driven prediction and origin identification of epidemics in population networks,” Royal Society Open Science, 2021. Publisher's VersionAbstract
Effective intervention strategies for epidemics rely on the identification of their origin and on the robustness of the predictions made by network disease models. We introduce a Bayesian uncertainty quantification framework to infer model parameters for a disease spreading on a network of communities from limited, noisy observations; the state-of-the-art computational framework compensates for the model complexity by exploiting massively parallel computing architectures. Using noisy, synthetic data, we show the potential of the approach to perform robust model fitting and additionally demonstrate that we can effectively identify the disease origin via Bayesian model selection. As disease-related data are increasingly available, the proposed framework has broad practical relevance for the prediction and management of epidemics.
A. Khosronejad, S. Kang, F. Wermelinger, P. Koumoutsakos, and F. Sotiropoulos, “A computational study of expiratory particle transport and vortex dynamics during breathing with and without face masks,” Physics of Fluids, vol. 33, no. 6, pp. 066605, 2021. Publisher's VersionAbstract

We present high-fidelity numerical simulations of expiratory biosol transport during normal breathing under indoor, stagnant air conditions with and without a facile mask. We investigate mask efficacy to suppress the spread of saliva particles that is underpinnings existing social distancing recommendations. The present simulations incorporate the effect of human anatomy and consider a spectrum of saliva particulate sizes that range from 0.1 to 10 μm while also accounting for their evaporation. The simulations elucidate the vorticity dynamics of human breathing and show that without a facile mask, saliva particulates could travel over 2.2 m away from the person. However, a non-medical grade face mask can drastically reduce saliva particulate propagation to 0.72 m away from the person. This study provides new quantitative evidence that facile masks can successfully suppress the spreading of saliva particulates due to normal breathing in indoor environments.

A. Economides, et al., “Hierarchical Bayesian Uncertainty Quantification for a Model of the Red Blood Cell,” American Physical Society (APS), vol. 15, no. 13, 2021. Publisher's VersionAbstract
Simulations of blood flows in microfluidic devices and physiological systems are gaining importance in complementing experimental and clinical studies. The predictive capabilities of these simulations hinge on the parameters of the red blood cell (RBC) model that are usually calibrated from experimental data. However, these parameter values may vary drastically when calibrated using different experimental quan- tities or experimental settings. In turn, the results of existing blood flow simulations largely depend on the utilized parameters that have been chosen to validate a particular experiment. We suggest a revision to this type of model calibration to properly integrate experimental data in the computational models and accordingly inform their predictions. In this context, we introduce the calibration of a popular RBC model using data-driven, hierarchical Bayesian inference. We employ data from classical experiments of RBC stretching by optical tweezers and tank treading in shear flows, and distinguish the calibration of the model parameters through single-level and hierarchical Bayesian uncertainty quantification. We find that the optimal model parameters depend not only on the data used for the inference but also on the way the data are used in the inference process. Single-level Bayesian models predict well the data used in their calibration, but are inferior to the hierarchical Bayesian model at predicting previously unseen data. This work demonstrates that the proper integration of experimental data is essential for the development of a robust and transferable RBC model. We believe that the present study can serve as a prototype across scientific fields, in revising the integration of computational models and heterogeneous experimental data.
G. Novati, H. L. de Laroussilhe, and P. Koumoutsakos, “Automating turbulence modelling by multi-agent reinforcement learning,” Nature Machine Intelligence, vol. 3, pp. 87–96, 2021.Abstract
Turbulent flow models are critical for applications such as aircraft design, weather forecasting and climate prediction. Existing
models are largely based on physical insight and engineering intuition. More recently, machine learning has been contributing
to this endeavour with promising results. However, all efforts have focused on supervised learning, which is difficult to generalize beyond training data. Here we introduce multi-agent reinforcement learning as an automated discovery tool of turbulence
models. We demonstrate the potential of this approach on large-eddy simulations of isotropic turbulence, using the recovery of
statistical properties of direct numerical simulations as a reward. The closure model is a control policy enacted by cooperating
agents, which detect critical spatio-temporal patterns in the flow field to estimate the unresolved subgrid-scale physics. Results
obtained with multi-agent reinforcement learning algorithms based on experience replay compare favourably with established
modelling approaches. Moreover, we show that the learned turbulence models generalize across grid sizes and flow conditions.

2020

M. Chatzimanolakis, et al., “Optimal allocation of limited test resources for the quantification of COVID-19 infections,” Swiss Medical Weekly, vol. 150, no. w20445, 2020. Publisher's VersionAbstract
The systematic identification of infected individuals is critical for the containment of the COVID-19 pandemic. Currently, the spread of the disease is mostly quantified by the reported numbers of infections, hospitalisations, recoveries and deaths; these quantities inform epidemiology models that provide forecasts for the spread of the epidemic and guide policy making. The veracity of these forecasts depends on the discrepancy between the numbers of reported, and unreported yet infectious, individuals. We combine Bayesian experimental design with an epidemiology model and propose a methodology for the optimal allocation of limited testing resources in space and time, which maximises the information gain for such unreported infections. The proposed approach is applicable at the onset and spread of the epidemic and can forewarn of a possible recurrence of the disease after relaxation of interventions. We examine its application in Switzerland; the open source software is, however, readily adaptable to countries around the world. We find that following the proposed methodology can lead to vastly less uncertain predictions for the spread of the disease, thus improving estimates of the effective reproduction number and the future number of unreported infections. This information can provide timely and systematic guidance for the effective identification of infectious individuals and for decision-making regarding lockdown measures and the distribution of vaccines.
Z. Y. Wan, P. Karnakov, P. Koumoutsakos, and T. P. Sapsis, “Bubbles in turbulent flows: Data-driven, kinematic models with history terms,” Int. J. Multiphas. Flow, vol. 129, pp. 103286, 2020. Publisher's VersionAbstract
We present data driven kinematic models for the motion of bubbles in high-Re turbulent fluid flows based on recurrent neural networks with long-short term memory enhancements. The models extend empirical relations, such as Maxey-Riley (MR) and its variants, whose applicability is limited when either the bubble size is large or the flow is very complex. The recurrent neural networks are trained on the trajectories of bubbles obtained by Direct Numerical Simulations (DNS) of the Navier Stokes equations for a two-component incompressible flow model. Long short term memory components exploit the time history of the flow field that the bubbles have encountered along their trajectories and the networks are further augmented by imposing rotational invariance to their structure. We first train and validate the formulated model using DNS data for a turbulent Taylor-Green vortex. Then we examine the model pre- dictive capabilities and its generalization to Reynolds numbers that are different from those of the train- ing data on benchmark problems, including a steady (Hill’s spherical vortex) and an unsteady (Gaussian vortex ring) flow field. We find that the predictions of the developed model are significantly improved compared with those obtained by the MR equation. Our results indicate that data-driven models with history terms are well suited in capturing the trajectories of bubbles in turbulent flows.
P. R. Vlachas, et al., “Backpropagation algorithms and Reservoir Computing in Recurrent Neural Networks for the forecasting of complex spatiotemporal dynamics,” Neural Networks, vol. 126, pp. 191 - 217, 2020. Publisher's VersionAbstract
We examine the efficiency of Recurrent Neural Networks in forecasting the spatiotemporal dynamics of high dimensional and reduced order complex systems using Reservoir Computing (RC) and Backpropagation through time (BPTT) for gated network architectures. We highlight advantages and limitations of each method and discuss their implementation for parallel computing architectures. We quantify the relative prediction accuracy of these algorithms for the long-term forecasting of chaotic systems using as benchmarks the Lorenz-96 and the Kuramoto–Sivashinsky (KS) equations. We find that, when the full state dynamics are available for training, RC outperforms BPTT approaches in terms of predictive performance and in capturing of the long-term statistics, while at the same time requiring much less training time. However, in the case of reduced order data, large scale RC models can be unstable and more likely than the BPTT algorithms to diverge. In contrast, RNNs trained via BPTT show superior forecasting abilities and capture well the dynamics of reduced order systems. Furthermore, the present study quantifies for the first time the Lyapunov Spectrum of the KS equation with BPTT, achieving similar accuracy as RC. This study establishes that RNNs are a potent computational framework for the learning and forecasting of complex spatiotemporal systems.
P. Karnakov, S. Litvinov, and P. Koumoutsakos, “A hybrid particle volume-of-fluid method for curvature estimation in multiphase flows,” Int. J. Multiphas. Flow, vol. 125, pp. 103209, 2020. Publisher's VersionAbstract
We present a particle method for estimating the curvature of interfaces in volume-of-fluid simulations of multiphase flows. The method is well suited for under-resolved interfaces, and it is shown to be more accurate than the parabolic fitting that is employed in such cases. The curvature is computed from the equilibrium positions of particles constrained to circular arcs and attracted to the interface. The pro- posed particle method is combined with the method of height functions at higher resolutions, and it is shown to outperform the current combinations of height functions and parabolic fitting. The algorithm is conceptually simple and straightforward to implement on new and existing software frameworks for multiphase flow simulations thus enhancing their capabilities in challenging flow problems. We evaluate the proposed hybrid method on a number of two- and three-dimensional benchmark flow problems and illustrate its capabilities on simulations of flows involving bubble coalescence and turbulent multiphase flows.
P. Weber, G. Arampatzis, G. Novati, S. Verma, C. Papadimitriou, and P. Koumoutsakos, “Optimal Flow Sensing for Schooling Swimmers,” Biomimetics, vol. 5, no. 1, 2020. Publisher's VersionAbstract
Fish schooling implies an awareness of the swimmers for their companions. In flow mediated environments, in addition to visual cues, pressure and shear sensors on the fish body are critical for providing quantitative information that assists the quantification of proximity to other fish. Here we examine the distribution of sensors on the surface of an artificial swimmer so that it can optimally identify a leading group of swimmers. We employ Bayesian experimental design coupled with numerical simulations of the two-dimensional Navier Stokes equations for multiple self-propelled swimmers. The follower tracks the school using information from its own surface pressure and shear stress. We demonstrate that the optimal sensor distribution of the follower is qualitatively similar to the distribution of neuromasts on fish. Our results show that it is possible to identify accurately the center of mass and the number of the leading swimmers using surface only information.

2019

G. Novati, L. Mahadevan, and P. Koumoutsakos, “Controlled gliding and perching through deep-reinforcement-learning,” Phys. Rev. Fluids, vol. 4, no. 9, 2019. Publisher's VersionAbstract
Controlled gliding is one of the most energetically efficient modes of transportation for
natural and human powered fliers. Here we demonstrate that gliding and landing strategies
with different optimality criteria can be identified through deep-reinforcement-learning
without explicit knowledge of the underlying physics. We combine a two-dimensional
model of a controlled elliptical body with deep-reinforcement-learning (D-RL) to achieve
gliding with either minimum energy expenditure, or fastest time of arrival, at a predetermined location. In both cases the gliding trajectories are smooth, although energy/time
optimal strategies are distinguished by small/high frequency actuations. We examine the
effects of the ellipse’s shape and weight on the optimal policies for controlled gliding.
We find that the model-free reinforcement learning leads to more robust gliding than
model-based optimal control strategies with a modest additional computational cost. We
also demonstrate that the gliders with D-RL can generalize their strategies to reach
the target location from previously unseen starting positions. The model-free character
and robustness of D-RL suggests a promising framework for developing robotic devices
capable of exploiting complex flow environments.
S. M. H. Hashemi, et al., “A versatile and membrane-less electrochemical reactor for the electrolysis of water and brine,” Energ. Environ. Sci. 2019. Publisher's VersionAbstract
Renewables challenge the management of energy supply and demand due to their intermittency. A promising solution is the direct conversion of the excess electrical energy into valuable chemicals in electrochemical reactors that are inexpensive, scalable, and compatible with irregular availability of electrical power. Membrane-less electrolyzers, deployed on a microfluidic platform, were recently shown to hold great promise for efficient electrolysis and cost-effective operation. The elimination of the membrane increases the reactor lifetime, reduces fabrication costs, and enables the deployment of liquid electrolytes with ionic conductivities that surpass those allowed by solid membranes. Here, we demonstrate a membrane-less architecture that enables unprecedented throughput by 3D printing a device that combines components such as the flow plates and the fluidic ports in a monolithic part while at the same time providing tight tolerances and smooth surfaces for precise flow conditioning. We show that inertial fluidic forces are effective even in milifluidic regimes and, therefore, are utilized to control the two-phase flows inside the device and prevent cross-contamination of the products. Simulations provide insight on governing fluid dynamics of coalescing bubbles and their rapid jumps away from the electrodes and help identify three key mechanisms for their fast and intriguing return towards the electrodes. Experiments and simulations are used to demonstrate the efficiency of the inertial separation mechanism in milichannels and at higher flow rates than in microchannels. We analyze the performance of the present device for two reactions: water splitting and the Chlor-Alkali process and find product purities of more than 99% and Faradaic efficiencies of more than 90%. The present membrane-less reactor - containing more efficient catalysts - provides close to 40 times higher throughput than its microfluidic counterpart and paves the way for realization of cost-effective and scalable electrochemical stacks that meet the performance and price targets of the renewable energy sector.
W. Byeon, et al., “Automated identification and deep classification of cut marks on bones and its paleoanthropological implications,” J. Comput. Sci. vol. 32, pp. 36 - 43, 2019. Publisher's VersionAbstract
The identification of cut marks and other bone surface modifications (BSM) provides evidence for the emergence of meat-eating in human evolution. This most crucial part of taphonomic analysis of the archaeological human record has been controversial due to highly subjective interpretations of BSM. Here, we use a sample of 79 trampling and cut marks to compare the accuracy in mark identification on bones by human experts and computer trained algorithms. We demonstrate that deep convolutional neural networks (DCNN) and support vector machines (SVM) can recognize marks with accuracy that far exceeds that of human experts. Automated recognition and analysis of BSM using DCNN can achieve an accuracy of 91% of correct identification of cut and trampling marks versus a much lower accuracy rate (63%) obtained by trained human experts. This success underscores the capability of machine learning algorithms to help resolve controversies in taphonomic research and, more specifically, in the study of bone surface modifications. We envision that the proposed methods can help resolve on-going controversies on the earliest human meat-eating behaviors in Africa and other issues such as the earliest occupation of America.
U. Rasthofer, F. Wermelinger, P. Karnakov, J. Šukys, and P. Koumoutsakos, “Computational study of the collapse of a cloud with 12500 gas bubbles in a liquid,” Phys. Rev. Fluids, vol. 4, pp. 063602, 2019. Publisher's VersionAbstract
We investigate the collapse of a cloud composed of 12500 gas bubbles in a liquid through large-scale simulations. The gas bubbles are discretized by a diffuse interface method, and a finite volume scheme is used to solve on a structured Cartesian grid the Euler equations. We investigate the propagation of the collapse wave front through the cloud and provide comparisons to existing models such as Mørch’s ordinary differential equations and a homogeneous mixture approach. We analyze the flow field to examine the evolution of individual gas bubbles and in particular their associated microjet. We find that the velocity magnitude of the microjets depends on the local strength of the collapse wave and hence on the radial position of the bubbles in the cloud. At the same time, the direction of the microjets is influenced by the distribution of the bubbles in its vicinity. We envision that the present, state-of-the-art, large-scale simulations will serve the further development of low-order models for bubble collapse.

2018

F. Wermelinger, U. Rasthofer, P. E. Hadjidoukas, and P. Koumoutsakos, “Petascale simulations of compressible flows with interfaces,” J. Comput. Sci. vol. 26, pp. 217–225, 2018. Publisher's VersionAbstract
We demonstrate a high throughput software for the efficient simulation of compressible multicomponent flow on high performance computing platforms. The discrete problem is represented on structured three-dimensional grids with non-uniform resolution. Discontinuous flow features are captured using a diffuse interface method. A distinguishing characteristic of the method is the proper treatment of the interface zone as a mixing region of liquid and gas. The governing equations are discretized by a Godunov-type finite volume method with explicit time stepping using a low-storage Runge-Kutta scheme. The presented flow solver Cubism-MPCF is based on our Cubism library which enables a highly optimized framework for the efficient treatment of stencil based problems on multicore architectures. The framework is general and not limited to applications in fluid dynamics. We validate our solver by classical benchmark examples. Furthermore, we examine a highly-resolved shock-induced bubble collapse and a cloud of $O(10^3)$ collapsing bubbles, which demonstrate the high potential of the proposed framework and solver.
Z. Y. Wan, P. R. Vlachas, P. Koumoutsakos, and T. P. Sapsis, “Data-assisted reduced-order modeling of extreme events in complex dynamical systems,” PLoS ONE, vol. 13, no. 5, pp. 1-22, 2018. Publisher's VersionAbstract
The prediction of extreme events, from avalanches and droughts to tsunamis and epidemics, depends on the formulation and analysis of relevant, complex dynamical systems. Such dynamical systems are characterized by high intrinsic dimensionality with extreme events having the form of rare transitions that are several standard deviations away from the mean. Such systems are not amenable to classical order-reduction methods through projection of the governing equations due to the large intrinsic dimensionality of the underlying attractor as well as the complexity of the transient events. Alternatively, data-driven techniques aim to quantify the dynamics of specific, critical modes by utilizing data-streams and by expanding the dimensionality of the reduced-order model using delayed coordinates. In turn, these methods have major limitations in regions of the phase space with sparse data, which is the case for extreme events. In this work, we develop a novel hybrid framework that complements an imperfect reduced order model, with data-streams that are integrated though a recurrent neural network (RNN) architecture. The reduced order model has the form of projected equations into a low-dimensional subspace that still contains important dynamical information about the system and it is expanded by a long short-term memory (LSTM) regularization. The LSTM-RNN is trained by analyzing the mismatch between the imperfect model and the data-streams, projected to the reduced-order space. The data-driven model assists the imperfect model in regions where data is available, while for locations where data is sparse the imperfect model still provides a baseline for the prediction of the system state. We assess the developed framework on two challenging prototype systems exhibiting extreme events. We show that the blended approach has improved performance compared with methods that use either data streams or the imperfect model alone. Notably the improvement is more significant in regions associated with extreme events, where data is sparse.
P. R. Vlachas, W. Byeon, Z. Y. Wan, T. P. Sapsis, and P. Koumoutsakos, “Data-driven forecasting of high-dimensional chaotic systems with long short-term memory networks,” P. Roy. Soc. A-Math. Phy. vol. 474, no. 2213, pp. 20170844, 2018. Publisher's VersionAbstract
We introduce a data-driven forecasting method for high-dimensional chaotic systems using long shortterm memory (LSTM) recurrent neural networks. The proposed LSTM neural networks perform inference of high-dimensional dynamical systems in their reduced order space and are shown to be an effective set of nonlinear approximators of their attractor. We demonstrate the forecasting performance of the LSTM and compare it with Gaussian processes (GPs) in time series obtained from the Lorenz 96 system, the Kuramoto–Sivashinsky equation and a prototype climate model. The LSTM networks outperform the GPs in short-term forecasting accuracy in all applications considered. A hybrid architecture, extending the LSTM with a mean stochastic model (MSM–LSTM), is proposed to ensure convergence to the invariant measure. This novel hybrid method is fully data-driven and extends the forecasting capabilities of LSTM networks.

2017

G. Novati, S. Verma, D. Alexeev, D. Rossinelli, W. M. Van Rees, and P. Koumoutsakos, “Synchronisation through learning for two self-propelled swimmers,” Bioinspir. Biomim. vol. 12, no. 3, pp. 036001, 2017. Publisher's VersionAbstract
The coordinated motion by multiple swimmers is a fundamental component in fish schooling. The flow field induced by the motion of each self-propelled swimmer implies non-linear hydrodynamic interactions among the members of a group. How do swimmers compensate for such hydrodynamic interactions in coordinated patterns? We provide an answer to this riddle though simulations of two, self-propelled, fish-like bodies that employ a learning algorithm to synchronise their swimming patterns. We distinguish between learned motion patterns and the commonly used a-priori specified movements, that are imposed on the swimmers without feedback from their hydrodynamic interactions. First, we demonstrate that two rigid bodies executing pre-specified motions, with an alternating leader and follower, can result in substantial drag-reduction and intermittent thrust generation. In turn, we study two self-propelled swimmers arranged in a leader-follower configuration, with a-priori specified body-deformations. These two self-propelled swimmers do not sustain their tandem configuration. The follower experiences either an increase or decrease in swimming speed, depending on the initial conditions, while the swimming of the leader remains largely unaffected. This indicates that a-priori specified patterns are not sufficient to sustain synchronised swimming. We then examine a tandem of swimmers where the leader has a steady gait and the follower learns to synchronize its motion, to overcome the forces induced by the leader’s vortex wake. The follower employs reinforcement learning to adapt its swimming-kinematics so as to minimize its lateral deviations from the leader’s path. Swimming in such a sustained synchronised tandem yields up to 30% reduction in energy expenditure for the follower, in addition to a 20% increase in its swimming-efficiency. The present results show that two self-propelled swimmers can be synchronised by adapting their motion patterns to compensate for flow-structure interactions. Moreover, swimmers can exploit the vortical structures of their flow field so that synchronised swimming is energetically beneficial.
B. Mosimann, et al., “Reference Ranges for Fetal Atrioventricular and Ventriculoatrial Time Intervals and Their Ratios during Normal Pregnancy,” Fetal Diagn. Ther. 2017. Publisher's VersionAbstract
Background: The diagnostic assessment of fetal arrhythmias relies on the measurements of atrioventricular (AV) and ventriculoatrial (VA) time intervals. Pulsed Doppler over in- and outflow of the left ventricle and tissue Doppler imaging are well-described methods, while Doppler measurements between the left brachiocephalic vein and the aortic arch are less investigated. The aim of this study was to compare these methods of measurement, to find influencing factors on AV and VA times and their ratio, and to create reference ranges. Methods: Echocardiography was performed between 16 and 40 weeks of gestation in normal singleton pregnancies. Nomograms for the individual measurements were created using quantile regression with Matlab Data Analytics. Statistical analyses were performed with GraphPad version 5.0 for Windows. Results: A total of 329 pregnant women were en- rolled. A significant correlation exists between AV and VA times and gestational age (GA) (p = 0.0104 to <0.0001, σ = 0.1412 to 0.3632). No correlation was found between the AV:VA ratio and GA (p = 0.08 to 0.60). All measurements differed significantly amongst the studied methods (p < 0.0001). Conclusions: AV and VA intervals increase proportionally with GA; no other independent influencing factors could be identified. As significant differences exist between the three methods of assessment, it is crucial to use appropriate reference ranges to diagnose pathologies.
L. Kulakova, G. Arampatzis, P. Angelikopoulos, P. Hadjidoukas, C. Papadimitriou, and P. Koumoutsakos, “Data driven inference for the repulsive exponent of the Lennard-Jones potential in molecular dynamics simulations,” Sci. Rep.-UK, vol. 7, no. 1, pp. 16576, 2017. Publisher's VersionAbstract
The Lennard-Jones (LJ) potential is a cornerstone of Molecular Dynamics (MD) simulations and among the most widely used computational kernels in science. The LJ potential models atomistic attraction and repulsion with century old prescribed parameters ($q=6$, $p=12$ respectively), originally related by a factor of two for simplicity of calculations. We propose the inference of the repulsion exponent through Hierarchical Bayesian uncertainty quantification We use experimental data of the radial distribution function and dimer interaction energies from quantum mechanics simulations. We find that the repulsion exponent $p\approx6.5$ provides an excellent fit for the experimental data of liquid argon, for a range of thermodynamic conditions, as well as for saturated argon vapour. Calibration using the quantum simulation data did not provide a good fit in these cases. However, values $p\approx12.7$ obtained by dimer quantum simulations are preferred for the argon gas while lower values are promoted by experimental data. These results show that the proposed LJ 6-p potential applies to a wider range of thermodynamic conditions, than the classical LJ 6-12 potential. We suggest that calibration of the repulsive exponent in the LJ potential widens the range of applicability and accuracy of MD simulations.

2016

S. Wu, P. Angelikopoulos, G. Tauriello, C. Papadimitriou, and P. Koumoutsakos, “Fusing heterogeneous data for the calibration of molecular dynamics force fields using hierarchical Bayesian models,” J. Chem. Phys. vol. 145, no. 24, pp. 244112, 2016. Publisher's VersionAbstract
We propose a hierarchical Bayesian framework to systematically integrate heterogeneous data for the calibration of force fields in Molecular Dynamics (MD) simulations. Our approach enables the fusion of diverse experimental data sets of the physico-chemical properties of a system at different thermodynamic conditions. We demonstrate the value of this framework for the robust calibration of MD force-fields for water using experimental data of its diffusivity, radial distribution function, and density. In order to address the high computational cost associated with the hierarchical Bayesian models, we develop a novel surrogate model based on the empirical interpolation method. Further computational savings are achieved by implementing a highly parallel transitional Markov chain Monte Carlo technique. The present method bypasses possible subjective weightings of the experimental data in identifying MD force-field parameters.
J. Chen, J. H. Walther, and P. Koumoutsakos, “Ultrafast cooling by covalently bonded graphene-carbon nanotube hybrid immersed in water,” Nanotechnology, vol. 27, no. 46, pp. 465705, 2016. Publisher's VersionAbstract
The increasing power density and the decreasing dimensions of transistors present severe thermal challenges to the design of modern microprocessors. Furthermore, new technologies such as three-dimensional chip-stack architectures require novel cooling solutions for their thermal management. Here, we demonstrate, through transient heat-dissipation simulations, that a covalently bonded graphene-carbon nanotube (G-CNT) hybrid immersed in water is a promising solution for the ultrafast cooling of such high-temperature and high heat-flux surfaces. The G-CNT hybrid offers a unique platform to integrate the superior axial heat transfer capability of individual CNTs via their parallel arrangement. The immersion of the G-CNT in water enables an additional heat dissipation path via the solid–liquid interaction, allowing for the sustainable cooling of the hot surface under a constant power input of up to 10 000 W cm-2.

2015

M. M. Hejlesen, P. Koumoutsakos, A. Leonard, and J. H. Walther, “Iterative Brinkman penalization for remeshed vortex methods,” J. Comput. Phys. vol. 280, pp. 547–562, 2015. Publisher's VersionAbstract
We introduce an iterative Brinkman penalization method for the enforcement of the no-slip boundary condition in remeshed vortex methods. In the proposed method, the Brinkman penalization is applied iteratively only in the neighborhood of the body. This allows for using significantly larger time steps, than what is customary in the Brinkman penalization, thus reducing its computational cost while maintaining the capability of the method to handle complex geometries. We demonstrate the accuracy of our method by considering challenging benchmark problems such as flow past an impulsively started cylinder and normal to an impulsively started and accelerated flat plate. We find that the present method enhances significantly the accuracy of the Brinkman penalization technique for the simulations of highly unsteady flows past complex geometries.
G. Tauriello and P. Koumoutsakos, “A comparative study of penalization and phase field methods for the solution of the diffusion equation in complex geometries,” J. Comput. Phys. vol. 283, pp. 388–407, 2015. Publisher's VersionAbstract
We present a comparative study of penalization and phase field methods for the solution of the diffusion equation in complex geometries embedded using simple Cartesian meshes. The two methods have been widely employed to solve partial differential equations in complex and moving geometries for applications ranging from solid and fluid mechanics to biology and geophysics. Their popularity is largely due to their discretization on Cartesian meshes thus avoiding the need to create body-fitted grids. At the same time, there are questions regarding their accuracy and it appears that the use of each one is confined by disciplinary boundaries. Here, we compare penalization and phase field methods to handle problems with Neumann and Robin boundary conditions. We discuss extensions for Dirichlet boundary conditions and in turn compare with methods that have been explicitly designed to handle Dirichlet boundary conditions. The accuracy of all methods is analyzed using one and two dimensional benchmark problems such as the flow induced by an oscillating wall and by a cylinder performing rotary oscillations. This comparative study provides information to decide which methods to consider for a given application and their incorporation in broader computational frameworks. We demonstrate that phase field methods are more accurate than penalization methods on problems with Neumann boundary conditions and we present an error analysis explaining this result.
P. E. Hadjidoukas, P. Angelikopoulos, C. Papadimitriou, and P. Koumoutsakos, “Pi4U: A high performance computing framework for Bayesian uncertainty quantification of complex models,” J. Comput. Phys. vol. 284, pp. 1–21, 2015. Publisher's VersionAbstract
We present Pi4U, an extensible framework, for non-intrusive Bayesian Uncertainty Quantification and Propagation (UQ+P) of complex and computationally demanding physical models, that can exploit massively parallel computer architectures. The framework incorporates Laplace asymptotic approximations as well as stochastic algorithms, along with distributed numerical differentiation and task-based parallelism for heterogeneous clusters. Sampling is based on the Transitional Markov Chain Monte Carlo (TMCMC) algorithm and its variants. The optimization tasks associated with the asymptotic approximations are treated via the Covariance Matrix Adaptation Evolution Strategy (CMA-ES). A modified subset simulation method is used for posterior reliability measurements of rare events. The framework accommodates scheduling of multiple physical model evaluations based on an adaptive load balancing library and shows excellent scalability. In addition to the software framework, we also provide guidelines as to the applicability and efficiency of Bayesian tools when applied to computationally demanding physical models. Theoretical and computational developments are demonstrated with applications drawn from molecular dynamics, structural dynamics and granular flow.
D. Rossinelli, B. Hejazialhosseini, W. van Rees, M. Gazzola, M. Bergdorf, and P. Koumoutsakos, “MRAG-I2D: Multi-resolution adapted grids for remeshed vortex methods on multicore architectures,” J. Comput. Phys. vol. 288, pp. 1–18, 2015. Publisher's VersionAbstract
We present MRAG-I2D, an open source software framework, for multiresolution simulations of two-dimensional, incompressible, viscous flows on multicore architectures. The spatiotemporal scales of the flow field are captured by remeshed vortex methods enhanced by high order average-interpolating wavelets and local time-stepping. The multiresolution solver of the Poisson equation relies on the development of a novel, tree-based multipole method. MRAG-I2D implements a number of HPC strategies to map efficiently the irregular computational workload of wavelet-adapted grids on multicore nodes. The capabilities of the present software are compared to the current state-of-the-art in terms of accuracy, compression rates and time-to-solution. Benchmarks include the inviscid evolution of an elliptical vortex, flow past an impulsively started cylinder at Re = 40 - 40 000 and simulations of self-propelled anguilliform swimmers. The results indicate that the present software has the same or better accuracy than state-of-the-art solvers while it exhibits unprecedented performance in terms of time-to-solution.
S. D. Finley, P. Angelikopoulos, P. KOUMOUTSAKOS, and A. S. Popel, “Pharmacokinetics of Anti-VEGF Agent Aflibercept in Cancer Predicted by Data-Driven, Molecular-Detailed Model,” CPT: PSP, vol. 4, no. 11, pp. 641–649, 2015. Publisher's VersionAbstract
Mathematical models can support the drug development process by predicting the pharmacokinetic (PK) properties of the drug and optimal dosing regimens. We have developed a pharmacokinetic model that includes a biochemical molecular interaction network linked to a whole-body compartment model. We applied the model to study the PK of the anti-vascular endothelial growth factor (VEGF) cancer therapeutic agent, aflibercept. Clinical data is used to infer model parameters using a Bayesian approach, enabling a quantitative estimation of the contributions of specific transport processes and molecular interactions of the drug that cannot be examined in other PK modeling, and insight into the mechanisms of aflibercept's antiangiogenic action. Additionally, we predict the plasma and tissue concentrations of unbound and VEGF-bound aflibercept. Thus, we present a computational framework that can serve as a valuable tool for drug development efforts.
P. B. de Reuille, et al., “MorphoGraphX: A platform for quantifying morphogenesis in 4D,” eLife, vol. 4, 2015. Publisher's VersionAbstract
Morphogenesis emerges from complex multiscale interactions between genetic and mechanical processes. To understand these processes, the evolution of cell shape, proliferation and gene expression must be quantified. This quantification is usually performed either in full 3D, which is computationally expensive and technically challenging, or on 2D planar projections, which introduces geometrical artifacts on highly curved organs. Here we present MorphoGraphX (www.MorphoGraphX.org), a software that bridges this gap by working directly with curved surface images extracted from 3D data. In addition to traditional 3D image analysis, we have developed algorithms to operate on curved surfaces, such as cell segmentation, lineage tracking and fluorescence signal quantification. The software's modular design makes it easy to include existing libraries, or to implement new algorithms. Cell geometries extracted with MorphoGraphX can be exported and used as templates for simulation models, providing a powerful platform to investigate the interactions between shape, genes and growth.
J. Chen, J. H. Walther, and P. Koumoutsakos, “Covalently Bonded Graphene-Carbon Nanotube Hybrid for High-Performance Thermal Interfaces,” Adv. Funct. Mater. vol. 25, no. 48, pp. 7539–7545, 2015. Publisher's VersionAbstract
The remarkable thermal properties of graphene and carbon nanotubes (CNTs) have been the subject of intensive investigations for the thermal management of integrated circuits. However, the small contact area of CNTs and the large anisotropic heat conduction of graphene have hindered their applications as effective thermal interface materials (TIMs). Here, a covalently bonded graphene–CNT (G-CNT) hybrid is presented that multiplies the axial heat transfer capability of individual CNTs through their parallel arrangement, while at the same time it provides a large contact area for efficient heat extraction. Through computer simulations, it is demonstrated that the G-CNT outperforms few-layer graphene by more than 2 orders of magnitude for the c-axis heat transfer, while its thermal resistance is 3 orders of magnitude lower than the state-of-the-art TIMs. We show that heat can be removed from the G-CNT by immersing it in a liquid. The heat transfer characteristics of G-CNT suggest that it has the potential to revolutionize the design of high-performance TIMs.
W. M. van Rees, M. Gazzola, and P. Koumoutsakos, “Optimal morphokinematics for undulatory swimmers at intermediate Reynolds numbers,” J. Fluid Mech. vol. 775, pp. 178–188, 2015. Publisher's VersionAbstract
Undulatory locomotion is an archetypal mode of propulsion for natural swimmers across scales. Undulatory swimmers convert transverse body oscillations into forward velocity by a complex interplay between their flexural movements, morphological features and the fluid environment. Natural evolution has produced a wide range of morphokinematic examples of undulatory swimmers that often serve as inspiration for engineering devices. It is, however, unknown to what extent natural swimmers are optimized for hydrodynamic performance. In this work, we reverse-engineer the morphology and gait for fast and efficient swimmers by coupling an evolution strategy to three-dimensional direct numerical simulations of flows at intermediate Reynolds numbers. The fastest swimmer is slender with a narrow tail fin and performs a sequence of C-starts to maximize its average velocity. The most efficient swimmer combines moderate transverse movements with a voluminous head, tapering into a streamlined profile via a pronounced inflection point. These optimal solutions outperform anguilliform swimming zebrafish in both efficiency and speed. We investigate the transition between morphokinematic solutions in the speed–energy space, laying the foundations for the design of high-performance artificial swimming devices.
M. U. Baumann, et al., “Placental plasticity in monochorionic twins: Impact on birth weight and placental weight,” Placenta, vol. 36, no. 9, pp. 1018–1023, 2015. Publisher's VersionAbstract
Introduction The knowledge about adaptive mechanisms of monochorionic placentas to fulfill the demands of two instead of one fetus is largely speculative. The aim of our study was to investigate the impact of chorionicity on birth weight and placental weight in twin pregnancies. Methods Forty Monochorionic (MC) and 43 dichorionic (DC) twin pregnancies were included in this retrospective study. Individual and total (sum of both twins) birth weights, placental weights ratios between placental and birth weights and observed-to-expected (O/E)-ratios were calculated and analyzed. Additionally, we investigated whether in twin pregnancies placental and birth weights follow the law of allometric metabolic scaling. Results MC pregnancies showed higher placental O/E-ratios than DC ones (2.25 \pm 0.85 versus 1.66 \pm 0.61; p < 0.05), whereas the total neonatal birth weight O/E-ratios were not different. In DC twins total placental weights correlated significantly with gestational age (r = 0.74, p < 0.001), but not in MC twins. Analysis of deliveries łeq}32 weeks revealed that the placenta to birth weight ratio in MC twins was higher than in matched DC twins (0.49 \pm 0.3 versus 0.24 \pm 0.03; p = 0.03). Allometric metabolic scaling revealed that dichorionic twin placentas scale with birth weight, while the monochorionic ones do not. Discussion The weight of MC placentas compared to that of DC is not gestational age dependent in the third trimester. Therefore an early accelerated placental growth pattern has to be postulated which leads to an excess placental mass particularly below 32 weeks of gestation. The monochorionic twins do not follow allometric metabolic scaling principle making them more vulnerable to placental compromise. Keywords Twin pregnancy Monochorionic Dichorionic Birth weight Placental weight Fetal:placental ratio Allometric metabolic scaling
W. M. van Rees, G. Novati, and P. Koumoutsakos, “Self-propulsion of a counter-rotating cylinder pair in a viscous fluid,” Phys. Fluids, vol. 27, no. 6, pp. 063102, 2015. Publisher's VersionAbstract
We study a self-propelling pair of steadily counter-rotating cylinders in simulations of a two-dimensional viscous fluid. We find two strikingly, opposite directions for the motion of the pair that is characterized by its width and rotational Reynolds number. At low Reynolds numbers and large widths, the cylinder pair moves similarly to an inviscid point vortex pair, while at higher Reynolds numbers and smaller widths, the pair moves in the opposite direction through a jet-like propulsion mechanism. Increasing further the Reynolds number, or decreasing the width, gives rise to non-polarised motion governed by the shedding direction and frequency of the boundary-layer vorticity. We discuss the fundamental physical mechanisms for these two types of motion and the transitions in the corresponding phase diagram. We discuss the fluid dynamics of each regime based on streamline plots, tracer particles, and the vorticity field. The counter rotating cylinder pair serves as a prototype for self-propelled bodies and suggests possible engineering devices composed of simple components and tunable by the rotation and width of the cylinder pair.
P. Angelikopoulos, C. Papadimitriou, and P. Koumoutsakos, “X-TMCMC: Adaptive kriging for Bayesian inverse modeling,” Comput. Method. Appl. M. vol. 289, pp. 409–428, 2015. Publisher's VersionAbstract
The Bayesian inference of models associated with large-scale simulations is prohibitively expensive even for massively parallel architectures. We demonstrate that we can drastically reduce this cost by combining adaptive kriging with the population-based Transitional Markov Chain Monte Carlo (TMCMC) techniques. For uni-modal posterior probability distribution functions (PDF), the proposed hybrid method can reduce the computational cost by an order of magnitude with the same computational resources. For complex posterior PDF landscapes we show that it is necessary to further extend the TMCMC by Langevin adjusted proposals. The proposed hybrid method exhibits high parallel efficiency. We demonstrate the capabilities of our method on test bed problems and on high fidelity simulations in structural dynamics.
A. Popadić, M. Praprotnik, P. Koumoutsakos, and J. H. Walther, “Continuum simulations of water flow past fullerene molecules,” Eur. Phys. J.-Spec. Top. vol. 224, no. 12, pp. 2321–2330, 2015. Publisher's VersionAbstract
We present continuum simulations of water flow past fullerene molecules. The governing Navier-Stokes equations are complemented with the Navier slip boundary condition with a slip length that is extracted from related molecular dynamics simulations. We find that several quantities of interest as computed by the present model are in good agreement with results from atomistic and atomistic-continuum simulations at a fraction of the cost. We simulate the flow past a single fullerene and an array of fullerenes and demonstrate that such nanoscale flows can be computed efficiently by continuum flow solvers, allowing for investigations into spatiotemporal scales inaccessible to atomistic simulations.
D. Alexeev, J. Chen, J. H. Walther, K. P. Giapis, P. Angelikopoulos, and P. Koumoutsakos, “Kapitza Resistance between Few-Layer Graphene and Water: Liquid Layering Effects,” Nano Lett. vol. 15, no. 9, pp. 5744–5749, 2015. Publisher's VersionAbstract
The Kapitza resistance (R_K) between few-layer graphene (FLG) and water was studied using molecular dynamics simulations. The R_K was found to depend on the number of the layers in the FLG though, surprisingly, not on the water block thickness. This distinct size dependence is attributed to the large difference in the phonon mean free path between the FLG and water. Remarkably, R_K is strongly dependent on the layering of water adjacent to the FLG, exhibiting an inverse proportionality relationship to the peak density of the first water layer, which is consistent with better acoustic phonon matching between FLG and water. These findings suggest novel ways to engineer the thermal transport properties of solid–liquid interfaces by controlling and regulating the liquid layering at the interface.
P. R. Jones, et al., “Sustaining dry surfaces under water,” Sci. Rep.-UK, vol. 5, no. 1, 2015. Publisher's VersionAbstract
Rough surfaces immersed under water remain practically dry if the liquid-solid contact is on roughness peaks, while the roughness valleys are filled with gas. Mechanisms that prevent water from invading the valleys are well studied. However, to remain practically dry under water, additional mechanisms need consideration. This is because trapped gas (e.g. air) in the roughness valleys can dissolve into the water pool, leading to invasion. Additionally, water vapor can also occupy the roughness valleys of immersed surfaces. If water vapor condenses, that too leads to invasion. These effects have not been investigated, and are critically important to maintain surfaces dry under water. In this work, we identify the critical roughness scale, below which it is possible to sustain the vapor phase of water and/or trapped gases in roughness valleys – thus keeping the immersed surface dry. Theoretical predictions are consistent with molecular dynamics simulations and experiments.
S. Wu, P. Angelikopoulos, C. Papadimitriou, R. Moser, and P. Koumoutsakos, “A hierarchical Bayesian framework for force field selection in molecular dynamics simulations,” Philos. T. Roy. Soc. A, vol. 374, no. 2060, pp. 20150032, 2015. Publisher's VersionAbstract
We present a hierarchical Bayesian framework for the selection of force fields in molecular dynamics (MD) simulations. The framework associates the variability of the optimal parameters of the MD potentials under different environmental conditions with the corresponding variability in experimental data. The high computational cost associated with the hierarchical Bayesian framework is reduced by orders of magnitude through a parallelized Transitional Markov Chain Monte Carlo method combined with the Laplace Asymptotic Approximation. The suitability of the hierarchical approach is demonstrated by performing MD simulations with prescribed parameters to obtain data for transport coefficients under different conditions, which are then used to infer and evaluate the parameters of the MD model. We demonstrate the selection of MD models based on experimental data and verify that the hierarchical model can accurately quantify the uncertainty across experiments; improve the posterior probability density function estimation of the parameters, thus, improve predictions on future experiments; identify the most plausible force field to describe the underlying structure of a given dataset. The framework and associated software are applicable to a wide range of nanoscale simulations associated with experimental data with a hierarchical structure.
G. Tauriello, H. M. Meyer, R. S. Smith, P. Koumoutsakos, and A. H. K. Roeder, “Variability and constancy in cellular growth of Arabidopsis sepals,” Plant Physiol. pp. 2342–2358, 2015. Publisher's VersionAbstract
Growth of tissues is highly reproducible; yet, growth of individual cells in a tissue is highly variable, and neighboring cells can grow at different rates. We analyzed the growth of epidermal cell lineages in the Arabidopsis (Arabidopsis thaliana) sepal to determine how the growth curves of individual cell lineages relate to one another in a developing tissue. To identify underlying growth trends, we developed a continuous displacement field to predict spatially averaged growth rates. We showed that this displacement field accurately describes the growth of sepal cell lineages and reveals underlying trends within the variability of in vivo cellular growth. We found that the tissue, individual cell lineages, and cell walls all exhibit growth rates that are initially low, accelerate to a maximum, and decrease again. Accordingly, these growth curves can be represented by sigmoid functions. We examined the relationships among the cell lineage growth curves and surprisingly found that all lineages reach the same maximum growth rate relative to their size. However, the cell lineages are not synchronized; each cell lineage reaches this same maximum relative growth rate but at different times. The heterogeneity in observed growth results from shifting the same underlying sigmoid curve in time and scaling by size. Thus, despite the variability in growth observed in our study and others, individual cell lineages in the developing sepal follow similarly shaped growth curves.

2014

W. M. Van Rees, D. Rossinelli, P. Hadjidoukas, and P. Koumoutsakos, “High performance CPU/GPU multiresolution Poisson solver,” Adv. Par. Com. vol. 1, no. 1, pp. 481–490, 2014. Publisher's VersionAbstract
We present a multipole-based N-body solver for 3D multiresolution, block-structured grids. The solver is designed for a single heterogeneous CPU/GPU compute node, and evaluates the multipole expansions on the CPU while offloading the compute-heavy particle-particle interactions to the GPU. The regular structure of the destination points is exploited for data parallelism on the CPU, to reduce data transfer to the GPU and to minimize memory accesses during evaluation of the direct and indirect interactions. The algorithmic improvements together with HPC techniques lead to 81% and 96% of the upper bound performance for the CPU and GPU parts, respectively.
A. Popadić, J. H. Walther, P. KOUMOUTSAKOS, and M. Praprotnik, “Continuum simulations of water flow in carbon nanotube membranes,” New J. Phys. vol. 16, no. 8, pp. 082001, 2014. Publisher's VersionAbstract
We propose the use of the Navier–Stokes equations subject to partial-slip boundary conditions to simulate water flows in Carbon NanoTube (CNT) membranes. The finite volume discretizations of the Navier–Stokes equations are combined with slip lengths extracted from molecular dynamics (MD) simulations to predict the pressure losses at the CNT entrance as well as the enhancement of the flow rate in the CNT. The flow quantities calculated from the present hybrid approach are in excellent agreement with pure MD results while they are obtained at a fraction of the computational cost. The method enables simulations of system sizes and times well beyond the present capabilities of MD simulations. Our simulations provide an asymptotic flow rate enhancement and indicate that the pressure losses at the CNT ends can be reduced by reducing their curvature. More importantly, our results suggest that flows at nanoscale channels can be described by continuum solvers with proper boundary conditions that reflect the molecular interactions of the liquid with the walls of the nanochannel.
F. Milde, G. Tauriello, H. Haberkern, and P. Koumoutsakos, “SEM++: A particle model of cellular growth, signaling and migration,” Computational particle mechanics, vol. 1, no. 2, pp. 211–227, 2014. Publisher's VersionAbstract
We present a discrete particle method to model biological processes from the sub-cellular to the inter-cellular level. Particles interact through a parametrized force field to model cell mechanical properties, cytoskeleton remodeling, growth and proliferation as well as signaling between cells. We discuss the guiding design principles for the selection of the force field and the validation of the particle model using experimental data. The proposed method is integrated into a multiscale particle framework for the simulation of biological systems.
P. E. Hadjidoukas, P. Angelikopoulos, D. Rossinelli, D. Alexeev, C. Papadimitriou, and P. Koumoutsakos, “Bayesian uncertainty quantification and propagation for discrete element simulations of granular materials,” Comput. Method. Appl. M. vol. 282, pp. 218–238, 2014. Publisher's VersionAbstract
Predictions in the behavior of granular materials using Discrete Element Methods (DEM) hinge on the employed interaction potentials. Here we introduce a data driven, Bayesian framework to quantify DEM predictions. Our approach relies on experimentally measured coefficients of restitution for single steel particle–wall collisions. The calibration data entail both tangential and normal coefficients of restitution, for varying impact angles and speeds of the bouncing particle. The parametric uncertainty in multiple Force–Displacement models is estimated using an enhanced Transitional Markov Chain Monte Carlo implemented efficiently on parallel computer architectures. In turn, the parametric model uncertainties are propagated to predict Quantities of Interest (QoI) for two testbed applications: silo discharge and vibration induced mass-segregation. This work demonstrates that the classical way of calibrating DEM potentials, through parameter optimization, is insufficient and it fails to provide robust predictions. The present Bayesian framework provides robust predictions for the behavior of granular materials using DEM simulations. Most importantly the results demonstrate the importance of including parametric and modeling uncertainties in the potentials employed in Discrete Element Methods.
M. Gazzola, B. Hejazialhosseini, and P. Koumoutsakos, “Reinforcement Learning and Wavelet Adapted Vortex Methods for Simulations of Self-propelled Swimmers,” SIAM J. Sci. Comput. vol. 36, no. 3, pp. B622–B639, 2014. Publisher's VersionAbstract
We present a numerical method for the simulation of collective hydrodynamics in self-propelled swimmers. Swimmers in a viscous incompressible flow are simulated with a remeshed vortex method coupled with Brinkman penalization and projection approach. The remeshed vortex methods are enhanced via wavelet based adaptivity in space and time. The method is validated on benchmark swimming problems. Furthermore the flow solver is integrated with a reinforcement learning algorithm, such that swimmers can learn to adapt their motion so as to optimally achieve a specified goal, such as fish schooling. The computational efficiency of the wavelet adapted remeshed vortex method is a key aspect for the effective coupling with learning algorithms. The suitability of this approach for the identification of swimming behaviors is assessed on a set of learning tasks.
E. R. Cruz-Chu, A. Malafeev, T. Pajarskas, I. V. Pivkin, and P. Koumoutsakos, “Structure and Response to Flow of the Glycocalyx Layer,” Biophys. J. vol. 106, no. 1, pp. 232–243, 2014. Publisher's VersionAbstract
The glycocalyx is a sugar-rich layer located at the luminal part of the endothelial cells. It is involved in key metabolic processes and its malfunction is related to several diseases. To understand the function of the glycocalyx, a molecular level characterization is necessary. In this article, we present large-scale molecular-dynamics simulations that provide a comprehensive description of the structure and dynamics of the glycocalyx. We introduce the most detailed, to-date, all-atom glycocalyx model, composed of lipid bilayer, proteoglycan dimers, and heparan sulfate chains with realistic sequences. Our results reveal the folding of proteoglycan ectodomain and the extended conformation of heparan sulfate chains. Furthermore, we study the glycocalyx response under shear flow and its role as a flypaper for binding fibroblast growth factors (FGFs), which are involved in diverse functions related to cellular differentiation, including angiogenesis, morphogenesis, and wound healing. The simulations show that the glycocalyx increases the effective concentration of FGFs, leading to FGF oligomerization, and acts as a lever to transfer mechanical stimulus into the cytoplasmic side of endothelial cells.
J. Chen, J. H. Walther, and P. Koumoutsakos, “Strain Engineering of Kapitza Resistance in Few-Layer Graphene,” Nano Lett. vol. 14, no. 2, pp. 819–825, 2014. Publisher's VersionAbstract
We demonstrate through molecular dynamics simulations that the Kapitza resistance in few-layer graphene (FLG) can be controlled by applying mechanical strain. For unstrained FLG, the Kapitza resistance decreases with the increase of thickness and reaches an asymptotic value of 6 × 10–10 m\^2 K/W at a thickness about 16 nm. Uniaxial cross-plane strain is found to increase the Kapitza resistance in FLG monotonically, when the applied strain varies from compressive to tensile. Moreover, uniaxial strain couples the in-plane and out-of-plane strain/stress when the surface of FLG is buckled. We find that with a compressive cross-plane stress of 2 GPa, the Kapitza resistance is reduced by about 50%. On the other hand it is almost tripled with a tensile cross-plane stress of 1 GPa. Remarkably, compressive in-plane strain can either increase or reduce the Kapitza resistance, depending on the specific way it is applied. Our study suggests that graphene can be exploited for both heat dissipation and insulation through strain engineering.

2013

W. M. van Rees, M. Gazzola, and P. Koumoutsakos, “Optimal shapes for anguilliform swimmers at intermediate Reynolds numbers,” J. Fluid Mech. vol. 722, 2013. Publisher's VersionAbstract
We investigate the optimal morphologies for fast and efficient anguilliform swimmers at intermediate Reynolds numbers, by combining an evolution strategy with three-dimensional viscous vortex methods. We show that anguilliform swimmer shapes enable the trapping and subsequent acceleration of regions of fluid transported along the entire body by the midline travelling wave. A sensitivity analysis of the optimal morphological traits identifies that the width thickness in the anterior of the body and the height of the caudal fin are critical factors for both speed and efficiency. The fastest swimmer without a caudal fin, however, still retains 80 % of its speed, showing that the entire body is used to generate thrust. The optimal shapes share several features with naturally occurring morphologies, but their overall appearances differ. This demonstrates that engineered swimmers can outperform biomimetic swimmers for the criteria considered here.
F. Milde, S. Lauw, P. Koumoutsakos, and M. L. Iruela-Arispe, “The mouse retina in 3D: quantification of vascular growth and remodeling,” Integrative Biology, vol. 5, no. 12, pp. 1426, 2013. Publisher's VersionAbstract
The mouse retina has become a prominent model for studying angiogenesis. The easy access and well-known developmental progression have significantly propelled our ability to examine and manipulate blood vessels in vivo. Nonetheless, most studies have restricted their evaluations to the superficial plexus (an upper vascular layer in contact with the vitreous). Here we present experimental data and quantification for the developmental progression of the full retina including the intermediate and deeper plexus that sprouts from the superficial layer. We analyze the origin and advancement of vertical sprouting and present the progression of vascular perfusion within the tissue. Furthermore, we introduce the use of Minkowsky functionals to quantify remodeling in the superficial and deeper plexus. The work expands information on the retina towards a 3D structure. This is of particular interest, as recent data have demonstrated differential effects of gene deletion on the upper and deeper plexus, highlighting the concept of distinct operational pathways during sprouting angiogenesis.
P. Angelikopoulos, C. Papadimitriou, and P. Koumoutsakos, “Data driven, predictive molecular dynamics for nanoscale flow simulations under uncertainty,” The Journal of Physical Chemistry B, vol. 117, no. 47, pp. 14808–14816, 2013. Publisher's VersionAbstract
For over five decades, molecular dynamics (MD) simulations have helped to elucidate critical mechanisms in a broad range of physiological systems and technological innovations. MD simulations are synergetic with experiments, relying on measurements to calibrate their parameters and probing {\textquotedblleft}what if scenarios{\textquotedblright} for systems that are difficult to investigate experimentally. However, in certain systems, such as nanofluidics, the results of experiments and MD simulations differ by several orders of magnitude. This discrepancy may be attributed to the spatiotemporal scales and structural information accessible by experiments and simulations. Furthermore, MD simulations rely on parameters that are often calibrated semiempirically, while the effects of their computational implementation on their predictive capabilities have only been sporadically probed. In this work, we show that experimental and MD investigations can be consolidated through a rigorous uncertainty quantification framework. We employ a Bayesian probabilistic framework for large scale MD simulations of graphitic nanostructures in aqueous environments. We assess the uncertainties in the MD predictions for quantities of interest regarding wetting behavior and hydrophobicity. We focus on three representative systems: water wetting of graphene, the aggregation of fullerenes in aqueous solution, and the water transport across carbon nanotubes. We demonstrate that the dominant mode of calibrating MD potentials in nanoscale fluid mechanics, through single values of water contact angle on graphene, leads to large uncertainties and fallible quantitative predictions. We demonstrate that the use of additional experimental data reduces uncertainty, improves the predictive accuracy of MD models, and consolidates the results of experiments and simulations.
B. Hejazialhosseini, D. Rossinelli, and P. Koumoutsakos, “3d shock-bubble interaction,” Physics of Fluids, vol. 25, no. 11, pp. 110816, 2013. Publisher's VersionAbstract
We present detailed visualizations of the interactions of a normal shock wave at Mach 3, with a spherical helium bubble immersed in air,1 with an interface Atwood number of −0.76 (Figure 1). The governing 3D Euler equations for two-phase compressible flows are solved using a finite volume solver with uniform resolution. We employ the 5th order WENO reconstruction of the primitive quantities, an HLL-type numerical flux, and the 3rd order TVD Runge-Kutta time stepping scheme. The software achieves 30% of the peak performance on a Cray XE6, using 4 × 109 cells. Extended simulations reveal that the shock passage compresses the bubble and generates baroclinic vorticity on the density interface. Initial distribution of the vorticity and compressions lead to the formation of an air jet, interface roll-ups, and the formation of a long lasting vortical core.
P. Koumoutsakos and J. Feigelman, “Multiscale stochastic simulations of chemical reactions with regulated scale separation,” Journal of Computational Physics, vol. 244, pp. 290–297, 2013. Publisher's VersionAbstract
We present a coupling of multiscale frameworks with accelerated stochastic simulation algorithms for systems of chemical reactions with disparate propensities. The algorithms regulate the propensities of the fast and slow reactions of the system, using alternating micro and macro sub-steps simulated with accelerated algorithms such as τ and R-leaping. The proposed algorithms are shown to provide significant speedups in simulations of stiff systems of chemical reactions with a trade-off in accuracy as controlled by a regulating parameter. More importantly, the error of the methods exhibits a cutoff phenomenon that allows for optimal parameter choices. Numerical experiments demonstrate that hybrid algorithms involving accelerated stochastic simulations can be, in certain cases, more accurate while faster, than their corresponding stochastic simulation algorithm counterparts.
J. H. Walther, K. Ritos, E. R. Cruz-Chu, C. M. Megaridis, and P. Koumoutsakos, “Barriers to superfast water transport in carbon nanotube membranes,” Nano Letters, vol. 13, no. 5, pp. 1910–1914, 2013. Publisher's VersionAbstract
Carbon nanotube (CNT) membranes hold the promise of extraordinary fast water transport for applications such as energy efficient filtration and molecular level drug delivery. However, experiments and computations have reported flow rate enhancements over continuum hydrodynamics that contradict each other by orders of magnitude. We perform large scale molecular dynamics simulations emulating for the first time the micrometer thick CNTs membranes used in experiments. We find transport enhancement rates that are length dependent due to entrance and exit losses but asymptote to 2 orders of magnitude over the continuum predictions. These rates are far below those reported experimentally. The results suggest that the reported superfast water transport rates cannot be attributed to interactions of water with pristine CNTs alone.
D. Franco, et al., “Accelerated endothelial wound healing on microstructured substrates under flow,” Biomaterials, vol. 34, no. 5, pp. 1488–1497, 2013. Publisher's VersionAbstract
Understanding and accelerating the mechanisms of endothelial wound healing is of fundamental interest for biotechnology and of significant medical utility in repairing pathologic changes to the vasculature induced by invasive medical interventions. We report the fundamental mechanisms that determine the influence of substrate topography and flow on the efficiency of endothelial regeneration. We exposed endothelial monolayers, grown on topographically engineered substrates (gratings), to controlled levels of flow-induced shear stress. The wound healing dynamics were recorded and analyzed in various configurations, defined by the relative orientation of an inflicted wound, the topography and the flow direction. Under flow perpendicular to the wound, the speed of endothelial regeneration was significantly increased on substrates with gratings oriented in the direction of the flow when compared to flat substrates. This behavior is linked to the dynamic state of cell-to-cell adhesions in the monolayer. In particular, interactions with the substrate topography counteract Vascular Endothelial Cadherin phosphorylation induced by the flow and the wounding. This effect contributes to modulating the mechanical connection between migrating cells to an optimal level, increasing their coordination and resulting in coherent cell motility and preservation of the monolayer integrity, thus accelerating wound healing. We further demonstrate that the reduction of vascular endothelial cadherin phosphorylation, through specific inhibition of Src activity, enhances endothelial wound healing in flows over flat substrates.
B. Hejazialhosseini, D. Rossinelli, and P. Koumoutsakos, “Vortex dynamics in 3d shock-bubble interaction,” Publisher, vol. 25, no. 11, pp. 110816, 2013. Publisher's VersionAbstract
The dynamics of shock-bubble interaction involve an interplay of vortex stretching, dilation, and baroclinic vorticity generation. Here, we quantify the interplay of these contributions through high resolution 3D simulations for several Mach and Atwood numbers. We present a volume rendering of density and vorticity magnitude fields of shock-bubble interaction at M = 3 and air/helium density ratio {η} = 7.25 to elucidate the evolution of the flow structures. We distinguish the vorticity growth rates due to baroclinicity, stretching, and dilatation at low and high Mach numbers as well as the late time evolution of the circulation. The results demonstrate that a number of analytical models need to be revised in order to predict the late time circulation of shock-bubble interactions at high Mach numbers. To this effect, we propose a simple model for the dependence of the circulation to Mach number and ambient to bubble density ratio for air/helium shock-bubble interactions.
P. Koumoutsakos, I. Pivkin, and F. Milde, “The fluid mechanics of cancer and its therapy,” Annual Review of Fluid Mechanics, vol. 45, no. 1, pp. 325–355, 2013. Publisher's VersionAbstract
Fluid mechanics is involved in the growth, progression, metastasis, and therapy of cancer. Blood vessels transport oxygen and nutrients to cancerous tissues, provide a route for metastasizing cancer cells to distant organs, and deliver drugs to tumors. The irregular and leaky tumor vasculature is responsible for increased interstitial pressure in the tumor microenvironment, whereas multiscale flow-structure interaction processes control tumor growth, metastasis, and nanoparticle-mediated drug delivery. We outline these flow-mediated processes, along with related experimental and computational methods for the diagnosis, predictive modeling, and therapy of cancer.
G. Tauriello and P. Koumoutsakos, “Coupling remeshed particle and phase field methods for the simulation of reaction-diffusion on the surface and the interior of deforming geometries,” SIAM Journal on Scientific Computing, vol. 35, no. 6, pp. B1285–B1303, 2013. Publisher's VersionAbstract
Reaction-diffusion systems on the surface and the interior of complex domains are potent models of growth in living organisms. The simulation of these models requires numerical methods capable of handling large deformations and the accurate coupling of the evolution of substances in the lumen and on the surface of the deforming geometries. Here, we develop a novel computational method to handle such problems by combining a remeshed particle method with a phase field method. Remeshed particle methods are well suited to discretizing deforming geometries, while the phase field method is used to impose boundary conditions that effectuate the coupling of substances evolving in their lumen and on their surfaces. We demonstrate that this hybrid method enables for the first time the accurate coupling of reaction-diffusion on a deformable surface and its interior. The method is validated on benchmark problems and the effect of lumen diffusion to a pattern forming reaction-diffusion system on a deforming surface is discussed.

2012

W. M. van Rees, F. Hussain, and P. Koumoutsakos, “Vortex tube reconnection at Re = 10^4,” Phys. Fluids, vol. 24, no. 7, pp. 075105, 2012. Publisher's VersionAbstract
We present simulations of the long-time dynamics of two anti-parallel vortex tubes with and without initial axial flow, at Reynolds number Re = \Gamma/ν = 104. Simulations were performed in a periodic domain with a remeshed vortex method using 785 × 106 particles. We quantify the vortex dynamics of the primary vortex reconnection that leads to the formation of elliptical rings with axial flow and report for the first time a subsequent collision of these rings. In the absence of initial axial flow, a -5/3 slope of the energy spectrum is observed during the first reconnection of the tubes. The resulting elliptical vortex rings experience a coiling of their vortex lines imparting an axial flow inside their cores. These rings eventually collide, exhibiting a -7/3 slope of the energy spectrum. Studies of vortex reconnection with an initial axial flow exhibit also the -7/3 slope during the initial collision as well as in the subsequent collision of the ensuing elliptical vortex rings. We quantify the detailed vortex dynamics of these collisions and examine the role of axial flow in the breakup of vortex structures.
M. Paolucci, et al., “Towards a living earth simulator,” Eur. Phys. J.-Spec. Top. vol. 214, no. 1, pp. 77–108, 2012. Publisher's VersionAbstract
The Living Earth Simulator (LES) is one of the core components of the FuturICT architecture. It will work as a federation of methods, tools, techniques and facilities supporting all of the FuturICT simulation-related activities to allow and encourage interactive exploration and understanding of societal issues. Society-relevant problems will be targeted by leaning on approaches based on complex systems theories and data science in tight interaction with the other components of FuturICT. The LES will evaluate and provide answers to real-world questions by taking into account multiple scenarios. It will build on present approaches such as agent-based simulation and modeling, multiscale modelling, statistical inference, and data mining, moving beyond disciplinary borders to achieve a new perspective on complex social systems.
F. Milde, D. Franco, A. Ferrari, V. Kurtcuoglu, D. Poulikakos, and P. Koumoutsakos, “Cell Image Velocimetry (CIV): boosting the automated quantification of cell migration in wound healing assays,” Integr. Biol. vol. 4, no. 11, pp. 1437, 2012. Publisher's VersionAbstract
Cell migration is commonly quantified by tracking the speed of the cell layer interface in wound healing assays. This quantification is often hampered by low signal to noise ratio, in particular when complex substrates are employed to emulate in vivo cell migration in geometrically complex environments. Moreover, information about the cell motion, readily available inside the migrating cell layers, is not usually harvested. We introduce Cell Image Velocimetry (CIV), a combination of cell layer segmentation and image velocimetry algorithms, to drastically enhance the quantification of cell migration by wound healing assays. The resulting software analyses the speed of the interface as well as the detailed velocity field inside the cell layers in an automated fashion. CIV is shown to be highly robust for images with low signal to noise ratio, low contrast and frame shifting and it is portable across various experimental settings. The modular design and parametrization of CIV is not restricted to wound healing assays and allows for the exploration and quantification of flow phenomena in any optical microscopy dataset. Here, we demonstrate the capabilities of CIV in wound healing assays over topographically engineered surfaces and quantify the relative merits of differently aligned gratings on cell migration.
J. H. Walther, M. Praprotnik, E. M. Kotsalis, and P. Koumoutsakos, “Multiscale simulation of water flow past a C540 fullerene,” J. Comput. Phys. vol. 231, no. 7, pp. 2677–2681, 2012. Publisher's VersionAbstract
We present a novel, three-dimensional, multiscale algorithm for simulations of water flow past a fullerene. We employ the Schwarz alternating overlapping domain method to couple molecular dynamics (MD) of liquid water around the C540 buckyball with a Lattice–Boltzmann (LB) description for the Navier–Stokes equations. The proposed method links the MD and LB domains using a fully three-dimensional interface and coupling of velocity gradients. The present overlapping domain method implicitly preserves the flux of mass and momentum and bridges flux-based and Schwarz domain decomposition algorithms. We use this method to determine the slip length and hydrodynamic radius for water flow past a buckyball.
M. Gazzola, C. Mimeau, A. A. Tchieu, and P. Koumoutsakos, “Flow mediated interactions between two cylinders at finite Re numbers,” Phys. Fluids, vol. 24, no. 4, pp. 043103, 2012. Publisher's VersionAbstract
We present simulations of two interacting moving cylinders immersed in a twodimensional incompressible, viscous flow. Simulations are performed by coupling a wavelet-adapted, remeshed vortex method with the Brinkman penalization and projection approach. This method is validated on benchmark problems and applied to simulations of a master-slave pair of cylinders. The master cylinder's motion is imposed and the slave cylinder is let free to respond to the flow.We study the relative role of viscous and inertia effects in the cylinders interactions and identify related sharp transitions in the response of the slave. The observed differences in the behavior of cylinders with respect to corresponding potential flow simulations are discussed. In addition, it is observed that in certain situations the finite size of the slave cylinders enhances the transport so that the cylinders are advected more effectively than passive tracers placed, respectively, at the same starting position.
M. Gazzola, W. M. van Rees, and P. Koumoutsakos, “C-start: optimal start of larval fish,” J. Fluid Mech. vol. 698, pp. 5–18, 2012. Publisher's VersionAbstract
We investigate the C-start escape response of larval fish by combining flow simulations using remeshed vortex methods with an evolutionary optimization. We test the hypothesis of the optimality of C-start of larval fish by simulations of larval-shaped, two- and three-dimensional self-propelled swimmers. We optimize for the distance travelled by the swimmer during its initial bout, bounding the shape deformation based on the larval mid-line curvature values observed experimentally. The best motions identified within these bounds are in good agreement with in vivo experiments and show that C-starts do indeed maximize escape distances. Furthermore we found that motions with curvatures beyond the ones experimentally observed for larval fish may result in even larger escape distances. We analyse the flow field and find that the effectiveness of the C-start escape relies on the ability of pronounced C-bent body configurations to trap and accelerate large volumes of fluid, which in turn correlates with large accelerations of the swimmer.
C. Conti, D. Rossinelli, and P. Koumoutsakos, “GPU and APU computations of Finite Time Lyapunov Exponent fields,” J. Comput. Phys. vol. 231, no. 5, pp. 2229–2244, 2012. Publisher's VersionAbstract
We present GPU and APU accelerated computations of Finite-Time Lyapunov Exponent (FTLE) fields. The calculation of FTLEs is a computationally intensive process, as in order to obtain the sharp ridges associated with the Lagrangian Coherent Structures an extensive resampling of the flow field is required. The computational performance of this resampling is limited by the memory bandwidth of the underlying computer architecture. The present technique harnesses data-parallel execution of many-core architectures and relies on fast and accurate evaluations of moment conserving functions for the mesh to particle interpolations. We demonstrate how the computation of FTLEs can be efficiently performed on a GPU and on an APU through OpenCL and we report over one order of magnitude improvements over multi-threaded executions in FTLE computations of bluff body flows.
P. Angelikopoulos, C. Papadimitriou, and P. Koumoutsakos, “Bayesian uncertainty quantification and propagation in molecular dynamics simulations: A high performance computing framework,” J. Chem. Phys. vol. 137, no. 14, pp. 144103, 2012. Publisher's VersionAbstract
ent a Bayesian probabilistic framework for quantifying and propagating the uncertainties in the parameters of force fields employed in molecular dynamics (MD) simulations. We propose a highly parallel implementation of the transitional Markov chain Monte Carlo for populating the posterior probability distribution of the MD force-field parameters. Efficient scheduling algorithms are proposed to handle the MD model runs and to distribute the computations in clusters with heterogeneous architectures. Furthermore, adaptive surrogate models are proposed in order to reduce the computational cost associated with the large number of MD model runs. The effectiveness and computational efficiency of the proposed Bayesian framework is demonstrated in MD simulations of liquid and gaseous argon.

2011

D. Rossinelli, B. Hejazialhosseini, D. G. Spampinato, and P. Koumoutsakos, “Multicore/Multi-GPU Accelerated Simulations of Multiphase Compressible Flows Using Wavelet Adapted Grids,” SIAM J. Sci. Comput. vol. 33, no. 2, pp. 512–540, 2011. Publisher's VersionAbstract
We present a computational method of coupling average interpolating wavelets with high-order finite volume schemes and its implementation on heterogeneous computer architectures for the simulation of multiphase compressible flows. The method is implemented to take advantage of the parallel computing capabilities of emerging heterogeneous multicore/multi-GPU architectures. A highly efficient parallel implementation is achieved by introducing the concept of wavelet blocks, exploiting the task-based parallelism for CPU cores, and by managing asynchronously an array of GPUs by means of OpenCL. We investigate the comparative accuracy of the GPU and CPU based simulations and analyze their discrepancy for two-dimensional simulations of shock-bubble interaction and Richtmeyer–Meshkov instability. The results indicate that the accuracy of the GPU/CPU heterogeneous solver is competitive with the one that uses exclusively the CPU cores. We report the performance improvements by employing up to 12 cores and 6 GPUs compared to the single-core execution. For the simulation of the shock-bubble interaction at Mach 3 with two million grid points, we observe a 100-fold speedup for the heterogeneous part and an overall speedup of 34.
D. Rossinelli, C. Conti, and P. Koumoutsakos, “Mesh-particle interpolations on graphics processing units and multicore central processing units,” Philos. T. Roy. Soc. A, vol. 369, no. 1944, pp. 2164–2175, 2011. Publisher's VersionAbstract
Particle–mesh interpolations are fundamental operations for particle-in-cell codes, as implemented in vortex methods, plasma dynamics and electrostatics simulations. In these simulations, the mesh is used to solve the field equations and the gradients of the fields are used in order to advance the particles. The time integration of particle trajectories is performed through an extensive resampling of the flow field at the particle locations. The computational performance of this resampling turns out to be limited by the memory bandwidth of the underlying computer architecture. We investigate how mesh–particle interpolation can be efficiently performed on graphics processing units (GPUs) and multicore central processing units (CPUs), and we present two implementation techniques. The single-precision results for the multicore CPU implementation show an acceleration of 45–70×, depending on system size, and an acceleration of 85–155× for the GPU implementation over an efficient single-threaded C++ implementation. In double precision, we observe a performance improvement of 30–40× for the multicore CPU implementation and 20–45× for the GPU implementation. With respect to the 16-threaded standard C++ implementation, the present CPU technique leads to a performance increase of roughly 2.8–3.7× in single precision and 1.7–2.4× in double precision, whereas the GPU technique leads to an improvement of 9× in single precision and 2.2–2.8× in double precision.
W. M. van Rees, A. Leonard, D. I. Pullin, and P. Koumoutsakos, “A comparison of vortex and pseudo-spectral methods for the simulation of periodic vortical flows at high Reynolds numbers,” J. Comput. Phys. vol. 230, no. 8, pp. 2794–2805, 2011. Publisher's VersionAbstract
We present a validation study for the hybrid particle-mesh vortex method against a pseudo-spectral method for the Taylor–Green vortex at Re=1600 as well as in the collision of two antiparallel vortex tubes at Re=10,000. In this study we present diagnostics such as energy spectra and enstrophy as computed by both methods as well as point-wise comparisons of the vorticity field. Using a fourth order accurate kernel for interpolation between the particles and the mesh, the results of the hybrid vortex method and of the pseudo-spectral method agree well in both flow cases. For the Taylor–Green vortex, the vorticity contours computed by both methods around the time of the energy dissipation peak overlap. The energy spectrum shows that only the smallest length scales in the flow are not captured by the vortex method. In the second flow case, where we compute the collision of two anti-parallel vortex tubes at Reynolds number 10,000, the vortex method results and the pseudo-spectral method results are in very good agreement up to and including the first reconnection of the tubes. The maximum error in the effective viscosity is about 2.5% for the vortex method and about 1% for the pseudo-spectral method. At later times the flows computed with the different methods show the same qualitative features, but the quantitative agreement on vortical structures is lost.
P. D. Koumoutsakos, B. Bayati, F. Milde, and G. Tauriello, “Particle Simulations of Morphogenesis,” Math. Mod. Meth. Appl. S. vol. 21, no. supp01, pp. 955–1006, 2011. Publisher's VersionAbstract
The simulation of the creation and evolution of biological forms requires the development of computational methods that are capable of resolving their hierarchical, spatial and temporal complexity. Computations based on interacting particles, provide a unique computational tool for discrete and continuous descriptions of morphogenesis of systems ranging from the molecular to the organismal level. The capabilities of particle methods hinge on the simplicity of their formulation which enables the formulation of a unifying computational framework encom- passing deterministic and stochastic models. In this paper, we discuss recent advances in particle methods for the simulation of biological systems at the mesoscopic and the macroscale level. We present results from applications of particle methods including reaction-diffusion on deforming surfaces, deterministic and stochastic descriptions of tumor growth and angiogenesis and discuss successes and challenges of this approach.
M. Gazzola, P. Chatelain, W. M. van Rees, and P. Koumoutsakos, “Simulations of single and multiple swimmers with non-divergence free deforming geometries,” J. Comput. Phys. vol. 230, no. 19, pp. 7093–7114, 2011. Publisher's VersionAbstract
We present a vortex particle method coupled with a penalization technique to simulate single and multiple swimmers in an incompressible, viscous flow in two and three dimensions. The proposed algorithm can handle arbitrarily deforming bodies and their corresponding non-divergence free deformation velocity fields. The method is validated on a number of benchmark problems with stationary and moving boundaries. Results include flows of tumbling objects and single and multiple self-propelled swimmers.
M. Gazzola, O. V. Vasilyev, and P. Koumoutsakos, “Shape optimization for drag reduction in linked bodies using evolution strategies,” Comput. Struct. vol. 89, no. 11-12, pp. 1224–1231, 2011. Publisher's VersionAbstract
We present results from the shape optimization of linked bodies for drag reduction in simulations of incompressible flow at moderate Reynolds numbers. The optimization relies on the covariance matrix adaptation evolution strategy (CMA-ES) and the flow simulations use vortex methods with the Brinkman penalization to enforce boundary conditions in complex bodies. We exploit the inherent parallelism of CMA-ES, by implementing a multi-host framework which allows for the distribution of the expensive cost function evaluations across parallel architectures, without being limited to one computing facility. This study repeats in silico for the first time Ingo Rechenberg's pioneering wind tunnel experiments for drag reduction that led to the inception of evolution strategies. The simulations confirm that the results of these experimental studies indicate a flat plate is not the optimal solution for drag reduction in linked bodies. We present the vorticity field of the flow and identify the governing mechanisms for this drag reduction by the slightly corrugated linked plate configuration.
B. Bayati, P. Chatelain, and P. Koumoutsakos, “Adaptive mesh refinement for stochastic reaction–diffusion processes,” J. Comput. Phys. vol. 230, no. 1, pp. 13–26, 2011. Publisher's VersionAbstract
We present an algorithm for Adaptive Mesh Refinement applied to mesoscopic stochastic simulations of spatially evolving reaction-diffusion processes. The transition rates for the diffusion process are derived on adaptive, locally refined structured meshes. Convergence of the diffusion process is presented and the fluctuations of the stochastic process are verified. Furthermore, a refinement criterion is proposed for the evolution of the adaptive mesh. The method is validated in simulations of reaction-diffusion processes as described by the Fisher-Kolmogorov and Gray-Scott equations.
G. Schwank, G. Tauriello, R. Yagi, E. Kranz, P. Koumoutsakos, and K. Basler, “Antagonistic Growth Regulation by Dpp and Fat Drives Uniform Cell Proliferation,” Dev. Cell, vol. 20, no. 1, pp. 123–130, 2011. Publisher's VersionAbstract
We use the Dpp morphogen gradient in the Drosophila wing disc as a model to address the fundamental question of how a gradient of a growth factor can produce uniform growth. We first show that proper expression and subcellular localization of components in the Fat tumor-suppressor pathway, which have been argued to depend on Dpp activity differences, are not reliant on the Dpp gradient. We next analyzed cell proliferation in discs with uniformly high Dpp or uniformly low Fat signaling activity and found that these pathways regulate growth in a complementary manner. While the Dpp mediator Brinker inhibits growth in the primordium primarily in the lateral regions, Fat represses growth mostly in the medial region. Together, our results indicate that the activities of both signaling pathways are regulated in a parallel rather than sequential manner and that uniform proliferation is achieved by their complementary action on growth.

2010

D. Rossinelli, B. Hejazialhosseini, M. Bergdorf, and P. Koumoutsakos, “Wavelet-adaptive solvers on multi-core architectures for the simulation of complex systems,” Concurr. Comp.-Prat. E. vol. 23, no. 2, pp. 172–186, 2010. Publisher's VersionAbstract
We build wavelet-based adaptive numerical methods for the simulation of advection-dominated flows that develop multiple spatial scales, with an emphasis on fluid mechanics problems. Wavelet-based adaptivity is inherently sequential and in this work we demonstrate that these numerical methods can be implemented in software that is capable of harnessing the capabilities of multi-core architectures while maintaining their computational efficiency. Recent designs in frameworks for multi-core software development allow us to rethink parallelism as task-based, where parallel tasks are specified and automatically mapped onto physical threads. This way of exposing parallelism enables the parallelization of algorithms that were considered inherently sequential, such as wavelet-based adaptive simulations. In this paper we present a framework that combines wavelet-based adaptivity with the task-based parallelism. We demonstrate the promising performance obtained by simulating various physical systems on different multi-core architectures using up to 16 cores.
D. Rossinelli, M. Bergdorf, G. - H. Cottet, and P. Koumoutsakos, “GPU accelerated simulations of bluff body flows using vortex particle methods,” J. Comput. Phys. vol. 229, no. 9, pp. 3316–3333, 2010. Publisher's VersionAbstract
We present a GPU accelerated solver for simulations of bluff body flows in 2D using a remeshed vortex particle method and the vorticity formulation of the Brinkman penalization technique to enforce boundary conditions. The efficiency of the method relies on fast and accurate particle-grid interpolations on GPUs for the remeshing of the particles and the computation of the field operators. The GPU implementation uses OpenGL so as to perform efficient particle-grid operations and a CUFFT-based solver for the Poisson Equation with unbounded boundary conditions. The accuracy and performance of the GPU simulations and their relative advantages/drawbacks over CPU based computations are reported in simulations of flows past an impulsively started circular cylinder from Reynolds numbers between 40 and 9,500. The results indicate up to two orders of magnitude speed up of the GPU implementation over the respective CPU implementations. The accuracy of the GPU computations depends on the Re number of the flow. For Re up to 1000 there is little difference between GPU and CPU calculations but this agreement deteriorates (albeit remaining to within 5% in drag calculations) for higher Re numbers as the single precision of the GPU adversely affects the accuracy of the simulations.
G. Morra, D. A. Yuen, L. Boschi, P. Chatelain, P. Koumoutsakos, and P. J. Tackley, “The fate of the slabs interacting with a density/viscosity hill in the mid-mantle,” Phys. Earth Planet. In. vol. 180, no. 3-4, pp. 271–282, 2010. Publisher's VersionAbstract
In the last two decades it has been proposed several times that a non-monotonic profile might fit the average lower mantle radial viscosity. Most proposed profiles consist in a more or less broad viscosity hill in the middle of the mantle, at a depth roughly between 1200 km and 2000 km. Also many tomographic models display strong signals of the presence of “fast” material lying at mid mantle depths and a recent spectral analysis of seismic tomography shows a very clear transition for degree up to around 16 at a less than 1500 km depth. Finally latest works, both theoretical and experimental, on the high-to-low-spin transition for periclase, have suggested that the high-spin to low-spin transition of Fe++ might lie at the heart of all these observations. To verify the dynamical compatibility between possible mantle profile and observed tomographic images and compare them with possible mineral physics scenarios, such as the spin transition, we employ here a recently developed Fast Multipole-accelerated Boundary Element Method (FMM-BEM), a numerical approach for solving the viscous momentum equation in a global spherical setting, for simulating the interaction of an individual slab with a mid mantle smooth discontinuity in density and viscosity. We have focused on the complexities induced to the behaviour of average and very large plates O (2000–10,000 km), characteristic of the Farallon, Tethys and Pacific plate subducting during the Cenozoic, demonstrating that the a mid mantle density and/or viscosity discontinuity produces a strong alteration of the sinking velocity and an intricate set of slab morphologies. We also employ the Kula–Farallon plate system subducting at 60 Ma as a paradigmatic case, which reveals the best high resolution tomography models and clearly suggests an interaction with a strong and/or denser layer in the mantle. Our 38 models show that a plate might or might not penetrate into the lowest mantle and might stall in the mid lower mantle for long periods, depending on the radial profiles of density and viscosity, within a realistic range (viscosity 1, 10 or 100 times more viscous of the rest of the mantle, and a change of differential density in the range -2% to 2%), of a transitional layer of 200 km or 500 km. We conclude that a layer with high viscosity or negative density would naturally trigger the observed geodynamic snapshot. We finally propose a scenario in which the long time accumulation of depleted slabs in the mid mantle would give rise to a partially chemically stratified mantle, starting from the less prominent high-spin to low-spin contribution on the basis of mantle density and rheology.
B. Hejazialhosseini, D. Rossinelli, M. Bergdorf, and P. Koumoutsakos, “High order finite volume methods on wavelet-adapted grids with local time-stepping on multicore architectures for the simulation of shock-bubble interactions,” J. Comput. Phys. vol. 229, no. 22, pp. 8364–8383, 2010. Publisher's VersionAbstract
We present a space–time adaptive solver for single- and multi-phase compressible flows that couples average interpolating wavelets with high-order finite volume schemes. The solver introduces the concept of wavelet blocks, handles large jumps in resolution and employs local time-stepping for efficient time integration. We demonstrate that the inherently sequential wavelet-based adaptivity can be implemented efficiently in multicore computer architectures using task-based parallelism and introducing the concept of wavelet blocks. We validate our computational method on a number of benchmark problems and we present simulations of shock-bubble interaction at different Mach numbers, demonstrating the accuracy and computational performance of the method.
I. Hanasaki, J. H. Walther, S. Kawano, and P. Koumoutsakos, “Coarse-grained molecular dynamics simulations of shear-induced instabilities of lipid bilayer membranes in water,” Phys. Rev. E, vol. 82, no. 5, 2010. Publisher's VersionAbstract
We study shear-induced instabilities of lipid bilayers immersed in water using coarse-grained molecular dynamics simulations. The shear imposed by the flow of the water induces initially microscopic structural changes of the membrane, starting with tilting of the molecules in the direction of the shear. The tilting propagates in the spanwise direction when the shear rate exceeds a critical value and the membrane undergoes a bucklinglike deformation in the direction perpendicular to the shear. The bucklinglike undulation continues until a localized Kelvin-Helmholtz-like instability leads to membrane rupture. We study the different modes of membrane undulation using membranes of different geometries and quantify the relative importance of the bucklinglike bending and the Kelvin-Helmholtz-like instability of the membrane.
P. Chatelain and P. Koumoutsakos, “A Fourier-based elliptic solver for vortical flows with periodic and unbounded directions,” J. Comput. Phys. vol. 229, no. 7, pp. 2425–2431, 2010. Publisher's VersionAbstract
We present a computationally efficient, adaptive solver for the solution of the Poisson and Helmholtz equation used in flow simulations in domains with combinations of unbounded and periodic directions. The method relies on using FFTs on an extended domain and it is based on the method proposed by Hockney and Eastwood for plasma simulations. The method is well-suited to problems with dynamically growing domains and in particular flow simulations using vortex particle methods. The efficiency of the method is demonstrated in simulations of trailing vortices.
B. Bayati, H. Owhadi, and P. Koumoutsakos, “A cutoff phenomenon in accelerated stochastic simulations of chemical kinetics via flow averaging (FLAVOR-SSA),” J. Chem. Phys. vol. 133, no. 24, pp. 244117, 2010. Publisher's VersionAbstract
We present a simple algorithm for the simulation of stiff, discrete-space, continuous-time Markov processes. The algorithm is based on the concept of flow averaging for the integration of stiff ordinary and stochastic differential equations and ultimately leads to a straightforward variation of the well-known stochastic simulation algorithm (SSA). The speedup that can be achieved by the present algorithm [flow averaging integrator SSA (FLAVOR-SSA)] over the classical SSA comes naturally at the expense of its accuracy. The error of the proposed method exhibits a cutoff phenomenon as a function of its speed-up, allowing for optimal tuning. Two numerical examples from chemical kinetics are provided to illustrate the efficiency of the method.

2009

P. Gonnet, J. H. Walther, and P. Koumoutsakos, “vartheta-SHAKE: An extension to SHAKE for the explicit treatment of angular constraints,” Comput. Phys. Commun. vol. 180, no. 3, pp. 360–364, 2009. Publisher's VersionAbstract
This paper presents \textvartheta-SHAKE, an extension to SHAKE, an algorithm for the resolution of holonomic constraints in molecular dynamics simulations, which allows for the explicit treatment of angular constraints. We show that this treatment is more efficient than the use of fictitious bonds, significantly reducing the overlap between the individual constraints and thus accelerating convergence. The new algorithm is compared with SHAKE, M-SHAKE, the matrix-based approach described by Ciccotti and Ryckaert and P-SHAKE for rigid water and octane.
T. Gebäck, M. M. Schulz, P. Koumoutsakos, and M. Detmar, “TScratch: a novel and simple software tool for automated analysis of monolayer wound healing assays.Biotechniques, vol. 46, no. 4, pp. 265–274, 2009. Publisher's VersionAbstract
Cell migration plays a major role in development, physiology, and disease, and is frequently evaluated in vitro by the monolayer wound healing assay. The assay analysis, however, is a time-consuming task that is often performed manually. In order to accelerate this analysis, we have developed TScratch, a new, freely available image analysis technique and associated software tool that uses the fast discrete curvelet transform to automate the measurement of the area occupied by cells in the images. This tool helps to significantly reduce the time needed for analysis and enables objective and reproducible quantification of assays. The software also offers a graphical user interface which allows easy inspection of analysis results and, if desired, manual modification of analysis parameters. The automated analysis was validated by comparing its results with manual-analysis results for a range of different cell lines. The comparisons demonstrate a close agreement for the vast majority of images that were examined and indicate that the present computational tool can reproduce statistically significant results in experiments with well-known cell migration inhibitors and enhancers.
T. Gebäck and P. Koumoutsakos, “Edge detection in microscopy images using curvelets,” BMC Bioinformatics, vol. 10, no. 1, pp. 75, 2009. Publisher's VersionAbstract
Despite significant progress in imaging technologies, the efficient detection of edges and elongated features in images of intracellular and multicellular structures acquired using light or electron microscopy is a challenging and time consuming task in many laboratories. We present a novel method, based on the discrete curvelet transform, to extract a directional field from the image that indicates the location and direction of the edges. This directional field is then processed using the non-maximal suppression and thresholding steps of the Canny algorithm to trace along the edges and mark them. Optionally, the edges may then be extended along the directions given by the curvelets to provide a more connected edge map. We compare our scheme to the Canny edge detector and an edge detector based on Gabor filters, and show that our scheme performs better in detecting larger, elongated structures possibly composed of several step or ridge edges. The proposed curvelet based edge detection is a novel and competitive approach for imaging problems. We expect that the methodology and the accompanying software will facilitate and improve edge detection in images available using light or electron microscopy.
M. Gazzola, C. J. Burckhardt, B. Bayati, M. Engelke, U. F. Greber, and P. Koumoutsakos, “A Stochastic Model for Microtubule Motors Describes the In Vivo Cytoplasmic Transport of Human Adenovirus,” PLoS Comput. Biol. vol. 5, no. 12, pp. e1000623, 2009. Publisher's VersionAbstract
Cytoplasmic transport of organelles, nucleic acids and proteins on microtubules is usually bidirectional with dynein and kinesin motors mediating the delivery of cargoes in the cytoplasm. Here we combine live cell microscopy, single virus tracking and trajectory segmentation to systematically identify the parameters of a stochastic computational model of cargo transport by molecular motors on microtubules. The model parameters are identified using an evolutionary optimization algorithm to minimize the Kullback-Leibler divergence between the in silico and the in vivo run length and velocity distributions of the viruses on microtubules. The present stochastic model suggests that bidirectional transport of human adenoviruses can be explained without explicit motor coordination. The model enables the prediction of the number of motors active on the viral cargo during microtubule-dependent motions as well as the number of motor binding sites, with the protein hexon as the binding site for the motors.
B. L. Falcón, et al., “Contrasting Actions of Selective Inhibitors of Angiopoietin-1 and Angiopoietin-2 on the Normalization of Tumor Blood Vessels,” Am. J. Pathol. vol. 175, no. 5, pp. 2159–2170, 2009. Publisher's VersionAbstract
Angiopoietin-1 (Ang1) and angiopoietin-2 (Ang2) have complex actions in angiogenesis and vascular remodeling due to their effects on Tie2 receptor signaling. Ang2 blocks Ang1-mediated activation of Tie2 in endothelial cells under certain conditions but is a Tie2 receptor agonist in others. We examined the effects of selective inhibitors of Ang1 (mL4-3) or Ang2 (L1-7[N]), alone or in combination, on the vasculature of human Colo205 tumors in mice. The Ang2 inhibitor decreased the overall abundance of tumor blood vessels by reducing tumor growth and keeping vascular density constant. After inhibition of Ang2, tumor vessels had many features of normal blood vessels (normalization), as evidenced by junctional accumulation of vascular endothelial-cadherin, junctional adhesion molecule-A, and platelet/endothelial cell adhesion molecule-1 in endothelial cells, increased pericyte coverage, reduced endothelial sprouting, and remodeling into smaller, more uniform vessels. The Ang1 inhibitor by itself had little noticeable effect on the tumor vasculature. However, when administered with the Ang2 inhibitor, the Ang1 inhibitor prevented tumor vessel normalization, but not the reduction in tumor vascularity produced by the Ang2 inhibitor. These findings are consistent with a model whereby inhibition of Ang2 leads to normalization of tumor blood vessels by permitting the unopposed action of Ang1, but decreases tumor vascularity primarily by blocking Ang2 actions.
H. A. Zambrano, J. H. Walther, P. Koumoutsakos, and I. F. Sbalzarini, “Thermophoretic Motion of Water Nanodroplets Confined inside Carbon Nanotubes,” Nano Lett. vol. 9, no. 1, pp. 66–71, 2009. Publisher's VersionAbstract
We study the thermophoretic motion of water nanodroplets confined inside carbon nanotubes using molecular dynamics simulations. We find that the nanodroplets move in the direction opposite the imposed thermal gradient with a terminal velocity that is linearly proportional to the gradient. The translational motion is associated with a solid body rotation of the water nanodroplet coinciding with the helical symmetry of the carbon nanotube. The thermal diffusion displays a weak dependence on the wetting of the water-carbon nanotube interface. We introduce the use of the moment scaling spectrum (MSS) in order to determine the characteristics of the motion of the nanoparticles inside the carbon nanotube. The MSS indicates that affinity of the nanodroplet with the walls of the carbon nanotubes is important for the isothermal diffusion and hence for the Soret coefficient of the system.
M. Bergdorf, I. F. Sbalzarini, and P. Koumoutsakos, “A Lagrangian particle method for reaction–diffusion systems on deforming surfaces,” J. Math Biol. vol. 61, no. 5, pp. 649–663, 2009. Publisher's VersionAbstract
Reaction-diffusion processes on complex deforming surfaces are fundamental to a number of biological processes ranging from embryonic development to cancer tumor growth and angiogenesis. The simulation of these processes using continuum reaction-diffusion models requires computational methods capable of accurately tracking the geometric deformations and discretizing on them the governing equations. We employ a Lagrangian level-set formulation to capture the deformation of the geometry and use an embedding formulation and an adaptive particle method to discretize both the level-set equations and the corresponding reaction-diffusion. We validate the proposed method and discuss its advantages and drawbacks through simulations of reaction-diffusion equations on complex and deforming geometries.
E. Mjolsness, D. Orendorff, P. Chatelain, and P. Koumoutsakos, “An exact accelerated stochastic simulation algorithm,” J. Chem. Phys. vol. 130, no. 14, pp. 144110, 2009. Publisher's VersionAbstract
An exact method for stochastic simulation of chemical reaction networks, which accelerates the stochastic simulation algorithm SSA , is proposed. The present “ER-leap” algorithm is derived from analytic upper and lower bounds on the multireaction probabilities sampled by SSA, together with rejection sampling and an adaptive multiplicity for reactions. The algorithm is tested on a number of well-quantified reaction networks and is found experimentally to be very accurate on test problems including a chaotic reaction network. At the same time ER-leap offers a substantial speedup over SSA with a simulation time proportional to the 2/3 power of the number of reaction events in a Galton–Watson process.
B. Bayati, P. Chatelain, and P. Koumoutsakos, “D-leaping: Accelerating stochastic simulation algorithms for reactions with delays,” J. Comput. Phys. vol. 228, no. 16, pp. 5908–5916, 2009. Publisher's VersionAbstract
We propose a novel, accelerated algorithm for the approximate stochastic simulation of biochemical systems with delays. The present work extends existing accelerated algorithms by distributing, in a time adaptive fashion, the delayed reactions so as to minimize the computational effort while preserving their accuracy. The accuracy of the present algorithm is assessed by comparing its results to those of the corresponding delay differential equations for a representative biochemical system. In addition, the fluctuations produced from the present algorithm are comparable to those from an exact stochastic simulation with delays. The algorithm is used to simulate biochemical systems that model oscillatory gene expression. The results indicate that the present algorithm is competitive with existing works for several benchmark problems while it is orders of magnitude faster for certain systems of biochemical reactions.
E. M. Kotsalis, J. H. Walther, E. Kaxiras, and P. Koumoutsakos, “Control algorithm for multiscale flow simulations of water,” Phys. Rev. E, vol. 79, no. 4, 2009. Publisher's VersionAbstract
We present a multiscale algorithm to couple atomistic water models with continuum incompressible flow simulations via a Schwarz domain decomposition approach. The coupling introduces an inhomogeneity in the description of the atomistic domain and prevents the use of periodic boundary conditions. The use of a mass conserving specular wall results in turn to spurious oscillations in the density profile of the atomistic description of water. These oscillations can be eliminated by using an external boundary force that effectively accounts for the virial component of the pressure. In this Rapid Communication, we extend a control algorithm, previously introduced for monatomic molecules, to the case of atomistic water and demonstrate the effectiveness of this approach. The proposed computational method is validated for the cases of equilibrium and Couette flow of water.
N. Hansen, A. S. P. Niederberger, L. Guzzella, and P. Koumoutsakos, “A Method for Handling Uncertainty in Evolutionary Optimization With an Application to Feedback Control of Combustion,” IEEE T. Evolut. Comput. vol. 13, no. 1, pp. 180–197, 2009. Publisher's VersionAbstract
We present a novel method for handling uncertainty in evolutionary optimization. The method entails quantification and treatment of uncertainty and relies on the rank based selection operator of evolutionary algorithms. The proposed uncertainty handling is implemented in the context of the covariance matrix adaptation evolution strategy (CMA-ES) and verified on test functions. The present method is independent of the uncertainty distribution, prevents premature convergence of the evolution strategy and is well suited for online optimization as it requires only a small number of additional function evaluations. The algorithm is applied in an experimental setup to the online optimization of feedback controllers of thermoacoustic instabilities of gas turbine combustors. In order to mitigate these instabilities, gain-delay or model-based controllers sense the pressure and command secondary fuel injectors. The parameters of these controllers are usually specified via a trial and error procedure. We demonstrate that their online optimization with the proposed methodology enhances, in an automated fashion, the online performance of the controllers, even under highly unsteady operating conditions, and it also compensates for uncertainties in the model-building and design process.

2009

F. Milde, M. Bergdorf, and P. Koumoutsakos, “A Hybrid Model of Sprouting Angiogenesis,” in Computational Science – ICCS 2008, Springer, 2008, pp. 167–176. Publisher's VersionAbstract

We present a computational model of tumor induced sprouting angiogenesis that involves a novel coupling of particle-continuum descriptions. The present 3D model of sprouting angiogenesis accounts for the effect of the extracellular matrix on capillary growth and considers both soluble and matrix-bound growth factors. The results of the simulations emphasize the role of the extracellular matrix and the different VEGF isoforms on branching behavior and the morphology of generated vascular networks.

F. Milde, M. Bergdorf, and P. Koumoutsakos, “A Hybrid Model for Three-Dimensional Simulations of Sprouting Angiogenesis,” Biophys. J. vol. 95, no. 7, pp. 3146–3160, 2008. Publisher's VersionAbstract
Recent advances in cancer research have identified critical angiogenic signaling pathways and the influence of the extracellular matrix on endothelial cell migration. These findings provide us with insight into the process of angiogenesis that can facilitate the development of effective computational models of sprouting angiogenesis. In this work, we present the first 3D model of sprouting angiogenesis that consider explicitly the effect of the extracellular matrix and of the soluble as well as matrix bound growth factors on capillary growth. The computational model relies on a hybrid particle-mesh representation of the blood vessels and it introduces an implicit representation of the vasculature that can accommodate detailed descriptions of nutrient transport. Extensive parametric studies reveal the role of the extracellular matrix structure and the distribution of the different VEGF isoforms on the dynamics and the morphology of the generated vascular networks.
S. E. Hieber and P. Koumoutsakos, “A Lagrangian particle method for the simulation of linear and nonlinear elastic models of soft tissue,” J. Comput. Phys. vol. 227, no. 21, pp. 9195–9215, 2008. Publisher's VersionAbstract
We present a novel Lagrangian particle method for the simulation of linear and nonlinear elastic models of soft tissue.Linear solids are represented by the Lagrangian formulation of the stress–strain relationship that is extended to nonlinear solids by using the Lagrangian evolution of the deformation gradient described in a moving framework. The present method introduces a level set description, along with the particles, to capture the body deformations and to enforce theboundary conditions. Furthermore, the accuracy of the method in cases of large deformations is ensured by implementinga particle remeshing procedure. The method is validated in several benchmark problems, in two and three dimensions andthe results compare well with the results of respective finite elements simulations. In simulations of large solid deformationunder plane strain compression, the finite element solver exhibits spurious structures that are not present in the Lagrangianparticle simulations. The particle simulations are compared with experimental results in an aspiration test of liver tissue.
S. E. Hieber and P. Koumoutsakos, “An immersed boundary method for smoothed particle hydrodynamics of self-propelled swimmers,” J. Comput. Phys. vol. 227, no. 19, pp. 8636–8654, 2008. Publisher's VersionAbstract
We present a novel particle method, combining remeshed Smoothed Particle Hydrodynamics with Immersed Boundary and Level Set techniques for the simulation of flows past complex deforming geometries. The present method retains the Lagrangian adaptivity of particle methods and relies on the remeshing of particle locations in order to ensure the accuracy of the method. In fact this remeshing step enables the introduction of Immersed Boundary Techniques used in grid based methods. The method is applied to simulations of flows of isothermal and compressible fluids past steady and unsteady solid boundaries that are described using a particle Level Set formulation. The method is validated with two and three-dimensional benchmark problems of flows past cylinders and spheres and it is shown to be well suited to simulations of large scale simulations using tens of millions of particles, on flow-structure interaction problems as they pertain to self-propelled anguilliform swimmers.
K. Fukagata, S. Kern, P. Chatelain, P. Koumoutsakos, and N. Kasagi, “Evolutionary optimization of an anisotropic compliant surface for turbulent friction drag reduction,” J. Turbul. vol. 9, 2008. Publisher's VersionAbstract
Direct numerical simulation (DNS) of the channel flow with an anisotropic compliant surface is performed in order to investigate its drag reduction effect in a fully developed turbulent flow. The computational domain is set to be 3δ × 2δ × 3δ , where δ is the channel half-width. The surface is passively driven by the pressure and wall-shear stress fluctuations, and the surface velocity provides a boundary condition for the fluid velocity field. An evolutionary optimization method (CMA-ES) is used to optimize the parameters of the anisotropic compliant surface. The optimization identifies several sets of parameters that result in a reduction of the friction drag with a maximum reduction rate of 8%. The primary mechanism for drag reduction is attributed to the decrease of the Reynolds shear stress (RSS) near the wall induced by the kinematics of the surface. The resultant wall motion is a uniform wave traveling downstream. The compliant wall, with the parameters found in the optimization study, is also tested in a computational domain that is doubled in the streamwise direction. The drag, however, is found to increase in the larger computational domain due to excessively large wall-normal velocity fluctuations.
U. Zimmerli and P. Koumoutsakos, “Simulations of Electrophoretic RNA Transport Through Transmembrane Carbon Nanotubes,” Biophys. J. vol. 94, no. 7, pp. 2546–2557, 2008. Publisher's VersionAbstract
The study of interactions between carbon nanotubes and cellular components, such as membranes and biomolecules, is fundamental for the rational design of nanodevices interfacing with biological systems. In this work, we use molecular dynamics simulations to study the electrophoretic transport of RNA through carbon nanotubes, embedded in membranes. Decorated and naked carbon nanotubes are inserted into a dodecane membrane and a dimyristoylphosphatidylcholine lipid bilayer, and the system is subjected to electrostatic potential differences. The transport properties of this artificial pore are determined by the structural modifications of the membrane in the vicinity of the nanotube openings and they are quantified by the nonuniform electrostatic potential maps at the entrance and inside the nanotube. The pore is used to transport electrophoretically a short RNA segment and we find that the speed of translocation exhibits an exponential dependence on the applied potential differences. The RNA is transported while undergoing a repeated stacking and unstacking process, affected by steric interactions with the membrane headgroups and by hydrophobic interaction with the walls of the nanotube. The RNA is structurally reorganized inside the nanotube, with its backbone solvated by water molecules near the axis of the tube and its bases aligned with the nanotube walls. Upon exiting the pore, the RNA interacts with the membrane headgroups and remains attached to the dodecane membrane while it is expelled into the solvent in the case of the lipid bilayer. The results of the simulations detail processes of molecular transport into cellular compartments through manufactured nanopores and they are discussed in the context of applications in biotechnology and nanomedicine.
A. Dupuis, P. Chatelain, and P. Koumoutsakos, “An immersed boundary–lattice-Boltzmann method for the simulation of the flow past an impulsively started cylinder,” J. Comput. Phys. vol. 227, no. 9, pp. 4486–4498, 2008. Publisher's VersionAbstract
We present a lattice-Boltzmann method coupled with an immersed boundary technique for the simulation of bluff body lows. The lattice-Boltzmann method for the modeling of the Navier–Stokes equations, is enhanced by a forcing term to ccount for the no-slip boundary condition on a non-grid conforming boundary. We investigate two alternatives of coupling he boundary forcing term with the grid nodes, namely the direct and the interpolated forcing techniques. The present B–IB methods are validated in simulations of the incompressible flow past an impulsively started cylinder at low and oderate Reynolds numbers. We present diagnostics such as the near wall vorticity field and the drag coefficient and comparisons ith previous computational and experimental works and assess the advantages and drawbacks of the two techniques. 2008 Elsevier Inc. All rights reserved.
D. Rossinelli, B. Bayati, and P. Koumoutsakos, “Accelerated stochastic and hybrid methods for spatial simulations of reaction–diffusion systems,” Chemical Phys. Lett. vol. 451, no. 1-3, pp. 136–140, 2008. Publisher's VersionAbstract
Spatial distributions characterize the evolution of reaction-diffusion models of several physical, chemical, and biological systems. We present two novel algorithms for the efficient simulation of these models: Spatial tau-Leaping (S tau-Leaping), employing a unified acceleration of the stochastic simulation of reaction and diffusion, and Hybrid tau-Leaping (H tau-Leaping), combining a deterministic diffusion approximation with a tau-Leaping acceleration of the stochastic reactions. The algorithms are validated by solving Fisher's equation and used to explore the role of the number of particles in pattern formation. The results indicate that the present algorithms have a nearly constant time complexity with respect to the number of events (reaction and diffusion), unlike the exact stochastic simulation algorithm which scales linearly.
P. Chatelain, A. Curioni, M. Bergdorf, D. Rossinelli, W. Andreoni, and P. Koumoutsakos, “Billion vortex particle direct numerical simulations of aircraft wakes,” Comput. Method. Appl. M. vol. 197, no. 13-16, pp. 1296–1304, 2008. Publisher's VersionAbstract
We present the Direct Numerical Simulations of high Reynolds numbers aircraft wakes employing vortex particle methods. The simulations involve a highly efficient implementation of vortex methods on massively parallel computers, enabling unprecedented simulations using billions of particles. The method relies on the Lagrangian discretization of the Navier-Stokes equations in vorticity-velocity form and relies on remeshing of the particles in order to ensure the convergence of the method. The remeshed particle locations are utilized for the computation of the field quantities, the discretization of the differential operators for diffusion and vortex stretching, and the solution of the Poisson equation for the advection velocity field. The method exhibits excellent scalability up to 16k BG/L nodes. The results include unprecedented Direct Numerical Simulations of the onset and the evolution of multiple wavelength instabilities induced by ambient noise in aircraft vortex wakes at Re = 6000. (c) 2007 Elsevier B.V. All rights reserved.
D. Rossinelli and P. Koumoutsakos, “Vortex methods for incompressible flow simulations on the GPU,” Visual Comput. vol. 24, no. 7-9, pp. 699–708, 2008. Publisher's VersionAbstract
We present a remeshed vortex particle method for incompressible flow simulations on GPUs. The particles are convected in a Lagrangian frame and are periodically reinitialized on a regular grid. The grid is used in addition to solve for the velocity–vorticity Poisson equation and for the computation of the diffusion operators. In the present GPU implementation of particle methods, the remeshing and the solution of the Poisson equation rely on fast and efficient mesh-particle interpolations. We demonstrate that particle remeshing introduces minimal artificial dissipation, enables a faster computation of differential operators on particles over grid-free techniques and can be efficiently implemented on GPUs. The results demonstrate that, contrary to common practice in particle simulations, it is necessary to remesh the (vortex) particle locations in order to solve accurately the equations they discretize, without compromising the speed of the method. The present method leads to simulations of incompressible vortical flows on GPUs with unprecedented accuracy and efficiency.
B. Bayati, P. Chatelain, and P. Koumoutsakos, “Multiresolution stochastic simulations of reaction–diffusion processes,” Phys. Chem. Chem. Phys. vol. 10, no. 39, pp. 5963, 2008. Publisher's VersionAbstract
Stochastic simulations of reaction–diffusion processes are used extensively for the modeling of complex systems in areas ranging from biology and social sciences to ecosystems and materials processing. These processes often exhibit disparate scales that render their simulation prohibitive even for massive computational resources. The problem is resolved by introducing a novel stochastic multiresolution method that enables the efficient simulation of reaction–diffusion processes as modeled by many-particle systems. The proposed method quantifies and efficiently handles the associated stiffness in simulating the system dynamics and its computational efficiency and accuracy are demonstrated in simulations of a model problem described by the Fisher–Kolmogorov equation. The method is general and can be applied to other many-particle models of physical processes.
P. Poncet, R. Hildebrand, G. - H. Cottet, and P. Koumoutsakos, “Spatially distributed control for optimal drag reduction of the flow past a circular cylinder,” J. Fluid Mech. vol. 599, 2008. Publisher's VersionAbstract
We report high drag reduction in direct numerical simulations of controlled flows past ircular cylinders at Reynolds numbers of 300 and 1000. The flow is controlled by the zimuthal component of the tangential velocity of the cylinder surface. Starting from spanwise-uniform velocity profile that leads to high drag reduction, the optimization rocedure identifies, for the same energy input, spanwise-varying velocity profiles that ead to higher drag reduction. The three-dimensional variations of the velocity field, orresponding to modes A and B of three-dimensional wake instabilities, are largely esponsible for this drag reduction. The spanwise wall velocity variations introduce treamwise vortex braids in the wake that are responsible for reducing the drag nduced by the primary spanwise vortices shed by the cylinder. The results demonstrate hat extending two-dimensional controllers to three-dimensional flows is not optimal s three-dimensional control strategies can lead efficiently to higher drag reduction.
G. Morra, P. Chatelain, P. Tackley, and P. Koumoutsakos, “Earth curvature effects on subduction morphology: Modeling subduction in a spherical setting,” Acta Geotech. vol. 4, no. 2, pp. 95–105, 2008. Publisher's VersionAbstract
{We present here the first application in Geodynamics of a (Fast ultipole) Accelerated Boundary Element Method (Accelerated-BEM) for tokes flow. The approach offers the advantages of a reduced number of computational elements and linear scaling with the problem size. We show that this numerical method can be fruitfully applied to the simulation of several geodynamic systems at the planetary scale in spherical coordinates and we suggest a general approach for modeling combined mantle convection and plate tectonics. The first part of the paper is devoted to the technical exposition of the new approach, while the second part focuses on the effect layed by Earth curvature on the subduction of a very wide oceanic ithosphere (W = 6000km and W = 9000km), comparing the effects of two different planetary radiuses (ER = 6371km

2007

S. Kern, P. Koumoutsakos, and K. Eschler, “Optimization of anguilliform swimming,” Phys. Fluids, vol. 19, no. 9, pp. 091102, 2007. Publisher's VersionAbstract
The European eel Anguilla anguilla migrates from the coasts of Europe to its spawning grounds in the Sargasso Sea. As the eels cover this 6000 km distance without feeding, anguilliform swimming has been regarded as a prime example of highly efficient aquatic propulsion.1 We investigate the hydrodynamics of anguilliform swimming motions using three-dimensional simulations of the fluid flow past a self-propelled body. An evolutionary optimization algorithm2 is used to determine the motion of the body for different objectives, linking swimming motion and biological function in a systematic fashion. The objectives are the swimming efficiency and the burst swimming speed of the swimmer as they pertain to migration and hunt/escape behavior, respectively. The kinematics of burst swimming is characterized by the large amplitude undulations of the tail and the straightness of the anterior part of the body. In contrast, during efficient swimming, significant lateral undulations are present along the entire length of the body. In burst swimming, the majority of the thrust is generated at the tail, whereas in efficient swimming, in addition to the tail, the central part of the body contributes significantly to the thrust.3 The wake, for both swimming modes, consists largely of a double row of vortex rings and corresponding lateral jets with an axis aligned with the swimming direction Fig. 1 and is consistent with experimental results.4
J. A. Helmuth, C. J. Burckhardt, P. Koumoutsakos, U. F. Greber, and I. F. Sbalzarini, “A novel supervised trajectory segmentation algorithm identifies distinct types of human adenovirus motion in host cells,” J. Struct. Biol. vol. 159, no. 3, pp. 347–358, 2007. Publisher's VersionAbstract
Biological trajectories can be characterized by transient patterns that may provide insight into the interactions of the moving object with its immediate environment. The accurate and automated identification of trajectory motifs is important for the understanding of the underlying mechanisms. In this work, we develop a novel trajectory segmentation algorithm based on supervised support vector classification. The algorithm is validated on synthetic data and applied to the identification of trajectory fingerprints of fluorescently tagged human adenovirus particles in live cells. In virus trajectories on the cell surface, periods of confined motion, slow drift, and fast drift are efficiently detected. Additionally, directed motion is found for viruses in the cytoplasm. The algorithm enables the linking of microscopic observations to molecular phenomena that are critical in many biological processes, including infectious pathogen entry and signal transduction. (c) 2007 Elsevier Inc. All rights reserved.
A. Dupuis, E. M. Kotsalis, and P. Koumoutsakos, “Coupling lattice Boltzmann and molecular dynamics models for dense fluids,” Phys. Rev. E, vol. 75, no. 4, 2007. Publisher's VersionAbstract
We propose a hybrid model, coupling lattice Boltzmann (LB) and molecular dynamics (MD) models, for the simulation of dense fluids. Time and length scales are decoupled by using an iterative Schwarz domain decomposition algorithm. The MD and LB formulations communicate via the exchange of velocities and velocity gradients at the interface. We validate the present LB-MD model in simulations of two- and three-dimensional flows of liquid argon past and through a carbon nanotube. Comparisons with existing hybrid algorithms and with reference MD solutions demonstrate the validity of the present approach.
A. Dupuis and P. Koumoutsakos, “Effects of Atomistic Domain Size on Hybrid Lattice Boltzmann–Molecular Dynamics Simulations of Dense Fluids,” Int. J. Mod. Phys. C, vol. 18, no. 04, pp. 644–651, 2007. Publisher's VersionAbstract
We present a convergence study for a hybrid Lattice Boltzmann-Molecular Dynamics model for the simulation of dense liquids. Time and length scales are decoupled by using an iterative Schwarz domain decomposition algorithm. The velocity field from the atomistic domain is introduced as forcing terms to the Lattice Boltzmann model of the continuum while the mean field of the continuum imposes mean field conditions for the atomistic domain. In the present paper we investigate the effect of varying the size of the atomistic subdomain in simulations of two dimensional flows of liquid argon past carbon nanotubes and assess the efficiency of the method.
P. Chatelain, G. - H. Cottet, and P. Koumoutsakos, “Particle Mesh Hydrodynamics for Astrophysics Simulations,” Int. J. Mod. Phys. C, vol. 18, no. 04, pp. 610–618, 2007. Publisher's VersionAbstract
We present a particle method for the simulation of three dimensional compressible hydrodynamics based on a hybrid Particle-Mesh discretization of the governing equations. The method is rooted on the regularization of particle locations as in remeshed Smoothed Particle Hydrodynamics (rSPH). The rSPH method was recently introduced to remedy problems associated with the distortion of computational elements in SPH, by periodically re-initializing the particle positions and by using high order interpolation kernels. In the PMH formulation, the particles solely handle the convective part of the compressible Euler equations. The particle quantities are then interpolated onto a mesh, where the pressure terms are computed. PMH, like SPH, is free of the convection CFL condition while at the same time it is more efficient as derivatives are computed on a mesh rather than particle-particle interactions. PMH does not detract from the adaptive character of SPH and allows for control of its accuracy. We present simulations of a benchmark astrophysics problem demonstrating the capabilities of this approach.
M. Bergdorf, P. Koumoutsakos, and A. Leonard, “Direct numerical simulations of vortex rings at Re_Gamma= 7500,” J. Fluid Mech. vol. 581, pp. 495–505, 2007. Publisher's VersionAbstract
We present direct numerical simulations of the turbulent decay of vortex rings with Re-Gamma = 7500. We analyse the vortex dynamics during the nonlinear stage of the instability along with the structure of the vortex wake during the turbulent stage. These simulations enable the quantification of vorticity dynamics and their correlation with structures from dye visualization and the observations of circulation decay that have been reported in related experimental works. Movies are available with the online version of the paper.
A. M. Altenhoff, J. H. Walther, and P. Koumoutsakos, “A stochastic boundary forcing for dissipative particle dynamics,” J. Comput. Phys. vol. 225, no. 1, pp. 1125–1136, 2007. Publisher's VersionAbstract
The method of dissipative particle dynamics (DPD) is an effective, coarse grained model of the hydrodynamics of complex fluids. DPD simulations of wall-bounded flows are however often associated with spurious fluctuations of the fluid properties near the wall. We present a novel stochastic boundary forcing for DPD simulations of wall-bounded flows, based on the identification of fluctuations in simulations of the corresponding homogeneous system at equilibrium. The present method is shown to enforce accurately the no-slip boundary condition, while minimizing spurious fluctuations of material properties, in a number of benchmark problems. (c) 2007 Elsevier Inc. All rights reserved.
P. A. E. Schoen, J. H. Walther, D. Poulikakos, and P. Koumoutsakos, “Phonon assisted thermophoretic motion of gold nanoparticles inside carbon nanotubes,” Appl. Phys. Lett. vol. 90, no. 25, pp. 253116, 2007. Publisher's VersionAbstract
The authors investigate the thermally driven mass transport of gold nanoparticles confined inside carbon nanotubes using molecular dynamics simulations. The observed thermophoretic motion of the gold nanoparticles correlates with the phonon dispersion exhibited by a standard carbon nanotube and, in particular, with the breathing mode of the tube. Additionally, the results show an increased static friction for gold nanoparticles confines inside a zig-zag carbon nanotube when increasing the size (length) of the nanoparticles. However, an unexpected, opposite trend is observed for the same nanoparticles inside armchair tubes. (c) 2007 American Institute of Physics.
E. M. Kotsalis, J. H. Walther, and P. Koumoutsakos, “Control of density fluctuations in atomistic-continuum simulations of dense liquids,” Phys. Rev. E, vol. 76, no. 1, 2007. Publisher's VersionAbstract
We present a control algorithm to eliminate spurious density fluctuations associated with the coupling of atomistic and continuum descriptions for dense liquids. A Schwartz domain decomposition algorithm is employed to couple molecular dynamics for the simulation of the atomistic system with a continuum solver for the simulation of the Navier-Stokes equations. The lack of periodic boundary conditions in the molecular dynamics simulations hinders the proper accounting for the virial pressure leading to spurious density fluctuations at the continuum-atomistic interface. An ad hoc boundary force is usually employed to remedy this situation. We propose the calculation of this boundary force using a control algorithm that explicitly cancels the density fluctuations. The results demonstrate that the present approach outperforms state-of-the-art algorithms. The conceptual and algorithmic simplicity of the method makes it suitable for any type of coupling between atomistic and continuum descriptions of dense fluids.

2006

N. Hansen, “An analysis of mutative Σ-self-adaptation on linear fitness functions,” Envol. Comput. vol. 14, no. 3, pp. 255–275, 2006.Abstract

This paper investigates sigma-self-adaptation for real valued evolutionary algorithms on linear fitness functions. We identify the step-size logarithm log sigma as a key quantity to understand strategy behavior. Knowing the bias of mutation, recombination, and selection on log sigma is sufficient to explain sigma-dynamics and strategy behavior in many cases, even from previously reported results on non-linear and/or noisy fitness functions. On a linear fitness function, if intermediate multi-recombination is applied on the object parameters, the i-th best and the i-th worst individual have the same sigma-distribution. Consequently, the correlation between fitness and step-size sigma is zero. Assuming additionally that sigma-changes due to mutation and recombination are unbiased, then sigma-self-adaptation enlarges sigma if and only if mu < lambda/2, given (mu, lambda)-truncation selection. Experiments show the relevance of the given assumptions.

P. A. E. Schoen, J. H. Walther, S. Arcidiacono, D. Poulikakos, and P. Koumoutsakos, “Nanoparticle Traffic on Helical Tracks: Thermophoretic Mass Transport through Carbon Nanotubes,” Nano Lett. vol. 6, no. 9, pp. 1910–1917, 2006. Publisher's VersionAbstract
Using molecular dynamics simulations, we demonstrate and quantify thermophoretic motion of solid gold nanoparticles inside carbon nanotubes subject to wall temperature gradients ranging from 0.4 to 25 K/nm. For temperature gradients below 1 K/nm, we find that the particles move "on tracks" in a predictable fashion as they follow unique helical orbits depending on the geometry of the carbon nanotubes. These findings markedly advance our knowledge of mass transport mechanisms relevant to nanoscale applications.
I. F. Sbalzarini, J. H. Walther, M. Bergdorf, S. E. Hieber, E. M. Kotsalis, and P. Koumoutsakos, “PPM – A highly efficient parallel particle–mesh library for the simulation of continuum systems,” J. Comput. Phys. vol. 215, no. 2, pp. 566–588, 2006. Publisher's VersionAbstract
This paper presents a highly efficient parallel particle-mesh (PPM) library, based on a unifying particle formulation for the simulation of continuous systems. In this formulation, the grid-free character of particle methods is relaxed by the introduction of a mesh for the reinitialization of the particles, the computation of the field equations, and the discretization of differential operators. The present utilization of the mesh does not detract from the adaptivity, the efficient handling of complex geometries, the minimal dissipation, and the good stability properties of particle methods. The coexistence of meshes and particles, allows for the development of a consistent and adaptive numerical method, but it presents a set of challenging parallelization issues that have hindered in the past the broader use of particle methods. The present library solves the key parallelization issues involving particle-mesh interpolations and the balancing of processor particle loading, using a novel adaptive tree for mixed domain decompositions along with a coloring scheme for the particle-mesh interpolation. The high parallel efficiency of the library is demonstrated in a series of benchmark tests on distributed memory and on a shared-memory vector architecture. The modularity of the method is shown by a range of simulations, from compressible vortex rings using a novel formulation of smooth particle hydrodynamics, to simulations of diffusion in real biological cell organelles. The present library enables large scale simulations of diverse physical problems using adaptive particle methods and provides a computational tool that is a viable alternative to mesh-based methods.
I. F. Sbalzarini, A. Hayer, A. Helenius, and P. Koumoutsakos, “Simulations of (An)Isotropic Diffusion on Curved Biological Surfaces,” Biophys. J. vol. 90, no. 3, pp. 878–885, 2006. Publisher's VersionAbstract
We present a computational particle method for the simulation of isotropic and anisotropic diffusion on curved biological surfaces that have been reconstructed from image data. The method is capable of handling surfaces of high curvature and complex shape, which are often encountered in biology. The method is validated on simple benchmark problems and is shown to be second-order accurate in space and time and of high parallel efficiency. It is applied to simulations of diffusion on the membrane of endoplasmic reticula ( ER) in live cells. Diffusion simulations are conducted on geometries reconstructed from real ER samples and are compared to fluorescence recovery after photobleaching experiments in the same ER samples using the transmembrane protein tsO45-VSV-G, C-terminally tagged with green fluorescent protein. Such comparisons allow derivation of geometry-corrected molecular diffusion constants for membrane components from fluorescence recovery after photobleaching data. The results of the simulations indicate that the diffusion behavior of molecules in the ER membrane differs significantly from the volumetric diffusion of soluble molecules in the lumen of the same ER. The apparent speed of recovery differs by a factor of similar to 4, even when the molecular diffusion constants of the two molecules are identical. In addition, the specific shape of the membrane affects the recovery half-time, which is found to vary by a factor of similar to 2 in different ER samples.
S. Kern and P. Koumoutsakos, “Simulations of optimized anguilliform swimming,” J. Exp. Biol. vol. 209, no. 24, pp. 4841–4857, 2006. Publisher's VersionAbstract
The hydrodynamics of anguilliform swimming motions was investigated using three-dimensional simulations of the fluid flow past a self-propelled body. The motion of the body is not specified a priori, but is instead obtained through an evolutionary algorithm used to optimize the swimming efficiency and the burst swimming speed. The results of the present simulations support the hypothesis that anguilliform swimmers modify their kinematics according to different objectives and provide a quantitative analysis of the swimming motion and the forces experienced by the body. The kinematics of burst swimming is characterized by the large amplitude of the tail undulations while the anterior part of the body remains straight. In contrast, during efficient swimming behavior significant lateral undulation occurs along the entire length of the body. In turn, during burst swimming, the majority of the thrust is generated at the tail, whereas in the efficient swimming mode, in addition to the tail, the middle of the body contributes significantly to the thrust. The burst swimming velocity is 42% higher and the propulsive efficiency is 15% lower than the respective values during efficient swimming. The wake, for both swimming modes, consists largely of a double row of vortex rings with an axis aligned with the swimming direction. The vortex rings are responsible for producing lateral jets of fluid, which has been documented in prior experimental studies. We note that the primary wake vortices are qualitatively similar in both swimming modes except that the wake vortex rings are stronger and relatively more elongated in the fast swimming mode. The present results provide quantitative information of three-dimensional fluid-body interactions that may complement related experimental studies. In addition they enable a detailed quantitative analysis, which may be difficult to obtain experimentally, of the different swimming modes linking the kinematics of the motion with the forces acting on the self-propelled body. Finally, the optimization procedure helps to identify, in a systematic fashion, links between swimming motion and biological function.
K. Fukagata, N. Kasagi, and P. Koumoutsakos, “A theoretical prediction of friction drag reduction in turbulent flow by superhydrophobic surfaces,” Phys. Fluids, vol. 18, no. 5, pp. 051703, 2006. Publisher's VersionAbstract
We present a theoretical prediction for the drag reduction rate achieved by superhydrophobic surfaces in a turbulent channel flow. The predicted drag reduction rate is in good agreement with results obtained from direct numerical simulations at Re-tau similar or equal to 180 and 400. The present theory suggests that large drag reduction is possible also at Reynolds numbers of practical interest (Re-tau similar to 10(5)-10(6)) by employing a hydrophobic surface, which induces a slip length on the order of ten wall units or more.
M. Bergdorf and P. Koumoutsakos, “A Lagrangian Particle-Wavelet Method,” Multiscale Model. Sim. vol. 5, no. 3, pp. 980–995, 2006. Publisher's VersionAbstract
This paper presents a novel, multiresolution Lagrangian particle method with enhanced, wavelet-based adaptivity. The method is formulated for transport problems and combines the efficiency of wavelet collocation with the inherent numerical stability of particle methods. The accuracy and efficiency of the present method is assessed on a number of benchmark problems pertaining to interface capturing and transport. The method is compared with existing techniques demonstrating its advantages and limitations. The present approach leads to a new generation of particle methods with multiresolution capabilities.
A. Auger, P. Chatelain, and P. Koumoutsakos, “R-leaping: Accelerating the stochastic simulation algorithm by reaction leaps,” J. Chem. Phys. vol. 125, no. 8, pp. 084103, 2006. Publisher's VersionAbstract
A novel algorithm is proposed for the acceleration of the exact stochastic simulation algorithm by a predefined number of reaction firings (R-leaping) that may occur across several reaction channels. In the present approach, the numbers of reaction firings are correlated binomial distributions and the sampling procedure is independent of any permutation of the reaction channels. This enables the algorithm to efficiently handle large systems with disparate rates, providing substantial computational savings in certain cases. Several mechanisms for controlling the accuracy and the appearance of negative species are described. The advantages and drawbacks of R-leaping are assessed by simulations on a number of benchmark problems and the results are discussed in comparison with established methods.

2005

P. Koumoutsakos, “Multicscale Flow Simulations using Particles,” Annu. Rev. Fluid Mech. vol. 37, no. 1, pp. 457–487, 2005. Publisher's VersionAbstract
Flow simulations are one of the archetypal multiscale problems. Simulations of turbulent and unsteady separated flows have to resolve a multitude of interacting scales, whereas molecular phenomena determine the structure of shocks and the validity of the no-slip boundary condition. Particle simulations of continuum and molecular phenomena can be formulated by following the motion of interacting particles that carry the physical properties of the flow. In this article we review Lagrangian, multiresolution, particle methods such as vortex methods and smooth particle hydrodynamics for the simulation of continuous flows and molecular dynamics for the simulation of flows at the atomistic scale. We review hybrid molecular-continuum simulations with an emphasis on the computational aspects of the problem. We identify the common computational characteristics of particle methods and discuss their properties that enable the formulation of a systematic framework for multiscale flow simulations.
E. M. Kotsalis, E. Demosthenous, J. H. Walther, S. C. Kassinos, and P. Koumoutsakos, “Wetting of doped carbon nanotubes by water droplets,” Chem. Phys. Lett. vol. 412, no. 4-6, pp. 250–254, 2005. Publisher's VersionAbstract
We study the wetting of doped single- and multi-walled carbon nanotubes by water droplets using molecular dynamics simulations. Chemisorbed hydrogen is considered as a model of surface impurities. We study systems with varying densities of surface impurities and we observe increased wetting, as compared to the pristine nanotube case, attributed to the surface dipole moment that changes the orientation of the interfacial water. We demonstrate that the nature of the impurity is important as here hydrogen induces the formation of an extended hydrogen bond network between the water molecules and the doping sites leading to enhanced wetting.
S. E. Hieber and P. Koumoutsakos, “A Lagrangian particle level set method,” J. Comput. Phys. vol. 210, no. 1, pp. 342–367, 2005. Publisher's VersionAbstract
We present a novel particle level set method for capturing interfaces. The level set equation is solved in a Lagrangian frame using particles that carry the level set information. A key aspect of the method involves a consistent remeshing procedure for the regularization of the particle locations when the particle map gets distorted by the advection field. The Lagrangian description of the level set method is inherently adaptive and exact in the case of solid body motions. The efficiency and accuracy of the method is demonstrated in several benchmark problems in two and three dimensions involving pure advection and curvature induced motion of the interface. The simplicity of the particle description is shown to be well suited for real time simulations of surfaces involving cutting and reconnection as in virtual surgery environments.
U. Zimmerli, P. G. Gonnet, J. H. Walther, and P. Koumoutsakos, “Curvature Induced L-Defects in Water Conduction in Carbon Nanotubes,” Nano Lett. vol. 5, no. 6, pp. 1017–1022, 2005. Publisher's VersionAbstract
We conduct molecular dynamics simulations to study the effect of the curvature induced static dipole moment of small open-ended single-walled carbon nanotubes (CNTs) immersed in water. This dipole moment generates a nonuniform electric field, changing the energy landscape in the CNT and altering the water conduction process. The CNT remains practically filled with water at all times, whereas intermittent filling is observed when the dipole term is not included. In addition, the dipole moment induces a preferential orientation of the water molecules near the end regions of the nanotube, which in turn causes a reorientation of the water chain in the middle of the nanotube. The most prominent feature of this reorientation is an L-defect in the chain of water molecules inside the CNT. The analysis of the water energetics and structural characteristics inside and in the vicinity of the CNT helps to identify the role of the dipole moment and to suggest possible mechanisms for controlled water and proton transport at the nanoscale.
H. Ewers, A. E. Smith, I. F. Sbalzarini, H. Lilie, P. Koumoutsakos, and A. Helenius, “Single-particle tracking of murine polyoma virus-like particles on live cells and artificial membranes,” P. Natl. Acad. Sci. vol. 102, no. 42, pp. 15110–15115, 2005. Publisher's VersionAbstract
The lateral mobility of individual murine polyoma virus-like particles (VLPs) bound to live cells and artificial lipid bilayers was studied by single fluorescent particle tracking using total internal reflection fluorescence microscopy. The particle trajectories were analyzed in terms of diffusion rates and modes of motion as described by the moment scaling spectrum. Although VLPs bound to their ganglioside receptor in lipid bilayers exhibited only free diffusion, analysis of trajectories on live 3T6 mouse fibroblasts revealed three distinct modes of mobility: rapid random motion, confined movement in small zones (30-60 nm in diameter), and confined movement in zones with a slow drift. After binding to the cell surface, particles typically underwent free diffusion for 5-10 s, and then they were confined in an actin filament-dependent manner without involvement of clathrin-coated pits or caveolae. Depletion of cholesterol dramatically reduced mobility of VLPs independently of actin, whereas inhibition of tyrosine kinases had no effect on confinement. The results suggested that clustering of ganglioside molecules by the multivalent VLPs induced transmembrane coupling that led to confinement of the virus/receptor complex by cortical actin filaments.
T. Werder, J. H. Walther, and P. Koumoutsakos, “Hybrid atomistic–continuum method for the simulation of dense fluid flows,” J. Comput. Phys. vol. 205, no. 1, pp. 373–390, 2005. Download PDFAbstract
We present a hybrid atomistic-continuum method for multiscale simulations of dense fluids. In this method, the atomistic part is described using a molecular dynamics description, while the continuum flow is described by a finite volume discretization of the incompressible Navier-Stokes equations. The two descriptions are combined in a domain decomposition formulation using the Schwarz alternating method. A novel method has been proposed in order to impose non-periodic velocity boundary conditions from the continuum to the atomistic domain, based on an effective boundary potential, consistent body forces, a particle insertion algorithm and specular walls. The extraction of velocity boundary conditions for the continuum from the atomistic domain is formulated by taking into account the associated statistical errors. The advantages and drawbacks of the proposed Schwarz decomposition method as compared to related flux-based schemes are discussed. The efficiency and applicability of the method is demonstrated by considering hybrid and full molecular dynamics simulations of the flow of a Lennard-Jones fluid past a carbon nanotube.
D. Bueche, N. N. Schraudolph, and P. Koumoutsakos, “Accelerating Evolutionary Algorithms With Gaussian Process Fitness Function Models,” IEEE T. Syst. Man Cy. C, vol. 35, no. 2, pp. 183–194, 2005. Publisher's VersionAbstract
We present an overview of evolutionary algorithms that use empirical models of the fitness function to accelerate convergence, distinguishing between evolution control and the surrogate approach. We describe the Gaussian process model and propose using it as an inexpensive fitness function surrogate. Implementation issues such as efficient and numerically stable computation, exploration versus exploitation, local modeling, multiple objectives and constraints, and failed evaluations are addressed. Our resulting Gaussian process optimization procedure clearly outperforms other evolutionary strategies on standard test functions as well as on a real-world problem: the optimization of stationary gas turbine compressor profiles.
I. F. Sbalzarini and P. Koumoutsakos, “Feature point tracking and trajectory analysis for video imaging in cell biology,” J. Struct. Biol. vol. 151, no. 2, pp. 182–195, 2005. Publisher's VersionAbstract
This paper presents a computationally efficient, two-dimensional, feature point tracking algorithm for the automated detection and quantitative analysis of particle trajectories as recorded by video imaging in cell biology. The tracking process requires no a priori mathematical modeling of the motion, it is self-initializing, it discriminates spurious detections, and it can handle temporary occlusion as well as particle appearance and disappearance from the image region. The efficiency of the algorithm is validated on synthetic video data where it is compared to existing methods and its accuracy and precision are assessed for a wide range of signal-to-noise ratios. The algorithm is well suited for video imaging in cell biology relying on low-intensity fluorescence microscopy. Its applicability is demonstrated in three case studies involving transport of low-density lipoproteins in endosomes, motion of fluorescently labeled Adenovirus-2 particles along microtubules, and tracking of quantum dots on the plasma membrane of live cells. The present automated tracking process enables the quantification of dispersive processes in cell biology using techniques such as moment scaling spectra.
M. Bergdorf, G. - H. Cottet, and P. Koumoutsakos, “Multilevel Adaptive Particle Methods for Convection-Diffusion Equations,” Multiscale Model. Sim. vol. 4, no. 1, pp. 328–357, 2005. Publisher's VersionAbstract
We present novel multilevel particle methods with extended adaptivity in areas where increased resolution is required. We present two complementary approaches as inspired by r-adaptivity and adaptive mesh refinement (AMR) concepts introduced infinite difference and finite element schemes. For the r-adaptivity a new class of particle-based mapping functions is introduced, while the particle AMR method uses particle remeshing in overlapping domains as a key element. The advantages and drawbacks of the proposed particle methods are illustrated based on results from the solution of one-dimensional convection-diffusion equations, while the extension of the method to higher dimensions is demonstrated in simulations of the inviscid evolution of an elliptical vortex.
I. F. Sbalzarini, A. Mezzacasa, A. Helenius, and P. Koumoutsakos, “Effects of Organelle Shape on Fluorescence Recovery after Photobleaching,” Biophys. J. vol. 89, no. 3, pp. 1482–1492, 2005. Publisher's VersionAbstract
The determination of diffusion coefficients from fluorescence recovery data is often complicated by geometric constraints imposed by the complex shapes of intracellular compartments. To address this issue, diffusion of proteins in the lumen of the endoplasmic reticulum (ER) is studied using cell biological and computational methods. Fluorescence recovery after photobleaching (FRAP) experiments are performed in tissue culture cells expressing GFP - KDEL, a soluble, fluorescent protein, in the ER lumen. The three-dimensional (3D) shape of the ER is determined by confocal microscopy and computationally reconstructed. Within these ER geometries diffusion of solutes is simulated using the method of particle strength exchange. The simulations are compared to experimental FRAP curves of GFP - KDEL in the same ER region. Comparisons of simulations in the 3D ER shapes to simulations in open 3D space show that the constraints imposed by the spatial confinement result in two-to fourfold underestimation of the molecular diffusion constant in the ER if the geometry is not taken into account. Using the same molecular diffusion constant in different simulations, the observed speed of fluorescence recovery varies by a factor of 2.5, depending on the particular ER geometry and the location of the bleached area. Organelle shape considerably influences diffusive transport and must be taken into account when relating experimental photobleaching data to molecular diffusion coefficients. This novel methodology combines experimental FRAP curves with high accuracy computer simulations of diffusion in the same ER geometry to determine the molecular diffusion constant of the solute in the particular ER lumen.
S. Arcidiacono, J. H. Walther, D. Poulikakos, D. Passerone, and P. Koumoutsakos, “Solidification of Gold Nanoparticles in Carbon Nanotubes,” Phys. Rev. Lett. vol. 94, no. 10, 2005. Publisher's VersionAbstract
The structure and the solidification of gold nanoparticles in a carbon nanotube are investigated using molecular dynamics simulations. The simulations indicate that the predicted solidification temperature of the enclosed particle is lower than its bulk counterpart, but higher than that observed for clusters placed in vacuum. A comparison with a phenomenological model indicates that, in the considered range of tube radii (R-CNT) of 0.5 < R-CNT < 1.6 nm, the solidification temperature depends mainly on the length of the particle with a minor dependence on R-CNT.
P. Poncet, G. - H. Cottet, and P. Koumoutsakos, “Control of three-dimensional wakes using evolution strategies,” C.R. Mecanique, vol. 333, no. 1, pp. 65–77, 2005. Publisher's VersionAbstract
We investigate three-dimensional cylinder wakes of incompressible fully developed flows at Re = 300, resulting from control induced by tangential motions of the cylinder surface. The motion of the cylinder surface, in two dimensions, is optimized using evolution strategies, resulting in significant drag reduction and drastic modification of the wake as compared to the uncontrolled flow. The quasi-optimal velocity profile obtained in 2D is modified by spanwise harmonics and applied to 3D flows. The results indicate important differences in the flow physics induced by two and three dimensional control strategies.
P. Poncet and P. Koumoutsakos, “Optimization of Vortex Shedding In 3-D Wakes Using Belt Actuators,” Int. J. Offshore Polar, vol. 15, no. 1, pp. 7–13, 2005. Publisher's VersionAbstract
This paper discusses the control of cylinder wakes via tangential wall velocity modifications. The wall velocity is piecewise constant (corresponding to belt actuators), and its amplitude is optimized using a clustering real coded genetic algorithm. This type of control significantly affects the vortical structures being shed in the wake, and it is shown that the flow gets significantly modified, resulting in a 3-dimensional body shedding 2-dimensional vortical structures in the near wake of the body. Depending on the energy involved in the control, one can also obtain a dramatic decrease of shedding amplitude.

2004

J. H. Walther, et al., “Water–carbon interactions III: The influence of surface and fluid impurities,” Phys. Chem. Chem. Phys. vol. 6, no. 8, pp. 1988–1995, 2004. Publisher's VersionAbstract
Molecular dynamics. simulations are performed to study the influence of surface and fluid impurities on water-carbon interactions. In order to quantify these interactions we consider the canonical problem of wetting of a doped flat graphitic surface by a water system with impurities. As model fluid impurities we consider aqueous solutions of potassium-chloride with molar concentrations up to 1.8 M. Quantum chemistry calculations are performed to derive pair potentials for the ion-graphite interactions. The contact angle is found to decrease weakly with increasing ionic concentration, from 90degrees at 0 M to 81degrees at 1.8 M concentration. The influence of solid impurities is found to be more significant. Thus, 10, 15, and 20% coverages of chemisorbed hydrogen result in contact angles of 90degrees, 74degrees and 60degrees, respectively.
S. D. Müller, I. Mezić, J. H. Walther, and P. Koumoutsakos, “Transverse momentum micromixer optimization with evolution strategies,” Comput. Fluids, vol. 33, no. 4, pp. 521–531, 2004. Publisher's VersionAbstract
We conduct a numerical study of mixing in a transverse momentum micromixer. Good values for actuation frequencies can be determined using simple kinematic arguments, and evolution strategies are introduced for the optimization of mixing by adjusting the control parameters in micromixer devices. It is shown that the chosen optimization algorithm can identify, in an automated fashion, effective actuation parameters. We find that optimal frequencies for increasing number of transverse channels are superposable despite the non-linear nature of the mixing process.
M. Milano, P. Koumoutsakos, and J. Schmidhuber, “Self-Organizing Nets for Optimization,” IEEE T. Neural Networ. vol. 15, no. 3, pp. 758–765, 2004. Publisher's VersionAbstract
Given some optimization problem and a series of typically expensive trials of solution candidates sampled from a search space, how can we efficiently select the next candidate? We address this fundamental problem by embedding simple optimization strategies in learning algorithms inspired by Kohonen's self-organizing maps and neural gas networks. Our adaptive nets or grids are used to identify and exploit search space regions that maximize the probability of generating points closer to the optima. Net nodes are attracted by candidates that lead to improved evaluations, thus, quickly biasing the active data selection process toward promising regions, without loss of ability to escape from local optima. On standard benchmark functions, our techniques perform more reliably than the widely used covariance matrix adaptation evolution strategy. The proposed algorithm is also applied to the problem of drag reduction in a flow past an actively controlled circular cylinder, leading to unprecedented drag reduction.
U. Zimmerli, M. Parrinello, and P. Koumoutsakos, “Dispersion corrections to density functionals for water aromatic interactions,” J. Chem. Phys. vol. 120, no. 6, pp. 2693–2699, 2004. Publisher's VersionAbstract
We investigate recently published methods for extending density functional theory to the description of long-range dispersive interactions. In all schemes an empirical correction consisting of a C(6)r(-6) term is introduced that is damped at short range. The coefficient C-6 is calculated either from average molecular or atomic polarizabilities. We calculate geometry-dependent interaction energy profiles for the water benzene cluster and compare the results with second-order Moller-Plesset calculations. Our results indicate that the use of the B3LYP functional in combination with an appropriate mixing rule and damping function is recommended for the interaction of water with aromatics. (C) 2004 American Institute of Physics.
J. H. Walther, R. L. Jaffe, E. M. Kotsalis, T. Werder, T. Halicioglu, and P. Koumoutsakos, “Hydrophobic hydration of C60 and carbon nanotubes in water,” Carbon, vol. 42, no. 5-6, pp. 1185–1194, 2004. Publisher's VersionAbstract
We perform molecular dynamics (MD) simulations to study the hydrophobic-hydrophilic behavior of pairs Of C-60 fullerene molecules and single wall carbon nanotubes in water. The interaction potentials involve a fully atomistic description of the fullerenes or carbon nanotubes and the water is modeled using the flexible SPC model. Both unconstrained and constrained MD simulations are carried out. We find that these systems display drying, as evidenced by expulsion of the interstitial water, when the C-60 and carbon nanotubes are separated by less than 12, and 9-10 Angstrom, respectively. From the constrained simulations, the computed mean force between two carbon nanotubes in water exhibits a maximum at a tube spacing of 5.0 Angstrom which corresponds to approximately one unstable layer of interstitial water molecules. The main contribution to the force stems from the van der Waals attraction between the carbon surfaces. The minimum in the potential of mean force has a value of - 17 kJ mol(-1) Angstrom(-1) at a tube spacing of 3.5 Angstrom.
J. H. Walther, T. Werder, R. L. Jaffe, and P. Koumoutsakos, “Hydrodynamic properties of carbon nanotubes,” Phys. Rev. E, vol. 69, no. 6, 2004. Publisher's VersionAbstract
We study water flowing past an array of single walled carbon nanotubes using nonequilibrium molecular dynamics simulations. For carbon nanotubes mounted with a tube spacing of 16.4x16.4 nm and diameters of 1.25 and 2.50 nm, respectively, we find drag coefficients in reasonable agreement with the macroscopic, Stokes-Oseen solution. The slip length is -0.11 nm for the 1.25 nm carbon nanotube, and 0.49 for the 2.50 nm tube for a flow speed of 50 m/s, respectively, and 0.28 nm for the 2.50 nm tube at 200 m/s. A slanted flow configuration with a stream- and spanwise velocity component of 100 ms(-1) recovers the two-dimensional results, but exhibits a significant 88 nm slip along the axis of the tube. These results indicate that slip depends on the particular flow configuration.
E. M. Kotsalis, J. H. Walther, and P. Koumoutsakos, “Multiphase water flow inside carbon nanotubes,” Int. J. Multiphas. Flow, vol. 30, no. 7-8, pp. 995–1010, 2004. Publisher's VersionAbstract
We present nonequilibrium molecular dynamics simulations of the flow of liquid-vapour water mixtures and mixtures of water and nitrogen inside carbon nanotubes. A new adaptive forcing scheme is proposed to impose a mean flow through the system. The flow of liquid water is characterised by a distinct layering of the water molecules in the vicinity of the boundary and a slip length that is found to increase with the radius of the carbon nanotube. Increasing the temperature and pressure of the system furthermore results in a decrease in the slip length. For the flow of mixtures of nitrogen and water we find that the slip length is reduced as compared to the slip for the pure water. The shorter slip length is attributed to the fact that nitrogen forms droplets at the carbon surface, thus partially shielding the bulk flow from the hydrophobic carbon surface.
S. Kern, S. D. Müller, N. Hansen, D. Büche, J. Ocenasek, and P. Koumoutsakos, “Learning probability distributions in continuous evolutionary algorithms – a comparative review,” Nat. Comput. vol. 3, no. 1, pp. 77–112, 2004. Publisher's VersionAbstract
We present a comparative review of Evolutionary Algorithms that generate new population members by sampling a probability distribution constructed during the optimization process. We present a unifying formulation for five such algorithms that enables us to characterize them based on the parametrization of the probability distribution, the learning methodology, and the use of historical information. The algorithms are evaluated on a number of test functions in order to assess their relative strengths and weaknesses. This comparative review helps to identify areas of applicability for the algorithms and to guide future algorithmic developments.
R. L. Jaffe, P. Gonnet, T. Werder, J. H. Walther, and P. Koumoutsakos, “Water–Carbon Interactions 2: Calibration of Potentials using Contact Angle Data for Different Interaction Models,” Mol. Simulat. vol. 30, no. 4, pp. 205–216, 2004. Publisher's VersionAbstract
Molecular dynamics simulations of water droplets on graphite are carried out to determine the contact angle for different water-carbon potential functions. Following the procedure of Werder et al. [ J. Phys. Chem. B , 107 (2003) 1345], the C-O Lennard-Jones well depth is varied to recover the experimental value for the contact angle (84-86degrees) using a 2000-molecule water droplet and compensating for the line tension effect that lowers the contact angle for increasing droplet size. For the discrete graphite surface model studied by Werder et al. , the effects of adding C-H Lennard-Jones interactions and changing the long-range cut-off distance are considered. In addition, a continuum graphite surface model is studied for which the water-graphite interaction energy depends only on the normal distance ( z ) from the water oxygen to the surface. This new model, with z(-10) repulsion and z(-4) attraction, is formulated in terms of the standard Lennard-Jones parameters, for which the recommended values are sigma(CO) =3.19 Angstrom and epsilon(CO) =0.3651 kJ/mol.
S. E. Hieber, J. H. Walther, and P. Koumoutsakos, “Remeshed smoothed particle hydrodynamics simulation of the mechanical behavior of human organs,” Technol. Health Care, vol. 12, no. 4, pp. 305–314, 2004. Publisher's VersionAbstract
In computer aided surgery the accurate simulation of the mechanical behavior of human organs is essential for the development of surgical simulators. In this paper we introduce particle based simulations of two different human organ materials modeled as linear viscoelastic solids. The constitutive equations for the material behavior are discretized using a particle approach based on the smoothed particle hydrodynamics (SPH) method while the body surface is tracked using level sets. A key aspect of this approach is its flexibility which allows the simulation of complex time varying topologies with large deformations. The accuracy of the original formulation is significantly enhanced by using a particle reinitialization technique resulting in remeshed smoothed particle hydrodynamics (rSPH). The mechanical parameters of the systems used in the simulations are derived from experimental measurements on human cadaver organs. We compare the mechanical behavior of liver- and kidney-like materials based on the dynamic simulations of a tensile test case. Moreover, we present a particle based reconstruction of the liver topology and its strain distribution under a small local load. Finally, we demonstrate a unified formulation of fluid structure interaction based on particle method.

2003

T. Werder, J. H. Walther, R. L. Jaffe, T. Halicioglu, and P. Koumoutsakos, “On the Water-Carbon Interaction for Use in Molecular Dynamics Simulations of Graphite and Carbon Nanotubes,” J. Phys. Chem. B, vol. 112, no. 44, pp. 14090–14090, 2003. Publisher's VersionAbstract
A systematic molecular dynamics study shows that the contact angle of a water droplet on graphite changes significantly as a function of the water-carbon interaction energy. Together with the observation that a linear relationship can be established between the contact angle and the water monomer binding energy on graphite, a new route to calibrate interaction potential parameters is presented. Through a variation of the droplet size in the range from 1000 to 17'500 water molecules, we determine the line tension to be positive and on the order of 2 × 10\^-10 J/m. To recover a macroscopic contact angle of 86°, a water monomer binding energy of -6.33 kJ mol\^-1 is required, which is obtained by applying a carbon-oxygen Lennard-Jones potential with the parameters eps_CO = 0.392 kJ mol\^-1 and sigma_CO = 3.19 \AA. For this new water-carbon interaction potential, we present density profiles and hydrogen bond distributions for a water droplet on graphite.
N. Hansen, S. D. Müller, and P. Koumoutsakos, “Reducing the Time Complexity of the Derandomized Evolution Strategy with Covariance Matrix Adaptation (CMA-ES),” Evol. Comput. vol. 11, no. 1, pp. 1–18, 2003. Publisher's VersionAbstract
This paper presents a novel evolutionary optimization strategy based on the derandomized evolution strategy with covariance matrix adaptation (CMA-ES). This new approach is intended to reduce the number of generations required for convergence to the optimum. Reducing the number of generations, i.e., the time complexity of the algorithm, is important if a large population size is desired: (1) to reduce the effect of noise; (2) to improve global search properties; and (3) to implement the algorithm on (highly) parallel machines. Our method results in a highly parallel algorithm which scales favorably with large numbers of processors. This is accomplished by efficiently incorporating the available information from a large population, thus significantly reducing the number of generations needed to adapt the covariance matrix. The original version of the CMA-ES was designed to reliably adapt the covariance matrix in small populations but it cannot exploit large populations efficiently. Our modifications scale up the efficiency to population sizes of up to 10n, where n is the problem dimension. This method has been applied to a large number of test problems, demonstrating that in many cases the CMA-ES can be advanced from quadratic to linear time complexity.

2001

N. Hansen and A. Ostermeier, “Completely derandomized self-adaptation in evolution strategies,” Evol. Comput. vol. 9, no. 2, pp. 159–195, 2001.Abstract

This paper puts forward two useful methods for self-adaptation of the mutation distribution – the concepts of derandomization and cumulation. Principle shortcomings of the concept of mutative strategy parameter control and two levels of derandomization are reviewed. Basic demands on the self-adaptation of arbitrary (normal) mutation distributions are developed. Applying arbitrary, normal Mutation distributions is equivalent to applying a general, linear problem encoding. The underlying objective of mutative strategy parameter control is roughly to favor previously selected mutation steps in the future. If this objective is pursued rigorously, a completely derandomized self-adaptation scheme results, which adapts arbitrary normal mutation distributions. This scheme, called covariance matrix adaptation (CMA), meets the previously stated demands. It can still be considerably improved by cumulation – utilizing an evolution path rather than single search steps. Simulations on various test functions reveal local and global search properties of the evolution strategy with and without covariance matrix adaptation. Their performances are comparable only on perfectly scaled functions. On badly scaled, nonseparable functions usually a speed up factor of several orders of magnitude is observed. On moderately mis-scaled functions a speed up factor of three to ten can be expected.

S. Zimmermann, P. Koumoutsakos, and W. Kinzelbach, “Simulation of Pollutant Transport Using a Particle Method,” J. Comput. Phys. vol. 173, no. 1, pp. 322–347, 2001. Publisher's VersionAbstract
A Lagrangian scheme is presented for the two-dimensional simulation of pollutant transport in a porous medium. The anisotropic extension of the particle strength exchange method is implemented to describe the diffusive–dispersive process. By applying the scheme to a benchmark problem with an analytical solution, the method is shown to be accurate and stable even in the limiting cases of vanishing diffusion– dispersion coefficients and high anisotropy ratios. The scheme is shown to perform well with spatially variable velocity fields also.
T. Werder, J. H. Walther, R. L. Jaffe, T. Halicioglu, F. Noca, and P. Koumoutsakos, “Molecular Dynamics Simulation of Contact Angles of Water Droplets in Carbon Nanotubes,” Nano Lett. vol. 1, no. 12, pp. 697–702, 2001. Publisher's VersionAbstract
We examine the equilibrium shape of the water droplet and deduce its contact angle with the CNT surface. This is done using isochore line plots of the water following the algorithm suggested by Nijmeijer. Additionally, radial water density profiles and radial hydrogen bond distributions are reported. The systems considered include zigzag CNTs with diameters of 25.0, 50.0, and 75.0 \AA respectively and initially flat 23.7 \AA thick drops (this corresponds to 12 layers of water molecules). The simulations indicate that pure water does not wet pristine single wall carbon nanotubes.
J. H. Walther and P. Koumoutsakos, “Three-Dimensional Vortex Methods for Particle-Laden Flows with Two-Way Coupling,” J. Comput. Phys. vol. 167, no. 1, pp. 39–71, 2001. Publisher's VersionAbstract
This paper presents a three-dimensional viscous vortex method for the simulation of particulate flows with two-way coupling. The flow is computed using Lagrangian vortex elements advected with the local velocity, while their strength is modified to account for viscous diffusion, vortex stretching, and generating vorticity induced by the particles. The solid particles move according to viscous drag and gravity, creating vorticity, which is discretised using vortex elements. This method adaptively tracks the evolution of the vorticity field and the generation of new computational elements to account for the vorticity source term. A key aspect of the present scheme is the remeshing of the computational elements to adaptively accommodate the production of vorticity induced by the solid particles, and to ensure sufficient support for the proper resolution of the diffusion equation. High-order moment-conserving formulas are implemented to maintain the adaptive character of the method while they remain local to minimize the computational cost. These formulas are also implemented in the particle–mesh interpolation of the field and particle quantities in the context of a Vortex-in-Cell algorithm. The method is validated against the results of a related finite-difference study for an axisymmetric swirling flow with particles. The method is then applied to the study of a three-dimensional particle blob falling under the effect of gravity. It is shown that drastically different behaviours are found depending on the presence of an initial vorticity field.
J. H. Walther, R. Jaffe, T. Halicioglu, and P. Koumoutsakos, “Carbon Nanotubes in Water: Structural Characteristics and Energetics,” J. Phys. Chem. B, vol. 105, no. 41, pp. 9980–9987, 2001. Publisher's VersionAbstract
We study the structural properties of water surrounding a carbon nanotube using molecular dynamics simulations. The interaction potentials involve a description of the carbon nanotube using Morse, harmonic bending, torsion, and Lennard-Jones potentials. The water is described by the flexible Simple Point Charge (SPC) model by Teleman et al., and the carbon-water interactions include a carbon-oxygen Lennard-Jones potential, and an electrostatic quadrupole moment acting between the carbon atoms and the charge sites of the water. Vibration of the breathing mode of the carbon nanotube in water is inferred from the oscillations in carbon-carbon van der Waals energy, and the inverse proportionality between the radius of the carbon nanotube and the breathing frequency is in good agreement with experimental values. The results indicate, that under the present conditions, the presence of the water has a negligible influence on the breathing frequency. The water at the carbon-water interface is found to have a HOH plane nearly tangential to the interface, and the water radial density profile exhibits the characteristic layering also found in the graphite-water system. The average number of hydrogen bonds decreases from a value of 3.73 in the bulk phase to a value of 2.89 at the carbon-water interface. The inclusion of the carbon quadrupole moment is found to have a negligible influence on the structural properties of the water. Energy changes that occur by the process of introducing a carbon nanotube in water are calculated. The creation of a cavity in the bulk water to accommodate the nanotube constitutes the largest energy contribution. Results include an analysis of surface structure and energy values for planar and for concave cylindrical surfaces of water.
J. H. Walther and P. Koumoutsakos, “Molecular Dynamics Simulation of Nanodroplet Evaporation,” J. Heat Transf. vol. 123, no. 4, pp. 741, 2001. Publisher's VersionAbstract
We conduct a series of high resolution simulations using several tens of thousands of computational molecules by efficiently employing tree data structures. We present a systematic convergence study to examine the influence of the different parameters imposed by the numerical method, including the size of the computational domain, cutoff radius, droplet diameter, heating frequency and size of heating region, initialisation period, and time step size. It is shown that molecular dynamics simulations converge to the D\^2 evaporation law to desired accuracy by using large number of computational molecules.
S. D. Mueller, J. H. Walther, and P. D. Koumoutsakos, “Evolution strategies for film cooling optimization,” AIAA J. vol. 39, pp. 537–539, 2001. Publisher's VersionAbstract
An evolutionary algorithm is implemented in a realistic automated design cycle of turbine blade film cooling. This design cycle involves the use of empirical formulas for the flow calculation and the use of a commercial software package for the simulation of the heat transfer problem. The overall process consists of an engineering multiobjective optimization problem with constraints. In this work we consider the solution of this problem in the context of an automated optimization cycle using evolution strategies.
P. Koumoutsakos, J. Freund, and D. Parekh, “Evolution strategies for automatic optimization of jet mixing,” AIAA J. vol. 39, pp. 967–969, 2001. Publisher's VersionAbstract
Evolution strategies (ES) are introduced for the optimization of active control parameters for enhancing jet mixing. It is shown that the evolution algorithms can identify, in an automated fashion, not only previously known effective actuations but also find good but previously unidentified parameters. In this study, simulations of model jets are used to demonstrate the feasibility of the methods. ES are robust, highly parallel, and portable algorithms that may be most useful in an experimental setting at realistic Reynolds numbers. Simulations of inviscid incompressible flows using vortex models, as well as direct numerical simulations (DNS) of very low-Reynolds-number compressible flows, are used in this study to evaluate different forcing parameters.