Volume 1, Issue 2013 -

Fractional derivative or fractional calculus plays a significant role in theoretical modeling of scientific and engineering problems. However, only relatively low order fractional derivatives are used at present. In general, it is not obvious what role a high fractional derivative can play and how to make use of arbitrarily high-order fractional derivatives. This work introduces arbitrarily high-order fractional partial differential equations (PDEs) to describe fractional hyperdiffusions. The fractional PDEs are constructed via fractional variational principle. A fast fractional Fourier transform (FFFT) is proposed to numerically integrate the high-order fractional PDEs so as to avoid stringent stability constraints in solving high-order evolution PDEs. The proposed high-order fractional PDEs are applied to the surface generation of proteins. We first validate the proposed method with a variety of test examples in two and three-dimensional settings. The impact of high-order fractional derivatives to surface analysis is examined. We also construct fractional PDE transform based on arbitrarily high-order fractional PDEs. We demonstrate that the use of arbitrarily high-order derivatives gives rise to time-frequency localization, the control of the spectral distribution, and the regulation of the spatial resolution in the fractional PDE transform. Consequently, the fractional PDE transform enables the mode decomposition of images, signals, and surfaces. The effect of the propagation time on the quality of resulting molecular surfaces is also studied. Computational efficiency of the present surface generation method is compared with the MSMS approach in Cartesian representation. We further validate the present method by examining some benchmark indicators of macromolecular surfaces, i.e., surface area, surface enclosed volume, surface electrostatic potential and solvation free energy. Extensive numerical experiments and comparison with an established surface model indicate that the proposed high-order fractional PDEs are robust, stable and efficient for biomolecular surface generation.

A general approach was proposed in this article to develop high-order exponentially fitted basis functions for finite element approximations of multi-dimensional drift-diffusion equations for modeling biomolecular electrodiffusion processes. Such methods are highly desirable for achieving numerical stability and efficiency. We found that by utilizing the one-to-one correspondence between the continuous piecewise polynomial space of degree k + 1 and the divergencefree vector space of degree k , one can construct high-order two-dimensional exponentially fitted basis functions that are strictly interpolative at a selected node set but are discontinuous on edges in general, spanning nonconforming finite element spaces. First order convergence was proved for the methods constructed from divergence-free Raviart-Thomas space RT 0 0 at two different node sets.

This review outlines the recent progress made in developing more accurate and efficient solutions to model electrostatics in systems comprised of bio-macromolecules and nanoobjects, the last one referring to objects that do not have biological function themselves but nowadays are frequently used in biophysical and medical approaches in conjunction with bio-macromolecules. The problem of modeling macromolecular electrostatics is reviewed from two different angles: as a mathematical task provided the specific definition of the system to be modeled and as a physical problem aiming to better capture the phenomena occurring in the real experiments. In addition, specific attention is paid to methods to extend the capabilities of the existing solvers to model large systems toward applications of calculations of the electrostatic potential and energies in molecular motors, mitochondria complex, photosynthetic machinery and systems involving large nanoobjects.

A stochastic generalized Born (GB) solver is presented which can give predictions of energies arbitrarily close to those that would be given by exact effective GB radii, and, unlike analytical GB solvers, these errors are Gaussian with estimates that can be easily obtained from the algorithm. This method was tested by computing the electrostatic solvation energies (Δ G solv ) and the electrostatic binding energies (Δ G bind ) of a set of DNA-drug complexes, a set of protein-drug complexes, a set of protein-protein complexes, and a set of RNA-peptide complexes. Its predictions of Δ G solv agree with those of the linearized Poisson-Boltzmann equation, but it does not predict Δ G bind well, although these predictions of Δ G bind may be marginally better than those of traditional analytical GB solvers. Apparently, the GB model itself must be improved before accurate estimates of Δ G bind can be obtained.

Zernike polynomials have been widely used in the description and shape retrieval of 3D objects. These orthonormal polynomials allow for efficient description and reconstruction of objects that can be scaled to fit within the unit ball. However, maps defined within box-shaped regions ¶ for example, rectangular prisms or cubes ¶ are not well suited to representation by Zernike polynomials, because these functions are not orthogonal over such regions. In particular, the representations require many expansion terms to describe object features along the edges and corners of the region. We overcome this problem by applying a Gram-Schmidt process to re-orthogonalize the Zernike polynomials so that they recover the orthonormality property over a specified box-shaped domain. We compare the shape retrieval performance of these new polynomial bases to that of the classical Zernike unit-ball polynomials.

rithms for solving the 3D electro-diffusion equations such as the Poisson-Nernst-Planck equations and the size-modified Poisson-Nernst-Planck equations in simulations of biomolecular systems in ionic liquid. A set of transformation methods based on the generalized Slotboom variables is used to solve the coupled equations. Calculations of the diffusion-reaction rate coefficients, electrostatic potential and ion concentrations for various systems verify the method’s validity and stability. The iterations between the Poisson equation and the Nernst- Planck equations in the primitive method and in the transformation method are compared to illustrate how the new method accelerates the convergence of the solution. To speed up the convergence, we introduce the DIIS (direct inversion of the iterative subspace) method including Simple Mixing and Anderson Mixing as under-relaxation techniques, the effectiveness of which on acceleration is shown by numerical tests. It is worth noting that the primitive method fails to solve the size-modified Poisson-Nernst-Planck equations for real protein systems but the transformation method succeeds in the simulations of the ACh-AChE reaction system and the DNA fragment. To improve the accuracy of the solution, we introduce high order elements and mesh adaptation based on an a posteriori error estimator. Numerical results indicate that our mesh adaptation process leads to quasi-optimal convergence. We implement our algorithms using the parallel adaptive finite element package PHG [53] and high parallel efficiency is obtained.

