max/iphreeqc
1
0
Fork 0
mirror of https://git.gfz-potsdam.de/naaice/iphreeqc.git synced 2025-12-16 08:38:23 +01:00
Code Issues Packages Projects Releases Wiki Activity

Moving class_main, phreeqc.cpp and phreeqc.h into phreeqcpp.sln

Browse Source 
Making a new configuration for compilation.

git-svn-id: svn://136.177.114.72/svn_GW/phreeqcpp/trunk@3850 1feff8c3-07ed-0310-ac33-dd36852eb9cd
This commit is contained in:
David L Parkhurst 2009-12-03 19:24:14 +00:00
parent 07d5fa7413
commit b33830cc5c
4 changed files with 2273 additions and 680 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

387
Phreeqc.cpp Normal file
View File

@ -0,0 +1,387 @@
#include "Phreeqc.h"
//class Phreeqc
//{
// Phreeqc(void);
// ~Phreeqc(void);
//};
#include <algorithm> // std::replace
Phreeqc::Phreeqc(void)
{
s_pTail = NULL;
user_database = NULL;
output = NULL;
log_file = NULL;
punch_file = NULL;
dump_file = NULL;
error_file = NULL;
rates = NULL;
tally_table = NULL;
spec = NULL;
cations = NULL;
anions = NULL;
neutrals = NULL;
IPRSNT = NULL;
M = NULL;
LGAMMA = NULL;
sit_params = NULL;
sit_LGAMMA = NULL;
sit_IPRSNT = NULL;
sit_M = NULL;
struct const_iso
{
const char *name;
LDBLE value;
LDBLE uncertainty;
};
struct const_iso temp_iso_defaults[] = {
{"13C", -10, 1},
{"13C(4)", -10, 1},
{"13C(-4)", -50, 5},
{"34S", 10, 1},
{"34S(6)", 10, 1},
{"34S(-2)", -30, 5},
{"2H", -28, 1},
{"18O", -5, .1},
{"87Sr", .71, .01},
{"11B", 20, 5}
};
int temp_count_iso_defaults =
(sizeof(temp_iso_defaults) / sizeof(struct const_iso));
count_iso_defaults = temp_count_iso_defaults;
iso_defaults = new iso[count_iso_defaults];
int i;
for (i = 0; i < temp_count_iso_defaults; i++)
{
iso_defaults[i].name = string_duplicate(temp_iso_defaults[i].name);
iso_defaults[i].value = temp_iso_defaults[i].value;
iso_defaults[i].uncertainty = temp_iso_defaults[i].uncertainty;
}
struct const_key keyword_temp[] = {
{"eof", 0},
{"end", 0},
{"solution_species", 0},
{"solution_master_species", 0},
{"solution", 0},
{"phases", 0},
{"pure_phases", 0},
{"reaction", 0},
{"mix", 0},
{"use", 0},
{"save", 0},
{"exchange_species", 0},
{"exchange_master_species", 0},
{"exchange", 0},
{"surface_species", 0},
{"surface_master_species", 0},
{"surface", 0},
{"reaction_temperature", 0},
{"inverse_modeling", 0},
{"gas_phase", 0},
{"transport", 0},
{"debug", 0},
{"selected_output", 0},
{"select_output", 0},
{"knobs", 0},
{"print", 0},
{"equilibrium_phases", 0},
{"equilibria", 0},
{"equilibrium", 0},
{"pure", 0},
{"title", 0},
{"comment", 0},
{"advection", 0},
{"kinetics", 0},
{"incremental_reactions", 0},
{"incremental", 0},
{"rates", 0},
{"solution_s", 0},
{"user_print", 0},
{"user_punch", 0},
{"solid_solutions", 0},
{"solid_solution", 0},
{"solution_spread", 0},
{"spread_solution", 0},
{"selected_out", 0},
{"select_out", 0},
{"user_graph", 0},
{"llnl_aqueous_model_parameters", 0},
{"llnl_aqueous_model", 0},
{"database", 0},
{"named_analytical_expression", 0},
{"named_analytical_expressions", 0},
{"named_expressions", 0},
{"named_log_k", 0},
{"isotopes", 0},
{"calculate_values", 0},
{"isotope_ratios", 0},
{"isotope_alphas", 0},
{"copy", 0},
{"pitzer", 0},
{"sit", 0}
#ifdef PHREEQC_CPP
,
{"solution_raw", 0},
{"exchange_raw", 0},
{"surface_raw", 0},
{"equilibrium_phases_raw", 0},
{"kinetics_raw", 0},
{"solid_solutions_raw", 0},
{"gas_phase_raw", 0},
{"reaction_raw", 0},
{"mix_raw", 0},
{"reaction_temperature_raw", 0},
{"dump", 0},
{"solution_modify", 0},
{"equilibrium_phases_modify", 0},
{"exchange_modify", 0},
{"surface_modify", 0},
{"solid_solutions_modify", 0},
{"gas_phase_modify", 0},
{"kinetics_modify", 0},
{"delete", 0},
{"run_cells", 0}
#endif
};
NKEYS = (sizeof(keyword_temp) / sizeof(struct const_key)); /* Number of valid keywords */
//keyword = (struct const_key *) PHRQ_malloc((size_t) (NKEYS * sizeof(const_key)));
keyword = new const_key[NKEYS];
for (i = 0; i < NKEYS; i++)
{
keyword[i].name = string_duplicate(keyword_temp[i].name);
keyword[i].keycount = 0;
}
// basic.c
struct const_key command_temp[] = {
{"and", tokand},
{"or", tokor},
{"xor", tokxor},
{"not", toknot},
{"mod", tokmod},
{"sqr", toksqr},
{"sqrt", toksqrt},
{"sin", toksin},
{"cos", tokcos},
{"tan", toktan},
{"arctan", tokarctan},
{"log", toklog},
{"exp", tokexp},
{"abs", tokabs},
{"sgn", toksgn},
{"str$", tokstr_},
{"val", tokval},
{"chr$", tokchr_},
{"asc", tokasc},
{"len", toklen},
{"mid$", tokmid_},
{"peek", tokpeek},
{"let", toklet},
{"print", tokprint},
{"punch", tokpunch},
#ifdef PHREEQ98
{"graph_x", tokgraph_x},
{"graph_y", tokgraph_y},
{"graph_sy", tokgraph_sy},
#endif
{"input", tokinput},
{"goto", tokgoto},
{"go to", tokgoto},
{"if", tokif},
{"end", tokend},
{"stop", tokstop},
{"for", tokfor},
{"next", toknext},
{"while", tokwhile},
{"wend", tokwend},
{"gosub", tokgosub},
{"return", tokreturn},
{"read", tokread},
{"data", tokdata},
{"restore", tokrestore},
{"gotoxy", tokgotoxy},
{"on", tokon},
{"dim", tokdim},
{"poke", tokpoke},
{"list", toklist},
{"run", tokrun},
{"new", toknew},
{"load", tokload},
{"merge", tokmerge},
{"save", toksave},
{"bye", tokbye},
{"quit", tokbye},
{"del", tokdel},
{"renum", tokrenum},
{"then", tokthen},
{"else", tokelse},
{"to", tokto},
{"step", tokstep},
{"tc", toktc},
{"tk", toktk},
{"time", toktime},
{"sim_time", toksim_time},
{"total_time", toktotal_time},
{"m0", tokm0},
{"m", tokm},
{"parm", tokparm},
{"act", tokact},
{"edl", tokedl},
{"surf", toksurf},
{"equi", tokequi},
{"kin", tokkin},
{"gas", tokgas},
{"s_s", toks_s},
{"misc1", tokmisc1},
{"misc2", tokmisc2},
{"mu", tokmu},
{"osmotic", tokosmotic},
{"alk", tokalk},
{"lk_species", toklk_species},
{"lk_named", toklk_named},
{"lk_phase", toklk_phase},
{"sum_species", toksum_species},
{"sum_gas", toksum_gas},
{"sum_s_s", toksum_s_s},
{"calc_value", tokcalc_value},
{"description", tokdescription},
{"sys", toksys},
{"instr", tokinstr},
{"ltrim", tokltrim},
{"rtrim", tokrtrim},
{"trim", toktrim},
{"pad", tokpad},
{"rxn", tokrxn},
{"dist", tokdist},
{"mol", tokmol},
{"la", tokla},
{"lm", toklm},
{"sr", toksr},
{"si", toksi},
{"step_no", tokstep_no},
{"cell_no", tokcell_no},
{"sim_no", toksim_no},
{"tot", toktot},
{"log10", toklog10},
{"charge_balance", tokcharge_balance},
{"percent_error", tokpercent_error},
{"put", tokput},
{"get", tokget},
{"exists", tokexists},
{"rem", tokrem},
{"change_por", tokchange_por},
{"get_por", tokget_por},
{"change_surf", tokchange_surf},
{"porevolume", tokporevolume},
{"sc", toksc},
{"gamma", tokgamma},
/* VP: Density Start */
{"lg", toklg},
{"rho", tokrho}
/* VP: Density End */
};
NCMDS = (sizeof(command_temp) / sizeof(struct const_key));
command = new const_key[NCMDS];
for (i = 0; i < NCMDS; i++)
{
command[i].name = string_duplicate(command_temp[i].name);
command[i].keycount = command_temp[i].keycount;
}
//cl1.c
x_arg = NULL, res_arg = NULL, scratch = NULL;
x_arg_max = 0, res_arg_max = 0, scratch_max = 0;
// dw.c
GASCON = 0.461522e0;
TZ = 647.073e0;
AA = 1.e0;
G1 = 11.e0;
G2 = 44.333333333333e0;
GF = 3.5e0;
// model.c
min_value = 1e-10;
normal = NULL;
ineq_array = NULL;
res = NULL;
cu = NULL;
zero = NULL;
delta1 = NULL;
iu = NULL;
is = NULL;
back_eq = NULL;
normal_max = 0;
ineq_array_max = 0;
res_max = 0;
cu_max = 0;
zero_max = 0;
delta1_max = 0;
iu_max = 0;
is_max = 0;
back_eq_max = 0;
// output.c
output_callbacks = new Phreeqc::output_callback[MAX_CALLBACKS];
count_output_callback = 0;
forward_output_to_log = 0;
// phqalloc.c
//s_pTail = NULL;
// transport.c
J_ij = NULL;
J_ij_il = NULL;
m_s = NULL;
default_data_base = string_duplicate("phreeqc.dat");
//FILE *input_file = NULL;
//FILE *database_file = NULL;
//FILE *output = NULL; /* OUTPUT_MESSAGE */
//FILE *log_file = NULL; /* OUTPUT_LOG */
//FILE *punch_file = NULL; /* OUTPUT_PUNCH */
//FILE *error_file = NULL; /* OUTPUT_ERROR */
//FILE *dump_file = NULL; /* OUTPUT_DUMP */
}
Phreeqc::~Phreeqc(void)
{
clean_up();
close_input_files();
close_output_files();
int i;
for (i = 0; i < count_iso_defaults; i++)
{
iso_defaults[i].name = (char *) free_check_null((void *) iso_defaults[i].name);
}
for (i = 0; i < NKEYS; i++)
{
keyword[i].name = (char *) free_check_null((void *) keyword[i].name);
}
delete[] keyword;
for (i = 0; i < NCMDS; i++)
{
command[i].name = (char *) free_check_null((void *) command[i].name);
}
delete[] command;
free_check_null(default_data_base);
PHRQ_free_all();
}
void Phreeqc::set_phast(int tf)
{
this->phast = tf;
}

