Electrostatics of a simple membrane model using Green''s functions formalism.

The electrostatics of a simple membrane model picturing a lipid bilayer as a low dielectric constant slab immersed in a homogeneous medium of high dielectric constant (water) can be accurately computed using the exact Green's functions obtainable for this geometry. We present an extensive discussion of the analysis and numerical aspects of the problem and apply the formalism and algorithms developed to the computation of the energy profiles of a test charge (e.g., ion) across the bilayer and a molecular model of the acetylcholine receptor channel embedded in it. The Green's function approach is a very convenient tool for the computer simulation of ionic transport across membrane channels and other membrane problems where a good and computationally efficient first-order treatment of dielectric polarization effects is crucial.     

Integral equation theories for predicting water structure around molecules

Water plays a crucial role in the structure and function of proteins and other biological macromolecules; thus, theories of aqueous solvation for these molecules are of great importance. However, water is a complex solvent whose properties are still not completely understood. Statistical mechanical integral equation theories predict the density distribution of water molecules around a solute so that all particles are fully represented and thus potentially both molecular and macroscopic properties are included. Here we discuss how several theoretical tools we have developed have been integrated into an integral equation theory designed for globular macromolecular solutes such as proteins. Our approach predicts the three-dimensional spatial and orientational distribution of water molecules around a solute. Beginning with a three-dimensional Ornstein-Zernike equation, a separation is made between a reference part dependent only on the spatial distribution of solvent and a perturbation part dependent also on the orientational distribution of solvent. The spatial part is treated at a molecular level by a modified hypernetted chain closure whereas the orientational part is treated as a Boltzmann prefactor using a quasi-continuum theory we developed for solvation of simple ions. A potential energy function for water molecules is also needed and the sticky dipole models of water, such as our recently developed soft-sticky dipole (SSD) model, are ideal for the proposed separation. Moreover, SSD water is as good as or better than three point models typically used for simulations of biological macromolecules in structural, dielectric and dynamics properties and yet is seven times faster in Monte Carlo and four times faster in molecular dynamics simulations. Since our integral equation theory accurately predicts results from Monte Carlo simulations for solvation of a variety of test cases from a single water or ion to ice-like clusters and ion pairs, the application of this theory to biological macromolecules is promising.     

Analysis and evaluation of channel models: simulations of alamethicin

Alamethicin is an antimicrobial peptide that forms stable channels with well-defined conductance levels. We have used extended molecular dynamics simulations of alamethicin bundles consisting of 4, 5, 6, 7, and 8 helices in a palmitoyl-oleolyl-phosphatidylcholine bilayer to evaluate and analyze channel models and to link the models to the experimentally measured conductance levels. Our results suggest that four helices do not form a stable water-filled channel and might not even form a stable intermediate. The lowest measurable conductance level is likely to correspond to the pentamer. At higher aggregation numbers the bundles become less symmetrical. Water properties inside the different-sized bundles are similar. The hexamer is the most stable model with a stability comparable with simulations based on crystal structures. The simulation was extended from 4 to 20 ns or several times the mean passage time of an ion. Essential dynamics analyses were used to test the hypothesis that correlated motions of the helical bundles account for high-frequency noise observed in open channel measurements. In a 20-ns simulation of a hexameric alamethicin bundle, the main motions are those of individual helices, not of the bundle as a whole. A detailed comparison of simulations using different methods to treat long-range electrostatic interactions (a twin range cutoff, Particle Mesh Ewald, and a twin range cutoff combined with a reaction field correction) shows that water orientation inside the alamethicin channels is sensitive to the algorithms used. In all cases, water ordering due to the protein structure is strong, although the exact profile changes somewhat. Adding an extra 4-nm layer of water only changes the water ordering slightly in the case of particle mesh Ewald, suggesting that periodicity artifacts for this system are not serious.     

Outflow boundary conditions for 3D simulations of non-periodic blood flow and pressure fields in deformable arteries

