00001 #include "globalt.h"
00002 #include "devriesmat.h"
00003 #include "kunzel.h"
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013 static double T0 = 273.15;
00014 static double Rvapd = 462.0;
00015 static double cwatd = 4180;
00016 static double cvapd = 1617;
00017 static double rold = 1000;
00018
00019 devriesmat::devriesmat ()
00020 {}
00021
00022 devriesmat::~devriesmat ()
00023 {}
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036 void devriesmat::matcond (matrix &d,long ri,long ci,long ipp)
00037 {
00038 long n;
00039 n = d.n;
00040
00041 switch (n){
00042 case 1:{
00043 matcond1d (d,ri,ci,ipp);
00044 break;
00045 }
00046 case 2:{
00047 matcond2d (d,ri,ci,ipp);
00048 break;
00049 }
00050 case 3:{
00051 matcond3d (d,ri,ci,ipp);
00052 break;
00053 }
00054 default:{
00055 fprintf (stderr,"\n unknown number of components of stress tensor is required");
00056 fprintf (stderr,"\n in function (%s, line %d).\n",__FILE__,__LINE__);
00057 }
00058 }
00059 }
00060
00061
00062
00063
00064
00065
00066
00067
00068
00069
00070 void devriesmat::matcond1d (matrix &d,long ri,long ci,long ipp)
00071 {
00072 double k;
00073 double x1,x2,x3;
00074 k = 0.0;
00075
00076 x1 = Tm->ip[ipp].av[0];
00077 x2 = Tm->ip[ipp].av[1];
00078 x3 = Tm->ip[ipp].av[2];
00079
00080 if((ri == 0) && (ci == 0))
00081 k = k11 (x1,x2,x3);
00082 if((ri == 0) && (ci == 1))
00083 k = k12 (x1,x2,x3);
00084 if((ri == 0) && (ci == 2))
00085 k = k13 (x1,x2,x3);
00086
00087 if((ri == 1) && (ci == 0))
00088 k = k21 (x1,x2,x3);
00089 if((ri == 1) && (ci == 1))
00090 k = k22 (x1,x2,x3);
00091 if((ri == 1) && (ci == 2))
00092 k = k23 (x1,x2,x3);
00093
00094 if((ri == 2) && (ci == 0))
00095 k = k31 (x1,x2,x3);
00096 if((ri == 2) && (ci == 1))
00097 k = k32 (x1,x2,x3);
00098 if((ri == 2) && (ci == 2))
00099 k = k33 (x1,x2,x3);
00100
00101 d[0][0] = k;
00102 }
00103
00104
00105
00106
00107
00108
00109
00110
00111
00112 void devriesmat::matcond2d (matrix &d,long ri,long ci,long ipp)
00113 {
00114 double k;
00115 double x1,x2,x3;
00116 k = 0.0;
00117
00118 x1 = Tm->ip[ipp].av[0];
00119 x2 = Tm->ip[ipp].av[1];
00120 x3 = Tm->ip[ipp].av[2];
00121
00122 if((ri == 0) && (ci == 0))
00123 k = k11 (x1,x2,x3);
00124 if((ri == 0) && (ci == 1))
00125 k = k12 (x1,x2,x3);
00126 if((ri == 0) && (ci == 2))
00127 k = k13 (x1,x2,x3);
00128
00129 if((ri == 1) && (ci == 0))
00130 k = k21 (x1,x2,x3);
00131 if((ri == 1) && (ci == 1))
00132 k = k22 (x1,x2,x3);
00133 if((ri == 1) && (ci == 2))
00134 k = k23 (x1,x2,x3);
00135
00136 if((ri == 2) && (ci == 0))
00137 k = k31 (x1,x2,x3);
00138 if((ri == 2) && (ci == 1))
00139 k = k32 (x1,x2,x3);
00140 if((ri == 2) && (ci == 2))
00141 k = k33 (x1,x2,x3);
00142
00143 fillm(0.0,d);
00144
00145 d[0][0] = k; d[0][1] = 0.0;
00146 d[1][0] = 0.0; d[1][1] = k;
00147 }
00148
00149
00150
00151
00152
00153
00154
00155
00156
00157
00158 void devriesmat::matcond3d (matrix &d,long ri,long ci,long ipp)
00159 {
00160 double k;
00161 double x1,x2,x3;
00162 k = 0.0;
00163
00164 x1 = Tm->ip[ipp].av[0];
00165 x2 = Tm->ip[ipp].av[1];
00166 x3 = Tm->ip[ipp].av[2];
00167
00168 if((ri == 0) && (ci == 0))
00169 k = k11 (x1,x2,x3);
00170 if((ri == 0) && (ci == 1))
00171 k = k12 (x1,x2,x3);
00172 if((ri == 0) && (ci == 2))
00173 k = k13 (x1,x2,x3);
00174
00175 if((ri == 1) && (ci == 0))
00176 k = k21 (x1,x2,x3);
00177 if((ri == 1) && (ci == 1))
00178 k = k22 (x1,x2,x3);
00179 if((ri == 1) && (ci == 2))
00180 k = k23 (x1,x2,x3);
00181
00182 if((ri == 2) && (ci == 0))
00183 k = k31 (x1,x2,x3);
00184 if((ri == 2) && (ci == 1))
00185 k = k32 (x1,x2,x3);
00186 if((ri == 2) && (ci == 2))
00187 k = k33 (x1,x2,x3);
00188
00189 fillm(0.0,d);
00190
00191 d[0][0]=k; d[0][1]=0.0; d[0][2]=0.0;
00192 d[1][0]=0.0; d[1][1]=k; d[1][2]=0.0;
00193 d[2][0]=0.0; d[2][1]=0.0; d[2][2]=k;
00194 }
00195
00196
00197
00198
00199
00200
00201
00202
00203
00204
00205 void devriesmat::matcap (double &c,long ri,long ci,long ipp)
00206 {
00207 double x1,x2,x3;
00208 c = 0.0;
00209
00210 x1 = Tm->ip[ipp].av[0];
00211 x2 = Tm->ip[ipp].av[1];
00212 x3 = Tm->ip[ipp].av[2];
00213
00214 if((ri == 0) && (ci == 0))
00215 c = c11 (x1,x2,x3);
00216 if((ri == 0) && (ci == 1))
00217 c = c12 (x1,x2,x3);
00218 if((ri == 0) && (ci == 2))
00219 c = c13 (x1,x2,x3);
00220
00221 if((ri == 1) && (ci == 0))
00222 c = c21 (x1,x2,x3);
00223 if((ri == 1) && (ci == 1))
00224 c = c22 (x1,x2,x3);
00225 if((ri == 1) && (ci == 2))
00226 c = c23 (x1,x2,x3);
00227
00228 if((ri == 2) && (ci == 0))
00229 c = c31 (x1,x2,x3);
00230 if((ri == 2) && (ci == 1))
00231 c = c32 (x1,x2,x3);
00232 if((ri == 2) && (ci == 2))
00233 c = c33 (x1,x2,x3);
00234 }
00235
00236
00237
00238
00239
00240
00241
00242
00243 void devriesmat::read(XFILE *in)
00244 {
00245 int i, j;
00246
00247 for (i = 2; i<19; i++)
00248 {
00249 xfscanf(in, "%d",&MatChar[i]);
00250 }
00251 xfscanf(in, "%d",&MatChar[19]);
00252 for (i = 2;i <= 19; i++)
00253 {
00254
00255 switch (MatChar[i]){
00256 case 0:{
00257 break;
00258 }
00259 case 1:{
00260 xfscanf(in, "%lf",&MatConst[i]);
00261
00262 break;
00263 }
00264 case 2:{
00265 if( MatChar[i-1] ==2)
00266 {
00267 xfscanf(in, "%lf %lf",&MatData[i][0][0],&MatData[i][1][0]);
00268 }
00269 else
00270 {
00271 xfscanf(in, "%lf %lf",&MatData[i][0][0],&MatData[i][1][0]);
00272 }
00273
00274 switch (int(MatData[i][1][0])){
00275 case 2:{
00276 for (j=1; j<=MatData[i][0][0];j++)
00277 {
00278 xfscanf(in, "%lf %lf",&MatData[i][0][j],&MatData[i][1][j] );
00279 }
00280 break;
00281 }
00282 case 3:{
00283 for (j=1; j<=MatData[i][0][0];j++)
00284 {
00285 xfscanf(in, "%lf %lf %lf",&MatData[i][0][j],&MatData[i][1][j],&MatData[i][2][j] );
00286 }
00287 break;
00288 }
00289 default:{
00290 fprintf (stderr,"\n\n unknown definition of Material KUNZEL-DATA is required");
00291 fprintf (stderr,"\n in function read (file %s, line %d)\n",__FILE__,__LINE__);
00292 }
00293 }
00294 xfscanf(in, "%d",&MatData[i][0][0]);
00295
00296 break;
00297 }
00298 case 3:{
00299
00300 break;
00301 }
00302 case 30:{
00303 if( MatChar[i-1] ==1)
00304 {
00305 xfscanf(in, "%lf %lf",&MatFunce[i][0],&MatFunce[i][1]);
00306 }
00307 else
00308 {
00309 xfscanf(in, "%lf %lf",&MatFunce[i][0],&MatFunce[i][1]);
00310 }
00311 break;
00312 }
00313 case 31:{
00314 if( MatChar[i-1] ==1)
00315 {
00316 if (i == 6)
00317 {
00318 xfscanf(in, "%lf %lf %lf",&MatFunce[i][0],&MatFunce[i][1], &MatFunce[i][2]);
00319 }
00320 else
00321 {
00322 xfscanf(in, "%lf %lf",&MatFunce[i][0],&MatFunce[i][1]);
00323 }
00324 }
00325 else
00326 {
00327 if (i == 6)
00328 {
00329 xfscanf(in, "%lf %lf %lf",&MatFunce[i][0],&MatFunce[i][1], &MatFunce[i][2]);
00330 }
00331 else
00332 {
00333 xfscanf(in, "%lf %lf",&MatFunce[i][0],&MatFunce[i][1]);
00334 }
00335 }
00336
00337 break;
00338 }
00339 case 32:{
00340 if( MatChar[i-1] ==1)
00341 {
00342 xfscanf(in, "%lf %lf",&MatFunce[i][0],&MatFunce[i][1]);
00343 }
00344 else
00345 {
00346 xfscanf(in, "%lf %lf",&MatFunce[i][0],&MatFunce[i][1]);
00347 }
00348 break;
00349 }
00350 default:{
00351 fprintf (stderr,"\n\n unknown definition of Material KUNZEL is required");
00352 fprintf (stderr,"\n in function read (file %s, line %d)\n",__FILE__,__LINE__);
00353 }
00354 }
00355 }
00356
00357 madripom = 0;
00358 }
00359
00360 void devriesmat::print(FILE *out)
00361 {
00362 fprintf (out,"\n devries material \n");
00363 }
00364
00365
00366
00367
00368
00369
00370
00371
00372
00373
00374
00375
00376
00377
00378
00379
00380
00381
00382
00383
00384
00385
00386 double devriesmat::k11(double x1,double x2,double x3)
00387 {
00388 double k11;
00389 double dwl, dwv;
00390 dwl = coeff_of_water_diffusion_by_moisture_grad(x1,x2,x3);
00391 dwv = coeff_dwv(x1,x2,x3);
00392 k11 = dwl + dwv;
00393
00394 return (k11*rold);
00395 }
00396
00397 double devriesmat::k12(double x1,double x2,double x3)
00398 {
00399 double k12;
00400
00401 double dtv, dtl;
00402 dtl = coeff_of_water_diffusion_by_temperature_grad(x1,x2,x3);
00403 dtv =coeff_dtv(x1,x2,x3);
00404
00405 k12 = dtv + dtl;
00406
00407 return (k12*rold);
00408 }
00409
00410 double devriesmat::k13(double x1,double x2,double x3)
00411 {
00412 double k13;
00413 k13 = 0.0;
00414 return (k13);
00415 }
00416
00417 double devriesmat::k21(double x1,double x2,double x3)
00418 {
00419 double k21;
00420 double dwl, dwv, lv;
00421 dwl =coeff_of_water_diffusion_by_moisture_grad(x1,x2,x3);
00422 dwv = coeff_dwv(x1,x2,x3);
00423 lv = latent_heat_of_evaporation_of_water(x1,x2,x3);
00424
00425
00426
00427 k21 = rold*cwatd*dwl *(x2-T0)+rold*dwv*(cvapd*(x2- T0)+lv);
00428
00429
00430 return (-k21);
00431 }
00432
00433 double devriesmat::k22(double x1,double x2,double x3)
00434 {
00435 double k22;
00436 double lambda, dtl, dtv, lv;
00437 dtl = coeff_of_water_diffusion_by_temperature_grad(x1,x2,x3);
00438 dtv = coeff_dtv(x1,x2,x3);
00439 CorD(10,kd,1.0,0,x1,lambda,a2,a3);
00440 lv = latent_heat_of_evaporation_of_water(x1,x2,x3);
00441
00442
00443 k22 = lambda-rold*cwatd*dtl*(x2-T0)-rold*dtv*(cvapd*(x2-T0)+lv);
00444
00445
00446 return (k22);
00447 }
00448
00449 double devriesmat::k23(double x1,double x2,double x3)
00450 {
00451 double k23;
00452
00453 k23 = 0.0;
00454
00455 return (k23);
00456 }
00457
00458 double devriesmat::k31(double x1,double x2,double x3)
00459 {
00460 double k31;
00461
00462 k31 = 0.0;
00463
00464 return (k31);
00465 }
00466
00467 double devriesmat::k32(double x1,double x2,double x3)
00468 {
00469 double k32;
00470
00471 k32 = 0.0;
00472
00473 return (k32);
00474 }
00475
00476 double devriesmat::k33(double x1,double x2,double x3)
00477 {
00478 double k33;
00479
00480 k33 = 0.0;
00481
00482 return (k33);
00483 }
00484
00485 double devriesmat::c11(double x1,double x2,double x3)
00486 {
00487 double c11;
00488 double a, drovdw,rov;
00489
00490 a = relative_volume_ration_a(x1,x2,x3);
00491 drovdw = derivation_drov_dw(x1,x2,x3);
00492 rov =partial_density_of_water_vapor(x1,x2,x3);
00493
00494 c11 = 1+a*drovdw/rold-rov/rold;
00495
00496 return(c11*rold);
00497
00498 }
00499
00500 double devriesmat::c12(double x1,double x2,double x3)
00501 {
00502 double c12;
00503 c12 = 0.0;
00504 return(c12);
00505 }
00506
00507 double devriesmat::c13(double x1,double x2,double x3)
00508 {
00509 double c13;
00510 c13 = 0.0;
00511
00512 return(c13);
00513 }
00514
00515 double devriesmat::c21(double x1,double x2,double x3)
00516 {
00517 double c21;
00518 double a,lv,rov,drovdw;
00519
00520 a = relative_volume_ration_a(x1,x2,x3);
00521 drovdw = derivation_drov_dw(x1,x2,x3);
00522 rov =partial_density_of_water_vapor(x1,x2,x3);
00523
00524
00525
00526 lv =latent_heat_of_evaporation_of_water(x1,x2,x3);
00527
00528
00529 c21 = rold*cwatd*(x2-T0)+a*drovdw*(cvapd+lv)-rov*(cvapd*(x2-T0)+lv);
00530
00531
00532
00533 return(c21);
00534 }
00535
00536 double devriesmat::c22(double x1,double x2,double x3)
00537 {
00538 double c22;
00539 double ros,cstar, cs,a,rov;
00540
00541 a = relative_volume_ration_a(x1,x2,x3);
00542 rov =partial_density_of_water_vapor(x1,x2,x3);
00543 CorD(2,kd,0,0,x1,ros,a2,a3);
00544 CorD(9,kd,0,0,x1,cs,a2,a3);
00545 cstar = specific_heat_capacities_star(x1,x2,x3);
00546
00547
00548 c22 =ros*cs+rold*x1*cwatd+a*rov*cvapd;
00549
00550
00551 return(c22);
00552 }
00553
00554 double devriesmat::c23(double x1,double x2,double x3)
00555 {
00556 double c23;
00557 c23 = 0.0;
00558
00559 return(c23);
00560 }
00561
00562 double devriesmat::c31(double x1,double x2,double x3)
00563 {
00564 double c31;
00565 c31 = 0.0;
00566 return(c31);
00567 }
00568
00569 double devriesmat::c32(double x1,double x2,double x3)
00570 {
00571 double c32;
00572 c32 = 0.0;
00573 return(c32);
00574 }
00575
00576 double devriesmat::c33(double x1,double x2,double x3)
00577 {
00578 double c33;
00579 c33 = 0.0;
00580
00581 return(c33);
00582 }
00583
00584
00585
00586 void devriesmat::auxiliarydata (double x1,double x2,double x3)
00587 {}
00588
00589
00590
00591
00592
00593
00594
00595
00596
00597
00598
00599
00600
00601
00602
00603
00604
00605 double devriesmat::transmission_transcoeff(double trc,long ri,long ci,long nn,long bc,long ipp)
00606 {
00607 long k;
00608 double new_trc,x1,x2,x3;
00609 new_trc = 0.0;
00610
00611 k=Gtt->give_dof(nn,0);
00612 if (k>0) {x1 = Lsrst->lhs[k-1]+Lsrst->lhsi[k-1];}
00613 if (k==0) {x1 = 0.0;}
00614 if (k<0) {x1 = Tb->lc[0].pv[0-k-1].getval();}
00615 k=Gtt->give_dof(nn,1);
00616 if (k>0) {x2 = Lsrst->lhs[k-1]+Lsrst->lhsi[k-1];}
00617 if (k==0) {x2 = 0.0;}
00618 if (k<0) {x2 = Tb->lc[0].pv[0-k-1].getval();}
00619 k=Gtt->give_dof(nn,2);
00620 if (k>0) {x3 = Lsrst->lhs[k-1]+Lsrst->lhsi[k-1];}
00621 if (k==0) {x3 = 0.0;}
00622 if (k<0) {x3 = Tb->lc[0].pv[0-k-1].getval();}
00623
00624 if((ri == 0) && (ci == 0))
00625 new_trc = get_transmission_transcoeff_11(x1,x2,x3,bc,ipp);
00626 if((ri == 0) && (ci == 1))
00627 new_trc = 0.0;
00628 if((ri == 0) && (ci == 2))
00629 new_trc = 0.0;
00630
00631 if((ri == 1) && (ci == 0))
00632 new_trc = 0.0;
00633 if((ri == 1) && (ci == 1))
00634 new_trc = 0.0;
00635 if((ri == 1) && (ci == 2))
00636 new_trc = 0.0;
00637
00638 if((ri == 2) && (ci == 0))
00639 new_trc = 0.0;
00640 if((ri == 2) && (ci == 1))
00641 new_trc = 0.0;
00642 if((ri == 2) && (ci == 2))
00643 new_trc = 0.0;
00644
00645 new_trc = new_trc*trc;
00646
00647 return (new_trc);
00648 }
00649
00650
00651
00652
00653
00654
00655
00656
00657
00658
00659
00660
00661 double devriesmat::transmission_nodval(double nodval,long ri,long ci,long nn,long bc,long ipp)
00662 {
00663 long k;
00664 double new_nodval,x1,x2,x3;
00665 new_nodval = 0.0;
00666
00667 k=Gtt->give_dof(nn,0);
00668 if (k>0) {x1 = Lsrst->lhs[k-1]+Lsrst->lhsi[k-1];}
00669 if (k==0) {x1 = 0.0;}
00670 if (k<0) {x1 = Tb->lc[0].pv[0-k-1].getval();}
00671 k=Gtt->give_dof(nn,1);
00672 if (k>0) {x2 = Lsrst->lhs[k-1]+Lsrst->lhsi[k-1];}
00673 if (k==0) {x2 = 0.0;}
00674 if (k<0) {x2 = Tb->lc[0].pv[0-k-1].getval();}
00675 k=Gtt->give_dof(nn,2);
00676 if (k>0) {x3 = Lsrst->lhs[k-1]+Lsrst->lhsi[k-1];}
00677 if (k==0) {x3 = 0.0;}
00678 if (k<0) {x3 = Tb->lc[0].pv[0-k-1].getval();}
00679
00680 if((ri == 0) && (ci == 0))
00681 new_nodval = get_transmission_nodval_11(nodval,x1,x2,x3,bc,ipp);
00682 if((ri == 0) && (ci == 1))
00683 new_nodval = 0.0;
00684 if((ri == 0) && (ci == 2))
00685 new_nodval = 0.0;
00686
00687 if((ri == 1) && (ci == 0))
00688 new_nodval = 0.0;
00689 if((ri == 1) && (ci == 1))
00690 new_nodval = 0.0;
00691 if((ri == 1) && (ci == 2))
00692 new_nodval = 0.0;
00693
00694 if((ri == 2) && (ci == 0))
00695 new_nodval = 0.0;
00696 if((ri == 2) && (ci == 1))
00697 new_nodval = 0.0;
00698 if((ri == 2) && (ci == 2))
00699 new_nodval = 0.0;
00700
00701 return (new_nodval);
00702 }
00703
00704
00705
00706
00707
00708
00709
00710
00711
00712
00713
00714
00715
00716 double devriesmat::transmission_flux(double nodval,long ri,long ci,long nn,long bc,long ipp)
00717 {
00718 long k;
00719 double flux,x1,x2,x3;
00720 flux = 0.0;
00721
00722 k=Gtt->give_dof(nn,0);
00723 if (k>0) {x1 = Lsrst->lhs[k-1]+Lsrst->lhsi[k-1];}
00724 if (k==0) {x1 = 0.0;}
00725 if (k<0) {x1 = Tb->lc[0].pv[0-k-1].getval();}
00726 k=Gtt->give_dof(nn,1);
00727 if (k>0) {x2 = Lsrst->lhs[k-1]+Lsrst->lhsi[k-1];}
00728 if (k==0) {x2 = 0.0;}
00729 if (k<0) {x2 = Tb->lc[0].pv[0-k-1].getval();}
00730 k=Gtt->give_dof(nn,2);
00731 if (k>0) {x3 = Lsrst->lhs[k-1]+Lsrst->lhsi[k-1];}
00732 if (k==0) {x3 = 0.0;}
00733 if (k<0) {x3 = Tb->lc[0].pv[0-k-1].getval();}
00734
00735 if((ri == 0) && (ci == 0))
00736 flux = get_transmission_flux_11(nodval,x1,x2,x3,bc,ipp);
00737 if((ri == 0) && (ci == 1))
00738 flux = 0.0;
00739 if((ri == 0) && (ci == 2))
00740 flux = 0.0;
00741
00742 if((ri == 1) && (ci == 0))
00743 flux = 0.0;
00744 if((ri == 1) && (ci == 1))
00745 flux = 0.0;
00746 if((ri == 1) && (ci == 2))
00747 flux = 0.0;
00748
00749 if((ri == 2) && (ci == 0))
00750 flux = 0.0;
00751 if((ri == 2) && (ci == 1))
00752 flux = 0.0;
00753 if((ri == 2) && (ci == 2))
00754 flux = 0.0;
00755
00756 return (flux);
00757 }
00758
00759
00760
00761
00762
00763
00764
00765
00766
00767
00768 double devriesmat::get_transmission_nodval_11(double bv,double x1,double x2,double x3,long bc,long ipp)
00769 {
00770 double nodval;
00771
00772 switch (bc){
00773 case 30:{
00774 nodval = bv;
00775 break;
00776 }
00777 default:{
00778 fprintf (stderr,"\n\n No real boundary condition is prescribed (%s, line %d).\n",__FILE__,__LINE__);
00779 exit(0);
00780 }
00781 }
00782
00783 return(nodval);
00784 }
00785
00786
00787
00788
00789
00790
00791
00792
00793
00794 double devriesmat::get_transmission_transcoeff_11(double x1,double x2,double x3,long bc,long ipp)
00795 {
00796 double trc;
00797
00798 switch (bc){
00799 case 30:{
00800 trc=1.0;
00801 break;
00802 }
00803 default:{
00804 fprintf (stderr,"\n\n No real boundary condition is prescribed (%s, line %d).\n",__FILE__,__LINE__);
00805 exit(0);
00806 }
00807 }
00808
00809 return(trc);
00810 }
00811
00812
00813
00814
00815
00816
00817
00818
00819
00820
00821 double devriesmat::get_transmission_flux_11(double bv,double x1,double x2,double x3,long bc,long ipp)
00822 {
00823 double flux;
00824
00825 switch (bc){
00826 case 30:{
00827 flux = (bv - x1);
00828 break;
00829 }
00830 default:{
00831 fprintf (stderr,"\n\n No real boundary condition is prescribed (%s, line %d).\n",__FILE__,__LINE__);
00832 exit(0);
00833 }
00834 }
00835
00836 return(flux);
00837 }
00838
00839
00840
00841
00842
00843
00844
00845
00846 double devriesmat::get_othervalue(long compother,long ipp, double x1,double x2,double x3)
00847 {
00848 double other;
00849
00850 switch (compother){
00851 case 0:{
00852 other = x1;
00853 break;
00854 }
00855 case 1:{
00856 other = x2;
00857 break;
00858 }
00859 case 2:{
00860
00861 double relhum, a1;
00862 sorption_izoterms_giva_data(0,x1,x2,x3,relhum,a1);
00863 other = relhum;
00864 break;
00865 }
00866 default:{
00867 fprintf (stderr,"\n\n unknown type of component is required in function (%s, line %d).\n",__FILE__,__LINE__);
00868
00869 }
00870 }
00871 return (other);
00872
00873 }
00874
00875
00876
00877
00878
00879
00880 void devriesmat::print_othervalue_name(FILE *out,long compother)
00881 {
00882 switch (compother){
00883 case 0:{
00884 fprintf (out,"Moisture (m3/m3) ");
00885 break;
00886 }
00887 case 1:{
00888 fprintf (out,"Temperature (K) ");
00889 break;
00890 }
00891 case 2:{
00892 fprintf (out,"Relative humidity (--) ");
00893 break;
00894 }
00895 default:{
00896 fprintf (stderr,"\n\n unknown type of component is required in function (%s, line %d).\n",__FILE__,__LINE__);
00897 }
00898 }
00899 }
00900
00901 void devriesmat::CorD(int cislochar, int &kvyhl,double in,int rhw, double x, double & y, double & z, double &z2)
00902 {
00903
00904 switch (MatChar[cislochar])
00905 {
00906 case 0:
00907 break;
00908 case 1: y = MatConst[cislochar];
00909 break;
00910 case 2: {
00911 int i = 1;
00912
00913 int psl, prad;
00914 prad = int(MatData[cislochar][0][0]);
00915 psl = int(MatData[cislochar][1][0]);
00916 if(rhw == 1)
00917 {
00918 if (x <=0)
00919 {
00920 x = 0;
00921 y = MatData[cislochar][0][1];
00922 if (psl == 3)
00923 {
00924 z = MatData[cislochar][2][1];
00925 }
00926
00927 }
00928 else {
00929 if (x > in)
00930 {
00931 x = in;
00932 y = MatData[cislochar][0][prad];
00933 if (psl == 3)
00934 {
00935 z = MatData[cislochar][2][prad];
00936 }
00937 }
00938 else
00939 {
00940 while (x > MatData[cislochar][1][i]){
00941 i++;
00942 }
00943
00944 y = MatData[cislochar][0][i-1] + ((x - MatData[cislochar][1][i-1])*(MatData[cislochar][0][i] - MatData[cislochar][0][i-1]))/(MatData[cislochar][1][i] - MatData[cislochar][1][i-1]);
00945 if (psl == 3)
00946 {
00947 z = MatData[cislochar][2][i-1] + ((x - MatData[cislochar][1][i-1])*(MatData[cislochar][2][i] - MatData[cislochar][2][i-1]))/(MatData[cislochar][1][i] - MatData[cislochar][1][i-1]);
00948 }
00949 }
00950 }
00951 }
00952 else {
00953 if (x <=0)
00954 {
00955 x = 0;
00956 y = MatData[cislochar][1][1];
00957 if (psl == 3)
00958 {
00959 z = MatData[cislochar][2][1];
00960 }
00961
00962 }
00963 else {
00964 if (x > in)
00965 {
00966 x = in;
00967 y = MatData[cislochar][1][prad];
00968 if (psl == 3)
00969 {
00970 z = MatData[cislochar][2][prad];
00971 }
00972 }
00973 else
00974 {
00975 while (x > MatData[cislochar][0][i]){
00976 i++;
00977 }
00978
00979 y = MatData[cislochar][1][i-1] + ((x - MatData[cislochar][0][i-1])*(MatData[cislochar][1][i] - MatData[cislochar][1][i-1]))/(MatData[cislochar][0][i] - MatData[cislochar][0][i-1]);
00980 if (psl == 3)
00981 {
00982 z = MatData[cislochar][2][i-1] + ((x - MatData[cislochar][0][i-1])*(MatData[cislochar][2][i] - MatData[cislochar][2][i-1]))/(MatData[cislochar][0][i] - MatData[cislochar][0][i-1]);
00983 }
00984 }
00985 }
00986 }
00987 kvyhl = 2;
00988 }
00989 break;
00990 case 30: y = MatFunce[cislochar][0];
00991 z = MatFunce[cislochar][1];
00992 kvyhl = 30;
00993 break;
00994 case 31: y = MatFunce[cislochar][0];
00995 z = MatFunce[cislochar][1];
00996 if (cislochar == 6)
00997 {
00998 z2 = MatFunce[cislochar][2];
00999 }
01000 kvyhl = 31;
01001 break;
01002 case 32: y = MatFunce[cislochar][0];
01003 z = MatFunce[cislochar][1];
01004 kvyhl = 32;
01005 break;
01006 case 33:
01007 break;
01008 }
01009 }
01010
01011 double devriesmat::saturation_water_vapour_pressure(double x1, double x2, double x3)
01012 {
01013 return (exp(23.5771 - 4042.9/(x2 - 37.58)));
01014 }
01015
01016 double devriesmat::derivation_saturation_water_vapour_pressure_temperature(double x1, double x2, double x3)
01017 {
01018 return ((4042.9 * exp(23.5771-4042.9/(x2 - 37.58)))/((x2-37.58)*(x2-37.58)));
01019 }
01020
01021 double devriesmat::partial_water_vapour_pressure_function(double x1, double x2, double x3)
01022 {
01023 double relhum, pvs;
01024
01025 pvs = saturation_water_vapour_pressure(x1,x2,x3);
01026 sorption_izoterms_giva_data(0,x1,x2,x3,relhum,a1);
01027
01028 return(pvs*relhum);
01029
01030 }
01031
01032
01033 void devriesmat::sorption_izotherm_derivation(double x1, double x2, double x3, double & derfi)
01034 {
01035 int i = 0;
01036 if ( x1 < 0)
01037 {
01038 x1 = 0;
01039 }
01040
01041 while (x1 > MatData[6][1][i]){
01042 i++;
01043 }
01044 derfi = (MatData [6][0][i] - MatData [6][0][i-1])/(MatData [6][1][i] - MatData [6][1][i-1]);
01045
01046 }
01047
01048 double devriesmat::derivation_specific_internal_energy_of_water_vapour(double x1, double x2, double x3)
01049 {
01050 return(cvapd );
01051 }
01052
01053 double devriesmat::get_rel_hum(double w)
01054 {
01055 double relhum;
01056 sorption_izoterms_giva_data(0,w,0,0,relhum,a1);
01057 return (relhum);
01058 }
01059
01060 double devriesmat::sortpion_isotherm_root_shifted(double x1, double w_hyg, double rh_hyg)
01061 {
01062 return(w_hyg/(1-sqrt(1-rh_hyg)));
01063 }
01064
01065 void devriesmat::give_data_si_root_dfidw(double x1, double x2, double x3, double rh_hyg, double w_sat, double w_hyg, double shift_w, double & dfdw)
01066 {
01067
01068 if ( x1< 0) x1 = 0;
01069 if (x1>w_sat) x1 = w_sat;
01070
01071 if(x1<w_hyg)
01072 {
01073 dfdw = 2*(shift_w-x1)/(shift_w*shift_w);
01074 }
01075 else
01076 {
01077 dfdw = (1-rh_hyg)/(w_sat-w_hyg);
01078 }
01079 }
01080
01081 double devriesmat::get_moisture(double rh)
01082 {
01083 int s;
01084 double moist,hmrh, mhmc, smc, u, a, n;
01085
01086 if (rh < 0) rh = 0;
01087
01088 CorD(6,s,0,0,0,a1,a2,a3);
01089
01090 switch (s)
01091 {
01092 case 2: CorD(6,s,1,0,rh,moist,a1,a3);
01093 break;
01094 case 30: CorD(6,s,1,0,rh,a1,hmrh,a3);
01095 CorD(7,s,0,0,rh,smc,a2,a3);
01096 CorD(6,s,0,0,rh,mhmc,hmrh,a3);
01097
01098
01099 break;
01100 case 31: CorD(6,s,0,0,0,u,a,n);
01101
01102 break ;
01103 }
01104 return(moist);
01105
01106 }
01107 void devriesmat::sorption_izoterms_giva_data(int kod,double x1, double x2, double x3, double & fi, double & dfdw)
01108 {
01109
01110
01111
01112 int s;
01113 double w_hyg,rh_hyg, smc, hvezda;
01114
01115 if (x1 < 0) x1 = 0;
01116
01117 CorD(6,s,0,0,0,a1,a2,a3);
01118
01119 switch (s)
01120 {
01121 case 2: CorD(6,s,1,1,x1,fi,a2,a3);
01122 sorption_izotherm_derivation(x1,x2,x3,dfdw);
01123 break;
01124 case 30:
01125 CorD(6,s,0,0,x1,w_hyg,rh_hyg,a3);
01126 hvezda =sortpion_isotherm_root_shifted(x1,w_hyg ,rh_hyg);
01127 CorD(7,s,0,0,x1,smc,a2,a3);
01128 switch (kod)
01129 {
01130 case 0:
01131 give_data_si_root_fi(x1,x2,x3,w_hyg, smc, rh_hyg,hvezda,fi);
01132 break;
01133 case 1: give_data_si_root_dfidw(x1,x2,x3,rh_hyg,smc,w_hyg,hvezda,dfdw);
01134 break;
01135 }
01136 break;
01137
01138
01139
01140
01141
01142
01143
01144
01145
01146 }
01147
01148
01149
01150 }
01151 void devriesmat::give_data_si_root_fi(double x1, double x2, double x3, double w_hyg, double w_sat, double rh_hyg, double shift_w, double & relh)
01152 {
01153
01154 if (x1 < w_hyg)
01155 {
01156 relh = 1 - (1-x1/shift_w)*(1-x1/shift_w);
01157 }
01158 else
01159 {
01160 if (x1 > w_sat)
01161 {
01162 relh = 1;
01163 }
01164 else
01165 {
01166 relh = rh_hyg+(1-rh_hyg)*(x1-w_hyg)/(w_sat-w_hyg);
01167 }
01168 }
01169
01170
01171
01172 }
01173
01174 double devriesmat::pressure_head(double x1, double x2, double x3)
01175 {
01176 double relhum;
01177 sorption_izoterms_giva_data(0,x1,0,0,relhum,a1);
01178 if (relhum<=0.1)
01179 {
01180 return(1e6);
01181 }
01182 else {
01183 return (Rvapd*x2*log(relhum)/9.80665);
01184 }
01185
01186
01187
01188
01189 }
01190
01191 double devriesmat::derivation_pressure_head(double x1, double x2, double x3)
01192 {
01193 double relhum, dfdw;
01194 sorption_izoterms_giva_data(0,x1,0,0,relhum,a1);
01195 sorption_izoterms_giva_data(1,x1,x2,x3,a1,dfdw);
01196 if (relhum<=0.1)
01197 {
01198 return(Rvapd*x2/(9.80665*0.1)*1/dfdw);
01199 }
01200 else {
01201 return (Rvapd*x2/(9.80665*relhum)*1/dfdw);
01202
01203 }
01204
01205 }
01206
01207
01208 double devriesmat::relative_volume_ration_a(double x1, double x2, double x3)
01209 {
01210 double pi;
01211 CorD(3,kd,0,0,x1,pi,a2,a3);
01212 return (pi-x1);
01213 }
01214
01215
01216 double devriesmat::surface_tension(double x1, double x2, double x3)
01217 {
01218 double teplota;
01219 teplota = x2-273.15;
01220
01221 return ((-0.1647*teplota + 75.882)/1000);
01222 }
01223
01224 double devriesmat::derivation_surface_tension_on_temperature(double x1, double x2, double x3)
01225 {
01226 return (-0.0001647);
01227 }
01228
01229
01230 double devriesmat::coeff_of_water_diffusion_by_moisture_grad(double x1, double x2, double x3)
01231 {
01232 double dhdw, kl;
01233
01234 CorD(8,kd,1,0,x1,kl,a2,a3);
01235 dhdw = derivation_pressure_head(x1,x2,x3);
01236 return(kl*dhdw);
01237
01238 }
01239
01240 double devriesmat::coeff_of_water_diffusion_by_temperature_grad(double x1, double x2, double x3)
01241 {
01242 double kl, sigm, h,dsigmdt;
01243 CorD(8,kd,1,0,x1,kl,a2,a3);
01244 h = pressure_head(x1,x2,x3);
01245 sigm = surface_tension(x1,x2,x3);
01246 dsigmdt = derivation_surface_tension_on_temperature(x1,x2,x3);
01247 return (kl*h*1/sigm*dsigmdt);
01248
01249
01250 }
01251
01252 double devriesmat::coeff_dwv(double x1, double x2, double x3)
01253 {
01254 double rol, mi, tor, a,d,p,pv,pvs,dfdw;
01255 rol = 1000;
01256 CorD(4,kd,0,0,x1,mi,a2,a3);
01257 tor = 1/mi;
01258 a =relative_volume_ration_a(x1,x2,x3);
01259 d = diffusion_coefficient_of_water_vapor_in_air(x1,x2,x3);
01260 p = 101325;
01261 pv = partial_water_vapour_pressure_function(x1,x2,x3);
01262 pvs = saturation_water_vapour_pressure(x1,x2,x3);
01263
01264 sorption_izoterms_giva_data(1,x1,x2,x3,a1,dfdw);
01265
01266 return (1/rol*tor*a*d*p/(p-pv)*1/(Rvapd*x2)*pvs *dfdw);
01267
01268
01269 }
01270
01271 double devriesmat::coeff_dtv(double x1, double x2, double x3)
01272 {
01273 double rold, mi, tor, a,d,p,pv,pvs,dfdw, relhum, dpvsdt, dfdt;
01274 rold = 1000;
01275 CorD(4,kd,0,0,x1,mi,a2,a3);
01276 tor = 1/mi;
01277 a =relative_volume_ration_a(x1,x2,x3);
01278 d = diffusion_coefficient_of_water_vapor_in_air(x1,x2,x3);
01279 p = 101325;
01280 pv = partial_water_vapour_pressure_function(x1,x2,x3);
01281 pvs = saturation_water_vapour_pressure(x1,x2,x3);
01282 dpvsdt = derivation_saturation_water_vapour_pressure_temperature(x1,x2,x3);
01283 sorption_izoterms_giva_data(0,x1,0,0,relhum,a1);
01284 sorption_izoterms_giva_data(1,x1,x2,x3,a1,dfdw);
01285 dfdt = 0.0;
01286
01287 return(1/rold*tor*a*d*p/(p-pv)*1/(Rvapd*x2)*(-relhum*pvs/x2+relhum*dpvsdt+pvs *dfdt));
01288
01289 }
01290
01291 double devriesmat::specific_heat_capacities_star(double x1, double x2, double x3)
01292 {
01293 double cs, cl, cv, ros, rov, a;
01294 CorD(9,kd,0,0,x1,cs,a2,a3);
01295 CorD(2,kd,0,0,x1,ros,a2,a3);
01296 a =relative_volume_ration_a(x1,x2,x3);
01297 cl = cwatd;
01298 cv = cvapd;
01299 rold = 1000;
01300 rov = partial_density_of_water_vapor(x1,x2,x3);
01301
01302 return (cs +rold/ros*cl*x1+rov/ros*a*cv);
01303
01304
01305 }
01306
01307 double devriesmat::partial_density_of_water_vapor(double x1, double x2, double x3)
01308 {
01309 double pv;
01310 pv = partial_water_vapour_pressure_function(x1,x2,x3);
01311 return (pv/(Rvapd*x2));
01312 }
01313
01314 double devriesmat::latent_heat_of_evaporation_of_water(double x1, double x2, double x3)
01315 {
01316 return (2.5008e6)*pow((273.15/x2),(0.167+x2*3.67e-4));
01317 }
01318
01319
01320 double devriesmat::derivation_drov_dw(double x1, double x2, double x3)
01321 {
01322 double pvs,dfdw;
01323 pvs = saturation_water_vapour_pressure(x1,x2,x3);
01324 sorption_izoterms_giva_data(1,x1,x2,x3,a1,dfdw);
01325
01326 return (1/(Rvapd*x2)*pvs*dfdw);
01327
01328 }
01329
01330 double devriesmat::diffusion_coefficient_of_water_vapor_in_air(double x1, double x2, double x3)
01331 {
01332 double mi;
01333 CorD(4,kd,0,0,x1,mi,a2,a3);
01334 return(2.3e-5 * pow(x2/273,1.81));
01335 }