Skip to content
Snippets Groups Projects
Commit d496c0fd authored by sjplimp's avatar sjplimp Committed by GitHub
Browse files

Merge pull request #426 from dstelter92/master

fix for temper_grem exchange probability
parents 4a90bca7 5c39dfd7
No related branches found
No related tags found
No related merge requests found
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
...@@ -283,8 +283,8 @@ void TemperGrem::command(int narg, char **arg) ...@@ -283,8 +283,8 @@ void TemperGrem::command(int narg, char **arg)
// compute weights // compute weights
volume = domain->xprd * domain->yprd * domain->zprd; volume = domain->xprd * domain->yprd * domain->zprd;
enth = pe + (pressref * volume); enth = pe + (pressref * volume);
weight = log(set_lambda[my_set_lambda] + (eta*(enth + h0))); weight = log(set_lambda[my_set_lambda] + (eta*(enth - h0)));
weight_cross = log(set_lambda[partner_set_lambda] + (eta*(enth + h0))); weight_cross = log(set_lambda[partner_set_lambda] + (eta*(enth - h0)));
if (me_universe > partner) { if (me_universe > partner) {
MPI_Send(&weight,1,MPI_DOUBLE,partner,0,universe->uworld); MPI_Send(&weight,1,MPI_DOUBLE,partner,0,universe->uworld);
...@@ -296,7 +296,7 @@ void TemperGrem::command(int narg, char **arg) ...@@ -296,7 +296,7 @@ void TemperGrem::command(int narg, char **arg)
} }
if (me_universe < partner) { if (me_universe < partner) {
boltz_factor = (weight - weight_partner + weight_cross - weight_cross_partner) * boltz_factor = (weight + weight_partner - weight_cross - weight_cross_partner) *
(1 / (boltz * eta)); (1 / (boltz * eta));
if (boltz_factor >= 0.0) swap = 1; if (boltz_factor >= 0.0) swap = 1;
else if (ranboltz->uniform() < exp(boltz_factor)) swap = 1; else if (ranboltz->uniform() < exp(boltz_factor)) swap = 1;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment