From 645dca23a831995bf723289497e3d531177c64d4 Mon Sep 17 00:00:00 2001
From: "Kurt A. O'Hearn" <ohearnk@msu.edu>
Date: Tue, 29 May 2018 15:45:47 -0400
Subject: [PATCH] PG-PuReMD: more function declaration / definition changes.
---
PG-PuReMD/src/allocate.c | 90 ++++++++++------------
PG-PuReMD/src/allocate.h | 29 ++++----
PG-PuReMD/src/box.c | 56 +++++++-------
PG-PuReMD/src/box.h | 35 ++++-----
PG-PuReMD/src/control.c | 4 +-
PG-PuReMD/src/control.h | 3 +-
PG-PuReMD/src/ffield.c | 131 +++++++++++++++++++--------------
PG-PuReMD/src/ffield.h | 4 +-
PG-PuReMD/src/geo_tools.c | 73 +++++++++---------
PG-PuReMD/src/geo_tools.h | 12 +--
PG-PuReMD/src/grid.c | 97 +++++++++++-------------
PG-PuReMD/src/grid.h | 10 +--
PG-PuReMD/src/index_utils.h | 12 ++-
PG-PuReMD/src/init_md.c | 112 ++++++++++++++--------------
PG-PuReMD/src/init_md.h | 20 ++---
PG-PuReMD/src/lin_alg.c | 18 +----
PG-PuReMD/src/list.c | 8 +-
PG-PuReMD/src/list.h | 8 +-
PG-PuReMD/src/lookup.c | 24 +++---
PG-PuReMD/src/lookup.h | 8 +-
PG-PuReMD/src/multi_body.c | 17 ++---
PG-PuReMD/src/puremd.c | 22 +++---
PG-PuReMD/src/reset_tools.c | 35 +++++----
PG-PuReMD/src/reset_tools.h | 20 ++---
PG-PuReMD/src/torsion_angles.c | 4 +-
PG-PuReMD/src/valence_angles.c | 3 -
26 files changed, 425 insertions(+), 430 deletions(-)
diff --git a/PG-PuReMD/src/allocate.c b/PG-PuReMD/src/allocate.c
index 3be198f2..f6829db1 100644
--- a/PG-PuReMD/src/allocate.c
+++ b/PG-PuReMD/src/allocate.c
@@ -38,7 +38,8 @@
#include "index_utils.h"
-void Init_Matrix_Row_Indices( sparse_matrix *H, int * max_row_entries )
+void Init_Matrix_Row_Indices( sparse_matrix * const H,
+ int * const max_row_entries )
{
int i;
@@ -58,8 +59,8 @@ void Init_Matrix_Row_Indices( sparse_matrix *H, int * max_row_entries )
* important: we cannot know the exact number of atoms that will fall into a
* process's box throughout the whole simulation. therefore
* we need to make upper bound estimates for various data structures */
-void PreAllocate_Space( reax_system *system, control_params *control,
- storage *workspace )
+void PreAllocate_Space( reax_system * const system, control_params * const control,
+ storage * const workspace )
{
/* determine capacity based on box vol & est atom volume */
system->local_cap = MAX( (int)(system->n * SAFE_ZONE), MIN_CAP );
@@ -85,7 +86,7 @@ void PreAllocate_Space( reax_system *system, control_params *control,
/************* system *************/
-void ReAllocate_System( reax_system *system, int local_cap, int total_cap )
+void ReAllocate_System( reax_system * const system, int local_cap, int total_cap )
{
system->my_atoms = srealloc( system->my_atoms, sizeof(reax_atom) * total_cap,
"ReAllocate_System::system->my_atoms" );
@@ -114,7 +115,7 @@ void ReAllocate_System( reax_system *system, int local_cap, int total_cap )
/************* workspace *************/
-void DeAllocate_Workspace( control_params *control, storage *workspace )
+void DeAllocate_Workspace( control_params * const control, storage * const workspace )
{
int i;
@@ -192,13 +193,6 @@ void DeAllocate_Workspace( control_params *control, storage *workspace )
// sfree( workspace->f_old, "DeAllocate_Workspace::f_old" );
sfree( workspace->v_const, "DeAllocate_Workspace::v_const" );
- /*workspace->realloc.far_nbrs = -1;
- workspace->realloc.Htop = -1;
- workspace->realloc.hbonds = -1;
- workspace->realloc.bonds = -1;
- workspace->realloc.num_3body = -1;
- workspace->realloc.gcell_atoms = -1;*/
-
/* storage for analysis */
if ( control->molecular_analysis || control->diffusion_coef )
{
@@ -240,8 +234,8 @@ void DeAllocate_Workspace( control_params *control, storage *workspace )
}
-void Allocate_Workspace( reax_system *system, control_params *control,
- storage *workspace, int local_cap, int total_cap )
+void Allocate_Workspace( reax_system * const system, control_params * const control,
+ storage * const workspace, int local_cap, int total_cap )
{
int i, total_real, total_rvec, local_rvec;
@@ -450,14 +444,15 @@ void Allocate_Workspace( reax_system *system, control_params *control,
}
-void Reallocate_Neighbor_List( reax_list *far_nbr_list, int n, int max_intrs )
+void Reallocate_Neighbor_List( reax_list * const far_nbr_list,
+ int n, int max_intrs )
{
Delete_List( far_nbr_list );
Make_List( n, max_intrs, TYP_FAR_NEIGHBOR, far_nbr_list );
}
-void Allocate_Matrix( sparse_matrix *H, int n, int n_max, int m )
+void Allocate_Matrix( sparse_matrix * const H, int n, int n_max, int m )
{
H->n = n;
H->n_max = n_max;
@@ -469,7 +464,7 @@ void Allocate_Matrix( sparse_matrix *H, int n, int n_max, int m )
}
-void Deallocate_Matrix( sparse_matrix *H )
+void Deallocate_Matrix( sparse_matrix * const H )
{
H->n = 0;
H->n_max = 0;
@@ -481,14 +476,14 @@ void Deallocate_Matrix( sparse_matrix *H )
}
-static void Reallocate_Matrix( sparse_matrix *H, int n, int n_max, int m )
+static void Reallocate_Matrix( sparse_matrix * const H, int n, int n_max, int m )
{
Deallocate_Matrix( H );
Allocate_Matrix( H, n, n_max, m );
}
-void Reallocate_HBonds_List( reax_system *system, reax_list *hbond_list )
+void Reallocate_HBonds_List( reax_system * const system, reax_list * const hbond_list )
{
Delete_List( hbond_list );
// Make_List( system->Hcap, system->total_hbonds, TYP_HBOND, hbond_list );
@@ -496,7 +491,7 @@ void Reallocate_HBonds_List( reax_system *system, reax_list *hbond_list )
}
-void Reallocate_Bonds_List( reax_system *system, reax_list *bond_list )
+void Reallocate_Bonds_List( reax_system * const system, reax_list * const bond_list )
{
Delete_List( bond_list );
Make_List( system->total_cap, system->total_bonds, TYP_BOND, bond_list );
@@ -504,18 +499,15 @@ void Reallocate_Bonds_List( reax_system *system, reax_list *bond_list )
/************* grid *************/
-int Estimate_GCell_Population( reax_system* system, MPI_Comm comm )
+int Estimate_GCell_Population( reax_system * const system, MPI_Comm comm )
{
int d, i, j, k, l, max_atoms, my_max, all_max;
ivec c;
- grid *g;
+ grid * const g = &system->my_grid;
grid_cell *gc;
- simulation_box *my_ext_box;
- reax_atom *atoms;
+ simulation_box * const my_ext_box = &system->my_ext_box;
+ reax_atom * const atoms = system->my_atoms;
- my_ext_box = &system->my_ext_box;
- g = &system->my_grid;
- atoms = system->my_atoms;
Reset_Grid( g );
for ( l = 0; l < system->n; l++ )
@@ -579,14 +571,13 @@ int Estimate_GCell_Population( reax_system* system, MPI_Comm comm )
}
-void Allocate_Grid( reax_system *system, MPI_Comm comm )
+void Allocate_Grid( reax_system * const system, MPI_Comm comm )
{
int i, j, k;
- grid *g;
+ grid * const g = &system->my_grid;
grid_cell *gc;
int total;
- g = &system->my_grid;
total = g->ncells[0] * g->ncells[1] * g->ncells[2];
/* allocate gcell reordering space */
@@ -650,7 +641,7 @@ void Allocate_Grid( reax_system *system, MPI_Comm comm )
}
-void Deallocate_Grid( grid *g )
+void Deallocate_Grid( grid * const g )
{
int i, j, k;
grid_cell *gc;
@@ -691,8 +682,8 @@ void Deallocate_Grid( grid *g )
*
* Note: buffers are (void *), type cast to the correct pointer type to access
* the allocated buffers */
-void Allocate_MPI_Buffers( mpi_datatypes *mpi_data, int est_recv,
- neighbor_proc *my_nbrs )
+void Allocate_MPI_Buffers( mpi_datatypes * const mpi_data, int est_recv,
+ neighbor_proc * const my_nbrs )
{
int i;
mpi_out_data *mpi_buf;
@@ -722,7 +713,7 @@ void Allocate_MPI_Buffers( mpi_datatypes *mpi_data, int est_recv,
}
-void Deallocate_MPI_Buffers( mpi_datatypes *mpi_data )
+void Deallocate_MPI_Buffers( mpi_datatypes * const mpi_data )
{
int i;
mpi_out_data *mpi_buf;
@@ -739,23 +730,19 @@ void Deallocate_MPI_Buffers( mpi_datatypes *mpi_data )
}
-void ReAllocate( reax_system *system, control_params *control,
- simulation_data *data, storage *workspace, reax_list **lists,
- mpi_datatypes *mpi_data )
+void ReAllocate( reax_system * const system, control_params * const control,
+ simulation_data * const data, storage * const workspace, reax_list ** const lists,
+ mpi_datatypes * const mpi_data )
{
int i, j, k;
int nflag, Nflag, Hflag, mpi_flag, total_send;
int renbr;
- reallocate_data *realloc;
- sparse_matrix *H;
- grid *g;
+ reallocate_data * const realloc = &workspace->realloc;
+ sparse_matrix * const H = &workspace->H;
+ grid * const g = &system->my_grid;
neighbor_proc *nbr_pr;
mpi_out_data *nbr_data;
- realloc = &workspace->realloc;
- g = &system->my_grid;
- H = &workspace->H;
-
/* IMPORTANT: LOOSE ZONES CHECKS ARE DISABLED FOR NOW BY &&'ing with FALSE!!! */
nflag = FALSE;
if ( system->n >= DANGER_ZONE * system->local_cap ||
@@ -903,9 +890,9 @@ void ReAllocate( reax_system *system, control_params *control,
realloc->gcell_atoms = -1;
}
- /* mpi buffers */
- // we have to be at a renbring step -
- // to ensure correct values at mpi_buffers for update_boundary_positions
+ /* mpi buffers:
+ * we have to be at a renbring step
+ * to ensure correct values at mpi_buffers for update_boundary_positions */
if ( renbr == FALSE )
{
mpi_flag = FALSE;
@@ -921,8 +908,9 @@ void ReAllocate( reax_system *system, control_params *control,
mpi_flag = FALSE;
for ( i = 0; i < MAX_NBRS; ++i )
{
- nbr_pr = &( system->my_nbrs[i] );
- nbr_data = &( mpi_data->out_buffers[i] );
+ nbr_pr = &system->my_nbrs[i];
+ nbr_data = &mpi_data->out_buffers[i];
+
if ( nbr_data->cnt >= nbr_pr->est_send * DANGER_ZONE )
{
mpi_flag = TRUE;
@@ -950,8 +938,8 @@ void ReAllocate( reax_system *system, control_params *control,
total_send = 0;
for ( i = 0; i < MAX_NBRS; ++i )
{
- nbr_pr = &( system->my_nbrs[i] );
- nbr_data = &( mpi_data->out_buffers[i] );
+ nbr_pr = &system->my_nbrs[i];
+ nbr_data = &mpi_data->out_buffers[i];
nbr_pr->est_send = MAX( nbr_data->cnt * SAFER_ZONE, MIN_SEND );
total_send += nbr_pr->est_send;
}
diff --git a/PG-PuReMD/src/allocate.h b/PG-PuReMD/src/allocate.h
index db9b865d..271ea91c 100644
--- a/PG-PuReMD/src/allocate.h
+++ b/PG-PuReMD/src/allocate.h
@@ -29,33 +29,34 @@
extern "C" {
#endif
-void Init_Matrix_Row_Indices( sparse_matrix *, int * );
+void Init_Matrix_Row_Indices( sparse_matrix * const, int * const );
-void PreAllocate_Space( reax_system*, control_params*, storage* );
+void PreAllocate_Space( reax_system * const, control_params * const, storage * const );
-void ReAllocate_System( reax_system*, int, int );
+void ReAllocate_System( reax_system * const, int, int );
-void Allocate_Workspace( reax_system*, control_params*, storage*,
+void Allocate_Workspace( reax_system * const, control_params * const, storage * const,
int, int );
-void Allocate_Grid( reax_system*, MPI_Comm );
+void Allocate_Grid( reax_system * const, MPI_Comm );
-void Deallocate_Grid( grid* );
+void Deallocate_Grid( grid * const );
-void Allocate_MPI_Buffers( mpi_datatypes*, int, neighbor_proc* );
+void Allocate_MPI_Buffers( mpi_datatypes * const, int, neighbor_proc * const );
-void Allocate_Matrix( sparse_matrix*, int, int, int );
+void Allocate_Matrix( sparse_matrix * const, int, int, int );
-void Deallocate_Matrix( sparse_matrix * );
+void Deallocate_Matrix( sparse_matrix * const );
-int Allocate_HBond_List( int, int, int*, int*, reax_list* );
+int Allocate_HBond_List( int, int, int * const, int * const, reax_list * const );
-int Allocate_Bond_List( int, int*, reax_list* );
+int Allocate_Bond_List( int, int * const, reax_list * const );
-void Deallocate_MPI_Buffers( mpi_datatypes * );
+void Deallocate_MPI_Buffers( mpi_datatypes * const );
-void ReAllocate( reax_system*, control_params*, simulation_data*, storage*,
- reax_list**, mpi_datatypes* );
+void ReAllocate( reax_system * const, control_params * const,
+ simulation_data * const, storage * const,
+ reax_list** const, mpi_datatypes * const );
#ifdef __cplusplus
}
diff --git a/PG-PuReMD/src/box.c b/PG-PuReMD/src/box.c
index 397d5e0d..5b7ee96f 100644
--- a/PG-PuReMD/src/box.c
+++ b/PG-PuReMD/src/box.c
@@ -29,7 +29,7 @@
#include "vector.h"
-void Make_Consistent( simulation_box* box )
+void Make_Consistent( simulation_box * const box )
{
real one_vol;
@@ -140,14 +140,14 @@ void Make_Consistent( simulation_box* box )
/* setup the entire simulation box */
void Setup_Big_Box( real a, real b, real c, real alpha, real beta, real gamma,
- simulation_box* box )
+ simulation_box * const box )
{
double c_alpha, c_beta, c_gamma, s_gamma, zi;
if ( IS_NAN_REAL(a) || IS_NAN_REAL(b) || IS_NAN_REAL(c)
|| IS_NAN_REAL(alpha) || IS_NAN_REAL(beta) || IS_NAN_REAL(gamma) )
{
- fprintf( stderr, "Invalid simulation box boundaries for big box (NaN). Terminating...\n" );
+ fprintf( stderr, "[ERROR] Invalid simulation box boundaries for big box (NaN). Terminating...\n" );
exit( INVALID_INPUT );
}
@@ -178,7 +178,7 @@ void Setup_Big_Box( real a, real b, real c, real alpha, real beta, real gamma,
}
-void Init_Box( rtensor box_tensor, simulation_box *box )
+void Init_Box( rtensor box_tensor, simulation_box * const box )
{
rvec_MakeZero( box->min );
rtensor_Copy( box->box, box_tensor );
@@ -192,13 +192,12 @@ void Init_Box( rtensor box_tensor, simulation_box *box )
/* setup my simulation box -- only the region that I own */
-void Setup_My_Box( reax_system *system, control_params *control )
+void Setup_My_Box( reax_system * const system, control_params * const control )
{
int d;
- simulation_box *big_box, *my_box;
+ simulation_box * const big_box = &system->big_box;
+ simulation_box * const my_box = &system->my_box;
- big_box = &(system->big_box);
- my_box = &(system->my_box);
rtensor_MakeZero( my_box->box );
for ( d = 0; d < 3; ++d )
@@ -216,18 +215,18 @@ void Setup_My_Box( reax_system *system, control_params *control )
/* setup my extended box -- my box together with the ghost regions */
-void Setup_My_Ext_Box( reax_system *system, control_params *control )
+void Setup_My_Ext_Box( reax_system * const system, control_params * const control )
{
int d;
ivec native_gcells, ghost_gcells;
rvec gcell_len;
- simulation_box *big_box, *my_box, *my_ext_box;
- boundary_cutoff *bc;
+ boundary_cutoff * const bc = &system->bndry_cuts;
+ simulation_box * const my_box = &system->my_box;
+ simulation_box * const my_ext_box = &system->my_ext_box;
+#if defined(DEBUG_FOCUS)
+ simulation_box * const big_box = &system->big_box;
+#endif
- big_box = &(system->big_box);
- my_box = &(system->my_box);
- my_ext_box = &(system->my_ext_box);
- bc = &(system->bndry_cuts);
rtensor_MakeZero( my_ext_box->box );
for ( d = 0; d < 3; ++d )
@@ -255,10 +254,9 @@ void Setup_My_Ext_Box( reax_system *system, control_params *control )
/******************** initialize parallel environment ***********************/
-void Setup_Boundary_Cutoffs( reax_system *system, control_params *control )
+void Setup_Boundary_Cutoffs( reax_system * const system, control_params * const control )
{
- boundary_cutoff *bc;
- bc = &(system->bndry_cuts);
+ boundary_cutoff * const bc = &system->bndry_cuts;
bc->ghost_nonb = control->nonb_cut;
bc->ghost_hbond = control->hbond_cut;
@@ -274,8 +272,8 @@ void Setup_Boundary_Cutoffs( reax_system *system, control_params *control )
}
-void Setup_Environment( reax_system *system, control_params *control,
- mpi_datatypes *mpi_data )
+void Setup_Environment( reax_system * const system, control_params * const control,
+ mpi_datatypes * const mpi_data )
{
ivec periodic = {1, 1, 1};
#if defined(DEBUG_FOCUS)
@@ -284,8 +282,8 @@ void Setup_Environment( reax_system *system, control_params *control,
/* initialize communicator - 3D mesh with wrap-arounds = 3D torus */
MPI_Cart_create( MPI_COMM_WORLD, 3, control->procs_by_dim, periodic, 1,
- &(mpi_data->comm_mesh3D) );
- MPI_Comm_rank( mpi_data->comm_mesh3D, &(system->my_rank) );
+ &mpi_data->comm_mesh3D );
+ MPI_Comm_rank( mpi_data->comm_mesh3D, &system->my_rank );
MPI_Cart_coords( mpi_data->comm_mesh3D, system->my_rank, 3,
system->my_coords );
@@ -299,11 +297,11 @@ void Setup_Environment( reax_system *system, control_params *control,
system->my_rank, system->my_coords[0],
system->my_coords[1], system->my_coords[2] );
sprintf( temp, "p%d big_box", system->my_rank );
- Print_Box( &(system->big_box), temp, stderr );
+ Print_Box( &system->big_box, temp, stderr );
sprintf( temp, "p%d my_box", system->my_rank );
- Print_Box( &(system->my_box), temp, stderr );
+ Print_Box( &system->my_box, temp, stderr );
sprintf( temp, "p%d ext_box", system->my_rank );
- Print_Box( &(system->my_ext_box), temp, stderr );
+ Print_Box( &system->my_ext_box, temp, stderr );
MPI_Barrier( MPI_COMM_WORLD );
fprintf( stderr, "p%d: parallel environment initialized\n",
@@ -312,8 +310,8 @@ void Setup_Environment( reax_system *system, control_params *control,
}
-void Scale_Box( reax_system *system, control_params *control,
- simulation_data *data, mpi_datatypes *mpi_data )
+void Scale_Box( reax_system * const system, control_params * const control,
+ simulation_data * const data, mpi_datatypes * const mpi_data )
{
int i, d;
real dt, lambda;
@@ -373,7 +371,7 @@ void Scale_Box( reax_system *system, control_params *control,
/* Scale velocities and positions at t+dt */
for ( i = 0; i < system->n; ++i )
{
- atom = &(system->my_atoms[i]);
+ atom = &system->my_atoms[i];
rvec_Scale( atom->v, lambda, atom->v );
atom->x[0] = mu[0] * atom->x[0];
atom->x[1] = mu[1] * atom->x[1];
@@ -401,7 +399,7 @@ void Scale_Box( reax_system *system, control_params *control,
}
-real Metric_Product( rvec x1, rvec x2, simulation_box* box )
+real Metric_Product( rvec x1, rvec x2, simulation_box * const box )
{
int i, j;
real dist = 0.0, tmp;
diff --git a/PG-PuReMD/src/box.h b/PG-PuReMD/src/box.h
index fa92c722..c210b5d1 100644
--- a/PG-PuReMD/src/box.h
+++ b/PG-PuReMD/src/box.h
@@ -29,44 +29,45 @@
extern "C" {
#endif
-void Make_Consistent( simulation_box* );
+void Make_Consistent( simulation_box * const );
/* initializes simulation boxes */
-void Setup_Big_Box( real, real, real, real, real, real, simulation_box* );
+void Setup_Big_Box( real, real, real, real, real, real, simulation_box * const );
-void Init_Box( rtensor, simulation_box* );
+void Init_Box( rtensor, simulation_box * const );
-void Setup_My_Box( reax_system*, control_params* );
+void Setup_My_Box( reax_system * const, control_params * const );
-void Setup_My_Ext_Box( reax_system*, control_params* );
+void Setup_My_Ext_Box( reax_system * const, control_params * const );
-void Setup_Environment( reax_system*, control_params*, mpi_datatypes* );
+void Setup_Environment( reax_system * const, control_params * const,
+ mpi_datatypes * const );
/* scales simulation box for NPT ensembles */
-void Scale_Box( reax_system*, control_params*,
- simulation_data*, mpi_datatypes* );
+void Scale_Box( reax_system * const, control_params * const,
+ simulation_data * const, mpi_datatypes * const );
/* applies transformation to/from Cartesian/ Triclinic coordinates */
/* use -1 flag for Cartesian -> Triclinic and +1 for otherway */
-//void Transform( rvec, simulation_box*, char, rvec );
+//void Transform( rvec, simulation_box * const, char, rvec );
-//void Distance_on_T3_Gen( rvec, rvec, simulation_box*, rvec );
+//void Distance_on_T3_Gen( rvec, rvec, simulation_box * const, rvec );
-//void Inc_on_T3_Gen( rvec, rvec, simulation_box* );
+//void Inc_on_T3_Gen( rvec, rvec, simulation_box * const );
-//int Get_Nbr_Box( simulation_box*, int, int, int );
+//int Get_Nbr_Box( simulation_box * const, int, int, int );
-//rvec Get_Nbr_Box_Press( simulation_box*, int, int, int );
+//rvec Get_Nbr_Box_Press( simulation_box * const, int, int, int );
-//void Inc_Nbr_Box_Press( simulation_box*, int, int, int, rvec );
+//void Inc_Nbr_Box_Press( simulation_box * const, int, int, int, rvec );
/* these functions assume that the coordinates are in triclinic system
this function returns cartesian norm but triclinic distance vector */
-//real Sq_Distance_on_T3( rvec, rvec, simulation_box*, rvec );
+//real Sq_Distance_on_T3( rvec, rvec, simulation_box * const, rvec );
-//void Inc_on_T3( rvec, rvec, simulation_box* );
+//void Inc_on_T3( rvec, rvec, simulation_box * const );
-//real Metric_Product( rvec, rvec, simulation_box* );
+//real Metric_Product( rvec, rvec, simulation_box * const );
#ifdef __cplusplus
}
diff --git a/PG-PuReMD/src/control.c b/PG-PuReMD/src/control.c
index b57fdaae..7232f004 100644
--- a/PG-PuReMD/src/control.c
+++ b/PG-PuReMD/src/control.c
@@ -30,8 +30,8 @@
#endif
-void Read_Control_File( const char * const control_file, control_params* control,
- output_controls *out_control )
+void Read_Control_File( const char * const control_file, control_params * const control,
+ output_controls * const out_control )
{
FILE *fp;
char *s, **tmp;
diff --git a/PG-PuReMD/src/control.h b/PG-PuReMD/src/control.h
index d7486f35..137b28a7 100644
--- a/PG-PuReMD/src/control.h
+++ b/PG-PuReMD/src/control.h
@@ -29,7 +29,8 @@
extern "C" {
#endif
-void Read_Control_File( const char * const, control_params*, output_controls* );
+void Read_Control_File( const char * const, control_params * const,
+ output_controls * const );
#ifdef __cplusplus
}
diff --git a/PG-PuReMD/src/ffield.c b/PG-PuReMD/src/ffield.c
index 26605995..e83142e8 100644
--- a/PG-PuReMD/src/ffield.c
+++ b/PG-PuReMD/src/ffield.c
@@ -30,8 +30,8 @@
#endif
-void Read_Force_Field_File( const char * const ffield_file, reax_interaction *reax,
- reax_system *system, control_params *control )
+void Read_Force_Field_File( const char * const ffield_file, reax_interaction * const reax,
+ reax_system * const system, control_params * const control )
{
FILE *fp;
char *s;
@@ -60,7 +60,11 @@ void Read_Force_Field_File( const char * const ffield_file, reax_interaction *re
c = Tokenize( s, &tmp );
/* reading the number of global parameters */
- n = atoi(tmp[0]);
+ if ( c > 0 )
+ {
+ n = atoi(tmp[0]);
+ }
+
if ( n < 1 )
{
fprintf( stderr, "[WARNING] p%d: number of globals in ffield file is 0!\n",
@@ -74,11 +78,14 @@ void Read_Force_Field_File( const char * const ffield_file, reax_interaction *re
/* see reax_types.h for mapping between l[i] and the lambdas used in ff */
for (i = 0; i < n; i++)
{
- fgets(s, MAX_LINE, fp);
- c = Tokenize(s, &tmp);
+ fgets( s, MAX_LINE, fp );
+ c = Tokenize( s, &tmp );
- val = (real) atof(tmp[0]);
- reax->gp.l[i] = val;
+ if ( c > 0 )
+ {
+ val = (real) atof(tmp[0]);
+ reax->gp.l[i] = val;
+ }
}
control->bo_cut = 0.01 * reax->gp.l[29];
@@ -88,7 +95,10 @@ void Read_Force_Field_File( const char * const ffield_file, reax_interaction *re
/* next line is number of atom types and some comments */
fgets( s, MAX_LINE, fp );
c = Tokenize( s, &tmp );
- reax->num_atom_types = atoi(tmp[0]);
+ if ( c > 0 )
+ {
+ reax->num_atom_types = atoi(tmp[0]);
+ }
/* 3 lines of comments */
fgets(s, MAX_LINE, fp);
@@ -134,7 +144,11 @@ void Read_Force_Field_File( const char * const ffield_file, reax_interaction *re
for ( j = 0; j < (int)(strlen(tmp[0])); ++j )
{
- reax->sbp[i].name[j] = toupper( tmp[0][j] );
+ if ( c > 0 )
+ {
+ reax->sbp[i].name[j] = toupper( tmp[0][j] );
+ --c;
+ }
}
#if defined(DEBUG_FOCUS)
@@ -164,59 +178,68 @@ void Read_Force_Field_File( const char * const ffield_file, reax_interaction *re
fgets( s, MAX_LINE, fp );
c = Tokenize( s, &tmp );
- val = atof(tmp[0]);
- reax->sbp[i].alpha = val;
- val = atof(tmp[1]);
- reax->sbp[i].gamma_w = val;
- val = atof(tmp[2]);
- reax->sbp[i].valency_boc = val;
- val = atof(tmp[3]);
- reax->sbp[i].p_ovun5 = val;
- val = atof(tmp[4]);
- val = atof(tmp[5]);
- reax->sbp[i].chi = val;
- val = atof(tmp[6]);
- reax->sbp[i].eta = 2.0 * val;
- val = atof(tmp[7]);
- reax->sbp[i].p_hbond = (int)val;
+ if ( c >= 8 )
+ {
+ val = atof(tmp[0]);
+ reax->sbp[i].alpha = val;
+ val = atof(tmp[1]);
+ reax->sbp[i].gamma_w = val;
+ val = atof(tmp[2]);
+ reax->sbp[i].valency_boc = val;
+ val = atof(tmp[3]);
+ reax->sbp[i].p_ovun5 = val;
+ val = atof(tmp[4]);
+ val = atof(tmp[5]);
+ reax->sbp[i].chi = val;
+ val = atof(tmp[6]);
+ reax->sbp[i].eta = 2.0 * val;
+ val = atof(tmp[7]);
+ reax->sbp[i].p_hbond = (int)val;
+ }
/* line 3 */
fgets( s, MAX_LINE, fp );
c = Tokenize( s, &tmp );
- val = atof(tmp[0]);
- reax->sbp[i].r_pi_pi = val;
- val = atof(tmp[1]);
- reax->sbp[i].p_lp2 = val;
- val = atof(tmp[2]);
- val = atof(tmp[3]);
- reax->sbp[i].b_o_131 = val;
- val = atof(tmp[4]);
- reax->sbp[i].b_o_132 = val;
- val = atof(tmp[5]);
- reax->sbp[i].b_o_133 = val;
- val = atof(tmp[6]);
- val = atof(tmp[7]);
+ if ( c >= 8 )
+ {
+ val = atof(tmp[0]);
+ reax->sbp[i].r_pi_pi = val;
+ val = atof(tmp[1]);
+ reax->sbp[i].p_lp2 = val;
+ val = atof(tmp[2]);
+ val = atof(tmp[3]);
+ reax->sbp[i].b_o_131 = val;
+ val = atof(tmp[4]);
+ reax->sbp[i].b_o_132 = val;
+ val = atof(tmp[5]);
+ reax->sbp[i].b_o_133 = val;
+ val = atof(tmp[6]);
+ val = atof(tmp[7]);
+ }
/* line 4 */
fgets( s, MAX_LINE, fp );
c = Tokenize( s, &tmp );
- val = atof(tmp[0]);
- reax->sbp[i].p_ovun2 = val;
- val = atof(tmp[1]);
- reax->sbp[i].p_val3 = val;
- val = atof(tmp[2]);
- val = atof(tmp[3]);
- reax->sbp[i].valency_val = val;
- val = atof(tmp[4]);
- reax->sbp[i].p_val5 = val;
- val = atof(tmp[5]);
- reax->sbp[i].rcore2 = val;
- val = atof(tmp[6]);
- reax->sbp[i].ecore2 = val;
- val = atof(tmp[7]);
- reax->sbp[i].acore2 = val;
+ if ( c >= 8 )
+ {
+ val = atof(tmp[0]);
+ reax->sbp[i].p_ovun2 = val;
+ val = atof(tmp[1]);
+ reax->sbp[i].p_val3 = val;
+ val = atof(tmp[2]);
+ val = atof(tmp[3]);
+ reax->sbp[i].valency_val = val;
+ val = atof(tmp[4]);
+ reax->sbp[i].p_val5 = val;
+ val = atof(tmp[5]);
+ reax->sbp[i].rcore2 = val;
+ val = atof(tmp[6]);
+ reax->sbp[i].ecore2 = val;
+ val = atof(tmp[7]);
+ reax->sbp[i].acore2 = val;
+ }
/* Inner-wall */
if ( reax->sbp[i].rcore2 > 0.01 && reax->sbp[i].acore2 > 0.01 )
@@ -321,7 +344,7 @@ void Read_Force_Field_File( const char * const ffield_file, reax_interaction *re
l = atoi(tmp[0]);
/* a line of comments */
- fgets(s, MAX_LINE, fp);
+ fgets( s, MAX_LINE, fp );
for (i = 0; i < l; i++)
{
@@ -334,7 +357,7 @@ void Read_Force_Field_File( const char * const ffield_file, reax_interaction *re
index1 = j * __N + k;
index2 = k * __N + j;
- if (j < reax->num_atom_types && k < reax->num_atom_types)
+ if ( j < reax->num_atom_types && k < reax->num_atom_types )
{
val = atof(tmp[2]);
reax->tbp[ index1 ].De_s = val;
diff --git a/PG-PuReMD/src/ffield.h b/PG-PuReMD/src/ffield.h
index 125464ba..ff069d52 100644
--- a/PG-PuReMD/src/ffield.h
+++ b/PG-PuReMD/src/ffield.h
@@ -29,8 +29,8 @@
extern "C" {
#endif
-void Read_Force_Field_File( const char * const, reax_interaction*,
- reax_system *, control_params* );
+void Read_Force_Field_File( const char * const, reax_interaction * const,
+ reax_system * const, control_params * const );
#ifdef __cplusplus
}
diff --git a/PG-PuReMD/src/geo_tools.c b/PG-PuReMD/src/geo_tools.c
index 6339d905..45b9cec9 100644
--- a/PG-PuReMD/src/geo_tools.c
+++ b/PG-PuReMD/src/geo_tools.c
@@ -30,7 +30,7 @@
/********************* geo format routines ******************/
-void Count_Geo_Atoms( FILE *geo, reax_system *system )
+static void Count_Geo_Atoms( FILE *geo, reax_system * const system )
{
int i, serial;
rvec x;
@@ -45,11 +45,14 @@ void Count_Geo_Atoms( FILE *geo, reax_system *system )
{
fscanf( geo, CUSTOM_ATOM_FORMAT,
&serial, element, name, &x[0], &x[1], &x[2] );
- Fit_to_Periodic_Box( &(system->big_box), &x );
+
+ Fit_to_Periodic_Box( &system->big_box, &x );
/* if the point is inside my_box, add it to my lists */
- if ( is_Inside_Box(&(system->my_box), x) )
+ if ( is_Inside_Box( &system->my_box, x ) == TRUE )
+ {
++system->n;
+ }
}
system->N = system->n;
@@ -67,8 +70,9 @@ void Count_Geo_Atoms( FILE *geo, reax_system *system )
}
-void Read_Geo_File( const char * const geo_file, reax_system* system, control_params *control,
- simulation_data *data, storage *workspace, mpi_datatypes *mpi_data )
+void Read_Geo_File( const char * const geo_file, reax_system * const system,
+ control_params * const control, simulation_data * const data,
+ storage * const workspace, mpi_datatypes * const mpi_data )
{
int i, j, serial, top;
char descriptor[9];
@@ -100,7 +104,7 @@ void Read_Geo_File( const char * const geo_file, reax_system* system, control_pa
fscanf( geo, CUSTOM_ATOM_FORMAT,
&serial, element, name, &x[0], &x[1], &x[2] );
- Fit_to_Periodic_Box( &(system->big_box), &x );
+ Fit_to_Periodic_Box( &system->big_box, &x );
#if defined(DEBUG)
fprintf( stderr, "atom%d: %s %s %f %f %f\n",
@@ -108,7 +112,7 @@ void Read_Geo_File( const char * const geo_file, reax_system* system, control_pa
#endif
/* if the point is inside my_box, add it to my list */
- if ( is_Inside_Box(&(system->my_box), x) )
+ if ( is_Inside_Box( &system->my_box, x ) == TRUE )
{
atom = &system->my_atoms[top];
atom->orig_id = serial;
@@ -145,7 +149,7 @@ void Read_Geo_File( const char * const geo_file, reax_system* system, control_pa
}
-int Read_Box_Info( reax_system *system, FILE *geo, int geo_format )
+static int Read_Box_Info( reax_system * const system, FILE *geo, int geo_format )
{
char *cryst;
char line[MAX_LINE + 1];
@@ -153,8 +157,7 @@ int Read_Box_Info( reax_system *system, FILE *geo, int geo_format )
char s_a[12], s_b[12], s_c[12], s_alpha[12], s_beta[12], s_gamma[12];
char s_group[12], s_zValue[12];
- /* initialize variables */
- fseek( geo, 0, SEEK_SET ); // set the pointer to the beginning of the file
+ fseek( geo, 0, SEEK_SET );
switch ( geo_format )
{
@@ -166,7 +169,7 @@ int Read_Box_Info( reax_system *system, FILE *geo, int geo_format )
}
/* locate the cryst line in the geo file, read it and
- initialize the big box */
+ * initialize the big box */
while ( fgets( line, MAX_LINE, geo ) )
{
if ( strncmp( line, cryst, 6 ) == 0 )
@@ -183,7 +186,7 @@ int Read_Box_Info( reax_system *system, FILE *geo, int geo_format )
/* compute full volume tensor from the angles */
Setup_Big_Box( atof(s_a), atof(s_b), atof(s_c),
atof(s_alpha), atof(s_beta), atof(s_gamma),
- &(system->big_box) );
+ &system->big_box );
return SUCCESS;
}
@@ -193,15 +196,14 @@ int Read_Box_Info( reax_system *system, FILE *geo, int geo_format )
}
-void Count_PDB_Atoms( FILE *geo, reax_system *system )
+static void Count_PDB_Atoms( FILE *geo, reax_system * const system )
{
char *endptr = NULL;
char line[MAX_LINE + 1];
char s_x[9], s_y[9], s_z[9];
rvec x;
- /* initialize variables */
- fseek( geo, 0, SEEK_SET ); /* set the pointer to the beginning of the file */
+ fseek( geo, 0, SEEK_SET );
system->bigN = 0;
system->n = 0;
system->N = 0;
@@ -222,17 +224,13 @@ void Count_PDB_Atoms( FILE *geo, reax_system *system )
s_z[8] = 0;
Make_Point( strtod( s_x, &endptr ), strtod( s_y, &endptr ),
strtod( s_z, &endptr ), &x );
- Fit_to_Periodic_Box( &(system->big_box), &x );
+ Fit_to_Periodic_Box( &system->big_box, &x );
/* if the point is inside my_box, add it to my lists */
- if ( is_Inside_Box(&(system->my_box), x) )
+ if ( is_Inside_Box( &system->my_box, x ) == TRUE )
{
++system->n;
}
-// if ( is_Inside_Box(&(system->my_ext_box), x) )
-// {
-// ++system->N;
-// }
}
}
@@ -247,8 +245,9 @@ void Count_PDB_Atoms( FILE *geo, reax_system *system )
}
-void Read_PDB_File( const char * const pdb_file, reax_system* system, control_params *control,
- simulation_data *data, storage *workspace, mpi_datatypes *mpi_data )
+void Read_PDB_File( const char * const pdb_file, reax_system * const system,
+ control_params * const control, simulation_data * const data,
+ storage * const workspace, mpi_datatypes * const mpi_data )
{
FILE *pdb;
char **tmp;
@@ -272,8 +271,8 @@ void Read_PDB_File( const char * const pdb_file, reax_system* system, control_pa
/* read box information */
if ( Read_Box_Info( system, pdb, PDB ) == FAILURE )
{
- fprintf( stderr, "Read_Box_Info: no CRYST line in the pdb file!" );
- fprintf( stderr, "terminating...\n" );
+ fprintf( stderr, "[ERROR] Read_Box_Info: no CRYST line in the pdb file!" );
+ fprintf( stderr, " Terminating...\n" );
MPI_Abort( MPI_COMM_WORLD, INVALID_GEO );
}
@@ -281,10 +280,11 @@ void Read_PDB_File( const char * const pdb_file, reax_system* system, control_pa
Count_PDB_Atoms( pdb, system );
PreAllocate_Space( system, control, workspace );
- /* start reading and processing the pdb file */
#if defined(DEBUG_FOCUS)
fprintf( stderr, "p%d: starting to read the pdb file\n", system->my_rank );
#endif
+
+ /* start reading and processing the pdb file */
fseek( pdb, 0, SEEK_SET );
c = 0;
c1 = 0;
@@ -378,17 +378,17 @@ void Read_PDB_File( const char * const pdb_file, reax_system* system, control_pa
strtod( &s_y[0], &endptr ),
strtod( &s_z[0], &endptr ), &x );
- Fit_to_Periodic_Box( &(system->big_box), &x );
+ Fit_to_Periodic_Box( &system->big_box, &x );
- if ( is_Inside_Box( &(system->my_box), x ) == TRUE )
+ if ( is_Inside_Box( &system->my_box, x ) == TRUE )
{
/* store orig_id, type, name and coord info of the new atom */
- atom = &(system->my_atoms[top]);
+ atom = &system->my_atoms[top];
pdb_serial = (int) strtod( &serial[0], &endptr );
atom->orig_id = pdb_serial;
Trim_Spaces( element );
- atom->type = Get_Atom_Type( &(system->reax_param), element );
+ atom->type = Get_Atom_Type( &system->reax_param, element );
strncpy( atom->name, atom_name, MAX_ATOM_NAME_LEN );
rvec_Copy( atom->x, x );
@@ -458,8 +458,9 @@ void Read_PDB_File( const char * const pdb_file, reax_system* system, control_pa
/* PDB serials are written without regard to the order, we'll see if this
* cause trouble, if so we'll have to rethink this approach
* Also, we do not write connect lines yet. */
-void Write_PDB_File( reax_system* system, reax_list* bonds, simulation_data *data,
- control_params *control, mpi_datatypes *mpi_data, output_controls *out_control )
+void Write_PDB_File( reax_system * const system, reax_list * const bond_list,
+ simulation_data * const data, control_params * const control,
+ mpi_datatypes * const mpi_data, output_controls * const out_control )
{
int i, cnt, me, np, buffer_req, buffer_len;
//int j, connect[4];
@@ -494,7 +495,7 @@ void Write_PDB_File( reax_system* system, reax_list* bonds, simulation_data *dat
buffer[0] = 0;
/* open pdb and write header */
- if (me == MASTER_NODE)
+ if ( me == MASTER_NODE )
{
/* Writing Box information */
gamma = ACOS( (system->big_box.box[0][0] * system->big_box.box[1][0] +
@@ -570,10 +571,10 @@ void Write_PDB_File( reax_system* system, reax_list* bonds, simulation_data *dat
/*
for(i=0; i < system->N; i++) {
count = 0;
- for(j = Start_Index(i, bonds); j < End_Index(i, bonds); ++j) {
- bo = bonds->bond_list[j].bo_data.BO;
+ for(j = Start_Index(i, bond_list); j < End_Index(i, bond_list); ++j) {
+ bo = bond_list->bond_list[j].bo_data.BO;
if (bo > 0.3) {
- connect[count] = bonds->bond_list[j].nbr+1;
+ connect[count] = bond_list->bond_list[j].nbr+1;
count++;
}
}
diff --git a/PG-PuReMD/src/geo_tools.h b/PG-PuReMD/src/geo_tools.h
index 1199c78a..d6b332d3 100644
--- a/PG-PuReMD/src/geo_tools.h
+++ b/PG-PuReMD/src/geo_tools.h
@@ -115,14 +115,14 @@ COLUMNS DATA TYPE FIELD DEFINITION
extern "C" {
#endif
-void Read_Geo_File( const char * const, reax_system*, control_params*,
- simulation_data*, storage*, mpi_datatypes* );
+void Read_Geo_File( const char * const, reax_system * const, control_params * const,
+ simulation_data * const, storage * const, mpi_datatypes * const );
-void Read_PDB_File( const char * const, reax_system*, control_params*,
- simulation_data*, storage*, mpi_datatypes* );
+void Read_PDB_File( const char * const, reax_system * const, control_params * const,
+ simulation_data * const, storage * const, mpi_datatypes * const );
-void Write_PDB_File( reax_system*, reax_list*, simulation_data*,
- control_params*, mpi_datatypes*, output_controls* );
+void Write_PDB_File( reax_system * const, reax_list * const, simulation_data * const,
+ control_params * const, mpi_datatypes * const, output_controls * const );
#ifdef __cplusplus
}
diff --git a/PG-PuReMD/src/grid.c b/PG-PuReMD/src/grid.c
index dddf9884..67570f13 100644
--- a/PG-PuReMD/src/grid.c
+++ b/PG-PuReMD/src/grid.c
@@ -32,7 +32,8 @@
/* determines the exchange boundaries with nbrs in terms of gcells */
-static void Mark_Grid_Cells( reax_system* system, grid *g, ivec procs, MPI_Comm comm )
+static void Mark_Grid_Cells( reax_system * const system, grid * const g,
+ ivec procs, MPI_Comm comm )
{
int x, y, z, d;
ivec r, nbr_coord, prdc;
@@ -136,7 +137,8 @@ static void Mark_Grid_Cells( reax_system* system, grid *g, ivec procs, MPI_Comm
/* finds the closest point between two grid cells denoted by c1 and c2.
* periodic boundary conditions are taken into consideration as well. */
-static void Find_Closest_Point( grid *g, ivec c1, ivec c2, rvec closest_point )
+static void Find_Closest_Point( grid * const g, ivec c1, ivec c2,
+ rvec closest_point )
{
int i, d;
@@ -161,7 +163,7 @@ static void Find_Closest_Point( grid *g, ivec c1, ivec c2, rvec closest_point )
/* mark gcells based on the kind of nbrs we will be looking for in them */
-static void Find_Neighbor_Grid_Cells( grid *g, control_params *control )
+static void Find_Neighbor_Grid_Cells( grid * const g, control_params * const control )
{
int d, top;
ivec ci, cj, cmin, cmax, span;
@@ -234,7 +236,7 @@ static void Find_Neighbor_Grid_Cells( grid *g, control_params *control )
}
-static void Reorder_Grid_Cells( grid *g )
+static void Reorder_Grid_Cells( grid * const g )
{
int i, j, k, x, y, z, top;
ivec dblock, nblocks;
@@ -284,23 +286,19 @@ static void Reorder_Grid_Cells( grid *g )
}
-void Setup_New_Grid( reax_system* system, control_params* control,
+void Setup_New_Grid( reax_system * const system, control_params * const control,
MPI_Comm comm )
{
int d, i, j, k;
- grid *g;
- simulation_box *my_box, *my_ext_box;
- boundary_cutoff *bc;
+ grid * const g = &system->my_grid;
+ simulation_box * const my_box = &system->my_box;
+ simulation_box * const my_ext_box = &system->my_ext_box;
+ boundary_cutoff * const bc = &system->bndry_cuts;
#if defined(DEBUG_FOCUS)
fprintf( stderr, "p%d: setup new grid\n", system->my_rank );
#endif
- g = &system->my_grid;
- my_box = &system->my_box;
- my_ext_box = &system->my_ext_box;
- bc = &system->bndry_cuts;
-
/* compute number of grid cells and props in each direction */
for ( d = 0; d < 3; ++d )
{
@@ -383,22 +381,18 @@ void Setup_New_Grid( reax_system* system, control_params* control,
}
-void Update_Grid( reax_system* system, control_params* control, MPI_Comm comm )
+void Update_Grid( reax_system * const system, control_params * const control,
+ MPI_Comm comm )
{
int d, i, j, k, itr;
ivec ci, native_cells, nonb_span, bond_span;
ivec ghost_span, ghost_nonb_span, ghost_bond_span, ghost_hbond_span;
rvec cell_len, inv_len;
- grid *g;
+ grid * const g = &system->my_grid;
grid_cell *gc;
- simulation_box *my_box;
- simulation_box *my_ext_box;
- boundary_cutoff *bc;
-
- g = &( system->my_grid );
- my_box = &( system->my_box );
- my_ext_box = &( system->my_ext_box );
- bc = &(system->bndry_cuts);
+ simulation_box * const my_box = &system->my_box;
+ simulation_box * const my_ext_box = &system->my_ext_box;
+ boundary_cutoff * const bc = &system->bndry_cuts;
/* compute number of grid cells and props in each direction */
for ( d = 0; d < 3; ++d )
@@ -499,20 +493,16 @@ void Update_Grid( reax_system* system, control_params* control, MPI_Comm comm )
/* Bin my (native) atoms into grid cells */
-void Bin_My_Atoms( reax_system *system, reallocate_data *realloc )
+void Bin_My_Atoms( reax_system * const system, reallocate_data * const realloc )
{
int i, j, k, l, d, max_atoms;
ivec c;
- simulation_box *big_box, *my_box, *my_ext_box;
- grid *g;
+// simulation_box * const big_box = &system->big_box;
+ simulation_box * const my_ext_box = &system->my_ext_box;
+ simulation_box * const my_box = &system->my_box;
+ grid * const g = &system->my_grid;
grid_cell *gc;
- reax_atom *atoms;
-
- big_box = &system->big_box;
- my_ext_box = &system->my_ext_box;
- my_box = &system->my_box;
- g = &system->my_grid;
- atoms = system->my_atoms;
+ reax_atom * const atoms = system->my_atoms;
#if defined(DEBUG)
fprintf( stderr, "p%d bin_my_atoms: entered\n", system->my_rank );
@@ -527,16 +517,21 @@ void Bin_My_Atoms( reax_system *system, reallocate_data *realloc )
{
for ( d = 0; d < 3; ++d )
{
- //if( atoms[l].x[d] < big_box->min[d] )
- // atoms[l].x[d] += big_box->box_norms[d];
- //else if( atoms[l].x[d] >= big_box->max[d] )
- // atoms[l].x[d] -= big_box->box_norms[d];
+// if( atoms[l].x[d] < big_box->min[d] )
+// {
+// atoms[l].x[d] += big_box->box_norms[d];
+// }
+// else if( atoms[l].x[d] >= big_box->max[d] )
+// {
+// atoms[l].x[d] -= big_box->box_norms[d];
+// }
+
if ( atoms[l].x[d] < my_box->min[d] || atoms[l].x[d] > my_box->max[d] )
{
- fprintf( stderr, "p%d: local atom%d [%f %f %f] is out of my box!\n",
+ fprintf( stderr, "[ERROR] p%d: local atom%d [%f %f %f] is out of my box!\n",
system->my_rank, l,
atoms[l].x[0], atoms[l].x[1], atoms[l].x[2] );
- fprintf( stderr, "p%d: my_box=[%f-%f, %f-%f, %f-%f]\n",
+ fprintf( stderr, "[ERROR] p%d: my_box=[%f-%f, %f-%f, %f-%f]\n",
system->my_rank, my_box->min[0], my_box->max[0],
my_box->min[1], my_box->max[1],
my_box->min[2], my_box->max[2] );
@@ -617,7 +612,7 @@ void Bin_My_Atoms( reax_system *system, reallocate_data *realloc )
/* Reorder atoms falling into the same gcell together in the atom list */
-void Reorder_My_Atoms( reax_system *system, storage *workspace )
+void Reorder_My_Atoms( reax_system * const system, storage * const workspace )
{
int i, l, x, y, z;
int top, old_id;
@@ -626,7 +621,7 @@ void Reorder_My_Atoms( reax_system *system, storage *workspace )
reax_atom *old_atom, *new_atoms;
/* allocate storage space for est_N */
- new_atoms = (reax_atom*) smalloc( sizeof(reax_atom) * system->total_cap,
+ new_atoms = smalloc( sizeof(reax_atom) * system->total_cap,
"Reorder_My_Atoms::new_atoms" );
top = 0;
g = &system->my_grid;
@@ -699,8 +694,8 @@ void Reorder_My_Atoms( reax_system *system, storage *workspace )
/* Determine the grid cell which a boundary atom falls within */
-static void Get_Boundary_Grid_Cell( grid *g, rvec base, rvec x, grid_cell **gc,
- rvec *cur_min, rvec *cur_max, ivec gcell_cood )
+static void Get_Boundary_Grid_Cell( grid * const g, rvec base, rvec x, grid_cell ** const gc,
+ rvec * const cur_min, rvec * const cur_max, ivec gcell_cood )
{
int d;
ivec c;
@@ -745,7 +740,7 @@ static void Get_Boundary_Grid_Cell( grid *g, rvec base, rvec x, grid_cell **gc,
/* Check if the current atom position falls within the
* boundaries of a grid cell */
-int is_Within_GCell( rvec x, rvec cur_min, rvec cur_max )
+static int is_Within_Grid_Cell( rvec x, rvec cur_min, rvec cur_max )
{
int d;
@@ -762,13 +757,13 @@ int is_Within_GCell( rvec x, rvec cur_min, rvec cur_max )
/* bin my boundary atoms into grid cells */
-void Bin_Boundary_Atoms( reax_system *system )
+void Bin_Boundary_Atoms( reax_system * const system )
{
int i, start, end;
rvec base, cur_min, cur_max;
- grid *g;
+ grid * const g = &system->my_grid;
grid_cell *gc;
- reax_atom *atoms;
+ reax_atom * const atoms = system->my_atoms;
simulation_box *ext_box;
ivec gcell_cood;
@@ -777,8 +772,6 @@ void Bin_Boundary_Atoms( reax_system *system )
system->my_rank, system->n, system->N );
#endif
- g = &system->my_grid;
- atoms = system->my_atoms;
start = system->n;
end = system->N;
if ( start == end )
@@ -794,7 +787,7 @@ void Bin_Boundary_Atoms( reax_system *system )
gc->top = 1;
/* error check */
- if ( is_Within_GCell( atoms[start].x, ext_box->min, ext_box->max ) == FALSE )
+ if ( is_Within_Grid_Cell( atoms[start].x, ext_box->min, ext_box->max ) == FALSE )
{
fprintf( stderr, "[ERROR] p%d: (start):ghost atom%d [%f %f %f] is out of my (grid cell) box!\n",
system->my_rank, start,
@@ -805,7 +798,7 @@ void Bin_Boundary_Atoms( reax_system *system )
for ( i = start + 1; i < end; i++ )
{
/* error check */
- if ( is_Within_GCell( atoms[i].x, ext_box->min, ext_box->max ) == FALSE )
+ if ( is_Within_Grid_Cell( atoms[i].x, ext_box->min, ext_box->max ) == FALSE )
{
fprintf( stderr, "[ERROR] p%d: (middle) ghost atom%d [%f %f %f] is out of my (grid cell) box!\n",
system->my_rank, i,
@@ -813,7 +806,7 @@ void Bin_Boundary_Atoms( reax_system *system )
MPI_Abort( MPI_COMM_WORLD, INVALID_INPUT );
}
- if ( is_Within_GCell( atoms[i].x, cur_min, cur_max ) == TRUE )
+ if ( is_Within_Grid_Cell( atoms[i].x, cur_min, cur_max ) == TRUE )
{
++gc->top;
}
diff --git a/PG-PuReMD/src/grid.h b/PG-PuReMD/src/grid.h
index cb124da7..eccfee1f 100644
--- a/PG-PuReMD/src/grid.h
+++ b/PG-PuReMD/src/grid.h
@@ -29,15 +29,15 @@
extern "C" {
#endif
-void Setup_New_Grid( reax_system*, control_params*, MPI_Comm );
+void Setup_New_Grid( reax_system * const, control_params * const, MPI_Comm );
-void Update_Grid( reax_system*, control_params*, MPI_Comm );
+void Update_Grid( reax_system * const, control_params * const, MPI_Comm );
-void Bin_My_Atoms( reax_system*, reallocate_data* );
+void Bin_My_Atoms( reax_system * const, reallocate_data * const );
-void Reorder_My_Atoms( reax_system*, storage* );
+void Reorder_My_Atoms( reax_system * const, storage * const );
-void Bin_Boundary_Atoms( reax_system* );
+void Bin_Boundary_Atoms( reax_system * const );
#ifdef __cplusplus
}
diff --git a/PG-PuReMD/src/index_utils.h b/PG-PuReMD/src/index_utils.h
index e6551b77..f4fb3fb4 100644
--- a/PG-PuReMD/src/index_utils.h
+++ b/PG-PuReMD/src/index_utils.h
@@ -5,26 +5,30 @@
/* Indexing routine for grid cells */
-static inline CUDA_HOST_DEVICE int index_grid_3d( int i, int j, int k, grid *g )
+static inline CUDA_HOST_DEVICE int index_grid_3d( int i, int j, int k,
+ const grid * const g )
{
return (i * g->ncells[1] * g->ncells[2]) + (j * g->ncells[2]) + k;
}
/* Indexing routine for grid cells, identical in functionality to above function */
-static inline CUDA_HOST_DEVICE int index_grid_3d_v( ivec x, grid *g )
+static inline CUDA_HOST_DEVICE int index_grid_3d_v( ivec x,
+ const grid * const g )
{
return (x[0] * g->ncells[1] * g->ncells[2]) + (x[1] * g->ncells[2]) + x[2];
}
/* Indexing routine for neighbors of binned atoms within grid cells */
-static inline CUDA_HOST_DEVICE int index_grid_nbrs( int i, int j, int k, int l, grid *g )
+static inline CUDA_HOST_DEVICE int index_grid_nbrs( int i, int j, int k, int l,
+ const grid * const g )
{
return (i * g->ncells[1] * g->ncells[2] * g->max_nbrs) +
(j * g->ncells[2] * g->max_nbrs) + (k * g->max_nbrs) + l;
}
/* Indexing routine for binned atoms within grid cells */
-static inline CUDA_HOST_DEVICE int index_grid_atoms( int i, int j, int k, int l, grid *g )
+static inline CUDA_HOST_DEVICE int index_grid_atoms( int i, int j, int k, int l,
+ const grid * const g )
{
return (i * g->ncells[1] * g->ncells[2] * g->max_atoms) +
(j * g->ncells[2] * g->max_atoms) + (k * g->max_atoms) + l;
diff --git a/PG-PuReMD/src/init_md.c b/PG-PuReMD/src/init_md.c
index 5e93f254..1364ab2f 100644
--- a/PG-PuReMD/src/init_md.c
+++ b/PG-PuReMD/src/init_md.c
@@ -56,8 +56,8 @@
#if defined(PURE_REAX)
/************************ initialize system ************************/
-int Reposition_Atoms( reax_system *system, control_params *control,
- simulation_data *data, mpi_datatypes *mpi_data, char *msg )
+static int Reposition_Atoms( reax_system * const system, control_params * const control,
+ simulation_data * const data, mpi_datatypes * const mpi_data, char * const msg )
{
int i;
rvec dx;
@@ -96,7 +96,7 @@ int Reposition_Atoms( reax_system *system, control_params *control,
-void Generate_Initial_Velocities( reax_system *system, real T )
+void Generate_Initial_Velocities( reax_system * const system, real T )
{
int i;
real m, scale, norm;
@@ -126,9 +126,9 @@ void Generate_Initial_Velocities( reax_system *system, real T )
}
-void Init_System( reax_system *system, control_params *control,
- simulation_data *data, storage *workspace,
- mpi_datatypes *mpi_data )
+void Init_System( reax_system * const system, control_params * const control,
+ simulation_data * const data, storage * const workspace,
+ mpi_datatypes * const mpi_data )
{
int i;
reax_atom *atom;
@@ -228,8 +228,8 @@ void Init_System( reax_system *system, control_params *control,
/************************ initialize simulation data ************************/
-void Init_Simulation_Data( reax_system *system, control_params *control,
- simulation_data *data )
+void Init_Simulation_Data( reax_system * const system, control_params * const control,
+ simulation_data * const data )
{
Reset_Simulation_Data( data );
@@ -333,7 +333,7 @@ void Init_Simulation_Data( reax_system *system, control_params *control,
#elif defined(LAMMPS_REAX)
-void Init_System( reax_system *system )
+void Init_System( reax_system * const system )
{
system->big_box.V = 0;
system->big_box.box_norms[0] = 0;
@@ -361,8 +361,8 @@ void Init_System( reax_system *system )
}
-void Init_Simulation_Data( reax_system *system, control_params *control,
- simulation_data *data )
+void Init_Simulation_Data( reax_system * const system, control_params * const control,
+ simulation_data * const data )
{
Reset_Simulation_Data( data );
@@ -378,26 +378,25 @@ void Init_Simulation_Data( reax_system *system, control_params *control,
/************************ initialize workspace ************************/
/* Initialize Taper params */
-void Init_Taper( control_params *control, storage *workspace )
+void Init_Taper( control_params * const control, storage * const workspace )
{
real d1, d7;
- real swa, swa2, swa3;
- real swb, swb2, swb3;
-
- swa = control->nonb_low;
- swb = control->nonb_cut;
+ const real swa = control->nonb_low;
+ const real swb = control->nonb_cut;
+ real swa2, swa3;
+ real swb2, swb3;
if ( FABS( swa ) > 0.01 )
{
fprintf( stderr, "[WARNING] non-zero lower Taper-radius cutoff in force field parameters\n" );
}
- if ( swb < 0 )
+ if ( swb < 0.0 )
{
fprintf( stderr, "[ERROR] negative upper Taper-radius cutoff in force field parameters\n" );
MPI_Abort( MPI_COMM_WORLD, INVALID_INPUT );
}
- else if ( swb < 5 )
+ else if ( swb < 5.0 )
{
fprintf( stderr, "[WARNING] very low Taper-radius cutoff in force field parameters (%f)\n", swb );
}
@@ -416,13 +415,13 @@ void Init_Taper( control_params *control, storage *workspace )
workspace->Tap[3] = 140.0 * (swa3 * swb + 3.0 * swa2 * swb2 + swa * swb3 ) / d7;
workspace->Tap[2] = -210.0 * (swa3 * swb2 + swa2 * swb3) / d7;
workspace->Tap[1] = 140.0 * swa3 * swb3 / d7;
- workspace->Tap[0] = (-35.0 * swa3 * swb2 * swb2 + 21.0 * swa2 * swb3 * swb2 +
- 7.0 * swa * swb3 * swb3 + swb3 * swb3 * swb ) / d7;
+ workspace->Tap[0] = (-35.0 * swa3 * swb2 * swb2 + 21.0 * swa2 * swb3 * swb2
+ + 7.0 * swa * swb3 * swb3 + swb3 * swb3 * swb ) / d7;
}
-void Init_Workspace( reax_system *system, control_params *control,
- storage *workspace )
+void Init_Workspace( reax_system * const system, control_params * const control,
+ storage * const workspace )
{
Allocate_Workspace( system, control, workspace, system->local_cap,
system->total_cap );
@@ -440,9 +439,10 @@ void Init_Workspace( reax_system *system, control_params *control,
}
-/************** setup communication data structures **************/
-void Init_MPI_Datatypes( reax_system *system, storage *workspace,
- mpi_datatypes *mpi_data )
+/* Setup communication data structures
+ * */
+void Init_MPI_Datatypes( reax_system * const system, storage * const workspace,
+ mpi_datatypes * const mpi_data )
{
int block[11];
int i;
@@ -619,10 +619,11 @@ void Init_MPI_Datatypes( reax_system *system, storage *workspace,
}
-/********************** allocate lists *************************/
-void Init_Lists( reax_system *system, control_params *control,
- simulation_data *data, storage *workspace, reax_list **lists,
- mpi_datatypes *mpi_data )
+/* Allocate and initialize lists
+ * */
+void Init_Lists( reax_system * const system, control_params * const control,
+ simulation_data * const data, storage * const workspace, reax_list ** const lists,
+ mpi_datatypes * const mpi_data )
{
int ret;
@@ -662,10 +663,10 @@ void Init_Lists( reax_system *system, control_params *control,
#if defined(PURE_REAX)
-void Initialize( reax_system *system, control_params *control,
- simulation_data *data, storage *workspace,
- reax_list **lists, output_controls *out_control,
- mpi_datatypes *mpi_data )
+void Initialize( reax_system * const system, control_params * const control,
+ simulation_data * const data, storage * const workspace,
+ reax_list ** const lists, output_controls * const out_control,
+ mpi_datatypes * const mpi_data )
{
Init_MPI_Datatypes( system, workspace, mpi_data );
@@ -693,10 +694,10 @@ void Initialize( reax_system *system, control_params *control,
}
-void Pure_Initialize( reax_system *system, control_params *control,
- simulation_data *data, storage *workspace,
- reax_list **lists, output_controls *out_control,
- mpi_datatypes *mpi_data )
+void Pure_Initialize( reax_system * const system, control_params * const control,
+ simulation_data * const data, storage * const workspace,
+ reax_list ** const lists, output_controls * const out_control,
+ mpi_datatypes * const mpi_data )
{
Init_Simulation_Data( system, control, data );
@@ -709,10 +710,10 @@ void Pure_Initialize( reax_system *system, control_params *control,
#elif defined(LAMMPS_REAX)
-void Initialize( reax_system *system, control_params *control,
- simulation_data *data, storage *workspace,
- reax_list **lists, output_controls *out_control,
- mpi_datatypes *mpi_data )
+void Initialize( reax_system * const system, control_params * const control,
+ simulation_data * const data, storage * const workspace,
+ reax_list ** const lists, output_controls * const out_control,
+ mpi_datatypes * const mpi_data )
{
Init_System( system );
@@ -736,12 +737,10 @@ void Initialize( reax_system *system, control_params *control,
#endif
-static void Finalize_System( reax_system *system, control_params *control,
- simulation_data *data )
+static void Finalize_System( reax_system * const system, control_params * const control,
+ simulation_data * const data )
{
- reax_interaction *reax;
-
- reax = &system->reax_param;
+ reax_interaction * const reax = &system->reax_param;
Deallocate_Grid( &system->my_grid );
@@ -766,14 +765,14 @@ static void Finalize_System( reax_system *system, control_params *control,
}
-static void Finalize_Simulation_Data( reax_system *system, control_params *control,
- simulation_data *data, output_controls *out_control )
+static void Finalize_Simulation_Data( reax_system * const system, control_params * const control,
+ simulation_data * const data, output_controls * const out_control )
{
}
-static void Finalize_Workspace( reax_system *system, control_params *control,
- storage *workspace )
+static void Finalize_Workspace( reax_system * const system, control_params * const control,
+ storage * const workspace )
{
int i;
@@ -920,7 +919,7 @@ static void Finalize_Workspace( reax_system *system, control_params *control,
}
-static void Finalize_Lists( control_params *control, reax_list **lists )
+static void Finalize_Lists( control_params * const control, reax_list ** const lists )
{
Delete_List( lists[FAR_NBRS] );
Delete_List( lists[BONDS] );
@@ -937,7 +936,7 @@ static void Finalize_Lists( control_params *control, reax_list **lists )
}
-static void Finalize_MPI_Datatypes( mpi_datatypes *mpi_data )
+static void Finalize_MPI_Datatypes( mpi_datatypes * const mpi_data )
{
int ret;
@@ -959,9 +958,10 @@ static void Finalize_MPI_Datatypes( mpi_datatypes *mpi_data )
/* Deallocate top-level data structures, close file handles, etc.
*
*/
-void Finalize( reax_system *system, control_params *control,
- simulation_data *data, storage *workspace, reax_list **lists,
- output_controls *out_control, mpi_datatypes *mpi_data, const int output_enabled )
+void Finalize( reax_system * const system, control_params * const control,
+ simulation_data * const data, storage * const workspace, reax_list ** const lists,
+ output_controls * const out_control, mpi_datatypes * const mpi_data,
+ const int output_enabled )
{
if ( control->tabulate )
{
diff --git a/PG-PuReMD/src/init_md.h b/PG-PuReMD/src/init_md.h
index 8ac1abff..80f6e5a9 100644
--- a/PG-PuReMD/src/init_md.h
+++ b/PG-PuReMD/src/init_md.h
@@ -29,21 +29,21 @@
extern "C" {
#endif
-void Generate_Initial_Velocities( reax_system *, real );
+void Generate_Initial_Velocities( reax_system * const, real );
-void Init_MPI_Datatypes( reax_system *, storage *, mpi_datatypes * );
+void Init_MPI_Datatypes( reax_system * const, storage * const, mpi_datatypes * const );
-void Initialize( reax_system*, control_params*, simulation_data*,
- storage*, reax_list**, output_controls*, mpi_datatypes* );
+void Initialize( reax_system * const, control_params * const, simulation_data * const,
+ storage * const, reax_list** const, output_controls * const, mpi_datatypes * const );
-void Pure_Initialize( reax_system*, control_params*, simulation_data*,
- storage*, reax_list**, output_controls*, mpi_datatypes* );
+void Pure_Initialize( reax_system * const, control_params * const, simulation_data * const,
+ storage * const, reax_list** const, output_controls * const, mpi_datatypes * const );
-void Init_Taper( control_params *, storage * );
+void Init_Taper( control_params * const, storage * const );
-void Finalize( reax_system *, control_params *,
- simulation_data *, storage *, reax_list **,
- output_controls *, mpi_datatypes *, const int );
+void Finalize( reax_system * const, control_params * const,
+ simulation_data * const, storage * const, reax_list ** const,
+ output_controls * const, mpi_datatypes * const, const int );
#ifdef __cplusplus
}
diff --git a/PG-PuReMD/src/lin_alg.c b/PG-PuReMD/src/lin_alg.c
index 178efee2..9e66451b 100644
--- a/PG-PuReMD/src/lin_alg.c
+++ b/PG-PuReMD/src/lin_alg.c
@@ -1445,11 +1445,9 @@ int CG( const reax_system * const system, const control_params * const control,
const sparse_matrix * const H, const real * const b,
const real tol, real * const x, const int fresh_pre )
{
- fprintf(stdout,"CG working\n");
- fflush(stdout);
- int i, j;
+ int i;
real tmp, alpha, beta, b_norm;
- real sig_old, sig_new, sig0;
+ real sig_old, sig_new;
//real *d, *r, *p, *z;
/*d = workspace->d;
@@ -1457,7 +1455,6 @@ int CG( const reax_system * const system, const control_params * const control,
p = workspace->q;
z = workspace->p;*/
-
Dist( system, mpi_data, x, REAL_PTR_TYPE, MPI_DOUBLE );
Sparse_MatVec( H, x, workspace->q, system->N );
@@ -1475,17 +1472,11 @@ int CG( const reax_system * const system, const control_params * const control,
Vector_Sum( workspace->r , 1.0, b, -1.0, workspace->q, system->n );
- // preconditioner
- /*for ( j = 0; j < system->n; ++j )
- {
- workspace->d[j] = workspace->r[j] * workspace->Hdia_inv[j];
- }*/
apply_preconditioner( system, workspace, control, workspace->r, workspace->d, fresh_pre, LEFT );
apply_preconditioner( system, workspace, control, workspace->d, workspace->p, fresh_pre, RIGHT );
b_norm = Parallel_Norm( b, system->n, mpi_data->world );
sig_new = Parallel_Dot( workspace->r, workspace->d, system->n, mpi_data->world );
- sig0 = sig_new;
#if defined(CG_PERFORMANCE)
if ( system->my_rank == MASTER_NODE )
@@ -1516,11 +1507,6 @@ int CG( const reax_system * const system, const control_params * const control,
Vector_Add( x, alpha, workspace->d, system->n );
Vector_Add( workspace->r, -alpha, workspace->q, system->n );
- /* preconditioner */
- /*for ( j = 0; j < system->n; ++j )
- {
- workspace->p[j] = workspace->r[j] * workspace->Hdia_inv[j];
- }*/
apply_preconditioner( system, workspace, control, workspace->r, workspace->d, FALSE, LEFT );
apply_preconditioner( system, workspace, control, workspace->d, workspace->p, FALSE, RIGHT );
diff --git a/PG-PuReMD/src/list.c b/PG-PuReMD/src/list.c
index 5ea13dfc..316b29ae 100644
--- a/PG-PuReMD/src/list.c
+++ b/PG-PuReMD/src/list.c
@@ -30,7 +30,7 @@
#endif
-void Print_List( reax_list* list )
+void Print_List( reax_list * const list )
{
int i;
@@ -52,7 +52,7 @@ void Print_List( reax_list* list )
* type: list interaction type
* l: pointer to list to be allocated
* */
-void Make_List( int n, int max_intrs, int type, reax_list *l )
+void Make_List( int n, int max_intrs, int type, reax_list * const l )
{
if ( l->allocated == TRUE )
{
@@ -109,7 +109,7 @@ void Make_List( int n, int max_intrs, int type, reax_list *l )
}
-void Delete_List( reax_list *l )
+void Delete_List( reax_list * const l )
{
if ( l->allocated == FALSE )
{
@@ -168,7 +168,7 @@ void Delete_List( reax_list *l )
* list: pointer to list
* max_intrs: max. num. of interactions for each list element
* */
-void Init_List_Indices( reax_list *list, int *max_intrs )
+void Init_List_Indices( reax_list * const list, int * const max_intrs )
{
int i;
diff --git a/PG-PuReMD/src/list.h b/PG-PuReMD/src/list.h
index aeeffce1..9bd9c8e1 100644
--- a/PG-PuReMD/src/list.h
+++ b/PG-PuReMD/src/list.h
@@ -29,13 +29,13 @@
extern "C" {
#endif
-void Print_List( reax_list* );
+void Print_List( reax_list * const );
-void Make_List( int, int, int, reax_list* );
+void Make_List( int, int, int, reax_list * const );
-void Delete_List( reax_list* );
+void Delete_List( reax_list * const );
-void Init_List_Indices( reax_list *, int * );
+void Init_List_Indices( reax_list * const, int * const );
#ifdef __cplusplus
}
diff --git a/PG-PuReMD/src/lookup.c b/PG-PuReMD/src/lookup.c
index 03a185e5..6c5f7927 100644
--- a/PG-PuReMD/src/lookup.c
+++ b/PG-PuReMD/src/lookup.c
@@ -39,8 +39,8 @@
/* Fills solution into x. Warning: will modify c and d! */
-static void Tridiagonal_Solve( const real *a, const real *b,
- real *c, real *d, real *x, unsigned int n )
+static void Tridiagonal_Solve( const real * const a, const real * const b,
+ real * const c, real * const d, real * const x, unsigned int n )
{
int i;
real id;
@@ -64,8 +64,8 @@ static void Tridiagonal_Solve( const real *a, const real *b,
}
-static void Natural_Cubic_Spline( const real *h, const real *f,
- cubic_spline_coef *coef, unsigned int n )
+static void Natural_Cubic_Spline( const real * const h, const real * const f,
+ cubic_spline_coef * const coef, unsigned int n )
{
int i;
real *a, *b, *c, *d, *v;
@@ -128,8 +128,8 @@ static void Natural_Cubic_Spline( const real *h, const real *f,
}
-static void Complete_Cubic_Spline( const real *h, const real *f, real v0, real vlast,
- cubic_spline_coef *coef, unsigned int n )
+static void Complete_Cubic_Spline( const real * const h, const real * const f,
+ real v0, real vlast, cubic_spline_coef * const coef, unsigned int n )
{
int i;
real *a, *b, *c, *d, *v;
@@ -184,7 +184,7 @@ static void Complete_Cubic_Spline( const real *h, const real *f, real v0, real v
}
-void LR_Lookup( LR_lookup_table *t, real r, LR_data *y )
+void LR_Lookup( LR_lookup_table * const t, real r, LR_data * const y )
{
int i;
real base, dif;
@@ -211,8 +211,8 @@ void LR_Lookup( LR_lookup_table *t, real r, LR_data *y )
}
-void Init_Lookup_Tables( reax_system *system, control_params *control,
- storage *workspace, mpi_datatypes *mpi_data )
+void Init_Lookup_Tables( reax_system * const system, control_params * const control,
+ storage * const workspace, mpi_datatypes * const mpi_data )
{
int i, j, r;
int num_atom_types;
@@ -353,8 +353,8 @@ void Init_Lookup_Tables( reax_system *system, control_params *control,
}
-void Finalize_LR_Lookup_Table( reax_system *system, control_params *control,
- storage *workspace, mpi_datatypes *mpi_data )
+void Finalize_LR_Lookup_Table( reax_system * const system, control_params * const control,
+ storage * const workspace, mpi_datatypes * const mpi_data )
{
int i, j;
int num_atom_types;
@@ -398,7 +398,7 @@ void Finalize_LR_Lookup_Table( reax_system *system, control_params *control,
/*
-void copy_LR_table_to_device( reax_system *system, control_params *control, int aggregated )
+void copy_LR_table_to_device( reax_system * const system, control_params * const control, int aggregated )
{
int i, j, r;
int num_atom_types;
diff --git a/PG-PuReMD/src/lookup.h b/PG-PuReMD/src/lookup.h
index 8ce37e0f..6a602163 100644
--- a/PG-PuReMD/src/lookup.h
+++ b/PG-PuReMD/src/lookup.h
@@ -29,11 +29,11 @@
extern "C" {
#endif
-void Init_Lookup_Tables( reax_system *, control_params *, storage *,
- mpi_datatypes * );
+void Init_Lookup_Tables( reax_system * const, control_params * const, storage * const,
+ mpi_datatypes * const );
-void Finalize_LR_Lookup_Table( reax_system *, control_params *,
- storage *, mpi_datatypes * );
+void Finalize_LR_Lookup_Table( reax_system * const, control_params * const,
+ storage * const, mpi_datatypes * const );
#ifdef _cplusplus
}
diff --git a/PG-PuReMD/src/multi_body.c b/PG-PuReMD/src/multi_body.c
index e9152f08..6c3e3886 100644
--- a/PG-PuReMD/src/multi_body.c
+++ b/PG-PuReMD/src/multi_body.c
@@ -49,7 +49,6 @@ void Atom_Energy( reax_system * const system, control_params * const control,
real exp_ovun2n, exp_ovun6, exp_ovun8;
real inv_exp_ovun1, inv_exp_ovun2, inv_exp_ovun2n, inv_exp_ovun8;
real e_un, CEunder1, CEunder2, CEunder3, CEunder4;
- const real p_lp1 = system->reax_param.gp.l[15];
real p_lp2;
const real p_lp3 = system->reax_param.gp.l[5];
real p_ovun2;
@@ -59,7 +58,7 @@ void Atom_Energy( reax_system * const system, control_params * const control,
const real p_ovun6 = system->reax_param.gp.l[6];
const real p_ovun7 = system->reax_param.gp.l[8];
const real p_ovun8 = system->reax_param.gp.l[9];
- single_body_parameters *sbp_i, *sbp_j;
+ single_body_parameters *sbp_i;
two_body_parameters *twbp;
bond_data *pbond;
bond_order_data *bo_ij;
@@ -80,8 +79,8 @@ void Atom_Energy( reax_system * const system, control_params * const control,
e_lp = p_lp2 * workspace->Delta_lp[i] * inv_expvd2;
data->my_en.e_lp += e_lp;
- dElp = p_lp2 * inv_expvd2 +
- 75.0 * p_lp2 * workspace->Delta_lp[i] * expvd2 * SQR(inv_expvd2);
+ dElp = p_lp2 * inv_expvd2 + 75.0 * p_lp2 * workspace->Delta_lp[i]
+ * expvd2 * SQR(inv_expvd2);
CElp = dElp * workspace->dDelta_lp[i];
workspace->CdDelta[i] += CElp; // lp - 1st term
@@ -124,8 +123,9 @@ void Atom_Energy( reax_system * const system, control_params * const control,
e_lph = p_lp3 * SQR( vov3 - 3.0 );
data->my_en.e_lp += e_lph;
- deahu2dbo = 2.*p_lp3 * (vov3 - 3.);
- deahu2dsbo = 2.*p_lp3 * (vov3 - 3.) * (-1. - 0.16 * POW(Di, 3.));
+ deahu2dbo = 2.0 * p_lp3 * (vov3 - 3.0);
+ deahu2dsbo = 2.0 * p_lp3 * (vov3 - 3.0)
+ * (-1.0 - 0.16 * POW(Di, 3.));
bo_ij->Cdbo += deahu2dbo;
workspace->CdDelta[i] += deahu2dsbo;
@@ -170,13 +170,12 @@ void Atom_Energy( reax_system * const system, control_params * const control,
j = bond_list->bond_list[pj].nbr;
type_j = system->my_atoms[j].type;
bo_ij = &bond_list->bond_list[pj].bo_data;
- sbp_j = &system->reax_param.sbp[ type_j ];
twbp = &system->reax_param.tbp[
index_tbp(type_i, type_j, system->reax_param.num_atom_types) ];
sum_ovun1 += twbp->p_ovun1 * twbp->De_s * bo_ij->BO;
- sum_ovun2 += (workspace->Delta[j] - dfvl * workspace->Delta_lp_temp[j]) *
- ( bo_ij->BO_pi + bo_ij->BO_pi2 );
+ sum_ovun2 += (workspace->Delta[j] - dfvl * workspace->Delta_lp_temp[j])
+ * ( bo_ij->BO_pi + bo_ij->BO_pi2 );
}
exp_ovun1 = p_ovun3 * EXP( p_ovun4 * sum_ovun2 );
diff --git a/PG-PuReMD/src/puremd.c b/PG-PuReMD/src/puremd.c
index ecbcdd1c..326d4296 100644
--- a/PG-PuReMD/src/puremd.c
+++ b/PG-PuReMD/src/puremd.c
@@ -54,7 +54,9 @@
#endif
#endif
+
/* CUDA-specific globals */
+#if defined(HAVE_CUDA)
reax_list **dev_lists;
storage *dev_workspace;
void *scratch;
@@ -62,12 +64,13 @@ void *host_scratch;
int BLOCKS, BLOCKS_POW_2, BLOCK_SIZE;
int BLOCKS_N, BLOCKS_POW_2_N;
int MATVEC_BLOCKS;
+#endif
static void Read_Config_Files( const char * const geo_file, const char * const ffield_file,
const char * const control_file,
- reax_system *system, control_params *control, simulation_data *data,
- storage *workspace, output_controls *out_control, mpi_datatypes *mpi_data )
+ reax_system * const system, control_params * const control, simulation_data * const data,
+ storage * const workspace, output_controls * const out_control, mpi_datatypes * const mpi_data )
{
Read_Force_Field_File( ffield_file, &system->reax_param, system, control );
@@ -99,9 +102,9 @@ static void Read_Config_Files( const char * const geo_file, const char * const f
}
-static void Post_Evolve( reax_system* system, control_params* control,
- simulation_data* data, storage* workspace, reax_list** lists,
- output_controls *out_control, mpi_datatypes *mpi_data )
+static void Post_Evolve( reax_system * const system, control_params * const control,
+ simulation_data * const data, storage * const workspace, reax_list ** const lists,
+ output_controls * const out_control, mpi_datatypes * const mpi_data )
{
int i;
rvec diff, cross;
@@ -139,9 +142,9 @@ static void Post_Evolve( reax_system* system, control_params* control,
#ifdef HAVE_CUDA
-static void Cuda_Post_Evolve( reax_system* system, control_params* control,
- simulation_data* data, storage* workspace, reax_list** lists,
- output_controls *out_control, mpi_datatypes *mpi_data )
+static void Cuda_Post_Evolve( reax_system * const system, control_params * const control,
+ simulation_data * const data, storage * const workspace, reax_list ** const lists,
+ output_controls * const out_control, mpi_datatypes * const mpi_data )
{
/* remove trans & rot velocity of the center of mass from system */
if ( control->ensemble != NVE && control->remove_CoM_vel &&
@@ -261,11 +264,12 @@ int simulate( const void * const handle )
output_controls *out_control;
mpi_datatypes *mpi_data;
puremd_handle *pmd_handle;
- real t_start = 0, t_elapsed;
+ real t_start, t_elapsed;
#if defined(DEBUG)
real t_begin, t_end;
#endif
+ t_start = 0;
ret = PUREMD_FAILURE;
if ( handle != NULL )
diff --git a/PG-PuReMD/src/reset_tools.c b/PG-PuReMD/src/reset_tools.c
index c2fb5a05..ca16d01d 100644
--- a/PG-PuReMD/src/reset_tools.c
+++ b/PG-PuReMD/src/reset_tools.c
@@ -36,7 +36,7 @@
#include "index_utils.h"
-void Reset_Atoms( reax_system* system, control_params *control )
+static void Reset_Atoms( reax_system * const system, control_params * const control )
{
int i;
reax_atom *atom;
@@ -61,7 +61,7 @@ void Reset_Atoms( reax_system* system, control_params *control )
}
-void Reset_Energies( energy_data *en )
+static void Reset_Energies( energy_data * const en )
{
en->e_bond = 0.0;
en->e_ov = 0.0;
@@ -82,13 +82,13 @@ void Reset_Energies( energy_data *en )
}
-void Reset_Temperatures( simulation_data *data )
+static void Reset_Temperatures( simulation_data * const data )
{
data->therm.T = 0.0;
}
-void Reset_Pressures( simulation_data *data )
+void Reset_Pressures( simulation_data * const data )
{
data->flex_bar.P_scalar = 0.0;
rtensor_MakeZero( data->flex_bar.P );
@@ -100,7 +100,7 @@ void Reset_Pressures( simulation_data *data )
}
-void Reset_Simulation_Data( simulation_data* data )
+void Reset_Simulation_Data( simulation_data * const data )
{
Reset_Energies( &data->my_en );
Reset_Energies( &data->sys_en );
@@ -109,7 +109,7 @@ void Reset_Simulation_Data( simulation_data* data )
}
-void Reset_Timing( reax_timing *rt )
+void Reset_Timing( reax_timing * const rt )
{
rt->total = Get_Time( );
rt->comm = 0.0;
@@ -124,7 +124,7 @@ void Reset_Timing( reax_timing *rt )
#ifdef TEST_FORCES
-void Reset_Test_Forces( reax_system *system, storage *workspace )
+void Reset_Test_Forces( reax_system * const system, storage * const workspace )
{
memset( workspace->f_ele, 0, system->total_cap * sizeof(rvec) );
memset( workspace->f_vdw, 0, system->total_cap * sizeof(rvec) );
@@ -143,7 +143,7 @@ void Reset_Test_Forces( reax_system *system, storage *workspace )
#endif
-void Reset_Workspace( reax_system *system, storage *workspace )
+void Reset_Workspace( reax_system * const system, storage * const workspace )
{
memset( workspace->total_bond_order, 0, system->total_cap * sizeof( real ) );
memset( workspace->dDeltap_self, 0, system->total_cap * sizeof( rvec ) );
@@ -157,7 +157,7 @@ void Reset_Workspace( reax_system *system, storage *workspace )
}
-void Reset_Grid( grid *g )
+void Reset_Grid( grid * const g )
{
int i, j, k;
@@ -176,7 +176,7 @@ void Reset_Grid( grid *g )
}
-void Reset_Out_Buffers( mpi_out_data *out_buf, int n )
+void Reset_Out_Buffers( mpi_out_data * const out_buf, int n )
{
int i;
@@ -187,14 +187,12 @@ void Reset_Out_Buffers( mpi_out_data *out_buf, int n )
}
-void Reset_Lists( reax_system *system, control_params *control,
- storage *workspace, reax_list **lists )
+void Reset_Lists( reax_system * const system, control_params * const control,
+ storage * const workspace, reax_list ** const lists )
{
int i;
- reax_list *bond_list, *hbond_list;
-
- bond_list = lists[BONDS];
- hbond_list = lists[HBONDS];
+ reax_list * const bond_list = lists[BONDS];
+ reax_list * const hbond_list = lists[HBONDS];
if ( system->N > 0 )
{
@@ -220,8 +218,9 @@ void Reset_Lists( reax_system *system, control_params *control,
}
-void Reset( reax_system *system, control_params *control,
- simulation_data *data, storage *workspace, reax_list **lists )
+void Reset( reax_system * const system, control_params * const control,
+ simulation_data * const data, storage * const workspace,
+ reax_list ** const lists )
{
Reset_Atoms( system, control );
diff --git a/PG-PuReMD/src/reset_tools.h b/PG-PuReMD/src/reset_tools.h
index 7ce7e124..f47ce364 100644
--- a/PG-PuReMD/src/reset_tools.h
+++ b/PG-PuReMD/src/reset_tools.h
@@ -29,24 +29,26 @@
extern "C" {
#endif
-void Reset_Pressures( simulation_data* );
+void Reset_Pressures( simulation_data * const );
-void Reset_Simulation_Data( simulation_data* );
+void Reset_Simulation_Data( simulation_data * const );
-void Reset_Timing( reax_timing* );
+void Reset_Timing( reax_timing * const );
-void Reset_Workspace( reax_system*, storage* );
+void Reset_Workspace( reax_system * const, storage * const );
-void Reset_Lists( reax_system*, control_params*, storage*, reax_list** );
+void Reset_Lists( reax_system * const, control_params * const,
+ storage * const, reax_list** const );
-void Reset_Grid( grid* );
+void Reset_Grid( grid * const );
-void Reset_Out_Buffers( mpi_out_data*, int );
+void Reset_Out_Buffers( mpi_out_data * const, int );
-void Reset( reax_system*, control_params*, simulation_data*, storage*, reax_list** );
+void Reset( reax_system * const, control_params * const,
+ simulation_data * const, storage * const, reax_list** const );
#ifdef TEST_FORCES
-void Reset_Test_Forces( reax_system*, storage* );
+void Reset_Test_Forces( reax_system * const, storage * const );
#endif
#ifdef __cplusplus
diff --git a/PG-PuReMD/src/torsion_angles.c b/PG-PuReMD/src/torsion_angles.c
index 4e3c90ca..a6722b21 100644
--- a/PG-PuReMD/src/torsion_angles.c
+++ b/PG-PuReMD/src/torsion_angles.c
@@ -152,7 +152,7 @@ void Torsion_Angles( reax_system * const system, control_params * const control,
{
int i, j, k, l, pi, pj, pk, pl, pij, plk;
int type_i, type_j, type_k, type_l;
- int start_j, end_j, start_k, end_k;
+ int start_j, end_j;
int start_pj, end_pj, start_pk, end_pk;
int num_frb_intrs;
real Delta_j, Delta_k;
@@ -213,8 +213,6 @@ void Torsion_Angles( reax_system * const system, control_params * const control,
&& bo_jk->BO > control->thb_cut
&& Num_Entries(pk, thb_list) )
{
- start_k = Start_Index( k, bond_list );
- end_k = End_Index( k, bond_list );
/* pj points to j on k's list */
pj = pbond_jk->sym_index;
diff --git a/PG-PuReMD/src/valence_angles.c b/PG-PuReMD/src/valence_angles.c
index 5fa133a9..7fd5f497 100644
--- a/PG-PuReMD/src/valence_angles.c
+++ b/PG-PuReMD/src/valence_angles.c
@@ -113,7 +113,6 @@ void Valence_Angles( reax_system * const system, control_params * const control,
real Cf7ij, Cf7jk, Cf8j, Cf9j;
real f7_ij, f7_jk, f8_Dj, f9_Dj;
real Ctheta_0, theta_0, theta_00, theta, cos_theta, sin_theta;
- real r_ij, r_jk;
real BOA_ij, BOA_jk;
rvec force, ext_press;
//rtensor temp_rtensor, total_rtensor;
@@ -198,7 +197,6 @@ void Valence_Angles( reax_system * const system, control_params * const control,
( j < system->n || pbond_ij->nbr < system->n ) )
{
i = pbond_ij->nbr;
- r_ij = pbond_ij->d;
type_i = system->my_atoms[i].type;
/* first copy 3-body intrs from previously computed ones where i>k.
@@ -261,7 +259,6 @@ void Valence_Angles( reax_system * const system, control_params * const control,
if ( j < system->n && BOA_jk > 0.0
&& (bo_ij->BO * bo_jk->BO > SQR(control->thb_cut)) )
{
- r_jk = pbond_jk->d;
thbh = &system->reax_param.thbp[
(type_i * system->reax_param.num_atom_types * system->reax_param.num_atom_types)
+ (type_j * system->reax_param.num_atom_types ) + type_k ];
--
GitLab