Skip to content
Snippets Groups Projects

Compare revisions

Changes are shown as if the source revision was being merged into the target revision. Learn more about comparing revisions.

Source

Select target project
No results found

Target

Select target project
  • SParTA/PuReMD
  • kaymakme/PuReMD
  • vermaasj/PuReMD
3 results
Show changes
Commits on Source (292)
Showing
with 4769 additions and 3452 deletions
......@@ -3,6 +3,7 @@
*.prs
*.pot
*.trj
*.res*
# TeX
*.aux
......
ACLOCAL_AMFLAGS = -I m4
AM_CFLAGS = -O2 -funroll-loops -Wall
SUBDIRS =
DIST_SUBDIRS =
EXTRA_DIST = data environ tools
......@@ -70,16 +72,16 @@ sPuReMD_lib_libspuremd_la_SOURCES = sPuReMD/src/allocate.c \
sPuReMD/src/valence_angles.h \
sPuReMD/src/vector.c \
sPuReMD/src/vector.h
sPuReMD_lib_libspuremd_la_CFLAGS = @DFT_CFLAGS@ @OMP_CFLAGS@
sPuReMD_lib_libspuremd_la_CPPFLAGS = -I sPuReMD/src
sPuReMD_lib_libspuremd_la_CFLAGS = @DFT_CFLAGS@ @CFLAGS_EXTRA@ @OMP_CFLAGS@
sPuReMD_lib_libspuremd_la_CPPFLAGS = -I sPuReMD/src @CPPFLAGS_EXTRA@
sPuReMD_lib_libspuremd_la_LDFLAGS = -version-info 1:0:0
sPuReMD_lib_libspuremd_la_LIBADD = @OMP_LIBS@ @L_LIBS@
sPuReMD_lib_libspuremd_la_LIBADD = @OMP_LIBS@ @Z_LIBS@ @L_LIBS@ @T_LIBS@ @G_LIBS@
bin_PROGRAMS = sPuReMD/bin/spuremd
sPuReMD_bin_spuremd_SOURCES = sPuReMD/src/driver.c
sPuReMD_bin_spuremd_CFLAGS = @DFT_CFLAGS@ @OMP_CFLAGS@
sPuReMD_bin_spuremd_CPPFLAGS = -I sPuReMD/src
sPuReMD_bin_spuremd_LDADD = sPuReMD/lib/libspuremd.la @OMP_LIBS@ @L_LIBS@
sPuReMD_bin_spuremd_CFLAGS = @DFT_CFLAGS@ @CFLAGS_EXTRA@ @OMP_CFLAGS@
sPuReMD_bin_spuremd_CPPFLAGS = -I sPuReMD/src @CPPFLAGS_EXTRA@
sPuReMD_bin_spuremd_LDADD = sPuReMD/lib/libspuremd.la @OMP_LIBS@ @Z_LIBS@ @L_LIBS@ @T_LIBS@ @G_LIBS@
check_PROGRAMS =
TESTS =
......@@ -180,8 +182,8 @@ PuReMD_bin_puremd_SOURCES = PuReMD/src/allocate.c \
PuReMD/src/vector.c \
PuReMD/src/vector.h
PuReMD_bin_puremd_CFLAGS = @DFT_CFLAGS@ @M_CFLAGS@
PuReMD_bin_puremd_CPPFLAGS = -I PuReMD/src
PuReMD_bin_puremd_CFLAGS = @DFT_CFLAGS@ @CFLAGS_EXTRA@ @M_CFLAGS@
PuReMD_bin_puremd_CPPFLAGS = -I PuReMD/src @CPPFLAGS_EXTRA@
PuReMD_bin_puremd_LDADD = @M_LIBS@ @L_LIBS@
endif
......@@ -190,15 +192,8 @@ if BUILD_GPU
if USE_CUDA
SUFFIXES = .cu
include cuda.am
endif
if USE_CUDA
# flags for CUDA compilation via NVCC (see cuda.am)
# Note: cc 13 for Tesla, cc 20 for Fermi, cc 30 for Kepler K10, cc 35 for Kepler K20/K40, etc.
#NVCCFLAGS += -gencode arch=compute_35,code=sm_35
NVCCFLAGS += $(NVCC_OPT_CODE)
NVCCFLAGS += --compiler-options "$(DEFS) $(NVCC_OPT_CODE_DEFS) @DFT_CFLAGS@"
#NVCCFLAGS += --ptxas-options -v
AM_NVCCFLAGS = --compiler-options="$(DEFS) $(AM_CFLAGS) @CFLAGS_EXTRA@ $(CFLAGS)" @NFLAGS@ @NFLAGS_EXTRA@
endif
bin_PROGRAMS = PuReMD-GPU/bin/spuremd
......@@ -315,8 +310,8 @@ PuReMD_GPU_bin_spuremd_SOURCES += PuReMD-GPU/src/cuda_allocate.cu \
nodist_EXTRA_PuReMD_GPU_bin_spuremd_SOURCES = PuReMD-GPU/src/dummy.c
endif
PuReMD_GPU_bin_spuremd_CFLAGS = @DFT_CFLAGS@
PuReMD_GPU_bin_spuremd_CPPFLAGS = -I PuReMD-GPU/src
PuReMD_GPU_bin_spuremd_CFLAGS = @DFT_CFLAGS@ @CFLAGS_EXTRA@
PuReMD_GPU_bin_spuremd_CPPFLAGS = -I PuReMD-GPU/src @CPPFLAGS_EXTRA@
if USE_CUDA
PuReMD_GPU_bin_spuremd_CFLAGS += @CUDA_CFLAGS@
......@@ -331,17 +326,18 @@ EXTRA_DIST += PG-PuReMD/src/cub
if USE_CUDA
SUFFIXES = .cu
include cuda.am
AM_NVCCFLAGS = --compiler-options="$(DEFS) $(AM_CFLAGS) -std=c++14 @CFLAGS_EXTRA@ $(CFLAGS) @OMP_CFLAGS@ @M_CFLAGS@" -std=c++14 @NFLAGS@ @NFLAGS_EXTRA@
endif
if USE_CUDA
# flags for CUDA compilation via NVCC (see cuda.am)
# Note: cc 13 for Tesla, cc 20 for Fermi, cc 30 for Kepler K10, cc 35 for Kepler K20/K40, etc.
#NVCCFLAGS += -gencode arch=compute_35,code=sm_35
NVCCFLAGS += $(NVCC_OPT_CODE)
NVCCFLAGS += --compiler-options "$(DEFS) $(NVCC_OPT_CODE_DEFS) @DFT_CFLAGS@ @M_CFLAGS@"
#NVCCFLAGS += --ptxas-options -v
if USE_HIP
SUFFIXES = .cu
include hip.am
AM_HIPCCFLAGS = $(DEFS) $(AM_CFLAGS) -std=c++14 @CFLAGS_EXTRA@ $(CFLAGS) @OMP_CFLAGS@ @M_CFLAGS@
endif
bin_PROGRAMS = PG-PuReMD/bin/pg-puremd
PG_PuReMD_bin_pg_puremd_SOURCES = PG-PuReMD/src/allocate.c \
PG-PuReMD/src/allocate.h \
......@@ -392,6 +388,7 @@ PG_PuReMD_bin_pg_puremd_SOURCES = PG-PuReMD/src/allocate.c \
PG-PuReMD/src/nonbonded.c \
PG-PuReMD/src/nonbonded.h \
PG-PuReMD/src/puremd.c \
PG-PuReMD/src/puremd.h \
PG-PuReMD/src/random.c \
PG-PuReMD/src/random.h \
PG-PuReMD/src/reax_types.h \
......@@ -413,6 +410,8 @@ PG_PuReMD_bin_pg_puremd_SOURCES = PG-PuReMD/src/allocate.c \
if USE_CUDA
PG_PuReMD_bin_pg_puremd_SOURCES += PG-PuReMD/src/cuda/cuda_allocate.cu \
PG-PuReMD/src/cuda/cuda_allocate.h \
PG-PuReMD/src/cuda/cuda_basic_comm.cu \
PG-PuReMD/src/cuda/cuda_basic_comm.h \
PG-PuReMD/src/cuda/cuda_bond_orders.cu \
PG-PuReMD/src/cuda/cuda_bond_orders.h \
PG-PuReMD/src/cuda/cuda_bonds.cu \
......@@ -423,6 +422,8 @@ PG_PuReMD_bin_pg_puremd_SOURCES += PG-PuReMD/src/cuda/cuda_allocate.cu \
PG-PuReMD/src/cuda/cuda_charges.h \
PG-PuReMD/src/cuda/cuda_copy.cu \
PG-PuReMD/src/cuda/cuda_copy.h \
PG-PuReMD/src/cuda/cuda_dense_lin_alg.cu \
PG-PuReMD/src/cuda/cuda_dense_lin_alg.h \
PG-PuReMD/src/cuda/cuda_environment.cu \
PG-PuReMD/src/cuda/cuda_environment.h \
PG-PuReMD/src/cuda/cuda_forces.cu \
......@@ -434,8 +435,6 @@ PG_PuReMD_bin_pg_puremd_SOURCES += PG-PuReMD/src/cuda/cuda_allocate.cu \
PG-PuReMD/src/cuda/cuda_init_md.h \
PG-PuReMD/src/cuda/cuda_integrate.cu \
PG-PuReMD/src/cuda/cuda_integrate.h \
PG-PuReMD/src/cuda/cuda_lin_alg.cu \
PG-PuReMD/src/cuda/cuda_lin_alg.h \
PG-PuReMD/src/cuda/cuda_list.cu \
PG-PuReMD/src/cuda/cuda_list.h \
PG-PuReMD/src/cuda/cuda_lookup.cu \
......@@ -454,7 +453,8 @@ PG_PuReMD_bin_pg_puremd_SOURCES += PG-PuReMD/src/cuda/cuda_allocate.cu \
PG-PuReMD/src/cuda/cuda_reduction.h \
PG-PuReMD/src/cuda/cuda_reset_tools.cu \
PG-PuReMD/src/cuda/cuda_reset_tools.h \
PG-PuReMD/src/cuda/cuda_shuffle.h \
PG-PuReMD/src/cuda/cuda_spar_lin_alg.cu \
PG-PuReMD/src/cuda/cuda_spar_lin_alg.h \
PG-PuReMD/src/cuda/cuda_system_props.cu \
PG-PuReMD/src/cuda/cuda_system_props.h \
PG-PuReMD/src/cuda/cuda_torsion_angles.cu \
......@@ -464,21 +464,81 @@ PG_PuReMD_bin_pg_puremd_SOURCES += PG-PuReMD/src/cuda/cuda_allocate.cu \
PG-PuReMD/src/cuda/cuda_valence_angles.cu \
PG-PuReMD/src/cuda/cuda_valence_angles.h \
PG-PuReMD/src/cuda/cuda_vector.h
if CUDA_DEBUG
PG_PuReMD_bin_pg_puremd_SOURCES += PG-PuReMD/src/cuda/cuda_validation.cu \
PG-PuReMD/src/cuda/cuda_validation.h
# dummy source to trigger C linking
nodist_EXTRA_PG_PuReMD_bin_pg_puremd_SOURCES = src/dummy.c
endif
if USE_HIP
PG_PuReMD_bin_pg_puremd_SOURCES += PG-PuReMD/src/hip/hip_allocate.cu \
PG-PuReMD/src/hip/hip_allocate.h \
PG-PuReMD/src/hip/hip_basic_comm.cu \
PG-PuReMD/src/hip/hip_basic_comm.h \
PG-PuReMD/src/hip/hip_bond_orders.cu \
PG-PuReMD/src/hip/hip_bond_orders.h \
PG-PuReMD/src/hip/hip_bonds.cu \
PG-PuReMD/src/hip/hip_bonds.h \
PG-PuReMD/src/hip/hip_box.cu \
PG-PuReMD/src/hip/hip_box.h \
PG-PuReMD/src/hip/hip_charges.cu \
PG-PuReMD/src/hip/hip_charges.h \
PG-PuReMD/src/hip/hip_copy.cu \
PG-PuReMD/src/hip/hip_copy.h \
PG-PuReMD/src/hip/hip_dense_lin_alg.cu \
PG-PuReMD/src/hip/hip_dense_lin_alg.h \
PG-PuReMD/src/hip/hip_environment.cu \
PG-PuReMD/src/hip/hip_environment.h \
PG-PuReMD/src/hip/hip_forces.cu \
PG-PuReMD/src/hip/hip_forces.h \
PG-PuReMD/src/hip/hip_helpers.h \
PG-PuReMD/src/hip/hip_hydrogen_bonds.cu \
PG-PuReMD/src/hip/hip_hydrogen_bonds.h \
PG-PuReMD/src/hip/hip_init_md.cu \
PG-PuReMD/src/hip/hip_init_md.h \
PG-PuReMD/src/hip/hip_integrate.cu \
PG-PuReMD/src/hip/hip_integrate.h \
PG-PuReMD/src/hip/hip_list.cu \
PG-PuReMD/src/hip/hip_list.h \
PG-PuReMD/src/hip/hip_lookup.cu \
PG-PuReMD/src/hip/hip_lookup.h \
PG-PuReMD/src/hip/hip_multi_body.cu \
PG-PuReMD/src/hip/hip_multi_body.h \
PG-PuReMD/src/hip/hip_neighbors.cu \
PG-PuReMD/src/hip/hip_neighbors.h \
PG-PuReMD/src/hip/hip_nonbonded.cu \
PG-PuReMD/src/hip/hip_nonbonded.h \
PG-PuReMD/src/hip/hip_post_evolve.cu \
PG-PuReMD/src/hip/hip_post_evolve.h \
PG-PuReMD/src/hip/hip_random.cu \
PG-PuReMD/src/hip/hip_random.h \
PG-PuReMD/src/hip/hip_reduction.cu \
PG-PuReMD/src/hip/hip_reduction.h \
PG-PuReMD/src/hip/hip_reset_tools.cu \
PG-PuReMD/src/hip/hip_reset_tools.h \
PG-PuReMD/src/hip/hip_spar_lin_alg.cu \
PG-PuReMD/src/hip/hip_spar_lin_alg.h \
PG-PuReMD/src/hip/hip_system_props.cu \
PG-PuReMD/src/hip/hip_system_props.h \
PG-PuReMD/src/hip/hip_torsion_angles.cu \
PG-PuReMD/src/hip/hip_torsion_angles.h \
PG-PuReMD/src/hip/hip_utils.cu \
PG-PuReMD/src/hip/hip_utils.h \
PG-PuReMD/src/hip/hip_valence_angles.cu \
PG-PuReMD/src/hip/hip_valence_angles.h \
PG-PuReMD/src/hip/hip_vector.h
# dummy source to trigger C linking
nodist_EXTRA_PG_PuReMD_bin_pg_puremd_SOURCES = src/dummy.c
endif
PG_PuReMD_bin_pg_puremd_CFLAGS = @DFT_CFLAGS@ @M_CFLAGS@
PG_PuReMD_bin_pg_puremd_CPPFLAGS = -I PG-PuReMD/src
PG_PuReMD_bin_pg_puremd_LDADD = @M_LIBS@ @L_LIBS@ -lstdc++
PG_PuReMD_bin_pg_puremd_CFLAGS = -std=c11 @CFLAGS_EXTRA@ @OMP_CFLAGS@ @M_CFLAGS@
PG_PuReMD_bin_pg_puremd_CPPFLAGS = -I PG-PuReMD/src @CPPFLAGS_EXTRA@
PG_PuReMD_bin_pg_puremd_LDADD = @OMP_LIBS@ @M_LIBS@ @L_LIBS@ -lstdc++
if USE_CUDA
PG_PuReMD_bin_pg_puremd_CFLAGS += @CUDA_CFLAGS@
PG_PuReMD_bin_pg_puremd_LDADD += @CUDA_LIBS@
endif
if USE_HIP
PG_PuReMD_bin_pg_puremd_CFLAGS += @HIP_CFLAGS@
PG_PuReMD_bin_pg_puremd_LDADD += @HIP_LIBS@
endif
endif
if BUILD_DOC
......
This diff is collapsed.
......@@ -33,18 +33,21 @@ void Init_Matrix_Row_Indices( sparse_matrix * const, int * const );
void PreAllocate_Space( reax_system * const, control_params * const, storage * const );
void ReAllocate_System( reax_system * const, int, int );
void Reallocate_System_Part1( reax_system * const, int );
void Allocate_Workspace( reax_system * const, control_params * const, storage * const,
int, int );
void Reallocate_System_Part2( reax_system * const, int, int );
void Allocate_Workspace_Part1( reax_system * const,
control_params * const, storage * const, int );
void Allocate_Workspace_Part2( reax_system * const,
control_params * const, storage * const, int );
void Allocate_Grid( reax_system * const, MPI_Comm );
void Deallocate_Grid( grid * const );
void Allocate_MPI_Buffers( mpi_datatypes * const, int, neighbor_proc * const );
void Allocate_Matrix( sparse_matrix * const, int, int, int );
void Allocate_Matrix( sparse_matrix * const, int, int, int, int );
void Deallocate_Matrix( sparse_matrix * const );
......@@ -54,7 +57,11 @@ int Allocate_Bond_List( int, int * const, reax_list * const );
void Deallocate_MPI_Buffers( mpi_datatypes * const );
void ReAllocate( reax_system * const, control_params * const,
void Reallocate_Part1( reax_system * const, control_params * const,
simulation_data * const, storage * const,
reax_list** const, mpi_datatypes * const );
void Reallocate_Part2( reax_system * const, control_params * const,
simulation_data * const, storage * const,
reax_list** const, mpi_datatypes * const );
......
......@@ -19,7 +19,10 @@
<http://www.gnu.org/licenses/>.
----------------------------------------------------------------------*/
#include "reax_types.h"
#if (defined(HAVE_CONFIG_H) && !defined(__CONFIG_H_))
#define __CONFIG_H_
#include "../../common/include/config.h"
#endif
#include "analyze.h"
......@@ -27,20 +30,25 @@
#include "list.h"
#include "vector.h"
#define MAX_FRAGMENT_TYPES 100
#define MAX_FRAGMENT_TYPES (100)
void Print_Molecule( reax_system *system, molecule *m, int mode, char *s )
void Print_Molecule( reax_system *system, molecule *m, int mode,
char * const s, size_t size )
{
int j, atom;
s[0] = 0;
if ( mode == 1 )
{
/* print molecule summary */
for ( j = 0; j < MAX_ATOM_TYPES; ++j )
{
if ( m->mtypes[j] )
sprintf( s, "%s%s%d", s, system->reax_param.sbp[j].name, m->mtypes[j] );
{
snprintf( s, size - 1, "%s%d",
system->reax_param.sbp[j].name, m->mtypes[j] );
s[size - 1] = '\0';
}
}
}
else if ( mode == 2 )
{
......@@ -48,10 +56,11 @@ void Print_Molecule( reax_system *system, molecule *m, int mode, char *s )
for ( j = 0; j < m->atom_count; ++j )
{
atom = m->atom_list[j];
sprintf( s, "%s%s(%d)",
s,
system->reax_param.sbp[system->my_atoms[atom].type].name,
atom );
snprintf( s, size - 1, "%s(%d)",
system->reax_param.sbp[system->my_atoms[atom].type].name,
atom );
s[size - 1] = '\0';
}
}
}
......@@ -91,15 +100,14 @@ void Analyze_Fragments( reax_system *system, control_params *control,
simulation_data *data, storage *workspace,
reax_list **lists, FILE *fout, int ignore )
{
int atom, i, flag;
int *mark = workspace->mark;
int num_fragments, num_fragment_types;
char fragment[MAX_ATOM_TYPES];
char fragments[MAX_FRAGMENT_TYPES][MAX_ATOM_TYPES];
int fragment_count[MAX_FRAGMENT_TYPES];
int atom, i, flag, *mark;
int num_fragments, num_fragment_types, fragment_count[MAX_FRAGMENT_TYPES];
char fragment[MAX_ATOM_TYPES], fragments[MAX_FRAGMENT_TYPES][MAX_ATOM_TYPES];
molecule m;
reax_list *new_bonds = lists[BONDS];
//list *old_bonds = lists[OLD_BONDS];
reax_list *new_bonds;
mark = workspace->mark;
new_bonds = lists[BONDS];
/* fragment analysis */
fprintf( fout, "step%d fragments\n", data->step );
......@@ -108,47 +116,51 @@ void Analyze_Fragments( reax_system *system, control_params *control,
memset( mark, 0, system->N * sizeof(int) );
for ( atom = 0; atom < system->N; ++atom )
{
if ( !mark[atom] )
{
/* discover a new fragment */
memset( m.mtypes, 0, MAX_ATOM_TYPES * sizeof(int) );
Visit_Bonds( atom, mark, m.mtypes, system, control, new_bonds, ignore );
++num_fragments;
Print_Molecule( system, &m, 1, fragment );
Print_Molecule( system, &m, 1, fragment, sizeof(fragment) );
/* check if a similar fragment already exists */
flag = 0;
for ( i = 0; i < num_fragment_types; ++i )
if ( !strcmp( fragments[i], fragment ) )
{
if ( strncmp( fragments[i], fragment, MAX_ATOM_TYPES ) == 0 )
{
fragment_count[i]++;
flag = 1;
break;
}
}
if ( flag == 0 )
{
/* it is a new one, add to the fragments list */
strcpy( fragments[num_fragment_types], fragment );
strncpy( fragments[num_fragment_types], fragment, sizeof(fragments[num_fragment_types]) - 1 );
fragments[num_fragment_types][sizeof(fragments[num_fragment_types]) - 1] = '\0';
fragment_count[num_fragment_types] = 1;
++num_fragment_types;
}
}
}
/* output the results of fragment analysis */
for ( i = 0; i < num_fragment_types; ++i )
{
/* strlen safe here as fragments[i] is NULL-terminated above */
if ( strlen(fragments[i]) )
{
fprintf( fout, "%d of %s\n", fragment_count[i], fragments[i] );
}
}
fprintf( fout, "\n" );
#if defined(DEBUG)
fflush( fout );
/* compare new bonds to old bonds */
//for( atom = 0; atom < system->N; ++atom ) {
// fprintf( fout, "atom: %d\n", atom ); fflush( fout );
// Compare_Bonding( atom, system, control, workspace,
// old_bonds, new_bonds, fout );
//}
//Copy_Bond_List( system, control, lists );
#endif
}
......@@ -175,8 +187,4 @@ void Analysis( reax_system *system, control_params *control,
Analyze_Fragments( system, control, data, workspace, lists,
out_control->ign, 1 );
}
#if defined(DEBUG)
fprintf( stderr, "analysis... done\n" );
#endif
}
This diff is collapsed.
......@@ -38,18 +38,14 @@ enum pointer_type
extern "C" {
#endif
void Dist( const reax_system * const, mpi_datatypes * const,
void*, int, MPI_Datatype );
void Dist( reax_system const * const, mpi_datatypes * const,
void const * const, int, MPI_Datatype );
void Coll( const reax_system * const, mpi_datatypes * const,
void*, int, MPI_Datatype );
real Parallel_Norm( const real * const, const int, MPI_Comm );
void Coll( reax_system const * const, mpi_datatypes * const,
void * const , int, MPI_Datatype );
real Parallel_Dot( const real * const, const real * const, const int, MPI_Comm );
real Parallel_Vector_Acc( const real * const, const int, MPI_Comm );
#if defined(TEST_FORCES)
void Coll_ids_at_Master( reax_system*, storage*, mpi_datatypes* );
......
This diff is collapsed.
......@@ -48,14 +48,16 @@ typedef struct
extern "C" {
#endif
void Add_dBond_to_Forces( int, int, storage * const, reax_list ** const );
void Add_dBond_to_Forces( int, int, reax_system*, simulation_data*,
storage*, reax_list** );
void Add_dBond_to_Forces_NPT( int, int, simulation_data * const, storage * const,
reax_list ** const );
void Add_dBond_to_Forces_NPT( int, int, reax_system*, simulation_data*,
storage*, reax_list** );
int BOp( storage * const, reax_list * const, real, int, int, const far_neighbor_data * const,
const single_body_parameters * const, const single_body_parameters * const,
const two_body_parameters * const );
int BOp( storage * const, reax_list * const, real, int, int, int,
ivec * const, real, rvec * const, int,
single_body_parameters const * const, single_body_parameters const * const,
two_body_parameters const * const );
void BO( reax_system * const, control_params * const, simulation_data * const,
storage * const, reax_list ** const, output_controls * const );
......
......@@ -19,16 +19,21 @@
<http://www.gnu.org/licenses/>.
----------------------------------------------------------------------*/
#include "reax_types.h"
#if (defined(HAVE_CONFIG_H) && !defined(__CONFIG_H_))
#define __CONFIG_H_
#include "../../common/include/config.h"
#endif
#if defined(PURE_REAX)
#include "bonds.h"
#include "bond_orders.h"
#include "list.h"
#include "tool_box.h"
#include "vector.h"
#elif defined(LAMMPS_REAX)
#include "reax_bonds.h"
#include "reax_bond_orders.h"
#include "reax_list.h"
#include "reax_tool_box.h"
......@@ -46,22 +51,24 @@ void Bonds( reax_system * const system, control_params * const control,
int start_i, end_i;
int type_i, type_j;
real ebond, pow_BOs_be2, exp_be12, CEbo;
const real gp3 = system->reax_param.gp.l[3];
const real gp4 = system->reax_param.gp.l[4];
const real gp7 = system->reax_param.gp.l[7];
const real gp10 = system->reax_param.gp.l[10];
const real gp37 = (int) system->reax_param.gp.l[37];
real gp3, gp4, gp7, gp10;
real exphu, exphua1, exphub1, exphuov, hulpov, estriph;
real decobdbo, decobdboua, decobdboub;
single_body_parameters *sbp_i, *sbp_j;
two_body_parameters *twbp;
bond_order_data *bo_ij;
reax_list * const bond_list = lists[BONDS];
reax_list *bond_list;
bond_list = lists[BONDS];
gp3 = system->reax_param.gp.l[3];
gp4 = system->reax_param.gp.l[4];
gp7 = system->reax_param.gp.l[7];
gp10 = system->reax_param.gp.l[10];
for ( i = 0; i < system->n; ++i )
{
start_i = Start_Index( i, bond_list );
end_i = End_Index( i, bond_list );
start_i = Start_Index(i, bond_list);
end_i = End_Index(i, bond_list);
for ( pj = start_i; pj < end_i; ++pj )
{
......@@ -77,74 +84,78 @@ void Bonds( reax_system * const system, control_params * const control,
index_tbp(type_i, type_j, system->reax_param.num_atom_types) ];
bo_ij = &bond_list->bond_list[pj].bo_data;
/* calculate the constants */
pow_BOs_be2 = POW( bo_ij->BO_s, twbp->p_be2 );
exp_be12 = EXP( twbp->p_be1 * ( 1.0 - pow_BOs_be2 ) );
CEbo = -twbp->De_s * exp_be12 *
( 1.0 - twbp->p_be1 * twbp->p_be2 * pow_BOs_be2 );
CEbo = -twbp->De_s * exp_be12
* (1.0 - twbp->p_be1 * twbp->p_be2 * pow_BOs_be2);
/* calculate the Bond Energy */
/* calculate bond energy */
ebond = -twbp->De_s * bo_ij->BO_s * exp_be12
- twbp->De_p * bo_ij->BO_pi
- twbp->De_pp * bo_ij->BO_pi2;
data->my_en.e_bond += ebond;
data->my_en->e_bond += ebond;
/* calculate derivatives of Bond Orders */
/* calculate derivatives of bond orders */
bo_ij->Cdbo += CEbo;
bo_ij->Cdbopi -= (CEbo + twbp->De_p);
bo_ij->Cdbopi2 -= (CEbo + twbp->De_pp);
bo_ij->Cdbopi -= CEbo + twbp->De_p;
bo_ij->Cdbopi2 -= CEbo + twbp->De_pp;
#ifdef TEST_ENERGY
#if defined(TEST_ENERGY)
//fprintf( out_control->ebond, "%6d%6d%24.15e%24.15e%24.15e\n",
fprintf( out_control->ebond, "%6d%6d%12.4f%12.4f%12.4f\n",
system->my_atoms[i].orig_id,
system->my_atoms[j].orig_id,
bo_ij->BO, ebond, data->my_en.e_bond );
bo_ij->BO, ebond, data->my_en->e_bond );
#endif
#ifdef TEST_FORCES
#if defined(TEST_FORCES)
Add_dBO( system, lists, i, pj, CEbo, workspace->f_be );
Add_dBOpinpi2( system, lists, i, pj,
-(CEbo + twbp->De_p), -(CEbo + twbp->De_pp),
workspace->f_be, workspace->f_be );
-(CEbo + twbp->De_p),
-(CEbo + twbp->De_pp),
workspace->f_be, workspace->f_be );
#endif
/* Stabilisation terminal triple bond */
/* Stabilisation terminal triple bond in C-O */
if ( bo_ij->BO >= 1.00 )
{
/* changed to avoid equality checks of doubles, but really need better approach... */
if ( gp37 == 2
|| (sbp_i->mass - 12.0000 < 1e-6 && sbp_j->mass - 15.9990 < 1e-6)
|| (sbp_j->mass - 12.0000 < 1e-6 && sbp_i->mass - 15.9990 < 1e-6) )
if ( (strncmp( sbp_i->name, "C", sizeof(sbp_i->name) ) == 0
&& strncmp( sbp_j->name, "O", sizeof(sbp_j->name) ) == 0)
|| (strncmp( sbp_i->name, "O", sizeof(sbp_i->name) ) == 0
&& strncmp( sbp_j->name, "C", sizeof(sbp_j->name) ) == 0) )
{
exphu = EXP( -gp7 * SQR(bo_ij->BO - 2.50) );
exphua1 = EXP( -gp3 * (workspace->total_bond_order[i] - bo_ij->BO) );
exphub1 = EXP( -gp3 * (workspace->total_bond_order[j] - bo_ij->BO) );
exphuov = EXP( gp4 * (workspace->Delta[i] + workspace->Delta[j]) );
//ba = SQR( bo_ij->BO - 2.5 );
exphu = EXP( -gp7 * SQR(bo_ij->BO - 2.5) );
//oboa = abo(j1) - boa;
//obob = abo(j2) - boa;
exphua1 = EXP(-gp3 * (workspace->total_bond_order[i] - bo_ij->BO));
exphub1 = EXP(-gp3 * (workspace->total_bond_order[j] - bo_ij->BO));
//ovoab = abo(j1) - aval(it1) + abo(j2) - aval(it2);
exphuov = EXP(gp4 * (workspace->Delta[i] + workspace->Delta[j]));
hulpov = 1.0 / (1.0 + 25.0 * exphuov);
estriph = gp10 * exphu * hulpov * (exphua1 + exphub1);
data->my_en.e_bond += estriph;
data->my_en->e_bond += estriph;
decobdbo = gp10 * exphu * hulpov * (exphua1 + exphub1)
* ( gp3 - 2.0 * gp7 * (bo_ij->BO - 2.50) );
decobdboua = -gp10 * exphu * hulpov * (gp3 * exphua1
+ 25.0 * gp4 * exphuov * hulpov * (exphua1 + exphub1));
decobdboub = -gp10 * exphu * hulpov * (gp3 * exphub1
+ 25.0 * gp4 * exphuov * hulpov * (exphua1 + exphub1));
* ( gp3 - 2.0 * gp7 * (bo_ij->BO - 2.5) );
decobdboua = -gp10 * exphu * hulpov
* (gp3 * exphua1 + 25.0 * gp4 * exphuov * hulpov * (exphua1 + exphub1));
decobdboub = -gp10 * exphu * hulpov
* (gp3 * exphub1 + 25.0 * gp4 * exphuov * hulpov * (exphua1 + exphub1));
bo_ij->Cdbo += decobdbo;
workspace->CdDelta[i] += decobdboua;
workspace->CdDelta[j] += decobdboub;
#ifdef TEST_ENERGY
#if defined(TEST_ENERGY)
//fprintf( out_control->ebond,
// "%6d%6d%24.15e%24.15e%24.15e%24.15e\n",
// system->my_atoms[i].orig_id, system->my_atoms[j].orig_id,
// estriph, decobdbo, decobdboua, decobdboub );
#endif
#ifdef TEST_FORCES
#if defined(TEST_FORCES)
Add_dBO( system, lists, i, pj, decobdbo, workspace->f_be );
Add_dDelta( system, lists, i, decobdboua, workspace->f_be );
Add_dDelta( system, lists, j, decobdboub, workspace->f_be );
......
......@@ -19,14 +19,26 @@
<http://www.gnu.org/licenses/>.
----------------------------------------------------------------------*/
#include "reax_types.h"
#include "box.h"
#if (defined(HAVE_CONFIG_H) && !defined(__CONFIG_H_))
#define __CONFIG_H_
#include "../../common/include/config.h"
#endif
#include "comm_tools.h"
#include "io_tools.h"
#include "system_props.h"
#include "vector.h"
#if defined(PURE_REAX)
#include "box.h"
#include "comm_tools.h"
#include "io_tools.h"
#include "system_props.h"
#include "vector.h"
#elif defined(LAMMPS_REAX)
#include "reax_box.h"
#include "reax_comm_tools.h"
#include "reax_io_tools.h"
#include "reax_system_props.h"
#include "reax_vector.h"
#endif
void Make_Consistent( simulation_box * const box )
......@@ -144,10 +156,10 @@ void Setup_Big_Box( real a, real b, real c, real alpha, real beta, real gamma,
{
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) )
if ( !isfinite(a) || !isfinite(b) || !isfinite(c)
|| !isfinite(alpha) || !isfinite(beta) || !isfinite(gamma) )
{
fprintf( stderr, "[ERROR] Invalid simulation box boundaries for big box (NaN). Terminating...\n" );
fprintf( stderr, "[ERROR] Invalid simulation box boundaries for big box (NaN or infinite). Terminating...\n" );
exit( INVALID_INPUT );
}
......@@ -169,8 +181,8 @@ void Setup_Big_Box( real a, real b, real c, real alpha, real beta, real gamma,
box->box[2][1] = c * zi;
box->box[2][2] = c * SQRT(1.0 - SQR(c_beta) - SQR(zi));
#if defined(DEBUG)
fprintf( stderr, "box is %8.2f x %8.2f x %8.2f\n",
#if defined(DEBUG_FOCUS)
fprintf( stderr, "[INFO] box = (%8.4f, %8.4f, %8.4f)\n",
box->box[0][0], box->box[1][1], box->box[2][2] );
#endif
......@@ -184,7 +196,7 @@ void Init_Box( rtensor box_tensor, simulation_box * const box )
rtensor_Copy( box->box, box_tensor );
Make_Consistent( box );
#if defined(DEBUG)
#if defined(DEBUG_FOCUS)
fprintf( stderr, "box is %8.2f x %8.2f x %8.2f\n",
box->box[0][0], box->box[1][1], box->box[2][2] );
#endif
......@@ -232,14 +244,14 @@ void Setup_My_Ext_Box( reax_system * const system, control_params * const contro
for ( d = 0; d < 3; ++d )
{
/* estimate the number of native cells */
native_gcells[d] = (int)(my_box->box_norms[d] / (control->vlist_cut / 2));
native_gcells[d] = (int) (my_box->box_norms[d] / (control->vlist_cut / 2.0));
if ( native_gcells[d] == 0 )
{
native_gcells[d] = 1;
}
gcell_len[d] = my_box->box_norms[d] / native_gcells[d];
ghost_gcells[d] = (int) CEIL(bc->ghost_cutoff / gcell_len[d]);
ghost_gcells[d] = (int) CEIL( bc->ghost_cutoff / gcell_len[d] );
/* extend my box with the ghost regions */
my_ext_box->min[d] = my_box->min[d] - ghost_gcells[d] * gcell_len[d];
......@@ -260,7 +272,7 @@ void Setup_Boundary_Cutoffs( reax_system * const system, control_params * const
bc->ghost_nonb = control->nonb_cut;
bc->ghost_hbond = control->hbond_cut;
bc->ghost_bond = 2 * control->bond_cut;
bc->ghost_bond = 2.0 * control->bond_cut;
bc->ghost_cutoff = MAX( control->vlist_cut, bc->ghost_bond );
#if defined(DEBUG_FOCUS)
......@@ -275,17 +287,21 @@ void Setup_Boundary_Cutoffs( reax_system * const system, control_params * const
void Setup_Environment( reax_system * const system, control_params * const control,
mpi_datatypes * const mpi_data )
{
int ret;
ivec periodic = {1, 1, 1};
#if defined(DEBUG_FOCUS)
char temp[100] = "";
char temp[100];
#endif
/* initialize communicator - 3D mesh with wrap-arounds = 3D torus */
MPI_Cart_create( MPI_COMM_WORLD, 3, control->procs_by_dim, periodic, 1,
ret = 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_Cart_coords( mpi_data->comm_mesh3D, system->my_rank, 3,
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Comm_rank( mpi_data->comm_mesh3D, &system->my_rank );
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Cart_coords( mpi_data->comm_mesh3D, system->my_rank, 3,
system->my_coords );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Setup_Boundary_Cutoffs( system, control );
Setup_My_Box( system, control );
......@@ -296,13 +312,18 @@ void Setup_Environment( reax_system * const system, control_params * const contr
fprintf( stderr, "p%d coord: %d %d %d\n",
system->my_rank, system->my_coords[0],
system->my_coords[1], system->my_coords[2] );
sprintf( temp, "p%d big_box", system->my_rank );
sprintf( temp, sizeof(temp) - 1, "p%d big_box", system->my_rank );
temp[sizeof(temp) - 1] = '\0';
Print_Box( &system->big_box, temp, stderr );
sprintf( temp, "p%d my_box", system->my_rank );
sprintf( temp, sizeof(temp) - 1, "p%d my_box", system->my_rank );
temp[sizeof(temp) - 1] = '\0';
Print_Box( &system->my_box, temp, stderr );
sprintf( temp, "p%d ext_box", system->my_rank );
sprintf( temp, sizeof(temp) - 1, "p%d ext_box", system->my_rank );
temp[sizeof(temp) - 1] = '\0';
Print_Box( &system->my_ext_box, temp, stderr );
MPI_Barrier( MPI_COMM_WORLD );
fprintf( stderr, "p%d: parallel environment initialized\n",
system->my_rank );
......@@ -324,7 +345,7 @@ void Scale_Box( reax_system * const system, control_params * const control,
if ( control->ensemble == iNPT )
{
mu[0] = POW( 1.0 + (dt / control->Tau_P[0]) * (data->iso_bar.P - control->P[0]),
1. / 3 );
1.0 / 3.0 );
if ( mu[0] < MIN_dV )
{
......@@ -342,8 +363,8 @@ void Scale_Box( reax_system * const system, control_params * const control,
{
for ( d = 0; d < 3; ++d )
{
mu[d] = POW(1.0 + (dt / control->Tau_P[d]) * (data->tot_press[d] - control->P[d]),
1. / 3 );
mu[d] = POW( 1.0 + (dt / control->Tau_P[d]) * (data->tot_press[d] - control->P[d]),
1.0 / 3.0 );
if ( mu[d] < MIN_dV )
{
......@@ -357,45 +378,45 @@ void Scale_Box( reax_system * const system, control_params * const control,
}
/* temperature scaler */
lambda = 1.0 + (dt / control->Tau_T) * (control->T / data->therm.T - 1.0);
lambda = 1.0 + ((dt * 1.0e-12) / control->Tau_T)
* (control->T / data->therm.T - 1.0);
if ( lambda < MIN_dT )
{
lambda = MIN_dT;
}
else if (lambda > MAX_dT )
lambda = SQRT( lambda );
if ( lambda > MAX_dT )
{
lambda = MAX_dT;
}
lambda = SQRT( lambda );
/* Scale velocities and positions at t+dt */
for ( i = 0; i < system->n; ++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];
atom->x[2] = mu[2] * atom->x[2];
}
Compute_Kinetic_Energy( system, data, mpi_data->comm_mesh3D );
#if defined(DEBUG)
fprintf( stderr, "damping - " );
#endif
/* update kinetic energy and temperature based on new positions and velocities */
Compute_Kinetic_Energy( system, data, mpi_data->comm_mesh3D );
/* update box & grid */
system->big_box.box[0][0] *= mu[0];
system->big_box.box[1][1] *= mu[1];
system->big_box.box[2][2] *= mu[2];
Make_Consistent( &(system->big_box) );
Make_Consistent( &system->big_box );
Setup_My_Box( system, control );
Setup_My_Ext_Box( system, control );
Update_Comm( system );
#if defined(DEBUG)
fprintf( stderr, "box & grid updated - " );
#endif
}
......
This diff is collapsed.
......@@ -29,9 +29,12 @@
extern "C" {
#endif
void Compute_Charges( reax_system * const, control_params * const,
simulation_data * const, storage * const,
const output_controls * const, mpi_datatypes * const );
int is_refactoring_step( control_params const * const, simulation_data * const );
void Compute_Charges( reax_system const * const, control_params const * const,
simulation_data * const,
storage * const, const output_controls * const,
mpi_datatypes * const );
#ifdef __cplusplus
}
......
This diff is collapsed.
......@@ -29,19 +29,30 @@
extern "C" {
#endif
void Check_MPI_Error( int, const char * const );
#if defined(NEUTRAL_TERRITORY)
void Setup_NT_Comm( reax_system * const, control_params * const,
mpi_datatypes * const );
void Estimate_NT_Atoms( reax_system * const, mpi_datatypes * const );
#endif
void Check_MPI_Error( int, const char * const, int );
void Setup_Comm( reax_system * const, control_params * const, mpi_datatypes * const );
void Update_Comm( reax_system * const );
void Estimate_Boundary_Atoms( reax_system * const, int, int, int, mpi_out_data * const );
void Count_Boundary_Atoms( reax_system const * const,
int, int, int, mpi_out_data * const, int * const );
void Sort_Boundary_Atoms( reax_system * const, int, int,
int, mpi_out_data * const, mpi_datatypes * const );
void Unpack_Estimate_Message( reax_system * const, int, void * const, int,
neighbor_proc * const, int );
void Unpack_Exchange_Message( reax_system * const, int, void * const,
int, neighbor_proc * const, int );
int SendRecv( reax_system * const, mpi_datatypes * const, MPI_Datatype, int * const,
message_sorter, unpacker, int );
int SendRecv( reax_system * const, mpi_datatypes * const, MPI_Datatype,
message_counter, message_sorter, unpacker, int );
void Comm_Atoms( reax_system * const, control_params * const, simulation_data * const, storage * const,
mpi_datatypes * const, int );
......
This diff is collapsed.
Subproject commit d622848f9fb62f13e5e064e1deb43b6bcbb12bad
Subproject commit c3cceac115c072fb63df1836ff46d8c60d9eb304
This diff is collapsed.
......@@ -3,39 +3,52 @@
#include "../reax_types.h"
#ifdef __cplusplus
extern "C" {
#endif
void * sCudaHostAllocWrapper( size_t, const char * const, int );
void dev_alloc_system( reax_system * );
void * sCudaHostReallocWrapper( void *, size_t, size_t, const char * const, int );
void dev_alloc_grid( reax_system * );
void * sCudaHostCallocWrapper( size_t, size_t, const char * const, int );
void dev_alloc_simulation_data( simulation_data * );
void sCudaFreeHostWrapper( void *, const char * const, int );
void dev_alloc_workspace( reax_system *, control_params *, storage *, int, int );
#ifdef __cplusplus
}
#endif
void dev_alloc_matrix( sparse_matrix *, int, int );
void Cuda_Allocate_System( reax_system *, control_params * );
void dev_alloc_control( control_params * );
void Cuda_Allocate_Grid( reax_system *, control_params * );
void dev_dealloc_grid_cell_atoms( reax_system * );
void Cuda_Allocate_Simulation_Data( simulation_data *, cudaStream_t );
void dev_alloc_grid_cell_atoms( reax_system *, int );
void Cuda_Allocate_Control( control_params * );
void dev_realloc_system( reax_system *, int , int , char * );
void Cuda_Allocate_Workspace_Part1( reax_system *, control_params *, storage *, int );
void dev_dealloc_workspace( control_params *, storage * );
void Cuda_Allocate_Workspace_Part2( reax_system *, control_params *, storage *, int );
void dev_dealloc_matrix( sparse_matrix * );
void Cuda_Allocate_Matrix( sparse_matrix * const, int, int, int, int, cudaStream_t );
void Cuda_ReAllocate( reax_system*, control_params*, simulation_data*, storage*,
reax_list**, mpi_datatypes* );
void Cuda_Deallocate_Grid_Cell_Atoms( reax_system * );
void Cuda_Allocate_Grid_Cell_Atoms( reax_system *, int );
void Cuda_Deallocate_Workspace_Part1( control_params *, storage * );
void Cuda_Deallocate_Workspace_Part2( control_params *, storage * );
void Cuda_Deallocate_Matrix( sparse_matrix * );
void Cuda_Reallocate_Part1( reax_system *, control_params *, simulation_data *, storage *,
reax_list **, mpi_datatypes * );
void Cuda_Reallocate_Part2( reax_system *, control_params *, simulation_data *, storage *,
reax_list **, mpi_datatypes * );
#ifdef __cplusplus
}
#endif
#endif
This diff is collapsed.