The simulation of blood flow and pressure in arteries requires outflow boundary conditions that incorporate models of downstream domains. We previously described a coupled multidomain method to couple analytical models of the downstream domains with 3D numerical models of the upstream vasculature. This prior work either included pure resistance boundary conditions or impedance boundary conditions based on assumed periodicity of the solution. However, flow and pressure in arteries are not necessarily periodic in time due to heart rate variability, respiration, complex transitional flow or acute physiological changes. We present herein an approach for prescribing lumped parameter outflow boundary conditions that accommodate transient phenomena. We have applied this method to compute haemodynamic quantities in different physiologically relevant cardiovascular models, including patient-specific examples, to study non-periodic flow phenomena often observed in normal subjects and in patients with acquired or congenital cardiovascular disease. The relevance of using boundary conditions that accommodate transient phenomena compared with boundary conditions that assume periodicity of the solution is discussed.     

Cation transport in nanopores

Using periodic boundary conditions and external electric potential field, we have simulated an ion current flow through a flexible nanopore using cations and an explicit extended simple point charge (SPC/E) water with molecular dynamics simulation. The simulation voltages range goes beyond the usual ionic channel measurements ( ± 1 V) and yields useful information about density profiles, current density distribution and current–voltage (I–V) characteristics.     

Orientational relaxation time of bottom-heavy squirmers in a semi-dilute suspension

One of the important quantities to characterize unsteady behaviour of a cell suspension is the orientational relaxation time, which is the time scale for a micro-organism to re-orientate to its preferred direction from disorientated conditions. In this paper, a swimming micro-organism is modelled as a squirming sphere with prescribed tangential surface velocity, in which the centre of mass of the sphere is displaced from the geometric centre (bottom-heaviness). The orientational relaxation time of bottom-heavy squirmers in a suspension is investigated both analytically and numerically. The three-dimensional movement of 64 identical squirmers in a fluid otherwise at rest, contained in a cube with periodic boundary conditions, is dynamically computed, for random initial positions and orientations. The effects of volume fraction of squirmers, the bottom-heaviness and the squirming mode on the relaxation time are discussed. The results for a semi-dilute suspension show that both the mean stresslet strength and the orientational relaxation time decrease from those for a dilute suspension. We also observe a stress overshoot in some cases. The mechanism for this is different from that for a visco-elastic fluid, and is explained by the change with time of the orientation of squirmers.     

Fluid-structure interaction of deformable aortic prostheses with a bileaflet mechanical valve

Two different aortic prostheses can be used for performing the Bentall procedure: a standard straight graft and the Valsalva graft that better reproduces the aortic root anatomy. The aim of the present work is to study the effect of the graft geometry on the blood flow when a bileaflet mechanical heart valve is used, as well as to evaluate the stress concentration near the suture line where the coronary arteries are connected to graft. An accurate three-dimensional numerical method is proposed, based on the immersed boundary technique. The method accounts for the interactions between the flow and the motion of the rigid leaflets and of the deformable aortic root, under physiological pulsatile conditions. The results show that the graft geometry only slightly influences the leaflets dynamics, while using the Valsalva graft the stress level near the coronary-root anastomoses is about half that obtained using the standard straight graft.     

Accuracy of femur reconstruction from sparse geometric data using a statistical shape model

Sparse geometric information from limited field-of-view medical images is often used to reconstruct the femur in biomechanical models of the hip and knee. However, the full femur geometry is needed to establish boundary conditions such as muscle attachment sites and joint axes which define the orientation of joint loads. Statistical shape models have been used to estimate the geometry of the full femur from varying amounts of sparse geometric information. However, the effect that different amounts of sparse data have on reconstruction accuracy has not been systematically assessed. In this study, we compared shape model and linear scaling reconstruction of the full femur surface from varying proportions of proximal and distal partial femur geometry in combination with morphometric and landmark data. We quantified reconstruction error in terms of surface-to-surface error as well as deviations in the reconstructed femur’s anatomical coordinate system which is important for biomechanical models. Using a partial proximal femur surface, mean shape model-based reconstruction surface error was 1.8 mm with 0.15° or less anatomic axis error, compared to 19.1 mm and 2.7–5.6° for linear scaling. Similar results were found when using a partial distal surface. However, varying amounts of proximal or distal partial surface data had a negligible effect on reconstruction accuracy. Our results show that given an appropriate set of sparse geometric data, a shape model can reconstruct full femur geometry with far greater accuracy than simple scaling.     

Controllers for imposing continuum-to-molecular boundary conditions in arbitrary fluid flow geometries

