jianmu
/
lammps
forked from lijiext/lammps
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@14911 f3b2605a-c512-4ea7-a41b-209d697bcdaa

This commit is contained in:
sjplimp 2016-04-30 18:07:04 +00:00
parent 00f38fdaf0
commit 997099253b
9 changed files with 113 additions and 112 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

1
src/BODY/fix_nh_body.cpp
View File

@ -104,7 +104,6 @@ void FixNHBody::nve_x()
AtomVecBody::Bonus *bonus = avec->bonus;
int *body = atom->body;
double **angmom = atom->angmom;
double *rmass = atom->rmass;
int *mask = atom->mask;
int nlocal = atom->nlocal;
if (igroup == atom->firstgroup) nlocal = atom->nfirst;

9
src/MANYBODY/pair_airebo.cpp
View File

@ -3184,12 +3184,11 @@ double PairAIREBO::piRCSpline(double Nij, double Nji, double Nijconj,
if (done==0) {
for (i=0; i<piCCdom[0][1]; i++)
if (Nij>=(double) i && Nij<=(double) i+1 || Nij==(double) i) x=i;
if (Nij>=(double) i && Nij<=(double) i+1) x=i;
for (i=0; i<piCCdom[1][1]; i++)
if (Nji>=(double) i && Nji<=(double) i+1 || Nji==(double) i) y=i;
if (Nji>=(double) i && Nji<=(double) i+1) y=i;
for (i=0; i<piCCdom[2][1]; i++)
if (Nijconj>=(double) i && Nijconj<=(double) i+1 ||
Nijconj==(double) i) z=i;
if (Nijconj>=(double) i && Nijconj<=(double) i+1) z=i;
for (i=0; i<64; i++) coeffs[i]=piCC[x][y][z][i];
piRC=Sptricubic(Nij,Nji,Nijconj,coeffs,dN3);
@ -3199,7 +3198,7 @@ double PairAIREBO::piRCSpline(double Nij, double Nji, double Nijconj,
// CH interaction
if (typei==0 && typej==1 || typei==1 && typej==0) {
if ((typei==0 && typej==1) || (typei==1 && typej==0)) {
// if the inputs are out of bounds set them back to a point in bounds
if (Nij<piCHdom[0][0] || Nij>piCHdom[0][1] ||

2
src/RIGID/fix_rigid_nh.cpp
View File

@ -84,7 +84,7 @@ FixRigidNH::FixRigidNH(LAMMPS *lmp, int narg, char **arg) :
error->all(FLERR,
"Cannot use fix rigid npt/nph on a non-periodic dimension");
if (pcouple == XYZ && dimension == 3 &&
if (pcouple == XYZ && dimension == 3 &&
(p_start[0] != p_start[1] || p_start[0] != p_start[2] ||
p_stop[0] != p_stop[1] || p_stop[0] != p_stop[2] ||
p_period[0] != p_period[1] || p_period[0] != p_period[2]))

190
src/USER-MISC/fix_ttm_mod.cpp
View File

@ -613,110 +613,110 @@ void FixTTMMod::end_of_step()
double stability_criterion = 0.0;
for (int ixnode = 0; ixnode < nxnodes; ixnode++)
for (int iynode = 0; iynode < nynodes; iynode++)
for (int iznode = 0; iznode < nznodes; iznode++)
T_electron_first[ixnode][iynode][iznode] =
T_electron[ixnode][iynode][iznode];
for (int iynode = 0; iynode < nynodes; iynode++)
for (int iznode = 0; iznode < nznodes; iznode++)
T_electron_first[ixnode][iynode][iznode] =
T_electron[ixnode][iynode][iznode];
do {
for (int ixnode = 0; ixnode < nxnodes; ixnode++)
for (int ixnode = 0; ixnode < nxnodes; ixnode++)
for (int iynode = 0; iynode < nynodes; iynode++)
for (int iznode = 0; iznode < nznodes; iznode++)
T_electron[ixnode][iynode][iznode] =
T_electron_first[ixnode][iynode][iznode];
stability_criterion = 1.0 -
2.0*inner_dt/el_specific_heat *
(el_thermal_conductivity*(1.0/dx/dx + 1.0/dy/dy + 1.0/dz/dz));
if (stability_criterion < 0.0) {
inner_dt = 0.25*el_specific_heat /
(el_thermal_conductivity*(1.0/dx/dx + 1.0/dy/dy + 1.0/dz/dz));
}
num_inner_timesteps = static_cast<unsigned int>(update->dt/inner_dt) + 1;
inner_dt = update->dt/double(num_inner_timesteps);
if (num_inner_timesteps > 1000000)
error->warning(FLERR,"Too many inner timesteps in fix ttm",0);
for (int ith_inner_timestep = 0; ith_inner_timestep < num_inner_timesteps;
ith_inner_timestep++) {
for (int ixnode = 0; ixnode < nxnodes; ixnode++)
for (int iynode = 0; iynode < nynodes; iynode++)
for (int iznode = 0; iznode < nznodes; iznode++)
T_electron_old[ixnode][iynode][iznode] =
T_electron[ixnode][iynode][iznode];
// compute new electron T profile
duration = duration + inner_dt;
for (int ixnode = 0; ixnode < nxnodes; ixnode++)
for (int iynode = 0; iynode < nynodes; iynode++)
for (int iznode = 0; iznode < nznodes; iznode++) {
int right_xnode = ixnode + 1;
int right_ynode = iynode + 1;
int right_znode = iznode + 1;
if (right_xnode == nxnodes) right_xnode = 0;
if (right_ynode == nynodes) right_ynode = 0;
if (right_znode == nznodes) right_znode = 0;
int left_xnode = ixnode - 1;
int left_ynode = iynode - 1;
int left_znode = iznode - 1;
if (left_xnode == -1) left_xnode = nxnodes - 1;
if (left_ynode == -1) left_ynode = nynodes - 1;
if (left_znode == -1) left_znode = nznodes - 1;
double skin_layer_d = double(skin_layer);
double ixnode_d = double(ixnode);
double surface_d = double(t_surface_l);
mult_factor = 0.0;
if (duration < width){
if (ixnode >= t_surface_l) mult_factor = (intensity/(dx*skin_layer_d))*exp((-1.0)*(ixnode_d - surface_d)/skin_layer_d);
}
if (ixnode < t_surface_l) net_energy_transfer_all[ixnode][iynode][iznode] = 0.0;
double cr_vac = 1;
if (T_electron_old[ixnode][iynode][iznode] == 0) cr_vac = 0;
double cr_v_l_x = 1;
if (T_electron_old[left_xnode][iynode][iznode] == 0) cr_v_l_x = 0;
double cr_v_r_x = 1;
if (T_electron_old[right_xnode][iynode][iznode] == 0) cr_v_r_x = 0;
double cr_v_l_y = 1;
if (T_electron_old[ixnode][left_ynode][iznode] == 0) cr_v_l_y = 0;
double cr_v_r_y = 1;
if (T_electron_old[ixnode][right_ynode][iznode] == 0) cr_v_r_y = 0;
double cr_v_l_z = 1;
if (T_electron_old[ixnode][iynode][left_znode] == 0) cr_v_l_z = 0;
double cr_v_r_z = 1;
if (T_electron_old[ixnode][iynode][right_znode] == 0) cr_v_r_z = 0;
if (cr_vac != 0) {
T_electron[ixnode][iynode][iznode] =
T_electron_old[ixnode][iynode][iznode] +
inner_dt/el_properties(T_electron_old[ixnode][iynode][iznode]).el_heat_capacity *
((cr_v_r_x*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[right_xnode][iynode][iznode]/2.0).el_thermal_conductivity*
(T_electron_old[right_xnode][iynode][iznode]-T_electron_old[ixnode][iynode][iznode])/dx -
cr_v_l_x*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[left_xnode][iynode][iznode]/2.0).el_thermal_conductivity*
(T_electron_old[ixnode][iynode][iznode]-T_electron_old[left_xnode][iynode][iznode])/dx)/dx +
(cr_v_r_y*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[ixnode][right_ynode][iznode]/2.0).el_thermal_conductivity*
(T_electron_old[ixnode][right_ynode][iznode]-T_electron_old[ixnode][iynode][iznode])/dy -
cr_v_l_y*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[ixnode][left_ynode][iznode]/2.0).el_thermal_conductivity*
(T_electron_old[ixnode][iynode][iznode]-T_electron_old[ixnode][left_ynode][iznode])/dy)/dy +
(cr_v_r_z*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[ixnode][iynode][right_znode]/2.0).el_thermal_conductivity*
(T_electron_old[ixnode][iynode][right_znode]-T_electron_old[ixnode][iynode][iznode])/dz -
cr_v_l_z*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[ixnode][iynode][left_znode]/2.0).el_thermal_conductivity*
(T_electron_old[ixnode][iynode][iznode]-T_electron_old[ixnode][iynode][left_znode])/dz)/dz);
T_electron[ixnode][iynode][iznode]+=inner_dt/el_properties(T_electron[ixnode][iynode][iznode]).el_heat_capacity*
(mult_factor -
net_energy_transfer_all[ixnode][iynode][iznode]/del_vol);
}
else T_electron[ixnode][iynode][iznode] =
T_electron_old[ixnode][iynode][iznode];
if ((T_electron[ixnode][iynode][iznode] > 0.0) && (T_electron[ixnode][iynode][iznode] < electron_temperature_min))
T_electron[ixnode][iynode][iznode] = T_electron[ixnode][iynode][iznode] + 0.5*(electron_temperature_min - T_electron[ixnode][iynode][iznode]);
stability_criterion = 1.0 -
2.0*inner_dt/el_specific_heat *
(el_thermal_conductivity*(1.0/dx/dx + 1.0/dy/dy + 1.0/dz/dz));
if (stability_criterion < 0.0) {
inner_dt = 0.25*el_specific_heat /
(el_thermal_conductivity*(1.0/dx/dx + 1.0/dy/dy + 1.0/dz/dz));
}
num_inner_timesteps = static_cast<unsigned int>(update->dt/inner_dt) + 1;
inner_dt = update->dt/double(num_inner_timesteps);
if (num_inner_timesteps > 1000000)
error->warning(FLERR,"Too many inner timesteps in fix ttm",0);
for (int ith_inner_timestep = 0; ith_inner_timestep < num_inner_timesteps;
ith_inner_timestep++) {
for (int ixnode = 0; ixnode < nxnodes; ixnode++)
for (int iynode = 0; iynode < nynodes; iynode++)
for (int iznode = 0; iznode < nznodes; iznode++)
T_electron_old[ixnode][iynode][iznode] =
T_electron[ixnode][iynode][iznode];
// compute new electron T profile
duration = duration + inner_dt;
for (int ixnode = 0; ixnode < nxnodes; ixnode++)
for (int iynode = 0; iynode < nynodes; iynode++)
for (int iznode = 0; iznode < nznodes; iznode++) {
int right_xnode = ixnode + 1;
int right_ynode = iynode + 1;
int right_znode = iznode + 1;
if (right_xnode == nxnodes) right_xnode = 0;
if (right_ynode == nynodes) right_ynode = 0;
if (right_znode == nznodes) right_znode = 0;
int left_xnode = ixnode - 1;
int left_ynode = iynode - 1;
int left_znode = iznode - 1;
if (left_xnode == -1) left_xnode = nxnodes - 1;
if (left_ynode == -1) left_ynode = nynodes - 1;
if (left_znode == -1) left_znode = nznodes - 1;
double skin_layer_d = double(skin_layer);
double ixnode_d = double(ixnode);
double surface_d = double(t_surface_l);
mult_factor = 0.0;
if (duration < width){
if (ixnode >= t_surface_l) mult_factor = (intensity/(dx*skin_layer_d))*exp((-1.0)*(ixnode_d - surface_d)/skin_layer_d);
}
if (ixnode < t_surface_l) net_energy_transfer_all[ixnode][iynode][iznode] = 0.0;
double cr_vac = 1;
if (T_electron_old[ixnode][iynode][iznode] == 0) cr_vac = 0;
double cr_v_l_x = 1;
if (T_electron_old[left_xnode][iynode][iznode] == 0) cr_v_l_x = 0;
double cr_v_r_x = 1;
if (T_electron_old[right_xnode][iynode][iznode] == 0) cr_v_r_x = 0;
double cr_v_l_y = 1;
if (T_electron_old[ixnode][left_ynode][iznode] == 0) cr_v_l_y = 0;
double cr_v_r_y = 1;
if (T_electron_old[ixnode][right_ynode][iznode] == 0) cr_v_r_y = 0;
double cr_v_l_z = 1;
if (T_electron_old[ixnode][iynode][left_znode] == 0) cr_v_l_z = 0;
double cr_v_r_z = 1;
if (T_electron_old[ixnode][iynode][right_znode] == 0) cr_v_r_z = 0;
if (cr_vac != 0) {
T_electron[ixnode][iynode][iznode] =
T_electron_old[ixnode][iynode][iznode] +
inner_dt/el_properties(T_electron_old[ixnode][iynode][iznode]).el_heat_capacity *
((cr_v_r_x*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[right_xnode][iynode][iznode]/2.0).el_thermal_conductivity*
(T_electron_old[right_xnode][iynode][iznode]-T_electron_old[ixnode][iynode][iznode])/dx -
cr_v_l_x*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[left_xnode][iynode][iznode]/2.0).el_thermal_conductivity*
(T_electron_old[ixnode][iynode][iznode]-T_electron_old[left_xnode][iynode][iznode])/dx)/dx +
(cr_v_r_y*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[ixnode][right_ynode][iznode]/2.0).el_thermal_conductivity*
(T_electron_old[ixnode][right_ynode][iznode]-T_electron_old[ixnode][iynode][iznode])/dy -
cr_v_l_y*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[ixnode][left_ynode][iznode]/2.0).el_thermal_conductivity*
(T_electron_old[ixnode][iynode][iznode]-T_electron_old[ixnode][left_ynode][iznode])/dy)/dy +
(cr_v_r_z*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[ixnode][iynode][right_znode]/2.0).el_thermal_conductivity*
(T_electron_old[ixnode][iynode][right_znode]-T_electron_old[ixnode][iynode][iznode])/dz -
cr_v_l_z*el_properties(T_electron_old[ixnode][iynode][iznode]/2.0+T_electron_old[ixnode][iynode][left_znode]/2.0).el_thermal_conductivity*
(T_electron_old[ixnode][iynode][iznode]-T_electron_old[ixnode][iynode][left_znode])/dz)/dz);
T_electron[ixnode][iynode][iznode]+=inner_dt/el_properties(T_electron[ixnode][iynode][iznode]).el_heat_capacity*
(mult_factor -
net_energy_transfer_all[ixnode][iynode][iznode]/del_vol);
}
else T_electron[ixnode][iynode][iznode] =
T_electron_old[ixnode][iynode][iznode];
if ((T_electron[ixnode][iynode][iznode] > 0.0) && (T_electron[ixnode][iynode][iznode] < electron_temperature_min))
T_electron[ixnode][iynode][iznode] = T_electron[ixnode][iynode][iznode] + 0.5*(electron_temperature_min - T_electron[ixnode][iynode][iznode]);
if (el_properties(T_electron[ixnode][iynode][iznode]).el_thermal_conductivity > el_thermal_conductivity)
el_thermal_conductivity = el_properties(T_electron[ixnode][iynode][iznode]).el_thermal_conductivity;
if ((T_electron[ixnode][iynode][iznode] > 0.0) && (el_properties(T_electron[ixnode][iynode][iznode]).el_heat_capacity < el_specific_heat))
if (el_properties(T_electron[ixnode][iynode][iznode]).el_thermal_conductivity > el_thermal_conductivity)
el_thermal_conductivity = el_properties(T_electron[ixnode][iynode][iznode]).el_thermal_conductivity;
if ((T_electron[ixnode][iynode][iznode] > 0.0) && (el_properties(T_electron[ixnode][iynode][iznode]).el_heat_capacity < el_specific_heat))
el_specific_heat = el_properties(T_electron[ixnode][iynode][iznode]).el_heat_capacity;
}
}
stability_criterion = 1.0 -
2.0*inner_dt/el_specific_heat *
(el_thermal_conductivity*(1.0/dx/dx + 1.0/dy/dy + 1.0/dz/dz));
}
}
stability_criterion = 1.0 -
2.0*inner_dt/el_specific_heat *
(el_thermal_conductivity*(1.0/dx/dx + 1.0/dy/dy + 1.0/dz/dz));
} while (stability_criterion < 0.0);
} while (stability_criterion < 0.0);
// output nodal temperatures for current timestep
if ((nfileevery) && !(update->ntimestep % nfileevery)) {
// compute atomic Ta for each grid point

10
src/USER-REAXC/reaxc_tool_box.cpp
View File

@ -54,7 +54,7 @@ int Tokenize( char* s, char*** tok )
}
/* safe malloc */
void *smalloc( long n, const char *name, MPI_Comm comm )
void *smalloc( rc_bigint n, const char *name, MPI_Comm comm )
{
void *ptr;
@ -77,19 +77,19 @@ void *smalloc( long n, const char *name, MPI_Comm comm )
/* safe calloc */
void *scalloc( int n, int size, const char *name, MPI_Comm comm )
void *scalloc( rc_bigint n, rc_bigint size, const char *name, MPI_Comm comm )
{
void *ptr;
if( n <= 0 ) {
fprintf( stderr, "WARNING: trying to allocate %d elements for array %s. ",
fprintf( stderr, "WARNING: trying to allocate %ld elements for array %s. ",
n, name );
fprintf( stderr, "returning NULL.\n" );
return NULL;
}
if( size <= 0 ) {
fprintf( stderr, "WARNING: elements size for array %s is %d. ",
fprintf( stderr, "WARNING: elements size for array %s is %ld. ",
name, size );
fprintf( stderr, "returning NULL.\n" );
return NULL;
@ -97,7 +97,7 @@ void *scalloc( int n, int size, const char *name, MPI_Comm comm )
ptr = calloc( n, size );
if( ptr == NULL ) {
fprintf( stderr, "ERROR: failed to allocate %d bytes for array %s",
fprintf( stderr, "ERROR: failed to allocate %ld bytes for array %s",
n*size, name );
MPI_Abort( comm, INSUFFICIENT_MEMORY );
}

4
src/USER-REAXC/reaxc_tool_box.h
View File

@ -37,7 +37,7 @@ double Get_Time( );
int Tokenize( char*, char*** );
/* from lammps */
void *smalloc( long, const char*, MPI_Comm );
void *scalloc( int, int, const char*, MPI_Comm );
void *smalloc( rc_bigint, const char*, MPI_Comm );
void *scalloc( rc_bigint, rc_bigint, const char*, MPI_Comm );
void sfree( void*, const char* );
#endif

2
src/USER-REAXC/reaxc_traj.cpp
View File

@ -101,7 +101,7 @@ int Write_Header( reax_system *system, control_params *control,
out_control->traj_title );
strncat( out_control->buffer, out_control->line, HEADER_LINE_LEN+1 );
sprintf( out_control->line, INT_LINE, "number_of_atoms:", system->bigN );
sprintf( out_control->line, BIGINT_LINE, "number_of_atoms:", system->bigN );
strncat( out_control->buffer, out_control->line, HEADER_LINE_LEN+1 );
sprintf( out_control->line, STR_LINE, "ensemble_type:",

1
src/USER-REAXC/reaxc_traj.h
View File

@ -36,6 +36,7 @@
#define HEADER_LINE_LEN 62
#define STR_LINE "%-37s%-24s\n"
#define INT_LINE "%-37s%-24d\n"
#define BIGINT_LINE "%-37s%-24ld\n"
#define INT2_LINE "%-36s%-12d,%-12d\n"
#define REAL_LINE "%-37s%-24.3f\n"
#define SCI_LINE "%-37s%-24g\n"

6
src/USER-REAXC/reaxc_types.h
View File

@ -68,8 +68,9 @@ typedef double rvec2[2];
typedef double rvec4[4];
// import LAMMPS' definition of tagint
// import LAMMPS' definition of tagint and bigint
typedef LAMMPS_NS::tagint rc_tagint;
typedef LAMMPS_NS::bigint rc_bigint;
typedef struct
{
@ -383,7 +384,8 @@ struct _reax_system
{
reax_interaction reax_param;
int n, N, bigN, numH;
rc_bigint bigN;
int n, N, numH;
int local_cap, total_cap, gcell_cap, Hcap;
int est_recv, est_trans, max_recved;
int wsize, my_rank, num_nbrs;