00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028 #include "kunzel.h"
00029 #include "globalt.h"
00030
00031
00032 kunmat::kunmat()
00033 {
00034 rho_m = 0.0;
00035 rho_w = 1000.0;
00036 }
00037
00038 kunmat::~kunmat()
00039 {
00040
00041
00042 }
00043
00044
00045
00046
00047
00048
00049
00050
00051
00052
00053
00054
00055
00056 void kunmat::matcond (matrix &d,long ri,long ci,long ipp)
00057 {
00058 long n;
00059 n = d.n;
00060
00061 switch (n){
00062 case 1:{
00063 matcond1d (d,ri,ci,ipp);
00064 break;
00065 }
00066 case 2:{
00067 matcond2d (d,ri,ci,ipp);
00068 break;
00069 }
00070 case 3:{
00071 matcond3d (d,ri,ci,ipp);
00072 break;
00073 }
00074 default:{
00075 print_err("unknown number of components of conductivity tensor is required",__FILE__,__LINE__,__func__);
00076 }
00077 }
00078 }
00079
00080
00081
00082
00083
00084
00085
00086
00087
00088
00089 void kunmat::matcond1d (matrix &d,long ri,long ci,long ipp)
00090 {
00091 double k;
00092 double phi,t;
00093 k = 0.0;
00094
00095 phi = Tm->ip[ipp].av[0];
00096 t = Tm->ip[ipp].av[1];
00097
00098 if((ri == 0) && (ci == 0))
00099 k = tokJ1(phi,t,ipp);
00100 if((ri == 0) && (ci == 1))
00101 k = tokJ2(phi,t,ipp);
00102 if((ri == 1) && (ci == 0))
00103 k = tokJ3(phi,t,ipp);
00104 if((ri == 1) && (ci == 1))
00105 k = tokJ4(phi,t,ipp);
00106
00107 d[0][0] = k;
00108 }
00109
00110
00111
00112
00113
00114
00115
00116
00117
00118 void kunmat::matcond2d (matrix &d,long ri,long ci,long ipp)
00119 {
00120 double k;
00121 double phi,t;
00122 k = 0.0;
00123
00124 phi = Tm->ip[ipp].av[0];
00125 t = Tm->ip[ipp].av[1];
00126
00127 if((ri == 0) && (ci == 0))
00128 k = tokJ1(phi,t,ipp);
00129 if((ri == 0) && (ci == 1))
00130 k = tokJ2(phi,t,ipp);
00131 if((ri == 1) && (ci == 0))
00132 k = tokJ3(phi,t,ipp);
00133 if((ri == 1) && (ci == 1))
00134 k = tokJ4(phi,t,ipp);
00135
00136 fillm(0.0,d);
00137
00138
00139
00140 d[0][0] = k; d[0][1] = 0.0;
00141 d[1][0] = 0.0; d[1][1] = k;
00142 }
00143
00144
00145
00146
00147
00148
00149
00150
00151
00152
00153 void kunmat::matcond3d (matrix &d,long ri,long ci,long ipp)
00154 {
00155 double k;
00156 double phi,t;
00157 k = 0.0;
00158
00159 phi = Tm->ip[ipp].av[0];
00160 t = Tm->ip[ipp].av[1];
00161
00162 if((ri == 0) && (ci == 0))
00163 k = tokJ1(phi,t,ipp);
00164 if((ri == 0) && (ci == 1))
00165 k = tokJ2(phi,t,ipp);
00166 if((ri == 1) && (ci == 0))
00167 k = tokJ3(phi,t,ipp);
00168 if((ri == 1) && (ci == 1))
00169 k = tokJ4(phi,t,ipp);
00170
00171 fillm(0.0,d);
00172
00173 d[0][0]=k; d[0][1]=0.0; d[0][2]=0.0;
00174 d[1][0]=0.0; d[1][1]=k; d[1][2]=0.0;
00175 d[2][0]=0.0; d[2][1]=0.0; d[2][2]=k;
00176 }
00177
00178
00179
00180
00181
00182
00183
00184
00185
00186
00187 void kunmat::matcap (double &c,long ri,long ci,long ipp)
00188 {
00189 double phi,t;
00190 c = 0.0;
00191
00192 phi = Tm->ip[ipp].av[0];
00193 t = Tm->ip[ipp].av[1];
00194
00195 if((ri == 0) && (ci == 0))
00196
00197 c = Tm->ip[ipp].eqother[1];
00198 if((ri == 0) && (ci == 1))
00199 c = 0.0;
00200 if((ri == 1) && (ci == 0))
00201 c = 0.0;
00202 if((ri == 1) && (ci == 1))
00203 c = DerivaceHustotyEntalpiePodleTeploty(phi,t,ipp);
00204 }
00205
00206
00207 void kunmat::values_correction (vector &nv, long ipp)
00208 {
00209
00210 relhum_check(nv[0],nv[1],ipp);
00211 }
00212
00213
00214
00215
00216
00217
00218
00219
00220
00221
00222
00223 void kunmat::relhum_check(double &rh,double &t,long ipp)
00224 {
00225 if (rh >= 1.0)
00226 rh = 1.0;
00227 if (rh <= 0.0)
00228 rh = 0.0;
00229
00230
00231 if (t >= 350.0)
00232 t = 350.0;
00233 if (t <= 240.0)
00234 t = 240.0;
00235
00236
00237
00238
00239 }
00240
00241
00242
00243
00244
00245
00246
00247
00248 double kunmat::DerivaceTlakuNasycenychParNaTeplote(double rh, double tk)
00249 {
00250
00251 return ((4042.9 * exp(23.5771-4042.9/(tk - 37.58)))/((tk-37.58)*(tk-37.58)));
00252 }
00253
00254
00255
00256
00257
00258
00259
00260
00261
00262
00263
00264
00265 double kunmat::PermeabilitaVodniPary(double rh, double tk,long ipp)
00266 {
00267 double P, Rv, Pa;
00268 double da, dp, mi, smc;
00269
00270 if (rh >1) rh = 1.0;
00271 if (rh < 0) rh = 0.0;
00272
00273
00274
00275 moist = Tm->ip[ipp].eqother[0];
00276 smc = Tm->ip[ipp].eqother[2];
00277
00278 CorD(4,kd,moist,mi,a2,a3);
00279
00280 Rv = 461.5;
00281 Pa = 101325;
00282 P = 101325;
00283
00284 da = (2.306e-5 * Pa)/(Rv * tk * P)*pow((tk/273.15),1.81);
00285
00286 if (mi==0.0)
00287 dp=0.0;
00288 else
00289 dp = da/mi;
00290
00291 return (dp);
00292 }
00293
00294
00295
00296
00297
00298
00299
00300
00301
00302
00303 double kunmat::tokJ2(double rh, double tk, long ipp)
00304 {
00305 double dp, dpvs;
00306
00307 if (rh>1.0) rh = 1.0;
00308 if (rh<0.0) rh = 0.0;
00309
00310 dpvs = DerivaceTlakuNasycenychParNaTeplote(rh,tk);
00311 dp = PermeabilitaVodniPary(rh, tk,ipp);
00312
00313 J2 = dp * rh * dpvs;
00314
00315
00316
00317 return (J2);
00318 }
00319
00320
00321
00322
00323
00324
00325
00326
00327
00328
00329
00330 double kunmat::tokJ4(double rh, double tk, long ipp)
00331 {
00332 double lv1, lambdaakt, dpvs, dp, moist,smc;
00333
00334 if (rh >1) rh = 1;
00335 if (rh < 0) rh = 0;
00336
00337
00338
00339 moist = Tm->ip[ipp].eqother[0];
00340 smc = Tm->ip[ipp].eqother[2];
00341
00342 CorD(10,kd,moist,lambdaakt,a2,a3);
00343
00344 dpvs = DerivaceTlakuNasycenychParNaTeplote(rh, tk);
00345 dp = PermeabilitaVodniPary(rh, tk,ipp);
00346
00347 lv1 = LatentHeatofEvaporationOfWater (tk);
00348
00349 J4 = lambdaakt + lv1 * dp * rh * dpvs;
00350
00351 return (J4);
00352 }
00353
00354
00355
00356
00357
00358
00359
00360
00361
00362
00363
00364
00365 double kunmat::tokJ3(double rh, double tk, long ipp)
00366 {
00367 double lv1;
00368 double dp, ps;
00369
00370 if (rh >1) rh = 1;
00371 if (rh < 0) rh = 0;
00372
00373 dp = PermeabilitaVodniPary(rh, tk,ipp);
00374 ps = TlakNasycenychVodnichParNaTeplote (rh, tk);
00375 lv1 = LatentHeatofEvaporationOfWater (tk);
00376
00377 J3 = lv1 * dp * ps;
00378
00379 return (J3);
00380 }
00381
00382
00383
00384
00385
00386
00387
00388 double kunmat::LatentHeatofEvaporationOfWater(double tk)
00389 {
00390 return (2.5008e6)*pow((273.15/tk),(0.167+tk*3.67e-4));
00391 }
00392
00393
00394
00395
00396
00397
00398
00399
00400
00401
00402
00403
00404
00405 double kunmat::TlakNasycenychVodnichParNaTeplote(double rh, double tk)
00406 {
00407 return (exp(23.5771 - 4042.9/(tk - 37.58)));
00408 }
00409
00410
00411
00412
00413
00414
00415
00416
00417 double kunmat::DerivaceHustotyEntalpiePodleTeploty(double rh, double tk, long ipp)
00418 {
00419 double rom, c, moist;
00420 if (rh >1) rh = 1;
00421 if (rh < 0) rh = 0;
00422
00423
00424 moist = Tm->ip[ipp].eqother[0];
00425
00426 CorD(2,kd,0,rom,a2,a3);
00427 CorD(9,kd,moist,c,a2,a3);
00428 rho_m = rom;
00429
00430 return (rom * c);
00431 }
00432
00433
00434
00435
00436
00437
00438
00439 double kunmat::HygroscopicMoisture(double rh, double tk)
00440 {
00441
00442 double mhmc, mhrh,a3;
00443 CorD(6,kd,rh,mhmc,mhrh,a3);
00444
00445 return (mhmc/(1-sqrt(1-mhrh)));
00446
00447 }
00448
00449
00450
00451
00452
00453
00454
00455
00456 double kunmat::DerivativeOfTheSorptionIsotherm(double rh, double tk)
00457 {
00458
00459 double mhv;
00460 if (rh >1) rh = 1;
00461 if (rh < 0) rh = 0;
00462
00463 mhv = HygroscopicMoisture (rh,tk);
00464
00465 return (1000 * mhv/(2*sqrt(1-rh)));
00466
00467 }
00468
00469
00470
00471
00472
00473
00474
00475 double kunmat::DerivativeOfTheRetentionCurve(double rh, double tk, long ipp)
00476 {
00477
00478 double smc, mhmc, mhrh;
00479
00480 smc = Tm->ip[ipp].eqother[2];
00481 CorD(6,kd,rh,mhmc,mhrh,a3);
00482
00483
00484 return (1000 * (smc - mhmc)/(1 - mhrh));
00485 }
00486
00487
00488
00489
00490
00491
00492
00493 double kunmat::DerivativeOfTheMoistureRetentionCharacteristik(double rh,double tk, long ipp)
00494 {
00495
00496 double drfdrh;
00497
00498 sorption_izoterms_giva_data(1,rh,tk,a1,drfdrh, ipp);
00499
00500 return (drfdrh);
00501
00502 }
00503
00504
00505
00506
00507
00508
00509
00510 double kunmat::tokJ1(double rh, double tk, long ipp)
00511 {
00512 double kapa1, mi, dp, dmretc, ps, drh, moist, smc;
00513
00514
00515
00516 if (rh < 0) rh = 0;
00517 if (rh >1) rh = 1;
00518
00519
00520 moist = Tm->ip[ipp].eqother[0];
00521
00522 smc = Tm->ip[ipp].eqother[2];
00523 CorD(4,kd,moist,mi,a2,a3);
00524 kapa1 = Tm->ip[ipp].eqother[3];
00525
00526
00527
00528
00529
00530
00531
00532
00533
00534
00535
00536
00537
00538
00539
00540
00541
00542
00543
00544
00545
00546
00547
00548
00549 dp = PermeabilitaVodniPary(rh, tk,ipp);
00550 ps = TlakNasycenychVodnichParNaTeplote (rh, tk);
00551
00552
00553
00554 dmretc = Tm->ip[ipp].eqother[1];
00555
00556 drh = kapa1 * dmretc;
00557
00558 J1 = drh + dp * ps;
00559
00560
00561
00562 return (J1);
00563 }
00564
00565
00566
00567
00568
00569
00570 void kunmat::read(XFILE *in)
00571 {
00572 int i,nr,nc;
00573
00574 for (i = 2; i<20; i++){
00575 xfscanf(in, "%d",&MatChar[i]);
00576 }
00577
00578 for (i = 2;i <= 19; i++)
00579 {
00580
00581 switch (MatChar[i]){
00582 case 0:{
00583
00584 break;
00585 }
00586 case 1:{
00587
00588
00589 data[i]=new gfunct ();
00590
00591
00592 data[i]->tfunc=stat;
00593
00594 xfscanf(in, "%lf",&data[i]->f);
00595
00596 break;
00597 }
00598 case 2:{
00599
00600
00601
00602
00603 xfscanf(in, "%ld %ld",&nr,&nc);
00604
00605
00606 data[i]=new gfunct (tab,nr);
00607
00608 data[i]->tabf[0].itype = piecewiselin;
00609
00610 data[i]->tabf[0].readval (in);
00611
00612
00613 xfscanf(in, "%d",&nr);
00614
00615 break;
00616 }
00617 case 3:{
00618
00619 break;
00620 }
00621 case 30:{
00622 if( MatChar[i-1] ==1)
00623 {
00624 xfscanf(in, "%lf %lf",&MatFunce[i][0],&MatFunce[i][1]);
00625 }
00626 else
00627 {
00628 xfscanf(in, "%lf %lf",&MatFunce[i][0],&MatFunce[i][1]);
00629 }
00630 break;
00631 }
00632 case 31:{
00633 if( MatChar[i-1] ==1)
00634 {
00635 if (i == 6)
00636 {
00637 xfscanf(in, "%lf %lf %lf",&MatFunce[i][0],&MatFunce[i][1], &MatFunce[i][2]);
00638 }
00639 else
00640 {
00641 xfscanf(in, "%lf %lf",&MatFunce[i][0],&MatFunce[i][1]);
00642 }
00643 }
00644 else
00645 {
00646 if (i == 6)
00647 {
00648 xfscanf(in, "%lf %lf %lf",&MatFunce[i][0],&MatFunce[i][1], &MatFunce[i][2]);
00649 }
00650 else
00651 {
00652 xfscanf(in, "%lf %lf",&MatFunce[i][0],&MatFunce[i][1]);
00653 }
00654 }
00655
00656 break;
00657 }
00658 case 32:{
00659 if( MatChar[i-1] ==1)
00660 {
00661 xfscanf(in, "%lf %lf",&MatFunce[i][0],&MatFunce[i][1]);
00662 }
00663 else
00664 {
00665 xfscanf(in, "%lf %lf",&MatFunce[i][0],&MatFunce[i][1]);
00666 }
00667 break;
00668 }
00669 default:{
00670 print_err("unknown definition of Material KUNZEL is required",__FILE__,__LINE__,__func__);
00671 }
00672 }
00673 }
00674
00675 madripom = 0;
00676 CorD(2,kd,0,rho_m,a2,a3);
00677 }
00678
00679
00680 void kunmat::print(FILE *out)
00681 {
00682
00683 int i, j;
00684
00685 fprintf (out,"\n");
00686 for (i = 2; i<19; i++)
00687 {
00688 fprintf(out, " %d ",MatChar[i]);
00689 }
00690 fprintf(out, " %d ",MatChar[19]);
00691 for (i = 2;i <= 19; i++)
00692 {
00693
00694 switch (MatChar[i]){
00695 case 0:{
00696
00697 break;
00698 }
00699 case 1:{
00700
00701 fprintf(out, "\n %e ",data[i]->f);
00702 break;
00703 }
00704 case 2:{
00705
00706
00707 fprintf(out, "\n %d %d ",data[i]->tabf[0].asize,2);
00708
00709 for (j=0;j<data[i]->tabf[0].asize;j++){
00710 fprintf(out, "\n %e %e ",data[i]->tabf[0].x[j],data[i]->tabf[0].y[j]);
00711 }
00712
00713 fprintf(out, "\n %d ",data[i]->tabf[0].asize);
00714 break;
00715 }
00716 case 3:{
00717
00718 break;
00719 }
00720 case 30:{
00721 if( MatChar[i-1] ==1)
00722 {
00723 fprintf(out, "\n %e %e ",MatFunce[i][0],MatFunce[i][1]);
00724 }
00725 else
00726 {
00727 fprintf(out, "\n %e %e ",MatFunce[i][0],MatFunce[i][1]);
00728 }
00729 break;
00730 }
00731 case 31:{
00732 if( MatChar[i-1] ==1)
00733 {
00734 if (i == 6)
00735 {
00736 fprintf(out, "\n %e %e %e ",MatFunce[i][0],MatFunce[i][1], MatFunce[i][2]);
00737 }
00738 else
00739 {
00740 fprintf(out, "\n %e %e ",MatFunce[i][0],MatFunce[i][1]);
00741 }
00742 }
00743 else
00744 {
00745 if (i == 6)
00746 {
00747 fprintf(out, "\n %e %e %e",MatFunce[i][0],MatFunce[i][1], MatFunce[i][2]);
00748 }
00749 else
00750 {
00751 fprintf(out, "\n %e %e ",MatFunce[i][0],MatFunce[i][1]);
00752 }
00753 }
00754
00755 break;
00756 }
00757 case 32:{
00758 if( MatChar[i-1] ==1)
00759 {
00760 fprintf(out, "\n %e %e ",MatFunce[i][0],MatFunce[i][1]);
00761 }
00762 else
00763 {
00764 fprintf(out, "\n %e %e ",MatFunce[i][0],MatFunce[i][1]);
00765 }
00766 break;
00767 }
00768 default:{
00769 print_err("unknown definition of Material KUNZEL is required",__FILE__,__LINE__,__func__);
00770 fprintf (out,"\n");
00771 }
00772 }
00773 }
00774 fprintf (out,"\n");
00775 }
00776
00777
00778
00779
00780
00781
00782
00783
00784
00785 double kunmat::sorptionizothermDerivation(double rh, double tk)
00786 {
00787 double a, b, n;
00788
00789 CorD(6,kd,rh,a,b,n);
00790
00791
00792 if (rh >1) rh = 1;
00793 if (rh < 0) rh = 0;
00794
00795
00796 return (1000 * (a/(b*rh*n))*pow((1-(log(rh))/b),(-1-(1/n))));
00797
00798 }
00799
00800
00801
00802
00803
00804
00805
00806
00807
00808
00809
00810
00811
00812
00813
00814
00815
00816 double kunmat::transmission_transcoeff(double trc,long ri,long ci,long nn,long bc,long ipp)
00817 {
00818 long k;
00819 double new_trc,h,t;
00820 new_trc = 0.0;
00821
00822 k=Gtt->give_dof(nn,0);
00823 if (k>0) {h = Lsrst->lhs[k-1]+Lsrst->lhsi[k-1];}
00824 if (k==0) {h = 0.0;}
00825 if (k<0) {h = Tb->lc[0].pv[0-k-1].getval();}
00826 k=Gtt->give_dof(nn,1);
00827 if (k>0) {t = Lsrst->lhs[k-1]+Lsrst->lhsi[k-1];}
00828 if (k==0) {t = 0.0;}
00829 if (k<0) {t = Tb->lc[0].pv[0-k-1].getval();}
00830
00831 if((ri == 0) && (ci == 0))
00832 new_trc = get_transmission_transcoeff_hh(h,t,bc,ipp);
00833 if((ri == 0) && (ci == 1))
00834 new_trc = 0.0;
00835
00836 if((ri == 1) && (ci == 0))
00837 new_trc = 0.0;
00838 if((ri == 1) && (ci == 1))
00839 new_trc = get_transmission_transcoeff_tt(h,t,bc,ipp);
00840
00841 new_trc = new_trc*trc;
00842
00843 return (new_trc);
00844 }
00845
00846
00847
00848
00849
00850
00851
00852
00853
00854
00855
00856
00857
00858
00859 double kunmat::transmission_nodval(double nodval,double trc2,long ri,long ci,long nn,long bc,long ipp)
00860 {
00861 long k;
00862 double new_nodval,h,t;
00863 new_nodval = 0.0;
00864
00865 k=Gtt->give_dof(nn,0);
00866 if (k>0) {h = Lsrst->lhs[k-1]+Lsrst->lhsi[k-1];}
00867 if (k==0) {h = 0.0;}
00868 if (k<0) {h = Tb->lc[0].pv[0-k-1].getval();}
00869 k=Gtt->give_dof(nn,1);
00870 if (k>0) {t = Lsrst->lhs[k-1]+Lsrst->lhsi[k-1];}
00871 if (k==0) {t = 0.0;}
00872 if (k<0) {t = Tb->lc[0].pv[0-k-1].getval();}
00873
00874 if((ri == 0) && (ci == 0))
00875 new_nodval = get_transmission_nodval_hh(nodval,h,t,bc,ipp);
00876 if((ri == 0) && (ci == 1))
00877 new_nodval = 0.0;
00878
00879 if((ri == 1) && (ci == 0))
00880 new_nodval = 0.0;
00881 if((ri == 1) && (ci == 1))
00882 new_nodval = get_transmission_nodval_tt(nodval,h,t,bc,ipp);
00883
00884 return (new_nodval);
00885 }
00886
00887
00888
00889
00890
00891
00892
00893
00894
00895
00896
00897
00898
00899
00900
00901 double kunmat::transmission_flux(double nodval,double trc2,long ri,long ci,long nn,long bc,long ipp)
00902 {
00903 long k;
00904 double flux,h,t;
00905 flux = 0.0;
00906
00907 k=Gtt->give_dof(nn,0);
00908 if (k>0) {h = Lsrst->lhs[k-1]+Lsrst->lhsi[k-1];}
00909 if (k==0) {h = 0.0;}
00910 if (k<0) {h = Tb->lc[0].pv[0-k-1].getval();}
00911 k=Gtt->give_dof(nn,1);
00912 if (k>0) {t = Lsrst->lhs[k-1]+Lsrst->lhsi[k-1];}
00913 if (k==0) {t = 0.0;}
00914 if (k<0) {t = Tb->lc[0].pv[0-k-1].getval();}
00915
00916 if((ri == 0) && (ci == 0))
00917 flux = get_transmission_flux_hh(nodval,h,t,bc,ipp);
00918 if((ri == 0) && (ci == 1))
00919 flux = 0.0;
00920
00921 if((ri == 1) && (ci == 0))
00922 flux = 0.0;
00923 if((ri == 1) && (ci == 1))
00924 flux = get_transmission_flux_tt(nodval,h,t,bc,ipp);
00925
00926 return (flux);
00927 }
00928
00929
00930
00931
00932
00933
00934
00935
00936
00937 double kunmat::get_transmission_transcoeff_hh(double rh,double t,long bc,long ipp)
00938 {
00939 double trc,pgws;
00940
00941 switch (bc){
00942 case 30:{
00943 pgws = TlakNasycenychVodnichParNaTeplote(rh,t);
00944 trc = pgws;
00945 break;
00946 }
00947 case 31:{
00948 pgws = TlakNasycenychVodnichParNaTeplote(rh,t);
00949 trc = pgws;
00950 break;
00951 }
00952 case 32:{
00953 trc = 0.0;
00954 break;
00955 }
00956 default:{
00957 print_err("no real boundary condition is prescribed",__FILE__,__LINE__,__func__);
00958 exit(0);
00959 }
00960 }
00961 return(trc);
00962 }
00963
00964
00965
00966
00967
00968
00969
00970
00971
00972 double kunmat::get_transmission_nodval_hh(double bv,double rh,double t,long bc,long ipp)
00973 {
00974 double nodval,pgws, pgwspred;
00975
00976 switch (bc){
00977 case 30:{
00978
00979 pgws = TlakNasycenychVodnichParNaTeplote(rh,t);
00980 pgwspred = TlakNasycenychVodnichParNaTeplote(rh,t);
00981
00982 nodval = bv*pgwspred;
00983
00984 break;
00985 }
00986 case 31:{
00987 nodval = bv;
00988 break;
00989 }
00990 case 32:{
00991 pgws = TlakNasycenychVodnichParNaTeplote(rh,t);
00992
00993 nodval = pgws*rh;
00994 bv = pgws*bv;
00995 nodval = bv - nodval;
00996
00997 break;
00998 }
00999 default:{
01000 print_err("no real boundary condition is prescribed",__FILE__,__LINE__,__func__);
01001 exit(0);
01002 }
01003 }
01004 return(nodval);
01005 }
01006
01007
01008
01009
01010
01011
01012
01013
01014
01015
01016 double kunmat::get_transmission_flux_hh(double bv,double rh,double t,long bc,long ipp)
01017 {
01018 double flux,pgws;
01019
01020 switch (bc){
01021 case 30:{
01022 pgws = TlakNasycenychVodnichParNaTeplote(rh, t);
01023
01024 flux = pgws*rh;
01025 bv = pgws*bv;
01026 flux = bv - flux;
01027
01028 break;
01029 }
01030 case 31:{
01031 pgws = TlakNasycenychVodnichParNaTeplote(rh, t);
01032 flux = pgws*rh;
01033
01034 flux = bv - flux;
01035 break;
01036 }
01037 case 32:{
01038 pgws = TlakNasycenychVodnichParNaTeplote(rh, t);
01039
01040 flux = pgws*rh;
01041 bv = pgws*bv;
01042 flux = bv - flux;
01043
01044 break;
01045 }
01046 default:{
01047 print_err("no real boundary condition is prescribed",__FILE__,__LINE__,__func__);
01048 exit(0);
01049 }
01050 }
01051 return(flux);
01052 }
01053
01054
01055
01056
01057
01058
01059
01060
01061
01062 double kunmat::get_transmission_transcoeff_tt(double h,double t,long bc,long ipp)
01063 {
01064 double trc;
01065
01066 switch (bc){
01067 case 30:{
01068 trc = 1.0;
01069 break;
01070 }
01071 default:{
01072 print_err("no real boundary condition is prescribed",__FILE__,__LINE__,__func__);
01073 exit(0);
01074 }
01075 }
01076
01077 return(trc);
01078 }
01079
01080
01081
01082
01083
01084
01085
01086
01087
01088 double kunmat::get_transmission_nodval_tt(double bv,double h,double t,long bc,long ipp)
01089 {
01090 double nodval;
01091
01092 switch (bc){
01093 case 30:{
01094 nodval = bv;
01095 break;
01096 }
01097 default:{
01098 print_err("no real boundary condition is prescribed",__FILE__,__LINE__,__func__);
01099 exit(0);
01100 }
01101 }
01102 return(nodval);
01103 }
01104
01105
01106
01107
01108
01109
01110
01111
01112
01113
01114 double kunmat::get_transmission_flux_tt(double bv,double h,double t,long bc,long ipp)
01115 {
01116 double flux;
01117
01118 switch (bc){
01119 case 30:{
01120 flux = (bv - t);
01121 break;
01122 }
01123 default:{
01124 print_err("no real boundary condition is prescribed",__FILE__,__LINE__,__func__);
01125 exit(0);
01126 }
01127 }
01128 return(flux);
01129 }
01130
01131
01132
01133
01134
01135
01136
01137
01138
01139
01140
01141 double kunmat::get_othervalue(long compother,double rh,double t, long ipp)
01142 {
01143 double other;
01144
01145 switch (compother){
01146 case 0:{
01147 other = rh;
01148 break;
01149 }
01150 case 1:{
01151 other = t;
01152 break;
01153 }
01154 case 2:{
01155 double w, a3;
01156 sorption_izoterms_giva_data(0,rh,t,w,a3, ipp);
01157 other = w;
01158 break;
01159 }
01160 case 3:{
01161 other = rh * TlakNasycenychVodnichParNaTeplote (rh, t);
01162 break;
01163 }
01164 case 4:{
01165 double w, a3;
01166 sorption_izoterms_giva_data(0,rh,t,w,a3, ipp);
01167 other = w*rho_w/rho_m;
01168 break;
01169 }
01170 case 5:{
01171 other = 0.0;
01172 break;
01173 }
01174 case 6:{
01175 other = 0.0;
01176 break;
01177 }
01178 case 7:{
01179 other = 0.0;
01180 break;
01181 }
01182 case 8:{
01183 other = 0.0;
01184 break;
01185 }
01186 case 9:{
01187 other = 0.0;
01188 break;
01189 }
01190 case 10:{
01191 other = 0.0;
01192 break;
01193 }
01194 default:{
01195 print_err("unknown type of component is required",__FILE__,__LINE__,__func__);
01196 }
01197 }
01198 return (other);
01199
01200 }
01201
01202
01203
01204
01205
01206
01207 void kunmat::print_othervalue_name(FILE *out,long compother)
01208 {
01209 switch (compother){
01210 case 0:{
01211 fprintf (out,"Relative humidity () ");
01212 break;
01213 }
01214 case 1:{
01215 fprintf (out,"Temperature (K) ");
01216 break;
01217 }
01218 case 2:{
01219 fprintf (out,"Moisture content w (m3/m3) ");
01220 break;
01221 }
01222 case 3:{
01223 fprintf (out,"Water vapour pressure (Pa) ");
01224 break;
01225 }
01226 case 4:{
01227 fprintf (out,"Moisture content u (kg/kg) ");
01228 break;
01229 }
01230 default:{
01231 print_err("unknown type of component is required",__FILE__,__LINE__,__func__);
01232 }
01233 }
01234 }
01235
01236
01237 void kunmat::give_data(double rh,double Mhmc, double Smc, double Mhrh, double & moistakt)
01238 {
01239 double hvezdicka;
01240 hvezdicka = Mhmc/(1-sqrt(1-Mhrh));
01241
01242 if (rh < Mhrh)
01243 {
01244 moistakt = (1-sqrt(1-rh)) * hvezdicka ;
01245 }
01246 else {
01247 if (rh > 1)
01248 {
01249 moistakt = Smc;
01250 }
01251 else
01252 {
01253 moistakt = Mhmc + (rh - Mhrh )/(1-Mhrh )* (Smc - Mhmc );
01254 }
01255 }
01256
01257 }
01258
01259
01260 double kunmat::kapa_exp(double a, double b, double rh, double tk, long ipp)
01261 {
01262 double moist;
01263
01264
01265
01266
01267
01268 return (a * exp(b * rh));
01269
01270 }
01271
01272
01273
01274
01275
01276
01277
01278
01279
01280
01281
01282
01283
01284
01285
01286
01287
01288
01289
01290
01291
01292
01293
01294
01295
01296
01297
01298
01299
01300
01301
01302
01303 void kunmat::CorD(long charid,long &kvyhl, double x, double &y, double &z, double &z2)
01304 {
01305
01306 switch (MatChar[charid]){
01307 case 0:{
01308 kvyhl = 0;
01309 break;
01310 }
01311 case 1:{
01312
01313
01314 y = data[charid]->getval(0.0);
01315 kvyhl = 1;
01316 break;
01317 }
01318 case 2:{
01319
01320
01321 y = data[charid]->getval(x);
01322 kvyhl = 2;
01323
01324
01325 long nr=data[charid]->tabf[0].asize;
01326
01327
01328 z2 = data[charid]->tabf[0].y[nr-1];
01329
01330 break;
01331 }
01332 case 30:{
01333 y = MatFunce[charid][0];
01334 z = MatFunce[charid][1];
01335 kvyhl = 30;
01336 break;
01337 }
01338 case 31:{
01339 y = MatFunce[charid][0];
01340 z = MatFunce[charid][1];
01341 if (charid == 6)
01342 {
01343 z2 = MatFunce[charid][2];
01344 }
01345 kvyhl = 31;
01346 break;
01347 }
01348 case 32:{
01349 y = MatFunce[charid][0];
01350 z = MatFunce[charid][1];
01351 kvyhl = 32;
01352 break;
01353 }
01354 case 33:{
01355 kvyhl = 33;
01356
01357 break;
01358 }
01359 default:{
01360 print_err("unknown definition of material parameter is required",__FILE__,__LINE__,__func__);
01361 }
01362 }
01363 }
01364
01365 void kunmat::sorption_izoterms_giva_data(long kod,double rh, double tk, double & moist, double & dmoistdrh, long ipp)
01366 {
01367
01368
01369
01370 long s;
01371 double hmrh, mhmc, smc, u, a, n;
01372
01373 if (rh < 0) rh = 0;
01374
01375 CorD(6,s,0,a1,a2,a3);
01376
01377 switch (s){
01378 case 2:{
01379 if(rh>0.976)
01380 {
01381
01382 CorD(7,kd,rh,smc,a2,a3);
01383 CorD(6,kd,rh,mhmc,a2,mhmc);
01384 give_data(rh,mhmc, smc, 0.976,moist);
01385 dmoistdrh =1000 * (smc - mhmc)/(1 - 0.976);
01386 }
01387 else
01388 {
01389 CorD(6,s,rh,moist,dmoistdrh,a3);
01390 dmoistdrh =1000* derivation_sorption_izoterm_data(rh,tk, 0);
01391 }
01392 break;
01393 }
01394 case 30:{
01395 CorD(6,s,rh,a1,hmrh,a3);
01396 CorD(7,s,rh,smc,a2,a3);
01397 CorD(6,s,rh,mhmc,hmrh,a3);
01398 give_data(rh,mhmc, smc, hmrh,moist);
01399 if (rh < hmrh)
01400 {
01401 dmoistdrh = DerivativeOfTheSorptionIsotherm (rh, tk);
01402 }
01403 else {
01404 dmoistdrh = DerivativeOfTheRetentionCurve(rh, tk, ipp);
01405 }
01406
01407 break;
01408 }
01409 case 31:{
01410 switch (kod){
01411 case 0:{
01412 CorD(6,kod,0,u,a,n);
01413 moist = si_kk_hansen(rh, 0,u,a,n);
01414 dmoistdrh = sorptionizothermDerivation(rh, tk);
01415 break;
01416 }
01417 case 1:{
01418 dmoistdrh = sorptionizothermDerivation(rh, tk);
01419 break;
01420 }
01421 default:{
01422 print_err("unknown definition of Sorption isotherm HANSEN is required",__FILE__,__LINE__,__func__);
01423 }
01424 }
01425 break ;
01426 }
01427 default:{
01428 print_err("unknown definition of Sorption isotherm is required",__FILE__,__LINE__,__func__);
01429 }
01430 }
01431
01432
01433 }
01434
01435 double kunmat::si_kk_hansen(double rh, double tk,double u, double a, double n)
01436 {
01437 return (u*pow((1-(log(rh))/a),(-1/n)));
01438 }
01439
01440
01441
01442
01443
01444
01445
01446
01447 double kunmat::derivation_sorption_izoterm_data(double x1)
01448 {
01449 double derbi;
01450
01451
01452 charid=6;
01453 derbi = data[charid]->getderiv (x1);
01454
01455 return derbi;
01456 }
01457
01458 void kunmat::aux_values (long ipp,double *in,double *inp, double *ineq,double *out)
01459 {
01460 double x1,x2,x1p, x2p, ypv1,ypv2,ypv3,ypv4,ypv5, ypv0;
01461
01462 double w,dwdf, kapaA;
01463 double t;
01464 t = Tp->time;
01465 x1 = in[0];
01466 x2 = in[1];
01467
01468 x1p = inp[0];
01469 x2p = inp[1];
01470
01471
01472 ypv0 =ineq[0];
01473 ypv1 =ineq[1];
01474 ypv2 =ineq[2];
01475 ypv3 =ineq[3];
01476 ypv4 =ineq[4];
01477 ypv5 =ineq[5];
01478
01479 sorption_izoterms_giva_data(0,x1,x2,w,dwdf, ipp);
01480
01481 kapa_values(MatChar[5],ipp,w,ypv0,ypv2, kapaA);
01482 out[3] = kapaA;
01483
01484
01485 out[5] = -10000.0;
01486
01487 out[0] = w;
01488 out[1] = dwdf;
01489
01490 out[4] = -10000.0;
01491
01492
01493
01494 out[2] = data[7]->getval(0.0);
01495
01496
01497
01498
01499 }
01500
01501
01502
01503
01504
01505
01506
01507
01508
01509 void kunmat::save_values (long ipp,double *out)
01510 {
01511 Tm->ip[ipp].eqother[0]=out[0];
01512 Tm->ip[ipp].eqother[1]=out[1];
01513 Tm->ip[ipp].eqother[2]=out[2];
01514 Tm->ip[ipp].eqother[3]=out[3];
01515 }
01516
01517
01518
01519
01520
01521
01522
01523
01524
01525
01526 void kunmat::give_values (long ipp,double *av,double *eq)
01527 {
01528 av[0] = Tm->ip[ipp].av[0];
01529 av[1] = Tm->ip[ipp].av[1];
01530
01531 eq[0] = Tm->ip[ipp].eqother[0];
01532 eq[1] = Tm->ip[ipp].eqother[1];
01533 eq[2] = Tm->ip[ipp].eqother[2];
01534 eq[3] = Tm->ip[ipp].eqother[3];
01535 }
01536
01537 void kunmat::kapa_values (long kod, long ipp,double x1, double xpv, double ineq1, double &kapa)
01538 {
01539 double kapa1, kapa2, kapa3, x1w;
01540 double smc, a, n, x2;
01541 long kd1;
01542
01543
01544 CorD(7,kd1,x1,smc,a,n);
01545
01546 if( x1 > smc)
01547 {
01548 x1 = smc;
01549 }
01550
01551 switch (kod){
01552 case 0: kapa = 0.0;
01553 break;
01554 case 1: CorD(5,kod,x1,kapa,kapa1,kapa2);
01555 break;
01556 case 2: CorD(5,kod,x1,kapa,kapa1,kapa2);
01557 break;
01558 case 30:{
01559 CorD(5,kod,x1,kapa1,kapa2,kapa3);
01560 kapa = kapa_exp(kapa1,kapa2,x1, x2, ipp);
01561 break;
01562 }
01563 default:{
01564 print_err("unknown kod is required",__FILE__,__LINE__,__func__);
01565 }
01566 }
01567 }
01568
01569 void kunmat::initvalues (long ipp,long ido)
01570 {
01571 double x1,x2;
01572
01573 x1 = Tm->ip[ipp].av[0];
01574 x2 = Tm->ip[ipp].av[1];
01575
01576 }
01577
01578