We present a new parallelised controller for steering an arbitrary geometric region of a molecular dynamics (MD) simulation towards a desired thermodynamic and hydrodynamic state. We show that the controllers may be applied anywhere in the domain to set accurately an initial MD state, or solely at boundary regions to prescribe non-periodic boundary conditions (PBCs) in MD simulations. The mean molecular structure and velocity autocorrelation function remain unchanged (when sampled a few molecular diameters away from the constrained region) when compared with those distributions measured using PBCs. To demonstrate the capability of our new controllers, we apply them as non-PBCs in parallel to a complex MD mixing nano-channel and in a hybrid MD continuum simulation with a complex coupling region. The controller methodology is easily extendable to polyatomic MD fluids.     

Conformation of nifedipine in hydrated 1,2-di-myristoyl-sn-glycero-3-phosphorylcholine bilayer molecular dynamics simulation

The conformation of nifedipine, a cardiac and smooth muscle calcium ion channel antagonist is studied in a hydrated bilayer of forty nine 1,2-di-myristoyl-sn-glycero-3-phosphorylcho-line (DMPC) molecules using molecular dynamics (MD) simulation technique. The simulation was carried out in conditions of constant number, volume and temperature (NVT) at 310 K, which is above the liquid crystalline (Lα) transition temperature of DMPC. The periodic boundary conditions were applied in three-dimensions. Thus the model represented an infinite bilayer. The important geometric parameters characteristic to DMPC and nifedipine molecules were calculated and compared with other theoretical and experimental results pertaining to nifedipine and other related dihydrophyridine (DHP) analogues. Our results suggest that conformational parameters required for antagonist activity are fairly conserved during the interaction of nifedipine with DMPC bilary and bilayer stabilizes the drug conformation in the bioactive form.     

Structure and dynamics of interfacial water in an Lalpha phase lipid bilayer from molecular dynamics simulations

Based on molecular dynamics simulations, an analysis of structure and dynamics is performed on interfacial water at a liquid crystalline dipalmitoylphosphatidycholine/water system. Water properties relevant for understanding NMR relaxation are emphasized. The first and second rank orientational order parameters of the water O-H bonds were calculated, where the second rank order parameter is in agreement with experimental determined quadrupolar splittings. Also, two different interfacial water regions (bound water regions) are revealed with respect to different signs of the second rank order parameter. The water reorientation correlation function reveals a mixture of fast and slow decaying parts. The fast (ps) part of the correlation function is due to local anisotropic water reorientation whereas the much slower part is due to more complicated processes including lateral diffusion along the interface and chemical exchange between free and bound water molecules. The 100-ns-long molecular dynamics simulation at constant pressure (1 atm) and at a temperature of 50 degrees C of 64 lipid molecules and 64 x 23 water molecules lack a slow water reorientation correlation component in the ns time scale. The (2)H(2)O powder spectrum of the dipalmitoylphosphatidycholine/water system is narrow and consequently, the NMR relaxation time T(2) is too short compared to experimental results.     

Model channel ion currents in NaCl-extended simple point charge water solution with applied-field molecular dynamics.

Using periodic boundary conditions and a constant applied field, we have simulated current flow through an 8.125-A internal diameter, rigid, atomistic channel with polar walls in a rigid membrane using explicit ions and extended simple point charge water. Channel and bath currents were computed from 10 10-ns trajectories for each of 10 different conditions of concentration and applied voltage. An electric field was applied uniformly throughout the system to all mobile atoms. On average, the resultant net electric field falls primarily across the membrane channel, as expected for two conductive baths separated by a membrane capacitance. The channel is rarely occupied by more than one ion. Current-voltage relations are concentration dependent and superlinear at high concentrations.     

RV functional imaging: 3-D echo-derived dynamic geometry and flow field simulations

We describe a novel functional imaging approach for quantitative analysis of right ventricular (RV) blood flow patterns in specific experimental animals (or humans) using real-time, three-dimensional (3-D) echocardiography (RT3D). The method is independent of the digital imaging modality used. It comprises three parts. First, a semiautomated segmentation aided by intraluminal contrast medium locates the RV endocardial surface. Second, a geometric scheme for dynamic RV chamber reconstruction applies a time interpolation procedure to the RT3D data to quantify wall geometry and motion at 400 Hz. A volumetric prism method validated the dynamic geometric reconstruction against simultaneous sonomicrometric canine measurements. Finally, the RV endocardial border motion information is used for mesh generation on a computational fluid dynamics solver to simulate development of the early RV diastolic inflow field. Boundary conditions (tessellated endocardial surface nodal velocities) for the solver are directly derived from the endocardial geometry and motion information. The new functional imaging approach may yield important kinematic information on the distribution of instantaneous velocities in the RV diastolic flow field of specific normal or diseased hearts.     

