diff --git a/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/PBasic.cpp b/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/PBasic.cpp index 1ccf4126..53dda855 100644 --- a/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/PBasic.cpp +++ b/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/PBasic.cpp @@ -4732,6 +4732,8 @@ factor(struct LOC_exec * LINK) case tokviscos: { + if (PhreeqcPtr->print_viscosity) + PhreeqcPtr->viscosity(nullptr); n.UU.val = (parse_all) ? 1 : PhreeqcPtr->viscos; } break; diff --git a/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/basicsubs.cpp b/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/basicsubs.cpp index c92e8c30..7724d062 100644 --- a/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/basicsubs.cpp +++ b/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/basicsubs.cpp @@ -205,8 +205,11 @@ diff_c(const char *species_name) Dw *= viscos_0_25 / viscos_0; } - if (s_ptr->dw_a_v_dif) + if (s_ptr->dw_a_v_dif && print_viscosity) + { + viscosity(nullptr); Dw *= pow(viscos_0 / viscos, s_ptr->dw_a_v_dif); + } return Dw; } @@ -235,8 +238,11 @@ setdiff_c(const char *species_name, double d, double d_v_d) Dw *= viscos_0_25 / viscos_0; } - if (d_v_d) - Dw *= pow(viscos_0 / viscos, d_v_d); + if (d_v_d && print_viscosity) + { + viscosity(nullptr); + Dw *= pow(viscos_0 / viscos, s_ptr->dw_a_v_dif); + } return Dw; } /* ---------------------------------------------------------------------- */ @@ -269,6 +275,8 @@ calc_SC(void) // } //} av = 0; + if (print_viscosity) + viscosity(nullptr); if (!Falk) { for (i = 0; i < (int)this->s_x.size(); i++) @@ -1171,6 +1179,29 @@ diff_layer_total(const char* total_name, const char* surface_name) return (0); } } + else if (strcmp_nocase("viscos_ddl", total_name) == 0) + { + if (dl_type_x != cxxSurface::NO_DL) + { + cxxSurfaceCharge* charge_ptr = use.Get_surface_ptr()->Find_charge(x[j]->surface_charge); + if (charge_ptr->Get_mass_water() > 0) + { + cxxSurface * surf_ptr = use.Get_surface_ptr(); + if (surf_ptr->Get_calc_viscosity()) + { + viscosity(surf_ptr); + viscosity(nullptr); + return charge_ptr->Get_DDL_viscosity(); + } + else + return charge_ptr->Get_DDL_viscosity() * viscos; + } + } + else + { + return (0); + } + } /* * find total moles of each element in diffuse layer... */ diff --git a/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/class_main.cpp b/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/class_main.cpp index 731f0227..79e6475f 100644 --- a/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/class_main.cpp +++ b/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/class_main.cpp @@ -1,6 +1,8 @@ +#ifndef NPP #ifdef DOS #include #endif +#endif #include "Phreeqc.h" #include "NameDouble.h" @@ -122,12 +124,18 @@ main_method(int argc, char *argv[]) { return errors; } +#ifdef NPP +#ifdef DOS + write_banner(); +#endif +#else #ifndef NO_UTF8_ENCODING #ifdef DOS SetConsoleOutputCP(CP_UTF8); #endif write_banner(); #endif +#endif /* * Initialize arrays @@ -214,11 +222,17 @@ main_method(int argc, char *argv[]) { return errors; } +#ifdef NPP +#ifdef DOS + write_banner(); +#endif +#else #ifndef NO_UTF8_ENCODING #ifdef DOS SetConsoleOutputCP(CP_UTF8); #endif write_banner(); +#endif #endif /* @@ -281,6 +295,7 @@ main_method(int argc, char *argv[]) return 0; } #endif //TEST_COPY +#ifdef NPP /* ---------------------------------------------------------------------- */ int Phreeqc:: write_banner(void) @@ -293,16 +308,69 @@ write_banner(void) char buffer[80]; int len, indent; screen_msg( - " █▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀█\n"); + " ÛßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßßÛ\n"); screen_msg( - " ║ ║\n"); + " º º\n"); /* version */ #ifdef NPP - len = sprintf(buffer, "* PHREEQC-%s *", "3.7.3"); + len = sprintf(buffer, "* PHREEQC-%s *", "3.8.0"); #else - len = snprintf(buffer, sizeof(buffer), "* PHREEQC-%s *", "@VERSION@"); + len = sprintf(buffer, "* PHREEQC-%s *", "@VERSION@"); +#endif + indent = (49 - len) / 2; + screen_msg(sformatf("%14cº%*c%s%*cº\n", ' ', indent, ' ', buffer, + 49 - indent - len, ' ')); + + screen_msg( + " º º\n"); + screen_msg( + " º A hydrogeochemical transport model º\n"); + screen_msg( + " º º\n"); + screen_msg( + " º by º\n"); + screen_msg( + " º D.L. Parkhurst and C.A.J. Appelo º\n"); + screen_msg( + " º º\n"); + + + /* date */ +#ifdef NPP + len = sprintf(buffer, "%s", "August 27, 2024, with bug-fixes and new items"); +#else + len = sprintf(buffer, "%s", "@VER_DATE@"); #endif + indent = (49 - len) / 2; + screen_msg(sformatf("%14cº%*c%s%*cº\n", ' ', indent, ' ', buffer, + 49 - indent - len, ' ')); + + screen_msg( + " ÛÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÜÛ\n\n"); + +#endif + return 0; +} +#else +/* ---------------------------------------------------------------------- */ +int Phreeqc:: +write_banner(void) +/* ---------------------------------------------------------------------- */ +{ +#ifdef TESTING + return OK; +#endif +#ifndef NO_UTF8_ENCODING + char buffer[80]; + int len, indent; + screen_msg( + " █▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀█\n"); + screen_msg( + " ║ ║\n"); + + /* version */ + len = snprintf(buffer, sizeof(buffer), "* PHREEQC-%s *", "@VERSION@"); indent = (44 - len) / 2; screen_msg(sformatf("%14c║%*c%s%*c║\n", ' ', indent, ' ', buffer, 44 - indent - len, ' ')); @@ -322,11 +390,7 @@ write_banner(void) /* date */ -#ifdef NPP - len = sprintf(buffer, "%s", "March 14, 2024, with bug-fixes and new items"); -#else len = snprintf(buffer, sizeof(buffer), "%s", "@VER_DATE@"); -#endif indent = (44 - len) / 2; screen_msg(sformatf("%14c║%*c%s%*c║\n", ' ', indent, ' ', buffer, 44 - indent - len, ' ')); @@ -336,6 +400,7 @@ write_banner(void) #endif return 0; } +#endif /* ---------------------------------------------------------------------- */ int Phreeqc:: @@ -507,7 +572,7 @@ process_file_names(int argc, char *argv[], std::istream **db_cookie, output_msg(sformatf(" Input file: %s\n", in_file.c_str())); output_msg(sformatf(" Output file: %s\n", out_file.c_str())); #ifdef NPP - output_msg(sformatf("Using PHREEQC: version 3.7.3, compiled March 14, 2024, with bug-fixes and new items\n")); + output_msg(sformatf("Using PHREEQC: version 3.8.2, compiled August 27, 2024, with bug-fixes and new items\n")); #endif output_msg(sformatf("Database file: %s\n\n", token.c_str())); #ifdef NPP @@ -516,7 +581,7 @@ process_file_names(int argc, char *argv[], std::istream **db_cookie, /* * local cleanup */ - + line = (char *) free_check_null(line); line_save = (char *) free_check_null(line_save); diff --git a/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/print.cpp b/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/print.cpp index 325e1243..033d4899 100644 --- a/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/print.cpp +++ b/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/print.cpp @@ -1,3 +1,4 @@ +//Note to encode in ANSI with NP++ #include "Utils.h" #include "Phreeqc.h" #include "phqalloc.h" @@ -270,14 +271,23 @@ print_diffuse_layer(cxxSurfaceCharge *charge_ptr) output_msg(sformatf( "\tWater in diffuse layer: %8.3e kg, %4.1f%% of total DDL-water.\n", (double) charge_ptr->Get_mass_water(), (double) d)); - if (charge_ptr->Get_DDL_viscosity()) + if (print_viscosity && d > 0) { - if (d == 100) - output_msg(sformatf( - "\t\t viscosity: %7.5f mPa s.\n", (double)charge_ptr->Get_DDL_viscosity())); + cxxSurface * surf_ptr = use.Get_surface_ptr(); + if (surf_ptr->Get_calc_viscosity()) + { + viscosity(surf_ptr); + viscosity(nullptr); + if (d == 100) + output_msg(sformatf( + "\t\t calculated viscosity: %7.5f mPa s.\n", (double)charge_ptr->Get_DDL_viscosity())); + else + output_msg(sformatf( + "\t\t calculated viscosity: %7.5f mPa s for this DDL water. (%7.5f mPa s for total DDL-water.)\n", (double)charge_ptr->Get_DDL_viscosity(), (double)use.Get_surface_ptr()->Get_DDL_viscosity())); + } else - output_msg(sformatf( - "\t\t viscosity: %7.5f mPa s for this DDL water. (%7.5f mPa s for total DDL-water.)\n", (double)charge_ptr->Get_DDL_viscosity(), (double)use.Get_surface_ptr()->Get_DDL_viscosity())); + output_msg(sformatf( + "\t\t viscosity: %7.5f mPa s for DDL water.\n", (double)charge_ptr->Get_DDL_viscosity() * viscos)); } if (use.Get_surface_ptr()->Get_debye_lengths() > 0 && d > 0) @@ -2256,8 +2266,8 @@ print_totals(void) //#ifdef NPP if (print_viscosity) { - output_msg(sformatf("%45s%9.5f", "Viscosity (mPa s) = ", - (double) viscos)); + viscosity(nullptr); + output_msg(sformatf("%45s%9.5f", "Viscosity (mPa s) = ", (double) viscos)); if (tc_x > 200 && !pure_water) { #ifdef NO_UTF8_ENCODING diff --git a/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/transport.cpp b/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/transport.cpp index 109fd662..c49ee26b 100644 --- a/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/transport.cpp +++ b/src/IPhreeqcPhast/IPhreeqc/phreeqcpp/transport.cpp @@ -962,6 +962,7 @@ transport(void) { snprintf(token, sizeof(token), "\nFor balancing negative concentrations in MCD, added in total to the system:"); + count_warnings = pr.warnings - 1; warning_msg(token); for (i = 0; i < count_moles_added; i++) { @@ -970,6 +971,7 @@ transport(void) snprintf(token, sizeof(token), "\t %.4e moles %s.", (double)moles_added[i].moles, moles_added[i].name); + count_warnings = pr.warnings - 1; warning_msg(token); } } @@ -3937,15 +3939,20 @@ find_J(int icell, int jcell, LDBLE mixf, LDBLE DDt, int stagnant) { if (s_ptr1->Get_dl_type() != cxxSurface::NO_DL) { + if (s_ptr1->Get_calc_viscosity()) + { + viscosity(s_ptr1); + ct[icell].visc1 = s_ptr1->Get_DDL_viscosity(); + } + else + ct[icell].visc1 = s_ptr1->Get_DDL_viscosity() * Utilities::Rxn_find(Rxn_solution_map, icell)->Get_viscosity(); + s_charge_p.assign(s_ptr1->Get_surface_charges().begin(), s_ptr1->Get_surface_charges().end()); s_com_p.assign(s_ptr1->Get_surface_comps().begin(), s_ptr1->Get_surface_comps().end()); if (s_ptr1->Get_only_counter_ions()) only_counter = TRUE; - ct[icell].visc1 = s_ptr1->Get_DDL_viscosity(); - if (s_ptr1->Get_calc_viscosity()) - ct[icell].visc1 /= Utilities::Rxn_find(Rxn_solution_map, icell)->Get_viscosity(); /* find the immobile surface charges with DL... */ for (i = 0; i < (int)s_charge_p.size(); i++) { @@ -3966,16 +3973,20 @@ find_J(int icell, int jcell, LDBLE mixf, LDBLE DDt, int stagnant) { if (s_ptr2->Get_dl_type() != cxxSurface::NO_DL) { + if (s_ptr2->Get_calc_viscosity()) + { + viscosity(s_ptr2); + ct[icell].visc2 = s_ptr2->Get_DDL_viscosity(); + } + else + ct[icell].visc2 = s_ptr2->Get_DDL_viscosity() * Utilities::Rxn_find(Rxn_solution_map, jcell)->Get_viscosity(); + s_charge_p.assign(s_ptr2->Get_surface_charges().begin(), s_ptr2->Get_surface_charges().end()); s_com_p.assign(s_ptr2->Get_surface_comps().begin(), s_ptr2->Get_surface_comps().end()); if (s_ptr2->Get_only_counter_ions()) only_counter = TRUE; - ct[icell].visc2 = s_ptr2->Get_DDL_viscosity(); - if (s_ptr2->Get_calc_viscosity()) - ct[icell].visc2 /= Utilities::Rxn_find(Rxn_solution_map, jcell)->Get_viscosity(); - for (i = 0; i < (int)s_charge_p.size(); i++) { for (i1 = 0; i1 < (int)s_com_p.size(); i1++) @@ -3990,6 +4001,7 @@ find_J(int icell, int jcell, LDBLE mixf, LDBLE DDt, int stagnant) } } } + viscosity(nullptr); if (!stagnant) { if (icell == 0) @@ -5975,8 +5987,11 @@ viscosity(cxxSurface *surf_ptr) /* ---------------------------------------------------------------------- */ { if (surf_ptr && !surf_ptr->Get_calc_viscosity()) + { + for (int i = 0; i < (int)surf_ptr->Get_surface_charges().size(); i++) + surf_ptr->Get_surface_charges()[i].Set_DDL_viscosity(surf_ptr->Get_DDL_viscosity()); return surf_ptr->Get_DDL_viscosity(); - + } /* from Atkins, 1994. Physical Chemistry, 5th ed. */ //viscos = @@ -6120,7 +6135,7 @@ viscosity(cxxSurface *surf_ptr) eq_plus += l_moles; l_mu_x += l_moles * s_x[i]->z; } - l_mu_x += fabs(eq_plus + eq_min); + l_mu_x += fabs(eq_plus + eq_min); // add surface charge l_mu_x /= (2 * l_water); eq_plus = eq_min = 0; } @@ -6140,6 +6155,8 @@ viscosity(cxxSurface *surf_ptr) { dw_t_visc = 0; t1 = l_moles / l_water; + if (t1 < 1e-9) + continue; l_z = fabs(s_x[i]->z); if (l_z) f_z = (l_z * l_z + l_z) / 2; @@ -6172,16 +6189,21 @@ viscosity(cxxSurface *surf_ptr) if (!surf_ptr) s_x[i]->dw_t_visc = dw_t_visc + t3; //output_msg(sformatf("\t%s\t%e\t%e\t%e\n", s_x[i]->name, t1, Bc, Dc )); } - // parms for A... + // parms for A and V_an. 7/26/24: added V_an calculation for gases z = 0 if ((l_z = s_x[i]->z) == 0) - continue; - Dw = s_x[i]->dw; - if (Dw) { - Dw *= (0.89 / viscos_0 * tk_x / 298.15); - if (s_x[i]->dw_t) - Dw *= exp(s_x[i]->dw_t / tk_x - s_x[i]->dw_t / 298.15); + if (s_x[i]->Jones_Dole[6]) + { + V_an += s_x[i]->logk[vm_tc] * s_x[i]->Jones_Dole[6] * l_moles; + m_an += l_moles; + } + continue; } + if ((Dw = s_x[i]->dw) == 0) + continue; + Dw *= (0.89 / viscos_0 * tk_x / 298.15); + if (s_x[i]->dw_t) + Dw *= exp(s_x[i]->dw_t / tk_x - s_x[i]->dw_t / 298.15); if (l_z < 0) { if (!strcmp(s_x[i]->name, "Cl-")) @@ -6191,21 +6213,18 @@ viscosity(cxxSurface *surf_ptr) V_an += V_Cl * l_moles; m_an += l_moles; } - else// if (s_x[i]->Jones_Dole[6]) + else if (s_x[i]->Jones_Dole[6]) { V_an += s_x[i]->logk[vm_tc] * s_x[i]->Jones_Dole[6] * l_moles; m_an += l_moles; } - if (Dw) - { - // anions for A... - m_min += l_moles; - t1 = l_moles * l_z; - eq_min -= t1; - eq_dw_min -= t1 / Dw; - } + // anions for A... + m_min += l_moles; + t1 = l_moles * l_z; + eq_min -= t1; + eq_dw_min -= t1 / Dw; } - else if (Dw) + else { // cations for A... m_plus += l_moles;