diff --git a/environ/param.gpu.water b/environ/param.gpu.water
index bed2200680dffb9cf0cf81d7a422a8b0dbd8fd36..3c38e57dd9a357e3d293dcedc2e0becf90f15e72 100644
--- a/environ/param.gpu.water
+++ b/environ/param.gpu.water
@@ -17,10 +17,12 @@ hbond_cutoff 7.50 ! cutoff distance for hydrogen b
charge_method 0 ! charge method: 0 = QEq, 1 = EEM, 2 = ACKS2
cm_q_net 0.0 ! net system charge
cm_solver_type 0 ! iterative linear solver for charge method: 0 = GMRES, 1 = GMRES_H, 2 = CG, 3 = SDM
+cm_solver_max_iters 20 ! max solver iterations
+cm_solver_restart 100 ! inner iterations of GMRES before restarting
cm_solver_q_err 1e-6 ! 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_solver_pre_comp_type 1 ! method used to compute preconditioner, if applicable
-cm_solver_pre_comp_refactor 100 ! number of steps before recomputing preconditioner
+cm_solver_pre_comp_refactor 1000 ! number of steps before recomputing preconditioner
cm_solver_pre_comp_droptol 0.0 ! threshold tolerance for dropping values in preconditioner computation, if applicable
cm_solver_pre_comp_sweeps 3 ! number of sweeps used to compute preconditioner (ILU_PAR)
cm_solver_pre_app_type 1 ! method used to apply preconditioner
diff --git a/sPuReMD/src/charges.c b/sPuReMD/src/charges.c
index ecd2827e18d42585568f723ad5371811dcbcfb3d..c2242449c0df822996beaf6d5cb7dbe34481540e 100644
--- a/sPuReMD/src/charges.c
+++ b/sPuReMD/src/charges.c
@@ -574,7 +574,7 @@ static real SuperLU_Factorize( const sparse_matrix * const A,
/* Diagonal (Jacobi) preconditioner computation */
-static real diag_pre_comp( const reax_system * const system, real * const Hdia_inv )
+static real diag_pre_comp( const sparse_matrix * const H, real * const Hdia_inv )
{
unsigned int i;
real start;
@@ -583,16 +583,9 @@ static real diag_pre_comp( const reax_system * const system, real * const Hdia_i
#pragma omp parallel for schedule(static) \
default(none) private(i)
- for ( i = 0; i < system->N; ++i )
- {
- Hdia_inv[i] = 1.0 / system->reaxprm.sbp[system->atoms[i].type].eta;
- }
-
- #pragma omp parallel for schedule(static) \
- default(none) private(i)
- for ( i = system->N; i < system->N_cm; ++i )
+ for ( i = 0; i < H->n; ++i )
{
- Hdia_inv[i] = 1.0;
+ Hdia_inv[i] = 1.0 / H->val[H->start[i + 1] - 1];
}
return Get_Timing_Info( start );
@@ -1485,7 +1478,7 @@ static void Compute_Preconditioner_QEq( const reax_system * const system,
{
case DIAG_PC:
data->timing.cm_solver_pre_comp +=
- diag_pre_comp( system, workspace->Hdia_inv );
+ diag_pre_comp( Hptr, workspace->Hdia_inv );
break;
case ICHOLT_PC:
@@ -1578,7 +1571,7 @@ static void Compute_Preconditioner_EEM( const reax_system * const system,
{
case DIAG_PC:
data->timing.cm_solver_pre_comp +=
- diag_pre_comp( system, workspace->Hdia_inv );
+ diag_pre_comp( Hptr, workspace->Hdia_inv );
break;
case ICHOLT_PC:
@@ -1819,7 +1812,7 @@ static void Compute_Preconditioner_ACKS2( const reax_system * const system,
{
case DIAG_PC:
data->timing.cm_solver_pre_comp +=
- diag_pre_comp( system, workspace->Hdia_inv );
+ diag_pre_comp( Hptr, workspace->Hdia_inv );
break;
case ICHOLT_PC:
@@ -2468,9 +2461,9 @@ static void QEq( reax_system * const system, control_params * const control,
break;
case CG_S:
- iters = CG( workspace, workspace->H, workspace->b_s, control->cm_solver_q_err,
+ iters = CG( workspace, control, workspace->H, workspace->b_s, control->cm_solver_q_err,
workspace->s[0], out_control->log ) + 1;
- iters += CG( workspace, workspace->H, workspace->b_t, control->cm_solver_q_err,
+ iters += CG( workspace, control, workspace->H, workspace->b_t, control->cm_solver_q_err,
workspace->t[0], out_control->log ) + 1;
// iters = CG( workspace, workspace->H, workspace->b_s, control->cm_solver_q_err,
// workspace->L, workspace->U, workspace->s[0], control->cm_solver_pre_app_type,
@@ -2481,9 +2474,9 @@ static void QEq( reax_system * const system, control_params * const control,
break;
case SDM_S:
- iters = SDM( workspace, workspace->H, workspace->b_s, control->cm_solver_q_err,
+ iters = SDM( workspace, control, workspace->H, workspace->b_s, control->cm_solver_q_err,
workspace->s[0], out_control->log ) + 1;
- iters += SDM( workspace, workspace->H, workspace->b_t, control->cm_solver_q_err,
+ iters += SDM( workspace,control, workspace->H, workspace->b_t, control->cm_solver_q_err,
workspace->t[0], out_control->log ) + 1;
// iters = SDM( workspace, workspace->H, workspace->b_s, control->cm_solver_q_err,
// workspace->L, workspace->U, workspace->s[0], control->cm_solver_pre_app_type,
@@ -2561,12 +2554,12 @@ static void EEM( reax_system * const system, control_params * const control,
break;
case CG_S:
- iters = CG( workspace, workspace->H, workspace->b_s, control->cm_solver_q_err,
+ iters = CG( workspace, control, workspace->H, workspace->b_s, control->cm_solver_q_err,
workspace->s[0], out_control->log ) + 1;
break;
case SDM_S:
- iters = SDM( workspace, workspace->H, workspace->b_s, control->cm_solver_q_err,
+ iters = SDM( workspace, control, workspace->H, workspace->b_s, control->cm_solver_q_err,
workspace->s[0], out_control->log ) + 1;
break;
@@ -2640,12 +2633,12 @@ static void ACKS2( reax_system * const system, control_params * const control,
break;
case CG_S:
- iters = CG( workspace, workspace->H, workspace->b_s, control->cm_solver_q_err,
+ iters = CG( workspace, control, workspace->H, workspace->b_s, control->cm_solver_q_err,
workspace->s[0], out_control->log ) + 1;
break;
case SDM_S:
- iters = SDM( workspace, workspace->H, workspace->b_s, control->cm_solver_q_err,
+ iters = SDM( workspace, control, workspace->H, workspace->b_s, control->cm_solver_q_err,
workspace->s[0], out_control->log ) + 1;
break;
diff --git a/sPuReMD/src/control.c b/sPuReMD/src/control.c
index 7f0134f3c69b25f13a0f12b7845ce0afe432dad0..e906fc85db5cda6898af5400242be5c23d0415b5 100644
--- a/sPuReMD/src/control.c
+++ b/sPuReMD/src/control.c
@@ -73,6 +73,8 @@ char Read_Control_File( FILE* fp, reax_system *system, control_params* control,
control->charge_method = QEQ_CM;
control->cm_q_net = 0.0;
control->cm_solver_type = GMRES_S;
+ control->cm_solver_max_iters = 100;
+ control->cm_solver_restart = 50;
control->cm_solver_q_err = 0.000001;
control->cm_domain_sparsify_enabled = FALSE;
control->cm_domain_sparsity = 1.0;
@@ -257,6 +259,16 @@ char Read_Control_File( FILE* fp, reax_system *system, control_params* control,
ival = atoi( tmp[1] );
control->cm_solver_type = ival;
}
+ else if ( strcmp(tmp[0], "cm_solver_max_iters") == 0 )
+ {
+ ival = atoi( tmp[1] );
+ control->cm_solver_max_iters = ival;
+ }
+ else if ( strcmp(tmp[0], "cm_solver_restart") == 0 )
+ {
+ ival = atoi( tmp[1] );
+ control->cm_solver_restart = ival;
+ }
else if ( strcmp(tmp[0], "cm_solver_q_err") == 0 )
{
val = atof( tmp[1] );
diff --git a/sPuReMD/src/grid.c b/sPuReMD/src/grid.c
index c45c1a2214d7589b3e4c1a647ebde1b9d3c7e56d..c422a48ae746f336cb0c3398c3990678d93ec35e 100644
--- a/sPuReMD/src/grid.c
+++ b/sPuReMD/src/grid.c
@@ -494,13 +494,12 @@ static inline void reax_atom_Copy( reax_atom *dest, reax_atom *src )
void Copy_Storage( reax_system *system, static_storage *workspace,
- int top, int old_id, int old_type,
- int *num_H, real **v, real **s, real **t,
- int *orig_id, rvec *f_old )
+ control_params *control, int top, int old_id, int old_type, int *num_H,
+ real **v, real **s, real **t, int *orig_id, rvec *f_old )
{
int i;
- for ( i = 0; i < RESTART + 1; ++i )
+ for ( i = 0; i < control->cm_solver_restart + 1; ++i )
{
v[i][top] = workspace->v[i][old_id];
}
@@ -529,11 +528,11 @@ void Copy_Storage( reax_system *system, static_storage *workspace,
}
-void Free_Storage( static_storage *workspace )
+void Free_Storage( static_storage *workspace, control_params * control )
{
int i;
- for ( i = 0; i < RESTART + 1; ++i )
+ for ( i = 0; i < control->cm_solver_restart + 1; ++i )
{
free( workspace->v[i] );
}
@@ -566,7 +565,8 @@ void Assign_New_Storage( static_storage *workspace,
}
-void Cluster_Atoms( reax_system *system, static_storage *workspace )
+void Cluster_Atoms( reax_system *system, static_storage *workspace,
+ control_params *control )
{
int i, j, k, l, top, old_id, num_H;
reax_atom *old_atom;
@@ -591,8 +591,8 @@ void Cluster_Atoms( reax_system *system, static_storage *workspace )
t[i] = (real *) calloc( system->N, sizeof( real ) );
}
- v = (real**) calloc( RESTART + 1, sizeof( real* ) );
- for ( i = 0; i < RESTART + 1; ++i )
+ v = (real**) calloc( control->cm_solver_restart + 1, sizeof( real* ) );
+ for ( i = 0; i < control->cm_solver_restart + 1; ++i )
{
v[i] = (real *) calloc( system->N, sizeof( real ) );
}
@@ -614,8 +614,8 @@ void Cluster_Atoms( reax_system *system, static_storage *workspace )
// fprintf( stderr, "%d <-- %d\n", top, old_id );
reax_atom_Copy( &(new_atoms[top]), old_atom );
- Copy_Storage( system, workspace, top, old_id, old_atom->type,
- &num_H, v, s, t, orig_id, f_old );
+ Copy_Storage( system, workspace, control, top, old_id, old_atom->type,
+ &num_H, v, s, t, orig_id, f_old );
++top;
}
@@ -626,7 +626,7 @@ void Cluster_Atoms( reax_system *system, static_storage *workspace )
free( system->atoms );
- Free_Storage( workspace );
+ Free_Storage( workspace, control );
system->atoms = new_atoms;
Assign_New_Storage( workspace, v, s, t, orig_id, f_old );
diff --git a/sPuReMD/src/grid.h b/sPuReMD/src/grid.h
index 41d7b57edb2ba271bd41f71c6a77d5f179370910..cab825018606dd56e9f2b057c89a7d43287ed1f7 100644
--- a/sPuReMD/src/grid.h
+++ b/sPuReMD/src/grid.h
@@ -30,7 +30,7 @@ void Update_Grid( reax_system* );
int Shift( int, int, int, grid* );
-void Cluster_Atoms( reax_system*, static_storage* );
+void Cluster_Atoms( reax_system *, static_storage *, control_params * );
void Bin_Atoms( reax_system*, static_storage* );
diff --git a/sPuReMD/src/init_md.c b/sPuReMD/src/init_md.c
index 2ddcbeb6ebfe05491d718b41ee390d40a6903575..4b93ec83c4d6dcfbb8c32ba2c7736889ba6fc635 100644
--- a/sPuReMD/src/init_md.c
+++ b/sPuReMD/src/init_md.c
@@ -266,7 +266,7 @@ void Init_Taper( control_params *control )
void Init_Workspace( reax_system *system, control_params *control,
- static_storage *workspace )
+ static_storage *workspace )
{
int i;
@@ -406,18 +406,18 @@ void Init_Workspace( reax_system *system, control_params *control,
/* GMRES storage */
case GMRES_S:
case GMRES_H_S:
- workspace->y = (real *) calloc( RESTART + 1, sizeof( real ) );
- workspace->z = (real *) calloc( RESTART + 1, sizeof( real ) );
- workspace->g = (real *) calloc( RESTART + 1, sizeof( real ) );
- workspace->h = (real **) calloc( RESTART + 1, sizeof( real*) );
- workspace->hs = (real *) calloc( RESTART + 1, sizeof( real ) );
- workspace->hc = (real *) calloc( RESTART + 1, sizeof( real ) );
- workspace->rn = (real **) calloc( RESTART + 1, sizeof( real*) );
- workspace->v = (real **) calloc( RESTART + 1, sizeof( real*) );
-
- for ( i = 0; i < RESTART + 1; ++i )
+ workspace->y = (real *) calloc( control->cm_solver_restart + 1, sizeof( real ) );
+ workspace->z = (real *) calloc( control->cm_solver_restart + 1, sizeof( real ) );
+ workspace->g = (real *) calloc( control->cm_solver_restart + 1, sizeof( real ) );
+ workspace->h = (real **) calloc( control->cm_solver_restart + 1, sizeof( real*) );
+ workspace->hs = (real *) calloc( control->cm_solver_restart + 1, sizeof( real ) );
+ workspace->hc = (real *) calloc( control->cm_solver_restart + 1, sizeof( real ) );
+ workspace->rn = (real **) calloc( control->cm_solver_restart + 1, sizeof( real*) );
+ workspace->v = (real **) calloc( control->cm_solver_restart + 1, sizeof( real*) );
+
+ for ( i = 0; i < control->cm_solver_restart + 1; ++i )
{
- workspace->h[i] = (real *) calloc( RESTART + 1, sizeof( real ) );
+ workspace->h[i] = (real *) calloc( control->cm_solver_restart + 1, sizeof( real ) );
workspace->rn[i] = (real *) calloc( system->N_cm * 2, sizeof( real ) );
workspace->v[i] = (real *) calloc( system->N_cm, sizeof( real ) );
}
diff --git a/sPuReMD/src/lin_alg.c b/sPuReMD/src/lin_alg.c
index 65cea774c8edfebba473f9a12555c7d26828604e..31bd16e1b1052cd302ffaed2b40ff465eab2dc35 100644
--- a/sPuReMD/src/lin_alg.c
+++ b/sPuReMD/src/lin_alg.c
@@ -1177,9 +1177,8 @@ static void apply_preconditioner( const static_storage * const workspace, const
/* generalized minimual residual iterative solver for sparse linear systems */
int GMRES( const static_storage * const workspace, const control_params * const control,
- simulation_data * const data, const sparse_matrix * const H,
- const real * const b, const real tol, real * const x,
- const FILE * const fout, const int fresh_pre )
+ simulation_data * const data, const sparse_matrix * const H, const real * const b,
+ const real tol, real * const x, const FILE * const fout, const int fresh_pre )
{
int i, j, k, itr, N, g_j, g_itr;
real cc, tmp1, tmp2, temp, ret_temp, bnorm, time_start;
@@ -1214,7 +1213,7 @@ int GMRES( const static_storage * const workspace, const control_params * const
}
/* GMRES outer-loop */
- for ( itr = 0; itr < MAX_ITR; ++itr )
+ for ( itr = 0; itr < control->cm_solver_max_iters; ++itr )
{
/* calculate r0 */
#pragma omp master
@@ -1295,7 +1294,7 @@ int GMRES( const static_storage * const workspace, const control_params * const
}
/* GMRES inner-loop */
- for ( j = 0; j < RESTART && FABS(workspace->g[j]) / bnorm > tol; j++ )
+ for ( j = 0; j < control->cm_solver_restart && FABS(workspace->g[j]) / bnorm > tol; j++ )
{
/* matvec */
#pragma omp master
@@ -1501,28 +1500,28 @@ int GMRES( const static_storage * const workspace, const control_params * const
// fprintf(fout,"GMRES outer:%d, inner:%d iters - residual norm: %25.20f\n",
// itr, j, fabs( workspace->g[j] ) / bnorm );
- // data->timing.solver_iters += itr * RESTART + j;
+ // data->timing.solver_iters += itr * control->cm_solver_restart + j;
- if ( g_itr >= MAX_ITR )
+ if ( g_itr >= control->cm_solver_max_iters )
{
fprintf( stderr, "GMRES convergence failed\n" );
// return -1;
- return g_itr * (RESTART + 1) + g_j + 1;
+ return g_itr * (control->cm_solver_restart + 1) + g_j + 1;
}
- return g_itr * (RESTART + 1) + g_j + 1;
+ return g_itr * (control->cm_solver_restart + 1) + g_j + 1;
}
-int GMRES_HouseHolder( const static_storage * const workspace, const control_params * const control,
- simulation_data * const data, const sparse_matrix * const H,
- const real * const b, real tol, real * const x,
- const FILE * const fout, const int fresh_pre )
+int GMRES_HouseHolder( const static_storage * const workspace,
+ const control_params * const control, simulation_data * const data,
+ const sparse_matrix * const H, const real * const b, real tol,
+ real * const x, const FILE * const fout, const int fresh_pre )
{
int i, j, k, itr, N;
real cc, tmp1, tmp2, temp, bnorm;
- real v[10000], z[RESTART + 2][10000], w[RESTART + 2];
- real u[RESTART + 2][10000];
+ real v[10000], z[control->cm_solver_restart + 2][10000], w[control->cm_solver_restart + 2];
+ real u[control->cm_solver_restart + 2][10000];
N = H->n;
bnorm = Norm( b, N );
@@ -1536,7 +1535,7 @@ int GMRES_HouseHolder( const static_storage * const workspace, const control_par
// memset( x, 0, sizeof(real) * N );
/* GMRES outer-loop */
- for ( itr = 0; itr < MAX_ITR; ++itr )
+ for ( itr = 0; itr < control->cm_solver_max_iters; ++itr )
{
/* compute z = r0 */
Sparse_MatVec( H, x, workspace->b_prm );
@@ -1546,7 +1545,7 @@ int GMRES_HouseHolder( const static_storage * const workspace, const control_par
}
Vector_Sum( z[0], 1., workspace->b_prc, -1., workspace->b_prm, N );
- Vector_MakeZero( w, RESTART + 1 );
+ Vector_MakeZero( w, control->cm_solver_restart + 1 );
w[0] = Norm( z[0], N );
Vector_Copy( u[0], z[0], N );
@@ -1557,7 +1556,7 @@ int GMRES_HouseHolder( const static_storage * const workspace, const control_par
// fprintf( stderr, "\n\n%12.6f\n", w[0] );
/* GMRES inner-loop */
- for ( j = 0; j < RESTART && fabs( w[j] ) / bnorm > tol; j++ )
+ for ( j = 0; j < control->cm_solver_restart && fabs( w[j] ) / bnorm > tol; j++ )
{
/* compute v_j */
Vector_Scale( z[j], -2 * u[j][j], u[j], N );
@@ -1661,7 +1660,7 @@ int GMRES_HouseHolder( const static_storage * const workspace, const control_par
}
// fprintf( stderr, "y: " );
- // for( i = 0; i < RESTART+1; ++i )
+ // for( i = 0; i < control->cm_solver_restart+1; ++i )
// fprintf( stderr, "%8.3f ", workspace->y[i] );
@@ -1712,20 +1711,21 @@ int GMRES_HouseHolder( const static_storage * const workspace, const control_par
//fprintf( fout,"GMRES outer:%d, inner:%d iters - residual norm: %15.10f\n",
// itr, j, fabs( workspace->g[j] ) / bnorm );
- if ( itr >= MAX_ITR )
+ if ( itr >= control->cm_solver_max_iters )
{
fprintf( stderr, "GMRES convergence failed\n" );
// return -1;
- return itr * (RESTART + 1) + j + 1;
+ return itr * (control->cm_solver_restart + 1) + j + 1;
}
- return itr * (RESTART + 1) + j + 1;
+ return itr * (control->cm_solver_restart + 1) + j + 1;
}
/* Preconditioned Conjugate Gradient */
-int PCG( static_storage *workspace, sparse_matrix *A, real *b, real tol,
- sparse_matrix *L, sparse_matrix *U, real *x, FILE *fout )
+int PCG( static_storage *workspace, const control_params * const control,
+ sparse_matrix *A, real *b, real tol,
+ sparse_matrix *L, sparse_matrix *U, real *x, FILE *fout )
{
int i, N;
real tmp, alpha, beta, b_norm, r_norm;
@@ -1780,8 +1780,9 @@ int PCG( static_storage *workspace, sparse_matrix *A, real *b, real tol,
/* Conjugate Gradient */
-int CG( static_storage *workspace, sparse_matrix *H,
- real *b, real tol, real *x, FILE *fout )
+int CG( const static_storage * const workspace, const control_params * const control,
+ const sparse_matrix * const H, const real * const b, const real tol, real * const x,
+ const FILE * const fout )
{
int i, j, N;
real tmp, alpha, beta, b_norm;
@@ -1805,7 +1806,7 @@ int CG( static_storage *workspace, sparse_matrix *H,
// sqrt(sig_new), Norm(workspace->d,N), Norm(workspace->q,N) );
//fprintf( stderr, "sig_new: %f\n", sig_new );
- for ( i = 0; i < 300 && SQRT(sig_new) / b_norm > tol; ++i )
+ for ( i = 0; i < control->cm_solver_max_iters && SQRT(sig_new) / b_norm > tol; ++i )
{
//for( i = 0; i < 300 && sig_new > SQR(tol)*sig0; ++i ) {
Sparse_MatVec( H, workspace->d, workspace->q );
@@ -1842,8 +1843,9 @@ int CG( static_storage *workspace, sparse_matrix *H,
/* Steepest Descent */
-int SDM( static_storage *workspace, sparse_matrix *H,
- real *b, real tol, real *x, FILE *fout )
+int SDM( const static_storage * const workspace, const control_params * const control,
+ const sparse_matrix * const H, const real * const b, const real tol, real * const x,
+ const FILE * const fout )
{
int i, j, N;
real tmp, alpha, beta, b_norm;
@@ -1863,7 +1865,7 @@ int SDM( static_storage *workspace, sparse_matrix *H,
sig = Dot( workspace->r, workspace->d, N );
sig0 = sig;
- for ( i = 0; i < 300 && SQRT(sig) / b_norm > tol; ++i )
+ for ( i = 0; i < control->cm_solver_max_iters && SQRT(sig) / b_norm > tol; ++i )
{
Sparse_MatVec( H, workspace->d, workspace->q );
diff --git a/sPuReMD/src/lin_alg.h b/sPuReMD/src/lin_alg.h
index 6678e1b4edc48f88944833e4927d5f54c59da292..e293a96864a0458e745bbf88fade8ec90ec5e861 100644
--- a/sPuReMD/src/lin_alg.h
+++ b/sPuReMD/src/lin_alg.h
@@ -56,11 +56,13 @@ int GMRES_HouseHolder( const static_storage * const, const control_params * cons
const real * const, const real, real * const,
const FILE * const, const int );
-int CG( static_storage*, sparse_matrix*,
- real*, real, real*, FILE* );
+int CG( const static_storage * const, const control_params * const,
+ const sparse_matrix * const, const real * const, const real, real * const,
+ const FILE * const );
-int SDM( static_storage*, sparse_matrix*,
- real*, real, real*, FILE* );
+int SDM( const static_storage * const, const control_params * const,
+ const sparse_matrix * const, const real * const, const real, real * const,
+ const FILE * const );
real condest( const sparse_matrix * const, const sparse_matrix * const );
diff --git a/sPuReMD/src/mytypes.h b/sPuReMD/src/mytypes.h
index 826dee383d23785b2fcd657b0e94f3745609042c..7487086ec260f15f8731834ca444c02638365d17 100644
--- a/sPuReMD/src/mytypes.h
+++ b/sPuReMD/src/mytypes.h
@@ -113,9 +113,6 @@
#define MAX_dT 4.00
#define MIN_dT 0.00
-#define MAX_ITR 10
-#define RESTART 50
-
#define ZERO 0.000000000000000e+00
#define ALMOST_ZERO 1e-10
#define NEG_INF -1e10
@@ -530,6 +527,8 @@ typedef struct
unsigned int charge_method;
unsigned int cm_solver_type;
real cm_q_net;
+ unsigned int cm_solver_max_iters;
+ unsigned int cm_solver_restart;
real cm_solver_q_err;
real cm_domain_sparsity;
unsigned int cm_domain_sparsify_enabled;
diff --git a/sPuReMD/src/neighbors.c b/sPuReMD/src/neighbors.c
index fc318ef2d7d7022fc97df29548017299cd1f5e63..1208c826606e61bf912d3623a9d4011eb95f99ee 100644
--- a/sPuReMD/src/neighbors.c
+++ b/sPuReMD/src/neighbors.c
@@ -390,7 +390,7 @@ void Generate_Neighbor_Lists( reax_system *system, control_params *control,
// fprintf( stderr, "atoms sorted - " );
#ifdef REORDER_ATOMS
- Cluster_Atoms( system, workspace );
+ Cluster_Atoms( system, workspace, control );
// fprintf( stderr, "atoms clustered - " );
#endif