changed print.cpp for Molar volume when not PR

Revised Tony organization of tidy_gas_phase.

git-svn-id: svn://136.177.114.72/svn_GW/phreeqc3/trunk@7935 1feff8c3-07ed-0310-ac33-dd36852eb9cd
This commit is contained in:
David L Parkhurst 2013-08-14 19:06:36 +00:00
parent 25a990751b
commit ec4fa77cba
2 changed files with 176 additions and 1 deletions

View File

@ -577,8 +577,18 @@ print_gas_phase(void)
output_msg(sformatf(" Gas volume: %10.2e liters\n",
(double) gas_phase_ptr->Get_volume()));
if(gas_phase_ptr->Get_total_moles() > 0)
output_msg(sformatf(" Molar volume: %10.2e liters/mole",
{
if (PR)
{
output_msg(sformatf(" Molar volume: %10.2e liters/mole",
(double) (gas_phase_ptr->Get_v_m())));
}
else
{
output_msg(sformatf(" Molar volume: %10.2e liters/mole",
(double) (gas_phase_ptr->Get_volume() / gas_phase_ptr->Get_total_moles())));
}
}
if (!numerical_fixed_volume && ((PR && gas_phase_ptr->Get_v_m() <= 0.016)))
output_msg(" WARNING: Program limit for Peng-Robinson.\n");
else

165
tidy.cpp
View File

@ -858,6 +858,170 @@ rewrite_eqn_to_primary(void)
int Phreeqc::
tidy_gas_phase(void)
/* ---------------------------------------------------------------------- */
{
int n_user, last;
LDBLE P, V_m;
bool PR;
/*
* Find all gases for each gas_phase in phase list
*/
for (std::set<int>::const_iterator nit = Rxn_new_gas_phase.begin(); nit != Rxn_new_gas_phase.end(); nit++)
{
std::map<int, cxxGasPhase>::iterator it = Rxn_gas_phase_map.find(*nit);
if (it == Rxn_gas_phase_map.end())
{
assert(false);
}
cxxGasPhase *gas_phase_ptr = &(it->second);
PR = false;
P = 0.0;
for (size_t j = 0; j < gas_phase_ptr->Get_gas_comps().size(); j++)
{
int k;
struct phase *phase_ptr = phase_bsearch(gas_phase_ptr->Get_gas_comps()[j].Get_phase_name().c_str(), &k, FALSE);
if (phase_ptr == NULL)
{
input_error++;
error_string = sformatf(
"Gas not found in PHASES database, %s.",
gas_phase_ptr->Get_gas_comps()[j].Get_phase_name().c_str());
error_msg(error_string, CONTINUE);
continue;
}
else
{
if (phase_ptr->t_c > 0 && phase_ptr->p_c > 0)
PR = true;
}
}
gas_phase_ptr->Set_pr_in(PR);
if (gas_phase_ptr->Get_new_def())
{
for (size_t j = 0; j < gas_phase_ptr->Get_gas_comps().size(); j++)
{
/*
* Fixed pressure
*/
if (gas_phase_ptr->Get_type() == cxxGasPhase::GP_PRESSURE)
{
if (gas_phase_ptr->Get_solution_equilibria())
{
input_error++;
error_string = sformatf(
"Gas phase %d: cannot use '-equilibrium' option with fixed pressure gas phase.",
gas_phase_ptr->Get_n_user());
error_msg(error_string, CONTINUE);
}
/* calculate moles */
if (gas_phase_ptr->Get_gas_comps()[j].Get_p_read() != NAN)
{
P += gas_phase_ptr->Get_gas_comps()[j].Get_p_read();
if (!PR)
gas_phase_ptr->Get_gas_comps()[j].Set_moles(
gas_phase_ptr->Get_gas_comps()[j].Get_p_read() * gas_phase_ptr->Get_volume() /
R_LITER_ATM / gas_phase_ptr->Get_temperature());
}
else
{
input_error++;
error_string = sformatf(
"Gas phase %d: partial pressure of gas component %s not defined.",
gas_phase_ptr->Get_n_user(), gas_phase_ptr->Get_gas_comps()[j].Get_phase_name().c_str());
error_msg(error_string, CONTINUE);
}
}
else
{
/*
* Fixed volume
*/
if (!gas_phase_ptr->Get_solution_equilibria())
{
if (gas_phase_ptr->Get_gas_comps()[j].Get_p_read() != NAN)
{
P += gas_phase_ptr->Get_gas_comps()[j].Get_p_read();
if (!PR)
gas_phase_ptr->Get_gas_comps()[j].Set_moles (
gas_phase_ptr->Get_gas_comps()[j].Get_p_read() *
gas_phase_ptr->Get_volume() / R_LITER_ATM /
gas_phase_ptr->Get_temperature());
}
else
{
input_error++;
error_string = sformatf(
"Gas phase %d: moles of gas component %s not defined.",
gas_phase_ptr->Get_n_user(),
gas_phase_ptr->Get_gas_comps()[j].Get_phase_name().c_str());
error_msg(error_string, CONTINUE);
}
}
}
}
if (PR && P > 0)
{
std::vector<struct phase *> phase_ptrs;
size_t j_PR;
std::vector<cxxGasComp> &gc = gas_phase_ptr->Get_gas_comps();
for (j_PR = 0; j_PR < gas_phase_ptr->Get_gas_comps().size(); j_PR++)
{
int k;
struct phase *phase_ptr = phase_bsearch(gas_phase_ptr->Get_gas_comps()[j_PR].Get_phase_name().c_str(), &k, FALSE);
if (gc[j_PR].Get_p_read() == 0)
continue;
phase_ptr->moles_x = gc[j_PR].Get_p_read() / P;
phase_ptrs.push_back(phase_ptr);
}
V_m = calc_PR(phase_ptrs, P, gas_phase_ptr->Get_temperature(), 0);
gas_phase_ptr->Set_v_m(V_m);
if (gas_phase_ptr->Get_type() == cxxGasPhase::GP_VOLUME)
{
gas_phase_ptr->Set_total_p(P);
}
for (j_PR = 0; j_PR < gas_phase_ptr->Get_gas_comps().size(); j_PR++)
{
int k;
struct phase *phase_ptr = phase_bsearch(gc[j_PR].Get_phase_name().c_str(), &k, FALSE);
if (gc[j_PR].Get_p_read() == 0)
{
gc[j_PR].Set_moles(0.0);
} else
{
gc[j_PR].Set_moles(phase_ptr->moles_x *
gas_phase_ptr->Get_volume() / V_m);
gas_phase_ptr->Set_total_moles(gas_phase_ptr->Get_total_moles() + gc[j_PR].Get_moles());
}
}
}
/*
* Duplicate gas phase, only if not solution equilibria
*/
if (!gas_phase_ptr->Get_solution_equilibria())
{
gas_phase_ptr->Set_new_def(false);
n_user = gas_phase_ptr->Get_n_user();
last = gas_phase_ptr->Get_n_user_end();
gas_phase_ptr->Set_n_user_end(n_user);
for (int j1 = n_user + 1; j1 <= last; j1++)
{
Utilities::Rxn_copy(Rxn_gas_phase_map, n_user, j1);
}
}
else
{
gas_phase_ptr->Set_new_def(true);
}
}
}
return (OK);
}
#ifdef SKIP
/* ---------------------------------------------------------------------- */
int Phreeqc::
tidy_gas_phase(void)
/* ---------------------------------------------------------------------- */
{
int n_user, last;
LDBLE P, V_m;
@ -1019,6 +1183,7 @@ tidy_gas_phase(void)
}
return (OK);
}
#endif
/* ---------------------------------------------------------------------- */
int Phreeqc::
tidy_inverse(void)