1340
Phreeqc.h Normal file
View File

File diff suppressed because it is too large Load Diff

737
System.cxx
View File

@ -1,143 +1,78 @@
#include "System.h"
#include <algorithm> // std::replace
extern void ORCH_write_chemistry_species(std::ostream & chemistry_dat);
cxxSystem::cxxSystem(void)
{
this->solution = NULL;
this->exchange = NULL;
this->ppassemblage = NULL;
this->gasphase = NULL;
this->ssassemblage = NULL;
this->kinetics = NULL;
this->surface = NULL;
this->mix = NULL;
this->reaction = NULL;
this->temperature = NULL;
this->solution = NULL;
this->exchange = NULL;
this->ppassemblage = NULL;
this->gasphase = NULL;
this->ssassemblage = NULL;
this->kinetics = NULL;
this->surface = NULL;
this->mix = NULL;
this->reaction = NULL;
this->temperature = NULL;
}
cxxSystem::~cxxSystem(void)
{
}
void
cxxSystem::Initialize(void)
{
this->solution = NULL;
this->exchange = NULL;
this->ppassemblage = NULL;
this->gasphase = NULL;
this->ssassemblage = NULL;
this->kinetics = NULL;
this->surface = NULL;
this->mix = NULL;
this->reaction = NULL;
this->temperature = NULL;
this->solution = NULL;
this->exchange = NULL;
this->ppassemblage = NULL;
this->gasphase = NULL;
this->ssassemblage = NULL;
this->kinetics = NULL;
this->surface = NULL;
this->mix = NULL;
this->reaction = NULL;
this->temperature = NULL;
}
void
cxxSystem::totalize(void)
{
//initialize
this->totals.clear();
//add solution
if (this->solution != NULL)
{
char token[MAX_LENGTH];
strcpy(token, "O");
this->totals[token] = this->solution->get_total_o();
strcpy(token, "H");
this->totals[token] = this->solution->get_total_h();
strcpy(token, "Charge");
this->totals[token] = this->solution->get_cb();
this->totals.add_extensive(this->solution->get_totals(), 1.0);
}
if (this->exchange != NULL)
char token[MAX_LENGTH];
strcpy(token, "O");
this->totals[token] = this->solution->get_total_o();
strcpy(token, "H");
this->totals[token] = this->solution->get_total_h();
strcpy(token, "Charge");
this->totals[token] = this->solution->get_cb();
this->totals.add_extensive(this->solution->get_totals(), 1.0);
}
if (this->exchange != NULL)
{
this->exchange->totalize();
this->totals.add_extensive(this->exchange->get_totals(), 1.0);
}
if (this->ppassemblage != NULL)
this->exchange->totalize();
this->totals.add_extensive(this->exchange->get_totals(), 1.0);
}
if (this->ppassemblage != NULL)
{
this->ppassemblage->totalize();
this->totals.add_extensive(this->ppassemblage->get_totals(), 1.0);
}
if (this->gasphase != NULL)
this->ppassemblage->totalize();
this->totals.add_extensive(this->ppassemblage->get_totals(), 1.0);
}
if (this->gasphase != NULL)
{
this->gasphase->totalize();
this->totals.add_extensive(this->gasphase->get_totals(), 1.0);
}
if (this->ssassemblage != NULL)
this->gasphase->totalize();
this->totals.add_extensive(this->gasphase->get_totals(), 1.0);
}
if (this->ssassemblage != NULL)
{
this->ssassemblage->totalize();
this->totals.add_extensive(this->ssassemblage->get_totals(), 1.0);
}
if (this->surface != NULL)
this->ssassemblage->totalize();
this->totals.add_extensive(this->ssassemblage->get_totals(), 1.0);
}
if (this->surface != NULL)
{
this->ssassemblage->totalize();
this->totals.add_extensive(this->surface->get_totals(), 1.0);
}
this->ssassemblage->totalize();
this->totals.add_extensive(this->surface->get_totals(), 1.0);
}
//Need to handle the following 3 reactions:
//
//this->kinetics = NULL;
@ -146,362 +81,206 @@ this->totals.add_extensive(this->surface->get_totals(), 1.0);
//this->temperature = NULL;
//this->totals.dump_raw(std::cerr, 1);
return;
}
#ifdef ORCHESTRA
void
void
cxxSystem::ORCH_write(std::ostream & chemistry_dat, std::ostream & input_dat,
std::ostream & output_dat)
{
//
// chemistry file
//
this->totalize();
// calculate orchestra components from PHREEQC totals
this->ORCH_components();
// write definitions for chemistry
chemistry_dat << std::
endl << "@class: species_reactions () {" << std::endl;
this->ORCH_write_chemistry_water(chemistry_dat);
this->ORCH_write_chemistry_primary(chemistry_dat);
this->ORCH_write_chemistry_total_O_H(chemistry_dat);
ORCH_write_chemistry_species(chemistry_dat);
this->ORCH_write_chemistry_water(chemistry_dat);
this->ORCH_write_chemistry_primary(chemistry_dat);
this->ORCH_write_chemistry_total_O_H(chemistry_dat);
ORCH_write_chemistry_species(chemistry_dat);
// add definitions for pure phases
if (this->ppassemblage != NULL)
this->ppassemblage->ORCH_write_chemistry(chemistry_dat);
// add definitions for other PHREEQC entities here
// finish up
chemistry_dat << std::endl << "}" << std::endl;
//
// input file
//
input_dat << std::endl << "@class: input_file_data () {" << std::endl;
this->ORCH_write_input(input_dat);
input_dat << std::endl << "}" << std::endl;
this->ORCH_write_input(input_dat);
input_dat << std::endl << "}" << std::endl;
//
// output file
//
output_dat << std::endl << "Output_every: 1" << std::endl;
output_dat << "Var:";
this->ORCH_write_output_vars(output_dat);
output_dat << "Var:";
this->ORCH_write_output_vars(output_dat);
// add definitions for pure phases
if (this->ppassemblage != NULL)
this->ppassemblage->ORCH_write_output_vars(output_dat);
//finish up
output_dat << std::endl;
};
void
cxxSystem::ORCH_write_chemistry_water(std::ostream & chemistry_dat)
{
//
// Write water entities
//
chemistry_dat << std::
endl << "//********* The water entities" << std::endl;
// e-, total hydrogen
chemistry_dat << "@entity(e-, diss, 0)" << std::endl;
chemistry_dat << "@Calc: (1, \"e-.act = 10.^(-pe)\")" << std::endl;
chemistry_dat << "@solve (pe, 1e-6, lin, 1, e-.liter, 1e-14)" << std::
chemistry_dat << "@Calc: (1, \"e-.act = 10.^(-pe)\")" << std::endl;
chemistry_dat << "@solve (pe, 1e-6, lin, 1, e-.liter, 1e-14)" << std::
endl;
// H+, charge balance
chemistry_dat << "@Calc: (1, \"H+.act = 10.^(-pH)\")" << std::endl;
chemistry_dat << "@solve (pH, 1e-6, lin, 1, H+.liter, 1e-14)" << std::
chemistry_dat << "@solve (pH, 1e-6, lin, 1, H+.liter, 1e-14)" << std::
endl;
// H2O
chemistry_dat << "@entity(" << s_h2o->
name << ", diss, 55.506)" << std::endl;
chemistry_dat << std::endl;
chemistry_dat << std::endl;
}
void
cxxSystem::ORCH_write_chemistry_primary(std::ostream & chemistry_dat)
{
chemistry_dat << std::endl << "//********* The primary species" << std::
chemistry_dat << std::endl << "//********* The primary species" << std::
endl;
//
// Write other master species definitions, i.e. primary entities
//
chemistry_dat << "@species(H+, 1)" << std::endl;
for (cxxNameDouble::iterator iter = this->totals.begin();
for (cxxNameDouble::iterator iter = this->totals.begin();
iter != this->totals.end(); ++iter)
{
std::string name(iter->first);
if (name == "H(1)" || name == "E" || name == "H" || name == "O"
std::string name(iter->first);
if (name == "H(1)" || name == "E" || name == "H" || name == "O"
|| name == "Charge")
continue;
struct element *elt;
char *element_name = string_hsave(name.c_str());
elt = element_store(element_name);
assert(elt != NULL);
struct species *s_ptr;
s_ptr = elt->master->s;
assert(s_ptr != NULL);
chemistry_dat << "@species(" << s_ptr->name << ", " << s_ptr->
struct element *elt;
char *element_name = string_hsave(name.c_str());
elt = element_store(element_name);
assert(elt != NULL);
struct species *s_ptr;
s_ptr = elt->master->s;
assert(s_ptr != NULL);
chemistry_dat << "@species(" << s_ptr->name << ", " << s_ptr->
z << ")" << std::endl;
// regular mole balance
chemistry_dat << "@primary_entity(" << s_ptr->
name << ", 1e-9, liter, 1.0e-9)" << std::endl;
}
chemistry_dat << std::endl;
}
chemistry_dat << std::endl;
}
void
cxxSystem::ORCH_write_chemistry_total_O_H(std::ostream & chemistry_dat)
{
chemistry_dat << std::
chemistry_dat << std::
endl << "//********* Define total hydrogen and oxygen" << std::endl;
// Write H total equation
chemistry_dat << "@var: total_hydrogen 0" << std::endl;
chemistry_dat << "@calc: (5, \"total_hydrogen = 0";
int i;
for (i = 0; i < count_s_x; i++)
chemistry_dat << "@calc: (5, \"total_hydrogen = 0";
int i;
for (i = 0; i < count_s_x; i++)
{
// write in terms of orchestra components
if (s_x[i]->primary != NULL
|| (s_x[i]->secondary != NULL
&& s_x[i]->secondary->in == TRUE))
{
if (s_x[i]->h != 0)
if (s_x[i]->h != 0)
{
chemistry_dat << "+";
if (s_x[i]->h != 1)
chemistry_dat << "+";
if (s_x[i]->h != 1)
{
chemistry_dat << s_x[i]->h << "*";
}
chemistry_dat << "{" << s_x[i]->name << ".liter}";
}
}
}
chemistry_dat << "\")" << std::endl;
chemistry_dat << s_x[i]->h << "*";
}
chemistry_dat << "{" << s_x[i]->name << ".liter}";
}
}
}
chemistry_dat << "\")" << std::endl;
// Write O total equation
chemistry_dat << "@var: total_oxygen 0" << std::endl;
chemistry_dat << "@calc: (5, \"total_oxygen = 0";
for (i = 0; i < count_s_x; i++)
chemistry_dat << "@calc: (5, \"total_oxygen = 0";
for (i = 0; i < count_s_x; i++)
{
if (s_x[i]->o != 0)
if (s_x[i]->o != 0)
{
// write in terms of orchestra components
if (s_x[i]->primary != NULL
|| (s_x[i]->secondary != NULL
&& s_x[i]->secondary->in == TRUE))
{
chemistry_dat << "+";
if (s_x[i]->o != 1)
chemistry_dat << "+";
if (s_x[i]->o != 1)
{
chemistry_dat << s_x[i]->o << "*";
chemistry_dat << s_x[i]->o << "*";
}
chemistry_dat << "{" << s_x[i]->name << ".liter}";
}
}
}
chemistry_dat << "\")" << std::endl;
chemistry_dat << std::endl;
}
chemistry_dat << "{" << s_x[i]->name << ".liter}";
}
}
}
chemistry_dat << "\")" << std::endl;
chemistry_dat << std::endl;
}
void
cxxSystem::ORCH_write_input(std::ostream & input_dat)
{
//
// Write orchestra input file info
//
//
// Write orchestra input file info
//
std::ostringstream headings, data;
data.precision(DBL_DIG - 1);
headings << "var: ";
data << "data: ";
headings << "tempc\t";
data << this->solution->get_tc() << "\t";
headings << "pH\t";
data << this->solution->get_ph() << "\t";
headings << "pe\t";
data << this->solution->get_pe() << "\t";
headings << "H2O.act\t";
data << 1 << "\t";
headings << "I\t";
data << this->solution->get_mu() << "\t";
for (cxxNameDouble::iterator iter = this->orch_totals.begin();
data.precision(DBL_DIG - 1);
headings << "var: ";
data << "data: ";
headings << "tempc\t";
data << this->solution->get_tc() << "\t";
headings << "pH\t";
data << this->solution->get_ph() << "\t";
headings << "pe\t";
data << this->solution->get_pe() << "\t";
headings << "H2O.act\t";
data << 1 << "\t";
headings << "I\t";
data << this->solution->get_mu() << "\t";
for (cxxNameDouble::iterator iter = this->orch_totals.begin();
iter != this->orch_totals.end(); ++iter)
{
headings << iter->first << ".liter" << "\t";
data << iter->second << "\t";
}
headings << iter->first << ".liter" << "\t";
data << iter->second << "\t";
}
// activity estimate
for (cxxNameDouble::iterator iter = this->totals.begin();
iter != this->totals.end(); ++iter)
{
if (iter->first == "O" || iter->first == "Charge")
if (iter->first == "O" || iter->first == "Charge")
continue;
std::string master_name;
struct master *master_ptr;
master_ptr = master_bsearch(iter->first);
assert(master_ptr != NULL);
cxxNameDouble ma = this->solution->get_master_activity();
cxxNameDouble::iterator it = ma.find(iter->first);
if (it == ma.end())
std::string master_name;
struct master *master_ptr;
master_ptr = master_bsearch(iter->first);
assert(master_ptr != NULL);
cxxNameDouble ma = this->solution->get_master_activity();
cxxNameDouble::iterator it = ma.find(iter->first);
if (it == ma.end())
{
it = ma.find(master_ptr->s->secondary->elt->name);
}
headings << master_ptr->s->name << ".act\t";
if (it != ma.end())
it = ma.find(master_ptr->s->secondary->elt->name);
}
headings << master_ptr->s->name << ".act\t";
if (it != ma.end())
{
data << pow(10., it->second) << "\t";
}
data << pow(10., it->second) << "\t";
}
else
{
data << 1e-9 << "\t";
}
}
data << 1e-9 << "\t";
}
}
// Isotopes
//s_oss << "-Isotopes" << std::endl;
/*
@ -511,260 +290,156 @@ data << 1e-9 << "\t";
}
}
*/
// Write data to string
input_dat << headings.str() << std::endl;
input_dat << data.str() << std::endl;
return;
input_dat << data.str() << std::endl;
return;
}
void
cxxSystem::ORCH_components(void)
{
// translate from H, O, charge to H+tot, e-tot, and H2Otot
cxxNameDouble::iterator it;
cxxNameDouble temp_totals;
cxxNameDouble temp_totals;
//
// Set names in orch_totals to master species names
//
this->orch_totals.clear();
temp_totals = this->totals;
temp_totals.erase("H");
temp_totals.erase("O");
temp_totals.erase("Charge");
temp_totals = this->totals;
temp_totals.erase("H");
temp_totals.erase("O");
temp_totals.erase("Charge");
//temp_totals.dump_raw(std::cerr, 1);
for (it = temp_totals.begin(); it != temp_totals.end(); it++)
{
struct element *elt_ptr;
elt_ptr = element_store(it->first);
assert(elt_ptr != NULL);
struct master *master_ptr;
master_ptr = elt_ptr->primary;
assert(master_ptr != NULL);
double coef = master_ptr->coef;
if (master_ptr->coef == 0)
struct element *elt_ptr;
elt_ptr = element_store(it->first);
assert(elt_ptr != NULL);
struct master *master_ptr;
master_ptr = elt_ptr->primary;
assert(master_ptr != NULL);
double coef = master_ptr->coef;
if (master_ptr->coef == 0)
{
coef = 1;
}
this->orch_totals[master_ptr->s->name] = it->second / coef;
}
coef = 1;
}
this->orch_totals[master_ptr->s->name] = it->second / coef;
}
//
// Calculate h2otot
//
it = this->totals.find("O");
assert(it != this->totals.end());
double h2otot = it->second;
assert(it != this->totals.end());
double h2otot = it->second;
// subtract O in master species
for (it = temp_totals.begin(); it != temp_totals.end(); it++)
{
struct element *elt_ptr;
elt_ptr = element_store(it->first);
struct master *master_ptr;
master_ptr = elt_ptr->primary;
h2otot -= it->second * master_ptr->s->o;
}
struct element *elt_ptr;
elt_ptr = element_store(it->first);
struct master *master_ptr;
master_ptr = elt_ptr->primary;
h2otot -= it->second * master_ptr->s->o;
}
//
// Calculate htot
//
it = this->totals.find("H");
assert(it != this->totals.end());
double htot = it->second - 2 * h2otot;
assert(it != this->totals.end());
double htot = it->second - 2 * h2otot;
// subtract O in master species
for (it = temp_totals.begin(); it != temp_totals.end(); it++)
{
struct element *elt_ptr;
elt_ptr = element_store(it->first);
struct master *master_ptr;
master_ptr = elt_ptr->primary;
htot -= it->second * master_ptr->s->h;
}
struct element *elt_ptr;
elt_ptr = element_store(it->first);
struct master *master_ptr;
master_ptr = elt_ptr->primary;
htot -= it->second * master_ptr->s->h;
}
//
// Calculate etot
//
it = this->totals.find("Charge");
assert(it != this->totals.end());
double etot = -it->second + htot;
assert(it != this->totals.end());
double etot = -it->second + htot;
// subtract O in master species
for (it = temp_totals.begin(); it != temp_totals.end(); it++)
{
struct element *elt_ptr;
elt_ptr = element_store(it->first);
struct master *master_ptr;
master_ptr = elt_ptr->primary;
etot += it->second * master_ptr->s->z;
}
struct element *elt_ptr;
elt_ptr = element_store(it->first);
struct master *master_ptr;
master_ptr = elt_ptr->primary;
etot += it->second * master_ptr->s->z;
}
//
// store h2otot, htot, etot in orch_totals
//
this->orch_totals["H2O"] = h2otot;
this->orch_totals["H+"] = htot;
this->orch_totals["e-"] = etot;
this->orch_totals.dump_raw(std::cerr, 1);
this->orch_totals["H+"] = htot;
this->orch_totals["e-"] = etot;
this->orch_totals.dump_raw(std::cerr, 1);
}
void
cxxSystem::ORCH_write_output_vars(std::ostream & outstream)
{
outstream << "\tnr_iter";
outstream << "\tnr_iter";
//
// Serialize solution
//
outstream << "\tstart_solution";
//tc
outstream << "\ttempc";
//pH
outstream << "\tpH";
//pe
outstream << "\tpe";
//mu
outstream << "\tI";
//ah2o
outstream << "\tH2O.act";
//total_h;
outstream << "\ttotal_hydrogen";
//total_o;
outstream << "\ttotal_oxygen";
//cb
outstream << "\tchargebalance";
//mass_water;
outstream << "\tH2O.con";
//total_alkalinity;
outstream << "\tAlkalinity";
//orchestra master variables
outstream << "\tH+.diss";
outstream << "\te-.diss";
outstream << "\tH2O.diss";
outstream << "\te-.diss";
outstream << "\tH2O.diss";
//
// Write totals
for (cxxNameDouble::iterator it = this->orch_totals.begin();
it != this->orch_totals.end(); it++)
{
if (it->first == "H+" || it->first == "H2O" || it->first == "e-")
if (it->first == "H+" || it->first == "H2O" || it->first == "e-")
continue;
outstream << "\t" << it->first << ".diss";
}
outstream << "\tend_totals";
for (cxxNameDouble::iterator it = this->orch_totals.begin();
outstream << "\t" << it->first << ".diss";
}
outstream << "\tend_totals";
for (cxxNameDouble::iterator it = this->orch_totals.begin();
it != this->orch_totals.end(); it++)
{
if (it->first == "H+" || it->first == "H2O" || it->first == "e-")
if (it->first == "H+" || it->first == "H2O" || it->first == "e-")
continue;
outstream << "\t" << it->first << ".act";
}
outstream << "\tend_master_activities";
outstream << "\t" << it->first << ".act";
}
outstream << "\tend_master_activities";
//
// Write all species activities and concentrations
//
int i;
for (i = 0; i < count_s_x; i++)
for (i = 0; i < count_s_x; i++)
{
std::string name(s_x[i]->name);
std::replace(name.begin(), name.end(), '(', '[');
std::replace(name.begin(), name.end(), ')', ']');
outstream << "\t" << name.c_str() << ".act" << "\t" << name.
std::string name(s_x[i]->name);
std::replace(name.begin(), name.end(), '(', '[');
std::replace(name.begin(), name.end(), ')', ']');
outstream << "\t" << name.c_str() << ".act" << "\t" << name.
c_str() << ".con";
}
outstream << "\tend_species";
}
outstream << "\tend_species";
}
#endif /*
*/
#endif /* */

