diff --git a/sPuReMD/src/allocate.c b/sPuReMD/src/allocate.c
index 15f425652cc274249f53e3ef6b5dbd8270d74fa7..6c07fe4f2e8ae5a25f4fc0108057d74a7b00b98e 100644
--- a/sPuReMD/src/allocate.c
+++ b/sPuReMD/src/allocate.c
@@ -269,7 +269,7 @@ void Reallocate( reax_system *system, control_params *control, static_storage *w
if ( realloc->num_far > 0 && nbr_flag )
{
- Reallocate_Neighbor_List( (*lists) + FAR_NBRS,
+ Reallocate_Neighbor_List( &(*lists)[FAR_NBRS],
system->N, realloc->num_far * SAFE_ZONE );
realloc->num_far = -1;
}
@@ -287,29 +287,31 @@ void Reallocate( reax_system *system, control_params *control, static_storage *w
if ( control->hb_cut > 0.0 && realloc->hbonds > 0 )
{
- Reallocate_HBonds_List(system->N, workspace->num_H, workspace->hbond_index,
- (*lists) + HBONDS );
+ Reallocate_HBonds_List( system->N, workspace->num_H, workspace->hbond_index,
+ &(*lists)[HBONDS] );
realloc->hbonds = -1;
}
num_bonds = est_3body = -1;
if ( realloc->bonds > 0 )
{
- Reallocate_Bonds_List( system->N, (*lists) + BONDS, &num_bonds, &est_3body );
+ Reallocate_Bonds_List( system->N, &(*lists)[BONDS], &num_bonds, &est_3body );
realloc->bonds = -1;
realloc->num_3body = MAX( realloc->num_3body, est_3body );
}
if ( realloc->num_3body > 0 )
{
- Delete_List( TYP_THREE_BODY, (*lists) + THREE_BODIES );
+ Delete_List( TYP_THREE_BODY, &(*lists)[THREE_BODIES] );
if ( num_bonds == -1 )
- num_bonds = ((*lists) + BONDS)->total_intrs;
+ {
+ num_bonds = (*lists)[BONDS].total_intrs;
+ }
realloc->num_3body *= SAFE_ZONE;
Make_List( num_bonds, realloc->num_3body,
- TYP_THREE_BODY, (*lists) + THREE_BODIES );
+ TYP_THREE_BODY, &(*lists)[THREE_BODIES] );
realloc->num_3body = -1;
#if defined(DEBUG_FOCUS)
fprintf( stderr, "reallocating 3 bodies\n" );
diff --git a/sPuReMD/src/analyze.c b/sPuReMD/src/analyze.c
index 675e900718cdf06c193cc1e7a7bfd6f1aad90ff1..1b260eea95d5e16ab443a8ba3d5ccae8a749e22a 100644
--- a/sPuReMD/src/analyze.c
+++ b/sPuReMD/src/analyze.c
@@ -64,8 +64,8 @@ void Copy_Bond_List( reax_system *system, control_params *control,
reax_list **lists )
{
int i, j, top_old;
- reax_list *new_bonds = (*lists) + BONDS;
- reax_list *old_bonds = (*lists) + OLD_BONDS;
+ reax_list *new_bonds = &(*lists)[BONDS];
+ reax_list *old_bonds = &(*lists)[OLD_BONDS];
for ( top_old = 0, i = 0; i < system->N; ++i )
{
@@ -94,8 +94,8 @@ void Copy_Bond_List( reax_system *system, control_params *control,
int Compare_Bond_Lists( int atom, control_params *control, reax_list **lists )
{
int oldp, newp;
- reax_list *new_bonds = (*lists) + BONDS;
- reax_list *old_bonds = (*lists) + OLD_BONDS;
+ reax_list *new_bonds = &(*lists)[BONDS];
+ reax_list *old_bonds = &(*lists)[OLD_BONDS];
/*fprintf( stdout, "\n%d\nold_bonds:", atom );
for( oldp = Start_Index( atom, old_bonds );
@@ -229,8 +229,8 @@ void Analyze_Molecules( reax_system *system, control_params *control,
int *old_mark = workspace->old_mark;
int num_old, num_new;
char s[MAX_MOLECULE_SIZE * 10];
- reax_list *new_bonds = (*lists) + BONDS;
- reax_list *old_bonds = (*lists) + OLD_BONDS;
+ reax_list *new_bonds = &(*lists)[BONDS];
+ reax_list *old_bonds = &(*lists)[OLD_BONDS];
molecule old_molecules[20], new_molecules[20];
fprintf( fout, "molecular analysis @ %d\n", data->step );
@@ -543,8 +543,8 @@ void Analyze_Fragments( reax_system *system, control_params *control,
char fragments[MAX_FRAGMENT_TYPES][MAX_ATOM_TYPES];
int fragment_count[MAX_FRAGMENT_TYPES];
molecule m;
- reax_list *new_bonds = (*lists) + BONDS;
- //reax_list *old_bonds = (*lists) + OLD_BONDS;
+ reax_list *new_bonds = &(*lists)[BONDS];
+// reax_list *old_bonds = &(*lists)[OLD_BONDS];
/* fragment analysis */
fprintf( fout, "step%d fragments\n", data->step );
@@ -609,14 +609,16 @@ void Analyze_Silica( reax_system *system, control_params *control,
int O_SI_O_count, SI_O_SI_count;
int si_coord[10], ox_coord[10];
real O_SI_O, SI_O_SI;
- reax_list *new_bonds = (*lists) + BONDS;
- reax_list *thb_intrs = (*lists) + THREE_BODIES;
+ reax_list *new_bonds = &(*lists)[BONDS];
+ reax_list *thb_intrs = &(*lists)[THREE_BODIES];
Analyze_Fragments( system, control, data, workspace, lists, fout, 0 );
/* analyze atom coordinations */
for ( i = 0; i < 10; ++i )
+ {
si_coord[i] = ox_coord[i] = 0;
+ }
for ( atom = 0; atom < system->N; ++atom )
{
@@ -624,8 +626,12 @@ void Analyze_Silica( reax_system *system, control_params *control,
for ( newp = Start_Index( atom, new_bonds );
newp < End_Index( atom, new_bonds ); ++newp )
+ {
if ( new_bonds->select.bond_list[newp].bo_data.BO >= control->bg_cut )
+ {
++coord;
+ }
+ }
if ( system->atoms[ atom ].type == SI_ATOM )
{
@@ -718,11 +724,12 @@ void Analyze_Silica( reax_system *system, control_params *control,
}
-
int Get_Type_of_Molecule( molecule *m )
{
if ( m->atom_count == 3 && m->mtypes[1] == 2 && m->mtypes[2] == 1 )
+ {
return WATER;
+ }
return UNKNOWN;
}
@@ -982,7 +989,7 @@ void Analysis( reax_system *system, control_params *control,
/****** Electric Dipole Moment ******/
if ( control->dipole_anal && steps % control->freq_dipole_anal == 0 )
Calculate_Dipole_Moment( system, control, data, workspace,
- (*lists) + BONDS, out_control->dpl );
+ &(*lists)[BONDS], out_control->dpl );
/****** Drift ******/
if ( control->diffusion_coef && steps % control->freq_diffusion_coef == 0 )
diff --git a/sPuReMD/src/bond_orders.c b/sPuReMD/src/bond_orders.c
index 52d5ad6162adbbdd2f5bb7e09f6259c5a1a57d51..d8f71ee014c4d6b072476a4e5fe1138fec31e425 100644
--- a/sPuReMD/src/bond_orders.c
+++ b/sPuReMD/src/bond_orders.c
@@ -37,10 +37,12 @@ static inline real Cf45( real p1, real p2 )
void Get_dBO( reax_system *system, reax_list **lists,
int i, int pj, real C, rvec *v )
{
- reax_list *bonds = (*lists) + BONDS;
- reax_list *dBOs = (*lists) + DBO;
+ reax_list *bonds;
+ reax_list *dBOs;
int start_pj, end_pj, k;
+ bonds = &(*lists)[BONDS];
+ dBOs = &(*lists)[DBO];
pj = bonds->select.bond_list[pj].dbond_index;
start_pj = Start_Index(pj, dBOs);
end_pj = End_Index(pj, dBOs);
@@ -56,11 +58,13 @@ void Get_dBO( reax_system *system, reax_list **lists,
void Get_dBOpinpi2( reax_system *system, reax_list **lists,
int i, int pj, real Cpi, real Cpi2, rvec *vpi, rvec *vpi2 )
{
- reax_list *bonds = (*lists) + BONDS;
- reax_list *dBOs = (*lists) + DBO;
+ reax_list *bonds;
+ reax_list *dBOs;
dbond_data *dbo_k;
int start_pj, end_pj, k;
+ bonds = &(*lists)[BONDS];
+ dBOs = &(*lists)[DBO];
pj = bonds->select.bond_list[pj].dbond_index;
start_pj = Start_Index(pj, dBOs);
end_pj = End_Index(pj, dBOs);
@@ -77,10 +81,12 @@ void Get_dBOpinpi2( reax_system *system, reax_list **lists,
void Add_dBO( reax_system *system, reax_list **lists,
int i, int pj, real C, rvec *v )
{
- reax_list *bonds = (*lists) + BONDS;
- reax_list *dBOs = (*lists) + DBO;
+ reax_list *bonds;
+ reax_list *dBOs;
int start_pj, end_pj, k;
+ bonds = &(*lists)[BONDS];
+ dBOs = &(*lists)[DBO];
pj = bonds->select.bond_list[pj].dbond_index;
start_pj = Start_Index(pj, dBOs);
end_pj = End_Index(pj, dBOs);
@@ -98,11 +104,13 @@ void Add_dBO( reax_system *system, reax_list **lists,
void Add_dBOpinpi2( reax_system *system, reax_list **lists,
int i, int pj, real Cpi, real Cpi2, rvec *vpi, rvec *vpi2 )
{
- reax_list *bonds = (*lists) + BONDS;
- reax_list *dBOs = (*lists) + DBO;
+ reax_list *bonds;
+ reax_list *dBOs;
dbond_data *dbo_k;
int start_pj, end_pj, k;
+ bonds = &(*lists)[BONDS];
+ dBOs = &(*lists)[DBO];
pj = bonds->select.bond_list[pj].dbond_index;
start_pj = Start_Index(pj, dBOs);
end_pj = End_Index(pj, dBOs);
@@ -119,10 +127,12 @@ void Add_dBOpinpi2( reax_system *system, reax_list **lists,
void Add_dBO_to_Forces( reax_system *system, reax_list **lists,
int i, int pj, real C )
{
- reax_list *bonds = (*lists) + BONDS;
- reax_list *dBOs = (*lists) + DBO;
+ reax_list *bonds;
+ reax_list *dBOs;
int start_pj, end_pj, k;
+ bonds = &(*lists)[BONDS];
+ dBOs = &(*lists)[DBO];
pj = bonds->select.bond_list[pj].dbond_index;
start_pj = Start_Index(pj, dBOs);
end_pj = End_Index(pj, dBOs);
@@ -138,11 +148,13 @@ void Add_dBO_to_Forces( reax_system *system, reax_list **lists,
void Add_dBOpinpi2_to_Forces( reax_system *system, reax_list **lists,
int i, int pj, real Cpi, real Cpi2 )
{
- reax_list *bonds = (*lists) + BONDS;
- reax_list *dBOs = (*lists) + DBO;
+ reax_list *bonds;
+ reax_list *dBOs;
dbond_data *dbo_k;
int start_pj, end_pj, k;
+ bonds = &(*lists)[BONDS];
+ dBOs = &(*lists)[DBO];
pj = bonds->select.bond_list[pj].dbond_index;
start_pj = Start_Index(pj, dBOs);
end_pj = End_Index(pj, dBOs);
@@ -201,8 +213,8 @@ void Calculate_dBO( int i, int pj, static_storage *workspace, reax_list **lists,
dbond_data *top_dbo;
/* Initializations */
- bonds = (*lists) + BONDS;
- dBOs = (*lists) + DBO;
+ bonds = &(*lists)[BONDS];
+ dBOs = &(*lists)[DBO];
j = bonds->select.bond_list[pj].nbr;
bo_ij = &(bonds->select.bond_list[pj].bo_data);
start_i = Start_Index( i, bonds );
@@ -333,7 +345,7 @@ void Calculate_dBO( int i, int pj, static_storage *workspace, reax_list **lists,
void Add_dBond_to_Forces_NPT( int i, int pj, reax_system *system,
simulation_data *data, static_storage *workspace, reax_list **lists )
{
- reax_list *bonds = (*lists) + BONDS;
+ reax_list *bonds;
bond_data *nbr_j, *nbr_k;
bond_order_data *bo_ij, *bo_ji;
dbond_coefficients coef;
@@ -346,6 +358,7 @@ void Add_dBond_to_Forces_NPT( int i, int pj, reax_system *system,
#endif
/* Initializations */
+ bonds = &(*lists)[BONDS];
nbr_j = &(bonds->select.bond_list[pj]);
j = nbr_j->nbr;
bo_ij = &(nbr_j->bo_data);
@@ -523,7 +536,7 @@ void Add_dBond_to_Forces_NPT( int i, int pj, reax_system *system,
void Add_dBond_to_Forces( int i, int pj, reax_system *system,
simulation_data *data, static_storage *workspace, reax_list **lists )
{
- reax_list *bonds = (*lists) + BONDS;
+ reax_list *bonds;
bond_data *nbr_j, *nbr_k;
bond_order_data *bo_ij, *bo_ji;
dbond_coefficients coef;
@@ -534,6 +547,7 @@ void Add_dBond_to_Forces( int i, int pj, reax_system *system,
#endif
/* Initializations */
+ bonds = &(*lists)[BONDS];
nbr_j = &(bonds->select.bond_list[pj]);
j = nbr_j->nbr;
bo_ij = &(nbr_j->bo_data);
@@ -728,7 +742,7 @@ void Calculate_Bond_Orders( reax_system *system, control_params *control,
p_lp1 = system->reaxprm.gp.l[15];
p_boc1 = system->reaxprm.gp.l[0];
p_boc2 = system->reaxprm.gp.l[1];
- bonds = (*lists) + BONDS;
+ bonds = &(*lists)[BONDS];
#ifdef _OPENMP
#pragma omp parallel default(shared)
@@ -758,8 +772,8 @@ void Calculate_Bond_Orders( reax_system *system, control_params *control,
top_dbo = 0;
top_dDelta = 0;
- dDeltas = (*lists) + DDELTA;
- dBOs = (*lists) + DBO;
+ dDeltas = &(*lists)[DDELTA];
+ dBOs = &(*lists)[DBO];
#endif
/* Calculate Deltaprime, Deltaprime_boc values */
diff --git a/sPuReMD/src/forces.c b/sPuReMD/src/forces.c
index 86b6ebce8310e4a8c044bae6aed709c863e565f8..24e5c2914a3e76b365a7dbf10acab4ee677f05b3 100644
--- a/sPuReMD/src/forces.c
+++ b/sPuReMD/src/forces.c
@@ -150,7 +150,7 @@ void Compute_NonBonded_Forces( reax_system *system, control_params *control,
fprintf( stderr, "qeq - " );
#endif
- if ( control->tabulate == 0)
+ if ( control->tabulate <= 0 )
{
vdW_Coulomb_Energy( system, control, data, workspace, lists, out_control );
}
@@ -179,7 +179,7 @@ void Compute_Total_Force( reax_system *system, control_params *control,
int i;
reax_list *bonds;
- bonds = (*lists) + BONDS;
+ bonds = &(*lists)[BONDS];
#ifdef _OPENMP
#pragma omp parallel default(shared)
@@ -233,8 +233,8 @@ void Validate_Lists( static_storage *workspace, reax_list **lists, int step, int
int i, flag;
reax_list *bonds, *hbonds;
- bonds = *lists + BONDS;
- hbonds = *lists + HBONDS;
+ bonds = &(*lists)[BONDS];
+ hbonds = &(*lists)[HBONDS];
/* far neighbors */
if ( Htop > Hmax * DANGER_ZONE )
@@ -625,23 +625,24 @@ void Init_Forces( reax_system *system, control_params *control,
bond_data *ibond, *jbond;
bond_order_data *bo_ij, *bo_ji;
- far_nbrs = *lists + FAR_NBRS;
- bonds = *lists + BONDS;
- hbonds = *lists + HBONDS;
+ far_nbrs = &(*lists)[FAR_NBRS];
+ bonds = &(*lists)[BONDS];
+ hbonds = &(*lists)[HBONDS];
H = workspace->H;
H_sp = workspace->H_sp;
Htop = 0;
H_sp_top = 0;
num_bonds = 0;
num_hbonds = 0;
- btop_i = btop_j = 0;
+ btop_i = 0;
+ btop_j = 0;
for ( i = 0; i < system->N; ++i )
{
atom_i = &(system->atoms[i]);
type_i = atom_i->type;
- start_i = Start_Index(i, far_nbrs);
- end_i = End_Index(i, far_nbrs);
+ start_i = Start_Index( i, far_nbrs );
+ end_i = End_Index( i, far_nbrs );
H->start[i] = Htop;
H_sp->start[i] = H_sp_top;
btop_i = End_Index( i, bonds );
@@ -661,7 +662,7 @@ void Init_Forces( reax_system *system, control_params *control,
flag = 0;
flag_sp = 0;
- if ((data->step - data->prev_steps) % control->reneighbor == 0)
+ if ( (data->step - data->prev_steps) % control->reneighbor == 0 )
{
if ( nbr_pj->d <= control->r_cut )
{
@@ -678,7 +679,7 @@ void Init_Forces( reax_system *system, control_params *control,
}
}
else if ((nbr_pj->d = Sq_Distance_on_T3(atom_i->x, atom_j->x, &(system->box),
- nbr_pj->dvec)) <= SQR(control->r_cut))
+ nbr_pj->dvec)) <= SQR(control->r_cut))
{
if ( nbr_pj->d <= SQR(control->r_sp_cut))
{
@@ -738,7 +739,7 @@ void Init_Forces( reax_system *system, control_params *control,
{
r2 = SQR( r_ij );
- if ( sbp_i->r_s > 0.0 && sbp_j->r_s > 0.0)
+ if ( sbp_i->r_s > 0.0 && sbp_j->r_s > 0.0 )
{
C12 = twbp->p_bo1 * POW( r_ij / twbp->r_s, twbp->p_bo2 );
BO_s = (1.0 + control->bo_cut) * EXP( C12 );
@@ -749,7 +750,7 @@ void Init_Forces( reax_system *system, control_params *control,
C12 = 0.0;
}
- if ( sbp_i->r_pi > 0.0 && sbp_j->r_pi > 0.0)
+ if ( sbp_i->r_pi > 0.0 && sbp_j->r_pi > 0.0 )
{
C34 = twbp->p_bo3 * POW( r_ij / twbp->r_p, twbp->p_bo4 );
BO_pi = EXP( C34 );
@@ -760,7 +761,7 @@ void Init_Forces( reax_system *system, control_params *control,
C34 = 0.0;
}
- if ( sbp_i->r_pi_pi > 0.0 && sbp_j->r_pi_pi > 0.0)
+ if ( sbp_i->r_pi_pi > 0.0 && sbp_j->r_pi_pi > 0.0 )
{
C56 = twbp->p_bo5 * POW( r_ij / twbp->r_pp, twbp->p_bo6 );
BO_pi2 = EXP( C56 );
@@ -937,31 +938,34 @@ void Init_Forces_Tab( reax_system *system, control_params *control,
bond_data *ibond, *jbond;
bond_order_data *bo_ij, *bo_ji;
- far_nbrs = *lists + FAR_NBRS;
- bonds = *lists + BONDS;
- hbonds = *lists + HBONDS;
-
+ far_nbrs = &(*lists)[FAR_NBRS];
+ bonds = &(*lists)[BONDS];
+ hbonds = &(*lists)[HBONDS];
H = workspace->H;
H_sp = workspace->H_sp;
Htop = 0;
H_sp_top = 0;
num_bonds = 0;
num_hbonds = 0;
- btop_i = btop_j = 0;
+ btop_i = 0;
+ btop_j = 0;
for ( i = 0; i < system->N; ++i )
{
atom_i = &(system->atoms[i]);
type_i = atom_i->type;
- start_i = Start_Index(i, far_nbrs);
- end_i = End_Index(i, far_nbrs);
+ start_i = Start_Index( i, far_nbrs );
+ end_i = End_Index( i, far_nbrs );
H->start[i] = Htop;
H_sp->start[i] = H_sp_top;
btop_i = End_Index( i, bonds );
sbp_i = &(system->reaxprm.sbp[type_i]);
ihb = ihb_top = -1;
+
if ( control->hb_cut > 0 && (ihb = sbp_i->p_hbond) == 1 )
+ {
ihb_top = End_Index( workspace->hbond_index[i], hbonds );
+ }
for ( pj = start_i; pj < end_i; ++pj )
{
@@ -971,7 +975,7 @@ void Init_Forces_Tab( reax_system *system, control_params *control,
flag = 0;
flag_sp = 0;
- if ((data->step - data->prev_steps) % control->reneighbor == 0)
+ if ( (data->step - data->prev_steps) % control->reneighbor == 0 )
{
if (nbr_pj->d <= control->r_cut)
{
@@ -988,7 +992,7 @@ void Init_Forces_Tab( reax_system *system, control_params *control,
}
}
else if ((nbr_pj->d = Sq_Distance_on_T3(atom_i->x, atom_j->x, &(system->box),
- nbr_pj->dvec)) <= SQR(control->r_cut))
+ nbr_pj->dvec)) <= SQR(control->r_cut))
{
if ( nbr_pj->d <= SQR(control->r_sp_cut))
{
@@ -1006,8 +1010,8 @@ void Init_Forces_Tab( reax_system *system, control_params *control,
twbp = &(system->reaxprm.tbp[type_i][type_j]);
H->j[Htop] = j;
- H->val[Htop] = Init_Charge_Matrix_Entry_Tab( system, control,
- workspace->LR, i, j, r_ij, OFF_DIAGONAL );
+ H->val[Htop] = Init_Charge_Matrix_Entry( system, control, i, j,
+ r_ij, OFF_DIAGONAL );
++Htop;
/* H_sp matrix entry */
@@ -1048,7 +1052,7 @@ void Init_Forces_Tab( reax_system *system, control_params *control,
{
r2 = SQR( r_ij );
- if ( sbp_i->r_s > 0.0 && sbp_j->r_s > 0.0)
+ if ( sbp_i->r_s > 0.0 && sbp_j->r_s > 0.0 )
{
C12 = twbp->p_bo1 * POW( r_ij / twbp->r_s, twbp->p_bo2 );
BO_s = (1.0 + control->bo_cut) * EXP( C12 );
@@ -1059,7 +1063,7 @@ void Init_Forces_Tab( reax_system *system, control_params *control,
C12 = 0.0;
}
- if ( sbp_i->r_pi > 0.0 && sbp_j->r_pi > 0.0)
+ if ( sbp_i->r_pi > 0.0 && sbp_j->r_pi > 0.0 )
{
C34 = twbp->p_bo3 * POW( r_ij / twbp->r_p, twbp->p_bo4 );
BO_pi = EXP( C34 );
@@ -1070,7 +1074,7 @@ void Init_Forces_Tab( reax_system *system, control_params *control,
C34 = 0.0;
}
- if ( sbp_i->r_pi_pi > 0.0 && sbp_j->r_pi_pi > 0.0)
+ if ( sbp_i->r_pi_pi > 0.0 && sbp_j->r_pi_pi > 0.0 )
{
C56 = twbp->p_bo5 * POW( r_ij / twbp->r_pp, twbp->p_bo6 );
BO_pi2 = EXP( C56 );
@@ -1097,7 +1101,6 @@ void Init_Forces_Tab( reax_system *system, control_params *control,
ibond->d = r_ij;
jbond->d = r_ij;
rvec_Copy( ibond->dvec, nbr_pj->dvec );
- //fprintf (stderr, " %f - %f - %f \n", nbr_pj->dvec[0], nbr_pj->dvec[1], nbr_pj->dvec[2]);
rvec_Scale( jbond->dvec, -1, nbr_pj->dvec );
ivec_Copy( ibond->rel_box, nbr_pj->rel_box );
ivec_Scale( jbond->rel_box, -1, nbr_pj->rel_box );
@@ -1162,8 +1165,8 @@ void Init_Forces_Tab( reax_system *system, control_params *control,
/* diagonal entry */
H->j[Htop] = i;
- H->val[Htop] = Init_Charge_Matrix_Entry_Tab( system, control, workspace->LR,
- i, j, r_ij, DIAGONAL );
+ H->val[Htop] = Init_Charge_Matrix_Entry( system, control, i, i,
+ r_ij, DIAGONAL );
++Htop;
H_sp->j[H_sp_top] = i;
@@ -1177,12 +1180,16 @@ void Init_Forces_Tab( reax_system *system, control_params *control,
}
}
+ Init_Charge_Matrix_Remaining_Entries( system, control, far_nbrs,
+ H, H_sp, &Htop, &H_sp_top );
+
// mark the end of j list
- H->start[i] = Htop;
- H_sp->start[i] = H_sp_top;
+ H->start[system->N_cm] = Htop;
+ H_sp->start[system->N_cm] = H_sp_top;
+
/* validate lists - decide if reallocation is required! */
- Validate_Lists( workspace, lists, data->step, system->N,
- H->m, Htop, num_bonds, num_hbonds );
+ Validate_Lists( workspace, lists,
+ data->step, system->N, H->m, Htop, num_bonds, num_hbonds );
#if defined(DEBUG_FOCUS)
fprintf( stderr, "step%d: Htop = %d, num_bonds = %d, num_hbonds = %d\n",
@@ -1304,7 +1311,7 @@ void Compute_Forces( reax_system *system, control_params *control,
real t_start, t_elapsed;
t_start = Get_Time( );
- if ( !control->tabulate )
+ if ( control->tabulate <= 0 )
{
Init_Forces( system, control, data, workspace, lists, out_control );
}
diff --git a/sPuReMD/src/four_body_interactions.c b/sPuReMD/src/four_body_interactions.c
index 26a15b7499a0da34f01651e2b7b36822008aa477..22d104763a7aaaca6c97dfed1619ff96e11672a3 100644
--- a/sPuReMD/src/four_body_interactions.c
+++ b/sPuReMD/src/four_body_interactions.c
@@ -164,8 +164,8 @@ void Four_Body_Interactions( reax_system *system, control_params *control,
p_tor3 = system->reaxprm.gp.l[24];
p_tor4 = system->reaxprm.gp.l[25];
p_cot2 = system->reaxprm.gp.l[27];
- bonds = (*lists) + BONDS;
- thb_intrs = (*lists) + THREE_BODIES;
+ bonds = &(*lists)[BONDS];
+ thb_intrs = &(*lists)[THREE_BODIES];
e_tor_total = 0.0;
e_con_total = 0.0;
diff --git a/sPuReMD/src/init_md.c b/sPuReMD/src/init_md.c
index 01b0fc5593bbdb562b374b86227764edbfbc553e..4e8ca7d8885e3335aded9c8ba8f66349757ecace 100644
--- a/sPuReMD/src/init_md.c
+++ b/sPuReMD/src/init_md.c
@@ -590,7 +590,7 @@ void Init_Lists( reax_system *system, control_params *control,
num_nbrs = Estimate_NumNeighbors( system, control, workspace, lists );
- Make_List( system->N, num_nbrs, TYP_FAR_NEIGHBOR, (*lists) + FAR_NBRS );
+ Make_List( system->N, num_nbrs, TYP_FAR_NEIGHBOR, &(*lists)[FAR_NBRS] );
#if defined(DEBUG_FOCUS)
fprintf( stderr, "memory allocated: far_nbrs = %ldMB\n",
@@ -669,7 +669,7 @@ void Init_Lists( reax_system *system, control_params *control,
else
{
Allocate_HBond_List( system->N, workspace->num_H, workspace->hbond_index,
- hb_top, (*lists) + HBONDS );
+ hb_top, &(*lists)[HBONDS] );
}
#if defined(DEBUG_FOCUS)
@@ -681,7 +681,7 @@ void Init_Lists( reax_system *system, control_params *control,
}
/* bonds list */
- Allocate_Bond_List( system->N, bond_top, (*lists) + BONDS );
+ Allocate_Bond_List( system->N, bond_top, &(*lists)[BONDS] );
num_bonds = bond_top[system->N - 1];
#if defined(DEBUG_FOCUS)
@@ -691,7 +691,7 @@ void Init_Lists( reax_system *system, control_params *control,
#endif
/* 3bodies list */
- Make_List( num_bonds, num_3body, TYP_THREE_BODY, (*lists) + THREE_BODIES );
+ Make_List( num_bonds, num_3body, TYP_THREE_BODY, &(*lists)[THREE_BODIES] );
#if defined(DEBUG_FOCUS)
fprintf( stderr, "estimated storage - num_3body: %d\n", num_3body );
@@ -700,9 +700,9 @@ void Init_Lists( reax_system *system, control_params *control,
#endif
#ifdef TEST_FORCES
- Make_List( system->N, num_bonds * 8, TYP_DDELTA, (*lists) + DDELTA );
+ Make_List( system->N, num_bonds * 8, TYP_DDELTA, &(*lists)[DDELTA] );
- Make_List( num_bonds, num_bonds * MAX_BONDS * 3, TYP_DBO, (*lists) + DBO );
+ Make_List( num_bonds, num_bonds * MAX_BONDS * 3, TYP_DBO, &(*lists)[DBO] );
#endif
sfree( hb_top, "Init_Lists::hb_top" );
@@ -1217,17 +1217,17 @@ void Finalize_Workspace( reax_system *system, control_params *control,
void Finalize_Lists( control_params *control, reax_list **lists )
{
- Delete_List( TYP_FAR_NEIGHBOR, (*lists) + FAR_NBRS );
+ Delete_List( TYP_FAR_NEIGHBOR, &(*lists)[FAR_NBRS] );
if ( control->hb_cut > 0.0 )
{
- Delete_List( TYP_HBOND, (*lists) + HBONDS );
+ Delete_List( TYP_HBOND, &(*lists)[HBONDS] );
}
- Delete_List( TYP_BOND, (*lists) + BONDS );
- Delete_List( TYP_THREE_BODY, (*lists) + THREE_BODIES );
+ Delete_List( TYP_BOND, &(*lists)[BONDS] );
+ Delete_List( TYP_THREE_BODY, &(*lists)[THREE_BODIES] );
#ifdef TEST_FORCES
- Delete_List( TYP_DDELTA, (*lists) + DDELTA );
- Delete_List( TYP_DBO, (*lists) + DBO );
+ Delete_List( TYP_DDELTA, &(*lists)[DDELTA] );
+ Delete_List( TYP_DBO, &(*lists)[DBO] );
#endif
}
diff --git a/sPuReMD/src/lin_alg.c b/sPuReMD/src/lin_alg.c
index 5ea019d003c2e42e8f2b20dc16cd8e429e633e3c..4f69152b3d1dc09484a6f3feaf04993f54b4bf88 100644
--- a/sPuReMD/src/lin_alg.c
+++ b/sPuReMD/src/lin_alg.c
@@ -3317,10 +3317,10 @@ int GMRES( const static_storage * const workspace, const control_params * const
if ( g_itr >= control->cm_solver_max_iters )
{
fprintf( stderr, "[WARNING] GMRES convergence failed (%d outer iters)\n", g_itr );
- return g_itr * control->cm_solver_restart + g_j;
+ return g_itr * (control->cm_solver_restart + 1) + g_j + 1;
}
- return g_itr * control->cm_solver_restart + g_j + 1;
+ return g_itr * (control->cm_solver_restart + 1) + g_j + 1;
}
diff --git a/sPuReMD/src/neighbors.c b/sPuReMD/src/neighbors.c
index c0351c857bc4709d729311922f3ecbef8192675a..b63a7fca09a8249c2ab2cab71ffc7741f9df9200 100644
--- a/sPuReMD/src/neighbors.c
+++ b/sPuReMD/src/neighbors.c
@@ -71,9 +71,11 @@ void Generate_Neighbor_Lists( reax_system *system, control_params *control,
t_start = Get_Time( );
// fprintf( stderr, "\n\tentered nbrs - " );
g = &( system->g );
- far_nbrs = (*lists) + FAR_NBRS;
+ far_nbrs = &(*lists)[FAR_NBRS];
+
Bin_Atoms( system, workspace );
// fprintf( stderr, "atoms sorted - " );
+
num_far = 0;
/* first pick up a cell in the grid */
@@ -365,25 +367,29 @@ static inline int is_Near_Neighbor( reax_list *near_nbrs, int atom1, int atom2 )
return TRUE;
}
+
void Generate_Neighbor_Lists( reax_system *system, control_params *control,
simulation_data *data, static_storage *workspace,
reax_list **lists, output_controls *out_control )
{
- int i, j, k;
- int x, y, z;
- int *nbr_atoms;
- int atom1, atom2, max;
- int num_far;
- int c, count;
- int grid_top;
- grid *g = &( system->g );
- reax_list *far_nbrs = (*lists) + FAR_NBRS;
+ int i, j, k;
+ int x, y, z;
+ int *nbr_atoms;
+ int atom1, atom2, max;
+ int num_far;
+ int c, count;
+ int grid_top;
+ grid *g;
+ reax_list *far_nbrs;
//int hb_type1, hb_type2;
- //reax_list *hbonds = (*lists) + HBOND;
+ //reax_list *hbonds = &(*lists)[HBOND];
//int top_hbond1, top_hbond2;
get_far_neighbors_function Get_Far_Neighbors;
far_neighbor_data new_nbrs[125];
+ g = &( system->g );
+ far_nbrs = &(*lists)[FAR_NBRS];
+
// fprintf( stderr, "\n\tentered nbrs - " );
if ( control->ensemble == iNPT || control->ensemble == sNPT ||
control->ensemble == NPT )
@@ -632,7 +638,7 @@ void Generate_Neighbor_Lists( reax_system *system, control_params *control,
far_neighbor_data new_nbrs[125];
g = &( system->g );
- far_nbrs = (*lists) + FAR_NBRS;
+ far_nbrs = &(*lists)[FAR_NBRS];
// fprintf( stderr, "\n\tentered nbrs - " );
if ( control->ensemble == iNPT ||
diff --git a/sPuReMD/src/print_utils.c b/sPuReMD/src/print_utils.c
index 30c071251bdfcfb9a55f147dbeab16b50994132b..8c3ea3095354baddad4c763e69357de35132c5a6 100644
--- a/sPuReMD/src/print_utils.c
+++ b/sPuReMD/src/print_utils.c
@@ -42,8 +42,8 @@ void Print_Bond_Orders( reax_system *system, control_params *control,
{
int i, pj, pk;
bond_order_data *bo_ij;
- reax_list *bonds = (*lists) + BONDS;
- reax_list *dBOs = (*lists) + DBO;
+ reax_list *bonds = &(*lists)[BONDS];
+ reax_list *dBOs = &(*lists)[DBO];
dbond_data *dbo_k;
/* bond orders */
@@ -473,7 +473,7 @@ void Print_Far_Neighbors( reax_system *system, control_params *control,
int i, j, id_i, id_j;
char fname[MAX_STR];
FILE *fout;
- reax_list *far_nbrs = &((*lists)[FAR_NBRS]);
+ reax_list *far_nbrs = &(*lists)[FAR_NBRS];
snprintf( fname, MAX_STR, "%.*s.far_nbrs", MAX_STR - 10, control->sim_name );
fout = fopen( fname, "w" );
@@ -518,7 +518,7 @@ void Print_Far_Neighbors2( reax_system *system, control_params *control,
int i, j, id_i, id_j;
char fname[MAX_STR];
FILE *fout;
- reax_list *far_nbrs = &((*lists)[FAR_NBRS]);
+ reax_list *far_nbrs = &(*lists)[FAR_NBRS];
snprintf( fname, MAX_STR, "%.*s.far_nbrs_lgj", MAX_STR - 14, control->sim_name );
fout = fopen( fname, "w" );
@@ -720,7 +720,7 @@ void Output_Results( reax_system *system, control_params *control,
out_control->append_traj_frame( system, control, data,
workspace, lists, out_control );
- //Write_PDB( system, *lists+BONDS, data, control, workspace, out_control );
+ //Write_PDB( system, &(*lists)[BONDS], data, control, workspace, out_control );
//t_elapsed = Get_Timing_Info( t_start );
//fprintf(stdout, "append_frame took %.6f seconds\n", t_elapsed );
}
diff --git a/sPuReMD/src/reset_utils.c b/sPuReMD/src/reset_utils.c
index 0aed2a953efe94d904eccc3d17c8725ee214f96e..3edfe9072c0b635f7d3cfb97e555dd37d19d4a37 100644
--- a/sPuReMD/src/reset_utils.c
+++ b/sPuReMD/src/reset_utils.c
@@ -118,8 +118,11 @@ void Reset_Neighbor_Lists( reax_system *system, control_params *control,
static_storage *workspace, reax_list **lists )
{
int i, tmp;
- reax_list *bonds = (*lists) + BONDS;
- reax_list *hbonds = (*lists) + HBONDS;
+ reax_list *bonds;
+ reax_list *hbonds;
+
+ bonds = &(*lists)[BONDS];
+ hbonds = &(*lists)[HBONDS];
for ( i = 0; i < system->N; ++i )
{
diff --git a/sPuReMD/src/single_body_interactions.c b/sPuReMD/src/single_body_interactions.c
index b6edce1db6b7becfd311334fec83568625276261..4a9ff3102d5d86e33db27d71c19ed101a465c54f 100644
--- a/sPuReMD/src/single_body_interactions.c
+++ b/sPuReMD/src/single_body_interactions.c
@@ -20,6 +20,7 @@
----------------------------------------------------------------------*/
#include "single_body_interactions.h"
+
#include "bond_orders.h"
#include "list.h"
#include "lookup.h"
@@ -45,8 +46,9 @@ void LonePair_OverUnder_Coordination_Energy( reax_system *system, control_params
two_body_parameters *twbp;
bond_data *pbond;
bond_order_data *bo_ij;
- reax_list *bonds = (*lists) + BONDS;
+ reax_list *bonds;
+ bonds = &(*lists)[BONDS];
/* Initialize parameters */
p_lp1 = system->reaxprm.gp.l[15];
p_lp3 = system->reaxprm.gp.l[5];
diff --git a/sPuReMD/src/three_body_interactions.c b/sPuReMD/src/three_body_interactions.c
index 159e0cc17f354c94172c01729a17f473f941b51b..485be5fdf934f5e0d113a35569a44ed44a5fb14f 100644
--- a/sPuReMD/src/three_body_interactions.c
+++ b/sPuReMD/src/three_body_interactions.c
@@ -93,9 +93,9 @@ void Three_Body_Interactions( reax_system *system, control_params *control,
real e_ang_total, e_pen_total, e_coa_total;
total_bo = workspace->total_bond_order;
- bonds = (*lists) + BONDS;
+ bonds = &(*lists)[BONDS];
bond_list = bonds->select.bond_list;
- thb_intrs = (*lists) + THREE_BODIES;
+ thb_intrs = &(*lists)[THREE_BODIES];
thb_list = thb_intrs->select.three_body_list;
/* global parameters used in these calculations */
p_pen2 = system->reaxprm.gp.l[19];
@@ -705,9 +705,9 @@ void Hydrogen_Bonds( reax_system *system, control_params *control,
#ifdef TEST_FORCES
num_hb_intrs = 0;
#endif
- bonds = (*lists) + BONDS;
+ bonds = &(*lists)[BONDS];
bond_list = bonds->select.bond_list;
- hbonds = (*lists) + HBONDS;
+ hbonds = &(*lists)[HBONDS];
hbond_list = hbonds->select.hbond_list;
/* loops below discover the Hydrogen bonds between i-j-k triplets.
diff --git a/sPuReMD/src/traj.c b/sPuReMD/src/traj.c
index 7972d0868c47fc992f909d9be457c3d454743e44..03cc5d74e1a644fdc10d84e8c6b57d9331c1a644 100644
--- a/sPuReMD/src/traj.c
+++ b/sPuReMD/src/traj.c
@@ -178,8 +178,8 @@ int Append_Custom_Frame( reax_system *system, control_params *control,
int frame_globals_len, num_bonds, num_thb_intrs;
real P;
char buffer[SIZE];
- reax_list *bonds = (*lists) + BONDS;
- reax_list *thb_intrs = (*lists) + THREE_BODIES;
+ reax_list *bonds = &(*lists)[BONDS];
+ reax_list *thb_intrs = &(*lists)[THREE_BODIES];
bond_data *bo_ij;
/* IMPORTANT: This whole part will go to init_trj after finalized! */
diff --git a/sPuReMD/src/two_body_interactions.c b/sPuReMD/src/two_body_interactions.c
index 664d7c0f3088d1b11beceb6e6a2647c49fb46cb6..d83db947b82936725eb37ad9ac9deb83029abbed 100644
--- a/sPuReMD/src/two_body_interactions.c
+++ b/sPuReMD/src/two_body_interactions.c
@@ -35,7 +35,7 @@ void Bond_Energy( reax_system *system, control_params *control,
real gp3, gp4, gp7, gp10, gp37, ebond_total;
reax_list *bonds;
- bonds = (*lists) + BONDS;
+ bonds = &(*lists)[BONDS];
gp3 = system->reaxprm.gp.l[3];
gp4 = system->reaxprm.gp.l[4];
gp7 = system->reaxprm.gp.l[7];
@@ -177,7 +177,7 @@ void vdW_Coulomb_Energy( reax_system *system, control_params *control,
p_vdW1 = system->reaxprm.gp.l[28];
p_vdW1i = 1.0 / p_vdW1;
- far_nbrs = (*lists) + FAR_NBRS;
+ far_nbrs = &(*lists)[FAR_NBRS];
e_vdW_total = 0.0;
e_ele_total = 0.0;
@@ -527,7 +527,7 @@ void Tabulated_vdW_Coulomb_Energy( reax_system *system, control_params *control,
reax_list *far_nbrs;
real e_vdW_total, e_ele_total;
- far_nbrs = (*lists) + FAR_NBRS;
+ far_nbrs = &(*lists)[FAR_NBRS];
steps = data->step - data->prev_steps;
update_freq = out_control->energy_update_freq;
update_energies = update_freq > 0 && steps % update_freq == 0;