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

git-svn-id: svn://136.177.114.72/svn_GW/phreeqcpp/trunk@1723 1feff8c3-07ed-0310-ac33-dd36852eb9cd

Browse Source
This commit is contained in:
David L Parkhurst 2007-01-30 01:34:18 +00:00
parent 93a7759e1b
commit 2b442eef30
9 changed files with 107 additions and 286 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
ExchComp.cxx
View File

@ -102,7 +102,7 @@ cxxExchComp::cxxExchComp(std::vector<cxxExchComp> &ec_vector, std::vector<double
this->la += it_ec->la*intensive;
this->charge_balance += it_ec->charge_balance*extensive;
this->phase_proportion += it_ec->phase_proportion*intensive;
this->totals.add_extensive(it_ec->totals, extensive);
this->totals.add(it_ec->totals, extensive);
it_ec++;
it_f++;
}

14
Exchange.h
View File

@ -1,15 +1,17 @@
#if !defined(EXCHANGE_H_INCLUDED)
#define EXCHANGE_H_INCLUDED
#include "NumKeyword.h"
#define EXTERNAL extern
#include "global.h"
#include <cassert> // assert
#include <map> // std::map
#include <string> // std::string
#include <list> // std::list
#include <list> // std::list
#include <vector> // std::vector
#include "char_star.h"
#include "ExchComp.h"
#include "cxxMix.h"
class cxxExchange : public cxxNumKeyword
{
@ -17,13 +19,14 @@ class cxxExchange : public cxxNumKeyword
public:
cxxExchange();
cxxExchange(struct exchange *);
cxxExchange(const std::map<int, cxxExchange> &exchange_map, cxxMix &mx);
~cxxExchange();
struct exchange *cxxExchange2exchange();
struct exch_comp *cxxExchComp2exch_comp();
void dump_xml(std::ostream& os, unsigned int indent = 0)const;
void dump_raw(std::ostream& s_oss, unsigned int indent)const;
void read_raw(CParser& parser);
@ -47,11 +50,6 @@ public:
void mpi_pack(std::vector<int>& ints, std::vector<double>& doubles);
void mpi_unpack(int *ints, int *ii, double *doubles, int *dd);
#endif
private:
void zero();
void add(const cxxExchange &ex, double intensive, double extensive);
// not written
void dump_xml(std::ostream& os, unsigned int indent = 0)const;
protected:
std::list<cxxExchComp> exchComps;

68
NameDouble.cxx
View File

@ -275,63 +275,21 @@ CParser::STATUS_TYPE cxxNameDouble::read_raw(CParser& parser, std::istream::pos_
return CParser::PARSER_OK;
}
void cxxNameDouble::add_extensive(const cxxNameDouble &addee, double factor)
//
// Sums two name doubles, this + factor*nd2
//
void cxxNameDouble::add(const cxxNameDouble &old, double factor)
//
// constructor for cxxNameDouble from list of elt_list
//
{
if (factor == 0) return;
assert (factor > 0.0);
for (cxxNameDouble::const_iterator it = addee.begin(); it != addee.end(); it++)
{
cxxNameDouble::iterator current = (*this).find(it->first);
if (current != (*this).end())
{
(*this)[it->first] = current->second + it->second * factor;
} else {
(*this)[it->first] = it->second * factor;
}
}
}
void cxxNameDouble::add_intensive(const cxxNameDouble &addee, double f1, double f2)
//
// Sums two name doubles, this*f1 + f2*nd2
//
{
assert(f1 > 0 && f2 > 0);
for (cxxNameDouble::const_iterator it = addee.begin(); it != addee.end(); it++)
{
cxxNameDouble::iterator current = (*this).find(it->first);
if (current != (*this).end())
{
(*this)[it->first] = f1*current->second + f2*it->second;
} else {
(*this)[it->first] = f2 * it->second;
}
}
}
void cxxNameDouble::add_log_activities(const cxxNameDouble &addee, double f1, double f2)
//
// Sums two name doubles, this*f1 + f2*nd2, assuming log values
//
{
assert (f1 > 0 && f2 > 0);
for (cxxNameDouble::const_iterator it = addee.begin(); it != addee.end(); it++)
{
cxxNameDouble::iterator current = (*this).find(it->first);
if (current != (*this).end())
{
double a1 = pow(10., current->second);
double a2 = pow(10., it->second);
(*this)[it->first] = log10(f1*a1 + f2*a2);
} else {
//double a2 = pow(10. it->second);
//(*this)[it->first] = log10(f2 * a2);
(*this)[it->first] = it->second + log10(f2);
}
}
}
for (cxxNameDouble::const_iterator it = old.begin(); it != old.end(); it++) {
cxxNameDouble::iterator current = (*this).find(it->first);
if (current != (*this).end()) {
(*this)[it->first] = current->second + it->second * factor;
} else {
(*this)[it->first] = it->second * factor;
}
}
}
void cxxNameDouble::add(char * key, double total)
//
// add to total for a specified element

