From f14c919a8f2b8e4ab759dbec7bfc81cedd497d8a Mon Sep 17 00:00:00 2001
From: "Kurt A. O'Hearn" <ohearnku@msu.edu>
Date: Mon, 11 Jan 2021 14:07:12 -0500
Subject: [PATCH] sPuReMD: add initial QM/MM support with work for integrating
with Amber.
---
configure.ac | 10 +
sPuReMD/src/bonds.c | 16 +
sPuReMD/src/charges.c | 25 +
sPuReMD/src/control.c | 141 +++--
sPuReMD/src/control.h | 5 +
sPuReMD/src/forces.c | 4 +-
sPuReMD/src/grid.c | 4 +
sPuReMD/src/hydrogen_bonds.c | 15 +-
sPuReMD/src/init_md.c | 9 +
sPuReMD/src/lin_alg.c | 12 +
sPuReMD/src/multi_body.c | 14 +
sPuReMD/src/nonbonded.c | 56 +-
sPuReMD/src/reax_types.h | 16 +-
sPuReMD/src/spuremd.c | 75 +--
sPuReMD/src/spuremd.h | 6 +
sPuReMD/src/torsion_angles.c | 31 +
sPuReMD/src/valence_angles.c | 26 +
sPuReMD/src/vector.h | 17 +
test/library/qmmm_amber/AVE/control | 75 +++
test/library/qmmm_amber/AVE/fort.3 | 781 ++++++++++++++++++++++++++
test/library/qmmm_amber/AVE/fort.4 | 387 +++++++++++++
test/library/qmmm_amber/control | 69 +++
test/library/qmmm_amber/tester_AVE.py | 520 +++++++++++++++++
23 files changed, 2203 insertions(+), 111 deletions(-)
create mode 100644 test/library/qmmm_amber/AVE/control
create mode 100644 test/library/qmmm_amber/AVE/fort.3
create mode 100644 test/library/qmmm_amber/AVE/fort.4
create mode 100644 test/library/qmmm_amber/control
create mode 100644 test/library/qmmm_amber/tester_AVE.py
diff --git a/configure.ac b/configure.ac
index ba8b3afa..18133ec5 100644
--- a/configure.ac
+++ b/configure.ac
@@ -61,6 +61,12 @@ AC_ARG_ENABLE([puremd-standalone],
[enable build for standalone PuReMD code @<:@default: yes@:>@])],
[puremd_standalone=${enableval}], [puremd_standalone=yes])
+# Build code in QM/MM mode.
+AC_ARG_ENABLE([qmmm],
+ [AS_HELP_STRING([--enable-qmmm],
+ [enable build for code in QM/MM mode @<:@default: no@:>@])],
+ [qmmm=${enableval}], [qmmm=no])
+
# Build LAMMPS/reaxc integration code.
AC_ARG_ENABLE([lammps-reaxc],
[AS_HELP_STRING([--enable-lammps-reaxc],
@@ -319,6 +325,10 @@ if test "x${pack_serial_enabled}" = "xyes" || test "x${pack_openmp_enabled}" = "
fi
AC_SUBST([G_LIBS], ["${GSL_LIBS}"])
+ # Build code in QM/MM mode
+ AS_IF([test "x${qmmm}" = "xyes"],
+ [AC_DEFINE([QMMM], [1], [Define to 1 to build PuReMD code in QMMM mode.])])
+
AC_LANG_POP([C])
fi
AM_CONDITIONAL([BUILD_S_OMP], [test "x${pack_serial_enabled}" = "xyes" || test "x${pack_openmp_enabled}" = "xyes"])
diff --git a/sPuReMD/src/bonds.c b/sPuReMD/src/bonds.c
index a257fc43..2e61e261 100644
--- a/sPuReMD/src/bonds.c
+++ b/sPuReMD/src/bonds.c
@@ -60,6 +60,10 @@ void Bonds( reax_system *system, control_params *control,
#endif
for ( i = 0; i < system->N; ++i )
{
+#if defined(QMMM)
+ if ( system->atoms[i].qmmm_mask == TRUE )
+ {
+#endif
start_i = Start_Index(i, bonds);
end_i = End_Index(i, bonds);
@@ -68,6 +72,12 @@ void Bonds( reax_system *system, control_params *control,
if ( i < bonds->bond_list[pj].nbr )
{
j = bonds->bond_list[pj].nbr;
+
+#if defined(QMMM)
+ if ( system->atoms[j].qmmm_mask == TRUE )
+ {
+#endif
+
type_i = system->atoms[i].type;
type_j = system->atoms[j].type;
sbp_i = &system->reax_param.sbp[type_i];
@@ -156,8 +166,14 @@ void Bonds( reax_system *system, control_params *control,
#endif
}
}
+#if defined(QMMM)
+ }
+#endif
}
}
+#if defined(QMMM)
+ }
+#endif
}
}
diff --git a/sPuReMD/src/charges.c b/sPuReMD/src/charges.c
index 4695942e..3d8b8021 100644
--- a/sPuReMD/src/charges.c
+++ b/sPuReMD/src/charges.c
@@ -1390,6 +1390,15 @@ static void QEq( reax_system * const system, control_params * const control,
break;
}
+#if defined(QMMM)
+ for ( int i = 0; i < system->N; ++i )
+ {
+ workspace->s[0][i] = system->atoms[i].q_init;
+ workspace->t[0][i] = system->atoms[i].q_init;
+ workspace->mask_qmmm[i] = system->atoms[i].qmmm_mask;
+ }
+#endif
+
switch ( control->cm_solver_type )
{
case GMRES_S:
@@ -1502,6 +1511,14 @@ static void EE( reax_system * const system, control_params * const control,
break;
}
+#if defined(QMMM)
+ for ( int i = 0; i < system->N; ++i )
+ {
+ workspace->s[0][i] = system->atoms[i].q_init;
+ workspace->mask_qmmm[i] = system->atoms[i].qmmm_mask;
+ }
+#endif
+
switch ( control->cm_solver_type )
{
case GMRES_S:
@@ -1613,6 +1630,14 @@ static void ACKS2( reax_system * const system, control_params * const control,
#undef SIZE
#endif
+#if defined(QMMM)
+ for ( int i = 0; i < system->N; ++i )
+ {
+ workspace->s[0][i] = system->atoms[i].q_init;
+ workspace->mask_qmmm[i] = system->atoms[i].qmmm_mask;
+ }
+#endif
+
switch ( control->cm_solver_type )
{
case GMRES_S:
diff --git a/sPuReMD/src/control.c b/sPuReMD/src/control.c
index 8019737a..db62f2a5 100644
--- a/sPuReMD/src/control.c
+++ b/sPuReMD/src/control.c
@@ -430,55 +430,68 @@ int Set_Control_Parameter( const char * const keyword,
}
-void Read_Control_File( const char * const control_file, reax_system * const system,
+/* Assign default values to control parameters
+ *
+ * NOTE: force field file parameters must be parsed before this
+ * */
+void Set_Control_Defaults( reax_system * const system,
control_params * const control, output_controls * const out_control )
{
- char *s, **tmp;
- int c, i, ret;
- FILE *fp;
-
- fp = sfopen( control_file, "r" );
-
- assert( fp != NULL );
+ int i;
- /* assign default values */
strncpy( control->sim_name, "default.sim", sizeof(control->sim_name) - 1 );
control->sim_name[sizeof(control->sim_name) - 1] = '\0';
+ strncpy( control->restart_from, "default.res", sizeof(control->restart_from) - 1 );
+ control->restart_from[sizeof(control->restart_from) - 1] = '\0';
control->restart = FALSE;
out_control->restart_format = WRITE_BINARY;
out_control->restart_freq = 0;
- strncpy( control->restart_from, "default.res", sizeof(control->restart_from) - 1 );
- control->restart_from[sizeof(control->restart_from) - 1] = '\0';
out_control->restart_freq = 0;
- control->random_vel = FALSE;
+ control->random_vel = FALSE;
control->reposition_atoms = 0;
-
control->ensemble = NVE;
control->nsteps = 0;
- control->dt = 0.25;
-
- control->geo_format = PDB;
- control->restrict_bonds = FALSE;
-
control->periodic_boundaries = TRUE;
-
+ control->restrict_bonds = FALSE;
+ control->tabulate = 0;
+ control->dt = 0.25;
control->reneighbor = 1;
- /* interaction cutoffs from force field global paramters */
- control->bo_cut = 0.01 * system->reax_param.gp.l[29];
- control->nonb_low = system->reax_param.gp.l[11];
- control->nonb_cut = system->reax_param.gp.l[12];
-
/* defaults values for other cutoffs */
- control->vlist_cut = control->nonb_cut + 2.5;
+ control->vlist_cut = system->reax_param.gp.l[12] + 2.5;
control->bond_cut = 5.0;
- control->bg_cut = 0.3;
+ control->nonb_low = system->reax_param.gp.l[11];
+ control->nonb_cut = system->reax_param.gp.l[12];
+ control->bo_cut = 0.01 * system->reax_param.gp.l[29];
control->thb_cut = 0.001;
control->hbond_cut = 0.0;
- control->tabulate = 0;
+ control->T_init = 0.0;
+ control->T_final = 300.0;
+ control->Tau_T = 1.0;
+ control->T_mode = 0;
+ control->T_rate = 1.0;
+ control->T_freq = 1.0;
+ control->P[0] = 0.000101325;
+ control->P[1] = 0.000101325;
+ control->P[2] = 0.000101325;
+ control->Tau_P[0] = 500.0;
+ control->Tau_P[1] = 500.0;
+ control->Tau_P[2] = 500.0;
+ control->Tau_PT = 500.0;
+ control->compressibility = 1.0;
+ control->press_mode = 0;
+ control->compute_pressure = FALSE;
+
+ control->remove_CoM_vel = 25;
+ control->geo_format = PDB;
+ control->dipole_anal = 0;
+ control->freq_dipole_anal = 0;
+ control->diffusion_coef = 0;
+ control->freq_diffusion_coef = 0;
+ control->restrict_type = 0;
control->charge_method = QEQ_CM;
control->cm_q_net = 0.0;
@@ -492,7 +505,8 @@ void Read_Control_File( const char * const control_file, reax_system * const sys
control->cm_init_guess_extrap1 = 3;
control->cm_init_guess_extrap2 = 2;
/* assign default values */
- strncpy( control->cm_init_guess_gd_model, "frozen_model.pb", sizeof(control->cm_init_guess_gd_model) - 1 );
+ strncpy( control->cm_init_guess_gd_model, "frozen_model.pb",
+ sizeof(control->cm_init_guess_gd_model) - 1 );
control->cm_init_guess_gd_model[sizeof(control->cm_init_guess_gd_model) - 1] = '\0';
control->cm_init_guess_win_size = 5;
control->cm_solver_pre_comp_type = JACOBI_PC;
@@ -503,55 +517,40 @@ void Read_Control_File( const char * const control_file, reax_system * const sys
control->cm_solver_pre_app_type = TRI_SOLVE_PA;
control->cm_solver_pre_app_jacobi_iters = 50;
- control->T_init = 0.0;
- control->T_final = 300.0;
- control->Tau_T = 1.0;
- control->T_mode = 0.0;
- control->T_rate = 1.0;
- control->T_freq = 1.0;
-
- control->P[0] = 0.000101325;
- control->P[1] = 0.000101325;
- control->P[2] = 0.000101325;
- control->Tau_P[0] = 500.0;
- control->Tau_P[1] = 500.0;
- control->Tau_P[2] = 500.0;
- control->Tau_PT = 500.0;
- control->compressibility = 1.0;
- control->press_mode = 0;
- control->compute_pressure = FALSE;
-
- control->remove_CoM_vel = 25;
+ control->molec_anal = NO_ANALYSIS;
+ control->freq_molec_anal = 0;
+ control->bg_cut = 0.3;
+ control->num_ignored = 0;
+ for ( i = 0; i < MAX_ATOM_TYPES; i++ )
+ {
+ control->ignore[i] = 0;
+ }
out_control->log_update_freq = 0;
-
out_control->write_steps = 0;
out_control->traj_compress = 0;
out_control->write = &fprintf;
out_control->traj_format = 0;
out_control->write_header = &Write_Custom_Header;
out_control->append_traj_frame = &Append_Custom_Frame;
-
strncpy( out_control->traj_title, "default_title", sizeof(out_control->traj_title) - 1 );
out_control->traj_title[sizeof(out_control->traj_title) - 1] = '\0';
out_control->atom_format = 0;
out_control->bond_info = 0;
out_control->angle_info = 0;
+}
- control->molec_anal = NO_ANALYSIS;
- control->freq_molec_anal = 0;
- control->num_ignored = 0;
- for ( i = 0; i < MAX_ATOM_TYPES; i++ )
- {
- control->ignore[i] = 0;
- }
- control->dipole_anal = 0;
- control->freq_dipole_anal = 0;
+void Read_Control_File( const char * const control_file, reax_system * const system,
+ control_params * const control, output_controls * const out_control )
+{
+ char *s, **tmp;
+ int c, i, ret;
+ FILE *fp;
- control->diffusion_coef = 0;
- control->freq_diffusion_coef = 0;
- control->restrict_type = 0;
+ fp = sfopen( control_file, "r" );
+
+ assert( fp != NULL );
if ( fp != NULL )
{
@@ -594,6 +593,15 @@ void Read_Control_File( const char * const control_file, reax_system * const sys
sfree( s, "Read_Control_File::s" );
}
+ sfclose( fp, "Read_Control_File::fp" );
+}
+
+
+void Set_Control_Derived_Values( reax_system * const system,
+ control_params * const control )
+{
+ int i;
+
/* determine target T */
if ( control->T_mode == 0 )
{
@@ -625,15 +633,4 @@ void Read_Control_File( const char * const control_file, reax_system * const sys
{
system->g.spread[i] = 2;
}
-
-#if defined(DEBUG_FOCUS)
- fprintf( stderr,
- "en=%d steps=%d dt=%.5f opt=%d T=%.5f P=%.5f %.5f %.5f\n",
- control->ensemble, control->nsteps, control->dt, control->tabulate,
- control->T, control->P[0], control->P[1], control->P[2] );
-
- fprintf( stderr, "control file read\n" );
-#endif
-
- sfclose( fp, "Read_Control_File::fp" );
}
diff --git a/sPuReMD/src/control.h b/sPuReMD/src/control.h
index 13ae7da0..75d74562 100644
--- a/sPuReMD/src/control.h
+++ b/sPuReMD/src/control.h
@@ -27,7 +27,12 @@
int Set_Control_Parameter( const char * const, const char ** const,
control_params * const, output_controls * const );
+void Set_Control_Defaults( reax_system * const, control_params * const,
+ output_controls * const );
+
void Read_Control_File( const char * const, reax_system * const, control_params * const,
output_controls * const );
+void Set_Control_Derived_Values( reax_system * const, control_params * const );
+
#endif
diff --git a/sPuReMD/src/forces.c b/sPuReMD/src/forces.c
index 3b6afe5b..73179c2b 100644
--- a/sPuReMD/src/forces.c
+++ b/sPuReMD/src/forces.c
@@ -1209,8 +1209,8 @@ static void Init_Forces_Tab( reax_system *system, control_params *control,
}
/* hydrogen bond lists */
- if ( control->hbond_cut > 0
- &&(ihb == H_ATOM || ihb == H_BONDING_ATOM)
+ if ( control->hbond_cut > 0.0
+ && (ihb == H_ATOM || ihb == H_BONDING_ATOM)
&& nbr_pj->d <= control->hbond_cut )
{
jhb = sbp_j->p_hbond;
diff --git a/sPuReMD/src/grid.c b/sPuReMD/src/grid.c
index 72790ddc..25f9591c 100644
--- a/sPuReMD/src/grid.c
+++ b/sPuReMD/src/grid.c
@@ -684,6 +684,10 @@ static inline void reax_atom_Copy( reax_atom *dest, reax_atom *src )
rvec_Copy( dest->v, src->v );
rvec_Copy( dest->f, src->f );
dest->q = src->q;
+#if defined(QMMM)
+ dest->qmmm_mask = src->qmmm_mask;
+ dest->q_init = src->q_init;
+#endif
}
diff --git a/sPuReMD/src/hydrogen_bonds.c b/sPuReMD/src/hydrogen_bonds.c
index dd9a65cd..ccbff7f6 100644
--- a/sPuReMD/src/hydrogen_bonds.c
+++ b/sPuReMD/src/hydrogen_bonds.c
@@ -86,7 +86,11 @@ void Hydrogen_Bonds( reax_system *system, control_params *control,
for ( j = 0; j < system->N; ++j )
{
/* j must be a hydrogen atom */
- if ( system->reax_param.sbp[system->atoms[j].type].p_hbond == H_ATOM )
+ if ( system->reax_param.sbp[system->atoms[j].type].p_hbond == H_ATOM
+#if defined(QMMM)
+ && system->atoms[j].qmmm_mask == TRUE
+#endif
+ )
{
type_j = system->atoms[j].type;
start_j = Start_Index( j, bonds );
@@ -123,6 +127,12 @@ void Hydrogen_Bonds( reax_system *system, control_params *control,
for ( pk = hb_start_j; pk < hb_end_j; ++pk )
{
k = hbond_list[pk].nbr;
+
+#if defined(QMMM)
+ if ( system->atoms[k].qmmm_mask == TRUE )
+ {
+#endif
+
type_k = system->atoms[k].type;
nbr_jk = hbond_list[pk].ptr;
r_jk = nbr_jk->d;
@@ -304,6 +314,9 @@ void Hydrogen_Bonds( reax_system *system, control_params *control,
#endif
}
}
+#if defined(QMMM)
+ }
+#endif
}
}
}
diff --git a/sPuReMD/src/init_md.c b/sPuReMD/src/init_md.c
index 0b13e1bc..e850ce32 100644
--- a/sPuReMD/src/init_md.c
+++ b/sPuReMD/src/init_md.c
@@ -683,6 +683,11 @@ static void Init_Workspace( reax_system *system, control_params *control,
workspace->perm_ilutp = NULL;
}
+#if defined(QMMM)
+ workspace->mask_qmmm = smalloc( system->N_max * sizeof( int ),
+ "Init_Workspace::workspace->mask_qmmm" );
+#endif
+
/* integrator storage */
workspace->a = smalloc( system->N_max * sizeof( rvec ),
"Init_Workspace::workspace->a" );
@@ -1558,6 +1563,10 @@ static void Finalize_Workspace( reax_system *system, control_params *control,
sfree( workspace->perm_ilutp, "Finalize_Workspace::workspace->perm_ilutp" );
}
+#if defined(QMMM)
+ sfree( workspace->mask_qmmm, "Init_Workspace::workspace->mask_qmmm" );
+#endif
+
/* integrator storage */
sfree( workspace->a, "Finalize_Workspace::workspace->a" );
sfree( workspace->f_old, "Finalize_Workspace::workspace->f_old" );
diff --git a/sPuReMD/src/lin_alg.c b/sPuReMD/src/lin_alg.c
index 53f736e3..84ab600f 100644
--- a/sPuReMD/src/lin_alg.c
+++ b/sPuReMD/src/lin_alg.c
@@ -3527,6 +3527,9 @@ int CG( const static_storage * const workspace, const control_params * const con
t_start = Get_Time( );
Vector_Sum( r, 1.0, b, -1.0, d, N );
+#if defined(QMMM)
+ Vector_Mask_qmmm( r, workspace->mask_qmmm, N );
+#endif
rnorm = Norm( r, N );
t_vops += Get_Timing_Info( t_start );
@@ -3537,6 +3540,9 @@ int CG( const static_storage * const workspace, const control_params * const con
t_start = Get_Time( );
Vector_Copy( p, z, N );
+#if defined(QMMM)
+ Vector_Mask_qmmm( p, workspace->mask_qmmm, N );
+#endif
sig_new = Dot( r, p, N );
t_vops += Get_Timing_Info( t_start );
@@ -3551,6 +3557,9 @@ int CG( const static_storage * const workspace, const control_params * const con
alpha = sig_new / tmp;
Vector_Add( x, alpha, p, N );
Vector_Add( r, -1.0 * alpha, d, N );
+#if defined(QMMM)
+ Vector_Mask_qmmm( r, workspace->mask_qmmm, N );
+#endif
rnorm = Norm( r, N );
t_vops += Get_Timing_Info( t_start );
@@ -3564,6 +3573,9 @@ int CG( const static_storage * const workspace, const control_params * const con
sig_new = Dot( r, z, N );
beta = sig_new / sig_old;
Vector_Sum( p, 1.0, z, beta, p, N );
+#if defined(QMMM)
+ Vector_Mask_qmmm( p, workspace->mask_qmmm, N );
+#endif
t_vops += Get_Timing_Info( t_start );
}
diff --git a/sPuReMD/src/multi_body.c b/sPuReMD/src/multi_body.c
index 1dd047df..b9b1e67a 100644
--- a/sPuReMD/src/multi_body.c
+++ b/sPuReMD/src/multi_body.c
@@ -59,6 +59,10 @@ void Atom_Energy( reax_system *system, control_params *control,
for ( i = 0; i < system->N; ++i )
{
+#if defined(QMMM)
+ if ( system->atoms[i].qmmm_mask == TRUE )
+ {
+#endif
/* set the parameter pointer */
type_i = system->atoms[i].type;
sbp_i = &system->reax_param.sbp[ type_i ];
@@ -137,10 +141,17 @@ void Atom_Energy( reax_system *system, control_params *control,
}
}
}
+#if defined(QMMM)
+ }
+#endif
}
for ( i = 0; i < system->N; ++i )
{
+#if defined(QMMM)
+ if ( system->atoms[i].qmmm_mask == TRUE )
+ {
+#endif
type_i = system->atoms[i].type;
sbp_i = &system->reax_param.sbp[ type_i ];
@@ -315,6 +326,9 @@ void Atom_Energy( reax_system *system, control_params *control,
fprintf( out_control->eov, "%6d%15.8f%15.8f\n",
i + 1/*workspace->orig_id[i]+1*/, e_un, data->E_Ov + data->E_Un );
+#endif
+#if defined(QMMM)
+ }
#endif
}
}
diff --git a/sPuReMD/src/nonbonded.c b/sPuReMD/src/nonbonded.c
index 2c3ae783..74383c8c 100644
--- a/sPuReMD/src/nonbonded.c
+++ b/sPuReMD/src/nonbonded.c
@@ -45,11 +45,18 @@ static void Compute_Polarization_Energy( reax_system* system,
#endif
for ( i = 0; i < system->N; i++ )
{
+#if defined(QMMM)
+ if ( system->atoms[i].qmmm_mask == TRUE )
+ {
+#endif
q = system->atoms[i].q;
type_i = system->atoms[i].type;
e_pol += KCALpMOL_to_EV * (system->reax_param.sbp[ type_i ].chi * q
+ (system->reax_param.sbp[ type_i ].eta / 2.0) * SQR( q ));
+#if defined(QMMM)
+ }
+#endif
}
}
else if ( control->charge_method == ACKS2_CM )
@@ -60,6 +67,10 @@ static void Compute_Polarization_Energy( reax_system* system,
#endif
for ( i = 0; i < system->N; i++ )
{
+#if defined(QMMM)
+ if ( system->atoms[i].qmmm_mask == TRUE )
+ {
+#endif
q = system->atoms[i].q;
type_i = system->atoms[i].type;
@@ -69,6 +80,9 @@ static void Compute_Polarization_Energy( reax_system* system,
/* energy due to coupling with kinetic energy potential */
e_pol += KCALpMOL_to_EV * system->atoms[i].q * workspace->s[0][ system->N + i ];
+#if defined(QMMM)
+ }
+#endif
}
}
@@ -138,6 +152,11 @@ void vdW_Coulomb_Energy( reax_system *system, control_params *control,
{
nbr_pj = &far_nbrs->far_nbr_list[pj];
j = nbr_pj->nbr;
+
+#if defined(QMMM)
+ if ( system->atoms[i].qmmm_mask == TRUE || system->atoms[j].qmmm_mask == TRUE )
+ {
+#endif
r_ij = nbr_pj->d;
twbp = &system->reax_param.tbp[ system->atoms[i].type ]
[ system->atoms[j].type ];
@@ -163,6 +182,10 @@ void vdW_Coulomb_Energy( reax_system *system, control_params *control,
dTap = dTap * r_ij + 2.0 * workspace->Tap[2];
dTap = dTap * r_ij + workspace->Tap[1];
+#if defined(QMMM)
+ if ( system->atoms[i].qmmm_mask == TRUE && system->atoms[j].qmmm_mask == TRUE )
+ {
+#endif
/* vdWaals Calculations */
if ( system->reax_param.gp.vdw_type == 1
|| system->reax_param.gp.vdw_type == 3 )
@@ -214,6 +237,18 @@ void vdW_Coulomb_Energy( reax_system *system, control_params *control,
CEvd = self_coef * ( (de_base + de_core) * Tap
+ (e_base + e_core) * dTap );
+#if defined(QMMM)
+ }
+ else
+ {
+ e_core = 0.0;
+ de_core = 0.0;
+ e_vdW = 0.0;
+ e_base = 0.0;
+ de_base = 0.0;
+ CEvd = 0.0;
+ }
+#endif
#if defined(DEBUG_FOCUS)
fprintf( stderr, "%6d%6d%24.12f%24.12f%24.12f%24.12f\n",
@@ -221,6 +256,10 @@ void vdW_Coulomb_Energy( reax_system *system, control_params *control,
e_base, de_base, e_core, de_core ); fflush( stderr );
#endif
+#if defined(QMMM)
+ if ( system->atoms[i].qmmm_mask == TRUE && system->atoms[j].qmmm_mask == TRUE )
+ {
+#endif
/* Coulomb Calculations */
dr3gamij_1 = r_ij * r_ij * r_ij + POW( twbp->gamma, -3.0 );
dr3gamij_3 = POW( dr3gamij_1 , 1.0 / 3.0 );
@@ -231,6 +270,16 @@ void vdW_Coulomb_Energy( reax_system *system, control_params *control,
de_clb = -C_ELE * (system->atoms[i].q * system->atoms[j].q)
* (r_ij * r_ij) / POW( dr3gamij_1, 4.0 / 3.0);
CEclmb = self_coef * (de_clb * Tap + e_clb * dTap);
+#if defined(QMMM)
+ }
+ else
+ {
+ e_clb = 0.0;
+ e_ele = 0.0;
+ de_clb = 0.0;
+ CEclmb = 0.0;
+ }
+#endif
#if defined(DEBUG_FOCUS)
fprintf( stderr, "%6d%6d%24.12f%24.12f\n",
@@ -343,6 +392,9 @@ void vdW_Coulomb_Energy( reax_system *system, control_params *control,
rvec_ScaledAdd( workspace->f_vdw[j], +CEvd, nbr_pj->dvec );
rvec_ScaledAdd( workspace->f_ele[i], -CEclmb, nbr_pj->dvec );
rvec_ScaledAdd( workspace->f_ele[j], +CEclmb, nbr_pj->dvec );
+#endif
+#if defined(QMMM)
+ }
#endif
}
}
@@ -801,7 +853,7 @@ void LR_vdW_Coulomb( reax_system *system, control_params *control,
/* Coulomb calculations */
dr3gamij_1 = r_ij * r_ij * r_ij
+ POW( twbp->gamma, -3.0 );
- dr3gamij_3 = POW( dr3gamij_1 , 1.0 / 3.0 );
+ dr3gamij_3 = POW( dr3gamij_1, 1.0 / 3.0 );
tmp = Tap / dr3gamij_3;
lr->H = EV_to_KCALpMOL * tmp;
@@ -813,7 +865,7 @@ void LR_vdW_Coulomb( reax_system *system, control_params *control,
twbp->gamma, Tap, dr3gamij_3,
system->atoms[i].q, system->atoms[j].q ); */
- lr->CEclmb = C_ELE * ( dTap - Tap * r_ij / dr3gamij_1 ) / dr3gamij_3;
+ lr->CEclmb = C_ELE * (dTap - Tap * r_ij / dr3gamij_1) / dr3gamij_3;
/* fprintf( stdout, "%d %d\t%g\t%g %g\t%g %g\t%g %g\n",
i+1, j+1, r_ij, e_vdW, CEvd * r_ij,
diff --git a/sPuReMD/src/reax_types.h b/sPuReMD/src/reax_types.h
index 8ba64f9a..5a57b49a 100644
--- a/sPuReMD/src/reax_types.h
+++ b/sPuReMD/src/reax_types.h
@@ -773,6 +773,12 @@ struct reax_atom
rvec f;
/* charge on this atom, in Coulombs */
real q;
+#if defined(QMMM)
+ /* TRUE if the atom is in the QM region, FALSE otherwise (atom in MM region) */
+ int qmmm_mask;
+ /* inital charge on this atom, in Coulombs */
+ real q_init;
+#endif
};
@@ -859,10 +865,12 @@ struct reax_system
int ffield_params_allocated;
/* number of local (non-periodic image) atoms for the current simulation */
int N;
+#if defined(QMMM)
/* number of local (non-periodic image) QM atoms for the current simulation in QMMM mode */
int N_qm;
/* number of local (non-periodic image) MM atoms for the current simulation in QMMM mode */
int N_mm;
+#endif
/* max. number of local (non-periodic image) atoms across all simulations */
int N_max;
/* dimension of the N x N sparse charge method matrix H */
@@ -1588,9 +1596,13 @@ struct static_storage
/* for hydrogen bonds */
int *hbond_index;
- rvec *v_const;
- rvec *f_old;
+#if defined(QMMM)
+ /* TRUE if the atom is in the QM region, FALSE otherwise (atom in MM region) */
+ int *mask_qmmm;
+#endif
rvec *a; // used in integrators
+ rvec *f_old;
+ rvec *v_const;
/* coefficient of dDelta for force calculations */
real *CdDelta;
diff --git a/sPuReMD/src/spuremd.c b/sPuReMD/src/spuremd.c
index 364b105b..eb70ead6 100644
--- a/sPuReMD/src/spuremd.c
+++ b/sPuReMD/src/spuremd.c
@@ -105,11 +105,15 @@ static void Read_Input_Files( const char * const geo_file,
Read_Force_Field( ffield_file, system, &system->reax_param );
}
+ Set_Control_Defaults( system, control, out_control );
+
if ( control_file != NULL )
{
Read_Control_File( control_file, system, control, out_control );
}
+ Set_Control_Derived_Values( system, control );
+
if ( geo_file != NULL )
{
if ( control->geo_format == CUSTOM )
@@ -147,6 +151,7 @@ static void Read_Input_Files( const char * const geo_file,
}
+#if defined(QMMM)
/* Allocate top-level data structures and parse input files
* for the first simulation
*
@@ -161,7 +166,7 @@ static void Read_Input_Files( const char * const geo_file,
* mm_pos_y: y-coordinate of MM atom positions, in Angstroms
* mm_pos_z: z-coordinate of MM atom positions, in Angstroms
* mm_q: charge of MM atom, in Coulombs
- * sim_box: simulation box information, where the entries are
+ * sim_box_info: simulation box information, where the entries are
* - box length per dimension (3 entries)
* - angles per dimension (3 entries)
* ffield_file: file containing force field parameters
@@ -169,12 +174,10 @@ static void Read_Input_Files( const char * const geo_file,
*/
void * setup_qmmm_( int qm_num_atoms, const int * const qm_types,
const double * const qm_pos_x, const double * const qm_pos_y,
- const double * const qm_pos_z,
- int mm_num_atoms, const int * const mm_types,
+ const double * const qm_pos_z, int mm_num_atoms, const int * const mm_types,
const double * const mm_pos_x, const double * const mm_pos_y,
const double * const mm_pos_z, const double * const mm_q,
- const double * const sim_box,
- const char * const ffield_file,
+ const double * const sim_box_info, const char * const ffield_file,
const char * const control_file )
{
int i;
@@ -229,16 +232,16 @@ void * setup_qmmm_( int qm_num_atoms, const int * const qm_types,
spmd_handle->system->N_qm = qm_num_atoms;
spmd_handle->system->N_mm = mm_num_atoms;
- spmd_handle->system->N = qm_num_atoms + mm_num_atoms;
+ spmd_handle->system->N = spmd_handle->system->N_qm + spmd_handle->system->N_mm;
PreAllocate_Space( spmd_handle->system, spmd_handle->control,
spmd_handle->workspace, spmd_handle->system->N );
- Setup_Box( sim_box[0], sim_box[1], sim_box[2],
- sim_box[3], sim_box[4], sim_box[5],
+ Setup_Box( sim_box_info[0], sim_box_info[1], sim_box_info[2],
+ sim_box_info[3], sim_box_info[4], sim_box_info[5],
&spmd_handle->system->box );
- for ( i = 0; i < qm_num_atoms; ++i )
+ for ( i = 0; i < spmd_handle->system->N_qm; ++i )
{
x[0] = qm_pos_x[i];
x[1] = qm_pos_y[i];
@@ -257,31 +260,33 @@ void * setup_qmmm_( int qm_num_atoms, const int * const qm_types,
rvec_MakeZero( spmd_handle->system->atoms[i].v );
rvec_MakeZero( spmd_handle->system->atoms[i].f );
spmd_handle->system->atoms[i].q = 0.0;
+ spmd_handle->system->atoms[i].q_init = 0.0;
-// mask[i] = 1;
+ spmd_handle->system->atoms[i].qmmm_mask = TRUE;
}
- for ( i = qm_num_atoms; i < qm_num_atoms + mm_num_atoms; ++i )
+ for ( i = spmd_handle->system->N_qm; i < spmd_handle->system->N; ++i )
{
- x[0] = mm_pos_x[i - qm_num_atoms];
- x[1] = mm_pos_y[i - qm_num_atoms];
- x[2] = mm_pos_z[i - qm_num_atoms];
+ x[0] = mm_pos_x[i - spmd_handle->system->N_qm];
+ x[1] = mm_pos_y[i - spmd_handle->system->N_qm];
+ x[2] = mm_pos_z[i - spmd_handle->system->N_qm];
Fit_to_Periodic_Box( &spmd_handle->system->box, x );
spmd_handle->workspace->orig_id[i] = i + 1;
// spmd_handle->system->atoms[i].type = Get_Atom_Type( &system->reax_param,
// element, sizeof(element) );
- spmd_handle->system->atoms[i].type = mm_types[i - qm_num_atoms];
+ spmd_handle->system->atoms[i].type = mm_types[i - spmd_handle->system->N_qm];
// strncpy( spmd_handle->system->atoms[i].name, atom_name,
// sizeof(spmd_handle->system->atoms[i].name) - 1 );
// spmd_handle->system->atoms[i].name[sizeof(spmd_handle->system->atoms[i].name) - 1] = '\0';
rvec_Copy( spmd_handle->system->atoms[i].x, x );
rvec_MakeZero( spmd_handle->system->atoms[i].v );
rvec_MakeZero( spmd_handle->system->atoms[i].f );
- spmd_handle->system->atoms[i].q = mm_q[i - qm_num_atoms];
+ spmd_handle->system->atoms[i].q = mm_q[i - spmd_handle->system->N_qm];
+ spmd_handle->system->atoms[i].q_init = mm_q[i - spmd_handle->system->N_qm];
-// mask[i] = 0;
+ spmd_handle->system->atoms[i].qmmm_mask = FALSE;
}
Read_Input_Files( NULL, ffield_file, control_file,
@@ -293,6 +298,7 @@ void * setup_qmmm_( int qm_num_atoms, const int * const qm_types,
return (void *) spmd_handle;
}
+#endif
/* Allocate top-level data structures and parse input files
@@ -528,7 +534,8 @@ int simulate( const void * const handle )
spmd_handle->data->timing.end = Get_Time( );
spmd_handle->data->timing.elapsed = Get_Timing_Info( spmd_handle->data->timing.start );
- if ( spmd_handle->output_enabled == TRUE )
+ if ( spmd_handle->output_enabled == TRUE
+ && spmd_handle->out_control->log_update_freq > 0 )
{
fprintf( spmd_handle->out_control->log, "total: %.2f secs\n", spmd_handle->data->timing.elapsed );
}
@@ -581,6 +588,7 @@ int cleanup( const void * const handle )
}
+#if defined(QMMM)
/* Reset for the next simulation by parsing input files and triggering
* reallocation if more space is needed
*
@@ -596,7 +604,7 @@ int cleanup( const void * const handle )
* mm_pos_y: y-coordinate of MM atom positions, in Angstroms
* mm_pos_z: z-coordinate of MM atom positions, in Angstroms
* mm_q: charge of MM atom, in Coulombs
- * sim_box: simulation box information, where the entries are
+ * sim_box_info: simulation box information, where the entries are
* - box length per dimension (3 entries)
* - angles per dimension (3 entries)
* ffield_file: file containing force field parameters
@@ -611,7 +619,7 @@ int reset_qmmm_( const void * const handle,
int mm_num_atoms, const int * const mm_types,
const double * const mm_pos_x, const double * const mm_pos_y,
const double * const mm_pos_z, const double * const mm_q,
- const double * const sim_box,
+ const double * const sim_box_info,
const char * const ffield_file, const char * const control_file )
{
int i, ret;
@@ -636,16 +644,16 @@ int reset_qmmm_( const void * const handle,
spmd_handle->system->N_qm = qm_num_atoms;
spmd_handle->system->N_mm = mm_num_atoms;
- spmd_handle->system->N = qm_num_atoms + mm_num_atoms;
+ spmd_handle->system->N = spmd_handle->system->N_qm + spmd_handle->system->N_mm;
PreAllocate_Space( spmd_handle->system, spmd_handle->control,
spmd_handle->workspace, spmd_handle->system->N );
- Setup_Box( sim_box[0], sim_box[1], sim_box[2],
- sim_box[3], sim_box[4], sim_box[5],
+ Setup_Box( sim_box_info[0], sim_box_info[1], sim_box_info[2],
+ sim_box_info[3], sim_box_info[4], sim_box_info[5],
&spmd_handle->system->box );
- for ( i = 0; i < qm_num_atoms; ++i )
+ for ( i = 0; i < spmd_handle->system->N_qm; ++i )
{
x[0] = qm_pos_x[i];
x[1] = qm_pos_y[i];
@@ -665,30 +673,30 @@ int reset_qmmm_( const void * const handle,
rvec_MakeZero( spmd_handle->system->atoms[i].f );
spmd_handle->system->atoms[i].q = 0.0;
-// mask[i] = 1;
+ spmd_handle->system->atoms[i].qmmm_mask = TRUE;
}
- for ( i = qm_num_atoms; i < qm_num_atoms + mm_num_atoms; ++i )
+ for ( i = spmd_handle->system->N_qm; i < spmd_handle->system->N; ++i )
{
- x[0] = mm_pos_x[i - qm_num_atoms];
- x[1] = mm_pos_y[i - qm_num_atoms];
- x[2] = mm_pos_z[i - qm_num_atoms];
+ x[0] = mm_pos_x[i - spmd_handle->system->N_qm];
+ x[1] = mm_pos_y[i - spmd_handle->system->N_qm];
+ x[2] = mm_pos_z[i - spmd_handle->system->N_qm];
Fit_to_Periodic_Box( &spmd_handle->system->box, x );
spmd_handle->workspace->orig_id[i] = i + 1;
// spmd_handle->system->atoms[i].type = Get_Atom_Type( &system->reax_param,
// element, sizeof(element) );
- spmd_handle->system->atoms[i].type = mm_types[i - qm_num_atoms];
+ spmd_handle->system->atoms[i].type = mm_types[i - spmd_handle->system->N_qm];
// strncpy( spmd_handle->system->atoms[i].name, atom_name,
// sizeof(spmd_handle->system->atoms[i].name) - 1 );
// spmd_handle->system->atoms[i].name[sizeof(spmd_handle->system->atoms[i].name) - 1] = '\0';
rvec_Copy( spmd_handle->system->atoms[i].x, x );
rvec_MakeZero( spmd_handle->system->atoms[i].v );
rvec_MakeZero( spmd_handle->system->atoms[i].f );
- spmd_handle->system->atoms[i].q = mm_q[i - qm_num_atoms];
+ spmd_handle->system->atoms[i].q = mm_q[i - spmd_handle->system->N_qm];
-// mask[i] = 0;
+ spmd_handle->system->atoms[i].qmmm_mask = FALSE;
}
Read_Input_Files( NULL, ffield_file, control_file,
@@ -714,6 +722,7 @@ int reset_qmmm_( const void * const handle,
return ret;
}
+#endif
/* Reset for the next simulation by parsing input files and triggering
@@ -773,6 +782,7 @@ int reset( const void * const handle, const char * const geo_file,
}
+#if defined(QMMM)
/* Getter for atom positions in QMMM mode
*
* handle: pointer to wrapper struct with top-level data structures
@@ -946,6 +956,7 @@ int get_atom_charges_qmmm_( const void * const handle, double * const qm_q,
return ret;
}
+#endif
/* Getter for atom positions
diff --git a/sPuReMD/src/spuremd.h b/sPuReMD/src/spuremd.h
index 74803b5a..ea021f8a 100644
--- a/sPuReMD/src/spuremd.h
+++ b/sPuReMD/src/spuremd.h
@@ -33,6 +33,7 @@
extern "C" {
#endif
+#if defined(QMMM)
void * setup_qmmm_( int, const int * const,
const double * const, const double * const,
const double * const, int, const int * const,
@@ -40,6 +41,7 @@ void * setup_qmmm_( int, const int * const,
const double * const, const double * const,
const double * const,
const char * const, const char * const );
+#endif
void * setup( const char * const, const char * const,
const char * const );
@@ -50,6 +52,7 @@ int simulate( const void * const );
int cleanup( const void * const );
+#if defined(QMMM)
int reset_qmmm_( const void * const, int,
const int * const,
const double * const, const double * const,
@@ -58,10 +61,12 @@ int reset_qmmm_( const void * const, int,
const double * const, const double * const,
const double * const,
const char * const, const char * const );
+#endif
int reset( const void * const, const char * const,
const char * const, const char * const );
+#if defined(QMMM)
int get_atom_positions_qmmm_( const void * const, double * const,
double * const, double * const, double * const,
double * const, double * const );
@@ -75,6 +80,7 @@ int get_atom_forces_qmmm_( const void * const, double * const,
double * const, double * const );
int get_atom_charges_qmmm_( const void * const, double * const, double * const );
+#endif
int get_atom_positions( const void * const, double * const,
double * const, double * const );
diff --git a/sPuReMD/src/torsion_angles.c b/sPuReMD/src/torsion_angles.c
index 56541a78..b6213632 100644
--- a/sPuReMD/src/torsion_angles.c
+++ b/sPuReMD/src/torsion_angles.c
@@ -206,6 +206,10 @@ void Torsion_Angles( reax_system *system, control_params *control,
#endif
for ( j = 0; j < system->N; ++j )
{
+#if defined(QMMM)
+ if ( system->atoms[j].qmmm_mask == TRUE )
+ {
+#endif
type_j = system->atoms[j].type;
Delta_j = workspace->Delta_boc[j];
start_j = Start_Index(j, bonds);
@@ -226,6 +230,11 @@ void Torsion_Angles( reax_system *system, control_params *control,
{
pbond_jk = &bonds->bond_list[pk];
k = pbond_jk->nbr;
+
+#if defined(QMMM)
+ if ( system->atoms[k].qmmm_mask == TRUE )
+ {
+#endif
bo_jk = &pbond_jk->bo_data;
BOA_jk = bo_jk->BO - control->thb_cut;
#if defined(_OPENMP)
@@ -280,6 +289,11 @@ void Torsion_Angles( reax_system *system, control_params *control,
if ( bo_ij->BO > control->thb_cut )
{
i = p_ijk->thb;
+
+#if defined(QMMM)
+ if ( system->atoms[i].qmmm_mask == TRUE )
+ {
+#endif
type_i = system->atoms[i].type;
r_ij = pbond_ij->d;
BOA_ij = bo_ij->BO - control->thb_cut;
@@ -320,6 +334,11 @@ void Torsion_Angles( reax_system *system, control_params *control,
{
p_jkl = &thb_intrs->three_body_list[pl];
l = p_jkl->thb;
+
+#if defined(QMMM)
+ if ( system->atoms[l].qmmm_mask == TRUE )
+ {
+#endif
plk = p_jkl->pthb; //pointer to l on k's bond_list!
pbond_kl = &bonds->bond_list[plk];
bo_kl = &pbond_kl->bo_data;
@@ -749,12 +768,24 @@ void Torsion_Angles( reax_system *system, control_params *control,
rvec_ScaledAdd( workspace->f_con[l], CEconj6, dcos_omega_dl );
#endif
} // pl check ends
+#if defined(QMMM)
+ }
+#endif
} // pl loop ends
+#if defined(QMMM)
+ }
+#endif
} // pi check ends
} // pi loop ends
} // k-j neighbor check ends
} // j<k && j-k neighbor check ends
+#if defined(QMMM)
+ }
+#endif
} // pk loop ends
+#if defined(QMMM)
+ }
+#endif
} // j loop
}
diff --git a/sPuReMD/src/valence_angles.c b/sPuReMD/src/valence_angles.c
index 3f9a9c8a..5f094610 100644
--- a/sPuReMD/src/valence_angles.c
+++ b/sPuReMD/src/valence_angles.c
@@ -250,6 +250,7 @@ void Valence_Angles( reax_system *system, control_params *control,
if ( BOA_ij >= 0.0 )
{
i = pbond_ij->nbr;
+
type_i = system->atoms[i].type;
//#if defined(_OPENMP)
// f_i = &workspace->f_local[tid * system->N + i];
@@ -306,6 +307,7 @@ void Valence_Angles( reax_system *system, control_params *control,
}
k = pbond_jk->nbr;
+
type_k = system->atoms[k].type;
p_ijk = &thb_list[num_thb_intrs];
//#if defined(_OPENMP)
@@ -408,6 +410,11 @@ void Valence_Angles( reax_system *system, control_params *control,
CEval8 = CEval4 / sin_theta;
e_ang = f7_ij * f7_jk * f8_Dj * expval12theta;
+#if defined(QMMM)
+ if ( system->atoms[i].qmmm_mask == TRUE
+ && system->atoms[j].qmmm_mask == TRUE
+ && system->atoms[k].qmmm_mask == TRUE )
+#endif
e_ang_total += e_ang;
#if defined(DEBUG_FOCUS)
@@ -436,6 +443,11 @@ void Valence_Angles( reax_system *system, control_params *control,
+ p_pen4 * exp_pen4 )) / SQR( trm_pen34 );
e_pen = p_pen1 * f9_Dj * exp_pen2ij * exp_pen2jk;
+#if defined(QMMM)
+ if ( system->atoms[i].qmmm_mask == TRUE
+ && system->atoms[j].qmmm_mask == TRUE
+ && system->atoms[k].qmmm_mask == TRUE )
+#endif
e_pen_total += e_pen;
#if defined(DEBUG_FOCUS)
@@ -465,6 +477,11 @@ void Valence_Angles( reax_system *system, control_params *control,
* EXP( -p_coa3 * SQR(workspace->total_bond_order[i] - BOA_ij) )
* EXP( -p_coa3 * SQR(workspace->total_bond_order[k] - BOA_jk) )
/ (1.0 + exp_coa2);
+#if defined(QMMM)
+ if ( system->atoms[i].qmmm_mask == TRUE
+ && system->atoms[j].qmmm_mask == TRUE
+ && system->atoms[k].qmmm_mask == TRUE )
+#endif
e_coa_total += e_coa;
CEcoa1 = -2.0 * p_coa4 * (BOA_ij - 1.5) * e_coa;
@@ -472,6 +489,12 @@ void Valence_Angles( reax_system *system, control_params *control,
CEcoa3 = -p_coa2 * exp_coa2 * e_coa / (1.0 + exp_coa2);
CEcoa4 = -2.0 * p_coa3 * (workspace->total_bond_order[i] - BOA_ij) * e_coa;
CEcoa5 = -2.0 * p_coa3 * (workspace->total_bond_order[k] - BOA_jk) * e_coa;
+#if defined(QMMM)
+ if ( system->atoms[i].qmmm_mask == TRUE
+ && system->atoms[j].qmmm_mask == TRUE
+ && system->atoms[k].qmmm_mask == TRUE )
+ {
+#endif
/* calculate force contributions */
#if defined(_OPENMP)
@@ -676,6 +699,9 @@ void Valence_Angles( reax_system *system, control_params *control,
Add_dDelta( system, lists, i, CEcoa4, workspace->f_coa );
Add_dDelta( system, lists, k, CEcoa5, workspace->f_coa );
/* end coalition forces */
+#endif
+#if defined(QMMM)
+ }
#endif
}
}
diff --git a/sPuReMD/src/vector.h b/sPuReMD/src/vector.h
index 2f3dd52e..6ca360ed 100644
--- a/sPuReMD/src/vector.h
+++ b/sPuReMD/src/vector.h
@@ -78,6 +78,23 @@ static inline void Vector_MakeZero( real * const v, const unsigned int k )
}
+#if defined(QMMM)
+static inline void Vector_Mask_qmmm( real * const v, int const * const mask,
+ const unsigned int k )
+{
+ unsigned int i;
+
+#if defined(_OPENMP)
+ #pragma omp for simd schedule(simd:static)
+#endif
+ for ( i = 0; i < k; ++i )
+ {
+ v[i] = (mask[i] == TRUE ? v[i] : 0.0);
+ }
+}
+#endif
+
+
static inline void Vector_Copy( real * const dest, const real * const v, const unsigned int k )
{
unsigned int i;
diff --git a/test/library/qmmm_amber/AVE/control b/test/library/qmmm_amber/AVE/control
new file mode 100644
index 00000000..0b785422
--- /dev/null
+++ b/test/library/qmmm_amber/AVE/control
@@ -0,0 +1,75 @@
+# General parameters
+ 0 ianaly 1: perform post-run analysis
+ 1 icentr
+ 1 itrans
+ 0 imetho 0: Normal MD-run 1: Energy minimisation 2:MD-energy minimisation
+ 1 igeofo 0:xyz-input geometry 1: Biograf input geometry 2: xmol-input geometry
+ 80.000 axis1 a (for non-periodical systems)
+ 80.000 axis2 b (for non-periodical systems)
+ 80.000 axis3 c (for non-periodical systems)
+ 0.0010 cutof2 BO-cutoff for valency angles and torsion angles
+ 0.300 cutof3 BO-cutoff for bond order for graphs
+ 4 icharg Charges. 1:EEM 2:- 3: Shielded EEM (default for crystals) 4: Full system EEM 5: Fixed (unit 26) 6: Fragment EEM
+ 1 ichaen Charges. 1:include charge energy 0: Do not include charge energy
+ 0 iappen 1: Append fort.7 and fort.8
+ 0 isurpr 1: Surpress lots of output 2: Read in all geometries at the same time
+ 25 irecon Frequency of reading control-file
+ 0 icheck 0: Normal run 1:Check first derivatives;2: Single run
+ 0 idebug 0: normal run 1: debug run
+ 0 ixmolo 0: only x,y,z-coordinates in xmolout 1: x,y,z + velocities + molnr. in xmolout
+10.0000 volcha volume change (%) with 'S' and 'B' labels
+ 0 iconne 0: Normal run 1: Run with fixed connection table 2: Read in from cnt.in
+ 0 imolde 0: Normal run 1: Run with fixed molecule definition (moldef.in)
+# MD-parameters
+ 3 imdmet MD-method. 1:Velocity Verlet+Berendsen 2:Hoover-Nose;3:NVE 4: NPT
+ 0.250 tstep MD-time step (fs)
+0001.00 mdtemp 1st MD-temperature
+0001.00 mdtem2 2nd MD-temperature
+0000.00 tincr Increase/decrease temperature
+ 2 itdmet 0: T-damp atoms 1: Energy cons 2:System 3: Mols 4: Anderson 5: Mols+2 types of damping 6: System with 2 types of damping
+ 1.0 tdamp1 1st Berendsen/Anderson temperature damping constant (fs)
+ 1.0 tdamp2 2nd Berendsen/Anderson temperature damping constant (fs)
+ 192 ntdamp Nr. of atoms with 1st Berendsen damping constant and 1st MD-temperature
+0000.00 mdpres MD-pressure (MPa)
+00100.0 pdamp1 Berendsen pressure damping constant (fs)
+ 0 inpt 0: Change all cell parameters in NPT-run 1: fixed x 2: fixed y 3: fixed z
+0000000 nmdit Number of MD-iterations
+ 000000 nmdeqi Number of MD-equilibrium iterations
+ 00001 ichupd Charge update frequency
+ 010 iout1 Output to unit 71 and unit 73
+ 0100 iout2 Save coordinates
+ 0 ivels 1:Set vels and accels from moldyn.vel to zero
+ 00025 itrafr Frequency of trarot-calls
+ 1 iout3 0: create moldyn.xxxx-files 1: do not create moldyn.xxxx-files
+ 0 iravel 1: Random initial velocities
+0.00001 endmd End point criterium for MD energy minimisation
+ 000000 iout6 Save velocity file
+ 000025 irten Frequency of removal of rotational and translational energy
+ 0 npreit Nr. of iterations in previous runs
+ 02.50 range Range for back-translation of atoms
+ 1 fixedc fixed charge param
+ 1e-14 convg Tol. of the solver
+# MM-parameters
+2.00000 endmm End point criterium for MM energy minimisation
+ 00001 imaxmo Maximum movement (1/1D6 A) during minimisation 0: Conjugate gradient
+ 00000 imaxit Maximum number of iterations
+ 100 iout4 Frequency of structure output during minimisation
+ 0 iout5 1:Remove fort.57 and fort.58 files
+1.00010 celopt Cell parameter change
+ 0 icelo2 Change all cell parameters (0) or only x/y/z axis (1/2/3)
+# FF-optimisation parameters
+ 1.00 parsca Parameter optimization: parameter step scaling
+ 0.0250 parext Parameter optimization: extrapolation
+ 0 icelop 0: No cell parameter optimisation 1:Cell parameter optimisation
+ 1 igeopt 0: Always use same start gemetries 1:Use latest geometries in optimisation
+ 0 iincop heat increment optimisation 1: yes 0: no
+1.00000 accerr Accepted increase in error force field
+ 251 nmoset Nr. of molecules in training set
+#Outdated parameters
+ 0 ideve1 0: Normal run 1:Check for radical/double bond distances
+ 2000 ideve2 Frequency of radical/double bond check
+ 0 nreac 0: reactive; 1: non-reactive; 2: Place default atoms
+ 0 ibiola 0: Use old Biograf-labels 1: Assign Biograf-labels
+100.000 tdhoov Hoover-Noose temperature damping constant (fs)
+ 01.000 achoov 100*Accuracy Hoover-Noose
+ 0 itfix 1:Keep temperature fixed at exactly tset
diff --git a/test/library/qmmm_amber/AVE/fort.3 b/test/library/qmmm_amber/AVE/fort.3
new file mode 100644
index 00000000..aef79db5
--- /dev/null
+++ b/test/library/qmmm_amber/AVE/fort.3
@@ -0,0 +1,781 @@
+BIOGRF 200
+DESCRP AMBER-REAX bgf
+REMARK .bgf-file generated by amber
+HETATM 1 C -1.65460 0.98104 0.08681 C 0 0 0.11473
+HETATM 2 O -1.42542 0.31336 -1.21380 O 0 0 -0.54053
+HETATM 3 C -0.67167 -0.83192 -1.03518 C 0 0 0.09807
+HETATM 4 C 0.88286 -0.45943 -0.91203 C 0 0 -0.05189
+HETATM 5 C -0.68256 1.31997 0.93751 C 0 0 -0.03019
+HETATM 6 C 1.74254 -0.98094 0.00111 C 0 0 -0.08214
+HETATM 7 H -0.92503 -1.51330 -1.84567 H 0 0 0.07882
+HETATM 8 H -0.72447 -1.29225 -0.01108 H 0 0 0.07563
+HETATM 9 H 0.96493 0.22311 -1.72226 H 0 0 0.07363
+HETATM 10 H -1.10310 1.79396 1.84567 H 0 0 0.04146
+HETATM 11 H 0.35955 1.07205 0.83497 H 0 0 0.04636
+HETATM 12 H 1.42234 -1.79396 0.64896 H 0 0 0.05007
+HETATM 13 H 2.75145 -0.55465 -0.06272 H 0 0 0.03641
+HETATM 14 H -2.75145 1.09686 0.08169 H 0 0 0.08958
+HETATM 15 H 12.47437 -2.46007 -0.86770 H 0 0 0.41700
+HETATM 16 H 1.51043 -9.55854 -6.96047 H 0 0 0.41700
+HETATM 17 O -1.52245 -4.17882 -0.98859 O 0 0 -0.83400
+HETATM 18 H -2.05818 -3.70939 -1.62801 H 0 0 0.41700
+HETATM 19 H -2.07281 -4.90875 -0.70481 H 0 0 0.41700
+HETATM 20 O 1.64531 7.58528 7.47500 O 0 0 -0.83400
+HETATM 21 H 1.28772 7.72415 8.35197 H 0 0 0.41700
+HETATM 22 H 2.10499 8.39997 7.27200 H 0 0 0.41700
+HETATM 23 O -4.80137 1.40579 7.54104 O 0 0 -0.83400
+HETATM 24 H -4.49650 2.23505 7.90928 H 0 0 0.41700
+HETATM 25 H -5.75577 1.47887 7.54405 H 0 0 0.41700
+HETATM 26 O -5.86344 0.78091 2.53354 O 0 0 -0.83400
+HETATM 27 H -5.38331 0.88531 3.35501 H 0 0 0.41700
+HETATM 28 H -5.85173 -0.16169 2.36740 H 0 0 0.41700
+HETATM 29 O 1.18766 6.14728 -5.29738 O 0 0 -0.83400
+HETATM 30 H 2.00788 6.13799 -5.79071 H 0 0 0.41700
+HETATM 31 H 1.45909 6.10316 -4.38053 H 0 0 0.41700
+HETATM 32 O 6.23125 -4.40479 6.82307 O 0 0 -0.83400
+HETATM 33 H 6.91478 -4.69426 6.21874 H 0 0 0.41700
+HETATM 34 H 5.42190 -4.46421 6.31547 H 0 0 0.41700
+HETATM 35 O -9.78227 0.91352 3.04722 O 0 0 -0.83400
+HETATM 36 H -9.09905 1.28231 3.60707 H 0 0 0.41700
+HETATM 37 H -9.47702 1.08017 2.15543 H 0 0 0.41700
+HETATM 38 O 0.27563 2.25473 11.64407 O 0 0 -0.83400
+HETATM 39 H -0.37657 1.56635 11.77448 H 0 0 0.41700
+HETATM 40 H 0.15242 2.53255 10.73640 H 0 0 0.41700
+HETATM 41 O 7.78028 5.27491 -3.31086 O 0 0 -0.83400
+HETATM 42 H 8.59723 4.77724 -3.27699 H 0 0 0.41700
+HETATM 43 H 7.28263 4.97029 -2.55205 H 0 0 0.41700
+HETATM 44 O -9.32232 4.78234 -3.62344 O 0 0 -0.83400
+HETATM 45 H -10.08586 5.22696 -3.25525 H 0 0 0.41700
+HETATM 46 H -9.14676 5.24619 -4.44213 H 0 0 0.41700
+HETATM 47 O 4.07673 -4.96939 5.17839 O 0 0 -0.83400
+HETATM 48 H 3.68738 -5.77251 5.52428 H 0 0 0.41700
+HETATM 49 H 3.92524 -5.01703 4.23445 H 0 0 0.41700
+HETATM 50 O 0.90425 8.50466 -3.48230 O 0 0 -0.83400
+HETATM 51 H 1.43251 9.13865 -3.96730 H 0 0 0.41700
+HETATM 52 H 0.75174 7.79381 -4.10493 H 0 0 0.41700
+HETATM 53 O 4.71819 -2.92704 2.37050 O 0 0 -0.83400
+HETATM 54 H 4.88041 -2.48863 1.53521 H 0 0 0.41700
+HETATM 55 H 5.58766 -3.03999 2.75454 H 0 0 0.41700
+HETATM 56 O -2.91781 -5.62346 8.14339 O 0 0 -0.83400
+HETATM 57 H -2.98733 -6.44008 7.64890 H 0 0 0.41700
+HETATM 58 H -2.40702 -5.85577 8.91886 H 0 0 0.41700
+HETATM 59 O -7.81577 -4.63150 5.13329 O 0 0 -0.83400
+HETATM 60 H -8.11606 -3.85362 4.66323 H 0 0 0.41700
+HETATM 61 H -7.99848 -4.44079 6.05333 H 0 0 0.41700
+HETATM 62 O 2.44516 2.28769 5.27087 O 0 0 -0.83400
+HETATM 63 H 2.25147 1.60422 4.62933 H 0 0 0.41700
+HETATM 64 H 1.92723 3.03963 4.98353 H 0 0 0.41700
+HETATM 65 H 3.86232 -11.05629 2.78783 H 0 0 0.41700
+HETATM 66 O -4.19147 6.54452 7.17018 O 0 0 -0.83400
+HETATM 67 H -3.57169 7.01302 7.72930 H 0 0 0.41700
+HETATM 68 H -5.03804 6.94911 7.35955 H 0 0 0.41700
+HETATM 69 O -2.39190 -5.34344 -8.05266 O 0 0 -0.83400
+HETATM 70 H -3.24081 -5.00668 -7.76599 H 0 0 0.41700
+HETATM 71 H -2.07678 -5.86076 -7.31151 H 0 0 0.41700
+HETATM 72 O 0.16193 -8.48154 5.96619 O 0 0 -0.83400
+HETATM 73 H -0.63200 -8.27232 6.45826 H 0 0 0.41700
+HETATM 74 H -0.06298 -9.27032 5.47276 H 0 0 0.41700
+HETATM 75 O 10.61884 4.48455 3.06688 O 0 0 -0.83400
+HETATM 76 O -7.43674 0.48914 -2.65133 O 0 0 -0.83400
+HETATM 77 H -8.17450 -0.11987 -2.68366 H 0 0 0.41700
+HETATM 78 H -7.26820 0.61183 -1.71711 H 0 0 0.41700
+HETATM 79 O -3.00211 -2.03050 2.74258 O 0 0 -0.83400
+HETATM 80 H -2.45061 -1.44224 3.25835 H 0 0 0.41700
+HETATM 81 H -2.48346 -2.83031 2.65573 H 0 0 0.41700
+HETATM 82 O 0.26145 -6.47543 -8.22080 O 0 0 -0.83400
+HETATM 83 H -0.46195 -5.88762 -8.43854 H 0 0 0.41700
+HETATM 84 H 0.67277 -6.67112 -9.06267 H 0 0 0.41700
+HETATM 85 O -5.53231 -2.29276 2.41034 O 0 0 -0.83400
+HETATM 86 H -4.57734 -2.22945 2.42632 H 0 0 0.41700
+HETATM 87 H -5.76805 -2.12769 1.49743 H 0 0 0.41700
+HETATM 88 O 5.69668 3.80971 6.90579 O 0 0 -0.83400
+HETATM 89 H 4.84654 4.18824 6.68171 H 0 0 0.41700
+HETATM 90 H 6.28034 4.09997 6.20485 H 0 0 0.41700
+HETATM 91 O 5.05535 2.90612 3.12351 O 0 0 -0.83400
+HETATM 92 H 4.61555 2.40467 2.43696 H 0 0 0.41700
+HETATM 93 H 5.64669 2.27762 3.53768 H 0 0 0.41700
+HETATM 94 O 12.10441 -1.86721 2.60602 O 0 0 -0.83400
+HETATM 95 H 11.84128 -2.70607 2.98458 H 0 0 0.41700
+HETATM 96 O -0.33891 -2.59062 -5.56421 O 0 0 -0.83400
+HETATM 97 H 0.08556 -2.83587 -6.38636 H 0 0 0.41700
+HETATM 98 H -1.00742 -1.95512 -5.82007 H 0 0 0.41700
+HETATM 99 O 1.18815 -10.80660 2.75850 O 0 0 -0.83400
+HETATM 100 H 1.26126 -11.22728 1.90181 H 0 0 0.41700
+HETATM 101 H 1.31700 -9.87512 2.57970 H 0 0 0.41700
+HETATM 102 O 10.28026 -6.27764 1.25187 O 0 0 -0.83400
+HETATM 103 H 10.39772 -6.82808 0.47762 H 0 0 0.41700
+HETATM 104 O 7.46402 -0.63599 1.22262 O 0 0 -0.83400
+HETATM 105 H 6.87255 -1.12202 0.64802 H 0 0 0.41700
+HETATM 106 H 7.74581 0.11326 0.69778 H 0 0 0.41700
+HETATM 107 H -11.23556 -0.64867 -3.79961 H 0 0 0.41700
+HETATM 108 O -5.16204 -9.26949 3.87501 O 0 0 -0.83400
+HETATM 109 H -4.76341 -9.13614 3.01504 H 0 0 0.41700
+HETATM 110 O -0.77280 3.62158 -2.08117 O 0 0 -0.83400
+HETATM 111 H -1.37041 3.56446 -1.33563 H 0 0 0.41700
+HETATM 112 H 0.10056 3.57501 -1.69218 H 0 0 0.41700
+HETATM 113 O -7.33073 -7.23562 2.96805 O 0 0 -0.83400
+HETATM 114 H -7.82138 -7.15014 3.78548 H 0 0 0.41700
+HETATM 115 H -6.64493 -6.57070 3.02970 H 0 0 0.41700
+HETATM 116 O -8.95003 -2.56815 -7.08701 O 0 0 -0.83400
+HETATM 117 H -8.30912 -2.02875 -7.55017 H 0 0 0.41700
+HETATM 118 H -8.47456 -2.90675 -6.32840 H 0 0 0.41700
+HETATM 119 O 3.29681 -1.12525 7.73678 O 0 0 -0.83400
+HETATM 120 H 4.21086 -1.39029 7.83919 H 0 0 0.41700
+HETATM 121 H 3.34306 -0.22377 7.41830 H 0 0 0.41700
+HETATM 122 H 6.77244 -5.67568 7.76257 H 0 0 0.41700
+HETATM 123 O -5.80451 3.53698 2.89728 O 0 0 -0.83400
+HETATM 124 H -6.43209 3.65801 3.60983 H 0 0 0.41700
+HETATM 125 H -5.69966 2.58799 2.82907 H 0 0 0.41700
+HETATM 126 O 2.37009 -8.90098 0.65776 O 0 0 -0.83400
+HETATM 127 H 2.67859 -9.56774 0.04418 H 0 0 0.41700
+HETATM 128 H 3.11538 -8.74645 1.23820 H 0 0 0.41700
+HETATM 129 O 1.86096 3.37904 -1.12249 O 0 0 -0.83400
+HETATM 130 H 1.96782 3.64661 -0.20968 H 0 0 0.41700
+HETATM 131 H 2.56926 2.75305 -1.27307 H 0 0 0.41700
+HETATM 132 O 3.64397 5.61553 -7.06748 O 0 0 -0.83400
+HETATM 133 H 4.06592 6.25380 -7.64263 H 0 0 0.41700
+HETATM 134 H 3.27777 4.96175 -7.66305 H 0 0 0.41700
+HETATM 135 O 5.72388 -2.06404 8.09509 O 0 0 -0.83400
+HETATM 136 H 6.44578 -2.04864 8.72347 H 0 0 0.41700
+HETATM 137 H 5.91652 -2.81014 7.52723 H 0 0 0.41700
+HETATM 138 O -4.06434 -6.19809 -4.93422 O 0 0 -0.83400
+HETATM 139 H -4.62090 -5.67414 -5.51037 H 0 0 0.41700
+HETATM 140 H -3.32550 -6.45663 -5.48513 H 0 0 0.41700
+HETATM 141 H 4.98272 -0.41895 9.32879 H 0 0 0.41700
+HETATM 142 O -2.92512 10.10889 5.04619 O 0 0 -0.83400
+HETATM 143 H -3.32716 9.34681 5.46311 H 0 0 0.41700
+HETATM 144 H -3.46291 10.84836 5.32938 H 0 0 0.41700
+HETATM 145 O -3.44407 -6.00847 -0.21847 O 0 0 -0.83400
+HETATM 146 H -3.21323 -6.18917 0.69273 H 0 0 0.41700
+HETATM 147 H -4.39621 -6.10233 -0.24748 H 0 0 0.41700
+HETATM 148 H -4.02684 -2.71273 -10.76308 H 0 0 0.41700
+HETATM 149 O -4.54862 8.53014 1.73190 O 0 0 -0.83400
+HETATM 150 H -4.99310 7.76190 2.09031 H 0 0 0.41700
+HETATM 151 H -4.62900 8.43242 0.78310 H 0 0 0.41700
+HETATM 152 O 10.78063 -0.96494 -4.59914 O 0 0 -0.83400
+HETATM 153 H 11.04495 -0.47626 -3.81969 H 0 0 0.41700
+HETATM 154 H 9.93990 -0.58135 -4.84867 H 0 0 0.41700
+HETATM 155 O 2.89254 5.17030 6.65149 O 0 0 -0.83400
+HETATM 156 H 2.51351 6.04876 6.68126 H 0 0 0.41700
+HETATM 157 H 3.27138 5.10248 5.77507 H 0 0 0.41700
+HETATM 158 O 7.84173 -6.55780 -5.33986 O 0 0 -0.83400
+HETATM 159 H 8.24917 -5.99135 -5.99512 H 0 0 0.41700
+HETATM 160 H 7.27221 -7.13770 -5.84542 H 0 0 0.41700
+HETATM 161 O 4.66649 2.51840 -4.25366 O 0 0 -0.83400
+HETATM 162 H 5.06330 2.37465 -5.11279 H 0 0 0.41700
+HETATM 163 H 3.72739 2.57505 -4.43003 H 0 0 0.41700
+HETATM 164 O 0.55960 4.33625 4.39654 O 0 0 -0.83400
+HETATM 165 H 0.74411 4.98633 3.71862 H 0 0 0.41700
+HETATM 166 H 0.24107 4.84881 5.13954 H 0 0 0.41700
+HETATM 167 O -2.68137 3.79140 -0.05413 O 0 0 -0.83400
+HETATM 168 H -3.61402 3.83339 -0.26537 H 0 0 0.41700
+HETATM 169 H -2.65213 3.79414 0.90262 H 0 0 0.41700
+HETATM 170 H 6.22380 2.51887 8.41741 H 0 0 0.41700
+HETATM 171 O -4.96717 8.54163 -4.69896 O 0 0 -0.83400
+HETATM 172 H -4.98358 7.62252 -4.43215 H 0 0 0.41700
+HETATM 173 H -4.06862 8.68940 -4.99392 H 0 0 0.41700
+HETATM 174 O -2.68488 -8.74332 1.83822 O 0 0 -0.83400
+HETATM 175 H -1.84251 -8.86831 2.27528 H 0 0 0.41700
+HETATM 176 H -3.08849 -8.00934 2.30147 H 0 0 0.41700
+HETATM 177 O -5.51579 -9.91124 -1.46756 O 0 0 -0.83400
+HETATM 178 H -5.82413 -9.48484 -0.66797 H 0 0 0.41700
+HETATM 179 H -5.23631 -9.18985 -2.03121 H 0 0 0.41700
+HETATM 180 O 8.79414 -4.47557 5.75492 O 0 0 -0.83400
+HETATM 181 H 9.18631 -3.66200 6.07198 H 0 0 0.41700
+HETATM 182 H 9.54081 -5.01997 5.50523 H 0 0 0.41700
+HETATM 183 O -8.17368 -2.44728 3.54684 O 0 0 -0.83400
+HETATM 184 H -7.25929 -2.19197 3.42455 H 0 0 0.41700
+HETATM 185 H -8.64632 -1.98358 2.85559 H 0 0 0.41700
+HETATM 186 O 10.24812 -1.08944 -0.91798 O 0 0 -0.83400
+HETATM 187 H 9.39204 -0.89976 -1.30188 H 0 0 0.41700
+HETATM 188 H 10.21798 -2.02451 -0.71559 H 0 0 0.41700
+HETATM 189 O -6.12070 -8.81526 0.99087 O 0 0 -0.83400
+HETATM 190 H -6.48226 -8.47893 1.81086 H 0 0 0.41700
+HETATM 191 H -6.07446 -8.04998 0.41778 H 0 0 0.41700
+HETATM 192 O -2.30669 7.56848 8.95256 O 0 0 -0.83400
+HETATM 193 H -1.57344 8.17235 9.07054 H 0 0 0.41700
+HETATM 194 H -2.77427 7.59348 9.78740 H 0 0 0.41700
+HETATM 195 O -10.00365 -0.90118 5.19505 O 0 0 -0.83400
+HETATM 196 H -10.31338 -0.29249 4.52438 H 0 0 0.41700
+HETATM 197 H -9.29096 -1.38062 4.77262 H 0 0 0.41700
+HETATM 198 O 2.52803 5.44788 -3.06166 O 0 0 -0.83400
+HETATM 199 H 3.47255 5.37164 -3.19692 H 0 0 0.41700
+HETATM 200 H 2.31444 4.73688 -2.45744 H 0 0 0.41700
+HETATM 201 H -9.02098 6.87753 5.11797 H 0 0 0.41700
+HETATM 202 O -4.32138 0.74601 -5.84398 O 0 0 -0.83400
+HETATM 203 H -4.06427 1.01751 -4.96284 H 0 0 0.41700
+HETATM 204 H -5.18643 0.35228 -5.73044 H 0 0 0.41700
+HETATM 205 O 6.31913 -9.35841 -2.07819 O 0 0 -0.83400
+HETATM 206 H 5.90959 -9.81600 -1.34395 H 0 0 0.41700
+HETATM 207 H 5.94046 -8.47963 -2.05397 H 0 0 0.41700
+HETATM 208 O 4.67133 0.94878 1.38187 O 0 0 -0.83400
+HETATM 209 H 5.47143 0.42481 1.34266 H 0 0 0.41700
+HETATM 210 H 4.40664 1.05055 0.46764 H 0 0 0.41700
+HETATM 211 O 8.94578 -6.76656 -2.89672 O 0 0 -0.83400
+HETATM 212 H 8.30537 -6.83658 -3.60468 H 0 0 0.41700
+HETATM 213 H 8.77622 -7.52947 -2.34405 H 0 0 0.41700
+HETATM 214 O -2.14844 -8.19920 7.23757 O 0 0 -0.83400
+HETATM 215 O 0.73695 -6.22008 8.23993 O 0 0 -0.83400
+HETATM 216 H 0.31143 -6.37994 7.39755 H 0 0 0.41700
+HETATM 217 H 1.64662 -6.02430 8.01542 H 0 0 0.41700
+HETATM 218 O -5.50895 -5.88655 5.98084 O 0 0 -0.83400
+HETATM 219 H -4.92626 -5.15222 5.78728 H 0 0 0.41700
+HETATM 220 H -6.36683 -5.60081 5.66678 H 0 0 0.41700
+HETATM 221 O -5.92957 -4.78518 3.28455 O 0 0 -0.83400
+HETATM 222 H -6.46862 -4.60858 4.05557 H 0 0 0.41700
+HETATM 223 H -5.68182 -3.91871 2.96194 H 0 0 0.41700
+HETATM 224 O -2.60894 -3.58537 -3.52131 O 0 0 -0.83400
+HETATM 225 H -3.19939 -4.14359 -4.02726 H 0 0 0.41700
+HETATM 226 H -2.02298 -3.20262 -4.17429 H 0 0 0.41700
+HETATM 227 O -5.81896 4.06062 -2.16226 O 0 0 -0.83400
+HETATM 228 H -6.21046 4.00545 -1.29053 H 0 0 0.41700
+HETATM 229 H -5.80883 3.15701 -2.47789 H 0 0 0.41700
+HETATM 230 O -7.89340 2.45586 -4.45248 O 0 0 -0.83400
+HETATM 231 H -8.46115 3.09429 -4.02088 H 0 0 0.41700
+HETATM 232 H -7.76776 1.76717 -3.79967 H 0 0 0.41700
+HETATM 233 O 0.18533 0.09029 -9.70858 O 0 0 -0.83400
+HETATM 234 H -0.43518 -0.61086 -9.90756 H 0 0 0.41700
+HETATM 235 H -0.32892 0.73415 -9.22152 H 0 0 0.41700
+HETATM 236 O -0.64980 5.91174 6.36486 O 0 0 -0.83400
+HETATM 237 H 0.02469 6.45545 6.77189 H 0 0 0.41700
+HETATM 238 H -1.11162 5.51148 7.10156 H 0 0 0.41700
+HETATM 239 O 2.05252 -3.22299 -4.03573 O 0 0 -0.83400
+HETATM 240 H 2.48145 -3.57481 -4.81577 H 0 0 0.41700
+HETATM 241 H 1.17720 -2.97631 -4.33440 H 0 0 0.41700
+HETATM 242 O 7.95718 -0.49263 -5.32217 O 0 0 -0.83400
+HETATM 243 H 7.84864 -1.40822 -5.06498 H 0 0 0.41700
+HETATM 244 H 7.06933 -0.19101 -5.51451 H 0 0 0.41700
+HETATM 245 H -6.19234 8.20464 -6.04733 H 0 0 0.41700
+HETATM 246 H 5.45757 8.98645 -2.59250 H 0 0 0.41700
+HETATM 247 H 4.45199 9.15837 -1.47436 H 0 0 0.41700
+HETATM 248 O 4.11316 6.41054 0.02398 O 0 0 -0.83400
+HETATM 249 H 3.70797 5.69074 0.50766 H 0 0 0.41700
+HETATM 250 H 3.43763 6.69594 -0.59119 H 0 0 0.41700
+HETATM 251 O -6.39003 -6.38728 -0.34118 O 0 0 -0.83400
+HETATM 252 H -6.99805 -5.75842 0.04751 H 0 0 0.41700
+HETATM 253 H -6.39050 -6.17535 -1.27462 H 0 0 0.41700
+HETATM 254 O 5.85842 -2.45347 -7.81959 O 0 0 -0.83400
+HETATM 255 H 5.74919 -3.26134 -7.31796 H 0 0 0.41700
+HETATM 256 H 6.79797 -2.40654 -7.99646 H 0 0 0.41700
+HETATM 257 O -11.19635 -2.65011 1.87006 O 0 0 -0.83400
+HETATM 258 H -11.98574 -2.23232 1.52576 H 0 0 0.41700
+HETATM 259 H -10.48426 -2.28785 1.34288 H 0 0 0.41700
+HETATM 260 H 0.89347 1.56175 -11.59387 H 0 0 0.41700
+HETATM 261 O 2.37851 4.04829 1.41053 O 0 0 -0.83400
+HETATM 262 H 2.63714 3.29698 1.94426 H 0 0 0.41700
+HETATM 263 H 1.97605 4.65580 2.03117 H 0 0 0.41700
+HETATM 264 O 5.02088 8.51370 -4.39435 O 0 0 -0.83400
+HETATM 265 H 4.89613 7.60373 -4.66383 H 0 0 0.41700
+HETATM 266 O 3.90683 1.39796 7.20617 O 0 0 -0.83400
+HETATM 267 H 3.36912 1.66992 6.46244 H 0 0 0.41700
+HETATM 268 H 3.82396 2.11386 7.83612 H 0 0 0.41700
+HETATM 269 O -4.11246 -6.60061 2.71738 O 0 0 -0.83400
+HETATM 270 H -4.75111 -5.91060 2.89697 H 0 0 0.41700
+HETATM 271 H -3.38899 -6.41942 3.31738 H 0 0 0.41700
+HETATM 272 O 1.67252 -7.35881 4.00741 O 0 0 -0.83400
+HETATM 273 H 1.17999 -7.64254 4.77757 H 0 0 0.41700
+HETATM 274 H 1.09573 -7.55414 3.26890 H 0 0 0.41700
+HETATM 275 O 11.64243 1.06247 -2.93088 O 0 0 -0.83400
+HETATM 276 H 11.10474 0.79668 -2.18491 H 0 0 0.41700
+HETATM 277 O 1.60401 8.94506 1.12675 O 0 0 -0.83400
+HETATM 278 H 1.66631 9.68842 1.72655 H 0 0 0.41700
+HETATM 279 H 0.66765 8.75161 1.08153 H 0 0 0.41700
+HETATM 280 O 5.53385 -1.98063 -0.22380 O 0 0 -0.83400
+HETATM 281 H 4.91081 -1.87857 -0.94327 H 0 0 0.41700
+HETATM 282 H 6.23470 -2.52089 -0.58869 H 0 0 0.41700
+HETATM 283 O -2.73895 5.99648 4.77908 O 0 0 -0.83400
+HETATM 284 H -3.37896 6.18026 5.46672 H 0 0 0.41700
+HETATM 285 H -1.92589 5.81567 5.25076 H 0 0 0.41700
+HETATM 286 O 8.94324 2.14955 3.10427 O 0 0 -0.83400
+HETATM 287 H 9.15310 1.45100 2.48441 H 0 0 0.41700
+HETATM 288 H 9.38540 2.92238 2.75290 H 0 0 0.41700
+HETATM 289 O 9.03025 5.19097 -0.06371 O 0 0 -0.83400
+HETATM 290 H 8.94117 6.11182 0.18190 H 0 0 0.41700
+HETATM 291 H 9.96961 5.01695 -0.00401 H 0 0 0.41700
+HETATM 292 O -0.90043 9.79934 -2.05315 O 0 0 -0.83400
+HETATM 293 H -1.50838 9.99666 -2.76567 H 0 0 0.41700
+HETATM 294 H -0.21802 9.26641 -2.46123 H 0 0 0.41700
+HETATM 295 O 4.12597 8.48114 2.82768 O 0 0 -0.83400
+HETATM 296 H 3.56805 8.68350 2.07667 H 0 0 0.41700
+HETATM 297 H 4.74748 7.83283 2.49655 H 0 0 0.41700
+HETATM 298 O -6.75930 -0.01303 -5.38251 O 0 0 -0.83400
+HETATM 299 H -7.49837 0.50653 -5.69884 H 0 0 0.41700
+HETATM 300 H -6.89058 -0.06462 -4.43576 H 0 0 0.41700
+HETATM 301 O -3.44196 -2.88438 8.51624 O 0 0 -0.83400
+HETATM 302 H -3.81926 -2.38937 7.78903 H 0 0 0.41700
+HETATM 303 H -3.42725 -3.79022 8.20727 H 0 0 0.41700
+HETATM 304 O 3.16078 -5.19646 2.29728 O 0 0 -0.83400
+HETATM 305 H 3.61289 -4.36555 2.44368 H 0 0 0.41700
+HETATM 306 H 2.23324 -4.96213 2.26574 H 0 0 0.41700
+HETATM 307 O 0.99214 -2.48416 6.68060 O 0 0 -0.83400
+HETATM 308 H 1.73823 -2.05894 7.10342 H 0 0 0.41700
+HETATM 309 H 0.48099 -2.84617 7.40440 H 0 0 0.41700
+HETATM 310 H 2.11044 -11.48747 -0.46525 H 0 0 0.41700
+HETATM 311 O 1.43034 11.28060 3.67448 O 0 0 -0.83400
+HETATM 312 H 0.86793 11.19041 4.44376 H 0 0 0.41700
+HETATM 313 O 2.42532 -4.88436 -1.69139 O 0 0 -0.83400
+HETATM 314 H 3.36089 -4.87298 -1.89337 H 0 0 0.41700
+HETATM 315 H 2.05966 -4.16748 -2.20966 H 0 0 0.41700
+HETATM 316 O -4.37211 3.98196 8.67165 O 0 0 -0.83400
+HETATM 317 H -4.96545 4.33616 9.33402 H 0 0 0.41700
+HETATM 318 H -4.29210 4.68145 8.02315 H 0 0 0.41700
+HETATM 319 O -4.46063 6.24538 -7.75163 O 0 0 -0.83400
+HETATM 320 H -3.61909 5.86985 -7.49275 H 0 0 0.41700
+HETATM 321 H -4.74592 5.70744 -8.49018 H 0 0 0.41700
+HETATM 322 O -5.07232 6.43449 3.43407 O 0 0 -0.83400
+HETATM 323 H -5.35113 5.60376 3.04888 H 0 0 0.41700
+HETATM 324 H -4.20656 6.25081 3.79871 H 0 0 0.41700
+HETATM 325 O 6.44127 -0.12335 6.01936 O 0 0 -0.83400
+HETATM 326 H 6.06509 0.75678 6.02918 H 0 0 0.41700
+HETATM 327 H 5.98490 -0.58674 6.72166 H 0 0 0.41700
+HETATM 328 O 3.75562 -9.12315 3.15419 O 0 0 -0.83400
+HETATM 329 H 3.26343 -8.30270 3.18300 H 0 0 0.41700
+HETATM 330 H 3.82236 -9.39469 4.06964 H 0 0 0.41700
+HETATM 331 O -10.61228 3.94743 5.15007 O 0 0 -0.83400
+HETATM 332 H -10.98648 3.93281 4.26917 H 0 0 0.41700
+HETATM 333 O 0.57695 -4.34367 0.73611 O 0 0 -0.83400
+HETATM 334 H 0.10938 -3.98990 -0.02050 H 0 0 0.41700
+HETATM 335 H 1.47039 -4.48270 0.42199 H 0 0 0.41700
+HETATM 336 O 5.19852 -0.25246 -6.45625 O 0 0 -0.83400
+HETATM 337 H 5.32935 -0.88512 -7.16254 H 0 0 0.41700
+HETATM 338 H 4.36422 -0.50448 -6.06045 H 0 0 0.41700
+HETATM 339 O 8.30016 -1.13794 4.38521 O 0 0 -0.83400
+HETATM 340 H 7.59024 -0.66781 4.82249 H 0 0 0.41700
+HETATM 341 H 8.33851 -0.75324 3.50956 H 0 0 0.41700
+HETATM 342 O -10.56875 -3.57863 -0.85652 O 0 0 -0.83400
+HETATM 343 H -9.91300 -4.03847 -1.38070 H 0 0 0.41700
+HETATM 344 O -9.91555 1.53454 0.24262 O 0 0 -0.83400
+HETATM 345 H -10.50565 1.15137 -0.40637 H 0 0 0.41700
+HETATM 346 H -10.20444 2.44367 0.32179 H 0 0 0.41700
+HETATM 347 H 12.00913 -0.70024 -0.17835 H 0 0 0.41700
+HETATM 348 H 12.47897 -1.20158 1.17031 H 0 0 0.41700
+HETATM 349 O 6.69664 7.50824 -0.07334 O 0 0 -0.83400
+HETATM 350 H 6.14664 6.72964 -0.16003 H 0 0 0.41700
+HETATM 351 H 6.31013 8.14290 -0.67671 H 0 0 0.41700
+HETATM 352 O -2.26700 8.44526 -5.81695 O 0 0 -0.83400
+HETATM 353 H -2.13618 9.17509 -5.21159 H 0 0 0.41700
+HETATM 354 H -1.70912 7.74695 -5.47434 H 0 0 0.41700
+HETATM 355 O 4.05205 -1.57307 -2.65767 O 0 0 -0.83400
+HETATM 356 H 4.10140 -2.52036 -2.78593 H 0 0 0.41700
+HETATM 357 H 3.62708 -1.24584 -3.45048 H 0 0 0.41700
+HETATM 358 H -11.49627 3.87188 -3.88160 H 0 0 0.41700
+HETATM 359 H 8.09367 6.23373 -5.12648 H 0 0 0.41700
+HETATM 360 O -4.94353 6.03694 -3.63876 O 0 0 -0.83400
+HETATM 361 H -5.21707 5.30994 -3.07941 H 0 0 0.41700
+HETATM 362 H -3.99849 5.92099 -3.73719 H 0 0 0.41700
+HETATM 363 O 1.85661 -7.42402 -1.74087 O 0 0 -0.83400
+HETATM 364 H 1.17884 -7.67424 -1.11297 H 0 0 0.41700
+HETATM 365 H 2.01666 -6.49705 -1.56382 H 0 0 0.41700
+HETATM 366 O 6.94195 6.83397 3.93913 O 0 0 -0.83400
+HETATM 367 H 6.36795 7.26666 4.57122 H 0 0 0.41700
+HETATM 368 H 7.75542 6.68099 4.41986 H 0 0 0.41700
+HETATM 369 O -1.36148 -11.05704 3.15380 O 0 0 -0.83400
+HETATM 370 H -0.41651 -10.96411 3.03289 H 0 0 0.41700
+HETATM 371 O 7.45869 3.88688 4.86320 O 0 0 -0.83400
+HETATM 372 H 7.86924 3.41197 4.14060 H 0 0 0.41700
+HETATM 373 H 7.04237 4.64201 4.44765 H 0 0 0.41700
+HETATM 374 O -1.50201 -7.70158 -1.95787 O 0 0 -0.83400
+HETATM 375 H -1.03907 -7.18269 -2.61566 H 0 0 0.41700
+HETATM 376 H -2.08443 -7.07655 -1.52619 H 0 0 0.41700
+HETATM 377 O 0.77453 -3.91483 9.28125 O 0 0 -0.83400
+HETATM 378 H 0.33533 -3.78157 10.12123 H 0 0 0.41700
+HETATM 379 H 0.53433 -4.80496 9.02393 H 0 0 0.41700
+HETATM 380 O 1.43188 6.23063 2.66750 O 0 0 -0.83400
+HETATM 381 H 0.69311 6.79579 2.89341 H 0 0 0.41700
+HETATM 382 H 2.16391 6.83362 2.53802 H 0 0 0.41700
+HETATM 383 O -8.70994 -1.20990 1.08300 O 0 0 -0.83400
+HETATM 384 H -8.21803 -1.86587 0.58909 H 0 0 0.41700
+HETATM 385 H -8.68323 -0.42661 0.53347 H 0 0 0.41700
+HETATM 386 O -7.82822 8.89435 2.05490 O 0 0 -0.83400
+HETATM 387 H -7.53469 8.59506 2.91542 H 0 0 0.41700
+HETATM 388 O -8.58987 6.22957 0.70702 O 0 0 -0.83400
+HETATM 389 H -8.20377 6.68768 -0.03950 H 0 0 0.41700
+HETATM 390 H -8.70956 6.90985 1.36969 H 0 0 0.41700
+HETATM 391 O 2.90643 3.26410 9.16992 O 0 0 -0.83400
+HETATM 392 H 1.96571 3.14889 9.03573 H 0 0 0.41700
+HETATM 393 H 3.13064 4.03330 8.64618 H 0 0 0.41700
+HETATM 394 O 4.81196 -7.07100 -1.94947 O 0 0 -0.83400
+HETATM 395 H 3.95665 -7.48614 -2.06053 H 0 0 0.41700
+HETATM 396 H 4.98223 -7.12117 -1.00888 H 0 0 0.41700
+HETATM 397 O 5.08592 5.72725 -4.07195 O 0 0 -0.83400
+HETATM 398 H 6.01514 5.71415 -3.84260 H 0 0 0.41700
+HETATM 399 H 5.03677 5.24916 -4.89975 H 0 0 0.41700
+HETATM 400 O 9.63937 1.36204 -6.87457 O 0 0 -0.83400
+HETATM 401 H 9.01771 0.70564 -6.56004 H 0 0 0.41700
+HETATM 402 H 9.17807 2.19503 -6.77682 H 0 0 0.41700
+HETATM 403 O 4.76634 -6.64142 0.64817 O 0 0 -0.83400
+HETATM 404 H 4.22087 -6.11798 1.23528 H 0 0 0.41700
+HETATM 405 H 5.65159 -6.30059 0.77628 H 0 0 0.41700
+HETATM 406 O -0.89941 -0.42265 4.04954 O 0 0 -0.83400
+HETATM 407 H -0.92545 0.51418 4.24424 H 0 0 0.41700
+HETATM 408 H -0.03459 -0.56322 3.66410 H 0 0 0.41700
+HETATM 409 H -11.61596 0.45658 3.20305 H 0 0 0.41700
+HETATM 410 O 2.44528 -1.73397 10.39924 O 0 0 -0.83400
+HETATM 411 H 2.61209 -1.67533 9.45851 H 0 0 0.41700
+HETATM 412 H 1.56176 -1.38219 10.50821 H 0 0 0.41700
+HETATM 413 O 2.47066 -2.16897 4.41589 O 0 0 -0.83400
+HETATM 414 H 1.93589 -2.27179 5.20309 H 0 0 0.41700
+HETATM 415 H 3.25412 -2.69106 4.58870 H 0 0 0.41700
+HETATM 416 O -1.08487 2.15551 4.72745 O 0 0 -0.83400
+HETATM 417 H -0.93568 1.69959 5.55576 H 0 0 0.41700
+HETATM 418 H -0.29623 2.68248 4.59868 H 0 0 0.41700
+HETATM 419 O -1.14714 -4.34202 2.68458 O 0 0 -0.83400
+HETATM 420 H -0.95067 -5.17838 3.10664 H 0 0 0.41700
+HETATM 421 H -0.45469 -4.23429 2.03256 H 0 0 0.41700
+HETATM 422 O 10.26020 1.30876 0.53961 O 0 0 -0.83400
+HETATM 423 H 10.43411 0.39822 0.30105 H 0 0 0.41700
+HETATM 424 H 9.53604 1.57113 -0.02872 H 0 0 0.41700
+HETATM 425 O 6.71490 4.12624 -0.77932 O 0 0 -0.83400
+HETATM 426 H 7.43941 4.69718 -0.52369 H 0 0 0.41700
+HETATM 427 H 6.04583 4.27451 -0.11104 H 0 0 0.41700
+HETATM 428 O -6.56002 1.13141 -0.03053 O 0 0 -0.83400
+HETATM 429 H -6.38997 0.86494 0.87298 H 0 0 0.41700
+HETATM 430 H -6.81566 2.05149 0.03520 H 0 0 0.41700
+HETATM 431 O -0.88532 7.82051 1.73006 O 0 0 -0.83400
+HETATM 432 H -1.69416 7.30914 1.70712 H 0 0 0.41700
+HETATM 433 H -1.17282 8.71545 1.91077 H 0 0 0.41700
+HETATM 434 O -1.23848 -10.31599 -2.64164 O 0 0 -0.83400
+HETATM 435 H -1.86511 -10.76301 -2.07266 H 0 0 0.41700
+HETATM 436 H -1.18580 -9.42778 -2.28873 H 0 0 0.41700
+HETATM 437 H -6.76024 -7.52965 -7.16640 H 0 0 0.41700
+HETATM 438 O -6.36403 -1.82293 -0.31996 O 0 0 -0.83400
+HETATM 439 H -6.20856 -0.88029 -0.37903 H 0 0 0.41700
+HETATM 440 H -5.61846 -2.21746 -0.77240 H 0 0 0.41700
+HETATM 441 O -5.00210 -4.59571 -6.68907 O 0 0 -0.83400
+HETATM 442 H -5.33979 -3.80354 -6.27116 H 0 0 0.41700
+HETATM 443 H -5.30541 -4.53648 -7.59501 H 0 0 0.41700
+HETATM 444 O -1.99926 4.11878 2.72816 O 0 0 -0.83400
+HETATM 445 H -1.75381 3.29018 3.13976 H 0 0 0.41700
+HETATM 446 H -2.30607 4.66433 3.45237 H 0 0 0.41700
+HETATM 447 O -3.99586 1.10287 4.58511 O 0 0 -0.83400
+HETATM 448 H -3.15962 1.56832 4.60171 H 0 0 0.41700
+HETATM 449 H -4.33291 1.18407 5.47732 H 0 0 0.41700
+HETATM 450 H -5.27196 -7.06546 -8.89856 H 0 0 0.41700
+HETATM 451 O 2.71222 8.28329 -1.46545 O 0 0 -0.83400
+HETATM 452 H 1.94144 8.40867 -2.01900 H 0 0 0.41700
+HETATM 453 H 2.42582 8.55035 -0.59202 H 0 0 0.41700
+HETATM 454 O -9.57763 -1.26796 -2.64866 O 0 0 -0.83400
+HETATM 455 H -9.41151 -2.11143 -3.06961 H 0 0 0.41700
+HETATM 456 H -9.98888 -1.49357 -1.81427 H 0 0 0.41700
+HETATM 457 O -2.72189 0.78078 9.66725 O 0 0 -0.83400
+HETATM 458 H -3.29611 1.21251 10.29979 H 0 0 0.41700
+HETATM 459 H -3.29253 0.57605 8.92652 H 0 0 0.41700
+HETATM 460 O -1.79943 4.19360 -10.02272 O 0 0 -0.83400
+HETATM 461 H -2.50283 3.65606 -10.38673 H 0 0 0.41700
+HETATM 462 H -2.22649 4.72269 -9.34898 H 0 0 0.41700
+HETATM 463 O -3.07333 -10.58769 -0.17281 O 0 0 -0.83400
+HETATM 464 H -2.87459 -10.09203 0.62158 H 0 0 0.41700
+HETATM 465 H -3.87822 -10.19186 -0.50703 H 0 0 0.41700
+HETATM 466 O -4.14459 -8.30711 -3.04888 O 0 0 -0.83400
+HETATM 467 H -4.20541 -7.53280 -3.60832 H 0 0 0.41700
+HETATM 468 H -3.22957 -8.32891 -2.76872 H 0 0 0.41700
+HETATM 469 O -6.40244 -2.97748 -5.15540 O 0 0 -0.83400
+HETATM 470 H -7.11989 -3.05153 -4.52610 H 0 0 0.41700
+HETATM 471 H -6.35367 -2.04239 -5.35408 H 0 0 0.41700
+HETATM 472 O -6.91490 3.85033 0.32613 O 0 0 -0.83400
+HETATM 473 H -6.42495 3.91710 1.14571 H 0 0 0.41700
+HETATM 474 H -7.70511 4.36916 0.47652 H 0 0 0.41700
+HETATM 475 O 6.62840 2.83802 -6.41474 O 0 0 -0.83400
+HETATM 476 H 6.90006 3.55205 -6.99144 H 0 0 0.41700
+HETATM 477 H 6.48407 2.09678 -7.00293 H 0 0 0.41700
+HETATM 478 O 0.16963 2.53329 8.83607 O 0 0 -0.83400
+HETATM 479 H 0.06891 1.84964 8.17371 H 0 0 0.41700
+HETATM 480 H -0.40838 3.23773 8.54301 H 0 0 0.41700
+HETATM 481 O -7.73082 3.27607 -7.04417 O 0 0 -0.83400
+HETATM 482 H -7.92855 3.17011 -6.11363 H 0 0 0.41700
+HETATM 483 H -8.06033 4.14815 -7.26126 H 0 0 0.41700
+HETATM 484 O 3.01412 -0.25348 -9.16032 O 0 0 -0.83400
+HETATM 485 H 3.35067 -1.12495 -9.36892 H 0 0 0.41700
+HETATM 486 H 2.07416 -0.31249 -9.33129 H 0 0 0.41700
+HETATM 487 O -1.60144 4.77965 8.60771 O 0 0 -0.83400
+HETATM 488 H -2.49082 4.44343 8.71809 H 0 0 0.41700
+HETATM 489 H -1.65233 5.69088 8.89634 H 0 0 0.41700
+HETATM 490 O 2.63212 -4.06327 -10.34548 O 0 0 -0.83400
+HETATM 491 H 1.91677 -4.00557 -10.97887 H 0 0 0.41700
+HETATM 492 H 2.53411 -4.93025 -9.95181 H 0 0 0.41700
+HETATM 493 H -8.16260 8.74724 -2.29250 H 0 0 0.41700
+HETATM 494 O 1.64386 1.88207 -4.00667 O 0 0 -0.83400
+HETATM 495 H 1.24591 2.25208 -4.79468 H 0 0 0.41700
+HETATM 496 H 1.22021 2.34549 -3.28418 H 0 0 0.41700
+HETATM 497 O 10.04210 -1.54664 6.39276 O 0 0 -0.83400
+HETATM 498 H 9.66714 -1.21099 5.57853 H 0 0 0.41700
+HETATM 499 O 7.37859 -5.62019 1.25693 O 0 0 -0.83400
+HETATM 500 H 7.62628 -4.89311 1.82811 H 0 0 0.41700
+HETATM 501 H 8.09455 -6.24926 1.34585 H 0 0 0.41700
+HETATM 502 O 5.08168 -5.41299 -7.57551 O 0 0 -0.83400
+HETATM 503 H 5.86262 -5.42927 -8.12877 H 0 0 0.41700
+HETATM 504 H 5.29108 -5.99913 -6.84831 H 0 0 0.41700
+HETATM 505 O 4.25060 5.34546 4.22885 O 0 0 -0.83400
+HETATM 506 H 5.05976 5.85379 4.28453 H 0 0 0.41700
+HETATM 507 H 4.51114 4.52572 3.80886 H 0 0 0.41700
+HETATM 508 O 7.65593 -3.68376 3.22965 O 0 0 -0.83400
+HETATM 509 H 7.88433 -4.18220 4.01426 H 0 0 0.41700
+HETATM 510 H 7.94435 -2.79122 3.42050 H 0 0 0.41700
+HETATM 511 O 7.69644 -2.97175 -4.00379 O 0 0 -0.83400
+HETATM 512 H 8.47729 -3.47487 -4.23485 H 0 0 0.41700
+HETATM 513 H 6.98653 -3.61375 -4.01271 H 0 0 0.41700
+HETATM 514 O 0.13722 -7.91205 1.76024 O 0 0 -0.83400
+HETATM 515 H -0.35299 -7.20036 1.34863 H 0 0 0.41700
+HETATM 516 H 0.84515 -8.10252 1.14479 H 0 0 0.41700
+HETATM 517 O -6.74347 -5.41264 -2.91010 O 0 0 -0.83400
+HETATM 518 H -5.94387 -5.13375 -3.35631 H 0 0 0.41700
+HETATM 519 H -7.24161 -5.87969 -3.58088 H 0 0 0.41700
+HETATM 520 O -7.44652 2.34018 7.67247 O 0 0 -0.83400
+HETATM 521 H -8.38900 2.23394 7.80160 H 0 0 0.41700
+HETATM 522 H -7.22924 3.14415 8.14432 H 0 0 0.41700
+HETATM 523 O 0.49814 -7.41774 -4.00280 O 0 0 -0.83400
+HETATM 524 H 1.10585 -7.11930 -4.67944 H 0 0 0.41700
+HETATM 525 H 1.03805 -7.50405 -3.21713 H 0 0 0.41700
+HETATM 526 H 8.36694 2.49895 6.98528 H 0 0 0.41700
+HETATM 527 O -6.84253 7.80913 -1.65907 O 0 0 -0.83400
+HETATM 528 H -5.89549 7.93662 -1.71477 H 0 0 0.41700
+HETATM 529 H -6.99620 6.96610 -2.08559 H 0 0 0.41700
+HETATM 530 O 0.98204 3.40426 -6.45759 O 0 0 -0.83400
+HETATM 531 H 0.90327 4.20053 -5.93225 H 0 0 0.41700
+HETATM 532 H 1.50468 3.66775 -7.21499 H 0 0 0.41700
+HETATM 533 H -2.32230 -10.04266 4.42941 H 0 0 0.41700
+HETATM 534 O -8.42746 0.57146 -7.93668 O 0 0 -0.83400
+HETATM 535 H -8.30027 1.51720 -7.86157 H 0 0 0.41700
+HETATM 536 O 4.29675 1.37781 -1.70238 O 0 0 -0.83400
+HETATM 537 H 4.57926 1.79786 -2.51477 H 0 0 0.41700
+HETATM 538 H 4.09831 0.47622 -1.95535 H 0 0 0.41700
+HETATM 539 O -1.94330 -0.43676 -6.70157 O 0 0 -0.83400
+HETATM 540 H -1.70408 0.39212 -7.11626 H 0 0 0.41700
+HETATM 541 H -2.76782 -0.25065 -6.25238 H 0 0 0.41700
+HETATM 542 O -5.46158 -1.47354 6.99725 O 0 0 -0.83400
+HETATM 543 H -5.34751 -0.59052 7.34868 H 0 0 0.41700
+HETATM 544 H -5.92977 -1.94650 7.68524 H 0 0 0.41700
+HETATM 545 H 4.07930 10.04324 3.77952 H 0 0 0.41700
+HETATM 546 O 0.06660 -0.20928 9.87424 O 0 0 -0.83400
+HETATM 547 H -0.15014 0.63944 9.48832 H 0 0 0.41700
+HETATM 548 H -0.73657 -0.72305 9.78948 H 0 0 0.41700
+HETATM 549 O -2.51352 0.37280 -10.77570 O 0 0 -0.83400
+HETATM 550 H -3.42753 0.34026 -10.49328 H 0 0 0.41700
+HETATM 551 O -2.48532 4.26680 -7.09615 O 0 0 -0.83400
+HETATM 552 H -3.33279 4.03514 -6.71620 H 0 0 0.41700
+HETATM 553 H -1.97758 4.60715 -6.35953 H 0 0 0.41700
+HETATM 554 H 5.45245 -2.90517 -9.75409 H 0 0 0.41700
+HETATM 555 H 4.09032 -3.06208 -10.39515 H 0 0 0.41700
+HETATM 556 O 7.21371 -4.00255 -0.95789 O 0 0 -0.83400
+HETATM 557 H 8.15888 -3.91484 -1.08115 H 0 0 0.41700
+HETATM 558 H 7.12265 -4.46163 -0.12291 H 0 0 0.41700
+HETATM 559 O 2.40899 3.83778 -8.89353 O 0 0 -0.83400
+HETATM 560 H 2.98360 3.20411 -9.32309 H 0 0 0.41700
+HETATM 561 H 1.82625 4.14263 -9.58903 H 0 0 0.41700
+HETATM 562 O 3.53253 -5.93545 7.86977 O 0 0 -0.83400
+HETATM 563 H 3.55944 -6.78674 8.30660 H 0 0 0.41700
+HETATM 564 H 4.29733 -5.46895 8.20694 H 0 0 0.41700
+HETATM 565 O 0.97513 -3.08556 -7.91344 O 0 0 -0.83400
+HETATM 566 H 1.77654 -3.43182 -7.52090 H 0 0 0.41700
+HETATM 567 H 1.15217 -3.08008 -8.85410 H 0 0 0.41700
+HETATM 568 O -3.41275 -2.25176 -8.59787 O 0 0 -0.83400
+HETATM 569 H -2.87308 -1.49037 -8.81065 H 0 0 0.41700
+HETATM 570 H -2.78268 -2.95239 -8.42948 H 0 0 0.41700
+HETATM 571 O 7.74809 -0.84979 -2.06751 O 0 0 -0.83400
+HETATM 572 H 7.28414 -0.35701 -2.74438 H 0 0 0.41700
+HETATM 573 H 7.68484 -1.76225 -2.34975 H 0 0 0.41700
+HETATM 574 O 6.57784 -8.62468 3.40201 O 0 0 -0.83400
+HETATM 575 H 5.69038 -8.94250 3.23576 H 0 0 0.41700
+HETATM 576 H 6.72520 -7.96187 2.72733 H 0 0 0.41700
+HETATM 577 H 11.63402 2.01999 1.20276 H 0 0 0.41700
+HETATM 578 H 11.55270 2.88015 2.44554 H 0 0 0.41700
+HETATM 579 O -2.69224 -1.75675 10.80550 O 0 0 -0.83400
+HETATM 580 H -2.99347 -2.15910 9.99088 H 0 0 0.41700
+HETATM 581 H -3.16243 -0.92416 10.84957 H 0 0 0.41700
+HETATM 582 O 2.91726 -4.37748 -6.57098 O 0 0 -0.83400
+HETATM 583 H 3.76697 -4.70784 -6.86270 H 0 0 0.41700
+HETATM 584 H 2.32801 -5.12743 -6.65214 H 0 0 0.41700
+HETATM 585 O -1.08831 1.89265 -8.09227 O 0 0 -0.83400
+HETATM 586 H -0.44589 2.36778 -7.56522 H 0 0 0.41700
+HETATM 587 H -1.61287 2.57821 -8.50589 H 0 0 0.41700
+HETATM 588 O 2.86029 -10.56198 -1.85038 O 0 0 -0.83400
+HETATM 589 H 2.54234 -10.06227 -2.60234 H 0 0 0.41700
+HETATM 590 H 3.61700 -11.04314 -2.18521 H 0 0 0.41700
+HETATM 591 O -0.02782 10.52005 6.01351 O 0 0 -0.83400
+HETATM 592 H -0.09484 9.57535 5.87465 H 0 0 0.41700
+HETATM 593 H -0.87242 10.86520 5.72410 H 0 0 0.41700
+HETATM 594 H 2.03877 -9.80840 5.58481 H 0 0 0.41700
+HETATM 595 H 1.76343 -10.88056 4.55251 H 0 0 0.41700
+HETATM 596 O -3.09165 1.36434 -3.40809 O 0 0 -0.83400
+HETATM 597 H -2.53659 1.12919 -2.66455 H 0 0 0.41700
+HETATM 598 H -2.79345 2.23749 -3.66289 H 0 0 0.41700
+HETATM 599 H 7.20645 -8.09502 5.02739 H 0 0 0.41700
+HETATM 600 O 4.07141 2.11995 -9.97116 O 0 0 -0.83400
+HETATM 601 H 4.82830 1.91723 -10.52094 H 0 0 0.41700
+HETATM 602 H 3.73104 1.26505 -9.70750 H 0 0 0.41700
+HETATM 603 O 5.06216 7.36447 6.59556 O 0 0 -0.83400
+HETATM 604 H 4.31120 7.62770 6.06359 H 0 0 0.41700
+HETATM 605 O -0.25130 0.52989 6.85064 O 0 0 -0.83400
+HETATM 606 H -1.03228 0.06657 6.54790 H 0 0 0.41700
+HETATM 607 H 0.27877 -0.14542 7.27398 H 0 0 0.41700
+HETATM 608 O -2.52058 10.11709 2.31315 O 0 0 -0.83400
+HETATM 609 H -3.26835 9.53915 2.16129 H 0 0 0.41700
+HETATM 610 H -2.54035 10.29213 3.25401 H 0 0 0.41700
+HETATM 611 O -11.14957 4.65746 0.40068 O 0 0 -0.83400
+HETATM 612 H -11.48478 4.56472 1.29246 H 0 0 0.41700
+HETATM 613 H -10.40364 5.25125 0.48574 H 0 0 0.41700
+HETATM 614 H -6.35737 -9.12296 -4.90288 H 0 0 0.41700
+HETATM 615 O -1.62373 5.47765 -3.92219 O 0 0 -0.83400
+HETATM 616 H -1.24599 4.92136 -3.24095 H 0 0 0.41700
+HETATM 617 H -0.87654 5.95421 -4.28389 H 0 0 0.41700
+HETATM 618 O -4.09146 -2.01346 -1.92537 O 0 0 -0.83400
+HETATM 619 H -3.76319 -1.16883 -2.23370 H 0 0 0.41700
+HETATM 620 H -3.51760 -2.65869 -2.33840 H 0 0 0.41700
+HETATM 621 O 7.91131 1.61057 -0.50124 O 0 0 -0.83400
+HETATM 622 H 7.33648 2.36557 -0.62683 H 0 0 0.41700
+HETATM 623 H 7.67867 1.01175 -1.21084 H 0 0 0.41700
+HETATM 624 O -4.06295 8.29877 -0.82034 O 0 0 -0.83400
+HETATM 625 H -4.23127 9.20561 -1.07636 H 0 0 0.41700
+HETATM 626 H -3.18175 8.11728 -1.14709 H 0 0 0.41700
+HETATM 627 O 7.05644 1.94118 -3.21256 O 0 0 -0.83400
+HETATM 628 H 7.72038 2.33119 -3.78116 H 0 0 0.41700
+HETATM 629 H 6.23065 2.32017 -3.51368 H 0 0 0.41700
+HETATM 630 O 9.85020 -3.99178 -0.39739 O 0 0 -0.83400
+HETATM 631 H 9.86541 -4.75973 0.17379 H 0 0 0.41700
+HETATM 632 H 10.71884 -3.97487 -0.79916 H 0 0 0.41700
+HETATM 633 O -2.25637 -6.27641 4.68283 O 0 0 -0.83400
+HETATM 634 H -1.67437 -5.62071 5.06697 H 0 0 0.41700
+HETATM 635 H -2.43028 -6.88990 5.39670 H 0 0 0.41700
+HETATM 636 H 3.99235 -7.34834 -8.02710 H 0 0 0.41700
+HETATM 637 O -1.78086 -7.10843 -6.06497 O 0 0 -0.83400
+HETATM 638 H -1.66434 -7.91445 -6.56795 H 0 0 0.41700
+HETATM 639 H -1.06427 -7.11639 -5.43041 H 0 0 0.41700
+HETATM 640 O -7.31236 7.59482 4.46127 O 0 0 -0.83400
+HETATM 641 H -7.28909 7.61187 5.41804 H 0 0 0.41700
+HETATM 642 H -6.44318 7.28565 4.20602 H 0 0 0.41700
+HETATM 643 H -5.43092 10.31183 1.82974 H 0 0 0.41700
+HETATM 644 O 5.01348 -4.45553 -2.38305 O 0 0 -0.83400
+HETATM 645 H 5.13301 -5.39295 -2.53534 H 0 0 0.41700
+HETATM 646 H 5.78979 -4.19227 -1.88881 H 0 0 0.41700
+HETATM 647 H 2.93270 -3.53082 10.32034 H 0 0 0.41700
+HETATM 648 H 3.27931 -4.84398 9.65215 H 0 0 0.41700
+HETATM 649 O -4.05359 -3.60361 5.52138 O 0 0 -0.83400
+HETATM 650 H -3.46705 -3.38966 4.79582 H 0 0 0.41700
+HETATM 651 H -4.53903 -2.79471 5.68345 H 0 0 0.41700
+HETATM 652 O -0.82696 -4.75665 6.06102 O 0 0 -0.83400
+HETATM 653 H 0.05580 -4.40574 5.94346 H 0 0 0.41700
+HETATM 654 H -1.23132 -4.17664 6.70624 H 0 0 0.41700
+HETATM 655 O 2.05725 0.26234 3.50811 O 0 0 -0.83400
+HETATM 656 H 2.81916 0.32976 2.93263 H 0 0 0.41700
+HETATM 657 H 2.20698 -0.53827 4.01094 H 0 0 0.41700
+HETATM 658 O -4.98307 3.26701 -6.45794 O 0 0 -0.83400
+HETATM 659 H -4.73367 2.36362 -6.26320 H 0 0 0.41700
+HETATM 660 H -5.93378 3.23534 -6.56458 H 0 0 0.41700
+HETATM 661 O 10.62945 1.59100 5.08413 O 0 0 -0.83400
+HETATM 662 H 10.27696 2.01969 5.86400 H 0 0 0.41700
+HETATM 663 H 10.09126 1.92186 4.36502 H 0 0 0.41700
+HETATM 664 O 1.75613 -6.90244 -6.13554 O 0 0 -0.83400
+HETATM 665 H 1.07061 -6.74561 -6.78493 H 0 0 0.41700
+HETATM 666 H 2.23627 -7.65903 -6.47207 H 0 0 0.41700
+HETATM 667 O -4.31764 3.20194 -8.96267 O 0 0 -0.83400
+HETATM 668 H -4.58167 2.28488 -9.03699 H 0 0 0.41700
+HETATM 669 H -4.62044 3.46878 -8.09472 H 0 0 0.41700
+HETATM 670 O -5.90723 -4.00919 -9.27177 O 0 0 -0.83400
+HETATM 671 H -5.57042 -4.28363 -10.12469 H 0 0 0.41700
+HETATM 672 H -6.25401 -3.12990 -9.42283 H 0 0 0.41700
+HETATM 673 O -9.07584 -3.58333 -3.91105 O 0 0 -0.83400
+HETATM 674 H -8.47711 -4.25368 -3.58182 H 0 0 0.41700
+HETATM 675 H -9.87724 -4.06159 -4.12378 H 0 0 0.41700
+HETATM 676 O 10.15990 3.44728 -3.62674 O 0 0 -0.83400
+HETATM 677 H 10.67066 2.67056 -3.85491 H 0 0 0.41700
+HETATM 678 O -0.40127 10.97249 0.69380 O 0 0 -0.83400
+HETATM 679 H -1.16182 10.85306 1.26260 H 0 0 0.41700
+HETATM 680 H -0.70052 10.68681 -0.16937 H 0 0 0.41700
+HETATM 681 O -5.43324 0.72930 -9.49308 O 0 0 -0.83400
+HETATM 682 H -5.35959 -0.14576 -9.11219 H 0 0 0.41700
+HETATM 683 H -6.32926 0.77040 -9.82729 H 0 0 0.41700
+HETATM 684 O -8.38182 -6.54287 -4.92307 O 0 0 -0.83400
+HETATM 685 H -8.24341 -7.47240 -5.10483 H 0 0 0.41700
+HETATM 686 O 2.66164 -0.73941 -5.11390 O 0 0 -0.83400
+HETATM 687 H 1.92784 -1.35258 -5.15630 H 0 0 0.41700
+HETATM 688 H 2.28878 0.05416 -4.72991 H 0 0 0.41700
+HETATM 689 O -8.17434 2.36772 5.01513 O 0 0 -0.83400
+HETATM 690 H -8.96595 2.88879 5.14954 H 0 0 0.41700
+HETATM 691 H -7.82156 2.23665 5.89525 H 0 0 0.41700
+HETATM 692 O -8.30784 -4.61066 1.18657 O 0 0 -0.83400
+HETATM 693 H -7.73024 -4.96191 1.86423 H 0 0 0.41700
+HETATM 694 H -8.89660 -4.01993 1.65629 H 0 0 0.41700
+HETATM 695 O -2.15914 -10.49054 -5.19114 O 0 0 -0.83400
+HETATM 696 H -2.04942 -10.61688 -4.24868 H 0 0 0.41700
+HETATM 697 H -3.05762 -10.17482 -5.28748 H 0 0 0.41700
+HETATM 698 O 10.23097 -4.71126 -4.30465 O 0 0 -0.83400
+HETATM 699 H 9.81293 -5.36884 -3.74872 H 0 0 0.41700
+HETATM 700 O 12.30451 -3.32117 -1.24965 O 0 0 -0.83400
+HETATM 701 H 12.35406 -3.17797 -2.19478 H 0 0 0.41700
+HETATM 702 O 2.44518 -9.45385 -7.13800 O 0 0 -0.83400
+HETATM 703 H 2.86778 -10.12197 -6.59832 H 0 0 0.41700
+HETATM 704 O 4.23979 -11.92627 2.91775 O 0 0 -0.83400
+HETATM 705 H 3.80636 -12.47542 2.26445 H 0 0 0.41700
+HETATM 706 H 11.33635 5.09870 2.91120 H 0 0 0.41700
+HETATM 707 H 10.46595 4.52913 4.01074 H 0 0 0.41700
+HETATM 708 H 12.54584 -1.40951 3.32147 H 0 0 0.41700
+HETATM 709 H 11.11996 -6.32486 1.70893 H 0 0 0.41700
+HETATM 710 O -12.07484 -0.35388 -4.15307 O 0 0 -0.83400
+HETATM 711 H -12.11808 -0.74534 -5.02550 H 0 0 0.41700
+HETATM 712 H -5.84813 -9.91937 3.72280 H 0 0 0.41700
+HETATM 713 O 7.00555 -6.50894 8.17194 O 0 0 -0.83400
+HETATM 714 H 7.04667 -7.12988 7.44464 H 0 0 0.41700
+HETATM 715 O 4.80864 -0.36754 10.26863 O 0 0 -0.83400
+HETATM 716 H 4.04338 -0.92698 10.40142 H 0 0 0.41700
+HETATM 717 O -3.99341 -2.55537 -11.70666 O 0 0 -0.83400
+HETATM 718 H -4.71205 -1.94546 -11.87338 H 0 0 0.41700
+HETATM 719 O 6.55326 2.21838 9.26440 O 0 0 -0.83400
+HETATM 720 H 7.43364 1.89269 9.07707 H 0 0 0.41700
+HETATM 721 O -9.94106 7.10134 5.25800 O 0 0 -0.83400
+HETATM 722 H -10.41822 6.30088 5.03933 H 0 0 0.41700
+HETATM 723 H -2.61561 -8.80200 6.65911 H 0 0 0.41700
+HETATM 724 H -2.26078 -8.57124 8.11233 H 0 0 0.41700
+HETATM 725 O -6.62490 7.96884 -6.86801 O 0 0 -0.83400
+HETATM 726 H -6.03650 7.33278 -7.27475 H 0 0 0.41700
+HETATM 727 O 5.35221 9.36969 -1.72172 O 0 0 -0.83400
+HETATM 728 O 0.54965 2.33862 -12.03489 O 0 0 -0.83400
+HETATM 729 H 0.42274 2.97633 -11.33243 H 0 0 0.41700
+HETATM 730 H 5.80872 8.79534 -4.85934 H 0 0 0.41700
+HETATM 731 H 12.21569 1.74495 -2.58184 H 0 0 0.41700
+HETATM 732 O 1.71726 -11.80876 0.34618 O 0 0 -0.83400
+HETATM 733 H 1.15252 -12.53034 0.06937 H 0 0 0.41700
+HETATM 734 H 0.86411 11.66996 3.00814 H 0 0 0.41700
+HETATM 735 H -11.18873 3.38246 5.66459 H 0 0 0.41700
+HETATM 736 H -11.37174 -4.08490 -0.97957 H 0 0 0.41700
+HETATM 737 O 12.79812 -0.94527 0.30504 O 0 0 -0.83400
+HETATM 738 O -12.33719 3.70429 -3.45617 O 0 0 -0.83400
+HETATM 739 H -12.32895 4.26469 -2.68021 H 0 0 0.41700
+HETATM 740 O 8.41208 6.69216 -5.90410 O 0 0 -0.83400
+HETATM 741 H 8.77391 6.00072 -6.45838 H 0 0 0.41700
+HETATM 742 H -1.54690 -11.96252 2.90487 H 0 0 0.41700
+HETATM 743 H -8.62042 9.40060 2.23478 H 0 0 0.41700
+HETATM 744 O -12.47600 0.25085 3.56942 O 0 0 -0.83400
+HETATM 745 H -13.05893 0.91354 3.19891 H 0 0 0.41700
+HETATM 746 O -7.35361 -7.06648 -7.75769 O 0 0 -0.83400
+HETATM 747 H -7.88122 -6.51271 -7.18218 H 0 0 0.41700
+HETATM 748 O -4.94059 -7.64439 -9.58505 O 0 0 -0.83400
+HETATM 749 H -5.09730 -7.16524 -10.39874 H 0 0 0.41700
+HETATM 750 O -8.93029 9.31180 -2.38273 O 0 0 -0.83400
+HETATM 751 H -9.54124 8.80207 -2.91484 H 0 0 0.41700
+HETATM 752 H 10.86679 -1.07070 6.49075 H 0 0 0.41700
+HETATM 753 O 9.28927 2.74153 7.06714 O 0 0 -0.83400
+HETATM 754 H 9.35845 3.57854 6.60795 H 0 0 0.41700
+HETATM 755 O -2.84091 -9.82108 5.20283 O 0 0 -0.83400
+HETATM 756 H -3.70069 -9.58447 4.85494 H 0 0 0.41700
+HETATM 757 H -9.29703 0.47687 -8.32544 H 0 0 0.41700
+HETATM 758 O 4.07550 10.96583 4.03457 O 0 0 -0.83400
+HETATM 759 H 3.18618 11.26360 3.84307 H 0 0 0.41700
+HETATM 760 H -2.47257 1.12325 -11.36847 H 0 0 0.41700
+HETATM 761 O 4.96964 -2.72172 -10.55999 O 0 0 -0.83400
+HETATM 762 O 12.16507 2.57057 1.77816 O 0 0 -0.83400
+HETATM 763 O 1.76471 -10.72125 5.49636 O 0 0 -0.83400
+HETATM 764 O 7.64671 -8.12743 5.87671 O 0 0 -0.83400
+HETATM 765 H 7.61827 -9.05112 6.12614 H 0 0 0.41700
+HETATM 766 H 4.81797 7.60372 7.48963 H 0 0 0.41700
+HETATM 767 O -6.96430 -8.98333 -5.62978 O 0 0 -0.83400
+HETATM 768 H -6.97823 -9.82016 -6.09428 H 0 0 0.41700
+HETATM 769 O 3.90419 -7.97132 -8.74845 O 0 0 -0.83400
+HETATM 770 H 3.36189 -8.67506 -8.39223 H 0 0 0.41700
+HETATM 771 O -5.87456 10.86223 1.18439 O 0 0 -0.83400
+HETATM 772 H -6.72929 10.44852 1.06389 H 0 0 0.41700
+HETATM 773 O 3.63662 -4.17379 10.23475 O 0 0 -0.83400
+HETATM 774 H 10.70592 4.18295 -3.90401 H 0 0 0.41700
+HETATM 775 H -9.29960 -6.38906 -5.14724 H 0 0 0.41700
+HETATM 776 H 10.31776 -5.13830 -5.15690 H 0 0 0.41700
+END
+
diff --git a/test/library/qmmm_amber/AVE/fort.4 b/test/library/qmmm_amber/AVE/fort.4
new file mode 100644
index 00000000..0f955730
--- /dev/null
+++ b/test/library/qmmm_amber/AVE/fort.4
@@ -0,0 +1,387 @@
+Reactive CHON-2017_weak force field: Dec13 2016, 2017(CH-retrain) Source: https:
+ 41 ! Number of general parameters
+ 50.0000 !Overcoordination parameter
+ 9.5469 !Overcoordination parameter
+ 1.6725 !Valency angle conjugation parameter
+ 1.7224 !Triple bond stabilisation parameter
+ 6.8702 !Triple bond stabilisation parameter
+ 60.4850 !C2-correction
+ 1.0588 !Povun6 Undercoordination Eq(12)
+ 4.6000 !Triple bond stabilisation parameter
+ 12.1176 !Povun7 Undercoordination Eq(12)
+ 13.3056 !Povun8 Undercoordination Eq(12)
+ -20.0000 !Triple bond stabilization energy
+ 0.0000 !Lower Taper-radius
+ 10.0000 !Upper Taper-radius
+ 2.8793 !Fe dimer correction
+ 33.8667 !Valency undercoordination
+ 6.0891 !Plp1 Lone pair param in Eq(8)
+ 1.0563 !Valency angle
+ 2.0384 !Valency angle parameter
+ 6.1431 !Fe dimer correction
+ 6.9290 !Double bond/angle parameter
+ 0.3989 !Double bond/angle parameter: overcoord
+ 3.9954 !Double bond/angle parameter: overcoord
+ -2.4837 !Fe dimer correction
+ 5.7796 !Torsion/BO parameter
+ 10.0000 !Torsion overcoordination
+ 1.9487 !Torsion overcoordination
+ -1.2327 !Reserved
+ 2.1645 !Conjugation
+ 1.5591 !vdWaals shielding
+ 0.1000 !Cutoff for bond order (*100)
+ 1.7602 !Valency angle conjugation parameter
+ 0.6991 !Povun4 Over/Undercoordination Eq(11b)
+ 50.0000 !Povun3 Over/Undercoordination Eq(11b)
+ 1.8512 !Valency/lone pair parameter
+ 0.5000 !ACKS2 softness parameter
+ 20.0000 !Scale factor (d) in dispersion (LJ)
+ 5.0000 !Gauss exponent for explicit electrons
+ 0.0000 !1: disable undecoord term in val angle
+ 0.7903 !Valency angle conjugation parameter
+ 0.0000 !1: triple bond stabilization for all b
+ 0.0000 !eReax taper radius for electrons and h
+ 11 ! Nr of atoms; cov.r; valency;a.m;Rvdw;Evdw;gammaEEM;cov.r2;#
+ alfa;gammavdW;valency;Eunder;Eover;chiEEM;etaEEM;n.u.
+ cov r3;Elp;Heat inc.;bo131;bo132;bo133;softcut;n.u.
+ ov/un;val1;n.u.;val3,vval4
+ C 1.3727 4.0000 12.0000 2.0270 0.1113 0.5516 1.1706 4.0000
+ 9.2293 4.5389 4.0000 30.0000 79.5548 4.4087 7.0601 0.0000
+ 1.1168 0.0000 181.0000 14.5210 24.9431 6.7313 0.8563 0.0000
+ -6.7437 5.6329 1.0564 4.0000 2.9663 0.0000 0.0000 0.0000
+ H 0.8924 1.0000 1.0080 1.6791 0.0709 0.7390 -0.1000 1.0000
+ 8.3519 39.1732 1.0000 0.0000 121.1250 3.5442 9.3848 1.0000
+ -0.1000 0.0000 61.6606 2.8222 2.1441 0.0003 1.0698 0.0000
+ -18.1423 5.3143 1.0338 1.0000 2.8793 0.0000 0.0000 0.0000
+ O 1.2450 2.0000 15.9990 2.3396 0.1000 1.1000 1.0548 6.0000
+ 9.3187 12.5083 4.0000 37.5000 116.0768 8.5000 8.4783 2.0000
+ 0.9049 0.1000 59.0626 3.4340 0.7722 0.0021 0.9745 0.0000
+ -3.5352 3.2703 1.0493 4.0000 2.9225 0.0000 0.0000 0.0000
+ N 1.2333 3.0000 14.0000 1.9690 0.1513 1.0000 1.1748 5.0000
+ 9.6326 9.1048 4.0000 32.4838 100.0000 6.6335 7.1473 2.0000
+ 1.0433 2.3175 119.9837 0.4648 9.7074 2.2521 0.9745 0.0000
+ -7.0000 4.0000 1.0183 4.0000 2.8793 0.0000 0.0000 0.0000
+ S 1.9673 2.0000 32.0600 2.1729 0.3000 1.0336 1.5359 6.0000
+ 10.3008 4.9055 4.0000 52.9998 112.1416 6.5000 8.2545 2.0000
+ 1.4601 9.7177 71.1843 5.7487 23.2859 12.7147 0.9745 0.0000
+ -11.0000 2.7466 1.0338 6.2998 2.8793 0.0000 0.0000 0.0000
+ Mg 1.8315 2.0000 24.3050 2.2464 0.1806 0.5020 1.0000 2.0000
+ 10.9186 27.1205 3.0000 38.0000 0.0000 0.9499 5.6130 0.0000
+ -1.3000 0.0000 220.0000 49.9248 0.3370 0.0000 0.0000 0.0000
+ -1.0823 2.3663 1.0564 6.0000 2.9663 0.0000 0.0000 0.0000
+ P 1.5994 3.0000 30.9738 1.7000 0.1743 1.0000 1.3000 5.0000
+ 9.1909 12.4054 5.0000 0.0000 0.0000 1.7836 6.9473 0.0000
+ -1.0000 22.0770 125.6300 0.2080 16.7535 24.5820 0.0000 0.0000
+ -8.6603 3.5000 1.0338 5.0000 2.8793 0.0000 0.0000 0.0000
+ Na 2.0300 1.0000 22.9898 2.3334 0.1481 0.8765 -1.0000 1.0000
+ 11.0000 9.8000 1.0000 0.0000 0.0000 -3.8501 5.9459 0.0000
+ -1.0000 0.0000 67.5458 100.0000 10.0000 0.2500 0.8563 0.0000
+ -2.5766 2.5000 1.0338 6.0000 2.5791 0.0000 0.0000 0.0000
+ Cu 1.9202 2.0000 63.5460 1.9221 0.2826 1.0000 0.1000 1.0000
+ 10.9889 100.0000 1.0000 0.0000 0.0000 2.7875 6.0000 0.0000
+ -1.0000 0.0000 80.7000 34.9555 0.4988 0.0000 0.8563 0.0000
+ -5.1872 3.1491 1.0000 4.0000 2.5791 0.0000 0.0000 0.0000
+ Cl 1.7140 1.0000 35.4500 1.9139 0.2000 0.3837 -1.0000 7.0000
+ 11.5345 10.1330 1.0000 0.0000 0.0000 9.9614 6.5316 0.0000
+ -1.0000 3.5750 143.1770 6.2293 5.2294 0.1542 0.8563 0.0000
+ -10.2080 2.9867 1.0338 6.2998 2.5791 0.0000 0.0000 0.0000
+ X -0.0998 2.0000 1.0080 2.0000 0.0000 1.0000 -0.1000 6.0000
+ 10.0000 2.5000 4.0000 0.0000 0.0000 8.5000 1.5000 0.0000
+ -0.1000 0.0000 -2.3700 8.7410 13.3640 0.6690 0.9745 0.0000
+ -11.0000 2.7466 1.0338 4.0000 2.8793 0.0000 0.0000 0.0000
+ 43 ! Nr of bonds; Edis1;LPpen;n.u.;pbe1;pbo5;13corr;pbo6
+ pbe2;pbo3;pbo4;n.u.;pbo1;pbo2;ovcorr
+ 1 1 91.8588 120.2415 50.6788 0.3723 -0.1000 1.0000 34.6637 0.7015
+ 6.8558 -0.1605 7.9139 1.0000 -0.0719 6.9788 1.0000 0.0000
+ 1 2 185.8889 0.0000 0.0000 -0.5074 0.0000 1.0000 6.0000 0.6970
+ 8.2974 1.0000 0.0000 1.0000 -0.0809 6.3908 0.0000 0.0000
+ 2 2 165.6186 0.0000 0.0000 -0.4194 0.0000 1.0000 6.0000 0.6197
+ 6.3245 1.0000 0.0000 1.0000 -0.0879 6.0634 0.0000 0.0000
+ 1 3 161.2413 119.3000 98.9000 -0.4044 -0.6000 1.0000 21.9000 0.8787
+ 1.4350 -0.3347 7.9000 0.9812 -0.2559 6.3865 0.0000 0.0000
+ 2 3 159.4722 0.0000 0.0000 -0.5725 0.0000 1.0000 6.0000 0.5505
+ 1.1150 1.0000 0.0000 0.0000 -0.0920 4.2790 0.0000 0.0000
+ 3 3 142.2858 85.0000 0.0000 0.2506 -0.1000 1.0000 29.7503 0.6051
+ 0.3451 -0.1055 9.0000 1.0000 -0.1225 5.5000 1.0000 0.0000
+ 1 4 194.7541 88.7334 62.9228 -1.6676 -0.1819 1.0000 28.4654 0.2106
+ 0.0228 -0.3690 6.5560 1.0000 -0.2971 4.6225 1.0000 0.0000
+ 3 4 128.8596 167.8643 40.0000 0.3819 -0.1539 1.0000 34.9972 0.1900
+ 1.0110 -0.3716 7.0805 1.0000 -0.1265 6.8843 1.0000 0.0000
+ 4 4 160.1592 82.5526 153.9884 0.4110 -0.0934 1.0000 12.4304 0.5899
+ 0.1538 -0.1473 11.9187 1.0000 -0.0753 5.4371 1.0000 0.0000
+ 2 4 174.3736 0.0000 0.0000 -0.3003 0.0000 1.0000 6.0000 0.5976
+ 0.5289 1.0000 0.0000 1.0000 -0.1960 5.2138 0.0000 0.0000
+ 1 5 150.8132 59.3363 55.2528 -0.0628 -0.5211 1.0000 18.9617 0.3219
+ 0.3317 -0.2289 7.5946 1.0000 -0.1946 5.9455 1.0000 0.0000
+ 2 5 143.4377 0.0000 0.0000 -0.2944 0.0000 1.0000 6.0000 0.6034
+ 9.5627 1.0000 0.0000 1.0000 -0.0516 7.0960 1.0000 0.0000
+ 3 5 0.0000 0.0000 0.0000 0.5563 -0.4038 1.0000 49.5611 0.6000
+ 0.4259 -0.4577 12.7569 1.0000 -0.1100 7.1145 1.0000 0.0000
+ 4 5 0.0000 0.0000 0.0000 0.4438 -0.2034 1.0000 40.3399 0.6000
+ 0.3296 -0.3153 9.1227 1.0000 -0.1805 5.6864 1.0000 0.0000
+ 5 5 140.8887 84.9350 68.6860 -0.4111 -0.4781 1.0000 17.8574 -0.1336
+ 0.2881 -0.2494 9.8436 1.0000 -0.1806 7.4732 1.0000 0.0000
+ 2 6 58.6896 0.0000 0.0000 -0.0203 -0.1418 1.0000 13.1260 0.0230
+ 8.2136 -0.1310 0.0000 1.0000 -0.2692 6.4254 0.0000 24.4461
+ 3 6 87.0227 0.0000 43.3991 0.0030 -0.3000 1.0000 36.0000 0.0250
+ 0.0087 -0.2500 12.0000 1.0000 -0.0439 6.6073 1.0000 24.4461
+ 6 6 32.3808 0.0000 0.0000 -0.0076 -0.2000 0.0000 16.0000 0.2641
+ 4.8726 -0.2000 10.0000 1.0000 -0.0729 4.6319 0.0000 0.0000
+ 1 7 0.0000 0.0000 0.0000 0.2171 -0.1418 1.0000 13.1260 0.6000
+ 0.3601 -0.2500 20.0000 1.0000 -0.2000 10.0000 1.0000 0.0000
+ 2 7 0.0000 0.0000 0.0000 0.2250 -0.1418 1.0000 13.1260 0.6000
+ 0.3912 -0.1310 0.0000 1.0000 -0.2000 10.0000 0.0000 0.0000
+ 3 7 250.0000 186.6645 0.0000 0.0894 -0.5000 1.0000 25.0000 0.4401
+ 0.7286 -0.2141 15.9041 1.0000 -0.1916 6.4133 1.0000 0.0000
+ 4 7 130.0000 0.0000 0.0000 0.2171 -0.1418 1.0000 13.1260 0.6000
+ 0.3601 -0.1310 10.7257 1.0000 -0.0869 5.3302 1.0000 0.0000
+ 6 7 0.1000 0.0000 0.0000 0.2500 -0.5000 1.0000 35.0000 0.6000
+ 0.5000 -0.5000 20.0000 1.0000 -0.2000 10.0000 1.0000 0.0000
+ 7 7 0.0000 0.0000 0.0000 0.2171 -0.5000 1.0000 35.0000 0.6000
+ 0.5000 -0.5000 20.0000 1.0000 -0.2000 10.0000 1.0000 0.0000
+ 2 8 0.0000 0.0000 0.0000 -1.0000 -0.3000 1.0000 36.0000 0.7000
+ 10.1151 -0.3500 25.0000 1.0000 -0.1053 8.2003 1.0000 0.0000
+ 3 8 76.0753 0.0000 0.0000 -0.4452 -0.3000 1.0000 36.0000 0.6433
+ 5.6834 -0.3500 25.0000 1.0000 -0.0539 8.0273 1.0000 0.0000
+ 4 8 0.0000 0.0000 0.0000 -1.0000 -0.3000 1.0000 36.0000 0.7000
+ 10.1151 -0.3500 25.0000 1.0000 -0.1053 8.2003 1.0000 0.0000
+ 6 8 0.1000 0.0000 0.0000 0.2500 -0.5000 1.0000 35.0000 0.6000
+ 0.5000 -0.5000 20.0000 1.0000 -0.2000 10.0000 1.0000 0.0000
+ 7 8 0.1000 0.0000 0.0000 0.2500 -0.5000 1.0000 35.0000 0.6000
+ 0.5000 -0.5000 20.0000 1.0000 -0.2000 10.0000 1.0000 0.0000
+ 8 8 27.8052 0.0000 0.0000 0.4022 0.3000 0.0000 25.0000 0.4894
+ 0.6222 -0.4000 12.0000 1.0000 -0.0500 5.3362 0.0000 0.0000
+ 4 6 0.0000 0.0000 0.0000 -1.0000 -0.3000 1.0000 36.0000 0.7000
+ 10.1151 -0.3500 25.0000 1.0000 -0.1053 8.2003 1.0000 0.0000
+ 1 9 0.0000 0.0000 0.0000 0.2000 -0.1418 1.0000 13.1260 0.5000
+ 0.5000 -0.2000 20.0000 1.0000 -0.1000 9.0000 0.0000 0.0000
+ 2 9 0.0000 0.0000 0.0000 0.2000 -0.1418 1.0000 13.1260 0.5000
+ 0.5000 -0.2000 20.0000 1.0000 -0.1000 9.0000 0.0000 0.0000
+ 3 9 81.4346 0.0000 0.0000 -0.1594 -0.3000 1.0000 36.0000 0.0025
+ 0.2904 -0.2500 12.0000 1.0000 -0.0742 9.3638 0.0000 0.0000
+ 4 9 96.5322 0.0000 0.0000 0.9970 -0.3000 1.0000 36.0000 0.5095
+ 0.7247 -0.2500 12.0000 1.0000 -0.1175 9.9985 0.0000 0.0000
+ 9 9 73.6263 0.0000 0.0000 0.0209 -0.2000 0.0000 16.0000 0.3414
+ 0.4703 -0.2000 15.0000 1.0000 -0.1319 5.9254 0.0000 0.0000
+ 2 10 109.1686 0.0000 0.0000 -0.1657 -0.2000 0.0000 16.0000 1.2500
+ 2.8463 -0.2000 15.0000 1.0000 -0.1111 5.2687 0.0000 0.0000
+ 3 10 0.0000 0.0000 0.0000 0.5000 -0.2000 0.0000 16.0000 0.5000
+ 1.0001 -0.2000 15.0000 1.0000 -0.1000 10.0000 0.0000 0.0000
+ 9 10 118.3052 0.0000 0.0000 -0.1168 -0.2000 0.0000 16.0000 0.0697
+ 2.9176 -0.2000 15.0000 1.0000 -0.1316 5.3624 0.0000 0.0000
+ 10 10 0.2500 0.0000 0.0000 0.1803 -0.2000 0.0000 16.0000 0.3356
+ 0.9228 -0.2000 15.0000 1.0000 -0.1178 5.6715 0.0000 0.0000
+ 1 8 0.0000 0.0000 0.0000 -1.0000 -0.3000 1.0000 36.0000 0.7000
+ 10.1151 -0.3500 25.0000 1.0000 -0.1053 8.2003 1.0000 0.0000
+ 1 10 0.0000 0.0000 0.0000 0.5000 -0.2000 0.0000 16.0000 0.5000
+ 1.0001 -0.2000 15.0000 1.0000 -0.1000 10.0000 0.0000 0.0000
+ 4 10 0.0000 0.0000 0.0000 0.5000 -0.2000 0.0000 16.0000 0.5000
+ 1.0001 -0.2000 15.0000 1.0000 -0.1000 10.0000 0.0000 0.0000
+ 24 ! Nr of off-diagonal terms; Ediss;Ro;gamma;rsigma;rpi;rpi2
+ 1 2 0.1165 1.3851 9.9415 1.1475 -1.0000 -1.0000
+ 1 3 0.0697 2.0716 9.6540 1.4310 1.1591 1.0781
+ 2 3 0.0283 1.2853 10.9173 0.9259 -1.0000 -1.0000
+ 1 4 0.1261 1.9631 9.7006 1.3387 1.2586 1.1688
+ 2 4 0.1154 1.4616 9.1198 1.0731 -1.0000 -1.0000
+ 3 4 0.1819 1.8000 9.4170 1.4262 1.0546 1.2607
+ 1 5 0.1618 1.7943 10.1042 1.7489 1.3150 1.4031
+ 2 5 0.0764 1.5838 10.1462 1.4206 -1.0000 -1.0000
+ 3 5 0.1022 1.9887 10.0605 1.5799 1.4000 -1.0000
+ 4 5 0.1505 1.9000 10.5104 1.8000 1.4000 -1.0000
+ 2 6 0.0100 1.6000 13.2979 1.8670 -1.0000 -1.0000
+ 3 6 0.0809 1.7000 11.4606 1.5177 -1.0000 -1.0000
+ 1 7 0.2221 1.7777 10.0383 -1.0000 -1.0000 -1.0000
+ 2 7 0.3000 1.9226 10.0193 -1.0000 -1.0000 -1.0000
+ 3 7 0.1541 1.8743 10.1246 1.7386 1.4412 -1.0000
+ 6 7 0.1801 1.8566 9.8498 0.1000 -1.0000 -1.0000
+ 3 8 0.1592 1.8283 11.7256 1.6655 -1.0000 -1.0000
+ 1 9 0.0500 1.7500 12.3500 0.1000 -1.0000 -1.0000
+ 2 9 0.0300 1.5200 12.5000 0.1000 -1.0000 -1.0000
+ 3 9 0.0348 1.7637 12.3562 1.7228 -1.0000 -1.0000
+ 4 9 0.0478 1.7704 12.8051 1.6100 -1.0000 -1.0000
+ 2 10 0.0568 1.6740 9.6297 1.2200 -1.0000 -1.0000
+ 3 10 0.1927 2.2551 11.2308 -1.0000 -1.0000 -1.0000
+ 9 10 0.1402 2.1604 10.9786 1.7505 -1.0000 -1.0000
+105 ! Nr of angles;at1;at2;at3;Thetao,o;ka;kb;pv1;pv2
+ 1 1 1 66.0347 45.0000 1.1248 0.0000 3.0000 70.0000 2.4142
+ 1 1 2 57.5203 11.1542 4.7426 -0.7339 -0.9176 0.0000 1.7607
+ 2 1 2 64.1492 18.6633 1.9552 0.0000 3.0000 0.0000 1.0010
+ 1 2 2 0.0000 0.0000 6.0000 0.0000 0.0000 0.0000 1.0400
+ 1 2 1 0.0000 7.5000 5.0000 0.0000 0.0000 0.0000 1.0400
+ 2 2 2 0.0000 27.9213 5.8635 0.0000 0.0000 0.0000 1.0400
+ 1 1 3 68.6071 46.1858 2.5058 -0.4029 -0.1784 58.6562 1.5863
+ 3 1 3 70.8589 10.0656 7.5000 -17.8180 1.1419 0.0000 1.4986
+ 1 1 4 54.9449 38.0899 0.5948 0.0000 1.4838 0.0000 1.5862
+ 3 1 4 86.4061 17.5133 2.2425 0.0000 0.4036 32.5000 1.3710
+ 4 1 4 69.6085 2.5000 0.1000 0.0000 0.9683 0.0000 1.5523
+ 2 1 3 54.7474 4.8106 9.4520 -30.0000 -1.0000 0.0000 4.2208
+ 2 1 4 75.3048 9.9497 2.4268 0.0000 0.1000 0.0000 1.5478
+ 1 2 4 0.0000 0.0019 6.3000 0.0000 0.0000 0.0000 1.0400
+ 1 3 1 -57.4049 4.5141 11.8498 2.0995 4.0000 0.0000 1.8452
+ 1 3 3 82.8564 45.0000 1.5233 0.0000 1.1034 68.1072 1.6192
+ 1 3 4 87.0300 18.3166 1.0893 0.0000 3.0479 0.0000 1.0825
+ 3 3 3 90.0000 21.1742 1.8914 0.0000 1.6310 50.0000 1.0000
+ 3 3 4 83.9925 27.3343 1.0888 0.0000 2.6572 0.0000 1.2246
+ 4 3 4 73.3858 43.0000 1.1230 0.0000 3.0072 0.0000 1.0000
+ 1 3 2 90.0000 9.5520 1.7468 0.0000 1.2146 0.0000 3.0000
+ 2 3 3 80.0453 50.0000 2.5336 0.0000 3.0000 0.0000 1.0821
+ 2 3 4 75.6278 14.6992 5.3049 0.0000 0.2502 0.0000 2.5848
+ 2 3 2 85.7876 9.3298 2.0747 0.0000 2.8632 0.0000 1.6905
+ 1 4 1 81.5674 13.5526 1.5408 0.0000 2.7249 0.0000 2.3448
+ 1 4 3 88.9691 8.9607 4.6383 0.0000 2.1900 0.0000 3.0000
+ 1 4 4 70.1079 13.1213 2.4301 0.0000 2.9290 0.0000 3.0000
+ 3 4 3 74.5618 45.0000 1.0990 -18.0069 3.0584 0.0000 1.0000
+ 3 4 4 79.9166 45.0000 0.8725 -0.9193 3.0760 0.0000 1.0000
+ 4 4 4 76.2179 30.5403 1.5752 0.0000 2.9239 0.0000 2.5095
+ 1 4 2 78.9698 8.3913 1.8560 0.0000 0.2894 0.0000 1.0108
+ 2 4 3 74.5764 10.5782 6.5751 0.0000 0.4346 0.0000 3.0000
+ 2 4 4 90.0000 41.1335 1.3912 0.0000 0.1000 0.0000 1.0596
+ 2 4 2 83.2464 11.7793 5.0092 0.0000 0.1777 0.0000 1.0711
+ 1 2 3 0.0000 16.8357 0.9686 0.0000 0.6225 0.0000 1.0384
+ 1 2 4 0.0000 16.8188 1.4585 0.0000 0.2078 0.0000 1.1171
+ 1 2 5 0.0000 15.0000 3.0000 0.0000 0.0000 0.0000 1.0400
+ 3 2 3 0.0000 1.0000 2.6859 0.0000 0.0000 0.0000 4.0000
+ 3 2 4 0.0000 18.3614 4.9976 0.0000 0.1224 0.0000 1.9411
+ 4 2 4 0.0000 20.0000 5.0000 0.0000 0.1460 0.0000 1.0000
+ 2 2 3 0.0000 5.8938 3.0000 0.0000 0.0000 0.0000 1.0629
+ 2 2 4 0.0000 0.0019 6.0000 0.0000 0.0000 0.0000 1.0400
+ 1 1 5 74.4180 33.4273 1.7018 0.1463 0.5000 0.0000 1.6178
+ 1 5 1 79.7037 28.2036 1.7073 0.1463 0.5000 0.0000 1.6453
+ 2 1 5 63.3289 29.4225 2.1326 0.0000 0.5000 0.0000 3.0000
+ 1 5 2 85.9449 38.3109 1.2492 0.0000 0.5000 0.0000 1.1000
+ 1 5 5 80.0000 25.0000 2.0000 0.0000 0.5000 0.0000 1.3830
+ 2 5 2 85.0000 15.1317 2.0000 0.0000 0.5000 0.0000 2.0000
+ 2 5 5 97.0064 32.1121 2.0242 0.0000 0.5000 0.0000 2.8568
+ 2 2 5 0.0000 0.0019 6.0000 0.0000 0.0000 0.0000 1.0400
+ 5 4 5 62.0000 33.4273 1.7018 0.1463 0.5000 0.0000 1.0500
+ 3 5 3 77.0699 39.4349 2.1313 -30.0000 0.9567 0.0000 1.1483
+ 1 5 3 70.0000 35.0000 3.4223 0.0000 1.3550 0.0000 1.2002
+ 1 5 4 70.0000 35.0000 3.4223 0.0000 1.3550 0.0000 1.2002
+ 4 1 5 60.0000 35.0000 3.0223 0.0000 2.3550 0.0000 1.2002
+ 3 5 4 70.0000 35.0000 3.4223 0.0000 1.3550 0.0000 1.2002
+ 5 1 7 70.0000 35.0000 3.4223 0.0000 1.3550 0.0000 1.2002
+ 1 3 5 73.0990 33.8942 1.2098 0.0000 0.8161 0.0000 1.1776
+ 3 3 5 83.9753 31.0715 3.5590 0.0000 0.8161 0.0000 1.1776
+ 2 3 5 76.9521 20.0000 2.0903 0.0000 1.0000 0.0000 1.0400
+ 2 6 2 0.0000 49.8261 0.2093 0.0000 2.0870 0.0000 2.2895
+ 2 2 6 0.0000 39.7818 3.1505 0.0000 1.1296 0.0000 1.1110
+ 6 2 6 0.0000 0.5047 0.8000 0.0000 0.8933 0.0000 4.6650
+ 2 6 6 0.0000 8.7037 0.0827 0.0000 3.5597 0.0000 1.1198
+ 3 6 3 0.0000 9.2317 0.1000 0.0000 1.0000 0.0000 1.0920
+ 6 3 6 0.0008 25.0000 8.0000 0.0000 1.0000 0.0000 3.0000
+ 2 3 6 66.0423 5.0000 1.0000 0.0000 1.0000 0.0000 1.2500
+ 2 6 3 0.0000 0.5000 0.1000 0.0000 1.0000 0.0000 3.0000
+ 3 3 6 70.0000 20.0000 1.0000 0.0000 1.0000 0.0000 1.2500
+ 3 7 3 90.0000 6.4822 1.2338 -24.7539 1.6651 0.0000 1.6853
+ 2 3 7 100.0000 13.2421 1.3110 0.0000 3.9965 0.0000 1.0400
+ 3 3 7 60.0000 40.0000 4.0000 0.0000 1.0000 0.0000 1.0400
+ 3 2 7 0.0000 10.0000 1.0000 0.0000 1.0000 0.0000 1.0400
+ 6 3 7 41.7961 3.4184 7.4578 0.0000 -0.2958 0.0000 1.0542
+ 7 3 7 89.7500 16.6529 7.5000 0.0000 -0.0234 0.0000 2.9500
+ 1 3 7 54.8627 5.2339 5.0296 0.0000 2.8402 0.0000 2.8529
+ 2 7 3 75.0000 25.0000 2.0000 0.0000 1.0000 0.0000 1.2500
+ 3 7 7 70.0000 25.0000 2.0000 0.0000 1.0000 0.0000 1.2500
+ 3 9 3 96.2265 4.5610 12.0000 0.0000 0.3211 0.0000 1.5204
+ 3 9 3 0.0000 9.1552 7.9919 0.0000 0.1660 0.0000 1.5386
+ 9 3 9 100.0000 10.1065 6.0000 0.0000 1.0000 0.0000 3.6601
+ 2 3 9 55.0417 3.5032 3.9979 0.0000 1.5171 0.0000 1.0400
+ 3 3 9 70.0000 30.0000 2.0000 0.0000 1.0000 0.0000 1.2500
+ 3 9 9 66.7783 14.3146 0.7911 0.0000 1.0000 0.0000 1.2333
+ 3 9 10 96.6924 9.4823 5.7883 0.0000 0.2248 0.0000 2.2640
+ 3 9 10 0.0000 3.8549 3.7230 0.0000 0.1482 0.0000 1.0400
+ 9 10 9 0.0000 11.2336 6.8851 0.0000 1.0000 0.0000 1.0893
+ 9 9 10 90.0000 5.0811 5.2147 0.0000 1.0000 0.0000 1.8538
+ 10 9 10 0.0100 21.1482 0.3506 0.0000 1.0000 0.0000 1.4361
+ 3 2 10 0.0000 0.0100 0.5211 0.0000 0.0000 0.0000 1.3859
+ 3 9 4 100.0000 28.1532 12.0000 0.0000 0.2932 0.0000 1.6489
+ 3 9 4 0.0000 22.7457 2.9039 0.0000 0.5593 0.0000 1.9764
+ 4 9 4 87.0081 27.6432 3.9735 0.0000 4.0000 0.0000 1.4578
+ 4 9 4 0.0000 22.8998 3.1077 0.0000 3.0000 0.0000 1.0696
+ 9 4 9 100.0000 10.1065 6.0000 0.0000 1.0000 0.0000 3.6601
+ 2 4 9 80.0000 3.5601 3.3645 0.0000 1.5171 0.0000 1.0400
+ 3 4 9 70.0000 30.0000 2.0000 0.0000 1.0000 0.0000 1.2500
+ 4 3 9 70.0000 30.0000 2.0000 0.0000 1.0000 0.0000 1.2500
+ 4 4 9 70.0000 30.0000 2.0000 0.0000 1.0000 0.0000 1.2500
+ 4 9 9 66.7783 14.3146 0.7911 0.0000 1.0000 0.0000 1.2333
+ 4 9 10 95.2122 5.7090 12.0000 0.0000 0.2248 0.0000 2.8936
+ 4 9 10 0.0000 9.0054 7.9511 0.0000 0.1482 0.0000 1.6245
+ 4 2 10 0.0000 15.0000 2.8900 0.0000 0.0000 0.0000 2.8774
+ 1 3 9 55.0000 15.0000 1.0000 0.0000 1.0000 0.0000 1.5000
+ 1 4 9 55.0000 15.0000 1.0000 0.0000 1.0000 0.0000 1.5000
+ 66 ! Nr of torsions;at1;at2;at3;at4;;V1;V2;V3;V2(BO);vconj;n.u;n
+ 1 1 1 1 0.3619 5.0000 0.2257 -9.0000 -2.4400 0.0000 0.0000
+ 1 1 1 2 -1.9878 78.0106 -0.8993 -6.0559 0.9721 0.0000 0.0000
+ 2 1 1 2 3.4027 97.0884 1.4919 -5.4079 3.5966 0.0000 0.0000
+ 1 1 1 3 2.6596 -21.6776 -0.8522 -4.1011 -4.0000 0.0000 0.0000
+ 2 1 1 3 -0.5174 94.7613 1.5000 -9.5447 -4.0000 0.0000 0.0000
+ 3 1 1 3 -2.3605 29.3576 -1.0000 -2.5000 -0.8614 0.0000 0.0000
+ 1 1 3 1 -2.5000 -21.9374 1.5000 -1.5000 -3.5912 0.0000 0.0000
+ 1 1 3 2 -2.4164 5.0000 -0.8158 -9.0000 -0.9000 0.0000 0.0000
+ 2 1 3 1 -2.5000 14.0094 0.1288 -1.1453 -3.5751 0.0000 0.0000
+ 2 1 3 2 -2.4308 80.0000 1.0000 -3.9796 -1.1000 0.0000 0.0000
+ 2 1 3 3 -0.7523 100.0000 1.0000 -4.2517 -2.8274 0.0000 0.0000
+ 3 1 3 1 -2.5000 5.2291 -1.0000 -7.9957 -3.0437 0.0000 0.0000
+ 3 1 3 2 2.5000 21.4033 -0.7653 -2.5000 -3.0476 0.0000 0.0000
+ 1 3 3 2 -2.5000 8.9548 -1.0000 -2.7194 -2.9498 0.0000 0.0000
+ 2 3 3 2 -2.5000 -25.0000 -1.0000 -2.5000 0.0000 0.0000 0.0000
+ 1 3 3 3 2.4980 -23.2133 0.9757 -2.5000 -0.9972 0.0000 0.0000
+ 2 3 3 3 -0.5841 58.4202 -0.6677 -6.7360 0.0000 0.0000 0.0000
+ 3 3 3 3 -2.5000 -24.6965 1.0000 -2.5551 -0.9000 0.0000 0.0000
+ 1 1 4 2 0.0244 75.4064 0.3846 -7.5031 -1.9825 0.0000 0.0000
+ 2 1 4 2 -0.7794 48.6352 0.2610 -5.5335 -2.1051 0.0000 0.0000
+ 3 1 4 2 1.0000 20.4656 1.0000 -3.8676 -2.5261 0.0000 0.0000
+ 3 1 1 4 -0.9880 8.6101 1.0000 -2.5000 -0.9511 0.0000 0.0000
+ 4 1 1 4 -0.3655 20.5499 -0.5000 -3.4432 -1.7241 0.0000 0.0000
+ 1 1 4 1 -0.2067 34.6291 0.2986 -5.5465 -1.6589 0.0000 0.0000
+ 3 1 4 1 -0.8383 20.7675 1.0000 -8.0000 -1.8038 0.0000 0.0000
+ 2 1 1 4 1.0000 0.0100 0.6954 -2.5000 -1.9000 0.0000 0.0000
+ 4 1 4 2 0.6604 53.3623 -0.3473 -4.0226 -2.0202 0.0000 0.0000
+ 2 1 4 1 0.7622 70.4937 0.2270 -5.3819 -0.0965 0.0000 0.0000
+ 4 1 3 1 -0.0727 18.9201 -0.5000 -2.5000 -3.0437 0.0000 0.0000
+ 4 1 3 2 -0.4860 -0.4199 -0.5000 -6.7056 -3.0476 0.0000 0.0000
+ 1 1 1 4 0.9135 -5.0000 1.0000 -6.9809 -2.0000 0.0000 0.0000
+ 4 1 4 1 -0.8907 11.2771 0.6391 -6.8467 -2.0051 0.0000 0.0000
+ 0 1 2 0 0.0000 0.0000 0.0000 -4.0000 0.0000 0.0000 0.0000
+ 0 2 2 0 0.0000 0.0000 0.0000 -4.0000 0.0000 0.0000 0.0000
+ 0 2 3 0 0.0000 0.1000 0.0200 -4.0000 0.0000 0.0000 0.0000
+ 0 1 1 0 0.0000 50.0000 0.3000 -4.0000 -2.0000 0.0000 0.0000
+ 0 3 3 0 0.5511 25.4150 1.1330 -5.1903 -1.0000 0.0000 0.0000
+ 0 1 4 0 0.2176 40.4126 0.3535 -3.9875 -2.0051 0.0000 0.0000
+ 0 2 4 0 0.0000 0.1000 0.0100 -5.0000 0.0000 0.0000 0.0000
+ 0 3 4 0 1.1397 61.3225 0.5139 -3.8507 -2.7831 0.0000 0.0000
+ 0 4 4 0 0.7265 44.3155 1.0000 -4.4046 -2.0000 0.0000 0.0000
+ 4 1 4 4 -0.0949 8.7582 0.3310 -7.9430 -2.0000 0.0000 0.0000
+ 0 1 5 0 0.8251 92.1468 0.7176 -4.2341 0.0000 0.0000 0.0000
+ 0 5 5 0 0.1291 -5.0000 0.9649 -5.0903 0.0000 0.0000 0.0000
+ 0 2 5 0 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
+ 0 6 6 0 0.0000 0.0000 0.1200 -2.4426 0.0000 0.0000 0.0000
+ 0 2 6 0 0.0000 0.0000 0.1200 -2.4847 0.0000 0.0000 0.0000
+ 0 3 6 0 0.0000 0.0000 0.1200 -2.4703 0.0000 0.0000 0.0000
+ 1 1 3 3 -0.0002 20.1851 0.1601 -9.0000 -2.0000 0.0000 0.0000
+ 1 3 3 1 0.0002 80.0000 -1.5000 -4.4848 -2.0000 0.0000 0.0000
+ 3 1 3 3 -0.1583 20.0000 1.5000 -9.0000 -2.0000 0.0000 0.0000
+ 2 3 9 3 -1.5000 6.8333 -0.1978 -1.4683 0.0000 0.0000 0.0000
+ 2 3 9 4 -0.6181 7.1542 -0.0047 -1.6577 0.0000 0.0000 0.0000
+ 2 4 9 3 -1.5000 1.7820 -1.0000 -5.4916 0.0000 0.0000 0.0000
+ 2 4 9 4 -0.1959 2.3626 -1.0000 -3.0702 0.0000 0.0000 0.0000
+ 2 1 4 9 0.0000 10.0000 0.3000 -6.0000 -1.0000 0.0000 0.0000
+ 2 3 9 10 0.1589 12.5000 0.4388 -1.5000 0.0000 0.0000 0.0000
+ 7 1 1 7 0.3063 77.0574 0.7057 -7.6050 -1.7255 0.0000 0.0000
+ 0 1 7 0 1.2722 34.8103 0.8518 -2.5118 0.0000 0.0000 0.0000
+ 0 7 7 0 -0.7934 35.6754 0.2832 -7.2724 0.0000 0.0000 0.0000
+ 1 1 1 7 1.4216 19.7387 0.2840 -10.4708 -1.7255 0.0000 0.0000
+ 2 1 3 7 -1.5000 99.6206 0.0516 -8.1129 0.0000 0.0000 0.0000
+ 2 3 7 3 -1.5000 -1.0000 0.2523 -2.7183 0.0000 0.0000 0.0000
+ 1 3 7 3 -1.5000 50.2166 0.9871 -3.6395 0.0000 0.0000 0.0000
+ 7 3 7 3 1.5000 -0.9128 0.0100 -6.2667 0.0000 0.0000 0.0000
+ 1 1 3 7 0.0000 0.1000 0.0100 -4.0000 0.0000 0.0000 0.0000
+ 9 ! Nr of hydrogen bonds;at1;at2;at3;Rhb;Dehb;vhb1
+ 3 2 3 2.1493 -3.8387 1.5230 28.0017
+ 3 2 4 2.0315 -6.5000 1.4500 19.5000
+ 4 2 3 1.7500 -2.3430 1.4500 19.5000
+ 4 2 4 1.8992 -2.6374 1.4500 19.5000
+ 3 2 5 1.5000 -2.0000 1.4500 19.5000
+ 4 2 5 1.5000 -2.0000 1.4500 19.5000
+ 5 2 3 1.5000 -2.0000 1.4500 19.5000
+ 5 2 4 1.5000 -2.0000 1.4500 19.5000
+ 5 2 5 1.5000 -2.0000 1.4500 19.5000
diff --git a/test/library/qmmm_amber/control b/test/library/qmmm_amber/control
new file mode 100644
index 00000000..c1631358
--- /dev/null
+++ b/test/library/qmmm_amber/control
@@ -0,0 +1,69 @@
+simulation_name AVE ! output files will carry this name + their specific extension
+ensemble_type 0 ! 0: NVE, 1: Berendsen NVT, 2: nose-Hoover NVT, 3: semi-isotropic NPT, 4: isotropic NPT, 5: anisotropic NPT
+nsteps 0 ! number of simulation steps
+dt 0.25 ! time step in fs
+periodic_boundaries 1 ! 0: no periodic boundaries, 1: periodic boundaries
+
+reposition_atoms 0 ! 0: just fit to periodic boundaries, 1: CoM to the center of box, 3: CoM to the origin
+restrict_bonds 0 ! enforce the bonds given in CONECT lines of pdb file for this many steps
+tabulate_long_range 0 ! denotes the granularity of long range tabulation, 0 means no tabulation
+energy_update_freq 1
+remove_CoM_vel 500 ! remove the translational and rotational vel around the center of mass at every 'this many' steps
+
+nbrhood_cutoff 5.0 ! near neighbors cutoff for bond calculations (Angstroms)
+bond_graph_cutoff 0.3 ! bond strength cutoff for bond graphs (Angstroms)
+thb_cutoff 0.001 ! cutoff value for three body interactions (Angstroms)
+hbond_cutoff 7.50 ! cutoff distance for hydrogen bond interactions (Angstroms)
+
+charge_method 0 ! charge method: 0 = QEq, 1 = EEM, 2 = ACKS2
+cm_q_net 0.0 ! net system charge
+cm_solver_type 2 ! iterative linear solver for charge method: 0 = GMRES(k), 1 = GMRES_H(k), 2 = CG, 3 = SDM
+cm_solver_max_iters 20 ! max solver iterations
+cm_solver_restart 100 ! inner iterations of before restarting (GMRES(k)/GMRES_H(k))
+cm_solver_q_err 1.0e-14 ! relative residual norm threshold used in solver
+cm_domain_sparsity 1.0 ! scalar for scaling cut-off distance, used to sparsify charge matrix (between 0.0 and 1.0)
+cm_init_guess_extrap1 3 ! order of spline extrapolation for initial guess (s)
+cm_init_guess_extrap2 2 ! order of spline extrapolation for initial guess (t)
+cm_solver_pre_comp_type 1 ! method used to compute preconditioner, if applicable
+cm_solver_pre_comp_refactor 1000 ! number of steps before recomputing preconditioner (-1 for dynamic refactoring)
+cm_solver_pre_comp_droptol 0.0 ! threshold tolerance for dropping values in preconditioner computation (ICHOLT/ILUT/FG-ILUT)
+cm_solver_pre_comp_sweeps 3 ! number of sweeps used to compute preconditioner (FG-ILUT)
+cm_solver_pre_comp_sai_thres 0.1 ! ratio of charge matrix NNZ's used to compute preconditioner (SAI)
+cm_solver_pre_app_type 1 ! method used to apply preconditioner (ICHOLT/ILUT/FG-ILUT)
+cm_solver_pre_app_jacobi_iters 50 ! num. Jacobi iterations used for applying precondition (ICHOLT/ILUT/FG-ILUT)
+
+temp_init 0.0 ! desired initial temperature of the simulated system
+temp_final 300.0 ! desired final temperature of the simulated system
+t_mass 0.16666 ! 0.16666 for Nose-Hoover nvt ! 100.0 for npt! in fs, thermal inertia parameter
+t_mode 0 ! 0: T-coupling only, 1: step-wise, 2: constant slope
+t_rate -100.0 ! in K
+t_freq 4.0 ! in ps
+
+pressure 0.000101325 ! desired pressure of the simulated system in GPa, 1atm = 0.000101325 GPa
+p_mass 5000.00 ! in fs, pressure inertia parameter
+compress 0.008134 ! in ps^2 * A / amu ( 4.5X10^(-5) bar^(-1) )
+press_mode 0 ! 0: internal + external pressure, 1: ext only, 2: int only
+
+geo_format 1 ! 0: custom, 1: pdb, 2: bgf 3: ASCII restart 3: binary restart
+write_freq 0 ! write trajectory after so many steps
+traj_compress 0 ! 0: no compression 1: uses zlib to compress trajectory output
+traj_format 0 ! 0: our own format (below options apply to this only), 1: xyz, 2: bgf, 3: pdb
+traj_title AVE ! (no white spaces)
+atom_info 1 ! 0: no atom info, 1: print basic atom info in the trajectory file
+atom_forces 0 ! 0: basic atom format, 1: print force on each atom in the trajectory file
+atom_velocities 0 ! 0: basic atom format, 1: print the velocity of each atom in the trajectory file
+bond_info 1 ! 0: do not print bonds, 1: print bonds in the trajectory file
+angle_info 1 ! 0: do not print angles, 1: print angles in the trajectory file
+test_forces 0 ! 0: normal run, 1: at every timestep print each force type into a different file
+
+molec_anal 0 ! 1: outputs newly formed molecules as the simulation progresses
+freq_molec_anal 0 ! perform molecular analysis at every 'this many' timesteps
+dipole_anal 0 ! 1: calculate a electric dipole moment of the system
+freq_dipole_anal 1 ! calculate electric dipole moment at every 'this many' steps
+diffusion_coef 0 ! 1: calculate diffusion coefficient of the system
+freq_diffusion_coef 1 ! calculate diffusion coefficient at every 'this many' steps
+restrict_type 2 ! -1: all types of atoms, 0 and up: only this type of atoms
+
+restart_format 1 ! 0: restarts in ASCII 1: restarts in binary
+restart_freq 0 ! 0: do not output any restart files. >0: output a restart file at every 'this many' steps
+
diff --git a/test/library/qmmm_amber/tester_AVE.py b/test/library/qmmm_amber/tester_AVE.py
new file mode 100644
index 00000000..d9cfacd8
--- /dev/null
+++ b/test/library/qmmm_amber/tester_AVE.py
@@ -0,0 +1,520 @@
+#!/bin/python3
+
+from ctypes import c_int, c_double, c_char, c_char_p, c_void_p, \
+ Structure, Union, POINTER, CFUNCTYPE, cdll
+import numpy as np
+import sqlite3 as sq3
+import pandas as pd
+from os import path
+
+
+class BondOrderData(Structure):
+ _fields_ = [
+ ("BO", c_double),
+ ("BO_s", c_double),
+ ("BO_pi", c_double),
+ ("BO_pi2", c_double),
+ ("Cdbo", c_double),
+ ("Cdbopi", c_double),
+ ("Cdbopi2", c_double),
+ ("C1dbo", c_double),
+ ("C2dbo", c_double),
+ ("C3dbo", c_double),
+ ("C1dbopi", c_double),
+ ("C2dbopi", c_double),
+ ("C3dbopi", c_double),
+ ("C4dbopi", c_double),
+ ("C1dbopi2", c_double),
+ ("C2dbopi2", c_double),
+ ("C3dbopi2", c_double),
+ ("C4dbopi2", c_double),
+ ("dBOp", c_double * 3),
+ ("dln_BOp_s", c_double * 3),
+ ("dln_BOp_pi", c_double * 3),
+ ("dln_BOp_pi2", c_double * 3),
+ ]
+
+
+class ThreeBodyData(Structure):
+ _fields_ = [
+ ("thb", c_int),
+ ("pthb", c_int),
+ ("theta", c_double),
+ ("cos_theta", c_double),
+ ("dcos_di", c_double * 3),
+ ("dcos_dj", c_double * 3),
+ ("dcos_dk", c_double * 3),
+ ]
+
+
+class BondData(Structure):
+ _fields_ = [
+ ("nbr", c_int),
+ ("sym_index", c_int),
+ ("dbond_index", c_int),
+ ("rel_box", c_int * 3),
+ ("d", c_double),
+ ("dvec", c_double * 3),
+ ("bo_data", BondOrderData),
+ ]
+
+
+class DBondData(Structure):
+ _fields_ = [
+ ("wrt", c_int),
+ ("dBO", c_double * 3),
+ ("dBOpi", c_double * 3),
+ ("dBOpi2", c_double * 3),
+ ]
+
+
+class DDeltaData(Structure):
+ _fields_ = [
+ ("wrt", c_int),
+ ("dVal", c_double * 3),
+ ]
+
+
+class FarNbrData(Structure):
+ _fields_ = [
+ ("nbr", c_int),
+ ("rel_box", c_int * 3),
+ ("d", c_double),
+ ("dvec", c_double * 3),
+ ]
+
+
+class NearNbrData(Structure):
+ _fields_ = [
+ ("nbr", c_int),
+ ("rel_box", c_int * 3),
+ ("d", c_double),
+ ("dvec", c_double * 3),
+ ]
+
+
+class HBondData(Structure):
+ _fields_ = [
+ ("nbr", c_int),
+ ("scl", c_int),
+ ("ptr", POINTER(FarNbrData)),
+ ]
+
+
+class Thermostat(Structure):
+ _fields_ = [
+ ("T", c_double),
+ ("xi", c_double),
+ ("v_xi", c_double),
+ ("v_xi_old", c_double),
+ ("G_xi", c_double),
+ ]
+
+
+class IsotropicBarostat(Structure):
+ _fields_ = [
+ ("P", c_double),
+ ("eps", c_double),
+ ("v_eps", c_double),
+ ("v_eps_old", c_double),
+ ("a_eps", c_double),
+ ]
+
+
+class FlexibleBarostat(Structure):
+ _fields_ = [
+ ("P", c_double * 9),
+ ("P_scalar", c_double),
+ ("eps", c_double),
+ ("v_eps", c_double),
+ ("v_eps_old", c_double),
+ ("a_eps", c_double),
+ ("h0", c_double * 9),
+ ("v_g0", c_double * 9),
+ ("v_g0_old", c_double * 9),
+ ("a_g0", c_double * 9),
+ ]
+
+
+class ReaxTiming(Structure):
+ _fields_ = [
+ ("start", c_double),
+ ("end", c_double),
+ ("elapsed", c_double),
+ ("total", c_double),
+ ("nbrs", c_double),
+ ("init_forces", c_double),
+ ("bonded", c_double),
+ ("nonb", c_double),
+ ("cm", c_double),
+ ("cm_sort_mat_rows", c_double),
+ ("cm_solver_pre_comp", c_double),
+ ("cm_solver_pre_app", c_double),
+ ("cm_solver_iters", c_int),
+ ("cm_solver_spmv", c_double),
+ ("cm_solver_vector_ops", c_double),
+ ("cm_solver_orthog", c_double),
+ ("cm_solver_tri_solve", c_double),
+ ("cm_last_pre_comp", c_double),
+ ("cm_total_loss", c_double),
+ ("cm_optimum", c_double),
+ ("num_retries", c_int),
+ ]
+
+
+class SimulationData(Structure):
+ _fields_ = [
+ ("sim_id", c_int),
+ ("step", c_int),
+ ("prev_step", c_int),
+ ("time", c_double),
+ ("M", c_double),
+ ("inv_M", c_double),
+ ("xcm", c_double * 3),
+ ("vcm", c_double * 3),
+ ("fcm", c_double * 3),
+ ("amcm", c_double * 3),
+ ("avcm", c_double * 3),
+ ("etran_cm", c_double),
+ ("erot_cm", c_double),
+ ("kinetic", c_double * 9),
+ ("virial", c_double * 9),
+ ("E_Tot", c_double),
+ ("E_Kin", c_double),
+ ("E_Pot", c_double),
+ ("E_BE", c_double),
+ ("E_Ov", c_double),
+ ("E_Un", c_double),
+ ("E_Lp", c_double),
+ ("E_Ang", c_double),
+ ("E_Pen", c_double),
+ ("E_Coa", c_double),
+ ("E_HB", c_double),
+ ("E_Tor", c_double),
+ ("E_Con", c_double),
+ ("E_vdW", c_double),
+ ("E_Ele", c_double),
+ ("E_Pol", c_double),
+ ("N_f", c_double),
+ ("t_scale", c_double * 3),
+ ("p_scale", c_double * 9),
+ ("therm", Thermostat),
+ ("iso_bar", IsotropicBarostat),
+ ("flex_bar", FlexibleBarostat),
+ ("inv_W", c_double),
+ ("int_press", c_double * 3),
+ ("ext_press", c_double * 3),
+ ("kin_press", c_double),
+ ("tot_press", c_double * 3),
+ ("timing", ReaxTiming),
+ ]
+
+
+class ReaxAtom(Structure):
+ _fields_ = [
+ ("type", c_int),
+ ("rel_map", c_int * 3),
+ ("name", c_char * 9),
+ ("x", c_double * 3),
+ ("v", c_double * 3),
+ ("f", c_double * 3),
+ ("q", c_double),
+ ]
+
+
+def create_db(name='spuremd.db'):
+ conn = sq3.connect(name)
+
+ conn.executescript("""
+ CREATE TABLE simulation(
+ id integer,
+ date text,
+ name text,
+ ensemble_type integer,
+ steps integer,
+ time_step integer,
+ restart_format integer,
+ random_velocity integer,
+ reposition_atoms integer,
+ peroidic_boundary integer,
+ geo_format integer,
+ restrict_bonds integer,
+ tabulate_long_range integer,
+ reneighbor integer,
+ vlist_cutoff real,
+ neighbor_cutoff real,
+ three_body_cutoff real,
+ hydrogen_bond_cutoff real,
+ bond_graph_cutoff real,
+ charge_method integer,
+ cm_q_net real,
+ cm_solver_type integer,
+ cm_solver_max_iters integer,
+ cm_solver_restart integer,
+ cm_solver_q_err real,
+ cm_domain_sparsity real,
+ cm_solver_pre_comp_type integer,
+ cm_solver_pre_comp_refactor integer,
+ cm_solver_pre_comp_droptol real,
+ cm_solver_pre_comp_sweeps integer,
+ cm_solver_pre_comp_sai_thres real,
+ cm_solver_pre_app_type integer,
+ cm_solver_pre_app_jacobi_iters integer,
+ temp_init real,
+ temp_final real,
+ temp_mass real,
+ temp_mode integer,
+ temp_rate real,
+ temp_freq integer,
+ pressure real,
+ pressure_mass real,
+ compress integer,
+ pressure_mode integer,
+ remove_center_of_mass integer,
+ debug_level integer,
+ write_freq integer,
+ traj_compress integer,
+ traj_format integer,
+ traj_title text,
+ atom_info integer,
+ atom_velocities integer,
+ atom_forces integer,
+ bond_info integer,
+ angle_info integer,
+ test_forces integer,
+ molecule_analysis integer,
+ freq_molecule_analysis integer,
+ ignore integer,
+ dipole_analysis integer,
+ freq_dipole_analysis integer,
+ diffusion_coefficient integer,
+ freq_diffusion_coefficient integer,
+ restrict_type integer,
+ PRIMARY KEY (id)
+ );
+
+ CREATE TABLE system_properties(
+ id integer,
+ step integer,
+ total_energy real,
+ potential_energy real,
+ kinetic_energy real,
+ temperature real,
+ volume real,
+ pressure real,
+ PRIMARY KEY (id, step)
+ );
+
+ CREATE TABLE potential(
+ id integer,
+ step integer,
+ bond_energy real,
+ atom_energy real,
+ lone_pair_energy real,
+ angle_energy real,
+ coa_energy real,
+ hydrogen_bond_energy real,
+ torsion_energy real,
+ conjugation_energy real,
+ van_der_waals_energy real,
+ coulombic_energy real,
+ polarization_energy real,
+ PRIMARY KEY (id, step)
+ );
+
+ CREATE TABLE trajectory(
+ id integer,
+ step integer,
+ atom_id integer,
+ position_x real,
+ position_y real,
+ position_z real,
+ charge real,
+ PRIMARY KEY (id, step, atom_id)
+ );
+
+ CREATE TABLE performance(
+ id integer,
+ step integer,
+ time_total real,
+ time_nbrs real,
+ time_init real,
+ time_bonded real,
+ time_nonbonded real,
+ time_cm real,
+ time_cm_sort real,
+ cm_solver_iters integer,
+ time_cm_pre_comp real,
+ time_cm_pre_app real,
+ time_cm_solver_spmv real,
+ time_cm_solver_vec_ops real,
+ time_cm_solver_orthog real,
+ time_cm_solver_tri_solve real,
+ PRIMARY KEY (id, step)
+ );
+ """)
+
+ conn.close()
+
+
+if __name__ == '__main__':
+ lib = cdll.LoadLibrary("libspuremd.so.1")
+
+ setup_qmmm = lib.setup_qmmm_
+ setup_qmmm.argtypes = [c_int, POINTER(c_int),
+ POINTER(c_double), POINTER(c_double), POINTER(c_double),
+ c_int, POINTER(c_int), POINTER(c_double), POINTER(c_double),
+ POINTER(c_double), POINTER(c_double), POINTER(c_double),
+ c_char_p, c_char_p]
+ setup_qmmm.restype = c_void_p
+
+ simulate = lib.simulate
+ simulate.argtypes = [c_void_p]
+ simulate.restype = c_int
+
+ cleanup = lib.cleanup
+ cleanup.argtypes = [c_void_p]
+ cleanup.restype = c_int
+
+ reset_qmmm = lib.reset_qmmm_
+ reset_qmmm.argtypes = [c_void_p, c_int, POINTER(c_int),
+ POINTER(c_double), POINTER(c_double), POINTER(c_double),
+ c_int, POINTER(c_int), POINTER(c_double), POINTER(c_double),
+ POINTER(c_double), POINTER(c_double), POINTER(c_double),
+ c_char_p, c_char_p]
+ reset_qmmm.restype = c_int
+
+ CALLBACKFUNC = CFUNCTYPE(None, c_int, POINTER(ReaxAtom),
+ POINTER(SimulationData))
+
+ setup_callback = lib.setup_callback
+ setup_callback.argtypes = [c_void_p, CALLBACKFUNC]
+ setup_callback.restype = c_int
+
+ set_output_enabled = lib.set_output_enabled
+ set_output_enabled.argtypes = [c_void_p, c_int]
+ set_output_enabled.restype = c_int
+
+ set_control_parameter = lib.set_control_parameter
+ set_control_parameter.argtypes = [c_void_p, c_char_p, POINTER(c_char_p)]
+ set_control_parameter.restype = c_int
+
+ get_atom_positions_qmmm = lib.get_atom_positions_qmmm_
+ get_atom_positions_qmmm.argtypes = [c_void_p,
+ POINTER(c_double), POINTER(c_double), POINTER(c_double),
+ POINTER(c_double), POINTER(c_double), POINTER(c_double)]
+ get_atom_positions_qmmm.restype = c_int
+
+ get_atom_forces_qmmm = lib.get_atom_forces_qmmm_
+ get_atom_forces_qmmm.argtypes = [c_void_p,
+ POINTER(c_double), POINTER(c_double), POINTER(c_double),
+ POINTER(c_double), POINTER(c_double), POINTER(c_double)]
+ get_atom_forces_qmmm.restype = c_int
+
+ get_atom_charges_qmmm = lib.get_atom_charges_qmmm_
+ get_atom_charges_qmmm.argtypes = [c_void_p, POINTER(c_double), POINTER(c_double)]
+ get_atom_charges_qmmm.restype = c_int
+
+ def get_simulation_step_results(num_atoms, atoms, data):
+ print("{0:24.15f} {1:24.15f} {2:24.15f}".format(
+ data[0].E_Tot, data[0].E_Kin, data[0].E_Pot))
+
+ # data from Amber
+ sim_box = (c_double * 6)(80.0, 80.0, 80.0, 90.0, 90.0, 90.0)
+ num_qm_atoms = 14
+ num_mm_atoms = 762
+ num_atoms = num_qm_atoms + num_mm_atoms
+
+ df = pd.read_csv('AVE/fort.3', sep='\s+', skiprows=[0,1,2,779,780],
+ names=['Keyword', 'Num', 'Tag', 'P_x', 'P_y', 'P_z', 'Elem',
+ 'Tag2', 'Tag3', 'Q'],
+ dtype={'Keyword': str, 'Num': np.int,
+ 'Tag': str, 'P_x': np.float64, 'P_y': np.float64,
+ 'P_z': np.float64, 'Elem': str,
+ 'Tag2': str, 'Tag3': str, 'Q': np.float64})
+
+ types_s = df['Elem'].to_list()
+ elems = {'C': 0, 'H': 1, 'O': 2}
+ types = [elems[t] for t in types_s]
+ p_x = df['P_x'].to_list()
+ p_y = df['P_y'].to_list()
+ p_z = df['P_z'].to_list()
+ q = df['Q'].to_list()
+
+ qm_types = (c_int * num_qm_atoms)(*types[0:num_qm_atoms])
+ qm_p_x = (c_double * num_qm_atoms)(*p_x[0:num_qm_atoms])
+ qm_p_y = (c_double * num_qm_atoms)(*p_y[0:num_qm_atoms])
+ qm_p_z = (c_double * num_qm_atoms)(*p_z[0:num_qm_atoms])
+ mm_types = (c_int * num_mm_atoms)(*types[num_qm_atoms:])
+ mm_p_x = (c_double * num_mm_atoms)(*p_x[num_qm_atoms:])
+ mm_p_y = (c_double * num_mm_atoms)(*p_y[num_qm_atoms:])
+ mm_p_z = (c_double * num_mm_atoms)(*p_z[num_qm_atoms:])
+ mm_q = (c_double * num_mm_atoms)(*q[num_qm_atoms:])
+
+ handle = setup_qmmm(c_int(num_qm_atoms), qm_types, qm_p_x, qm_p_y, qm_p_z,
+ c_int(num_mm_atoms), mm_types, mm_p_x, mm_p_y, mm_p_z, mm_q, sim_box,
+ b"AVE/fort.4", None)
+
+ ret = setup_callback(handle, CALLBACKFUNC(get_simulation_step_results))
+ if ret != 0:
+ print("[ERROR] setup_callback returned {0}".format(ret))
+
+ ret = set_output_enabled(handle, c_int(1))
+ if ret != 0:
+ print("[ERROR] set_output_enabled returned {0}".format(ret))
+
+ d = {
+ b"ensemble_type": (c_char_p)(b"0"),
+ b"nsteps": (c_char_p)(b"0"),
+ b"dt": (c_char_p)(b"0.25"),
+ b"periodic_boundaries": (c_char_p)(b"1"),
+ b"reposition_atoms": (c_char_p)(b"0"),
+ b"reneighbor": (c_char_p)(b"1"),
+ b"tabulate_long_range": (c_char_p)(b"0"),
+ b"vlist_buffer": (c_char_p)(b"2.5"),
+ b"nbrhood_cutoff": (c_char_p)(b"5.0"),
+ b"thb_cutoff": (c_char_p)(b"0.001"),
+ b"hbond_cutoff": (c_char_p)(b"7.5"),
+ b"bond_graph_cutoff": (c_char_p)(b"0.3"),
+ b"charge_method": (c_char_p)(b"0"),
+ b"cm_q_net": (c_char_p)(b"0.0"),
+ b"cm_solver_type": (c_char_p)(b"2"),
+ b"cm_solver_max_iters": (c_char_p)(b"200"),
+ b"cm_solver_q_err": (c_char_p)(b"1.0e-14"),
+ b"cm_solver_pre_comp_type": (c_char_p)(b"1"),
+ b"write_freq": (c_char_p)(b"1"),
+ }
+ for keyword, values in d.items():
+ ret = set_control_parameter(handle, keyword, values)
+ if ret != 0:
+ print("[ERROR] set_control_parameter returned {0}".format(ret))
+
+ print("{0:24}|{1:24}|{2:24}".format("Total Energy", "Kinetic Energy", "Potential Energy"))
+
+ ret = simulate(handle)
+ if ret != 0:
+ print("[ERROR] simulate returned {0}".format(ret))
+
+ qm_f_x = (c_double * num_qm_atoms)()
+ qm_f_y = (c_double * num_qm_atoms)()
+ qm_f_z = (c_double * num_qm_atoms)()
+ mm_f_x = (c_double * num_mm_atoms)()
+ mm_f_y = (c_double * num_mm_atoms)()
+ mm_f_z = (c_double * num_mm_atoms)()
+ ret = get_atom_forces_qmmm(handle, qm_f_x, qm_f_y, qm_f_z, mm_f_x, mm_f_y, mm_f_z)
+ if ret != 0:
+ print("[ERROR] get_atom_forces_qmmm returned {0}".format(ret))
+
+ qm_q = (c_double * num_qm_atoms)()
+ mm_q = (c_double * num_mm_atoms)()
+ ret = get_atom_charges_qmmm(handle, qm_q, mm_q)
+ if ret != 0:
+ print("[ERROR] get_atom_charges_qmmm returned {0}".format(ret))
+
+ print("\n{0:6}|{1:24}|{2:24}|{3:24}|{4:24}".format("i", "F_x", "F_y", "F_z", "Q"))
+ for i, (f_x, f_y, f_z, q) in enumerate(zip(qm_f_x, qm_f_y, qm_f_z, qm_q)):
+ print("{0:6d} {1:24.12f} {2:24.12f} {3:24.12f} {4:24.12f}".format(i + 1, f_x, f_y, f_z, q))
+ for i, (f_x, f_y, f_z, q) in enumerate(zip(mm_f_x, mm_f_y, mm_f_z, mm_q)):
+ print("{0:6d} {1:24.12f} {2:24.12f} {3:24.12f} {4:24.12f}".format(i + num_qm_atoms + 1, f_x, f_y, f_z, q))
+
+ cleanup(handle)
--
GitLab