From 058375c2bfae4c48ea90712095936716e64643ab Mon Sep 17 00:00:00 2001 From: David Parkhurst Date: Wed, 31 Oct 2018 14:01:22 -0600 Subject: [PATCH 1/5] removed check of ss when sum of components is small --- model.cpp | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/model.cpp b/model.cpp index f89990a2..7a60f4e0 100644 --- a/model.cpp +++ b/model.cpp @@ -532,7 +532,7 @@ check_residuals(void) { if (x[i]->ss_in == FALSE) continue; - if (x[i]->moles <= MIN_TOTAL_SS) + if (x[i]->moles <= 1e2*MIN_TOTAL) continue; if (residual[i] >= epsilon || residual[i] <= -epsilon /* || stop_program == TRUE */ ) @@ -4207,6 +4207,7 @@ residuals(void) else if (x[i]->type == SS_MOLES) { residual[i] = x[i]->f * LOG_10; + if (x[i]->moles <= 1e2*MIN_TOTAL) continue; if (fabs(residual[i]) > l_toler && x[i]->ss_in == TRUE) { if (print_fail) From 906cfd49c7e93d59ec0372896197c5ffc444defb Mon Sep 17 00:00:00 2001 From: David Parkhurst Date: Sat, 3 Nov 2018 15:30:19 -0600 Subject: [PATCH 2/5] Check value of nmix --- transport.cpp | 13 ++++++++++++- 1 file changed, 12 insertions(+), 1 deletion(-) diff --git a/transport.cpp b/transport.cpp index d4eb8df1..16e3233e 100644 --- a/transport.cpp +++ b/transport.cpp @@ -7,6 +7,7 @@ #include "SSassemblage.h" #include "cxxKinetics.h" #include "Solution.h" +#include LDBLE F_Re3 = F_C_MOL / (R_KJ_DEG_MOL * 1e3); LDBLE tk_x2; // average tx_x of icell and jcell @@ -1107,8 +1108,15 @@ init_mix(void) } else { + const int max_int = (int)(std::numeric_limits::max); + if ((int)round(2.25 * maxmix) > max_int) + { + char token[MAX_LENGTH]; + sprintf(token, "Calculated number of mixes %g, is beyond program limit,\nERROR: please decrease time_step.", 2.25 * maxmix); + error_msg(token, STOP); + } if (bcon_first == 1 || bcon_last == 1) - l_nmix = 1 + (int) floor(2.25 * maxmix); + l_nmix = 1 + (long int) floorl(2.25 * maxmix); else l_nmix = 1 + (int) floor(1.5 * maxmix); @@ -2046,6 +2054,7 @@ multi_D(LDBLE DDt, int mobile_cell, int stagnant) * 1. obtain J_ij... */ il_calcs = find_J(icell, jcell, mixf, DDt, stagnant); + if (find_current) { if (i < last_c) @@ -2078,6 +2087,8 @@ multi_D(LDBLE DDt, int mobile_cell, int stagnant) continue; } } + if (!ct[icell].J_ij_count_spec) + continue; /* * 2. sum up the primary or secondary master_species From d74c8ff5544c8a18eae9d9f9ba4690ead369b49e Mon Sep 17 00:00:00 2001 From: David Parkhurst Date: Sat, 3 Nov 2018 22:33:06 -0600 Subject: [PATCH 3/5] Corrected syntax of integer limit --- model.cpp | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/model.cpp b/model.cpp index 7a60f4e0..84cc1db0 100644 --- a/model.cpp +++ b/model.cpp @@ -1597,7 +1597,8 @@ ineq(int in_kode) (size_t) (count_unknowns + 1) * sizeof(LDBLE)); ineq_array[l_count_rows * max_column_count + i] = 1.0; ineq_array[l_count_rows * max_column_count + count_unknowns] = - 0.99 * x[i]->moles - MIN_TOTAL_SS; + x[i]->moles; + //0.99 * x[i]->moles - MIN_TOTAL_SS; back_eq[l_count_rows] = i; l_count_rows++; } From b10df16b03374366d5fe67e1cfb69c15ce55275c Mon Sep 17 00:00:00 2001 From: David Parkhurst Date: Sat, 3 Nov 2018 22:35:15 -0600 Subject: [PATCH 4/5] Corrected syntax of integer limit, previous commit actually changed ss convergence parameter, used to multiply by 0.99 --- transport.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/transport.cpp b/transport.cpp index 16e3233e..0f573de4 100644 --- a/transport.cpp +++ b/transport.cpp @@ -1108,7 +1108,7 @@ init_mix(void) } else { - const int max_int = (int)(std::numeric_limits::max); + const int max_int = (int) (std::numeric_limits::max)(); if ((int)round(2.25 * maxmix) > max_int) { char token[MAX_LENGTH]; From c43c9af38d99e809ffdb44b1f3b6c8db06c8e934 Mon Sep 17 00:00:00 2001 From: David Parkhurst Date: Tue, 6 Nov 2018 10:17:35 -0700 Subject: [PATCH 5/5] tweaked ss, changed surf function per Kinniburgh --- basicsubs.cpp | 79 ++++++++++++++++++++++++++++++--------------------- model.cpp | 15 ++++++++++ 2 files changed, 62 insertions(+), 32 deletions(-) diff --git a/basicsubs.cpp b/basicsubs.cpp index 003726aa..c48d56b7 100644 --- a/basicsubs.cpp +++ b/basicsubs.cpp @@ -2445,46 +2445,61 @@ surf_total(const char *total_name, const char *surface_name) // surface matches, now match element or redox state struct rxn_token *rxn_ptr; - for (rxn_ptr = s_x[j]->rxn_s->token + 1; rxn_ptr->s != NULL; rxn_ptr++) + if (s_x[j]->mole_balance == NULL) { - if (redox && rxn_ptr->s->secondary) + for (rxn_ptr = s_x[j]->rxn_s->token + 1; rxn_ptr->s != NULL; rxn_ptr++) { - token = rxn_ptr->s->secondary->elt->name; - } - else if (!redox && rxn_ptr->s->secondary) - { - token = rxn_ptr->s->secondary->elt->primary->elt->name; - } - else if (!redox && rxn_ptr->s->primary) - { - token = rxn_ptr->s->primary->elt->name; - } - else - { - continue; - } - if (strcmp(token.c_str(), total_name) == 0) - { - t += rxn_ptr->coef * s_x[j]->moles; - break; - } - else - // sum all sites in case total_name is a surface name without underscore surf ("Hfo_w", "Hfo") - { - if (rxn_ptr->s->type == SURF) + if (redox && rxn_ptr->s->secondary) { - if (token.find("_") != std::string::npos) + token = rxn_ptr->s->secondary->elt->name; + } + else if (!redox && rxn_ptr->s->secondary) + { + token = rxn_ptr->s->secondary->elt->primary->elt->name; + } + else if (!redox && rxn_ptr->s->primary) + { + token = rxn_ptr->s->primary->elt->name; + } + else + { + continue; + } + if (strcmp(token.c_str(), total_name) == 0) + { + t += rxn_ptr->coef * s_x[j]->moles; + break; + } + else + // sum all sites in case total_name is a surface name without underscore surf ("Hfo_w", "Hfo") + { + if (rxn_ptr->s->type == SURF) { - token = token.substr(0, token.find("_")); - } - if (strcmp(token.c_str(), total_name) == 0) - { - t += rxn_ptr->coef * s_x[j]->moles; - break; + if (token.find("_") != std::string::npos) + { + token = token.substr(0, token.find("_")); + } + if (strcmp(token.c_str(), total_name) == 0) + { + t += rxn_ptr->coef * s_x[j]->moles; + break; + } } } } } + else + { + for (int i = 0; s_x[j]->next_secondary[i].elt != NULL; i++) + { + token = s_x[j]->next_secondary[i].elt->name; + if (strcmp(token.c_str(), total_name) == 0) + { + t += s_x[j]->next_secondary[i].coef * s_x[j]->moles; + break; + } + } + } } return t; } diff --git a/model.cpp b/model.cpp index 84cc1db0..3086158c 100644 --- a/model.cpp +++ b/model.cpp @@ -1599,6 +1599,7 @@ ineq(int in_kode) ineq_array[l_count_rows * max_column_count + count_unknowns] = x[i]->moles; //0.99 * x[i]->moles - MIN_TOTAL_SS; + //0.9 * x[i]->moles; back_eq[l_count_rows] = i; l_count_rows++; } @@ -3265,6 +3266,20 @@ reset(void) } factor = f0; } + if (x[i]->moles > 0 && -delta[i] > 1e1*x[i]->moles) + { + f0 = (-delta[i])/(1e1*x[i]->moles); + if (f0 > factor) + { + if (debug_model == TRUE) + { + output_msg(sformatf( + "%-10.10s, Precipitating too much mineral.\t%f\n", + x[i]->description, (double)f0)); + } + factor = f0; + } + } } } }