The Poisson-Boltzmann (PB) model is an effective approach for the electrostatics analysis of solvated biomolecules. The nonlinearity associated with the PB equation is critical when the underlying electrostatic potential is strong, but is extremely difficult to solve numerically. In this paper, we construct two operator splitting alternating direction implicit (ADI) schemes to efficiently and stably solve the nonlinear PB equation in a pseudo-transient continuation approach. The operator splitting framework enables an analytical integration of the nonlinear term that suppresses the nonlinear instability. A standard finite difference scheme weighted by piecewise dielectric constants varying across the molecular surface is employed to discretize the nonhomogeneous diffusion term of the nonlinear PB equation, and yields tridiagonal matrices in the Douglas and Douglas-Rachford type ADI schemes. The proposed time splitting ADI schemes are different from all existing pseudo-transient continuation approaches for solving the classical nonlinear PB equation in the sense that they are fully implicit. In a numerical benchmark example, the steady state solutions of the fully-implicit ADI schemes based on different initial values all converge to the time invariant analytical solution, while those of the explicit Euler and semi-implicit ADI schemes blow up when the magnitude of the initial solution is large. For the solvation analysis in applications to real biomolecules with various sizes, the time stability of the proposed ADI schemes can be maintained even using very large time increments, demonstrating the efficiency and stability of the present methods for biomolecular simulation.

We analyze and suggest improvements to a recently developed approximate continuum-electrostatic model for proteins. The model, called BIBEE/I (boundary-integral based electrostatics estimation with interpolation), was able to estimate electrostatic solvation free energies to within a mean unsigned error of 4% on a test set of more than 600 proteins¶a significant improvement over previous BIBEE models. In this work, we tested the BIBEE/I model for its capability to predict residue-by-residue interactions in protein–protein binding, using the widely studied model system of trypsin and bovine pancreatic trypsin inhibitor (BPTI). Finding that the BIBEE/I model performs surprisingly less well in this task than simpler BIBEE models, we seek to explain this behavior in terms of the models’ differing spectral approximations of the exact boundary-integral operator. Calculations of analytically solvable systems (spheres and tri-axial ellipsoids) suggest two possibilities for improvement. The first is a modified BIBEE/I approach that captures the asymptotic eigenvalue limit correctly, and the second involves the dipole and quadrupole modes for ellipsoidal approximations of protein geometries. Our analysis suggests that fast, rigorous approximate models derived from reduced-basis approximation of boundaryintegral equations might reach unprecedented accuracy, if the dipole and quadrupole modes can be captured quickly for general shapes.

Poisson and Poisson-Boltzmann equations (PE and PBE) are widely used in molecular modeling to estimate the electrostatic contribution to the free energy of a system. In such applications, PE often needs to be solved multiple times for a large number of system configurations. This can rapidly become a highly demanding computational task. To accelerate such calculations we implemented a graphical processing unit (GPU) PE solver described in this work. The GPU solver performance is compared to that of our central processing unit (CPU) implementation of the solver. During the performance analysis the following three characteristics were studied: (1) precision associated with the modeled system discretization on the grid, (2) numeric precision associated with the floating point representation of real numbers (this is done via comparison of calculations with single precision (SP) and double precision (DP)), and (3) execution time. Two types of example calculations were carried out to evaluate the solver performance: (1) solvation energy of a single ion and a small protein (lysozyme), and (2) a single ion potential in a large ion-channel (α-hemolysin). In addition, influence of various boundary condition (BC) choices was analyzed, to determine the most appropriate BC for the systems that include a membrane, typically represented by a slab with the dielectric constant of low value. The implemented GPU PE solver is overall about 7 times faster than the CPU-based version (including all four cores). Therefore, a single computer equipped with multiple GPUs can offer a computational power comparable to that of a small cluster. Our calculations showed that DP versions of CPU and GPU solvers provide nearly identical results. SP versions of the solvers have very similar behavior: in the grid scale range of 1-4 grids/Å the difference between SP and DP versions is less than the difference stemming from the system discretization. We found that for the membrane protein, the use of a focusing technique with periodic boundary conditions in rough grid provides significantly better results than using a focusing technique with the electric potential set to zero at the boundaries.

Multi-scale modeling plays an important role in understanding the structure and biological functionalities of large biomolecular complexes. In this paper, we present an efficient computational framework to construct multi-scale models from atomic resolution data in the Protein Data Bank (PDB), which is accelerated by multi-core CPU and programmable Graphics Processing Units (GPU). A multi-level summation of Gaussian kernel functions is employed to generate implicit models for biomolecules. The coefficients in the summation are designed as functions of the structure indices, which specify the structures at a certain level and enable a local resolution control on the biomolecular surface. A method called neighboring search is adopted to locate the grid points close to the expected biomolecular surface, and reduce the number of grids to be analyzed. For a specific grid point, a KD-tree or bounding volume hierarchy is applied to search for the atoms contributing to its density computation, and faraway atoms are ignored due to the decay of Gaussian kernel functions. In addition to density map construction, three modes are also employed and compared during mesh generation and quality improvement to generate high quality tetrahedral meshes: CPU sequential, multi-core CPU parallel and GPU parallel. We have applied our algorithm to several large proteins and obtained good results.