diff --git a/sPuReMD/src/allocate.c b/sPuReMD/src/allocate.c
index 210b5b0fe20d1c2fdc1339e75d147841db138f99..80354102ef2f54cfd5a9ffeccd47782f1b8aafa6 100644
--- a/sPuReMD/src/allocate.c
+++ b/sPuReMD/src/allocate.c
@@ -59,15 +59,11 @@ void Reallocate_Neighbor_List( list *far_nbrs, int n, int num_intrs )
{
Delete_List( TYP_FAR_NEIGHBOR, far_nbrs );
- if (!Make_List( n, num_intrs, TYP_FAR_NEIGHBOR, far_nbrs ))
- {
- fprintf(stderr, "Problem in initializing far nbrs list. Terminating!\n");
- exit( CANNOT_INITIALIZE );
- }
+ Make_List( n, num_intrs, TYP_FAR_NEIGHBOR, far_nbrs );
#if defined(DEBUG_FOCUS)
fprintf( stderr, "num_far = %d, far_nbrs = %d -> reallocating!\n",
- num_intrs, far_nbrs->num_intrs );
+ num_intrs, far_nbrs->total_intrs );
fprintf( stderr, "memory allocated: far_nbrs = %ldMB\n",
num_intrs * sizeof(far_neighbor_data) / (1024 * 1024) );
#endif
@@ -79,7 +75,8 @@ int Allocate_Matrix( sparse_matrix **pH, int n, int m )
{
sparse_matrix *H;
- if ( (*pH = (sparse_matrix*) malloc(sizeof(sparse_matrix))) == NULL )
+ if ( (*pH = (sparse_matrix*) smalloc( sizeof(sparse_matrix),
+ "Allocate_Matrix::pH" )) == NULL )
{
return FAILURE;
}
@@ -88,9 +85,12 @@ int Allocate_Matrix( sparse_matrix **pH, int n, int m )
H->n = n;
H->m = m;
- if ( (H->start = (unsigned int*) malloc(sizeof(unsigned int) * (n + 1))) == NULL
- || (H->j = (unsigned int*) malloc(sizeof(unsigned int) * m)) == NULL
- || (H->val = (real*) malloc(sizeof(real) * m)) == NULL )
+ if ( (H->start = (unsigned int*) smalloc( sizeof(unsigned int) * (n + 1),
+ "Allocate_Matrix::H->start" )) == NULL
+ || (H->j = (unsigned int*) smalloc( sizeof(unsigned int) * m,
+ "Allocate_Matrix::H->j" )) == NULL
+ || (H->val = (real*) smalloc( sizeof(real) * m,
+ "Allocate_Matrix::H->val" )) == NULL )
{
return FAILURE;
}
@@ -143,11 +143,7 @@ int Allocate_HBond_List( int n, int num_h, int *h_index, int *hb_top,
}
num_hbonds = hb_top[n - 1];
- if ( !Make_List(num_h, num_hbonds, TYP_HBOND, hbonds ) )
- {
- fprintf( stderr, "not enough space for hbonds list. terminating!\n" );
- exit( CANNOT_INITIALIZE );
- }
+ Make_List( num_h, num_hbonds, TYP_HBOND, hbonds );
for ( i = 0; i < n; ++i )
{
@@ -181,7 +177,7 @@ int Reallocate_HBonds_List( int n, int num_h, int *h_index, list *hbonds )
fprintf( stderr, "reallocating hbonds\n" );
#endif
- hb_top = (int *) calloc( n, sizeof(int) );
+ hb_top = (int *) scalloc( n, sizeof(int), "Reallocate_HBonds_List::hb_top" );
for ( i = 0; i < n; ++i )
{
if ( h_index[i] >= 0 )
@@ -212,11 +208,7 @@ int Allocate_Bond_List( int n, int *bond_top, list *bonds )
}
num_bonds = bond_top[n - 1];
- if ( !Make_List(n, num_bonds, TYP_BOND, bonds ) )
- {
- fprintf( stderr, "not enough space for bonds list. terminating!\n" );
- exit( CANNOT_INITIALIZE );
- }
+ Make_List( n, num_bonds, TYP_BOND, bonds );
Set_Start_Index( 0, 0, bonds );
Set_End_Index( 0, 0, bonds );
@@ -244,7 +236,7 @@ int Reallocate_Bonds_List( int n, list *bonds, int *num_bonds, int *est_3body )
fprintf( stderr, "reallocating bonds\n" );
#endif
- bond_top = (int *) calloc( n, sizeof(int) );
+ bond_top = (int *) scalloc( n, sizeof(int), "Reallocate_Bonds_List::hb_top" );
*est_3body = 0;
for ( i = 0; i < n; ++i )
@@ -313,15 +305,11 @@ void Reallocate( reax_system *system, static_storage *workspace, list **lists,
Delete_List( TYP_THREE_BODY, (*lists) + THREE_BODIES );
if ( num_bonds == -1 )
- num_bonds = ((*lists) + BONDS)->num_intrs;
+ num_bonds = ((*lists) + BONDS)->total_intrs;
realloc->num_3body *= SAFE_ZONE;
- if ( !Make_List( num_bonds, realloc->num_3body,
- TYP_THREE_BODY, (*lists) + THREE_BODIES ) )
- {
- fprintf( stderr, "Problem in initializing angles list. Terminating!\n" );
- exit( CANNOT_INITIALIZE );
- }
+ Make_List( num_bonds, realloc->num_3body,
+ TYP_THREE_BODY, (*lists) + THREE_BODIES );
realloc->num_3body = -1;
#if defined(DEBUG_FOCUS)
fprintf( stderr, "reallocating 3 bodies\n" );
@@ -348,7 +336,8 @@ void Reallocate( reax_system *system, static_storage *workspace, list **lists,
// reallocate g->atoms
sfree( g->atoms[i][j][k], "Reallocate::g->atoms[i][j][k]" );
g->atoms[i][j][k] = (int*)
- calloc(workspace->realloc.gcell_atoms, sizeof(int));
+ scalloc( workspace->realloc.gcell_atoms, sizeof(int),
+ "Reallocate::g->atoms[i][j][k]" );
}
}
}
diff --git a/sPuReMD/src/analyze.c b/sPuReMD/src/analyze.c
index 4bb9cf4b6cbf850d0b2ba4abb25aa5b16e18a8df..d3582fea995288eb94560cbb5bb2be1b38abc53c 100644
--- a/sPuReMD/src/analyze.c
+++ b/sPuReMD/src/analyze.c
@@ -23,6 +23,7 @@
#include "box.h"
#include "list.h"
+#include "tool_box.h"
#include "vector.h"
@@ -848,14 +849,17 @@ void Calculate_Density_3DMesh( reax_system *system, simulation_data *data,
mesh_dims[0], mesh_dims[1], mesh_dims[2] );
/* allocate counter for each mesh cell */
- cell_counter = (int ***) calloc( mesh_dims[0], sizeof(int) );
+ cell_counter = (int ***) scalloc( mesh_dims[0], sizeof(int),
+ "Calculate_Density_3DMesh::cell_counter" );
for ( i = 0; i < mesh_dims[0]; ++i )
{
- cell_counter[i] = (int **) calloc( mesh_dims[1], sizeof(int) );
+ cell_counter[i] = (int **) scalloc( mesh_dims[1], sizeof(int),
+ "Calculate_Density_3DMesh::cell_counter[i]" );
for ( j = 0; j < mesh_dims[1]; ++j )
- cell_counter[i][j] = (int *) calloc( mesh_dims[2], sizeof(int) );
+ cell_counter[i][j] = (int *) scalloc( mesh_dims[2], sizeof(int),
+ "Calculate_Density_3DMesh::cell_counter[i][j]" );
}
@@ -902,7 +906,8 @@ void Calculate_Density_Slice( reax_system *system, simulation_data *data,
/* allocate counter */
num_slices = system->box.box_norms[2] / slice_thickness + 1.;
- slice_occ = (int*) calloc( num_slices, sizeof(int) );
+ slice_occ = (int*) scalloc( num_slices, sizeof(int),
+ "Calculate_Density_Slice::slice_occ" );
/* distribute atoms to slices */
for ( i = 0; i < system->N; ++i )
diff --git a/sPuReMD/src/control.c b/sPuReMD/src/control.c
index fc7a36a0cc2a637f4651c7dcc51849f900423cfe..127543e0b473e7cf7dacdb1af2d7856331a0da98 100644
--- a/sPuReMD/src/control.c
+++ b/sPuReMD/src/control.c
@@ -30,7 +30,7 @@ char Read_Control_File( FILE* fp, reax_system *system, control_params* control,
output_controls *out_control )
{
char *s, **tmp;
- int c, i;
+ int i;
real val;
int ival;
@@ -136,17 +136,18 @@ char Read_Control_File( FILE* fp, reax_system *system, control_params* control,
control->restrict_type = 0;
/* memory allocations */
- s = (char*) malloc(sizeof(char) * MAX_LINE);
- tmp = (char**) malloc(sizeof(char*)*MAX_TOKENS);
+ s = (char*) smalloc( sizeof(char) * MAX_LINE, "Read_Control_File::s" );
+ tmp = (char**) smalloc( sizeof(char*) * MAX_TOKENS, "Read_Control_File::tmp" );
for ( i = 0; i < MAX_TOKENS; i++ )
{
- tmp[i] = (char*) malloc(sizeof(char) * MAX_LINE);
+ tmp[i] = (char*) smalloc( sizeof(char) * MAX_LINE,
+ "Read_Control_File::tmp[i]" );
}
/* read control parameters file */
while (fgets(s, MAX_LINE, fp))
{
- c = Tokenize(s, &tmp);
+ Tokenize(s, &tmp);
if ( strcmp(tmp[0], "simulation_name") == 0 )
{
diff --git a/sPuReMD/src/ffield.c b/sPuReMD/src/ffield.c
index 226faa403968cacba3711d3376f7bc0c4c4b08a7..56d8a310439b94aaace6c11c2ab55d3f64bc09cd 100644
--- a/sPuReMD/src/ffield.c
+++ b/sPuReMD/src/ffield.c
@@ -30,7 +30,7 @@ char Read_Force_Field( FILE* fp, reax_interaction* reax )
char *s;
char **tmp;
char ****tor_flag;
- int c, i, j, k, l, m, n, o, p, cnt;
+ int i, j, k, l, m, n, o, p, cnt;
real val;
s = (char*) malloc( sizeof(char) * MAX_LINE );
@@ -46,7 +46,7 @@ char Read_Force_Field( FILE* fp, reax_interaction* reax )
/* line 2 is number of global parameters */
fgets( s, MAX_LINE, fp );
- c = Tokenize( s, &tmp );
+ Tokenize( s, &tmp );
/* reading the number of global parameters */
n = atoi(tmp[0]);
@@ -63,7 +63,7 @@ char Read_Force_Field( FILE* fp, reax_interaction* reax )
for (i = 0; i < n; i++)
{
fgets(s, MAX_LINE, fp);
- c = Tokenize(s, &tmp);
+ Tokenize(s, &tmp);
val = (real) atof(tmp[0]);
@@ -73,7 +73,7 @@ char Read_Force_Field( FILE* fp, reax_interaction* reax )
/* next line is number of atom types and some comments */
fgets( s, MAX_LINE, fp );
- c = Tokenize( s, &tmp );
+ Tokenize( s, &tmp );
reax->num_atom_types = atoi(tmp[0]);
@@ -144,7 +144,7 @@ char Read_Force_Field( FILE* fp, reax_interaction* reax )
{
/* line one */
fgets( s, MAX_LINE, fp );
- c = Tokenize( s, &tmp );
+ Tokenize( s, &tmp );
for ( j = 0; j < strlen( tmp[0] ); ++j )
reax->sbp[i].name[j] = toupper( tmp[0][j] );
@@ -169,7 +169,7 @@ char Read_Force_Field( FILE* fp, reax_interaction* reax )
/* line two */
fgets( s, MAX_LINE, fp );
- c = Tokenize( s, &tmp );
+ Tokenize( s, &tmp );
val = atof(tmp[0]);
reax->sbp[i].alpha = val;
@@ -193,7 +193,7 @@ char Read_Force_Field( FILE* fp, reax_interaction* reax )
/* line 3 */
fgets( s, MAX_LINE, fp );
- c = Tokenize( s, &tmp );
+ Tokenize( s, &tmp );
val = atof(tmp[0]);
reax->sbp[i].r_pi_pi = val;
@@ -212,7 +212,7 @@ char Read_Force_Field( FILE* fp, reax_interaction* reax )
/* line 4 */
fgets( s, MAX_LINE, fp );
- c = Tokenize( s, &tmp );
+ Tokenize( s, &tmp );
val = atof(tmp[0]);
reax->sbp[i].p_ovun2 = val;
@@ -305,7 +305,7 @@ char Read_Force_Field( FILE* fp, reax_interaction* reax )
/* next line is number of two body combination and some comments */
fgets(s, MAX_LINE, fp);
- c = Tokenize(s, &tmp);
+ Tokenize(s, &tmp);
l = atoi(tmp[0]);
/* a line of comments */
@@ -315,7 +315,7 @@ char Read_Force_Field( FILE* fp, reax_interaction* reax )
{
/* line 1 */
fgets(s, MAX_LINE, fp);
- c = Tokenize(s, &tmp);
+ Tokenize(s, &tmp);
j = atoi(tmp[0]) - 1;
k = atoi(tmp[1]) - 1;
@@ -351,7 +351,7 @@ char Read_Force_Field( FILE* fp, reax_interaction* reax )
/* line 2 */
fgets(s, MAX_LINE, fp);
- c = Tokenize(s, &tmp);
+ Tokenize(s, &tmp);
val = atof(tmp[0]);
reax->tbp[j][k].p_be2 = val;
@@ -486,13 +486,13 @@ char Read_Force_Field( FILE* fp, reax_interaction* reax )
/* these are two body offdiagonal terms that are different from the
combination rules defined above */
fgets(s, MAX_LINE, fp);
- c = Tokenize(s, &tmp);
+ Tokenize(s, &tmp);
l = atoi(tmp[0]);
for (i = 0; i < l; i++)
{
fgets(s, MAX_LINE, fp);
- c = Tokenize(s, &tmp);
+ Tokenize(s, &tmp);
j = atoi(tmp[0]) - 1;
k = atoi(tmp[1]) - 1;
@@ -553,13 +553,13 @@ char Read_Force_Field( FILE* fp, reax_interaction* reax )
/* next line is number of 3-body params and some comments */
fgets( s, MAX_LINE, fp );
- c = Tokenize( s, &tmp );
+ Tokenize( s, &tmp );
l = atoi( tmp[0] );
for ( i = 0; i < l; i++ )
{
fgets(s, MAX_LINE, fp);
- c = Tokenize(s, &tmp);
+ Tokenize(s, &tmp);
j = atoi(tmp[0]) - 1;
k = atoi(tmp[1]) - 1;
@@ -623,13 +623,13 @@ char Read_Force_Field( FILE* fp, reax_interaction* reax )
/* next line is number of 4-body params and some comments */
fgets( s, MAX_LINE, fp );
- c = Tokenize( s, &tmp );
+ Tokenize( s, &tmp );
l = atoi( tmp[0] );
for ( i = 0; i < l; i++ )
{
fgets( s, MAX_LINE, fp );
- c = Tokenize( s, &tmp );
+ Tokenize( s, &tmp );
j = atoi(tmp[0]) - 1;
k = atoi(tmp[1]) - 1;
@@ -710,13 +710,13 @@ char Read_Force_Field( FILE* fp, reax_interaction* reax )
/* next line is number of hydrogen bond params and some comments */
fgets( s, MAX_LINE, fp );
- c = Tokenize( s, &tmp );
+ Tokenize( s, &tmp );
l = atoi( tmp[0] );
for ( i = 0; i < l; i++ )
{
fgets( s, MAX_LINE, fp );
- c = Tokenize( s, &tmp );
+ Tokenize( s, &tmp );
j = atoi(tmp[0]) - 1;
k = atoi(tmp[1]) - 1;
diff --git a/sPuReMD/src/forces.c b/sPuReMD/src/forces.c
index 0d9ae4000d0246be0e1485ad3c022ab60d028d8f..9eff53b1eb4c22aa4eb7f96c2bb5253d333ce5ca 100644
--- a/sPuReMD/src/forces.c
+++ b/sPuReMD/src/forces.c
@@ -274,14 +274,14 @@ void Validate_Lists( static_storage *workspace, list **lists, int step, int n,
exit( INSUFFICIENT_MEMORY );
}
- if ( End_Index(i, bonds) >= bonds->num_intrs - 2 )
+ if ( End_Index(i, bonds) >= bonds->total_intrs - 2 )
{
workspace->realloc.bonds = 1;
- if ( End_Index(i, bonds) > bonds->num_intrs )
+ if ( End_Index(i, bonds) > bonds->total_intrs )
{
fprintf( stderr, "step%d-bondchk failed: i=%d end(i)=%d bond_end=%d\n",
- step, flag, End_Index(i, bonds), bonds->num_intrs );
+ step, flag, End_Index(i, bonds), bonds->total_intrs );
exit( INSUFFICIENT_MEMORY );
}
}
@@ -313,14 +313,14 @@ void Validate_Lists( static_storage *workspace, list **lists, int step, int n,
}
if ( Num_Entries(i, hbonds) >=
- (hbonds->num_intrs - Start_Index(i, hbonds)) * DANGER_ZONE )
+ (hbonds->total_intrs - Start_Index(i, hbonds)) * DANGER_ZONE )
{
workspace->realloc.hbonds = 1;
- if ( End_Index(i, hbonds) > hbonds->num_intrs )
+ if ( End_Index(i, hbonds) > hbonds->total_intrs )
{
fprintf( stderr, "step%d-hbondchk failed: i=%d end(i)=%d hbondend=%d\n",
- step, flag, End_Index(i, hbonds), hbonds->num_intrs );
+ step, flag, End_Index(i, hbonds), hbonds->total_intrs );
exit( INSUFFICIENT_MEMORY );
}
}
diff --git a/sPuReMD/src/init_md.c b/sPuReMD/src/init_md.c
index 44e8637180fe88bc56ded8d992844b376b7eda68..850c044878e044045acd5a67b1721d5301924254 100644
--- a/sPuReMD/src/init_md.c
+++ b/sPuReMD/src/init_md.c
@@ -517,11 +517,7 @@ void Init_Lists( reax_system *system, control_params *control,
num_nbrs = Estimate_NumNeighbors( system, control, workspace, lists );
- if ( !Make_List(system->N, num_nbrs, TYP_FAR_NEIGHBOR, (*lists) + FAR_NBRS) )
- {
- fprintf(stderr, "Problem in initializing far nbrs list. Terminating!\n");
- exit( CANNOT_INITIALIZE );
- }
+ Make_List( system->N, num_nbrs, TYP_FAR_NEIGHBOR, (*lists) + FAR_NBRS );
#if defined(DEBUG_FOCUS)
fprintf( stderr, "memory allocated: far_nbrs = %ldMB\n",
@@ -612,11 +608,7 @@ void Init_Lists( reax_system *system, control_params *control,
#endif
/* 3bodies list */
- if (!Make_List(num_bonds, num_3body, TYP_THREE_BODY, (*lists) + THREE_BODIES))
- {
- fprintf( stderr, "Problem in initializing angles list. Terminating!\n" );
- exit( CANNOT_INITIALIZE );
- }
+ Make_List( num_bonds, num_3body, TYP_THREE_BODY, (*lists) + THREE_BODIES );
#if defined(DEBUG_FOCUS)
fprintf( stderr, "estimated storage - num_3body: %d\n", num_3body );
@@ -625,17 +617,9 @@ void Init_Lists( reax_system *system, control_params *control,
#endif
#ifdef TEST_FORCES
- if (!Make_List( system->N, num_bonds * 8, TYP_DDELTA, (*lists) + DDELTA ))
- {
- fprintf( stderr, "Problem in initializing dDelta list. Terminating!\n" );
- exit( CANNOT_INITIALIZE );
- }
+ Make_List( system->N, num_bonds * 8, TYP_DDELTA, (*lists) + DDELTA );
- if ( !Make_List( num_bonds, num_bonds * MAX_BONDS * 3, TYP_DBO, (*lists) + DBO ) )
- {
- fprintf( stderr, "Problem in initializing dBO list. Terminating!\n" );
- exit( CANNOT_INITIALIZE );
- }
+ Make_List( num_bonds, num_bonds * MAX_BONDS * 3, TYP_DBO, (*lists) + DBO );
#endif
sfree( hb_top, "Init_Lists::hb_top" );
diff --git a/sPuReMD/src/list.c b/sPuReMD/src/list.c
index de5b0438610721380039d412d1f9b927260dec15..351ee8894392cf465e3305d57a9e30711761377c 100644
--- a/sPuReMD/src/list.c
+++ b/sPuReMD/src/list.c
@@ -24,108 +24,68 @@
#include "tool_box.h"
-char Make_List( int n, int num_intrs, int type, list* l )
+void Make_List( int n, int total_intrs, int type, list* l )
{
- char success = 1;
-
l->n = n;
- l->num_intrs = num_intrs;
+ l->total_intrs = total_intrs;
- l->index = (int*) malloc( n * sizeof(int) );
- l->end_index = (int*) malloc( n * sizeof(int) );
-
- if ( l->index == NULL )
- {
- success = 0;
- }
- if ( l->end_index == NULL )
- {
- success = 0;
- }
+ l->index = (int*) smalloc( n * sizeof(int), "Make_List::l->index" );
+ l->end_index = (int*) smalloc( n * sizeof(int), "Make_List::l->end_index" );
switch ( type )
{
case TYP_VOID:
- l->select.v = (void *) malloc( l->num_intrs * sizeof(void) );
- if ( l->select.v == NULL )
- {
- success = 0;
- }
+ l->select.v = (void *) smalloc( l->total_intrs * sizeof(void),
+ "Make_List::l->select.v" );
break;
case TYP_THREE_BODY:
l->select.three_body_list = (three_body_interaction_data*)
- malloc( l->num_intrs * sizeof(three_body_interaction_data) );
- if ( l->select.three_body_list == NULL )
- {
- success = 0;
- }
+ smalloc( l->total_intrs * sizeof(three_body_interaction_data),
+ "Make_List::l->select.three_body_list" );
break;
case TYP_BOND:
l->select.bond_list = (bond_data*)
- malloc(l->num_intrs * sizeof(bond_data));
- if ( l->select.bond_list == NULL )
- {
- success = 0;
- }
+ smalloc( l->total_intrs * sizeof(bond_data),
+ "Make_List::l->select.bond_list" );
break;
case TYP_DBO:
l->select.dbo_list = (dbond_data*)
- malloc( l->num_intrs * sizeof(dbond_data) );
- if ( l->select.dbo_list == NULL )
- {
- success = 0;
- }
+ smalloc( l->total_intrs * sizeof(dbond_data),
+ "Make_List::l->select.dbo_list" );
break;
case TYP_DDELTA:
l->select.dDelta_list = (dDelta_data*)
- malloc( l->num_intrs * sizeof(dDelta_data) );
- if ( l->select.dDelta_list == NULL )
- {
- success = 0;
- }
+ smalloc( l->total_intrs * sizeof(dDelta_data),
+ "Make_List::l->select.dDelta_list" );
break;
case TYP_FAR_NEIGHBOR:
l->select.far_nbr_list = (far_neighbor_data*)
- malloc( l->num_intrs * sizeof(far_neighbor_data) );
- if (l->select.far_nbr_list == NULL)
- {
- success = 0;
- }
+ smalloc( l->total_intrs * sizeof(far_neighbor_data),
+ "Make_List::l->select.far_nbr_list" );
break;
case TYP_NEAR_NEIGHBOR:
l->select.near_nbr_list = (near_neighbor_data*)
- malloc( l->num_intrs * sizeof(near_neighbor_data) );
- if ( l->select.near_nbr_list == NULL )
- {
- success = 0;
- }
+ smalloc( l->total_intrs * sizeof(near_neighbor_data),
+ "Make_List::l->select.near_nbr_list" );
break;
case TYP_HBOND:
l->select.hbond_list = (hbond_data*)
- malloc( l->num_intrs * sizeof(hbond_data) );
- if ( l->select.hbond_list == NULL )
- {
- success = 0;
- }
+ smalloc( l->total_intrs * sizeof(hbond_data),
+ "Make_List::l->select.hbond_list" );
break;
default:
- l->select.v = (void *) malloc( l->num_intrs * sizeof(void) );
- if ( l->select.v == NULL )
- {
- success = 0;
- }
+ l->select.v = (void *) smalloc( l->total_intrs * sizeof(void),
+ "Make_List::l->select.v" );
break;
}
-
- return success;
}
diff --git a/sPuReMD/src/list.h b/sPuReMD/src/list.h
index d57ad6f39ea088fbf30ee321e4fdcd27de0a54e4..a27015824eb4d3bfe9318d88f3caf2b5aa0b29cd 100644
--- a/sPuReMD/src/list.h
+++ b/sPuReMD/src/list.h
@@ -25,7 +25,7 @@
#include "mytypes.h"
-char Make_List( int, int, int, list* );
+void Make_List( int, int, int, list* );
void Delete_List( int, list* );
diff --git a/sPuReMD/src/mytypes.h b/sPuReMD/src/mytypes.h
index 7144082d2209062cf9e55e3dda244c5751e48107..ccc54b8fd5b0e304e6a120a8edfb1acfe45a5cc8 100644
--- a/sPuReMD/src/mytypes.h
+++ b/sPuReMD/src/mytypes.h
@@ -963,19 +963,34 @@ typedef struct
/* interaction lists */
typedef struct
{
+ /* num. entries in list */
int n;
- int num_intrs;
+ /* sum of max. interactions per atom */
+ int total_intrs;
+ /* starting position of atom's interactions */
int *index;
+ /* ending position of atom's interactions */
int *end_index;
+ /* max. num. of interactions per atom */
+ int *max_intrs;
+ /* interaction list (polymorphic via union dereference) */
union
{
+ /* typeless interaction list */
void *v;
+ /* three body interaction list */
three_body_interaction_data *three_body_list;
+ /* bond interaction list */
bond_data *bond_list;
+ /* bond interaction list */
dbond_data *dbo_list;
+ /* test forces interaction list */
dDelta_data *dDelta_list;
+ /* far neighbor interaction list */
far_neighbor_data *far_nbr_list;
+ /* near neighbor interaction list */
near_neighbor_data *near_nbr_list;
+ /* hydrogen bond interaction list */
hbond_data *hbond_list;
} select;
} list;
diff --git a/sPuReMD/src/neighbors.c b/sPuReMD/src/neighbors.c
index e3a19fe29bb354835594ab863c7b8e5443d11982..5e23558fb3e67a9248aa1bd4c9960995844df546 100644
--- a/sPuReMD/src/neighbors.c
+++ b/sPuReMD/src/neighbors.c
@@ -141,13 +141,13 @@ void Generate_Neighbor_Lists( reax_system *system, control_params *control,
}
}
- if ( num_far > far_nbrs->num_intrs * DANGER_ZONE )
+ if ( num_far > far_nbrs->total_intrs * DANGER_ZONE )
{
workspace->realloc.num_far = num_far;
- if ( num_far > far_nbrs->num_intrs )
+ if ( num_far > far_nbrs->total_intrs )
{
fprintf( stderr, "step%d-ran out of space on far_nbrs: top=%d, max=%d",
- data->step, num_far, far_nbrs->num_intrs );
+ data->step, num_far, far_nbrs->total_intrs );
exit( INSUFFICIENT_MEMORY );
}
}
@@ -559,9 +559,9 @@ void Generate_Neighbor_Lists( reax_system *system, control_params *control,
// fprintf( stderr, "restrictions applied-" );
- /* verify nbrlists, count num_intrs, sort nearnbrs */
- near_nbrs->num_intrs = 0;
- far_nbrs->num_intrs = 0;
+ /* verify nbrlists, count total_intrs, sort nearnbrs */
+ near_nbrs->total_intrs = 0;
+ far_nbrs->total_intrs = 0;
for ( i = 0; i < system->N - 1; ++i )
{
if ( End_Index(i, near_nbrs) > Start_Index(i + 1, near_nbrs) )
@@ -572,7 +572,7 @@ void Generate_Neighbor_Lists( reax_system *system, control_params *control,
exit( RUNTIME_ERROR );
}
- near_nbrs->num_intrs += Num_Entries(i, near_nbrs);
+ near_nbrs->total_intrs += Num_Entries(i, near_nbrs);
if ( End_Index(i, far_nbrs) > Start_Index(i + 1, far_nbrs) )
{
@@ -582,7 +582,7 @@ void Generate_Neighbor_Lists( reax_system *system, control_params *control,
exit( RUNTIME_ERROR );
}
- far_nbrs->num_intrs += Num_Entries(i, far_nbrs);
+ far_nbrs->total_intrs += Num_Entries(i, far_nbrs);
}
for ( i = 0; i < system->N; ++i )
@@ -723,7 +723,7 @@ void Generate_Neighbor_Lists( reax_system *system, control_params *control,
}
}
- far_nbrs->num_intrs = num_far;
+ far_nbrs->total_intrs = num_far;
fprintf( stderr, "nbrs done, num_far: %d\n", num_far );
#if defined(DEBUG)
diff --git a/sPuReMD/src/three_body_interactions.c b/sPuReMD/src/three_body_interactions.c
index c54def418822fc7c2fb325ecf571b7c7940bd8b6..533c6523158db2691a15e765df4e38d14c47b3a5 100644
--- a/sPuReMD/src/three_body_interactions.c
+++ b/sPuReMD/src/three_body_interactions.c
@@ -640,13 +640,13 @@ void Three_Body_Interactions( reax_system *system, control_params *control,
data->E_Pen += e_pen_total;
data->E_Coa += e_coa_total;
- if ( num_thb_intrs >= thb_intrs->num_intrs * DANGER_ZONE )
+ if ( num_thb_intrs >= thb_intrs->total_intrs * DANGER_ZONE )
{
workspace->realloc.num_3body = num_thb_intrs;
- if ( num_thb_intrs > thb_intrs->num_intrs )
+ if ( num_thb_intrs > thb_intrs->total_intrs )
{
fprintf( stderr, "step%d-ran out of space on angle_list: top=%d, max=%d",
- data->step, num_thb_intrs, thb_intrs->num_intrs );
+ data->step, num_thb_intrs, thb_intrs->total_intrs );
exit( INSUFFICIENT_MEMORY );
}
}
diff --git a/sPuReMD/src/traj.c b/sPuReMD/src/traj.c
index f679186d5e543e3e77a5aa9b49de4cb05338c462..017d6574facc0a5bddc13f3d56bca6b527e8c00f 100644
--- a/sPuReMD/src/traj.c
+++ b/sPuReMD/src/traj.c
@@ -20,14 +20,16 @@
----------------------------------------------------------------------*/
#include "traj.h"
+
#include "list.h"
+
/************************************************/
/* CUSTOM FORMAT ROUTINES */
/************************************************/
-int Write_Custom_Header(reax_system *system, control_params *control,
- static_storage *workspace, output_controls *out_control)
+int Write_Custom_Header( reax_system *system, control_params *control,
+ static_storage *workspace, output_controls *out_control )
{
int i, header_len, control_block_len, frame_format_len;
// char buffer[2048];
@@ -96,19 +98,27 @@ int Write_Custom_Header(reax_system *system, control_params *control,
case OPT_ATOM_FULL:
sprintf( atom_format, "Atom_Full: %s", ATOM_FULL );
break;
+ default:
+ break;
}
strcat( frame_format, atom_format );
bond_format[0] = OPT_NOBOND;
if ( out_control->bond_info == OPT_BOND_BASIC )
+ {
sprintf( bond_format, "Bond_Line: %s", BOND_BASIC );
+ }
else if ( out_control->bond_info == OPT_BOND_FULL )
+ {
sprintf( bond_format, "Bond_Line_Full: %s", BOND_FULL );
+ }
strcat( frame_format, bond_format );
angle_format[0] = OPT_NOANGLE;
if ( out_control->angle_info == OPT_ANGLE_BASIC )
+ {
sprintf( angle_format, "Angle_Line: %s", ANGLE_BASIC );
+ }
strcat( frame_format, angle_format );
frame_format_len = strlen( frame_format );
@@ -138,11 +148,13 @@ int Write_Custom_Header(reax_system *system, control_params *control,
ATOM_MAPPING_LEN * system->N );
for ( i = 0; i < system->N; ++i )
+ {
out_control->write( out_control->trj, ATOM_MAPPING,
workspace->orig_id[i],
system->atoms[i].type,
system->atoms[i].name,
system->reaxprm.sbp[ system->atoms[i].type ].mass );
+ }
fflush( out_control->trj );
@@ -151,8 +163,8 @@ int Write_Custom_Header(reax_system *system, control_params *control,
int Append_Custom_Frame( reax_system *system, control_params *control,
- simulation_data *data, static_storage *workspace,
- list **lists, output_controls *out_control )
+ simulation_data *data, static_storage *workspace,
+ list **lists, output_controls *out_control )
{
int i, j, pi, pk, pk_j;
int write_atoms, write_bonds, write_angles;
@@ -164,7 +176,6 @@ int Append_Custom_Frame( reax_system *system, control_params *control,
list *thb_intrs = (*lists) + THREE_BODIES;
bond_data *bo_ij;
-
/* IMPORTANT: This whole part will go to init_trj after finalized! */
switch ( out_control->atom_format )
{
@@ -187,6 +198,7 @@ int Append_Custom_Frame( reax_system *system, control_params *control,
default:
atom_line_len = 0;
write_atoms = 0;
+ break;
}
@@ -207,10 +219,16 @@ int Append_Custom_Frame( reax_system *system, control_params *control,
if ( write_bonds )
{
for ( i = 0; i < system->N; ++i )
+ {
for ( j = Start_Index( i, bonds ); j < End_Index( i, bonds ); ++j )
+ {
if ( i < bonds->select.bond_list[j].nbr &&
bonds->select.bond_list[j].bo_data.BO >= control->bg_cut )
+ {
++num_bonds;
+ }
+ }
+ }
}
@@ -230,11 +248,15 @@ int Append_Custom_Frame( reax_system *system, control_params *control,
if ( write_angles )
{
for ( j = 0; j < system->N; ++j )
+ {
for ( pi = Start_Index(j, bonds); pi < End_Index(j, bonds); ++pi )
+ {
if ( bonds->select.bond_list[pi].bo_data.BO >= control->bg_cut )
+ {
// physical j&i bond
for ( pk = Start_Index( pi, thb_intrs );
pk < End_Index( pi, thb_intrs ); ++pk )
+ {
if ( bonds->select.bond_list[pi].nbr <
thb_intrs->select.three_body_list[pk].thb )
{
@@ -242,20 +264,30 @@ int Append_Custom_Frame( reax_system *system, control_params *control,
pk_j = thb_intrs->select.three_body_list[pk].pthb;
if ( bonds->select.bond_list[pk_j].bo_data.BO >= control->bg_cut )
+ {
// physical j&k bond
++num_thb_intrs;
+ }
}
+ }
+ }
+ }
+ }
}
-
-
/* get correct pressure */
if ( control->ensemble == NPT || control->ensemble == sNPT )
+ {
P = data->flex_bar.P_scalar;
+ }
else if ( control->ensemble == iNPT )
+ {
P = data->iso_bar.P;
- else P = 0;
-
+ }
+ else
+ {
+ P = 0;
+ }
/* calculate total frame length*/
sprintf( buffer, FRAME_GLOBALS,
@@ -365,7 +397,9 @@ int Append_Custom_Frame( reax_system *system, control_params *control,
if ( out_control->bond_info == 1 )
{
for ( i = 0; i < system->N; ++i )
+ {
for ( j = Start_Index( i, bonds ); j < End_Index( i, bonds ); ++j )
+ {
if ( i < bonds->select.bond_list[j].nbr &&
bonds->select.bond_list[j].bo_data.BO >= control->bg_cut )
{
@@ -375,11 +409,15 @@ int Append_Custom_Frame( reax_system *system, control_params *control,
workspace->orig_id[bo_ij->nbr],
bo_ij->d, bo_ij->bo_data.BO );
}
+ }
+ }
}
else if ( out_control->bond_info == 2 )
{
for ( i = 0; i < system->N; ++i )
+ {
for ( j = Start_Index( i, bonds ); j < End_Index( i, bonds ); ++j )
+ {
if ( i < bonds->select.bond_list[j].nbr &&
bonds->select.bond_list[j].bo_data.BO >= control->bg_cut )
{
@@ -390,6 +428,8 @@ int Append_Custom_Frame( reax_system *system, control_params *control,
bo_ij->d, bo_ij->bo_data.BO, bo_ij->bo_data.BO_s,
bo_ij->bo_data.BO_pi, bo_ij->bo_data.BO_pi2 );
}
+ }
+ }
}
fflush( out_control->trj );
@@ -403,11 +443,15 @@ int Append_Custom_Frame( reax_system *system, control_params *control,
num_thb_intrs * angle_line_len, num_thb_intrs );
for ( j = 0; j < system->N; ++j )
+ {
for ( pi = Start_Index(j, bonds); pi < End_Index(j, bonds); ++pi )
+ {
if ( bonds->select.bond_list[pi].bo_data.BO >= control->bg_cut )
+ {
// physical j&i bond
for ( pk = Start_Index( pi, thb_intrs );
pk < End_Index( pi, thb_intrs ); ++pk )
+ {
if ( bonds->select.bond_list[pi].nbr <
thb_intrs->select.three_body_list[pk].thb )
{
@@ -415,13 +459,19 @@ int Append_Custom_Frame( reax_system *system, control_params *control,
// get k's pointer on j's bond list
if ( bonds->select.bond_list[pk_j].bo_data.BO >= control->bg_cut )
+ {
// physical j&k bond
out_control->write( out_control->trj, ANGLE_BASIC,
workspace->orig_id[bonds->select.bond_list[pi].nbr],
workspace->orig_id[j],
workspace->orig_id[thb_intrs->select.three_body_list[pk].thb],
RAD2DEG(thb_intrs->select.three_body_list[pk].theta) );
+ }
}
+ }
+ }
+ }
+ }
}
fflush( out_control->trj );
@@ -434,30 +484,36 @@ void Read_Traj( output_controls *out_control, char *traj_name )
{
int skip_all, skip_part, n;
char size_buffer[50];
- // char read_buffer[2048];
+ gzFile trj;
- out_control->trj = gzopen( traj_name, "r" );
+ trj = gzopen( traj_name, "r" );
fprintf( stderr, "file opened!\n" );
- while ( !gzeof( out_control->trj ) )
+ while ( !gzeof( trj ) )
{
- if ( gzgets( out_control->trj, size_buffer, 50 ) == Z_NULL )
+ if ( gzgets( trj, size_buffer, 50 ) == Z_NULL )
+ {
break;
+ }
fprintf( stderr, "read line\n" );
if ( strlen( size_buffer ) >= SIZE_INFO_LEN3 )
+ {
sscanf( size_buffer, "%d %d %d", &skip_all, &skip_part, &n );
+ }
else
+ {
sscanf( size_buffer, "%d %d", &skip_all, &skip_part );
+ }
fprintf( stderr, "%d %d\n", skip_all, skip_part );
- gzseek( out_control->trj, skip_part, SEEK_CUR );
+ gzseek( trj, skip_part, SEEK_CUR );
}
- gzclose( out_control->trj );
+ gzclose( trj );
}
@@ -468,35 +524,34 @@ void Read_Traj( output_controls *out_control, char *traj_name )
/********************************************************/
int Write_xyz_Header( reax_system *system, control_params *control,
- static_storage* workspace, output_controls *out_control )
+ static_storage* workspace, output_controls *out_control )
{
fflush( out_control->trj );
- return 1;
+ return SUCCESS;
}
int Append_xyz_Frame( reax_system *system, control_params *control,
- simulation_data *data, static_storage *workspace,
- list **lists, output_controls *out_control )
+ simulation_data *data, static_storage *workspace,
+ list **lists, output_controls *out_control )
{
int i;
out_control->write( out_control->trj, "%d\n", system->N );
out_control->write( out_control->trj, "%d\t%8.3f\t%8.3f\t%8.3f\t%8.3f\n",
- data->step,
- data->E_Tot, data->E_Pot,
- E_CONV * data->E_Kin, data->therm.T );
+ data->step, data->E_Tot, data->E_Pot,
+ E_CONV * data->E_Kin, data->therm.T );
for ( i = 0; i < system->N; ++i )
+ {
out_control->write( out_control->trj, "%3s %10.5f %10.5f %10.5f\n",
- system->reaxprm.sbp[ system->atoms[i].type ].name,
- system->atoms[i].x[0],
- system->atoms[i].x[1],
- system->atoms[i].x[2] );
+ system->reaxprm.sbp[ system->atoms[i].type ].name,
+ system->atoms[i].x[0], system->atoms[i].x[1], system->atoms[i].x[2] );
+ }
fflush( out_control->trj );
- return 1;
+ return SUCCESS;
}