@DuanChen:2012b] and multiscale analysis [@Wei:2009], in conjugation with other physical models, such as electrostatics, elasticity and molecular mechanics [@Wei:2013]. Geometry modeling and annotation of biomolecular surfaces together with physical features, such as electrostatics and lipophilicity, provide some of the best predictions of biomolecular solvation free energies [@BaoWang:2016FFTS; @BaoWang:2016HPK], protein-drug binding affinities [@BaoWang:2016FFTB], protein mutation energy changes [@ZXCang:2016a] and protein-protein interaction hot spots [@Darnell:2008; @Demerdash:2009].

Theoretical modeling of the structure-function relationship of biomolecules is usually based on fundamental laws of physics, i.e., quantum mechanics (QM), molecular mechanism (MM), continuum mechanics, statistical mechanics, thermodynamics, etc. QM methods are indispensable for chemical reactions and protein degradations [@Warshel:1976; @Cui:2002; @YZhang:2009a]. Molecular dynamics (MD) [@McCammon:1977] is a powerful tool for the understanding of the biomolecular conformational landscapes and elucidating collective motion and fluctuation. Currently, MD is a main workhorse in molecular biophysics for biomolecular modeling and simulation. However, all-electron or all-atom representations and long-time integrations lead to such an excessively large number of degrees of freedom that their application to real-time large-scale dynamics of large proteins or multiprotein complexes becomes prohibitively expensive. For instance, current computer simulations of protein folding take many months to come up with a very poor copy of what Nature administers 