Publications

2016
S. Litvinov, Q. Xie, X. Hu, N. Adams, and M. Ellero, “Simulation of individual polymer chains and polymer solutions with smoothed dissipative particle dynamics,” Fluids, vol. 1, no. 1, pp. 7, 2016. Publisher's VersionAbstract
In an earlier work (Litvinov et al., Phys.Rev.E 77, 066703 (2008)), a model for a polymer molecule in solution based on the smoothed dissipative particle dynamics method (SDPD) has been presented. In the present paper, we show that the model can be extended to three-dimensional situations and simulate effectively diluted and concentrated polymer solutions. For an isolated suspended polymer, calculated static and dynamic properties agree well with previous numerical studies and theoretical predictions based on the Zimm model. This implies that hydrodynamic interactions are fully developed and correctly reproduced under the current simulated conditions. Simulations of polymer solutions and melts are also performed using a reverse Poiseuille flow setup. The resulting steady rheological properties (viscosity, normal stress coefficients) are extracted from the simulations and the results are compared with the previous numerical studies, showing good results.
D. Azarnykh, S. Litvinov, X. Bian, and N. A. Adams, “Determination of macroscopic transport coefficients of a dissipative particle dynamics solvent,” Physical Review E, vol. 93, no. 1, 2016. Publisher's VersionAbstract
We present an approach to determine macroscopic transport coefficients of a dissipative particle dynamics (DPD) solvent. Shear viscosity, isothermal speed of sound, and bulk viscosity result from DPD-model input parameters and can be determined only a posteriori. For this reason approximate predictions of these quantities are desirable in order to set appropriate DPD input parameters. For the purpose of deriving an improved approximate prediction we analyze the autocorrelation of shear and longitudinal modes in Fourier space of a DPD solvent for Kolmogorov flow. We propose a fitting function with nonexponential properties which gives a good approximation to these autocorrelation functions. Given this fitting function we improve significantly the capability of a priori determination of macroscopic solvent transport coefficients in comparison to previously used exponential fitting functions.
N. Karathanasopoulos and P. Angelikopoulos, “Optimal structural arrangements of multilayer helical assemblies,” International Journal of Solids and Structures, vol. 78-79, pp. 1–8, 2016. Publisher's VersionAbstract

We report a quantitative framework to guide the braiding pattern design of multilayer helical assemblies. We optimize the structural pattern so as to maximize the construction’s resistance to axial loads and concurrently minimize its torsional propensity. To that extent, we consider helical assemblies comprised of up to five layers, for which we identify favorable structural patterns, providing a database that covers most practical applications.

2015
P. E. Hadjidoukas, P. Angelikopoulos, L. Kulakova, C. Papadimitriou, and P. Koumoutsakos, “Exploiting Task-Based Parallelism in Bayesian Uncertainty Quantification,” in Euro-Par 2015: Parallel Processing, Springer, 2015, pp. 532–544. Publisher's Version
D. Rossinelli, et al., “The In-Silico Lab-on-a-Chip: Petascale and High-Throughput Simulations of Microfluidics at Cell Resolution,” in Proceedings of the International Conference for High Performance Computing, Networking, Storage and Analysis – SC \textquotesingle15, 2015, no. Article 2. Publisher's VersionAbstract
We present simulations of blood and cancer cell separation in complex microfluidic channels with subcellular resolu- tion, demonstrating unprecedented time to solution, per- forming at 65.5% of the available 39.4 PetaInstructions/s in the 18, 688 nodes of the Titan supercomputer. These simulations outperform by one to three orders of magnitude the current state of the art in terms of numbers of simulated cells and computational elements. The com- putational setup emulates the conditions and the geometric complexity of microfluidic experiments and our results re- produce the experimental findings. These simulations pro- vide sub-micron resolution while accessing time scales rele- vant to engineering designs. We demonstrate an improvement of up to 45X over com- peting state-of-the-art solvers, thus establishing the frontiers of simulations by particle based methods. Our simulations redefine the role of computational science for the develop- ment of microfluidics – a technology that is becoming as important to medicine as integrated circuits have been to computers
The In-Silico Lab-on-a-Chip: Petascale and High-Throughput Simulations of Microfluidics at Cell Resolution
C. Conti, D. Rossinelli, D. Alexeev, K. Lykov, P. Hadjidoukas, and P. Koumoutsakos, “Video: The In-Silico Lab-on-a-Chip - Catching a Needle in a Flowing Haystack,” in 68th Annual Meeting of the APS Division of Fluid Dynamics - Gallery of Fluid Motion, 2015.
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.
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.
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.
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.
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.
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. 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.
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.
MRAG-i2d: multi-resolution adapted grids for remeshed vortex methods on multicore architectures
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.
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. 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.
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
F. Huhn, W. M. van Rees, M. Gazzola, D. Rossinelli, G. Haller, and P. Koumoutsakos, “Quantitative flow analysis of swimming dynamics with coherent Lagrangian vortices,” Chaos, vol. 25, no. 8, pp. 087405, 2015. Publisher's Version
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.

Pages