First fix of day.

git-svn-id: svn://136.177.114.72/svn_GW/phreeqc3/branches/concrete@10758 1feff8c3-07ed-0310-ac33-dd36852eb9cd
This commit is contained in:
David L Parkhurst 2016-01-28 21:24:32 +00:00
parent 7c56b9195e
commit 72991f1acc
2 changed files with 51 additions and 51 deletions

View File

@ -782,7 +782,7 @@ read_transport(void)
*/
if (count_por == 0)
{
if (old_cells < max_cells && multi_Dflag)
if (old_cells <= max_cells && multi_Dflag)
{
multi_Dpor = (multi_Dpor < 1e-10 ? 1e-10 : multi_Dpor);
if (multi_Dpor > 1e-10)
@ -792,23 +792,49 @@ read_transport(void)
error_string = sformatf(
"No porosities were read; used the minimal value %8.2e from -multi_D.", multi_Dpor);
warning_msg(error_string);
for (i = 1; i <= max_cells; i++)
for (i = 0; i < all_cells; i++)
cell_data[i].por = multi_Dpor;
}
}
else
{
for (i = 1; i <= count_por; i++)
cell_data[i].por = pors[i - 1];
if (max_cells > count_por)
if ((stag_data->exch_f > 0) && (stag_data->count_stag == 1))
{
error_string = sformatf(
"Porosities were read for %d cells. Last value is used till cell %d.",
count_por, max_cells);
"Mobile porosities were read, but mobile/immobile porosity was also defined in -stagnant. Using the values from -stagnant for mobile/immobile exchange and tortuosity factors.");
warning_msg(error_string);
for (i = count_por; i <= max_cells; i++)
cell_data[i + 1].por = pors[count_por - 1];
for (i = 1; i <= max_cells; i++)
cell_data[i].por = stag_data->th_m;
for (i++; i <= 2 * max_cells + 1; i++)
cell_data[i].por = stag_data->th_im;
}
else
{
for (i = 1; i <= count_por; i++)
cell_data[i].por = pors[i - 1];
if (max_cells > count_por)
{
error_string = sformatf(
"Porosities were read for %d cells. Last value is used till cell %d.",
count_por, max_cells);
warning_msg(error_string);
for (i = count_por; i <= max_cells; i++)
cell_data[i + 1].por = pors[count_por - 1];
}
}
}
if (interlayer_Dflag && !multi_Dflag)
{
input_error++;
error_string = sformatf(
"-multi_D must be defined, when -interlayer_D true.");
error_msg(error_string, CONTINUE);
}
for (i = 0; i < all_cells; i++)
{
interlayer_Dpor = (interlayer_Dpor < 1e-10 ? 1e-10 : interlayer_Dpor);
cell_data[i].por_il = interlayer_Dpor;
}
count_cells = max_cells;
/*
@ -870,24 +896,6 @@ read_transport(void)
else if (simul_tr == 1)
for (i = 1; i < all_cells; i++)
cell_data[i].print = TRUE;
/*
* Fill in porosities
*/
if (interlayer_Dflag && !multi_Dflag)
{
input_error++;
error_string = sformatf(
"-multi_D must be defined, when -interlayer_D true.");
error_msg(error_string, CONTINUE);
}
for (i = 0; i < all_cells; i++)
{
//multi_Dpor = (multi_Dpor < 1e-10 ? 1e-10 : multi_Dpor);
interlayer_Dpor = (interlayer_Dpor < 1e-10 ? 1e-10 : interlayer_Dpor);
//cell_data[i].por = multi_Dpor;
cell_data[i].por_il = interlayer_Dpor;
}
/*
* Calculate dump_modulus
*/
@ -976,8 +984,8 @@ read_transport(void)
* free storage for length, disp, punch
*/
length = (LDBLE *) free_check_null(length);
disp = (LDBLE *) free_check_null(disp);
pors = (LDBLE *) free_check_null(pors);
disp = (LDBLE *)free_check_null(disp);
pors = (LDBLE *)free_check_null(pors);
punch_temp = (int *) free_check_null(punch_temp);
print_temp = (int *) free_check_null(print_temp);
@ -1270,10 +1278,6 @@ dump_cpp(void)
fs << token;
sprintf(token, "\t-punch_cells");
fs << token;
//if (stag_data->count_stag > 0)
// j = 1 + (1 + stag_data->count_stag) * count_cells;
//else
// j = count_cells;
l = 0;
for (int i = 0; i < all_cells; i++)
{
@ -1292,10 +1296,6 @@ dump_cpp(void)
fs << token;
sprintf(token, "\t-print_cells");
fs << token;
//if (stag_data->count_stag > 0)
// j = 1 + (1 + stag_data->count_stag) * count_cells;
//else
// j = count_cells;
l = 0;
for (int i = 0; i < all_cells; i++)
{

View File

@ -337,8 +337,8 @@ transport(void)
/*
* Define stagnant/mobile mix structure, if not read explicitly.
*
* With count_stag = 1, mix factors are calculated from exchange factor <EFBFBD>
* (= exch_f), mobile <EFBFBD>_m (= th_m) and immobile <EFBFBD>_im (= th_im) porosity.
* With count_stag = 1, mix factors are calculated from exchange factor alpha
* (= exch_f), mobile th_m and immobile th_im porosity.
* These variables are read under keyword TRANSPORT, after stagnant, in
* structure stag_data.
* MIX 'cell_no' in input file can be an alternative for the calculation here.
@ -1588,7 +1588,7 @@ fill_spec(int l_cell_no)
struct master *master_ptr;
LDBLE dum, dum2;
LDBLE lm;
LDBLE por, por_il, temp_factor, temp_il_factor, viscos;
LDBLE por, por_il, viscos_f, viscos_il_f, viscos;
bool x_max_done = false;
s_ptr2 = NULL;
@ -1620,7 +1620,7 @@ fill_spec(int l_cell_no)
sol_D[l_cell_no].tk_x = tk_x;
temp_factor = temp_il_factor = 1.0;
viscos_f = viscos_il_f = 1.0;
if (l_cell_no == 0)
{
por = cell_data[1].por;
@ -1637,10 +1637,10 @@ fill_spec(int l_cell_no)
por_il = cell_data[l_cell_no].por_il;
}
if (por < multi_Dpor_lim)
por = temp_factor = 0.0;
por = viscos_f = 0.0;
if (por_il < interlayer_Dpor_lim)
por_il = temp_il_factor = 0.0;
por_il = viscos_il_f = 0.0;
/*
* correct diffusion coefficient for temperature and viscosity, D_T = D_298 * Tk * viscos_298 / (298 * viscos)
* modify viscosity effect: Dw(TK) = Dw(298.15) * exp(dw_t / TK - dw_t / 298.15), SC data from Robinson and Stokes, 1959
@ -1649,8 +1649,8 @@ fill_spec(int l_cell_no)
/*
* put temperature factor in por_factor which corrects for porous medium...
*/
temp_factor *= tk_x * viscos_0_25 / (298.15 * viscos);
temp_il_factor *= tk_x * viscos_0_25 / (298.15 * viscos);
viscos_f *= tk_x * viscos_0_25 / (298.15 * viscos);
viscos_il_f *= tk_x * viscos_0_25 / (298.15 * viscos);
sol_D[l_cell_no].viscos_f = tk_x * viscos_0_25 / (298.15 * viscos);
count_spec = count_exch_spec = 0;
@ -1735,17 +1735,17 @@ fill_spec(int l_cell_no)
sol_D[l_cell_no].spec[count_spec].z = s_ptr2->z;
if (s_ptr2->dw == 0)
sol_D[l_cell_no].spec[count_spec].Dwt =
default_Dw * temp_il_factor;
default_Dw * viscos_il_f;
else
{
if (s_ptr2->dw_t)
{
sol_D[l_cell_no].spec[count_spec].Dwt = s_ptr2->dw *
exp(s_ptr2->dw_t / 298.15 - s_ptr2->dw_t / tk_x) * temp_il_factor;
exp(s_ptr2->dw_t / 298.15 - s_ptr2->dw_t / tk_x) * viscos_il_f;
sol_D[l_cell_no].spec[count_spec].dw_t = s_ptr2->dw_t;
}
else
sol_D[l_cell_no].spec[count_spec].Dwt = s_ptr2->dw * temp_il_factor;
sol_D[l_cell_no].spec[count_spec].Dwt = s_ptr2->dw * viscos_il_f;
}
count_exch_spec++;
count_spec++;
@ -1765,17 +1765,17 @@ fill_spec(int l_cell_no)
sol_D[l_cell_no].spec[count_spec].lg = s_ptr->lg;
sol_D[l_cell_no].spec[count_spec].z = s_ptr->z;
if (s_ptr->dw == 0)
sol_D[l_cell_no].spec[count_spec].Dwt = default_Dw * temp_factor;
sol_D[l_cell_no].spec[count_spec].Dwt = default_Dw * viscos_f;
else
{
if (s_ptr->dw_t)
{
sol_D[l_cell_no].spec[count_spec].Dwt = s_ptr->dw *
exp(s_ptr->dw_t / tk_x - s_ptr->dw_t / 298.15) * temp_factor;
exp(s_ptr->dw_t / tk_x - s_ptr->dw_t / 298.15) * viscos_f;
sol_D[l_cell_no].spec[count_spec].dw_t = s_ptr->dw_t;
}
else
sol_D[l_cell_no].spec[count_spec].Dwt = s_ptr->dw * temp_factor;
sol_D[l_cell_no].spec[count_spec].Dwt = s_ptr->dw * viscos_f;
}
if (sol_D[l_cell_no].spec[count_spec].Dwt * pow(por, multi_Dn) > diffc_max)
diffc_max = sol_D[l_cell_no].spec[count_spec].Dwt * pow(por, multi_Dn);