4
NameDouble.h
View File

@ -44,9 +44,7 @@ public:
CParser::STATUS_TYPE read_raw(CParser& parser, std::istream::pos_type& pos);
void add_extensive(const cxxNameDouble &old, double factor);
void add_intensive(const cxxNameDouble &addee, double fthis, double f2);
void add_log_activities(const cxxNameDouble &addee, double fthis, double f2);
void add(const cxxNameDouble &old, double factor);
void add(char * key, double total);
void insert(char *str, double d) {

172
Solution.cxx
View File

@ -27,21 +27,20 @@ cxxSolution::cxxSolution()
//
: cxxNumKeyword()
{
this->tc = 25.0;
this->ph = 7.0;
this->pe = 4.0;
this->mu = 1e-7;
this->ah2o = 1.0;
this->total_h = 111.1;
this->total_o = 55.55;
this->cb = 0.0;
this->mass_water = 1.0;
this->total_alkalinity = 0.0;
this->totals.type = cxxNameDouble::ND_ELT_MOLES;
this->master_activity.type = cxxNameDouble::ND_SPECIES_LA;
this->species_gamma.type = cxxNameDouble::ND_SPECIES_GAMMA;
tc = 25.0;
ph = 7.0;
pe = 4.0;
mu = 1e-7;
ah2o = 1.0;
total_h = 111.1;
total_o = 55.55;
cb = 0.0;
mass_water = 1.0;
total_alkalinity = 0.0;
totals.type = cxxNameDouble::ND_ELT_MOLES;
master_activity.type = cxxNameDouble::ND_SPECIES_LA;
species_gamma.type = cxxNameDouble::ND_SPECIES_GAMMA;
}
/*
cxxSolution::cxxSolution(double zero)
//
// empty cxxSolution constructor
@ -63,7 +62,6 @@ cxxSolution::cxxSolution(double zero)
master_activity.type = cxxNameDouble::ND_SPECIES_LA;
species_gamma.type = cxxNameDouble::ND_SPECIES_GAMMA;
}
*/
cxxSolution::cxxSolution(struct solution *solution_ptr)
//
// constructor for cxxSolution from struct solution
@ -77,18 +75,18 @@ isotopes(solution_ptr)
{
this->set_description(solution_ptr->description);
this->n_user = solution_ptr->n_user;
this->n_user_end = solution_ptr->n_user_end;
this->tc = solution_ptr->tc;
this->ph = solution_ptr->ph;
this->pe = solution_ptr->solution_pe;
this->mu = solution_ptr->mu;
this->ah2o = solution_ptr->ah2o;
this->total_h = solution_ptr->total_h;
this->total_o = solution_ptr->total_o;
this->cb = solution_ptr->cb;
this->mass_water = solution_ptr->mass_water;
this->total_alkalinity = solution_ptr->total_alkalinity;
n_user = solution_ptr->n_user;
n_user_end = solution_ptr->n_user_end;
tc = solution_ptr->tc;
ph = solution_ptr->ph;
pe = solution_ptr->solution_pe;
mu = solution_ptr->mu;
ah2o = solution_ptr->ah2o;
total_h = solution_ptr->total_h;
total_o = solution_ptr->total_o;
cb = solution_ptr->cb;
mass_water = solution_ptr->mass_water;
total_alkalinity = solution_ptr->total_alkalinity;
// Totals filled in constructor, just save description and moles
@ -103,83 +101,6 @@ isotopes(solution_ptr)
// Species_gamma in constructor
}
#ifdef SKIP
cxxSolution::cxxSolution( const std::map<int, cxxSolution>& solutions, cxxMix &mix)
//
// constructor for cxxSolution from struct solution
//cxxSolution *cxxStorageBin::mix_cxxSolutions(cxxMix &mixmap)
{
/*
* mixes solutions based on cxxMix structure, returns new solution
* return solution must be freed by calling method
*/
double intensive, extensive;
/*
* Zero out global solution data
*/
this->zero();
/*
* Determine sum of mixing fractions
*/
extensive = 0.0;
std::map<int, double> *mixcomps = mix.comps();
std::map<int, double>::const_iterator it;
for (it = mixcomps->begin(); it != mixcomps->end(); it++) {
extensive += it->second;
}
/*
* Add solutions
*/
for (it = mixcomps->begin(); it != mixcomps->end(); it++)
{
const cxxSolution *cxxsoln_ptr1 = &(solutions.find(it->first)->second);
if (cxxsoln_ptr1 == NULL)
{
sprintf(error_string, "Solution %d not found in mix_cxxSolutions.", it->first);
error_msg(error_string, CONTINUE);
input_error++;
return;
}
intensive = it->second/extensive;
double ext = it->second;
this->add(*cxxsoln_ptr1, intensive, it->second);
}
}
#endif
cxxSolution::cxxSolution( const std::map<int, cxxSolution>& solutions, cxxMix &mix)
//
// constructor for cxxSolution from mixture of solutions
//
{
//
// Zero out solution data
//
this->zero();
//
// Mix solutions
//
std::map<int, double> *mixcomps = mix.comps();
std::map<int, double>::const_iterator it;
for (it = mixcomps->begin(); it != mixcomps->end(); it++)
{
const cxxSolution *cxxsoln_ptr1 = &(solutions.find(it->first)->second);
if (cxxsoln_ptr1 == NULL)
{
sprintf(error_string, "Solution %d not found in mix_cxxSolutions.", it->first);
error_msg(error_string, CONTINUE);
input_error++;
return;
}
this->add(*cxxsoln_ptr1, it->second);
}
}
#ifdef SKIP
cxxSolution::cxxSolution(const cxxSolution &old, double intensive, double extensive)
//
@ -748,23 +669,6 @@ void cxxSolution::read_raw(CParser& parser)
return;
}
void cxxSolution::zero()
{
this->tc = 0.0;
this->ph = 0.0;
this->pe = 0.0;
this->mu = 0.0;
this->ah2o = 0.0;
this->total_h = 0.0;
this->total_o = 0.0;
this->cb = 0.0;
this->mass_water = 0.0;
this->total_alkalinity = 0.0;
this->totals.type = cxxNameDouble::ND_ELT_MOLES;
this->master_activity.type = cxxNameDouble::ND_SPECIES_LA;
this->species_gamma.type = cxxNameDouble::ND_SPECIES_GAMMA;
}
#ifdef SKIP
void cxxSolution::add(const cxxSolution &old, double intensive, double extensive)
//
// Add existing solution to "this" solution
@ -818,31 +722,7 @@ void cxxSolution::add(const cxxSolution &old, double intensive, double extensive
cxxNameDouble species_gamma;
*/
}
#endif
void cxxSolution::add(const cxxSolution &addee, double extensive)
//
// Add existing solution to "this" solution
//
{
double ext1 = this->mass_water;
double ext2 = addee.mass_water * extensive;
double f1 = ext1/(ext1 + ext2);
double f2 = ext2/(ext1 + ext2);
this->tc = f1*this->tc + f2*addee.tc;
this->ph = f1*this->ph + f2*addee.ph;
this->pe = f1*this->pe + f2*addee.pe;
this->mu = f1*this->mu + f2*addee.mu;
this->ah2o = f1*this->mu + f2*addee.ah2o;
this->total_h += addee.total_h * extensive;
this->total_o += addee.total_o * extensive;
this->cb += addee.cb * extensive;
this->mass_water += addee.mass_water * extensive;
this->total_alkalinity += addee.total_alkalinity * extensive;
this->totals.add_extensive(addee.totals, extensive);
this->master_activity.add_log_activities(addee.master_activity, f1, f2);
this->species_gamma.add_intensive(addee.species_gamma, f1, f2);
this->isotopes.add(addee.isotopes, f2, extensive);
}
double cxxSolution::get_total(char *string)const
{
cxxNameDouble::const_iterator it = this->totals.find(string);

22
Solution.h
View File

@ -11,7 +11,7 @@
#include <cassert> // assert
#include <map> // std::map
#include <string> // std::string
#include <list> // std::list
#include <vector> // std::vector
#include "char_star.h"
@ -21,11 +21,12 @@ class cxxSolution : public cxxNumKeyword
public:
cxxSolution();
//cxxSolution(double zero);
cxxSolution(double zero);
cxxSolution(struct solution *);
cxxSolution(const std::map<int, cxxSolution> &solution_map, cxxMix &mx);
//cxxSolution(const cxxSolution &old, double intensive, double extensive);
//cxxSolution(const cxxSolution&);
cxxSolution(const cxxSolution &old, double intensive, double extensive);
//cxxSolution(const cxxSolution&);
~cxxSolution();
//static cxxSolution& read(CParser& parser);
@ -76,10 +77,13 @@ public:
struct solution *cxxSolution2solution();
void dump_raw(std::ostream& s_oss, unsigned int indent)const;
void dump_xml(std::ostream& os, unsigned int indent = 0)const;
void dump_raw(std::ostream& s_oss, unsigned int indent)const;
void read_raw(CParser& parser);
void add(const cxxSolution &sol, double intensive, double extensive);
#ifdef USE_MPI
void mpi_pack(std::vector<int>& ints, std::vector<double>& doubles);
@ -87,14 +91,6 @@ public:
void mpi_send(int task_number);
void mpi_recv(int task_number);
#endif
private:
void zero();
//void add(const cxxSolution &sol, double intensive, double extensive);
void add(const cxxSolution &addee, double extensive);
// not checked
void dump_xml(std::ostream& os, unsigned int indent = 0)const;
protected:
double tc;
double ph;

93
SolutionIsotopeList.cxx
View File

@ -32,55 +32,54 @@ cxxSolutionIsotopeList::cxxSolutionIsotopeList(struct solution *solution_ptr)
}
void cxxSolutionIsotopeList::add(cxxSolutionIsotopeList old, double intensive, double extensive)
{
for (cxxSolutionIsotopeList::const_iterator itold = old.begin(); itold != old.end(); ++itold) {
//for (std::list <cxxSolutionIsotope>::const_iterator itold = old.isotopes.begin(); itold != old.isotopes.end(); ++itold) {
bool found = false;
for (cxxSolutionIsotopeList::iterator it = this->begin(); it != this->end(); ++it)
{
//for (std::list <cxxSolutionIsotope>::iterator it = this->isotopes.begin(); it != this->isotopes.end(); ++it) {
if ((it->isotope_number == itold->isotope_number) &&
(it->elt_name == itold->elt_name) &&
(it->isotope_name == itold->isotope_name))
{
it->total += itold->total * extensive;
it->ratio += itold->ratio * intensive;
it->ratio_uncertainty += itold->ratio_uncertainty * intensive;
it->ratio_uncertainty_defined = (it->ratio_uncertainty_defined || itold->ratio_uncertainty_defined);
found = true;
break;
}
}
if (!found)
{
cxxSolutionIsotope iso;
iso.total = itold->total * extensive;
iso.ratio = itold->ratio * intensive;
iso.ratio_uncertainty = itold->ratio_uncertainty * intensive;
iso.ratio_uncertainty_defined = itold->ratio_uncertainty_defined;
this->push_back(iso);
}
}
for (cxxSolutionIsotopeList::const_iterator itold = old.begin(); itold != old.end(); ++itold) {
//for (std::list <cxxSolutionIsotope>::const_iterator itold = old.isotopes.begin(); itold != old.isotopes.end(); ++itold) {
bool found = false;
for (cxxSolutionIsotopeList::iterator it = this->begin(); it != this->end(); ++it) {
//for (std::list <cxxSolutionIsotope>::iterator it = this->isotopes.begin(); it != this->isotopes.end(); ++it) {
if ((it->isotope_number == itold->isotope_number) &&
(it->elt_name == itold->elt_name) &&
(it->isotope_name == itold->isotope_name)) {
it->total += itold->total * extensive;
it->ratio += itold->ratio * intensive;
it->ratio_uncertainty += itold->ratio_uncertainty * intensive;
it->ratio_uncertainty_defined = (it->ratio_uncertainty_defined || itold->ratio_uncertainty_defined);
found = true;
break;
}
}
if (!found) {
cxxSolutionIsotope iso;
iso.total = itold->total * extensive;
iso.ratio = itold->ratio * intensive;
iso.ratio_uncertainty = itold->ratio_uncertainty * intensive;
iso.ratio_uncertainty_defined = itold->ratio_uncertainty_defined;
this->push_back(iso);
}
}
}
struct isotope * cxxSolutionIsotopeList::cxxSolutionIsotopeList2isotope()
{
struct isotope *iso;
if (this->size() <= 0) {
return NULL;
} else {
iso = (struct isotope *) PHRQ_malloc((size_t) ((this->size()) * sizeof(struct isotope)));
if (iso == NULL) malloc_error();
int i = 0;
for (cxxSolutionIsotopeList::iterator it = this->begin(); it != this->end(); ++it) {
iso[i].isotope_number = it->isotope_number;
iso[i].elt_name = it->elt_name;
iso[i].total = it->total;
iso[i].ratio = it->ratio;
iso[i].ratio_uncertainty = it->ratio_uncertainty;
iso[i].master = it->master();
iso[i].primary = it->primary();
i++;
}
}
return(iso);
struct isotope *iso;
if (this->size() <= 0) {
return NULL;
} else {
iso = (struct isotope *) PHRQ_malloc((size_t) ((this->size()) * sizeof(struct isotope)));
if (iso == NULL) malloc_error();
int i = 0;
for (cxxSolutionIsotopeList::iterator it = this->begin(); it != this->end(); ++it) {
iso[i].isotope_number = it->isotope_number;
iso[i].elt_name = it->elt_name;
iso[i].total = it->total;
iso[i].ratio = it->ratio;
iso[i].ratio_uncertainty = it->ratio_uncertainty;
iso[i].master = it->master();
iso[i].primary = it->primary();
i++;
}
}
return(iso);
}