Molecular modeling of p-chlorophenoxyacetic acid binding to the CLC-0 channel

Molecular simulation techniques were applied to predict the interaction of the CLC-0 Cl(-) channel and the channel-blocking molecule p-chlorophenoxyacetic acid (CPA). A three-dimensional model of the CLC-0 channel was constructed on the basis of the homology with the bacterial Cl(-) channel StCLC, the structure of which has been solved by X-ray crystallography. Docking of the CPA molecule was obtained by using a geometric recognition algorithm, yielding 5000 possible conformations. By restraining the simulation to those conformations in which CPA is near the intracellular mouth of the channel, the CPA-protein complex models were reduced to three sets of conformations, which are interconvertible within 2 ns when molecular dynamics is applied to the system. Point mutations of CLC-0 at three different positions predicted to interact with CPA in these configurations did, however, not greatly alter CPA inhibition, suggesting a deeper final binding location. In the model, binding of CPA to a more internal position in the ionic pathway was obtained by applying a constant force vector to CPA, pushing it toward the center of the channel. This technique allowed us to outline the possible intrachannel pathway of CPA and to describe qualitatively the binding sites and energy barriers of this pathway. The consistency of the obtained models and the experimental data indicates that the CLC-0-CPA complex model is reasonable and can be used in further biological studies, such as rational design of blocking agents of and mutagenesis of CLC Cl(-) channels.     

Numerical simulation of local blood flow in the carotid and cerebral arteries under altered gravity

A computational fluid dynamics (CFD) approach was presented to model the blood flows in the carotid bifurcation and the brain arteries under altered gravity. Physical models required for CFD simulation were introduced including a model for arterial wall motion due to fluid-wall interactions, a shear thinning fluid model of blood, a vascular bed model for outflow boundary conditions, and a model for autoregulation mechanism. The three-dimensional unsteady incompressible Navier-Stokes equations coupled with these models were solved iteratively using the pseudocompressibility method and dual time stepping. Gravity source terms were added to the Navier-Stokes equations to take the effect of gravity into account. For the treatment of complex geometry, a chimera overset grid technique was adopted to obtain connectivity between arterial branches. For code validation, computed results were compared with experimental data for both steady-state and time-dependent flows. This computational approach was then applied to blood flows through a realistic carotid bifurcation and two Circle of Willis models, one using an idealized geometry and the other using an anatomical data set. A three-dimensional Circle of Willis configuration was reconstructed from subject-specific magnetic resonance images using an image segmentation method. Through the numerical simulation of blood flow in two model problems, namely, the carotid bifurcation and the brain arteries, it was observed that the altered gravity has considerable effects on arterial contraction/dilatation and consequent changes in flow conditions.     

The Kohonen self-organizing map: a tool for the clustering and alignment of single particles imaged using random conical tilt

An important step in determining the three-dimensional structure of single macromolecules is to bring common features in the images into register through alignment and classification. Here, we took advantage of the striking computational properties of the Kohonen self-organizing map (SOM) to align and classify images of channels obtained by random conical geometry into more homogeneous subsets. First, we used simulations with artificially created images to deduce simple geometrical rules governing the mapping of bounded (differing in size and shape) and unbounded (differing in in-plane orientation) variations in the output plane. Second, we measured the effect of noise on the accuracy of the algorithm to separate homogeneous subsets. Finally, we applied the rules ascertained in the previous steps to separate freeze-fracture images of the cytoplasmic and external domains of the small (approximately 118 kDa) aquaporin-0 water channel. Comparison with the results obtained from a similar input set using alignment-through-classification showed that both methods converged to stable classes exhibiting the same overall shapes (tetragonal and octagonal) for the cytoplasmic and external views of the channel. Processing with the SOM, however, was simplified by the utilization of the geometric rules governing the mapping of bounded and unbounded variations as well as the lack of subjectivity in selecting the reference images during alignment.     

Conformational investigation of a cyclic enterobacterial common antigen employing NMR spectroscopy and molecular dynamics simulations

