23,15 → 23,16 |
#include <stdlib.h> |
#include <math.h> |
#include "aliast.h" |
#include "constrel.h" |
#include "consol_hawf3.h" |
#include "globalt.h" |
#include "C60bazant.h" //testing |
#include "C60baroghel.h" //testing |
#include "baroghel_reten.h" |
#include "bazant_reten.h" |
#include "globmatt.h" |
|
con_hawf3mat::con_hawf3mat() |
{ |
compress = 0; //compressible grains: 0=no; 1=yes |
compress = 0; //compressible grains: 0=no; 1=yes |
kintr_type = 0; //intrinsic permability calculation type |
krw_type = 0; //relative permeability calculation type |
krg_type = 0; //relative permability calculation type |
137,6 → 138,10 |
|
pw_bc = 0.0; //free boundary pressure |
rel_gas_press = 0; |
|
//parameters for the artificial material: |
kww0 = kwg0 = kwt0 = kgw0 = kgg0 = kgt0 = ktw0 = ktg0 = ktt0 = 0.0; |
capww0 = capwg0 = capwt0 = capgw0 = capgg0 = capgt0 = captw0 = captg0 = captt0 =0.0; |
} |
con_hawf3mat::~con_hawf3mat() |
{} |
158,8 → 163,11 |
xfscanf (in,"%le %le %le %le %le %le %le %d %d %d %d %d %d %d %d", &alpha0, &ks0, &phi0, &rhos0, &pw_bc, &tau0, &kintr0, &kintr_type, &krw_type, &krg_type, &deff_type, &sr_type, &lambda_type, &cps_type, &betas_type); |
|
switch (model_type){ |
case artificial3:{//artificial isotropic material for non-isotherma air-water flow |
xfscanf (in,"%le %le %le %le %le %le %le %le %le", &kww0, &kwg0, &kwt0, &kgw0, &kgg0, &kgt0, &ktw0, &ktg0, &ktt0); |
xfscanf (in,"%le %le %le %le %le %le %le %le %le", &capww0, &capwg0, &capwt0, &capgw0, &capgg0, &capgt0, &captw0, &captg0, &captt0); |
} |
case lewis_and_schrefler3:{//Lewis and Schrefler's model approach |
lewis_ret.read(in); |
break; |
} |
case lewis_and_schrefler3_mefel:{//Lewis and Schrefler's model approach coupled with mechanics |
172,141 → 180,152 |
} |
} |
|
|
//intrinsic permebility calculation: |
switch (kintr_type){ |
case 0:{//constant |
break; |
if(model_type != artificial3){ |
//intrinsic permebility calculation: |
switch (kintr_type){ |
case 0:{//constant |
break; |
} |
case 1:{//dependent on porosity |
xfscanf (in,"%le %le", &bb1, &phi01); |
break; |
} |
default:{ |
print_err ("unknown model type is required", __FILE__, __LINE__, __func__); |
abort(); |
} |
} |
|
//relative permebility calculation: |
switch (krw_type){ |
case 0:{//constant |
break; |
} |
case 1: |
case 2:{//dependent on saturation degree: |
xfscanf (in,"%le %le", &sirr, &ssat); |
break; |
} |
case 3:{//exponential |
xfscanf (in,"%le", &lambda_krw); |
break; |
} |
case 9:{//bazant |
break; |
} |
default:{ |
print_err ("unknown model type is required", __FILE__, __LINE__, __func__); |
abort(); |
} |
} |
|
//gas relative permebility calculation: |
switch (krg_type){ |
case 0:{//constant |
break; |
} |
case 1: |
case 2:{//dependent on saturation degree: |
//xfscanf (in,"%le %le", &sirr, &ssat); //not finished |
break; |
} |
case 3:{//exponential |
//xfscanf (in,"%le", &lambda_krw); //not finished |
break; |
} |
default:{ |
print_err ("unknown model type is required", __FILE__, __LINE__, __func__); |
abort(); |
} |
} |
|
//retention curve type: |
switch (sr_type){ |
case baroghel_sr:{//Baroghle-Bouny approach |
baroghel_ret.read(in); |
break; |
} |
case bazant_sr:{//Bazant approach |
break; |
} |
case lewis_and_schrefler_sr:{//Lewis and Schrefler's book |
lewis_ret.read(in); |
break; |
} |
case van_genuchten_sr:{//partially saturated medium =Van Genuchten model |
van_genuchten_ret.read(in); |
break; |
} |
case van_genuchten2_sr:{//partially saturated medium =Van Genuchten model |
//van_genuchten_ret.read2(in); |
break; |
} |
case mefel_sr:{//from MEFEL |
break; |
} |
case table_sr:{//from table |
//reading of retention curve: |
data.read (in); |
break; |
} |
default:{ |
print_err ("unknown model type is required", __FILE__, __LINE__, __func__); |
abort(); |
} |
} |
|
//conductivity calculation: |
switch (lambda_type){ |
case 0:{//constant |
xfscanf (in,"%le", &lambda0); |
break; |
} |
case 1:{//dependent on moisture |
xfscanf (in,"%le %le", &lambda_dry, &lambda_wet); |
break; |
} |
case 2:{//dependent on moisture |
xfscanf (in,"%le %le %le %le", &lambda_dry, &lambda_wet, &sr_dry, &sr_wet); |
break; |
} |
default:{ |
print_err ("unknown model type is required", __FILE__, __LINE__, __func__); |
abort(); |
} |
} |
|
//specific heat calculation: |
switch (cps_type){ |
case 0:{//constant |
xfscanf (in,"%le", &cps0); |
break; |
} |
case 1:{//dependent on moisture |
xfscanf (in,"%le %le", &cps_dry, &cps_wet); |
break; |
} |
default:{ |
print_err ("unknown model type is required", __FILE__, __LINE__, __func__); |
abort(); |
} |
} |
|
//thermal expansion calculation: |
switch (betas_type){ |
case 0:{//constant |
xfscanf (in,"%le", &betas0); |
break; |
} |
case 1:{//dependent on moisture |
xfscanf (in,"%le %le", &betas_dry, &betas_wet); |
break; |
} |
default:{ |
print_err ("unknown model type is required", __FILE__, __LINE__, __func__); |
abort(); |
} |
} |
} |
case 1:{//dependent on porosity |
xfscanf (in,"%le %le", &bb1, &phi01); |
break; |
} |
default:{ |
print_err ("unknown model type is required", __FILE__, __LINE__, __func__); |
abort(); |
} |
} |
|
//relative permebility calculation: |
switch (krw_type){ |
case 0:{//constant |
break; |
} |
case 1: |
case 2:{//dependent on saturation degree: |
xfscanf (in,"%le %le", &sirr, &ssat); |
break; |
} |
case 3:{//exponential |
xfscanf (in,"%le", &lambda_krw); |
break; |
} |
default:{ |
print_err ("unknown model type is required", __FILE__, __LINE__, __func__); |
abort(); |
} |
} |
|
//gas relative permebility calculation: |
switch (krg_type){ |
case 0:{//constant |
break; |
} |
case 1: |
case 2:{//dependent on saturation degree: |
//xfscanf (in,"%le %le", &sirr, &ssat); //not finished |
break; |
} |
case 3:{//exponential |
//xfscanf (in,"%le", &lambda_krw); //not finished |
break; |
} |
default:{ |
print_err ("unknown model type is required", __FILE__, __LINE__, __func__); |
abort(); |
} |
} |
|
//retention curve type: |
switch (sr_type){ |
case lewis_and_schrefler_sr:{//Lewis and Schrefler's book |
lewis_ret.read(in); |
break; |
} |
case van_genuchten_sr:{//partially saturated medium =Van Genuchten model |
van_genuchten_ret.read(in); |
break; |
} |
case van_genuchten2_sr:{//partially saturated medium =Van Genuchten model |
//van_genuchten_ret.read2(in); |
break; |
} |
case mefel_sr:{//from MEFEL |
break; |
} |
case table_sr:{//from table |
//reading of retention curve: |
data.read (in); |
break; |
} |
default:{ |
print_err ("unknown model type is required", __FILE__, __LINE__, __func__); |
abort(); |
} |
} |
|
//conductivity calculation: |
switch (lambda_type){ |
case 0:{//constant |
xfscanf (in,"%le", &lambda0); |
break; |
} |
case 1:{//dependent on moisture |
xfscanf (in,"%le %le", &lambda_dry, &lambda_wet); |
break; |
} |
case 2:{//dependent on moisture |
xfscanf (in,"%le %le %le %le", &lambda_dry, &lambda_wet, &sr_dry, &sr_wet); |
break; |
} |
default:{ |
print_err ("unknown model type is required", __FILE__, __LINE__, __func__); |
abort(); |
} |
} |
|
//specific heat calculation: |
switch (cps_type){ |
case 0:{//constant |
xfscanf (in,"%le", &cps0); |
break; |
} |
case 1:{//dependent on moisture |
xfscanf (in,"%le %le", &cps_dry, &cps_wet); |
break; |
} |
default:{ |
print_err ("unknown model type is required", __FILE__, __LINE__, __func__); |
abort(); |
} |
} |
|
//thermal expansion calculation: |
switch (betas_type){ |
case 0:{//constant |
xfscanf (in,"%le", &betas0); |
break; |
} |
case 1:{//dependent on moisture |
xfscanf (in,"%le %le", &betas_dry, &betas_wet); |
break; |
} |
default:{ |
print_err ("unknown model type is required", __FILE__, __LINE__, __func__); |
abort(); |
} |
} |
} |
|
|
444,11 → 463,17 |
*/ |
double con_hawf3mat::get_cdiff(double pw,double pg,double t) |
{ |
double cdiff; |
double cdiff,pgw; |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
//actualized - check if pgw>pg |
pgw = 0.0; |
get_pgw(pw,pg,t); |
if (pgw > pg) |
pg = pgw; |
|
//conversion pow(a,b) -> exp(b*log(a)) |
//cdiff = dv0*p0/pg*pow((t/t0),bv0); |
if (t < tcr) |
512,10 → 537,16 |
*/ |
double con_hawf3mat::get_sw(double pw, double pg, double t, long ipp) |
{ |
double sw; |
double sw,pc; |
sw = 0.0; |
pc = 0.0; |
|
switch (sr_type){ |
case bazant_sr:{//Bazant |
pc = pg - pw; |
sw = bazant_ret.sat(pc,t); |
break; |
} |
case lewis_and_schrefler_sr:{//Lewis and Schrefler's book |
sw = lewis_ret.sw(pw); |
break; |
571,10 → 602,15 |
*/ |
double con_hawf3mat::get_ds_dpc(double pw, double pg, double t, long ipp) |
{ |
double dsw_dpc; |
dsw_dpc = 0.0; |
double dsw_dpc,pc; |
dsw_dpc = pc = 0.0; |
|
switch (sr_type){ |
case bazant_sr:{//Bazant |
pc = pg - pw; |
dsw_dpc = bazant_ret.dsat_dpc(pc,t); |
break; |
} |
case lewis_and_schrefler_sr:{//Lewis and Schrefler's book |
dsw_dpc = -lewis_ret.dsw_dpw(pw); |
break; |
620,10 → 656,15 |
*/ |
double con_hawf3mat::get_ds_dt(double pw, double pg, double t, long ipp) |
{ |
double dsw_dt; |
dsw_dt = 0.0; |
double dsw_dt,pc; |
dsw_dt = pc = 0.0; |
|
switch (sr_type){ |
case bazant_sr:{//Bazant |
pc = pg - pw; |
dsw_dt = bazant_ret.dsat_dt(pc,t); |
break; |
} |
case lewis_and_schrefler_sr:{//Lewis and Schrefler's book |
dsw_dt = lewis_ret.dsw_dt(pw); |
break; |
652,6 → 693,7 |
abort(); |
} |
} |
|
return(dsw_dt); |
} |
|
766,6 → 808,11 |
krw = pow(sw,lambda_krw); |
break; |
} |
case 9:{//bazant |
sw = get_sw(pw,pg,t,ipp); |
krw = exp(10.0*log(sw)); |
break; |
} |
|
default:{ |
print_err ("unknown model type is required", __FILE__, __LINE__, __func__); |
789,7 → 836,8 |
*/ |
double con_hawf3mat::get_krg(double pw, double pg, double t, long ipp) |
{ |
double sw,krg; |
double krg; |
//double sw; |
krg=1.0; |
|
switch (krg_type){ |
798,17 → 846,17 |
break; |
} |
case 1:{//dependent on saturation degree |
sw = get_sw(pw,pg,t,ipp); |
//sw = get_sw(pw,pg,t,ipp); |
krg = 1.0; //not finished |
break; |
} |
case 2:{//dependent on saturation degree |
sw = get_sw(pw,pg,t,ipp); |
//sw = get_sw(pw,pg,t,ipp); |
krg = 1.0; //not finished |
break; |
} |
case 3:{//exponential |
sw = get_sw(pw,pg,t,ipp); |
//sw = get_sw(pw,pg,t,ipp); |
krg = 1.0; //not finished |
break; |
} |
960,21 → 1008,21 |
double con_hawf3mat::get_dpgw_dpc(double pw,double pg,double t) |
{ |
double dpgw_dpc,pgws,rhow,tt,pc; |
|
|
//critical point of water check |
if(t >= tcr) |
tt = tcr; |
tt = tcr; |
else |
tt = t; |
|
tt = t; |
|
pgws = get_pgws(tt); |
rhow = get_rhow(tt); |
|
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
pg = pg + p0; |
|
pc = pg - pw; |
|
|
dpgw_dpc = -pgws*exp(-pc*mw/rhow/gasr/tt)*mw/rhow/gasr/tt; |
|
return(dpgw_dpc); |
991,7 → 1039,7 |
*/ |
double con_hawf3mat::get_dpgw_dt(double pw,double pg,double t) |
{ |
double dpgw_dt,dpgws_dt,tt,pgws,rhow,pc; |
double dpgw_dt,dpgws_dt,tt,pgws,rhow,pc,drhow_dt; |
|
//critical point of water check |
if(t >= tcr) |
1004,8 → 1052,12 |
pc = pg - pw; |
dpgws_dt = get_dpgws_dt(tt); |
|
drhow_dt = get_drhow_dt(tt); |
|
if(t < tcr) |
dpgw_dt = dpgws_dt*exp(-pc*mw/rhow/gasr/tt) + pgws*exp(-pc*mw/rhow/gasr/tt)*(pc*mw/rhow/gasr/tt/tt); |
//actualized form: |
//dpgw_dt = dpgws_dt*exp(-pc*mw/rhow/gasr/tt) + pgws*exp(-pc*mw/rhow/gasr/tt)*(pc*mw/rhow/gasr/tt)*(drhow_dt/rhow + 1/tt); |
else |
dpgw_dt = dpgws_dt*exp(-pc*mw/rhow/gasr/tt); |
|
1138,12 → 1190,13 |
//d2rhow_dtdpc = 0.0; |
} |
|
rhow = rhow0;//testing??!! |
//rhow = rhow0;//testing??!! |
|
return(rhow); |
} |
|
|
|
/** |
function computes compresibility coefficient of water |
|
1257,6 → 1310,8 |
drhow_dt = 0.0; |
} |
|
drhow_dt = 0.0;//testing??!! |
|
return(drhow_dt); |
} |
|
1487,7 → 1542,6 |
rhoga = (pg - pgw)*ma/gasr/t; |
else |
rhoga = 0.0; |
//rhoga = (pg - pgw)*ma/gasr/t; |
|
return(rhoga); |
} |
1614,6 → 1668,7 |
break; |
} |
case 3:{// |
sw = get_sw(pw,pg,t,ipp); |
lambdaeff = pow(lambda_wet,sw)*pow(lambda_dry,(1-sw)); //this is from Patek's thesis p. 18, but it seems to be strange |
break; |
} |
1755,7 → 1810,11 |
k = get_ktg(pw,pg,t,ipp);// *scale_g;//scaling |
if((ri == 2) && (ci == 2)) |
k = get_ktt1(pw,pg,t,ipp);// *scale_t;//scaling |
|
|
//if(Tp->time >= 4.345400e+05) |
//fprintf(Outt,"matcond k = %e, ri = %ld, ci = %ld, ipp = %ld\n",k,ri,ci,ipp); |
//fflush(Outt); |
|
fillm(0.0,d); |
|
d[0][0] = k; d[0][1] = 0.0; |
1847,7 → 1906,11 |
c = get_captg(pw,pg,t,ipp);// *scale_g;//scaling |
if((ri == 2) && (ci == 2)) |
c = get_captt(pw,pg,t,ipp);// *scale_t;//scaling |
|
|
//if(Tp->time >= 4.345400e+05) |
//fprintf(Outt,"matcap c = %e, ri = %ld, ci = %ld, ipp = %ld, pw = %e, pg = %e, t = %e\n",c,ri,ci,ipp,pw,pg,t); |
//fflush(Outt); |
|
check_math_errel(ipp); |
} |
|
2201,6 → 2264,7 |
fwu = 0.0; |
|
switch (model_type){ |
case artificial3: |
case lewis_and_schrefler3: |
break; |
case lewis_and_schrefler3_mefel:{//Lewis and Schrefler's book |
2255,6 → 2319,7 |
fgu = 0.0; |
|
switch (model_type){ |
case artificial3: |
case lewis_and_schrefler3: |
break; |
case lewis_and_schrefler3_mefel:{//Lewis and Schrefler's book |
2306,6 → 2371,7 |
ftu = 0.0; |
|
switch (model_type){ |
case artificial3: |
case lewis_and_schrefler3: |
break; |
case lewis_and_schrefler3_mefel:{//Lewis and Schrefler's book |
2534,44 → 2600,31 |
@param ipp - number of integration point |
|
TKr, 18.5.2005 |
TKr 06/12/2022 - actualized |
*/ |
void con_hawf3mat::gaspress_check(double pw,double &pg,double t,long ipp) |
{ |
switch (model_type){ |
case lewis_and_schrefler3://Lewis and Schrefler's book |
case van_genuchten3://partially saturated medium = Van Genuchten model |
case lewis_and_schrefler3_mefel:{//Lewis and Schrefler's book, saturation degree ind its derivative are obtained from mefel; |
break; |
} |
case lewis_and_schrefler3_2:{//Lewis and Schrefler's book p. 381 |
break; |
} |
default:{ |
print_err ("unknown model type is required", __FILE__, __LINE__, __func__); |
abort(); |
} |
} |
|
/* double pgw; |
|
double pgw; |
|
//general |
//gas pressure check |
pgw = get_pgw(pw,pg,t); |
if (pgw > pg) |
pg = pgw; |
|
|
if (pgw > pg){ |
pg = pgw; |
//fprintf(Outt,"gas pressure was modified in integration point No. %ld\n",ipp); |
} |
|
//general |
//fully saturated state |
if(pw < 100.0) |
pg = 101325.0; |
|
if(pgw < 100.0){ |
pg = 101325.0; |
} |
|
// M. Starnoni 24-11-2010 |
if (pg<0.0){ |
pg = 0.0; |
print_err("gas pressure was modified",__FILE__,__LINE__,__func__); |
if (pg < 0.0){ |
print_err("gas pressure was modified in integration point No. %ld",__FILE__,__LINE__,__func__,ipp); |
pg = 0.0; |
} |
*/ |
} |
|
|
2585,33 → 2638,10 |
|
|
TKr, 18.5.2005 |
TKr 06/12/2022 - actualized |
*/ |
void con_hawf3mat::waterpress_check(double &pw,double pg,double t,long ipp) |
{ |
switch (model_type){ |
case lewis_and_schrefler3:{//Lewis and Schrefler's book |
//this must be checked: |
//if (pw >= 0.0) |
//pw = 0.0; |
|
//storing into gauss point |
//Tm->ip[ipp].av[0] = pw; |
break; |
} |
case van_genuchten3:{//partially saturated medium = Van Genuchten model |
break; |
} |
case lewis_and_schrefler3_mefel:{//Lewis and Schrefler's book, saturation degree ind its derivative are obtained from mefel; |
break; |
} |
case lewis_and_schrefler3_2:{//Lewis and Schrefler's book p. 381 |
break; |
} |
default:{ |
print_err ("unknown model type is required", __FILE__, __LINE__, __func__); |
abort(); |
} |
} |
} |
|
/** |
2620,13 → 2650,13 |
@param nv - array with variables |
|
*/ |
void con_hawf3mat::values_correction (vector &nv) |
void con_hawf3mat::values_correction (vector &nv, long ipp) |
{ |
// pore water pressure |
waterpress_check(nv[0],nv[1],nv[2],0); |
waterpress_check(nv[0],nv[1],nv[2],ipp); |
|
// gas pressure |
gaspress_check(nv[0],nv[1],nv[2],0); |
gaspress_check(nv[0],nv[1],nv[2],ipp); |
} |
|
|
2644,35 → 2674,39 |
double capww; |
double alpha,n,ks,sw,rhogw,sg,rhow,kw,dpgw_dpc,ds_dpc; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
alpha = get_alpha(pw,pg,t,ipp); |
n = get_porosity(pw,pg,t,ipp); |
ks = get_ks(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
rhogw = get_rhogw(pw,pg,t); |
sg = 1.0 - sw; |
rhow = get_rhow(t); |
kw = get_kw(pw,pg,t,ipp); |
dpgw_dpc = get_dpgw_dpc(pw,pg,t); |
ds_dpc = get_ds_dpc(pw,pg,t,ipp);//this is equal to minus derivative with respect to suction |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
if(compress == 1){ |
//compressible grains: |
capww = (alpha - n)/ks*sw*(rhogw*sg + rhow*sw) + rhow*sw*n/kw; |
capww = capww - sg*n*mw/t/gasr*dpgw_dpc; |
capww = capww - ((alpha - n)/ks*(rhogw*sg*(pw-pg) + rhow*sw*(pw-pg)) + n*(rhow - rhogw))*ds_dpc; |
} |
if(model_type == artificial3) |
capww = capww0; |
else{ |
//incompressible grains: |
capww = -sg*n*mw/t/gasr*dpgw_dpc; |
capww = capww - n*(rhow - rhogw)*ds_dpc; |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
alpha = get_alpha(pw,pg,t,ipp); |
n = get_porosity(pw,pg,t,ipp); |
ks = get_ks(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
rhogw = get_rhogw(pw,pg,t); |
sg = 1.0 - sw; |
rhow = get_rhow(t); |
kw = get_kw(pw,pg,t,ipp); |
dpgw_dpc = get_dpgw_dpc(pw,pg,t); |
ds_dpc = get_ds_dpc(pw,pg,t,ipp);//this is equal to minus derivative with respect to suction |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
if(compress == 1){ |
//compressible grains: |
capww = (alpha - n)/ks*sw*(rhogw*sg + rhow*sw) + rhow*sw*n/kw; |
capww = capww - sg*n*mw/t/gasr*dpgw_dpc; |
capww = capww - ((alpha - n)/ks*(rhogw*sg*(pw-pg) + rhow*sw*(pw-pg)) + n*(rhow - rhogw))*ds_dpc; |
} |
else{ |
//incompressible grains: |
capww = -sg*n*mw/t/gasr*dpgw_dpc; |
capww = capww - n*(rhow - rhogw)*ds_dpc; |
} |
} |
return(capww); |
} |
2693,35 → 2727,39 |
double capwg; |
double alpha,n,ks,sw,rhogw,sg,rhow,kw,dpgw_dpc,ds_dpc; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
alpha = get_alpha(pw,pg,t,ipp); |
n = get_porosity(pw,pg,t,ipp); |
ks = get_ks(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
rhogw = get_rhogw(pw,pg,t); |
sg = 1.0 - sw; |
rhow = get_rhow(t); |
kw = get_kw(pw,pg,t,ipp); |
dpgw_dpc = get_dpgw_dpc(pw,pg,t); |
ds_dpc = get_ds_dpc(pw,pg,t,ipp);//this is equal to minus derivative with respect to suction |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
if(compress == 1){ |
//compressible grains: |
capwg = (alpha - n)/ks*sg*(rhogw*sg + rhow*sw); |
capwg = capwg + sg*n*mw/t/gasr*dpgw_dpc; |
capwg = capwg + ((alpha - n)/ks*(rhogw*sg*(pw-pg) + rhow*sw*(pw-pg)) + n*(rhow - rhogw))*ds_dpc; |
} |
if(model_type == artificial3) |
capwg = capwg0; |
else{ |
//incompressible grains: |
capwg = sg*n*mw/t/gasr*dpgw_dpc; |
capwg = capwg + n*(rhow - rhogw)*ds_dpc; |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
alpha = get_alpha(pw,pg,t,ipp); |
n = get_porosity(pw,pg,t,ipp); |
ks = get_ks(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
rhogw = get_rhogw(pw,pg,t); |
sg = 1.0 - sw; |
rhow = get_rhow(t); |
kw = get_kw(pw,pg,t,ipp); |
dpgw_dpc = get_dpgw_dpc(pw,pg,t); |
ds_dpc = get_ds_dpc(pw,pg,t,ipp);//this is equal to minus derivative with respect to suction |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
if(compress == 1){ |
//compressible grains: |
capwg = (alpha - n)/ks*sg*(rhogw*sg + rhow*sw); |
capwg = capwg + sg*n*mw/t/gasr*dpgw_dpc; |
capwg = capwg + ((alpha - n)/ks*(rhogw*sg*(pw-pg) + rhow*sw*(pw-pg)) + n*(rhow - rhogw))*ds_dpc; |
} |
else{ |
//incompressible grains: |
capwg = sg*n*mw/t/gasr*dpgw_dpc; |
capwg = capwg + n*(rhow - rhogw)*ds_dpc; |
} |
} |
return(capwg); |
} |
2741,39 → 2779,42 |
double capwt; |
double betasgw,sw,sg,dpgw_dt,pgw,alpha,n,ks,rhogw,rhow,ds_dt; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
betasgw = get_betasgw(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
sg = 1.0 - sw; |
dpgw_dt = get_dpgw_dt(pw,pg,t); |
pgw = get_pgw(pw,pg,t); |
alpha = get_alpha(pw,pg,t,ipp); |
n = get_porosity(pw,pg,t,ipp); |
ks = get_ks(pw,pg,t,ipp); |
rhogw = get_rhogw(pw,pg,t); |
rhow = get_rhow(t); |
ds_dt = get_ds_dt(pw,pg,t,ipp); |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
if(compress == 1){ |
//compressible grains: |
capwt = -betasgw; |
capwt = capwt + sg*n*mw/t/gasr*(dpgw_dt-pgw/t); |
capwt = capwt + ((alpha - n)/ks*(rhogw*sg*(pw-pg) + rhow*sw*(pw-pg)) + n*(rhow - rhogw))*ds_dt; |
} |
if(model_type == artificial3) |
capwt = capwt0; |
else{ |
//incompressible grains: |
capwt = -betasgw; |
capwt = capwt + sg*n*mw/t/gasr*(dpgw_dt-pgw/t); |
capwt = capwt + n*(rhow - rhogw)*ds_dt; |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
betasgw = get_betasgw(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
sg = 1.0 - sw; |
dpgw_dt = get_dpgw_dt(pw,pg,t); |
pgw = get_pgw(pw,pg,t); |
alpha = get_alpha(pw,pg,t,ipp); |
n = get_porosity(pw,pg,t,ipp); |
ks = get_ks(pw,pg,t,ipp); |
rhogw = get_rhogw(pw,pg,t); |
rhow = get_rhow(t); |
ds_dt = get_ds_dt(pw,pg,t,ipp); |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
if(compress == 1){ |
//compressible grains: |
capwt = -betasgw; |
capwt = capwt + sg*n*mw/t/gasr*(dpgw_dt-pgw/t); |
capwt = capwt + ((alpha - n)/ks*(rhogw*sg*(pw-pg) + rhow*sw*(pw-pg)) + n*(rhow - rhogw))*ds_dt; |
} |
else{ |
//incompressible grains: |
capwt = -betasgw; |
capwt = capwt + sg*n*mw/t/gasr*(dpgw_dt-pgw/t); |
capwt = capwt + n*(rhow - rhogw)*ds_dt; |
} |
} |
|
return(capwt); |
} |
|
2791,26 → 2832,29 |
double kww; |
double rhow,krw,muw,rhog,mg,dg,dpgw_dpc,kintr; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
rhow = get_rhow(t); |
krw = get_krw(pw,pg,t,ipp); |
muw = get_muw(t); |
rhog =get_rhog(pw,pg,t); |
mg = get_mg(pw,pg,t); |
dg = get_dg(pw,pg,t,ipp); |
dpgw_dpc = get_dpgw_dpc(pw,pg,t); |
kintr = get_kintr(pw,pg,t,ipp); |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
kww = rhow*kintr*krw/muw;//water |
kww = kww - rhog*ma*mw/mg/mg*dg/pg*dpgw_dpc;//water vapour |
|
if(model_type == artificial3) |
kww = kww0; |
else{ |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
rhow = get_rhow(t); |
krw = get_krw(pw,pg,t,ipp); |
muw = get_muw(t); |
rhog =get_rhog(pw,pg,t); |
mg = get_mg(pw,pg,t); |
dg = get_dg(pw,pg,t,ipp); |
dpgw_dpc = get_dpgw_dpc(pw,pg,t); |
kintr = get_kintr(pw,pg,t,ipp); |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
kww = rhow*kintr*krw/muw;//water |
kww = kww - rhog*ma*mw/mg/mg*dg/pg*dpgw_dpc;//water vapour |
} |
return(kww); |
} |
|
2829,30 → 2873,33 |
double kwg; |
double rhogw,krg,mug,dg,pgw,rhow,krw,muw,rhog,mg,dpgw_dpc,kintr; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
rhogw = get_rhogw(pw,pg,t); |
krg = get_krg(pw,pg,t,ipp); |
mug = get_mug(pw,pg,t); |
dg = get_dg(pw,pg,t,ipp); |
pgw = get_pgw(pw,pg,t); |
rhow = get_rhow(t); |
krw = get_krw(pw,pg,t,ipp); |
muw = get_muw(t); |
rhog = get_rhog(pw,pg,t); |
mg = get_mg(pw,pg,t); |
dpgw_dpc = get_dpgw_dpc(pw,pg,t); |
kintr = get_kintr(pw,pg,t,ipp); |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
kwg = rhogw*kintr*krg/mug;//air |
kwg = kwg + rhog*ma*mw/mg/mg*dg*(dpgw_dpc/pg - pgw/pg/pg);//water vapour |
|
if(model_type == artificial3) |
kwg = kwg0; |
else{ |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
rhogw = get_rhogw(pw,pg,t); |
krg = get_krg(pw,pg,t,ipp); |
mug = get_mug(pw,pg,t); |
dg = get_dg(pw,pg,t,ipp); |
pgw = get_pgw(pw,pg,t); |
rhow = get_rhow(t); |
krw = get_krw(pw,pg,t,ipp); |
muw = get_muw(t); |
rhog = get_rhog(pw,pg,t); |
mg = get_mg(pw,pg,t); |
dpgw_dpc = get_dpgw_dpc(pw,pg,t); |
kintr = get_kintr(pw,pg,t,ipp); |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
kwg = rhogw*kintr*krg/mug;//air |
kwg = kwg + rhog*ma*mw/mg/mg*dg*(dpgw_dpc/pg - pgw/pg/pg);//water vapour |
} |
return(kwg); |
} |
|
2870,34 → 2917,38 |
{ |
double capgw,alpha,n,ks,sw,sg,pc,rhoga,ds_dpc,dpgw_dpc; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
alpha = get_alpha(pw,pg,t,ipp); |
n = get_porosity(pw,pg,t,ipp); |
ks = get_ks(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
sg = 1.0 - sw; |
pc = pg - pw; |
rhoga = get_rhoga(pw,pg,t); |
ds_dpc = get_ds_dpc(pw,pg,t,ipp); |
dpgw_dpc = get_dpgw_dpc(pw,pg,t);//this is equal to minus derivative with respect to suction |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
if(compress == 1){ |
//compressible grains: |
capgw = (alpha - n)/ks*sw*sg*rhoga; |
capgw = capgw + ((alpha - n)/ks*sg*pc + n)*rhoga*ds_dpc; |
capgw = capgw + sg*n*mw/t/gasr*dpgw_dpc; |
} |
if(model_type == artificial3) |
capgw = capgw0; |
else{ |
//incompressible grains: |
capgw = n*rhoga*ds_dpc; |
capgw = capgw + sg*n*mw/t/gasr*dpgw_dpc; |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
alpha = get_alpha(pw,pg,t,ipp); |
n = get_porosity(pw,pg,t,ipp); |
ks = get_ks(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
sg = 1.0 - sw; |
pc = pg - pw; |
rhoga = get_rhoga(pw,pg,t); |
ds_dpc = get_ds_dpc(pw,pg,t,ipp); |
dpgw_dpc = get_dpgw_dpc(pw,pg,t);//this is equal to minus derivative with respect to suction |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
if(compress == 1){ |
//compressible grains: |
capgw = (alpha - n)/ks*sw*sg*rhoga; |
capgw = capgw + ((alpha - n)/ks*sg*pc + n)*rhoga*ds_dpc; |
capgw = capgw + sg*n*mw/t/gasr*dpgw_dpc; |
} |
else{ |
//incompressible grains: |
capgw = n*rhoga*ds_dpc; |
capgw = capgw + sg*n*mw/t/gasr*dpgw_dpc; |
} |
} |
return(capgw); |
} |
2917,36 → 2968,40 |
{ |
double capgg,alpha,n,ks,sw,sg,pc,rhoga,ds_dpc,dpgw_dpc; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
alpha = get_alpha(pw,pg,t,ipp); |
n = get_porosity(pw,pg,t,ipp); |
ks = get_ks(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
sg = 1.0 - sw; |
rhoga = get_rhoga(pw,pg,t); |
pc = pg - pw; |
ds_dpc = get_ds_dpc(pw,pg,t,ipp);//this is equal to minus derivative with respect to suction |
dpgw_dpc = get_dpgw_dpc(pw,pg,t); |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
if(compress == 1){ |
//compressible grains: |
capgg = (alpha - n)/ks*sg*sg*rhoga; |
capgg = capgg - ((alpha - n)/ks*sg*pc + n)*rhoga*ds_dpc; |
capgg = capgg + sg*n*ma/t/gasr; |
capgg = capgg - sg*n*mw/t/gasr*dpgw_dpc; |
} |
if(model_type == artificial3) |
capgg = capgg0; |
else{ |
//incompressible grains: |
capgg = -n*rhoga*ds_dpc; |
capgg = capgg + sg*n*ma/t/gasr; |
capgg = capgg - sg*n*mw/t/gasr*dpgw_dpc; |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
alpha = get_alpha(pw,pg,t,ipp); |
n = get_porosity(pw,pg,t,ipp); |
ks = get_ks(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
sg = 1.0 - sw; |
rhoga = get_rhoga(pw,pg,t); |
pc = pg - pw; |
ds_dpc = get_ds_dpc(pw,pg,t,ipp);//this is equal to minus derivative with respect to suction |
dpgw_dpc = get_dpgw_dpc(pw,pg,t); |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
if(compress == 1){ |
//compressible grains: |
capgg = (alpha - n)/ks*sg*sg*rhoga; |
capgg = capgg - ((alpha - n)/ks*sg*pc + n)*rhoga*ds_dpc; |
capgg = capgg + sg*n*ma/t/gasr; |
capgg = capgg - sg*n*mw/t/gasr*dpgw_dpc; |
} |
else{ |
//incompressible grains: |
capgg = -n*rhoga*ds_dpc; |
capgg = capgg + sg*n*ma/t/gasr; |
capgg = capgg - sg*n*mw/t/gasr*dpgw_dpc; |
} |
} |
return(capgg); |
} |
2965,39 → 3020,42 |
{ |
double capgt,rhoga,betasg,alpha,n,ks,sw,sg,pc,ds_dt,dpgw_dt,pgw; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
rhoga = get_rhoga(pw,pg,t); |
betasg = get_betasg(pw,pg,t,ipp); |
alpha = get_alpha(pw,pg,t,ipp); |
n = get_porosity(pw,pg,t,ipp); |
ks = get_ks(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
sg = 1.0 - sw; |
pc = pg - pw; |
ds_dt = get_ds_dt(pw,pg,t,ipp); |
dpgw_dt = get_dpgw_dt(pw,pg,t); |
pgw = get_pgw(pw,pg,t); |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
if(compress == 1){ |
//compressible grains: |
capgt = -1.0*rhoga*betasg; |
capgt = capgt - ((alpha - n)/ks*sg*pc + n)*rhoga*ds_dt; |
capgt = capgt - (sg*n*ma/t/t/gasr*pg + sg*n*ma/t/gasr*(dpgw_dt-pgw/t)); |
} |
if(model_type == artificial3) |
capgt = capgt0; |
else{ |
//incompressible grains: |
capgt = -1.0*rhoga*betasg; |
capgt = capgt - n*rhoga*ds_dt; |
capgt = capgt - (sg*n*ma/t/t/gasr*pg + sg*n*ma/t/gasr*(dpgw_dt-pgw/t)); |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
rhoga = get_rhoga(pw,pg,t); |
betasg = get_betasg(pw,pg,t,ipp); |
alpha = get_alpha(pw,pg,t,ipp); |
n = get_porosity(pw,pg,t,ipp); |
ks = get_ks(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
sg = 1.0 - sw; |
pc = pg - pw; |
ds_dt = get_ds_dt(pw,pg,t,ipp); |
dpgw_dt = get_dpgw_dt(pw,pg,t); |
pgw = get_pgw(pw,pg,t); |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
if(compress == 1){ |
//compressible grains: |
capgt = -1.0*rhoga*betasg; |
capgt = capgt - ((alpha - n)/ks*sg*pc + n)*rhoga*ds_dt; |
capgt = capgt - (sg*n*ma/t/t/gasr*pg + sg*n*ma/t/gasr*(dpgw_dt-pgw/t)); |
} |
else{ |
//incompressible grains: |
capgt = -1.0*rhoga*betasg; |
capgt = capgt - n*rhoga*ds_dt; |
capgt = capgt - (sg*n*ma/t/t/gasr*pg + sg*n*ma/t/gasr*(dpgw_dt-pgw/t)); |
} |
} |
|
return(capgt); |
} |
|
3016,21 → 3074,24 |
{ |
double kgw,rhog,dg,mg,dpgw_dpc; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
rhog = get_rhog(pw,pg,t); |
dg = get_dg(pw,pg,t,ipp); |
dpgw_dpc = get_dpgw_dpc(pw,pg,t); |
mg = get_mg(pw,pg,t); |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
kgw = rhog*ma*mw/mg/mg*dg/pg*dpgw_dpc;//water vapour |
|
if(model_type == artificial3) |
kgw = kgw0; |
else{ |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
rhog = get_rhog(pw,pg,t); |
dg = get_dg(pw,pg,t,ipp); |
dpgw_dpc = get_dpgw_dpc(pw,pg,t); |
mg = get_mg(pw,pg,t); |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
kgw = rhog*ma*mw/mg/mg*dg/pg*dpgw_dpc;//water vapour |
} |
return(kgw); |
} |
|
3049,27 → 3110,30 |
{ |
double kgg,rhoga,krg,mug,rhog,dg,mg,dpgw_dpc,pgw,kintr; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
rhoga = get_rhoga(pw,pg,t); |
krg = get_krg(pw,pg,t,ipp); |
mug = get_mug(pw,pg,t); |
rhog = get_rhog(pw,pg,t); |
dg = get_dg(pw,pg,t,ipp); |
mg = get_mg(pw,pg,t); |
dpgw_dpc = get_dpgw_dpc(pw,pg,t); |
pgw = get_pgw(pw,pg,t); |
kintr = get_kintr(pw,pg,t,ipp); |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
kgg = rhoga*kintr*krg/mug;//air |
kgg = kgg - rhog*ma*mw/mg/mg*dg*(dpgw_dpc/pg-pgw/pg/pg);//water vapour |
|
if(model_type == artificial3) |
kgg = kgg0; |
else{ |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
rhoga = get_rhoga(pw,pg,t); |
krg = get_krg(pw,pg,t,ipp); |
mug = get_mug(pw,pg,t); |
rhog = get_rhog(pw,pg,t); |
dg = get_dg(pw,pg,t,ipp); |
mg = get_mg(pw,pg,t); |
dpgw_dpc = get_dpgw_dpc(pw,pg,t); |
pgw = get_pgw(pw,pg,t); |
kintr = get_kintr(pw,pg,t,ipp); |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
kgg = rhoga*kintr*krg/mug;//air |
kgg = kgg - rhog*ma*mw/mg/mg*dg*(dpgw_dpc/pg-pgw/pg/pg);//water vapour |
} |
return(kgg); |
} |
|
3108,33 → 3172,38 |
double con_hawf3mat::get_captw(double pw,double pg,double t,long ipp) |
{ |
double captw,dhvap,rhow,alpha,n,ks,sw,kw,dsw_dpc; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
dhvap = get_dhvap(t); |
rhow = get_rhow(t); |
alpha = get_alpha(pw,pg,t,ipp); |
n = get_porosity(pw,pg,t,ipp); |
ks = get_ks(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
kw = get_kw(pw,pg,t,ipp); |
dsw_dpc = get_ds_dpc(pw,pg,t,ipp);//this is equal to minus derivative with respect to suction |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
if(compress == 1){ |
//compressible grains: |
captw = -1.0*dhvap*rhow*((alpha - n)/ks*sw*sw + sw*n/kw); |
captw = captw + dhvap*rhow*((alpha - n)/ks*sw*pw - (alpha - n)/ks*sw*pg + n)*dsw_dpc; |
} |
|
if(model_type == artificial3) |
captw = captw0; |
else{ |
//incompressible grains: |
captw = dhvap*rhow*n*dsw_dpc; |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
dhvap = get_dhvap(t); |
rhow = get_rhow(t); |
alpha = get_alpha(pw,pg,t,ipp); |
n = get_porosity(pw,pg,t,ipp); |
ks = get_ks(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
kw = get_kw(pw,pg,t,ipp); |
dsw_dpc = get_ds_dpc(pw,pg,t,ipp);//this is equal to minus derivative with respect to suction |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
if(compress == 1){ |
//compressible grains: |
captw = -1.0*dhvap*rhow*((alpha - n)/ks*sw*sw + sw*n/kw); |
captw = captw + dhvap*rhow*((alpha - n)/ks*sw*pw - (alpha - n)/ks*sw*pg + n)*dsw_dpc; |
} |
else{ |
//incompressible grains: |
captw = dhvap*rhow*n*dsw_dpc; |
} |
} |
|
return(captw); |
} |
|
3190,34 → 3259,44 |
{ |
double captg,dhvap,rhow,alpha,n,ks,sw,sg,kw,dsw_dpc; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
if(model_type == artificial3) |
captg = captg0; |
else{ |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
dhvap = get_dhvap(t); |
rhow = get_rhow(t); |
alpha = get_alpha(pw,pg,t,ipp); |
n = get_porosity(pw,pg,t,ipp); |
ks = get_ks(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
sg = 1.0 - sw; |
kw = get_kw(pw,pg,t,ipp); |
dsw_dpc = get_ds_dpc(pw,pg,t,ipp);//this is equal to minus derivative with respect to suction |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
//if(Tp->time >= 4.345400e+05 && ipp == 25188){ |
//fprintf(Outt,"dhvap = %e, rhow = %e, alpha = %e, n = %e, ks = %e, sw = %e, sg = %e, kw = %e, dsw_dpc = %e, pw = %e, pg = %e, t = %e\n", |
// dhvap,rhow,alpha,n,ks,sw,sg,kw,dsw_dpc,pw,pg,t); |
//fflush(Outt); |
//} |
|
dhvap = get_dhvap(t); |
rhow = get_rhow(t); |
alpha = get_alpha(pw,pg,t,ipp); |
n = get_porosity(pw,pg,t,ipp); |
ks = get_ks(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
sg = 1.0 - sw; |
kw = get_kw(pw,pg,t,ipp); |
dsw_dpc = get_ds_dpc(pw,pg,t,ipp);//this is equal to minus derivative with respect to suction |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
if(compress == 1){ |
//compressible grains: |
captg = -dhvap*rhow*(alpha - n)/ks*sg*sw;//water |
captg = captg - dhvap*rhow*((alpha - n)/ks*sw*pw - (alpha - n)/ks*sw*pg + n)*dsw_dpc; |
if(compress == 1){ |
//compressible grains: |
captg = -dhvap*rhow*(alpha - n)/ks*sg*sw;//water |
captg = captg - dhvap*rhow*((alpha - n)/ks*sw*pw - (alpha - n)/ks*sw*pg + n)*dsw_dpc; |
} |
else{ |
//incompressible grains: |
captg = -dhvap*rhow*n*dsw_dpc; |
} |
} |
else{ |
//incompressible grains: |
captg = -dhvap*rhow*n*dsw_dpc; |
} |
|
|
return(captg); |
} |
|
3235,35 → 3314,39 |
{ |
double captt,rhocp,dhvap,betasw,rhow,alpha,n,ks,sw,dsw_dt; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
rhocp = get_rhocp(pw,pg,t,ipp); |
dhvap = get_dhvap(t); |
betasw = get_betasw(pw,pg,t,ipp); |
rhow = get_rhow(t); |
alpha = get_alpha(pw,pg,t,ipp); |
n = get_porosity(pw,pg,t,ipp); |
ks = get_ks(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
dsw_dt = get_ds_dt(pw,pg,t,ipp); |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
if(compress == 1){ |
//compressible grains: |
captt = rhocp; |
captt = captt + dhvap*betasw*rhow; |
captt = captt - dhvap*(rhow*((alpha - n)/ks*sw*pw - (alpha - n)/ks*sw*pg + n)*dsw_dt); |
} |
if(model_type == artificial3) |
captt = captt0; |
else{ |
//incompressible grains: |
captt = rhocp; |
captt = captt + dhvap*betasw*rhow; |
captt = captt - dhvap*n*dsw_dt; |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
rhocp = get_rhocp(pw,pg,t,ipp); |
dhvap = get_dhvap(t); |
betasw = get_betasw(pw,pg,t,ipp); |
rhow = get_rhow(t); |
alpha = get_alpha(pw,pg,t,ipp); |
n = get_porosity(pw,pg,t,ipp); |
ks = get_ks(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
dsw_dt = get_ds_dt(pw,pg,t,ipp); |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
if(compress == 1){ |
//compressible grains: |
captt = rhocp; |
captt = captt + dhvap*betasw*rhow; |
captt = captt - dhvap*(rhow*((alpha - n)/ks*sw*pw - (alpha - n)/ks*sw*pg + n)*dsw_dt); |
} |
else{ |
//incompressible grains: |
captt = rhocp; |
captt = captt + dhvap*betasw*rhow; |
captt = captt - dhvap*n*dsw_dt; |
} |
} |
return(captt); |
} |
3283,19 → 3366,24 |
{ |
double dhvap,rhow,krw,muw,ktw,kintr; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
dhvap = get_dhvap(t); |
rhow = get_rhow(t); |
krw = get_krw(pw,pg,t,ipp); |
muw = get_muw(t); |
kintr = get_kintr(pw,pg,t,ipp); |
if(model_type == artificial3) |
ktw = ktw0; |
else{ |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
dhvap = get_dhvap(t); |
rhow = get_rhow(t); |
krw = get_krw(pw,pg,t,ipp); |
muw = get_muw(t); |
kintr = get_kintr(pw,pg,t,ipp); |
|
//correction 02/12/2022 |
ktw = -dhvap*(rhow*kintr*krw/muw);//water |
} |
|
ktw = dhvap*(rhow*kintr*krw/muw);//water |
|
return(ktw); |
} |
|
3315,15 → 3403,19 |
{ |
double lambdaeff,ktt1; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
if(model_type == artificial3) |
ktt1 = ktt0; |
else{ |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
lambdaeff = get_lambdaeff(pw,pg,t,ipp); |
|
ktt1 = lambdaeff; |
} |
|
lambdaeff = get_lambdaeff(pw,pg,t,ipp); |
|
ktt1 = lambdaeff; |
|
return(ktt1); |
} |
|
3341,21 → 3433,24 |
{ |
double n,sw,cpw,rhow,krw,muw,ktt2a,kintr; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
n = get_porosity(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
rhow = get_rhow(t); |
cpw = get_cpw(); |
krw = get_krw(pw,pg,t,ipp); |
muw = get_muw(t); |
kintr = get_kintr(pw,pg,t,ipp); |
|
ktt2a = n*sw*rhow*cpw*kintr*krw/muw; |
|
if(model_type == artificial3) |
ktt2a = 0.0; |
else{ |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
n = get_porosity(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
rhow = get_rhow(t); |
cpw = get_cpw(); |
krw = get_krw(pw,pg,t,ipp); |
muw = get_muw(t); |
kintr = get_kintr(pw,pg,t,ipp); |
|
ktt2a = n*sw*rhow*cpw*kintr*krw/muw; |
} |
return(ktt2a); |
} |
|
3373,22 → 3468,26 |
double n,sw,sg,rhog,cpg,krg,mug,kintr; |
double ktt2b; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
if(model_type == artificial3) |
ktt2b = 0.0; |
else{ |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
n = get_porosity(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
sg = 1.0 - sw; |
rhog = get_rhog(pw,pg,t); |
cpg = get_cpg(pw,pg,t); |
krg = get_krg(pw,pg,t,ipp); |
mug = get_mug(pw,pg,t); |
kintr = get_kintr(pw,pg,t,ipp); |
|
ktt2b = n*sg*rhog*cpg*kintr*krg/mug; |
} |
|
n = get_porosity(pw,pg,t,ipp); |
sw = get_sw(pw,pg,t,ipp); |
sg = 1.0 - sw; |
rhog = get_rhog(pw,pg,t); |
cpg = get_cpg(pw,pg,t); |
krg = get_krg(pw,pg,t,ipp); |
mug = get_mug(pw,pg,t); |
kintr = get_kintr(pw,pg,t,ipp); |
|
ktt2b = n*sg*rhog*cpg*kintr*krg/mug; |
|
return(ktt2b); |
} |
|
3406,14 → 3505,18 |
{ |
double rhow,ktt2c; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
rhow = get_rhow(t); |
|
ktt2c = rhow; |
if(model_type == artificial3) |
ktt2c = 0.0; |
else{ |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
rhow = get_rhow(t); |
|
ktt2c = rhow; |
} |
|
return(ktt2c); |
} |
3432,15 → 3535,19 |
{ |
double rhog,ktt2d; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
if(model_type == artificial3) |
ktt2d = 0.0; |
else{ |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
rhog = get_rhog(pw,pg,t); |
|
ktt2d = rhog; |
} |
|
rhog = get_rhog(pw,pg,t); |
|
ktt2d = rhog; |
|
return(ktt2d); |
} |
|
3459,22 → 3566,25 |
double fc1; |
double rhow,rhogw,rhog,krw,muw,krg,mug,kintr; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
rhow = get_rhow(t); |
rhogw = get_rhogw(pw,pg,t); |
rhog = get_rhog(pw,pg,t); |
krw = get_krw(pw,pg,t,ipp); |
muw = get_muw(t); |
krg = get_krg(pw,pg,t,ipp); |
mug = get_mug(pw,pg,t); |
kintr = get_kintr(pw,pg,t,ipp); |
|
fc1 = rhogw*kintr*krg*rhog/mug + rhow*kintr*krw*rhow/muw; |
|
if(model_type == artificial3) |
fc1 = 0.0; |
else{ |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
rhow = get_rhow(t); |
rhogw = get_rhogw(pw,pg,t); |
rhog = get_rhog(pw,pg,t); |
krw = get_krw(pw,pg,t,ipp); |
muw = get_muw(t); |
krg = get_krg(pw,pg,t,ipp); |
mug = get_mug(pw,pg,t); |
kintr = get_kintr(pw,pg,t,ipp); |
|
fc1 = rhogw*kintr*krg*rhog/mug + rhow*kintr*krw*rhow/muw; |
} |
return(fc1); |
} |
|
3492,19 → 3602,24 |
double fg; |
double rhoga,rhog,krg,mug,kintr; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
if(model_type == artificial3) |
fg = 0.0; |
else{ |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
rhoga = get_rhoga(pw,pg,t); |
rhog = get_rhog(pw,pg,t); |
krg = get_krg(pw,pg,t,ipp); |
mug = get_mug(pw,pg,t); |
kintr = get_kintr(pw,pg,t,ipp); |
|
fg = rhoga*kintr*krg*rhog/mug; |
fg = 0.0;//no gravity force is included for the moist air |
} |
|
rhoga = get_rhoga(pw,pg,t); |
rhog = get_rhog(pw,pg,t); |
krg = get_krg(pw,pg,t,ipp); |
mug = get_mug(pw,pg,t); |
kintr = get_kintr(pw,pg,t,ipp); |
|
fg = rhoga*kintr*krg*rhog/mug; |
fg = 0.0;//no gravity force is included for the moist air |
return(fg); |
} |
|
3522,19 → 3637,23 |
double ft1; |
double dhvap,rhow,krw,muw,kintr; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
if(model_type == artificial3) |
ft1 = 0.0; |
else{ |
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
waterpress_check(pw,pg,t,ipp); |
|
dhvap = get_dhvap(t); |
rhow = get_rhow(t); |
krw = get_krw(pw,pg,t,ipp); |
muw = get_muw(t); |
kintr = get_kintr(pw,pg,t,ipp); |
|
ft1 = -dhvap*rhow*kintr*krw*rhow/muw; |
} |
|
dhvap = get_dhvap(t); |
rhow = get_rhow(t); |
krw = get_krw(pw,pg,t,ipp); |
muw = get_muw(t); |
kintr = get_kintr(pw,pg,t,ipp); |
|
ft1 = -dhvap*rhow*kintr*krw*rhow/muw; |
|
return(ft1); |
} |
|
3567,7 → 3686,7 |
pgws = get_pgws(t); |
rhow = get_rhow(t); |
|
trc=-1.0*pgws*mw*mw/gasr/gasr/t/t/rhow;//1st member of T. serie |
trc=1.0*pgws*mw*mw/gasr/gasr/t/t/rhow;//1st member of T. serie |
|
break; |
} |
3651,7 → 3770,7 |
//pgws = get_pgws(t); |
//rhow = get_rhow(t); |
|
trc = 1.0;//get_pwrh(bv,t);//tady??!! |
trc = -1.0;//get_pwrh(bv,t);//tady??!! |
//trc=-1.0*pgws*mw*mw/gasr/gasr/t/t/rhow;//1st member of T. serie |
|
break; |
3728,7 → 3847,8 |
double con_hawf3mat::get_transmission_nodval_ww(double bv,double pw,double pg,double t,long bc,long ipp) |
{ |
double new_nodval; |
|
state_eq tt; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
//water pressure check |
3736,7 → 3856,8 |
|
switch (bc){//type of prescribed variable |
case 30:{ |
new_nodval = 0.0; |
new_nodval = tt.get_pcrh(bv,t); |
new_nodval = pg - new_nodval; |
break; |
} |
case 40:{//simulation of free soil surface |
3765,7 → 3886,8 |
*/ |
double con_hawf3mat::get_transmission_flux_ww(double bv,double pw,double pg,double t,long bc,long ipp) |
{ |
double flux,trc; |
double flux,trc,pc; |
state_eq tt; |
|
//gas pressure check |
gaspress_check(pw,pg,t,ipp); |
3774,7 → 3896,11 |
|
switch (bc){//type of prescribed variable |
case 30:{//for testing - boundary flux |
flux = 0.0; |
pc = pg - pw; |
|
flux = tt.get_pcrh(bv,t); |
flux = -1.0*(flux - pc); |
|
break; |
} |
case 40:{//simulation of free soil surface |
3938,8 → 4064,12 |
other = r[0]; |
break; |
} |
case 6:{//moisture content |
other = get_w(r[0],r[1],r[2],ipp); |
//case 6:{//moisture content |
//other = get_w(r[0],r[1],r[2],ipp); |
//break; |
//} |
case 6:{//relative humidity |
other = get_pgw(r[0],r[1],r[2])/get_pgws(r[2]); |
break; |
} |
case 7:{//dsw_dpw |
3991,8 → 4121,12 |
fprintf (out,"Pore water pressure (Pa) "); |
break; |
} |
case 6:{//moisture content |
fprintf (out,"Moisture content (kg/kg) "); |
//case 6:{//moisture content |
//fprintf (out,"Moisture content (kg/kg) "); |
//break; |
//} |
case 6:{//relative humidity |
fprintf (out,"Relative humidity () "); |
break; |
} |
case 7:{//dS_dpw |
4137,25 → 4271,6 |
|
pg = Tm->ip[ipp].av[1]; |
|
switch (model_type){ |
case lewis_and_schrefler3:{//Lewis and Schrefler's book |
//this must be checked: |
if (pg <= 0.0) |
pg = 0.0; |
|
//storing into gauss point |
Tm->ip[ipp].av[1] = pg; |
break; |
} |
case van_genuchten3://partially saturated medium = Van Genuchten model |
case lewis_and_schrefler3_mefel:{//Lewis and Schrefler's book, saturation degree ind its derivative are obtained from mefel |
break; |
} |
default:{ |
print_err ("unknown model type is required", __FILE__, __LINE__, __func__); |
abort(); |
} |
} |
return(pg); |
} |
|
4175,6 → 4290,9 |
pw = Tm->ip[ipp].av[0]; |
pg = Tm->ip[ipp].av[1]; |
|
if(rel_gas_press == 1) |
pg = pg + p0; |
|
switch (model_type){ |
case lewis_and_schrefler3:{//Lewis and Schrefler's book |
suction = -1.0*(pg - pw);//this is correct, because pore water pressure is negative |
4309,6 → 4427,7 |
void con_hawf3mat::give_reqntq(long *antq) |
{ |
switch (model_type){ |
case artificial3: |
case lewis_and_schrefler3:{//Lewis and Schrefler's book |
break; |
} |