191
class_main.cpp Normal file
View File

@ -0,0 +1,191 @@
#include "Phreeqc.h"
#include "output.h"
#include "phrqproto.h"
#include "input.h"
/* ----------------------------------------------------------------------
* MAIN
* ---------------------------------------------------------------------- */
int
main(int argc, char *argv[])
/*
* Main program for PHREEQC
*/
{
//int errors;
#if defined(WIN32_MEMORY_DEBUG)
int tmpDbgFlag;
/*
* Set the debug-heap flag to keep freed blocks in the
* heap's linked list - This will allow us to catch any
* inadvertent use of freed memory
*/
tmpDbgFlag = _CrtSetDbgFlag(_CRTDBG_REPORT_FLAG);
//tmpDbgFlag |= _CRTDBG_DELAY_FREE_MEM_DF;
tmpDbgFlag |= _CRTDBG_LEAK_CHECK_DF;
///tmpDbgFlag |= _CRTDBG_CHECK_ALWAYS_DF;
_CrtSetDbgFlag(tmpDbgFlag);
//_crtBreakAlloc = 9482;
#endif
/*if (svnid == NULL)
fprintf(stderr, " ");*/
Phreeqc phreeqc_instance;
phreeqc_instance.main_method(argc, argv);
}
int CLASS_QUALIFIER
main_method(int argc, char *argv[])
/*
* Main program for PHREEQC
*/
{
int errors;
void *db_cookie = NULL;
void *input_cookie = NULL;
#if defined(WIN32_MEMORY_DEBUG)
int tmpDbgFlag;
/*
* Set the debug-heap flag to keep freed blocks in the
* heap's linked list - This will allow us to catch any
* inadvertent use of freed memory
*/
tmpDbgFlag = _CrtSetDbgFlag(_CRTDBG_REPORT_FLAG);
//tmpDbgFlag |= _CRTDBG_DELAY_FREE_MEM_DF;
tmpDbgFlag |= _CRTDBG_LEAK_CHECK_DF;
///tmpDbgFlag |= _CRTDBG_CHECK_ALWAYS_DF;
_CrtSetDbgFlag(tmpDbgFlag);
//_crtBreakAlloc = 9482;
#endif
//if (svnid == NULL)
// fprintf(stderr, " ");
phast = FALSE;
/*
* Add callbacks for error_msg and warning_msg
*/
if (add_output_callback(phreeqc_handler, this) != OK)
{
fprintf(stderr, "ERROR: %s\n",
"NULL pointer returned from malloc or realloc.");
fprintf(stderr, "ERROR: %s\n", "Program terminating.");
return -1;
}
/*
* Open input/output files
*/
errors = process_file_names(argc, argv, &db_cookie, &input_cookie, TRUE);
if (errors != 0)
{
//clean_up();
return errors;
}
#ifdef DOS
write_banner();
#endif
/*
* Initialize arrays
*/
errors = do_initialize();
if (errors != 0)
{
//clean_up();
return errors;
}
/*
* Load database into memory
*/
errors = read_database(getc_callback, db_cookie);
if (errors != 0)
{
//clean_up();
return errors;
}
/*
* Read input data for simulation
*/
errors = run_simulations(getc_callback, input_cookie);
if (errors != 0)
{
//clean_up();
return errors;
}
/*
* Display successful status
*/
errors = do_status();
if (errors != 0)
{
//clean_up();
return errors;
}
//clean_up();
//close_input_files();
//close_output_files();
return 0;
}
#ifdef DOS
/* ---------------------------------------------------------------------- */
int CLASS_QUALIFIER
write_banner(void)
/* ---------------------------------------------------------------------- */
{
char buffer[80];
int len, indent;
output_msg(OUTPUT_SCREEN,
" ÛßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßÛ\n");
output_msg(OUTPUT_SCREEN,
" º º\n");
/* version */
len = sprintf(buffer, "* PHREEQC-%s *", "@VERSION@");
indent = (44 - len) / 2;
output_msg(OUTPUT_SCREEN, "%14cº%*c%s%*cº\n", ' ', indent, ' ', buffer,
44 - indent - len, ' ');
output_msg(OUTPUT_SCREEN,
" º º\n");
output_msg(OUTPUT_SCREEN,
" º A hydrogeochemical transport model º\n");
output_msg(OUTPUT_SCREEN,
" º º\n");
output_msg(OUTPUT_SCREEN,
" º by º\n");
output_msg(OUTPUT_SCREEN,
" º D.L. Parkhurst and C.A.J. Appelo º\n");
output_msg(OUTPUT_SCREEN,
" º º\n");
/* date */
len = sprintf(buffer, "%s", "@VER_DATE@");
indent = (44 - len) / 2;
output_msg(OUTPUT_SCREEN, "%14cº%*c%s%*cº\n", ' ', indent, ' ', buffer,
44 - indent - len, ' ');
output_msg(OUTPUT_SCREEN,
" ÛÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÛ\n\n");
return 0;
}
#endif