High performance computing can play a catalytic role in advancing physics, biology and medicine. We focus on the development of new computational models and corresponding numerical methods suitable for the next generation of supercomputers.Our research interests lie in the area of multiscale/multiphysics modeling and parallel large-scale simulations of physical and biological systems.
Research Projects and Areas
Our research involves fields requiring expertise in mathematics, biology, physics and computer science.
Some of the past and current projects and corresponding numerical methods are illustrated in the image above.
Selected projects are described in more details below.
On the image: Snapshot from the simulation of the whole RBC cell stretched under the flow with a line attachment to the microfluidic channel wall. The inset is the local view of the bilayer-cytoskeletal interaction. The bilayer is shown as gray surface while the cytoskeleton is shown as a triangular network with the contour. The contour is the magnitude of the bilayer-cytoskeletal interaction force on a junctional complex.
We developed a systematic coarse-graining procedure for modeling red blood cells (RBCs) using arguments based on mean-field theory.This general coarse-graining procedure, published in the PRL, does not employ any fitting parameters. Our spectrin-based RBC model lies between continuum and atomic scales and can be used for arbitrary levels of coarse-graining. The model takes into account bending and in-plane shear energy, viscous effects of the membrane, and constraints of total area and volume. We performed systematic fully three-dimensional computational simulations as well as microfluidic experiments to quantitatively study the flow dynamics of the RBCs at the smallest relevant scale and different physiological conditions. Recently, we have extended this model by representing the lipid bilayer and the cytoskeleton of spectrin network as two distinct components. We employed this new model to investigate the effects of the bilayer-cytoskeletal viscoelastic interactions. Currently there are no experimental techniques that directly measure the mechanical characteristics of these interactions. By applying this new two-component whole-cell model, we reconciled and resolved several controversies and issues in RBC mechanics. Our computational framework provides a broad general methodology for meso-scale simulations of the flow of cells. In addition, it can be used to explore important problems involving cell physiology and pathological states mediated by protein mutations, such as the bilayer loss in hereditary spherocytosis and the bilayer–cytoskeleton uncoupling in sickle cell anemia.
Large scale high performance simulations
On the image: Snapshot from flow simulations in microfluidic device.
In collaboration with the groups of Petros Koumoutsakos (ETH Zurich) and George Karniadakis (Brown University) we developed an efficient highly scalable framework for simulations of cells in flow based on coarse-grained red blood cell (RBC) model and Dissipative Particle Dynamics (DPD) method. We performed large-scale simulations of cell flow in several microfluidic devices which involved solvent, RBCs, white blood cells (WBCs) and circulating tumor cells (CTCs) in complex geometries with volumes of up to 150mm3 and containing up to 1.43 Billion deforming cells. To the best of our knowledge, these are the biggest simulations of flow with deformable cells to date, which were selected as the finalist of the 2015 Gordon Bell Prize Award in high performance computing.
Related publications: CS15_paper
On the image: Three organisms (Olga, Igor and Antoine) working on bioleaching of copper minerals hundreds of meters below the Earth's surface. This environmentally friendly method to extract metals from ores reduces the release of toxic compounds associated with traditional mining techniques. Our goal is to define the key factors affecting biofilm formation and allow optimisation of biofilm growth rates as a function of microbial composition and genetic profile, thereby making this technique more attractive to mining companies.
Bioleaching of copper minerals such as chalcopyrite (the largest copper resource in the word) is usually done in engineered heaps and accounts for approximately 15% of the present world copper production. Bioleaching is an environmentally friendly method to extract metals from ores as it reduces the release of toxic compounds associated with traditional mining techniques. Biofilm formation on the surface of copper containing sulphide minerals is a vital stage in the biotechnological process of biomining. Our goals is to study biofilm formation by some of the most important moderately thermophilic bioleaching bacteria, Acidithiobacillus caldus, Leptospirillum ferriphilum and Sulfobacillus thermosulfidooxidans. In collaboration with several experimental research groups and two companies actively carrying out biomining and other biotechnological applications, we plan to utilize large "omics" data sets for extensive bioinformatic analysis, aimed at the formulation of multi-species biological network models. Specifically, an ODE-based computational module is being developed to reverse engineer metatranscriptomics data (RNAseq) and quantitative metaproteomics data from different microbial mixtures employed in stirred reactor experiments to obtain interaction networks of the principal molecular players. The network model will be integrated into particle-based model of biofilm growth at the single-cell level. The computational framework will be validated using experimental data from confocal laser scanning microscopy and epifluorescence microscopy of biofilms. The in-silico system will be used to define the key factors affecting biofilm formation and will allow optimisation of biofilm growth rates as a function of microbial composition and genetic profile. This will allow to reduce industrial operational costs, thereby making this environmentally friendly technique more attractive to mining companies.
Coarse-grained protein model
On the image: Snapshot from the folding simulations of the WW-domain using coarse-grained polarizable protein model. Comparison of least square fitted backward-mapped structures from the simulations (cartoon) to the experimental structures (ribbon).
Computer simulations on an atomistic scale are suitable for studying biomolecular systems and can be used complementary to experiments for the exploration of free energy landscapes and binding free energies. However, at present all-atom approaches are limited in their timescale and can be computationally expensive. Coarse-graining of proteins ranges over different length and time-scales and several methods have been used for the description of these biomolecular systems, including Go-models in protein folding, where the native structure is biased using restraint potentials. Very recently, polarizability effects have been incorporated into some coarse-grained forcefields. In addition, charge-charge interactions have been added for coarse-grained protein simulations to optimize the stability of secondary and tertiary structure formation. We developed a new coarse-grained forcefield for the description of proteins on the basis of Dissipative Particle Dynamics method, thus further extending the scope of applications of DPD to simulations of biomolecular systems. Our new model is based on the electrostatic polarization of the protein backbone and a detailed representation of the sidechains in combination with a polarizable water model. We defined our model parameters using the experimental structures of two proteins, TrpZip2 and TrpCage. Backmapping and subsequent replica-exchange molecular dynamics runs verified our approach and showed convergence to the experimental structures on the atomistic level. We validated our model on five different proteins: GB1, the WW-domain, the B-domain of Protein A, the Peripheral binding subunit and Villin headpiece. We emphasize that our model is free from elastic network restraints needed for instance in the popular MARTINI coarse-grained model to keep the secondary structure stable, which we consider as a strong advantage of our new forcefield. We believe that our model is suitable for the coarse-grained description of proteins and has the potential to improve sampling of native states in coarse-grained protein simulations.
Related publications: PCCP_paper, JChemP_paper
Transport processes in tumor induced microcirculation
On the image: Atomistic model of the glycocalyx, snapshot from NAMD simulations. Sugar chains are shown in green licorice representation. Transmembrane proteins are pictured as red surfaces, lipid bilayer as van der Waals beads. Water and ions are not shown.
Over 85% of human cancers are solid tumors. In order to grow tumors require supply of oxygen and nutrients. Therefore the blood flow in tumor vasculature network is a key regulator of the tumor development in the vascular phase. The effectiveness of cancer therapy in solid tumors depends on the adequate delivery of the therapeutic agent to tumor cells. Chemotherapeutic drugs are often administered systemically. In this case, drug delivery to cells in tumor involves several processes, including transport within a vascular network and transport across vessel walls. In collaboration with the group of Petros Koumoutsakos from ETH Zurich, we study the transport processes in healthy and tumor induced microcirculation. We build on our combined expertise in modeling of RBC dynamics in flow using experimentally validated multiscale models, on the development of continuum, discrete and hybrid angiogenesis models, and on the capability to implement highly-scalable particle codes on emergent computer architectures. Within the project, we recently built an atomistic model of the glycocalyx. To the best of our knowledge this is the most detailed model to date. The glycocalyx is a polymeric layer that is the first barrier for solute exchange between the blood and tissues and as such it is involved in key metabolic processes. It is now well established that information about molecular motions within the glycocalyx is necessary to advance our knowledge of its function. Current experimental techniques are not yet sensitive enough to elucidate the detailed dynamic interplay among the different glycocalyx components. Here, MD simulations act as the unique exploratory tool, not only to complement experimental findings but also to propose testable predictions. Using computational code NAMD, we have studied glycocalyx conformational and dynamic properties, albumin binding to it and the mechanics of the glycocalyx under the shear flow. We envision that the present study will enhance our understanding of the relative importance of phenomena associated with transport processes in the vasculature. The results will help to quantify transport phenomena in healthy and tumor induced microcirculation thus contributing to the development of rational strategies for cancer therapy.
Related publications: AnnRev paper
Stochastic modeling of platelet aggregation
On the image: Thrombus growing on a blood vessel wall in our computer simulations. Flow structure interaction illustrated by streamlines close to developing thrombus. For clarity most platelets not interacting with thrombi are omitted.
We developed a computational stochastic model of platelet aggregation in 3D blood flows by accounting for the movements of all platelets individually involved. The model incorporates information about platelet adhesion molecule behavior, cell mechanics, and fibrinogen formation. To make simulations of thousands platelets tractable, we employed spectal/hp element discretization of the unsteady Stokes equations in combination with the Force Coupling Method (FCM), in which platelets are represented as force envelopes based on a spatial distribution of finite force multipoles. Our simulations demonstrated the dependence of thrombus growth rate on blood velocity as found experimentally by Begent & Born (Nature 1970).Thrombus growth rates are affected by the velocity of the blood flow, but do not simply increase with it. As the velocity increases, the growth rates exhibit a maximum, and subsequently decrease. In addition, we found that thrombus growth rate is enhanced by modest pulsatility but this effect is less pronounced when pulsations are amplified due in part to more embolization. These large-scale simulations required the development of efficient, highly optimized parallel algorithms for high performance computers.
Coarse-grained stochastic molecular dynamics simulations
On the image: Comparison of density, velocity, temperature and stress profiles for Poiseuille flow in MD and corresponding coarse-grained DPD simulations. Adaptive boundary conditions are used for controlling density fluctuations near solid walls.
Dissipative particle dynamics (DPD) is a relatively new, potentially very effective approach in simulating mesoscale hydrodynamics. It can simulate efficiently complex liquids and dense suspensions using only a few thousand virtual particles and at speed-up factors of more than one hundred thousand compared to molecular dynamics. The DPD model consists of particles which correspond to coarse-grained entities, thus representing molecular clusters rather than individual atoms. DPD can be thought of as a coarse-grained version of Molecular Dynamics (MD), but it employs dissipative and stochastic forces to account for the eliminated degrees of freedom. Similar to the effort that has been going on with the lattice Boltzmann method (LBM), another mesoscopic simulation technique, a systematic verification and validation of DPD is required to evaluate its accuracy, efficiency and robustness. In a series of papers we analyzed the fundamental modeling ideas of DPD. Unlike an MD simulation where the choice of potential is based on a theoretical model of the physical system to be simulated, a DPD simulation involves potentials of a form independent of the physical system. We proposed a process of choosing the DPD parameters and determining the DPD length and time scales for different levels of coarse-graining such that the DPD simulations correspond to an MD simulation of a Lennard-Jones (LJ) liquid. We investigated and identified the limits of coarse-graining procedure in DPD. Unlike the MD method, the soft repulsion between DPD particles cannot prevent fluid particles from penetrating solid boundaries, and thus extra effort is required to impose accurately the no-slip wall boundary condition. We developed an adaptive method to impose no-slip conditions and to control anomalous density fluctuations near the solid walls in DPD simulations. The formalism of this method is not restricted to DPD but it could remedy a similar pathology in other particle-based simulations.