4
StorageBin.cxx
View File

@ -495,7 +495,7 @@ void cxxStorageBin::remove(int n)
// Surface
this->Surfaces.erase(n);
}
#ifdef SKIP
cxxSolution *cxxStorageBin::mix_cxxSolutions(cxxMix &mixmap)
{
@ -536,7 +536,7 @@ cxxSolution *cxxStorageBin::mix_cxxSolutions(cxxMix &mixmap)
}
return(cxxsoln_ptr);
}
#endif
struct system *cxxStorageBin::cxxStorageBin2system(int n)
//
// make a system from storagebin

14
StorageBin.h
View File

@ -42,8 +42,8 @@ public:
}
return (NULL);
}
void setSolution(int n_user, cxxSolution &entity) {
Solutions[n_user] = entity;
void setSolution(int n_user, cxxSolution *entity) {
Solutions[n_user] = *entity;
}
void removeSolution(int n_user) {
Solutions.erase(n_user);
@ -137,17 +137,9 @@ public:
struct system *cxxStorageBin2system(int i);
//cxxSolution *mix_cxxSolutions(cxxMix &mixmap);
cxxSolution *mix_cxxSolutions(cxxMix &mixmap);
cxxExchange *mix_cxxExchange(cxxMix &mixmap);
const std::map<int, cxxSolution>& getSolutions()const {return this->Solutions;};
const std::map<int, cxxExchange>& getExchangers()const {return this->Exchangers;};
const std::map<int, cxxGasPhase>& getGasPhases()const {return this->GasPhases;};
const std::map<int, cxxKinetics>& getKinetics()const {return this->Kinetics;};
const std::map<int, cxxPPassemblage>& getPPassemblages()const {return this->PPassemblages;};
const std::map<int, cxxSSassemblage>& getSSassemblages()const {return this->SSassemblages;};
const std::map<int, cxxSurface>& getSurfaces()const {return this->Surfaces;};
#ifdef USE_MPI
void mpi_send(int n, int task_number);
void mpi_recv(int task_number);