The three-dimensional structure of a cyclic enterobacterial common antigen (ECA) having four trisaccharide repeating units has been investigated by NMR spectroscopy and molecular dynamics simulations. Three different NMR parameters were determined: (a) (1)H,(1)H cross-relaxation rates from NOE experiments were used for determination of proton-proton distances; (b) trans-glycosidic (3)J(C,H) scalar coupling constants analyzed via a Karplus-type relationship provided information on torsion angles; and (c) (1)H,(13)C one-bond dipolar couplings obtained in a dilute liquid-crystalline medium were interpreted in terms of the orientational order and molecular conformations. The molecular dynamics simulations of the dodecasaccharide were performed with explicit water and counterions, which are important factors that strongly influence molecular conformation. Subsequently, the results from computer simulation were used to generate a three-dimensional structure of the cyclic ECA which is consistent with the experimental NMR parameters.     

Concept and development of an orthotropic FE model of the proximal femur

PURPOSE: In contrast to many isotropic finite-element (FE) models of the femur in literature, it was the object of our study to develop an orthotropic FE "model femur" to realistically simulate three-dimensional bone remodelling. METHODS: The three-dimensional geometry of the proximal femur was reconstructed by CT scans of a pair of cadaveric femurs at equal distances of 2mm. These three-dimensional CT models were implemented into an FE simulation tool. Well-known "density-determined" bony material properties (Young's modulus; Poisson's ratio; ultimate strength in pressure, tension and torsion; shear modulus) were assigned to each FE of the same "CT-density-characterized" volumetric group.In order to fix the principal directions of stiffness in FE areas with the same "density characterization", the cadaveric femurs were cut in 2mm slices in frontal (left femur) and sagittal plane (right femur). Each femoral slice was scanned into a computer-based image processing system. On these images, the principal directions of stiffness of cancellous and cortical bone were determined manually using the orientation of the trabecular structures and the Haversian system. Finally, these geometric data were matched with the "CT-density characterized" three-dimensional femur model. In addition, the time and density-dependent adaptive behaviour of bone remodelling was taken into account by implementation of Carter's criterion. RESULTS: In the constructed "model femur", each FE is characterized by the principal directions of the stiffness and the "CT-density-determined" material properties of cortical and cancellous bone. Thus, on the basis of anatomic data a three-dimensional FE simulation reference model of the proximal femur was realized considering orthotropic conditions of bone behaviour. CONCLUSIONS: With the orthotropic "model femur", the fundamental basis has been formed to realize realistic simulations of the dynamical processes of bone remodelling under different loading conditions or operative procedures (osteotomies, total hip replacements, etc).     

Validation of CFD simulations of cerebral aneurysms with implication of geometric variations

BACKGROUND: Computational fluid dynamics (CFD) simulations using medical-image-based anatomical vascular geometry are now gaining clinical relevance. This study aimed at validating the CFD methodology for studying cerebral aneurysms by using particle image velocimetry (PIV) measurements, with a focus on the effects of small geometric variations in aneurysm models on the flow dynamics obtained with CFD. METHOD OF APPROACH: An experimental phantom was fabricated out of silicone elastomer to best mimic a spherical aneurysm model. PIV measurements were obtained from the phantom and compared with the CFD results from an ideal spherical aneurysm model (S1). These measurements were also compared with CFD results, based on the geometry reconstructed from three-dimensional images of the experimental phantom. We further performed CFD analysis on two geometric variations, S2 and S3, of the phantom to investigate the effects of small geometric variations on the aneurysmal flow field. Results. We found poor agreement between the CFD results from the ideal spherical aneurysm model and the PIV measurements from the phantom, including inconsistent secondary flow patterns. The CFD results based on the actual phantom geometry, however, matched well with the PIV measurements. CFD of models S2 and S3 produced qualitatively similar flow fields to that of the phantom but quantitatively significant changes in key hemodynamic parameters such as vorticity, positive circulation, and wall shear stress. CONCLUSION: CFD simulation results can closely match experimental measurements as long as both are performed on the same model geometry. Small geometric variations on the aneurysm model can significantly alter the flow-field and key hemodynamic parameters. Since medical images are subjected to geometric uncertainties, image-based patient-specific CFD results must be carefully scrutinized before providing clinical feedback.