diff --git a/KineticsComp.cxx b/KineticsComp.cxx index f777a0c9..0d015b8b 100644 --- a/KineticsComp.cxx +++ b/KineticsComp.cxx @@ -59,7 +59,7 @@ cxxKineticsComp::~cxxKineticsComp() } struct kinetics_comp * - cxxKineticsComp::cxxKineticsComp2kinetics_comp(PHREEQC_PTR_ARG_COMMA std::map < std::string, cxxKineticsComp > + cxxKineticsComp::cxxKineticsComp2kinetics_comp(PHREEQC_PTR_ARG_COMMA std::list < cxxKineticsComp > &el) // // Builds kinetics_comp structure from of cxxKineticsComp @@ -72,20 +72,20 @@ struct kinetics_comp * P_INSTANCE_POINTER malloc_error(); int i = 0; - for (std::map < std::string, cxxKineticsComp >::iterator it = el.begin(); + for (std::list < cxxKineticsComp >::iterator it = el.begin(); it != el.end(); ++it) { - if ((*it).second.rate_name.size() == 0) + if ((*it).rate_name.size() == 0) kinetics_comp_ptr[i].rate_name = NULL; else - kinetics_comp_ptr[i].rate_name = P_INSTANCE_POINTER string_hsave((*it).second.rate_name.c_str()); - kinetics_comp_ptr[i].list = (*it).second.namecoef.name_coef(P_INSTANCE); - kinetics_comp_ptr[i].count_list = (int) (*it).second.namecoef.size(); - kinetics_comp_ptr[i].tol = (*it).second.tol; - kinetics_comp_ptr[i].m = (*it).second.m; + kinetics_comp_ptr[i].rate_name = P_INSTANCE_POINTER string_hsave((*it).rate_name.c_str()); + kinetics_comp_ptr[i].list = (*it).namecoef.name_coef(P_INSTANCE); + kinetics_comp_ptr[i].count_list = (int) (*it).namecoef.size(); + kinetics_comp_ptr[i].tol = (*it).tol; + kinetics_comp_ptr[i].m = (*it).m; kinetics_comp_ptr[i].initial_moles = 0.; - kinetics_comp_ptr[i].m0 = (*it).second.m0; - kinetics_comp_ptr[i].moles = (*it).second.moles; + kinetics_comp_ptr[i].m0 = (*it).m0; + kinetics_comp_ptr[i].moles = (*it).moles; kinetics_comp_ptr[i].count_c_params = 0; kinetics_comp_ptr[i].c_params = NULL; /* @@ -93,15 +93,15 @@ struct kinetics_comp * kinetics_comp_ptr[i].d_params = NULL; */ - kinetics_comp_ptr[i].count_d_params = (int) (*it).second.d_params.size(); + kinetics_comp_ptr[i].count_d_params = (int) (*it).d_params.size(); kinetics_comp_ptr[i].d_params = NULL; - if ((*it).second.d_params.size() > 0) + if ((*it).d_params.size() > 0) { kinetics_comp_ptr[i].d_params = (double *) - P_INSTANCE_POINTER PHRQ_malloc((size_t) ((*it).second.d_params.size() * sizeof(double))); + P_INSTANCE_POINTER PHRQ_malloc((size_t) ((*it).d_params.size() * sizeof(double))); if (kinetics_comp_ptr[i].d_params == NULL) P_INSTANCE_POINTER malloc_error(); - std::copy((*it).second.d_params.begin(), (*it).second.d_params.end(), + std::copy((*it).d_params.begin(), (*it).d_params.end(), kinetics_comp_ptr[i].d_params); } i++; diff --git a/KineticsComp.h b/KineticsComp.h index 754034d7..acecc8c7 100644 --- a/KineticsComp.h +++ b/KineticsComp.h @@ -17,7 +17,7 @@ public: cxxKineticsComp(struct kinetics_comp *); ~cxxKineticsComp(); - static struct kinetics_comp *cxxKineticsComp2kinetics_comp(PHREEQC_PTR_ARG_COMMA std::map < std::string, cxxKineticsComp > &el); + static struct kinetics_comp *cxxKineticsComp2kinetics_comp(PHREEQC_PTR_ARG_COMMA std::list < cxxKineticsComp > &el); void dump_xml(std::ostream & os, unsigned int indent = 0) const; diff --git a/NameDouble.cxx b/NameDouble.cxx index c36fead0..2c16de98 100644 --- a/NameDouble.cxx +++ b/NameDouble.cxx @@ -114,6 +114,7 @@ cxxNameDouble::cxxNameDouble(struct master_activity *ma, int count, } this->type = type; } +/* cxxNameDouble::cxxNameDouble(struct name_coef *nc, int count) // // constructor for cxxNameDouble from list of elt_list @@ -128,6 +129,30 @@ cxxNameDouble::cxxNameDouble(struct name_coef *nc, int count) } this->type = ND_NAME_COEF; } +*/ +cxxNameDouble::cxxNameDouble(struct name_coef *nc, int count) + // + // constructor for cxxNameDouble from list of elt_list + // +{ + int i; + for (i = 0; i < count; i++) + { + if (nc[i].name == NULL) + continue; + + if ((*this).find(nc[i].name) == (*this).end()) + { + (*this)[nc[i].name] = nc[i].coef; + } + else + { + (*this)[nc[i].name] = (*this).find(nc[i].name)->second + nc[i].coef; + } + + } + this->type = ND_NAME_COEF; +} cxxNameDouble::~cxxNameDouble() { @@ -462,7 +487,70 @@ cxxNameDouble::multiply(double extensive) it->second *= extensive; } } +void +cxxNameDouble::merge_redox(const cxxNameDouble & source) +// +// Merges source into this +// Accounts for possible conflicts between redox state and +// totals +// +{ + for (cxxNameDouble::const_iterator sit = source.begin(); sit != source.end(); sit++) + { + + std::string redox_name = sit->first; + std::string elt_name; + size_t pos = redox_name.find("("); + bool redox; + if (pos != std::string::npos) + { + redox = true; + elt_name = redox_name.substr(0, pos - 1); + } + else + { + redox = false; + elt_name = redox_name; + } + if (redox) + { + // Remove elt_name, if present + if ((*this).find(elt_name) != (*this).end()) + { + (*this).erase((*this).find(elt_name)); + } + // Put in redox name + (*this)[redox_name] = sit->second; + } + else + { + std::string substring; + substring.append(elt_name); + substring.append("("); + + // Remove all redox + bool deleted = true; + while (deleted) + { + deleted = false; + cxxNameDouble::iterator current = (*this).begin(); + for ( ; current != (*this).end(); current++) + { + if (current->first.find(substring) == 0) + { + (*this).erase(current); + deleted = true; + break; + } + } + } + // Put in elt name + (*this)[elt_name] = sit->second; + } + + } +} #ifdef USE_MPI void cxxNameDouble::mpi_pack(std::vector < int >&ints, diff --git a/NameDouble.h b/NameDouble.h index 09e0f541..e94c533e 100644 --- a/NameDouble.h +++ b/NameDouble.h @@ -66,6 +66,9 @@ class cxxNameDouble:public void multiply(double factor); + void + merge_redox(const cxxNameDouble & source); + void insert(char *str, double d) { diff --git a/Solution.cxx b/Solution.cxx index 0eecc1c3..f5a15ba5 100644 --- a/Solution.cxx +++ b/Solution.cxx @@ -630,13 +630,27 @@ cxxSolution::read_raw(PHREEQC_PTR_ARG_COMMA CParser & parser, bool check) continue; case 0: // totals - if (this->totals.read_raw(P_INSTANCE_COMMA parser, next_char) != - CParser::PARSER_OK) + //if (this->totals.read_raw(P_INSTANCE_COMMA parser, next_char) != + // CParser::PARSER_OK) + //{ + // parser.incr_input_error(); + // parser. + // error_msg("Expected element name and moles for totals.", + // CParser::OT_CONTINUE); + //} { - parser.incr_input_error(); - parser. - error_msg("Expected element name and moles for totals.", - CParser::OT_CONTINUE); + cxxNameDouble temp_totals; + if (temp_totals.read_raw(P_INSTANCE_COMMA parser, next_char) != CParser::PARSER_OK) + { + parser.incr_input_error(); + parser. + error_msg("Expected element name and moles for totals.", + CParser::OT_CONTINUE); + } + else + { + this->totals.merge_redox(temp_totals); + } } opt_save = 0; break; diff --git a/cxxKinetics.cxx b/cxxKinetics.cxx index 02d50055..9856cd8b 100644 --- a/cxxKinetics.cxx +++ b/cxxKinetics.cxx @@ -66,7 +66,8 @@ totals(kinetics_ptr->totals) { cxxKineticsComp ec(&(kinetics_ptr->comps[i])); std::string str(ec.get_rate_name()); - this->kineticsComps[str] = ec; + //this->kineticsComps[str] = ec; + this->kineticsComps.push_back(ec); } // steps @@ -231,12 +232,12 @@ cxxKinetics::dump_raw(std::ostream & s_oss, unsigned int indent) const s_oss << "-cvode_order " << this->cvode_order << std::endl; // kineticsComps structures - for (std::map < std::string, cxxKineticsComp >::const_iterator it = + for (std::list < cxxKineticsComp >::const_iterator it = kineticsComps.begin(); it != kineticsComps.end(); ++it) { s_oss << indent1; s_oss << "-component" << std::endl; - (*it).second.dump_raw(s_oss, indent + 2); + (*it).dump_raw(s_oss, indent + 2); } // totals @@ -421,6 +422,28 @@ cxxKinetics::read_raw(PHREEQC_PTR_ARG_COMMA CParser & parser, bool check) reread.set_echo_file(CParser::EO_NONE); reread.set_echo_stream(CParser::EO_NONE); + std::list < cxxKineticsComp >::iterator kit; + bool found = false; + for (kit = this->kineticsComps.begin(); kit != this->kineticsComps.end(); kit++) + { + if (kit->get_rate_name() == ec.get_rate_name()) + { + found = true; + break; + } + } + if (found) + { + kit->read_raw(P_INSTANCE_COMMA reread, false); + } + else + { + cxxKineticsComp ec1; + ec1.read_raw(P_INSTANCE_COMMA reread, false); + std::string str(ec1.get_rate_name()); + this->kineticsComps.push_back(ec1); + } + /* if (this->kineticsComps.find(ec.get_rate_name()) != this->kineticsComps.end()) { cxxKineticsComp & comp = this->kineticsComps.find(ec.get_rate_name())->second; @@ -433,6 +456,7 @@ cxxKinetics::read_raw(PHREEQC_PTR_ARG_COMMA CParser & parser, bool check) std::string str(ec1.get_rate_name()); this->kineticsComps[str] = ec1; } + */ } useLastLine = true; break; @@ -554,11 +578,18 @@ cxxKinetics::mpi_pack(std::vector < int >&ints, { ints.push_back(this->n_user); ints.push_back((int) this->kineticsComps.size()); + std::list < cxxKineticsComp >::iterator it; + for (it = this->kineticsComps.begin(); it != this->kineticsComps.end(); it++) + { + (*it).mpi_pack(ints, doubles); + } + /* for (std::map < std::string, cxxKineticsComp >::iterator it = this->kineticsComps.begin(); it != this->kineticsComps.end(); it++) { (*it).second.mpi_pack(ints, doubles); } + */ ints.push_back((int) this->steps.size()); for (std::vector < double >::iterator it = this->steps.begin(); it != this->steps.end(); it++) @@ -588,8 +619,9 @@ cxxKinetics::mpi_unpack(int *ints, int *ii, double *doubles, int *dd) { cxxKineticsComp kc; kc.mpi_unpack(ints, &i, doubles, &d); - std::string str(kc.get_rate_name()); - this->kineticsComps[str] = kc; + //std::string str(kc.get_rate_name()); + //this->kineticsComps[str] = kc; + this->kineticsComps.push_back(kc); } n = ints[i++]; this->steps.clear(); @@ -664,21 +696,31 @@ cxxKinetics::add(const cxxKinetics & addee, double extensive) if (extensive == 0.0) return; //std::map < std::string, cxxKineticsComp> kineticsComps; - for (std::map < std::string, cxxKineticsComp >::const_iterator itadd = + for (std::list < cxxKineticsComp >::const_iterator itadd = addee.kineticsComps.begin(); itadd != addee.kineticsComps.end(); ++itadd) { - std::map < std::string, cxxKineticsComp >::iterator it = this->kineticsComps.find((*itadd).first); - if (it != this->kineticsComps.end()) + bool found(false); + std::list < cxxKineticsComp >::iterator it; + for (it = this->kineticsComps.begin(); it != this->kineticsComps.end(); it++) { - (*it).second.add((*itadd).second, extensive); + if ((*it).get_rate_name() == (*itadd).get_rate_name()) + { + found = true; + break; + } + } + + //std::map < std::string, cxxKineticsComp >::iterator it = this->kineticsComps.find((*itadd).first); + if (found) + { + (*it).add((*itadd), extensive); } else { - cxxKineticsComp entity = (*itadd).second; + cxxKineticsComp entity = (*itadd); entity.multiply(extensive); - std::string str(entity.get_rate_name()); - this->kineticsComps[str] = entity; + this->kineticsComps.push_back(entity); } } //std::vector steps; diff --git a/cxxKinetics.h b/cxxKinetics.h index 750b4d03..014b911d 100644 --- a/cxxKinetics.h +++ b/cxxKinetics.h @@ -43,7 +43,8 @@ class cxxKinetics:public cxxNumKeyword void add(const cxxKinetics & addee, double extensive); protected: - std::map < std::string, cxxKineticsComp > kineticsComps; + //std::map < std::string, cxxKineticsComp > kineticsComps; + std::list < cxxKineticsComp > kineticsComps; std::vector < double >steps; cxxNameDouble totals; double step_divide;