fixes bug with cd_psi in model.cpp.

Bug with printed saturation indices with pressure.

prep: rewrote calc_delta_v(reaction *r_ptr, bool phase) to handle phase and species separately.
print: calculate pressure term for rxn before printing value.

git-svn-id: svn://136.177.114.72/svn_GW/phreeqc3/trunk@7547 1feff8c3-07ed-0310-ac33-dd36852eb9cd
This commit is contained in:
David L Parkhurst 2013-03-12 18:07:31 +00:00
parent 9849ed1aee
commit 6899bc4ba7
3 changed files with 34 additions and 12 deletions

View File

@ -4138,6 +4138,10 @@ residuals(void)
if (charge_ptr->Get_grams() == 0)
{
residual[i] = 0.0;
cd_psi.clear();
cd_psi.push_back(0.0);
cd_psi.push_back(0.0);
cd_psi.push_back(0.0);
}
else
{

View File

@ -5521,22 +5521,36 @@ calc_delta_v(reaction *r_ptr, bool phase)
{
/* calculate delta_v from molar volumes */
int p = -1;
LDBLE d_v = 0.0;
if (phase)
p = 1; /* for phases: reactants have coef's < 0, products have coef's > 0, v.v. for species */
for (size_t i = 0; r_ptr->token[i].name; i++)
{
if (!r_ptr->token[i].s)
continue;
if (!strcmp(r_ptr->token[i].s->name, "H+"))
continue;
if (!strcmp(r_ptr->token[i].s->name, "e-"))
continue;
else if (r_ptr->token[i].s->logk[vm_tc])
d_v += p * r_ptr->token[i].coef * r_ptr->token[i].s->logk[vm_tc];
/* for phases: reactants have coef's < 0, products have coef's > 0, v.v. for species */
for (size_t i = 1; r_ptr->token[i].s /*|| r_ptr->token[i].s*/ ; i++)
{
if (!r_ptr->token[i].s)
continue;
if (!strcmp(r_ptr->token[i].s->name, "H+"))
continue;
if (!strcmp(r_ptr->token[i].s->name, "e-"))
continue;
else if (r_ptr->token[i].s->logk[vm_tc])
d_v += r_ptr->token[i].coef * r_ptr->token[i].s->logk[vm_tc];
}
}
else
{
for (size_t i = 0; r_ptr->token[i].name /*|| r_ptr->token[i].s*/ ; i++)
{
if (!r_ptr->token[i].s)
continue;
if (!strcmp(r_ptr->token[i].s->name, "H+"))
continue;
if (!strcmp(r_ptr->token[i].s->name, "e-"))
continue;
else if (r_ptr->token[i].s->logk[vm_tc])
d_v += - r_ptr->token[i].coef * r_ptr->token[i].s->logk[vm_tc];
}
}
return d_v;
}

View File

@ -1283,6 +1283,10 @@ print_pp_assemblage(void)
PR_inprint = true;
phase_ptr->pr_in = false;
}
phase_ptr->rxn->logk[delta_v] = calc_delta_v(phase_ptr->rxn, true) -
phase_ptr->logk[vm0];
if (phase_ptr->rxn->logk[delta_v])
mu_terms_in_logk = true;
lk = k_calc(phase_ptr->rxn->logk, tk_x, patm_x * PASCAL_PER_ATM);
if (PR_inprint)
phase_ptr->pr_in = true;