Simulate Bovine Serum Albumin (BSA) confined in a silica nanopore and compare dynamics to results from QENS experiments. The goal of the experiments is to investigate the dynamics of single protein cluster under ultra-confinement using BASIS



Eugene Mamontov
Kevin Slater
Panchao Yin


Details of the experiment

IPTS-18245 - proposal document

Two samples:

  • DBSA: encapsulated deuterated water and hydrogenated protein.
  • HBSA: dry protein powder, containing no water (or just crystallographic) and hydrogenated protein.

Susceptibility plots:


Comparison of MSD (label says 'Elastic Intensity' but it's wrong):


Maximum MSD corresponds to ~0.4A, very minimal flucuation. Thus, there is very little deuterated water surrounding the protein in the silica-encapsulated DBSA sample.

Slides that Panchao sent me providing more info on the samples:

The last slide shows the relative weights of BSA (66054.4 amu per monomer), silica (60.08 amu per SiO2 molecule), and water (18.015 amu per molecule), corresponding to 1545 SiO2 and 622 water molecules for each BSA monomer.




Some references regarding parameterization of force-fields describing amorphous silica

References for mean square displacements in bulk materials containing silica

Simulations of proteins interacting with silica


Gromacs tutorial for solvated Lyzosyme


. *.top: (ASCII) The topology contains all the information necessary to define the molecule within a simulation. Thi- information includes nonbonded parameters (atom types and charges) as well as bonded parameters (bonds, angles, and dihedrals). Does not contain coordinates for the atoms.

  • *.gro: (ASCII) GROMACS-formatted structure file that contains residue and atom names of our system with the nomenclature of the specified force field. Also coordinates.

  • *.itp: (ASCII) Indicates a harmonic restraint of specified force constant to be applied for certain atoms of the simulation.

  • *.mdp: (ASCII) Run configuration parameters.

  • *.tpr: (binary) Run input. Integrates information from the .top, .gro, and *.mdp files.

  • *.log: (ASCII) standard output of the Gromacs run.

  • *.edr: (binary) values for various thermodynamic properties of the system versus time.

  • *.trr: (binary) trajectory (coordinates, and optionally velocities) versus time, in full-precision.

  • *.cpt: (binary) Checkpoint containing all the state variables, including coordinates and velocities, that allow one to continue a simulation.

  • *.xvg: (ASCII) Output of gmx energy, a file to be read by the plotting tool Xmgrace.

  • *.ndx: (ASCII) index containing user-defined sets of atoms.


Solvated Monomeric BSA box

In subdirectory simulations/monomer/solvatedBSA.

Obtaining a working initial PDB file:
Structure 3v03.pdb is a dimer. We select the second chain:

We clean up a bit the PDB file before obtaining the preprocessed Gromacs file. - Removing crystallographic water Oxygen atoms. - Removing residue "ACT". - Submitting the PDB file to [pdb_hydro](http://lorentz.dynstr.pasteur.fr/pdb_hydro.php) server for generation of missing sidechains. We obtain structure [3v03B_noH.pdb](simulation/monomer/solvatedBSA/3v03B_noH.pdb) and log file [3v03B_PDB_Hydro.stdout](simulation/monomer/solvatedBSA/3v03B_PDB_Hydro.stdout). There are missing hydrogen atoms. - Removed residues "CA" with residue numbers 584 to 586. These are not part of the BSA protein. - Added hydrogens with command: gmx pdb2gmx -f 3v03B_noH.pdb -o 3v03B.pdb -p 3v03B.top -i 3v03B.itp -water spce -ff oplsaa &> 3v03B_pdb2gmx.log to obtain [3v03B.pdb](simulation/monomer/solvatedBSA/3v03B.pdb), among other output files. Solvating the protein with water: Use OPLSAA force field and SPCE water model: gmx pdb2gmx -f 3v03B.pdb -o 3v03B.gro -p 3v03B.top -i 3v03B.itp -water spce -ff oplsaa &> 3v03B_pdb2gmx.log outputs: - 3v03B_pdb2gmx.log states 579 residues, 9173 atoms, and a total charge of -16e. We need to place 16 Na atoms in order to neutralize the system. - 3v03B.gro processed binary gromacs file. - 3v03B.itp positional restraints file. - 3v03B.top topology of the isolated molecule. The below command centers the protein in the box (-c), and places it at least 1.0 nm from the box edge (-d 1.0). The box type is defined as a cube (-bt cubic): gmx editconf -f 3v03B.gro -o 3v03B_newbox.gro -c -d 1.0 -bt cubic Fill with solvent (water). Solvation is accomplished using solvate: cp 3v03B.top 3v03B_solv.top gmx solvate -cp 3v03B_newbox.gro -cs spc216.gro -o 3v03B_solv.gro -p 3v03B_solv.top -scale 0.7 Add ions ([ions.mdp](simulation/monomer/solvatedBSA/ions.mdp)): cp 3v03B_solv.top 3v03B_solv_ions.top gmx grompp -f ions.mdp -c 3v03B_solv.gro -p 3v03B_solv_ions.top -o ions.tpr gmx genion -s ions.tpr -o 3v03B_solv_ions.gro -p 3v03B_solv_ions.top -neutral # select group number 13. Minimize the protein: (In subdirectory simulation/monomer/solvatedBSA/minimize) Minimize options specified in file [minimize.mdp](simulation/monomer/solvatedBSA/minimize/minimize.mdp). gmx grompp -f minimize.mdp -c 3v03B_solv_ions.gro -p 3v03B_solv_ions.top -o em.tpr # creates binary input em.tpr. gmx mdrun -ntmpi 16 -ntomp 1 -pin on -deffnm em Relax the water around the protein: (In subdirectory simulation/monomer/solvatedBSA/minimize/relax) 50ps NVT run using restraints on the protein atoms. Coordinates defining the equilibrium positions for the restraints are in file 3v03B.itp. This file is encoded in the topology file 3v03B_solv_ions.top. Options for the run specified in file [relax_1.mdp](simulation/monomer/solvatedBSA/minimize/relax/relax_1.mdp). gmx grompp -f relax_1.mdp -c em.gro -p 3v03B_solv_ions.top -o relax_1.tpr gmx mdrun -ntmpi 16 -ntomp 1 -pin on -deffnm relax_1 # use 16 MPI threads, pin each thread to a diffeent core. Continue with a 500ps NPT run using restraints on the protein atoms. The goal is to shrink the simulation box to tipycal densities. Input file is [relax_2.mdp](simulation/monomer/solvatedBSA/minimize/relax/relax_2.mdp). gmx grompp -f relax_2.mdp -c relax_1.gro -t relax_1.cpt -p 3v03B_solv_ions.top -o relax_2.tpr time gmx mdrun -ntmpi 16 -ntomp 1 -pin on -deffnm relax_2 Recover some thermodynamics of the cell shrinkage: gmx energy -f relax_2.edr -o relax_2.xvg
[relax_2.agr](simulation/monomer/solvatedBSA/minimize/relax/relax_2.agr) [relax_2.gif](simulation/monomer/solvatedBSA/minimize/relax/relax_2.gif), [relax_2.vmd](simulation/monomer/solvatedBSA/minimize/relax/relax_2.vmd)

Center protein in the simulation box:
gmx trjconv -f relax_2.cpt -s relax_2.tpr -pbc mol -ur compact -center -o relax_2_centered.gro -boxcenter tric

  • Select group for centering: Group 1 ( Protein) has 9173 elements
  • Select group for output: Group 0 ( System) has 84738 elements
    Produces file relax_2_centered.gro.

Select closest N water molecules to protein:
The mass of the protein is 66054.497 a.m.u according to 3v03B_pdb2gmx.log. Thus a hydration level $h=0.1$ corresponds to 6605.4497 a.m.u of water, corresponding to 6605.4497/18.01528 = 366 water molecules.

In the experiment, the estimated hydration level is 0.17. Thus, we have to include 1.7*366 = 622 molecules.

Save the following selection to file selection_h0.17_waters.txt:

close_solvent_atoms = group "SOL" and within 0.1874 of group "Protein"; close_solvent_molecules = same resindex as close_solvent_atoms; group "Protein" or close_solvent_molecules;

We found a cutoff = 0.1874 to retain the 622 closest water molecules, simply by counting the number of lines containing "SOL" in file stripped_system.pdb (see below).

Create index file stripped_system_h0.17.ndx with protein plus close waters, and also PDB file stripped_system_h0.17.pdb gmx select -f relax_2_centered.gro -s relax_2_centered.gro -sf selection_h0.17_waters.txt -on stripped_system_h0.17.ndx -ofpdb stripped_system_h0.17.pdb -pdbatoms selected


Amorphous silica box from Materials Studio post

Issue #4

Silica $(SiO_2)_N$ has a sort-range order, whereby a Si atom is located at the center of a tetrahedron, whose vertices are defined by four O atoms.

In subdirectory silica/amorphous_from_md/ and also Materials Studio project silica in virtual machine LEIOAVIRTUAL.
I downloaded 21A box from a post in the Materials Studio community].

Since the maximum diameter of BSA is ~77A, we use SiO2_21A_3d_NoConnect.pdb to create a 5x5x5 supercell of size ~107A. This precludes the protein from "feeling itself" via the periodic boundary conditions.

Inserting the protein in the silica box:

We start with file relax_2_centered.gro from directory simulation/monomer/solvatedBSA/minimize/relax, which is contained in a cubic box of 9.41445nm. We have to insert the protein in the silica box SiO2_107A_3d_NoConnect.pdb, which is a cubic box of 10.6974nm. Thus, we resize the box of the protein plus water:

gmx editconf -f relax_2_centered.gro -o relax_2_centered_resized.gro -bt cubic -box 10.6974


Next we select the 622 water molecules closest to the protein as we did before:

gmx select -f relax_2_centered_resized.gro -s relax_2_centered_resized.gro -sf selection_h0.17_waters.txt -on stripped_system_h0.17.ndx -ofpdb stripped_system_h0.17.pdb -pdbatoms selected


We load stripped_system_h0.17.pdb and SiO2_107A_3d_NoConnect.pdb in VMD to check the protein is in the middle of the silica box:

We need to "clean up" files stripped_system_h0.17.pdb and SiO2_107A_3d_NoConnect.pdb before merging protein and silica into a single PDB. We use merge_silica_and_protein.py for this purpose

python merge_silica_and_protein.py creates the following files:

Removing $SiO_2$ molecules overlapping the protein and water


Simulations of small cristobalite crystal

In subdirectory simulation/silica/cristobalite/cristobalite_2_2_2

With Materials Studio and

Small cristobalite crystal containing 96 atoms

Create topology with x2top utitily

We have modified oplsa.ff/ force field in simulation/silica/cristobalite/cristobalite_2_2_2/oplsa.ff/ accounting for silicon type SIO and oxygen type OSI present in the crystal. We modified files ffnonbonded.itp, ffbonded.ipt, and atomtypes.atp.

We create topology file with gmx x2top. This utitility requires description of connectivity between atom types in file oplsa.f/atomname2type.n2t.

gmx x2top -v -f cristobalite_2_2_2.gro -o cristobalite_2_2_2.top -ff oplsaa -name SIL -pbc -noparam

After this, a minimization (minin) followed by equilibration at either constant volume (relax_1) or constant pressure (relax_2) are carried out. Run relax_2 is continued with another NPT run (relax_3) which can be used to gather data. Script gmx_commands.sh containst the relevant gmx submissions.

Relaxation at constant pressure relax_2.mdp:

  • 0.1fs time step.
  • verlet-buffer-tolerance=-1 to use nstlist=1 and rlist.
  • tau_p=100 for slow relaxation of the volume

Continuation at constant pressure relax_3.mdp:

  • gmx grompp -f relax_3.mdp -t relax_2.cpt -c relax_2.tpr -o relax_3.tpr # checkpoint required
  • 1fs time step
  • verlet-buffer-tolerance=-1 to use nstlist=10



From the Crystallography Open Database we search for the crystal structure of alpha-cristobalite, one of the morphologies of $SiO_2$:


In Materials studio, the supercell replicating the unit cell 21 times along X and Y and 15 times along Z is able to surround BSA.


A post in materials studio prescribes a recipe to build amorphous silica starting with the amorphous cell module. We can take parts of this recipe to create amorphous silica from a cristobalite supercell.

From: Jason DeJoannis
I think we built that structure in partnership with Corning. It is documented in this paper:
Vyas, Dickinson and Armstrong-Poston, Mol Sim, 32, 2, 135 (2006).

Build bulk glass with Amorphous Cell
Pack cell with SiO2 and Si2O4 fragments at density of 2g/mL
cvff_aug forcefield with Ewald
Anneal the bulk
NVT at 7000K for 20ps
NPT for 20ps
Reduce T in 1000K increments with 20ps at each temperature down to 300K
Create the surface
Cleave and add vacuum slab 50-60A
Extend surface supercell as needed
Anneal the surface
Start at 600K and cool to 300K in steps of 100


Protein confined in cristobalite

In Materials studio, the supercell replicating the unit cell 21 times along X and Y and 15 times along Z is able to surround BSA.


Inserting the protein in the silica box:

We start with file relax_2_centered.gro from directory simulation/monomer/solvatedBSA/minimize/relax, which is contained in a cubic box of 9.41445nm. We have to insert the protein in the cristobalite box cristobalite_21_21_15.pdb, which has dimensions 10.4406x10.4406x10.3835 nm. Thus, we resize the box of the protein plus water:

gmx editconf -f relax_2_centered.gro -o relax_2_centered_resized.gro -bt triclinic -box 10.4406 10.4406 10.3835


Next we select the 622 water molecules closest to the protein as we did before:

gmx select -f relax_2_centered_resized.gro -s relax_2_centered_resized.gro -sf selection_h0.17_waters.txt -on stripped_system_h0.17.ndx -ofpdb stripped_system_h0.17.pdb -pdbatoms selected


We need to "clean up" files stripped_system_h0.17.pdb and cristobalite_21_21_15.pdb before merging protein and silica into a single PDB. We use merge_silica_and_protein.py for this purpose. python merge_silica_and_protein.py creates the following files:


Removing $SiO_2$ molecules overlapping the protein and water

The cristobalite lattice projected onto the XY plane can be viewed as a distorted squared lattice with the Si atoms occupying the corners, and the oxygen atoms located in the middle of the square sides. See below:


With this arrangement, we can univoquely associate two oxygens to each Si, out of the four first-neighbors. We chose the oxygen with the largest X-component of its position relative to the Si, and the oxygen with the largest Y-component. This way we can paritition the cristobalite crystal in a unique way.

Script carve_silica.py uses silicanet.py takes the previous partition to create a hole in the cristobalite lattice. It takes as input file SiO2_protein_0.17.pdb and partitions the cristobalite lattice into $SiO_2$ molecules. Molecules having any of its atoms overlapping the protein+water are removed. The overlapp parameter is tried with different values, from 1.0 to 2.5 Angstroms, producing different carved files. Also the names of the silica atoms will match those of the force field in subdirectory simulation/silica/cristobalite/confineBSA/poretop/confinedBSA.ff/:


From the PDB files, it seems that overlapp > 1.5 leaves too many $SiO_2$ molecules and steric clashes that will not be resolved by a minimization run.

Simulating the confining cristobalite without bonds
As a first approximation, we assume each cristobalite atom is independent of the others, and it's subjected to a harmonic restraint. For cristobalite atoms in the interface of the protein+water, we can use as spring constant the value extracted from the experimental MSD of the dry protein (~2E-06 nm^2/K). For the bulk cristobalite atoms located away from the protein+water, we can use the spring constant derived from the simulations of the small cristobalite (~7E-07 nm^2/K).

in directory simulation/silica/cristobalite/confineBSA/poretop/nobonds we have modified force field oplsaa.ff/ containing nonbonded interactions for the cristobalite atom types SIO, OSI, and OA. We modified the following files:

We use SiO2carved_ovl1.5_protein_0.17.pdb and create a suitable .gro file with script adaptPDB.py.

python adaptPDB.py # produces confinedBSA_0.gro

Create topology and restraint files:

gmx pdb2gmx -f confinedBSA_0.gro -o confinedBSA.gro -p confinedBSA.top -i confinedBSA.itp -water spce -ff oplsaa &> confinedBSA.log

produces the following files:

The total charge of the system is reported to be -342.900. To make the system neutral, we increase the positive charge of the SIO atoms. There are 23694 SIO atoms in the silica topology file confinedBSA_Other2.itp. Thus, the new charge will be 1.1 + 342.9/23694 = 1.1144720182324641. We edit file confinedBSA_Other2.itp and make this change.

We also updated the SIO charge in the following files of the force field:

  • ions.itp
  • residues.rtp (this is the file from which the silica charges are read by pdb2gmx)
  • ffnonbonded.itp

We renamed Other2 to Silica, modifying confinedBSA.top and producing producing new files confinedBSA_Silica.itp and confinedBSA_Silica_pr.itp

Topology for the carved cristobalite

script evaluate_bonding.py uses silicanet.py to calculate average bond and angle distances in the cyrstobalite using the cristobalite in SiO2_protein_0.17.pdb. The results: bond_Si_O = 1.57833041748 angle_O_Si_O = 109.149737561 angle_Si_O_Si = 146.498743845

Working directory: simulation/silica/cristobalite/confineBSA/poretop/yesbonds

We use Zhang's force field with the previously found bond and angle average quantities. Resulting force field in subdirectory simulation/silica/cristobalite/confineBSA/poretop/yesbonds/oplsaa.ff/.

From Simulating the confining cristobalite without bonds, we now that the charge of SIO should be 1.1144720182324641 to render a neutral silica+protein+water system. The charges are thus:

  • SIO +1.1144720182324641
  • OSI -0.55
  • OA -0.90
  • HO +0.40

If we only want to simulate the silica, then a charge of SIO equal to 1.113796741791171 will render a neutral silica sytem.

To create a topology file for the carved cristobalite (without the protein) we edit SiO2carved_ovl1.5_protein_0.17.pdb, removing the protein and water atoms, and then saving as SiO2carved_ovl1.5.pdb. Then, script topology_silica.py will input this file and use module silicatop.py to generate an include topology file. we produce SiO2carved_ovl1.5_spw0.itp for a neutral silica+protein+water system, and SiO2carved_ovl1.5_s0.itp for a neutral silica system.

We simulate the silica alone in directory simulation/silica/cristobalite/confineBSA/poretop/yesbonds/onlysilica/ by copying SiO2carved_ovl1.5_s0.itp into SIL.itp. We had to tweak the last SIO with a slightly different charge (1.11505874179 instead of 1.11379674179) so that gmx grommp would report a neutral system. This is a numerical inaccuracy of gmx grommp. Also we had to comment the info in the [ defaults ] section because this info is previously set up in oplsaa.ff/forcefield.itp.

With script shift_indexes.py we shift the residue and atom indexes of SiO2carved_ovl1.5.pdb to create confiningSIL.pdb. This file starts with atom number one and has residue number equal to one.

We manually created confiningSIL.top and minimize.mdp

gmx editconf -f confiningSIL.pdb -o confiningSIL.gro gmx grommp -f minimize.mdp -c confiningSIL.gro -p confiningSIL.top -o minimize.tpr gmx mdrun -ntmpi 16 -ntomp 1 -pin on -deffnm minimize


The minimization movie shows an expanding hole, causing a deformation in the crystal structure close to the hole surface. The expansion is due to the net negative charge caused by the dangling oxygens. These oxygens repel each other, and they quickly shield from each other by "turning" towards the bulk crystal. There are some $SiO_2$ molecules in the hole flying around. The expansion due to the negative charge will be less pronounced once we include the protein because it will act as a shielding insulator, but the steric clashes between protein and surface with induce an expansion on the hole.

Inserting the protein in the carved crystobalite with bonds

In directory simulation/silica/cristobalite/confineBSA/poretop/yesbonds/wholesystem/ we copied SiO2carved_ovl1.5_spw0.itp into SIL.itp and modified this file as follows:

  • commented out the [ defaults ] section
  • changed slightly the charge of the last SIO atom (from 1.11447201823 to 1.11468801823). gmx grompp complained the system had no net zero charge.

We copied file SiO2carved_ovl1.5_protein_0.17.pdb into confinedBSA.pdb

We created topology file confinedBSA.top by hand, retrieving all info for the protein and water from previous topology file of the solvated BSA, 3v03B_solv_ions.top

Input run files minimize.mdp and relax_1.mdp

gmx editconf -f confinedBSA.pdb -o confinedBSA.gro gmx grompp -f minimize.mdp -c confinedBSA.gro -p confinedBSA.top -o minimize.tpr gmx mdrun -ntmpi 16 -ntomp 1 -pin on -deffnm minimize gmx grompp -f relax_1.mdp -c minimize.gro -p confinedBSA.top -o relax_1.tpr gmx mdrun -ntmpi 16 -ntomp 1 -pin on -deffnm relax_1 # 4ns in the NVT gmx grompp -f relax_2.mdp -c relax_1.gro -p confinedBSA.top -o relax_2.tpr gmx mdrun -ntmpi 16 -ntomp 1 -pin on -deffnm relax_2 # 2ns in the NPT to relax the box size gmx grompp -c relax_2.tpr -t relax_2.cpt -p confinedBSA.top -f equil_1.mdp -o equil_1.tpr # continuation but in the NVT gmx mdrun -ntmpi 16 -ntomp 1 -pin on -deffnm equil_1


The 2ns relaxation at NPT (relax_2) proved insufficiently long since the pressure became negative in the equil_1 run. Thus, we extend the simulation for 10ns in Titan by appending 1ns jobs. For this we need to transfer to Titan all files relax_2.* and the following

In Titan, we just run the substitution script by passing as argument the number of the firs run:
./relax_2.submitnext.sh 0
Will trigger submission of the first nanosecond, followed by 9 more runs up to a total of 10ns

After relaxation, we run an equilibrium run for 1ns in the NVT ensemble, and then append subsequent runs for a total of 12ns:


Comparing confinements and OA charges

We simulate four systems in order to find out which system provides us with the best comparison to experimental MSD at 300K:

  • $q_{OA}=0.90$ and silica-system overlap 1.5A
  • $q_{OA}=0.55$ and silica-system overlap 1.5A
  • $q_{OA}=0.55$ and silica-system overlap 2.5A
  • $q_{OA}=0.90$ and silica-system overlap 2.5A

We already simulated the first system. The other three systems we create and simulate in subdirectories simulation/confinedBSA/h017_q0.55_ovl1.5, simulation/confinedBSA/h017_q0.55_ovl2.5, and simulation/confinedBSA/h017_q0.90_ovl2.5, respectively.

With gmx msd we calculate MSD on the equil_1.xtc trajectories


The system $q_{OA}=0.90$ and silica-system overlap 2.5A approaches the most the experimental point. We select this system for the temperature scan.


Protein confined in cristobalite with $q_{OA}=0.90$ and overlap=2.5A

Files confinedBSA.gro and confinedBSA.pdb are same as for system $q_{OA}=0.55$ and overlap=2.5A. We use topology_silica.py to create SiO2carved_q0.90_ovl2.5.itp which we rename to SIL.itp.

As usual, we ran gmx grompp -f minimize.mdp -c confinedBSA.gro -p confinedBSA.top -o minimize.tpr and use the reported excess charge to edit file SIL.itp and there tune the charge of the SIO atoms to render a neutral system.

We run the following jobs at 300K in Titan:

Temperatures scan within subdirectory simulation/confinedBSA/h017_q0.90_ovl2.5/Tscan

We run the temperatures scan in Titan from $T_i=40$ to $T_e=420$ in increments of $20K$. Preparation of input files prior to submitting jobs to the queue done by running bash_commands.sh.

The following jobs:


The INTERFACE force field


Latest version 1.5 in subdirectory simulation/forcefields/INTERFACE_FF_1_5/FORCE_FIELDS. Contains the force-field in PCFF and CHARMM format, plus a series of bulk and surface silica systems.


Zhang's force field

Prof. Yingchun Liu sent me a force-field for silica after a personal request. This after reading Zheng16. This paper referenced Wensink00 paper.

Located in subdirectory simulation/forcefields/Zheng16. Contains bonded and nonbonded parameters.


Protocol for initial configuration for confinedBSA

Issue #6

A possible series of steps:

  • Solvate protein with water, and equilibrate the system. Center protein in the simulation box.
  • Strip all but the closest $N$ water molecules, thus selecting a particular hydration level. This is the hydratedBSA system.
  • Separately, simulate a system of silica to create amorphous silica. Dimensions of this system should be the same as the solvated equilibrated sytem. This annealing simulation entails melting the silica network followed by slow cooling.
  • Overlap the hydratedBSA system onto the silica system, then remove overlapping $SiO_2$ molecules. This is the initial confinedBSA system. This system has a broken silica network due to unpaired electrons of the $Si$ atoms close to the water+protein.
  • Fix the silica network with an annealing process in which the water+protein atoms are held frozen.
  • Freeze the silica network and equilibrate the water + protein. This is the initial and equilibrated confinedBSA system.
  • Run equilibrium simulations of the confinedBSA system. The silica network is always kept frozen so simulations are in the NVT ensemble.


HTML and Markdown Syntax Examples

local link: link</br> remote link: http://ambermd.org font face="courier new"
$$S_{model}(Q,E)=A(Q)\cdot S_{elastic}(E) + B(Q)\cdot S_{simulation}(Q,E)\otimes S_{elastic}(E) + C(Q)+D(Q)\cdot E$$

 Quoted text 

image caption
some text

In [ ]: