00001 #include <string.h>
00002 #include <math.h>
00003 #include <stdlib.h>
00004 #include "climatcond.h"
00005 #include "globalt.h"
00006 #include "aliast.h"
00007 #include "iotools.h"
00008
00009 climatcond::climatcond (void)
00010 {
00011
00012 T0 = 273.15;
00013
00014 Cwat = 4180.0;
00015
00016 hvap = 2256000.0;
00017
00018 cvap = 1617.0;
00019
00020
00021 last_time = -1.0;
00022
00023 zero = 1.0e-10;
00024
00025
00026
00027 mattyp = 0;
00028
00029 basicdirect = 0.0;
00030
00031 basicinclin = 0.0;
00032
00033
00034 KIND1 = KIND2 = KIND3 = KIND4 = KIND5 = KIND6 = WCorPH = RHorVP = 0;
00035
00036 cfmax = WindFact = WindFactVD = coefAlf = AbsCoef = Albedo = Latitude = EmisCoeff = ExchCoeff = ExchCoeffVD = RainCoeff = MinRainTemp = MinRainFlux = u = a = n = alfaK = Ax = Bx = Cx = 0.0;
00037
00038 SI = numberSI = d = 0;
00039 SIdata = NULL;
00040
00041 aa=0.0; b=0.0;
00042 eps=0.0;
00043 m=0.0;
00044 alpha=0.0;
00045 sig=5.67e-8;
00046 tn = 0;
00047
00048
00049
00050 TepVal=0.0;
00051
00052 RHVal=0.0;
00053
00054 DifSVal=0.0;
00055
00056 DirSVal=0.0;
00057
00058 WDVal=0.0;
00059
00060 WSVal=0.0;
00061
00062 VRVal=0.0;
00063
00064 LWVal=0.0;
00065
00066 SLWVal=0.0;
00067
00068 NorRainVal=0.0;
00069
00070 GasFluxVal=0.0;
00071
00072 HeatFluxVal=0.0;
00073
00074 ShWRadVal=0.0;
00075
00076 LoWRadVal=0.0;
00077
00078 RadVal=0.0;
00079
00080 TempbVal=0.0;
00081
00082 files0 = new char[1000];
00083 files1 = new char[1000];
00084 files2 = new char[1000];
00085 files3 = new char[1000];
00086 files4 = new char[1000];
00087 files5 = new char[1000];
00088 files6 = new char[1000];
00089 files7 = new char[1000];
00090 files8 = new char[1000];
00091 files14 = new char[1000];
00092 files15 = new char[1000];
00093
00094 }
00095
00096 climatcond::~climatcond (void)
00097 {
00098 if(SIdata != NULL)
00099 delete [] SIdata ;
00100
00101
00102
00103
00104 delete [] files0;
00105 delete [] files1;
00106 delete [] files2;
00107 delete [] files3;
00108 delete [] files4;
00109 delete [] files5;
00110 delete [] files6;
00111 delete [] files7;
00112 delete [] files8;
00113 delete [] files14;
00114 delete [] files15;
00115
00116 }
00117
00118
00119
00120
00121
00122
00123
00124
00125
00126
00127
00128
00129
00130 void climatcond::read (XFILE *in, long readclimfile)
00131 {
00132 long klier=0;
00133
00134
00135 xfscanf (in,"%d",(long*)&mattyp);
00136
00137 xfscanf (in,"%lf",&basicdirect);
00138
00139 xfscanf (in,"%lf",&basicinclin);
00140
00141
00142
00143
00144 xfscanf (in,"%ld",&KIND1);
00145
00146 switch (KIND1){
00147 case 0:{
00148
00149 break;
00150 }
00151 case 12:{
00152
00153 xfscanf (in,"%lf",&coefAlf);
00154 break;
00155 }
00156 case 13:{
00157
00158
00159 xfscanf (in,"%lf %lf",&coefAlf, &WindFact);
00160 break;
00161 }
00162 case 14:{
00163
00164
00165
00166
00167
00168
00169 xfscanf (in,"%lf %lf %lf %lf %lf",&aa,&b,&eps,&alpha,&m);
00170 break;
00171 }
00172 case 15:{
00173
00174
00175
00176
00177 xfscanf (in,"%lf %lf %lf",&b,&aa,&eps);
00178
00179 klier=1;
00180 break;
00181 }
00182 case 16:{
00183
00184
00185 xfscanf (in,"%lf %lf",&tint,&coefAlf);
00186 break;
00187 }
00188 default:{
00189 print_err("unknown definition of Read clim.data-KIND1 is required",__FILE__,__LINE__,__func__);
00190 }
00191 }
00192
00193
00194
00195
00196 xfscanf (in,"%ld",&KIND2);
00197
00198 switch (KIND2){
00199 case 0:{
00200
00201 break;
00202 }
00203 case 21:{
00204
00205 xfscanf (in,"%lf",&AbsCoef);
00206 break;
00207 }
00208 case 22:{
00209
00210
00211
00212 xfscanf (in,"%lf %lf %lf",&AbsCoef,&Albedo,&Latitude);
00213 break;
00214 }
00215 default:{
00216 print_err("unknown definition of Read clim.data-KIND2 is required",__FILE__,__LINE__,__func__);
00217 }
00218 }
00219
00220
00221
00222
00223 xfscanf (in,"%ld",&KIND3);
00224
00225 switch (KIND3){
00226 case 0:{
00227
00228 break;
00229 }
00230 case 31:{
00231
00232 xfscanf (in,"%lf",&EmisCoeff);
00233 break;
00234 }
00235 case 32:{
00236
00237 xfscanf (in,"%lf",&EmisCoeff);
00238 break;
00239 }
00240 default:{
00241 print_err("unknown definition of Read clim.data-KIND3 is required",__FILE__,__LINE__,__func__);
00242 }
00243 }
00244
00245
00246
00247 xfscanf (in,"%ld",&KIND4);
00248
00249 switch (KIND4){
00250 case 0:{
00251
00252 break;
00253 }
00254 default:{
00255 print_err("unknown definition of Read clim.data-KIND4 is required",__FILE__,__LINE__,__func__);
00256 }
00257 }
00258
00259
00260
00261
00262
00263 xfscanf (in,"%ld",&KIND5);
00264
00265 switch (KIND5){
00266 case 0:{
00267
00268 break;
00269 }
00270 case 51:{
00271
00272
00273
00274
00275 xfscanf (in,"%lf %lf %lf %lf",&ExchCoeff,&RainCoeff,&MinRainTemp,&MinRainFlux);
00276 break;
00277 }
00278 case 52:{
00279
00280
00281
00282
00283 xfscanf (in,"%lf %lf %lf %lf",&ExchCoeff,&RainCoeff,&MinRainTemp,&MinRainFlux);
00284 break;
00285 }
00286 default:{
00287 print_err("unknown definition of Read clim.data-KIND5 is required",__FILE__,__LINE__,__func__);
00288 }
00289 }
00290
00291
00292
00293
00294
00295 xfscanf (in,"%ld",&KIND6);
00296
00297 switch (KIND6){
00298 case 0:{
00299
00300 break;
00301 }
00302 case 61:{
00303
00304
00305
00306
00307 xfscanf (in,"%lf %ld",&ExchCoeffVD,&RHorVP);
00308 break;
00309 }
00310 case 62:{
00311
00312
00313
00314
00315
00316 xfscanf (in,"%lf %ld %lf",&ExchCoeffVD,&RHorVP,&WindFactVD);
00317 break;
00318 }
00319 case 63:{
00320
00321
00322
00323
00324
00325 xfscanf (in,"%lf %ld %lf",&ExchCoeffVD,&RHorVP,&alfaK);
00326 break;
00327 }
00328 case 64:{
00329
00330
00331
00332
00333
00334 xfscanf (in,"%lf %ld %lf",&ExchCoeffVD,&RHorVP,&alfaK);
00335 break;
00336 }
00337 case 65:{
00338
00339
00340
00341
00342 xfscanf (in,"%lf %ld %lf %lf",&Ax,&RHorVP,&Bx,&Cx);
00343 break;
00344 }
00345 case 66:{
00346
00347
00348
00349
00350 xfscanf (in,"%ld %d",&RHorVP,&d);
00351
00352 exchangecoeff.tfunc=tab;
00353 exchangecoeff.tabf = new tablefunct();
00354 exchangecoeff.tabf->asize = d;
00355 exchangecoeff.read (in);
00356
00357 break;
00358 }
00359
00360 case 67:{
00361
00362
00363
00364 xfscanf (in,"%lf %lf %lf",&rhint,&ExchCoeffVD,&rhtemp);
00365 break;
00366 }
00367 case 68:{
00368
00369
00370
00371 xfscanf (in,"%lf %lf %lf",&rhint,&ExchCoeffVD,&rhtemp);
00372 break;
00373 }
00374 case 70:{
00375 xfscanf (in,"%lf %lf %lf",&rhint, &rhalfa, &rhtemp);
00376 break;
00377 }
00378 default:{
00379 print_err ("unknown definition of Read clim.data-KIND6 is required",__FILE__,__LINE__,__func__);
00380 }
00381 }
00382
00383
00384
00385 xfscanf(in, " %a", files0);
00386 if (readclimfile)
00387 read_climatic_files(files0,0);
00388
00389
00390 xfscanf(in, " %a", files1);
00391 if (readclimfile)
00392 read_climatic_files(files1,1);
00393
00394
00395 xfscanf(in, " %a", files2);
00396 if (readclimfile)
00397 read_climatic_files(files2,2);
00398
00399
00400 xfscanf(in, " %a", files3);
00401 if (readclimfile)
00402 read_climatic_files(files3,3);
00403
00404
00405 xfscanf(in, " %a", files4);
00406 if (readclimfile)
00407 read_climatic_files(files4,4);
00408
00409
00410 xfscanf(in, " %a", files5);
00411 if (readclimfile)
00412 read_climatic_files(files5,5);
00413
00414
00415 xfscanf(in, " %a", files6);
00416 if (readclimfile)
00417 read_climatic_files(files6,6);
00418
00419
00420 xfscanf(in, " %a", files7);
00421 if (readclimfile)
00422 read_climatic_files(files7,7);
00423
00424
00425 xfscanf(in, " %a", files8);
00426 if (readclimfile)
00427 read_climatic_files(files8,8);
00428
00429
00430 xfscanf(in, " %a", files14);
00431 if (readclimfile)
00432 read_climatic_files(files14,14);
00433
00434 if (klier==1){
00435
00436 xfscanf(in, " %a", files15);
00437 if (readclimfile)
00438 read_climatic_files(files15,15);
00439 }
00440
00441 }
00442
00443
00444
00445
00446
00447
00448
00449
00450 void climatcond::print(FILE *out)
00451 {
00452 fprintf (out,"\n %ld ",mattyp);
00453 fprintf (out," %le ",basicdirect);
00454 fprintf (out," %le ",basicinclin);
00455
00456
00457 fprintf (out,"\n %ld ",KIND1);
00458
00459 switch (KIND1){
00460 case 0:{
00461
00462 break;
00463 }
00464 case 11:{
00465
00466 fprintf (out," %e ",coefAlf);
00467 break;
00468 }
00469 case 12:{
00470
00471 fprintf (out," %e ",coefAlf);
00472 break;
00473 }
00474 case 13:{
00475
00476 fprintf (out," %e %e ",coefAlf,WindFact);
00477 break;
00478 }
00479 case 14:{
00480
00481
00482
00483
00484
00485
00486 fprintf (out," %e %e %e %e %e ",aa,b,eps,alpha,m);
00487 break;
00488 }
00489 case 15:{
00490
00491
00492
00493
00494 fprintf (out," %ld %le %le",d,aa,eps);
00495 break;
00496 }
00497 case 16:{
00498
00499 fprintf (out," %e %e",tint,coefAlf);
00500 break;
00501 }
00502 default:{
00503 print_err("unknown definition of Read clim.data-KIND1 is required",__FILE__,__LINE__,__func__);
00504 }
00505 }
00506
00507
00508
00509 fprintf (out,"\n %ld ",KIND2 );
00510
00511 switch (KIND2){
00512 case 0:{
00513
00514 break;
00515 }
00516 case 21:{
00517
00518 fprintf (out," %e ",AbsCoef );
00519 break;
00520 }
00521 case 22:{
00522
00523 fprintf (out," %e %e %e ",AbsCoef , Albedo , Latitude );
00524 break;
00525 }
00526 default:{
00527 print_err ("unknown definition of clim.data-KIND2 is required",__FILE__,__LINE__,__func__);
00528 }
00529 }
00530
00531
00532 fprintf (out,"\n %ld ",KIND3);
00533
00534 switch (KIND3){
00535 case 0:{
00536
00537 break;
00538 }
00539 case 31:{
00540
00541 fprintf (out," %e ",EmisCoeff );
00542 break;
00543 }
00544 case 32:{
00545
00546 fprintf (out," %e ",EmisCoeff );
00547 break;
00548 }
00549 default:{
00550 print_err ("unknown definition of clim.data-KIND3 is required",__FILE__,__LINE__,__func__);
00551 }
00552 }
00553
00554
00555 fprintf (out,"\n %ld ",KIND4);
00556
00557 switch (KIND4){
00558 case 0:{
00559
00560 break;
00561 }
00562 case 41:{
00563
00564 fprintf (out," %ld ",WCorPH );
00565 break;
00566 }
00567 default:{
00568 print_err ("unknown definition of clim.data-KIND4 is required",__FILE__,__LINE__,__func__);
00569 }
00570 }
00571
00572
00573 fprintf (out,"\n %ld ",KIND5);
00574
00575
00576 switch (KIND5){
00577 case 0:{
00578
00579 break;
00580 }
00581 case 51:{
00582
00583
00584
00585
00586
00587 fprintf (out," %e %e %e %e ",ExchCoeff , RainCoeff , MinRainTemp, MinRainFlux);
00588 break;
00589 }
00590 case 52:{
00591
00592
00593
00594
00595
00596 fprintf (out," %e %e %e %e ",ExchCoeff , RainCoeff , MinRainTemp, MinRainFlux);
00597 break;
00598 }
00599 default:{
00600 print_err ("unknown definition of clim.data-KIND5 is required",__FILE__,__LINE__,__func__);
00601 }
00602 }
00603
00604
00605 fprintf (out,"\n %ld ",KIND6);
00606
00607
00608 switch (KIND6){
00609 case 0:{
00610 break;
00611 }
00612 case 61:{
00613
00614
00615 fprintf (out," %e %ld ",ExchCoeffVD, RHorVP);
00616 break;
00617 }
00618 case 62:{
00619
00620
00621
00622 fprintf (out," %e %ld %e ", ExchCoeffVD, RHorVP, WindFactVD);
00623 break;
00624 }
00625 case 63:{
00626
00627
00628 fprintf (out," %e %ld %e ", ExchCoeffVD, RHorVP, alfaK);
00629 break;
00630 }
00631 case 64:{
00632
00633
00634 fprintf (out," %e %ld %e ", ExchCoeffVD, RHorVP, alfaK);
00635 break;
00636 }
00637 case 65:{
00638 fprintf (out," %e %ld %e %e ", Ax, RHorVP, Bx, Cx);
00639 break;
00640 }
00641 case 66:{
00642 fprintf (out," %ld ",RHorVP);
00643 exchangecoeff.print (out);
00644 break;
00645 }
00646 case 67:{
00647
00648
00649
00650 fprintf (out," %e %e %e",rhint, ExchCoeffVD, rhtemp);
00651 break;
00652 }
00653 default:{
00654 print_err ("unknown definition of clim.data-KIND6 is required",__FILE__,__LINE__,__func__);
00655 }
00656 }
00657
00658
00659 fprintf (out,"\n%s",files0);
00660 fprintf (out,"\n%s",files1);
00661 fprintf (out,"\n%s",files2);
00662 fprintf (out,"\n%s",files3);
00663 fprintf (out,"\n%s",files4);
00664 fprintf (out,"\n%s",files5);
00665 fprintf (out,"\n%s",files6);
00666 fprintf (out,"\n%s",files7);
00667 fprintf (out,"\n%s",files8);
00668 fprintf (out,"\n%s",files14);
00669 if (KIND1==15)
00670 fprintf (out,"\n%s",files15);
00671
00672 fprintf (out,"\n");
00673 }
00674
00675
00676
00677
00678
00679
00680
00681
00682
00683
00684 void climatcond::read_climatic_files (char *Cfile, long cod)
00685 {
00686 long time;
00687 long i;
00688 double Val;
00689
00690 XFILE *ffile;
00691
00692
00693
00694
00695
00696
00697
00698
00699
00700
00701
00702
00703
00704
00705
00706
00707
00708
00709
00710
00711
00712 if (strcmp(Cfile,"0") == 0)
00713 return;
00714
00715 ffile = xfopen(Cfile,"r");
00716
00717 if (ffile==NULL){
00718 print_err ("Input file %s for climatic condition is not found. \n",__FILE__,__LINE__,__func__,Cfile);
00719 abort();
00720 }
00721
00722 xfscanf (ffile,"%ld",&tn);
00723
00724 switch (cod){
00725 case 0:{
00726
00727 TepValue.tfunc=tab;
00728 TepValue.tabf = new tablefunct (tn,1);
00729 break;
00730 }
00731 case 1:{
00732
00733 RHValue.tfunc=tab;
00734 RHValue.tabf = new tablefunct (tn,1);
00735 break;
00736 }
00737 case 2:{
00738
00739 DifSValue.tfunc=tab;
00740 DifSValue.tabf = new tablefunct (tn,1);
00741 break;
00742 }
00743 case 3:{
00744
00745 DirSValue.tfunc=tab;
00746 DirSValue.tabf = new tablefunct (tn,1);
00747 break;
00748 }
00749 case 4:{
00750
00751 WDValue.tfunc=tab;
00752 WDValue.tabf = new tablefunct (tn,1);
00753 break;
00754 }
00755 case 5:{
00756
00757 WSValue.tfunc=tab;
00758 WSValue.tabf = new tablefunct (tn,1);
00759 break;
00760 }
00761 case 6:{
00762
00763 VRValue.tfunc=tab;
00764 VRValue.tabf = new tablefunct (tn,1);
00765 break;
00766 }
00767 case 7:{
00768
00769 LWValue.tfunc=tab;
00770 LWValue.tabf = new tablefunct (tn,1);
00771 break;
00772 }
00773 case 8:{
00774
00775 SLWValue.tfunc=tab;
00776 SLWValue.tabf = new tablefunct (tn,1);
00777 break;
00778 }
00779 case 9:{
00780
00781 NorRainValue.tfunc=tab;
00782 NorRainValue.tabf = new tablefunct (tn,1);
00783 break;
00784 }
00785 case 10:{
00786
00787 GasFluxValue.tfunc=tab;
00788 GasFluxValue.tabf = new tablefunct (tn,1);
00789 break;
00790 }
00791 case 11:{
00792
00793 HeatFluxValue.tfunc=tab;
00794 HeatFluxValue.tabf = new tablefunct (tn,1);
00795 break;
00796 }
00797 case 12:{
00798
00799 ShWRadValue.tfunc=tab;
00800 ShWRadValue.tabf = new tablefunct (tn,1);
00801 break;
00802 }
00803 case 13:{
00804
00805 LoWRadValue.tfunc=tab;
00806 LoWRadValue.tabf = new tablefunct (tn,1);
00807 break;
00808 }
00809 case 14:{
00810
00811 RadValue.tfunc=tab;
00812 RadValue.tabf = new tablefunct (tn,1);
00813 break;
00814 }
00815 case 15:{
00816
00817 TempbValue.tfunc=tab;
00818 TempbValue.tabf = new tablefunct (tn,1);
00819 break;
00820 }
00821 default:{
00822 print_err ("unknown definition of ReadCllimaticFiles is required",__FILE__,__LINE__,__func__);
00823 }
00824 }
00825
00826
00827
00828
00829
00830 for (i=0;i<tn;i++){
00831
00832
00833
00834
00835
00836
00837
00838 xfscanf (ffile,"%ld %lf", &time, &Val);
00839 switch (cod){
00840 case 0:{
00841
00842 TepValue.tabf->x[i]=time;
00843 TepValue.tabf->y[i]=Val;
00844 break;
00845 }
00846 case 1:{
00847
00848 RHValue.tabf->x[i]=time;
00849 RHValue.tabf->y[i]=Val;
00850 break;
00851 }
00852 case 2:{
00853
00854 DifSValue.tabf->x[i]=time;
00855 DifSValue.tabf->y[i]=Val;
00856 break;
00857 }
00858 case 3:{
00859
00860 DirSValue.tabf->x[i]=time;
00861 DirSValue.tabf->y[i]=Val;
00862 break;
00863 }
00864 case 4:{
00865
00866 WDValue.tabf->x[i]=time;
00867 WDValue.tabf->y[i]=Val;
00868 break;
00869 }
00870 case 5:{
00871
00872 WSValue.tabf->x[i]=time;
00873 WSValue.tabf->y[i]=Val;
00874 break;
00875 }
00876 case 6:{
00877
00878 VRValue.tabf->x[i]=time;
00879 VRValue.tabf->y[i]=Val;
00880 break;
00881 }
00882 case 7:{
00883
00884 LWValue.tabf->x[i]=time;
00885 LWValue.tabf->y[i]=Val;
00886 break;
00887 }
00888 case 8:{
00889
00890 SLWValue.tabf->x[i]=time;
00891 SLWValue.tabf->y[i]=Val;
00892 break;
00893 }
00894 case 9:{
00895
00896 NorRainValue.tabf->x[i]=time;
00897 NorRainValue.tabf->y[i]=Val;
00898 break;
00899 }
00900 case 10:{
00901
00902 GasFluxValue.tabf->x[i]=time;
00903 GasFluxValue.tabf->y[i]=Val;
00904 break;
00905 }
00906 case 11:{
00907
00908 HeatFluxValue.tabf->x[i]=time;
00909 HeatFluxValue.tabf->y[i]=Val;
00910 break;
00911 }
00912 case 12:{
00913
00914 ShWRadValue.tabf->x[i]=time;
00915 ShWRadValue.tabf->y[i]=Val;
00916 break;
00917 }
00918 case 13:{
00919
00920 LoWRadValue.tabf->x[i]=time;
00921 LoWRadValue.tabf->y[i]=Val;
00922 break;
00923 }
00924 case 14:{
00925
00926 RadValue.tabf->x[i]=time;
00927 RadValue.tabf->y[i]=Val;
00928 break;
00929 }
00930 case 15:{
00931
00932 TempbValue.tabf->x[i]=time;
00933 TempbValue.tabf->y[i]=Val;
00934 break;
00935 }
00936 default:{
00937 print_err ("unknown definition of ReadCllimaticFiles is required",__FILE__,__LINE__,__func__);
00938 }
00939 }
00940
00941
00942 }
00943
00944
00945 xfclose (ffile);
00946
00947 }
00948
00949
00950
00951
00952
00953
00954
00955
00956
00957
00958
00959
00960
00961 double climatcond::give_flux (long flid, double *nval,long ipp,long nid)
00962 {
00963 double flux,time;
00964 mednamest mn;
00965
00966
00967 time = Tp->time;
00968
00969
00970 prepare_values (time);
00971
00972
00973 mn = Tp->mednam;
00974
00975 switch (mn){
00976 case 1:{
00977
00978 flux = give_heat_flux (0,nval[0],ipp, nid);
00979 break;
00980 }
00981 case 2:{
00982
00983 flux = give_moisture_flux (nval[0], -1000,ipp, nid);
00984 break;
00985 }
00986 case 10:{
00987
00988 switch (flid){
00989 case 0:{
00990
00991 flux = give_moisture_flux (nval[0], nval[1],ipp, nid);
00992 break;
00993 }
00994 case 1:{
00995
00996 flux = give_heat_flux (nval[0],nval[1],ipp, nid);
00997 break;
00998 }
00999 default:{
01000 print_err ("unknown flux is required",__FILE__,__LINE__,__func__);
01001 }
01002 }
01003 break;
01004 }
01005 case 25:{
01006
01007 switch (flid){
01008 case 0:{
01009
01010 flux = give_moisture_flux (nval[0], 293.15,ipp, nid);
01011 break;
01012 }
01013 case 1:{
01014
01015 flux = give_salt_flux (nval[0],nval[1], 293.15, nid);
01016 break;
01017 }
01018 default:{
01019 print_err ("unknown flux is required",__FILE__,__LINE__,__func__);
01020 }
01021 }
01022 break;
01023 }
01024 case 30:{
01025
01026 switch (flid){
01027 case 0:{
01028
01029 flux = give_moisture_flux (nval[0],nval[3],ipp,nid);
01030 break;
01031 }
01032 case 1:{
01033
01034 flux = give_salt_flux (nval[0],nval[1],nval[3],nid);
01035 break;
01036 }
01037 case 3:{
01038
01039 flux = give_heat_flux (nval[0],nval[3],ipp,nid);
01040 break;
01041 }
01042 default:{
01043 print_err ("unknown flux is required",__FILE__,__LINE__,__func__);
01044 }
01045 }
01046 break;
01047 }
01048 default:{
01049 print_err ("unknown medium/media flux is required",__FILE__,__LINE__,__func__);
01050 }
01051 }
01052
01053 return flux;
01054 }
01055
01056
01057
01058
01059
01060
01061
01062
01063
01064
01065
01066
01067
01068
01069
01070
01071
01072
01073
01074
01075
01076
01077
01078 double climatcond::give_heat_flux (double nval_0,double temp,long ipp, long nid)
01079 {
01080 double JH1, JH2, JH3, JH7, JH8;
01081 double time;
01082 mednamest mn;
01083
01084
01085 mn = Tp->mednam;
01086
01087
01088 time = Tp->time;
01089
01090 switch (KIND1){
01091 case 0:{
01092 JH1 = 0;
01093 break;
01094 }
01095 case 11:{
01096
01097 JH1 = condit_heat_flux (11,temp);
01098 break;
01099 }
01100 case 12:{
01101
01102 JH1 = condit_heat_flux (12,temp);
01103 break;
01104 }
01105 case 13:{
01106
01107 JH1 = condit_heat_flux (13,temp);
01108 break;
01109 }
01110 case 14:{
01111
01112 JH1 = condit_heat_flux (14,temp);
01113 break;
01114 }
01115 case 15:{
01116
01117 JH1 = condit_heat_flux (15,temp);
01118 break;
01119 }
01120 case 16:{
01121
01122 JH1 = condit_heat_flux (16,temp);
01123 break;
01124 }
01125 default:{
01126 print_err ("unknown definition of Get_Heat_Flux-KIND1 is required",__FILE__,__LINE__,__func__);
01127 }
01128 }
01129
01130
01131
01132
01133 switch (KIND2){
01134 case 0:{
01135
01136 JH2 = 0;
01137 break;
01138 }
01139 case 21:{
01140
01141 JH2 =-1* condit_short_wave_radiat_flux(21,time);
01142 break;
01143 }
01144 case 22:{
01145
01146 JH2 =-1* condit_short_wave_radiat_flux(22,time);
01147 break;
01148 }
01149 default:{
01150 print_err ("unknown definition of Get_Heat_Flux-KIND2 is required",__FILE__,__LINE__,__func__);
01151 }
01152 }
01153
01154
01155
01156 switch(KIND3){
01157 case 0:{
01158
01159 JH3 = 0;
01160 break;
01161 }
01162 case 31:{
01163
01164 JH3 =-1* condit_long_wave_radiat_flux (31);
01165 break;
01166 }
01167 case 32:{
01168
01169 JH3 =-1* condit_long_wave_radiat_flux (32);
01170 break;
01171 }
01172 default:{
01173 print_err ("unknown definition of Get_Heat_Flux-KIND3 is required",__FILE__,__LINE__,__func__);
01174 }
01175 }
01176
01177
01178 switch(mn){
01179 case 30:
01180 case 10:{
01181 switch (KIND5){
01182 case 0:{
01183 JH7 = 0;
01184 break;
01185 }
01186 case 51:{
01187
01188 JH7 = condit_rain_heat_flux(51,nval_0,nid);
01189 break;
01190 }
01191 case 52:{
01192
01193 JH7 = condit_rain_heat_flux(52,nval_0,nid);
01194 break;
01195 }
01196 default:{
01197 print_err ("unknown definition of Get_Heat_Flux-KIND5 is required",__FILE__,__LINE__,__func__);
01198 }
01199 }
01200 break;
01201 }
01202 default:{
01203 JH7 = 0;
01204 break;
01205 }
01206 }
01207
01208
01209 switch(mn){
01210 case 30:
01211 case 10:{
01212 switch (KIND6){
01213 case 0:{
01214 JH8 =0;
01215 break;
01216 }
01217 case 61:{
01218
01219 JH8 = condit_vapour_diffusion_heat_flux(61,nval_0,temp,ipp, nid);
01220 break;
01221 }
01222 case 62:{
01223
01224
01225
01226 JH8 = condit_vapour_diffusion_heat_flux(62,nval_0,temp,ipp, nid);
01227 break;
01228 }
01229 case 63:{
01230
01231 JH8 = condit_vapour_diffusion_heat_flux(63,nval_0,temp,ipp, nid);
01232 break;
01233 }
01234 case 64:{
01235
01236 JH8 = condit_vapour_diffusion_heat_flux(64,nval_0,temp,ipp, nid);
01237 break;
01238 }
01239 case 65:{
01240
01241 JH8 = condit_vapour_diffusion_heat_flux(65,nval_0,temp,ipp, nid);
01242 break;
01243 }
01244 case 66:{
01245
01246 JH8 = condit_vapour_diffusion_heat_flux(66,nval_0,temp,ipp, nid) ;
01247 break;
01248 }
01249 case 67:{
01250
01251 JH8 = condit_vapour_diffusion_heat_flux(67,nval_0,temp,ipp, nid) ;
01252 break;
01253 }
01254 default:{
01255 print_err ("unknown definition of Get_Heat_Flux-KIND6 is required",__FILE__,__LINE__,__func__);
01256 }
01257 }
01258 break;
01259 }
01260 default:{
01261 JH8 = 0;
01262 break;
01263 }
01264 }
01265
01266
01267 return(JH1 +JH2 +JH3 + JH7 + JH8);
01268 }
01269
01270
01271
01272
01273
01274
01275
01276
01277
01278
01279
01280
01281 double climatcond::condit_heat_flux (long tkind, double temperstr)
01282 {
01283 double betaWind,coefAlf0;
01284
01285
01286 switch (tkind){
01287 case 11:{
01288
01289
01290 Jkq = HeatFluxVal;
01291 break;
01292 }
01293 case 12:{
01294
01295
01296 Jkq = -1.0 * coefAlf * (TepVal - temperstr);
01297 break;
01298 }
01299 case 13:{
01300
01301
01302
01303
01304
01305
01306 if ((fabs(basicdirect - WDVal)) <= M_PI){
01307 betaWind = fabs (basicdirect - WDVal);
01308 }else{
01309 betaWind = 2.0*M_PI - (fabs (basicdirect - WDVal));
01310 }
01311
01312 if (betaWind < (M_PI/2.0)){
01313 coefAlf0 = coefAlf + WindFact*sqrt(WSVal);
01314 }else{
01315 coefAlf0 = coefAlf;
01316 }
01317
01318 Jkq = -1.0 *(coefAlf0 * (TepVal - temperstr));
01319
01320 break;
01321 }
01322 case 14:{
01323
01324
01325
01326
01327
01328
01329 coefAlf0 = b + aa * WSVal;
01330 Jkq = -1.0*(coefAlf0 *(TepVal - temperstr));
01331 Jkq += -1.0*eps*sig*(TepVal*TepVal*TepVal*TepVal - temperstr*temperstr*temperstr*temperstr);
01332 Jkq += -1.0*alpha*m*RadVal;
01333
01334 break;
01335 }
01336 case 15:{
01337
01338
01339
01340 Jkq = b *(temperstr - TepVal);
01341 Jkq += aa * (temperstr - TempbVal);
01342 Jkq += eps*sig*(temperstr*temperstr*temperstr*temperstr - TepVal*TepVal*TepVal*TepVal);
01343 Jkq += RadVal;
01344
01345 break;
01346 }
01347 case 16:{
01348
01349
01350
01351 Jkq = -1.0 * coefAlf *((tint+273.15) - temperstr);
01352
01353 break;
01354 }
01355
01356 default:{
01357 print_err("unknown definition of ConditHeatFlux is required",__FILE__,__LINE__,__func__);
01358 }
01359 }
01360
01361
01362 return Jkq;
01363 }
01364
01365
01366
01367
01368
01369
01370
01371
01372
01373
01374 double climatcond::condit_short_wave_radiat_flux (long swkind, double time)
01375 {
01376 double MAXSUNDECL = 0.4089306437;
01377 double HOURFREQU = 0.26179939;
01378 double YEARFREQU =0.017214206;
01379 double SPRINGBEGIN =80;
01380 double th, dt,day, sdt,cdt, hold, hnew, sold, snew, alf, r;
01381 double JkQdif, JkQdir,Jq;
01382
01383 switch (swkind){
01384 case 21:{
01385
01386
01387
01388 Jq=AbsCoef*ShWRadVal;
01389 break;
01390 }
01391 case 22:{
01392
01393
01394
01395
01396
01397
01398
01399
01400
01401 th = time/3600;
01402 day = time/86400;
01403 dt = MAXSUNDECL * sin((YEARFREQU * (day-SPRINGBEGIN )));
01404
01405 sdt = sin(Latitude)*sin(dt);
01406 cdt = cos(Latitude)*cos(dt);
01407
01408 hold = asin(sdt - cdt*cos(HOURFREQU * (th - 2)));
01409 hnew = asin(sdt - cdt*cos(HOURFREQU * (th - 1)));
01410
01411 sold = cos(dt)*(sin(HOURFREQU * (th-2)))/(cos(hold));
01412 snew = cos(dt)*(sin(HOURFREQU * (th-1)))/(cos(hnew));
01413
01414 if((th>12) && (snew > sold)){
01415 alf = 2*M_PI + asin(snew);
01416 }
01417 if ((th <= 12) || (snew <= sold)){
01418 alf = M_PI - asin(snew);
01419 }
01420
01421 r = cos(basicinclin) + sin(basicinclin)*(cos(alf -basicdirect))/(tan(hnew));
01422
01423 if (r<0.0){
01424 JkQdir = 0.0;
01425 }else{
01426 if (r>0.5){
01427 JkQdir = 0.5 *DirSVal;
01428 }else{
01429 JkQdir = r * DirSVal;
01430 }
01431 }
01432
01433 JkQdif = DifSVal*cos(basicinclin/2.0)*cos(basicinclin/2.0) + Albedo *(DifSVal +DirSVal)*sin(basicinclin/2.0)*sin(basicinclin/2.0);
01434
01435 Jq=AbsCoef*(JkQdif +JkQdir);
01436
01437 break;
01438 }
01439 default:{
01440 print_err("unknown definition of ConditHeatFlux is required",__FILE__,__LINE__,__func__);
01441 }
01442 }
01443
01444 return (Jq);
01445
01446 }
01447
01448
01449
01450
01451
01452
01453
01454 double climatcond::condit_long_wave_radiat_flux (long lwkind)
01455 {
01456 double flux;
01457
01458 switch (lwkind){
01459 case 31:{
01460
01461
01462
01463 flux = EmisCoeff*LoWRadVal;
01464 break;
01465 }
01466 case 32:{
01467
01468
01469
01470
01471
01472 flux = (EmisCoeff * (LWVal + SLWVal));
01473 break;
01474 }
01475 default:{
01476 print_err("unknown definition of ConditLongWaveRadiatFlux is required",__FILE__,__LINE__,__func__);
01477 }
01478 }
01479
01480 return (flux);
01481 }
01482
01483
01484
01485
01486
01487
01488 double climatcond::condit_rain_heat_flux (long rainkind, double nval_0,long nid)
01489 {
01490 double RainMoistureFlux, TemperRain;
01491 double uwT;
01492
01493
01494
01495
01496
01497 RainMoistureFlux =-1.0 * condit_rain_moisture_flux (rainkind, nval_0,nid);
01498
01499 TemperRain = pow (RHVal,0.1247) * (109.8 + TepVal) - 109.8;
01500
01501 uwT = Cwat *(TemperRain - T0);
01502
01503 return uwT*RainMoistureFlux;
01504 }
01505
01506
01507
01508
01509
01510
01511 double climatcond::condit_vapour_diffusion_heat_flux (long vdkind,double nval_0, double nval_1,long ipp, long nid)
01512 {
01513 double JkmvDiff, Tvap, hv;
01514
01515
01516
01517
01518
01519
01520
01521
01522
01523
01524
01525
01526
01527
01528 JkmvDiff = condit_vapour_diffusion_moisture_flux (vdkind,nval_0,nval_1,ipp, nid);
01529
01530 Tvap = (TepVal + nval_1 )/2;
01531
01532
01533
01534 hv = cvap * (Tvap -T0) + hvap ;
01535
01536 return (-1.0 * hv * JkmvDiff);
01537 }
01538
01539
01540
01541
01542
01543
01544
01545
01546
01547
01548
01549
01550
01551
01552
01553
01554
01555
01556 double climatcond::give_moisture_flux (double nval_0, double nval_1,long ipp, long nid)
01557 {
01558
01559
01560
01561
01562
01563
01564
01565 double jh4, jh5, jh6;
01566
01567
01568 switch (KIND4){
01569 case 0:{
01570
01571 jh4 = 0;
01572 break;
01573 }
01574 case 41:{
01575 jh4 = condit_water_contact ();
01576 break;
01577 }
01578 default:{
01579 print_err("unknown definition of Get_moisture_flux-KIND4 is required",__FILE__,__LINE__,__func__);
01580 }
01581 }
01582
01583 switch (KIND5){
01584 case 0:{
01585
01586 jh5 = 0;
01587 break;
01588 }
01589 case 51:{
01590
01591 jh5 = -1.0 * condit_rain_moisture_flux (51,nval_0,nid);
01592 break;
01593 }
01594 case 52:{
01595
01596 jh5 = -1.0 * condit_rain_moisture_flux (52,nval_0,nid);
01597 break;
01598 }
01599 default:{
01600 print_err("unknown definition of Get_moisture_flux-KIND5 is required",__FILE__,__LINE__,__func__);
01601 }
01602 }
01603
01604 if (nval_1 != -1000){
01605 switch(KIND6){
01606 case 0:{
01607 jh6 = 0;
01608 break;
01609 }
01610 case 61:{
01611 jh6 = condit_vapour_diffusion_moisture_flux (61,nval_0,nval_1,ipp, nid);
01612 break;
01613 }
01614 case 62:{
01615
01616
01617
01618 jh6 = condit_vapour_diffusion_moisture_flux (62,nval_0,nval_1,ipp, nid);
01619 break;
01620 }
01621 case 63:{
01622 jh6 = condit_vapour_diffusion_moisture_flux (63,nval_0,nval_1,ipp, nid);
01623 break;
01624 }
01625 case 64:{
01626 jh6 = condit_vapour_diffusion_moisture_flux (64,nval_0,nval_1,ipp, nid);
01627 break;
01628 }
01629 case 65:{
01630 jh6 = condit_vapour_diffusion_moisture_flux (65,nval_0,nval_1,ipp, nid);
01631 break;
01632 }
01633 case 66:{
01634 jh6 = condit_vapour_diffusion_moisture_flux (66,nval_0,nval_1,ipp, nid);
01635 break;
01636 }
01637 case 67:{
01638 jh6 = condit_vapour_diffusion_moisture_flux (67,nval_0,nval_1,ipp, nid);
01639 break;
01640 }
01641 case 68:{
01642 double betaExch, pv,VapPres;
01643 betaExch = ExchCoeffVD;
01644 pv = nval_0;
01645 VapPres = rhint;
01646 jh6 = (betaExch *(pv - VapPres));
01647 fprintf(Outt, " tok je %le\n",jh6);
01648 break;
01649 }
01650 default:{
01651 print_err("unknown definition of Get_moisture_flux-KIND6 is required",__FILE__,__LINE__,__func__);
01652 }
01653 }
01654 }
01655
01656 return(jh4 +jh5 +jh6);
01657 }
01658
01659
01660
01661
01662
01663
01664
01665
01666
01667
01668
01669
01670 double climatcond::condit_water_contact ()
01671 {
01672 return 0.0;
01673 }
01674
01675
01676
01677
01678
01679
01680
01681
01682
01683
01684 double climatcond::condit_rain_moisture_flux (long rainkind, double nval_0,long nid)
01685 {
01686 double NorRain, MoistureAkt, Smc;
01687 double kwind, betawind, krainEff, jkNorRain, jkmaxWat;
01688 double c1, c2, c3, c4, c5;
01689
01690
01691
01692 if (WCorPH == 1)
01693 abort();
01694
01695
01696 switch (mattyp){
01697 case 150:{
01698 MoistureAkt = nval_0;
01699 Smc = 10000;
01700 break;
01701 }
01702 case 154:{
01703 MoistureAkt = Tt->nodes[nid].eqother[1];
01704 Smc = Tt->nodes[nid].eqother[3];
01705 break;
01706 }
01707 case 155:{
01708 MoistureAkt = Tt->nodes[nid].eqother[0];
01709 Smc = Tt->nodes[nid].eqother[2];
01710
01711
01712 break;
01713 }
01714 case 158:{
01715 MoistureAkt = Tt->nodes[nid].eqother[0];
01716 Smc = Tt->nodes[nid].eqother[3];
01717 break;
01718 }
01719 case 203:{
01720 MoistureAkt = nval_0;
01721 Smc = Tt->nodes[nid].eqother[3];
01722 break;
01723 }
01724 case 202:{
01725 break;
01726 }
01727 case 157:{
01728 MoistureAkt = nval_0;
01729 break;
01730 }
01731 case 159:{
01732 MoistureAkt = nval_0;
01733 Smc = Tt->nodes[nid].eqother[3];
01734 break;
01735 }
01736 case 156:{
01737 MoistureAkt = nval_0;
01738 Smc = Tt->nodes[nid].eqother[2];
01739 break;
01740 }
01741 default:{
01742 print_err("unknown definition of Condit_Rain_moistureFlux is required",__FILE__,__LINE__,__func__);
01743 }
01744 }
01745
01746
01747
01748
01749
01750 if (Smc == 0.0)
01751 return 0.0;
01752
01753
01754
01755
01756
01757
01758 switch (rainkind){
01759 case 51:{
01760 kwind = 1;
01761 NorRain = VRVal;
01762 break;
01763 }
01764 case 52:{
01765
01766
01767
01768 if ((fabs(basicdirect - WDVal)) <= M_PI){
01769 betawind = fabs(basicdirect - WDVal);
01770 }else{
01771 betawind = 2*M_PI - fabs(basicdirect - WDVal);
01772 }
01773
01774 if ((betawind <(M_PI/2)) && (WSVal >0) && (VRVal >0)){
01775
01776 c1 = cos (betawind);
01777 c2 = 1141*sqrt(3600*VRVal/(WSVal *WSVal*WSVal*WSVal));
01778 c3 = sqrt (1 + c2);
01779 c5 = sqrt(3600*VRVal);
01780 c5 = sqrt(c5);
01781 c4 = exp(-12/(5*c5));
01782 kwind = c1/c3 * c4;
01783
01784 }else{
01785 kwind = 0;
01786 }
01787
01788 NorRain = VRVal ;
01789
01790 break;
01791 }
01792 default:{
01793 print_err("unknown definition of Condit_Rain_moistureFlux is required",__FILE__,__LINE__,__func__);
01794 }
01795 }
01796
01797
01798 if (NorRain <= MinRainFlux){
01799 krainEff = 1;
01800 }else{
01801 krainEff = RainCoeff ;
01802 }
01803
01804 if (TepVal <= MinRainTemp){
01805 jkNorRain = 0;
01806 }else{
01807 jkNorRain = kwind * krainEff * NorRain;
01808 }
01809
01810 jkmaxWat = ExchCoeff * (Smc - MoistureAkt);
01811
01812 if (jkNorRain <jkmaxWat){
01813 return jkNorRain;
01814 }else{
01815 return jkmaxWat;
01816 }
01817
01818 }
01819
01820
01821
01822
01823
01824
01825
01826
01827
01828
01829
01830
01831
01832
01833 double climatcond::condit_vapour_diffusion_moisture_flux (long VDKind,double nval_0, double nval_1,long ipp, long nid)
01834 {
01835
01836 long idm;
01837 double VapPres,VapPresF,MoistureAkt, RHAkt;
01838 double betaWind, psate, psati, betaExch, pv, smc;
01839 double rov,M,R;
01840
01841
01842 M = 0.01801528;
01843
01844 R = 8.314472;
01845
01846
01847
01848 smc = MatConstC[1];
01849
01850
01851
01852 switch (mattyp){
01853 case 150:{
01854
01855
01856
01857
01858 MoistureAkt = nval_0;
01859
01860 idm = Tm->ip[ipp].idm;
01861 RHAkt = Tm->bazped[idm].inverse_sorption_isotherm (MoistureAkt);
01862
01863
01864 psati = exp(23.5771 - 4042.9/(nval_1 - 37.58));
01865 pv = psati * RHAkt;
01866 break;
01867 }
01868
01869 case 154:{
01870
01871
01872 MoistureAkt = Tt->nodes[nid].eqother[1];
01873 smc = Tt->nodes[nid].eqother[3];
01874 RHAkt = Tt->nodes[nid].eqother[2];
01875
01876
01877 psati = exp(23.5771 - 4042.9/(nval_1 - 37.58));
01878 pv = psati * RHAkt ;
01879 break;
01880 }
01881
01882 case 155:{
01883
01884
01885 MoistureAkt = Tt->nodes[nid].eqother[0];
01886 smc = Tt->nodes[nid].eqother[2];
01887 switch (int(Tt->nodes[nid].eqother[4])){
01888 case 1:{
01889 RHAkt = nval_0;
01890
01891
01892 psati = exp(23.5771 - 4042.9/(nval_1 - 37.58));
01893 pv = psati * RHAkt ;
01894 break;
01895 }
01896
01897 case 2:
01898 case 3:
01899 case 4:{
01900 pv = nval_0;
01901 break;
01902 }
01903 case 5:
01904 case 6:{
01905
01906 pv = nval_0;
01907 MoistureAkt = 0.0;
01908 break;
01909 }
01910 case 7:{
01911
01912
01913 rov = nval_0;
01914 pv = rov*R*nval_1/M;
01915 MoistureAkt = 0.0001;
01916 break;
01917 }
01918 case 8:{
01919
01920
01921
01922 pv = nval_0;
01923 MoistureAkt = Tt->nodes[nid].eqother[0];
01924 TepVal = 298.15;
01925 break;
01926 }
01927 default:{
01928
01929 }
01930 }
01931 break;
01932 }
01933
01934 case 158:{
01935
01936
01937 smc = Tt->nodes[nid].eqother[3];
01938 MoistureAkt = Tt->nodes[nid].eqother[0];
01939
01940
01941 psati = exp(23.5771 - 4042.9/(nval_1 - 37.58));
01942 pv = psati * RHAkt ;
01943 break;
01944 }
01945
01946 case 159:{
01947
01948
01949 MoistureAkt = nval_0;
01950 RHAkt = Tt->nodes[nid].eqother[1];
01951 smc = Tt->nodes[nid].eqother[3];
01952
01953
01954 psati = exp(23.5771 - 4042.9/(nval_1 - 37.58));
01955 pv = psati * RHAkt ;
01956 break;
01957 }
01958
01959 case 203:{
01960
01961
01962 MoistureAkt = nval_0;
01963 RHAkt = Tt->nodes[nid].eqother[1];
01964 smc = Tt->nodes[nid].eqother[3];
01965
01966
01967 psati = exp(23.5771 - 4042.9/(nval_1 - 37.58));
01968 pv = psati * RHAkt ;
01969 break;
01970 }
01971
01972 case 202:
01973 case 201:
01974 case 200:{
01975
01976
01977 MoistureAkt = nval_0;
01978 RHAkt = Tt->nodes[nid].eqother[0];
01979 smc = Tt->nodes[nid].eqother[2];
01980
01981 psati = exp(23.5771 - 4042.9/(273.15 - 37.58));
01982 pv = psati * RHAkt ;
01983 break;
01984 }
01985
01986 case 157:
01987 case 156:{
01988
01989
01990 MoistureAkt = nval_0;
01991 RHAkt = Tt->nodes[nid].eqother[0];
01992 smc = Tt->nodes[nid].eqother[2];
01993
01994 psati = exp(23.5771 - 4042.9/(nval_1 - 37.58));
01995 pv = psati * RHAkt ;
01996 break;
01997 }
01998 default:{
01999 print_err("unknown material type is required",__FILE__,__LINE__,__func__);
02000 }
02001 }
02002
02003
02004
02005
02006 if (smc == 0)
02007 return(0);
02008
02009
02010
02011
02012
02013
02014
02015
02016
02017
02018
02019
02020
02021
02022
02023
02024
02025 switch (VDKind){
02026 case 61:{
02027
02028 betaExch = ExchCoeffVD;
02029 break;
02030 }
02031 case 62:{
02032
02033
02034
02035
02036 WDVal = (M_PI/180.0)*WDVal;
02037 if (fabs (basicdirect - WDVal) <= M_PI){
02038 betaWind = fabs (basicdirect-WDVal);
02039 }
02040 else{
02041 betaWind = 2*M_PI - fabs(basicdirect-WDVal);
02042 }
02043
02044 if (betaWind < (M_PI/2)){
02045 betaExch = ExchCoeffVD +WindFactVD * sqrt(WSVal);
02046 }
02047 else{
02048 betaExch = ExchCoeffVD ;
02049 }
02050 break;
02051 }
02052 case 63:{
02053
02054 betaExch = ExchCoeffVD*MoistureAkt + alfaK;
02055 break;
02056 }
02057 case 64:{
02058
02059 betaExch = ExchCoeffVD*exp(alfaK*MoistureAkt);
02060 break;
02061 }
02062 case 65:{
02063
02064 betaExch = Ax*MoistureAkt*MoistureAkt + Bx*MoistureAkt + Cx;
02065 break;
02066 }
02067 case 66:{
02068
02069 betaExch = exchangecoeff.getval (MoistureAkt);
02070 break;
02071 }
02072 case 67:{
02073
02074
02075
02076 betaExch = ExchCoeffVD;
02077
02078 TepVal = rhtemp + 273.15;
02079
02080 RHVal=rhint/100;
02081 break;
02082 }
02083 default:{
02084 print_err("unknown type of VDKind is required",__FILE__,__LINE__,__func__);
02085 }
02086 }
02087
02088
02089
02090
02091
02092 psate = exp(23.5771 - 4042.9/(TepVal - 37.58));
02093
02094
02095 if (RHorVP == 0){
02096 VapPres = psate * RHVal;
02097 }else{
02098 VapPres = VapPresF;
02099 }
02100
02101
02102 double tokp;
02103
02104 if ((VapPres - pv) <0)
02105 {
02106
02107
02108 tokp = (betaExch *(pv - VapPres));
02109 }
02110 else
02111 {
02112
02113 tokp = (((smc-MoistureAkt)/smc)*betaExch *(pv -VapPres));
02114
02115
02116 }
02117
02118 if(mattyp==155)
02119 {
02120 if (int(Tt->nodes[nid].eqother[4])==7) {
02121 double rove;
02122 rove = VapPres*M/(R*TepVal);
02123
02124 }
02125 }
02126
02127
02128
02129 return(tokp);
02130 }
02131
02132
02133
02134
02135
02136
02137
02138
02139
02140
02141
02142
02143
02144
02145
02146
02147
02148
02149
02150
02151 double climatcond::give_salt_flux(double w, double cf, double temp, long nid)
02152 {
02153 double jh6salt;
02154
02155
02156
02157
02158
02159
02160
02161
02162
02163
02164
02165
02166
02167
02168 if (temp != -1000){
02169 switch(KIND6){
02170 case 0:{
02171 jh6salt = 0;
02172 break;
02173 }
02174 case 61:{
02175 jh6salt = condit_vapour_diffusion_moisture_flux_salt (61,w,cf,temp, nid);
02176 break;
02177 }
02178 case 62:{
02179 jh6salt = condit_vapour_diffusion_moisture_flux_salt (62,w,cf,temp, nid);
02180 break;
02181 }
02182 case 63:{
02183 jh6salt = condit_vapour_diffusion_moisture_flux_salt (63,w,cf,temp, nid);
02184 break;
02185 }
02186 case 64:{
02187 jh6salt = condit_vapour_diffusion_moisture_flux_salt (64,w,cf,temp, nid);
02188 break;
02189 }
02190 case 65:{
02191 jh6salt = condit_vapour_diffusion_moisture_flux_salt (65,w,cf,temp, nid);
02192 break;
02193 }
02194 case 66:{
02195 jh6salt = condit_vapour_diffusion_moisture_flux_salt (66,w,cf,temp, nid);
02196 break;
02197 }
02198 default:{
02199 print_err("unknown definition of Get_salt_flux-KIND6 is required",__FILE__,__LINE__,__func__);
02200 }
02201 }
02202 }
02203
02204
02205 return jh6salt;
02206 }
02207
02208
02209
02210
02211
02212
02213
02214
02215
02216 void climatcond::compute_values (double time)
02217 {
02218
02219 if (TepValue.tabf != NULL){
02220 if (TepValue.tabf->asize>0){
02221 TepVal = TepValue.getval (time);
02222 TepVal = TepVal + 273.15;
02223 }
02224 }
02225
02226 if (RHValue.tabf != NULL){
02227 if (RHValue.tabf->asize>0){
02228 RHVal = RHValue.getval (time);
02229 RHVal = RHVal/100.0;
02230 }
02231 }
02232
02233 if (DifSValue.tabf != NULL){
02234 if (DifSValue.tabf->asize>0){
02235 DifSVal = DifSValue.getval (time);
02236 }
02237 }
02238
02239 if (DirSValue.tabf != NULL){
02240 if (DirSValue.tabf->asize>0){
02241 DirSVal = DirSValue.getval (time);
02242 }
02243 }
02244
02245 if (WDValue.tabf != NULL){
02246 if (WDValue.tabf->asize>0){
02247 WDVal = WDValue.getval (time);
02248 WDVal = WDVal * M_PI / 180.0;
02249 }
02250 }
02251
02252 if (WSValue.tabf != NULL){
02253 if (WSValue.tabf->asize>0){
02254 WSVal = WSValue.getval (time);
02255 }
02256 }
02257
02258 if (VRValue.tabf != NULL){
02259 if (VRValue.tabf->asize>0){
02260 VRVal = VRValue.getval (time);
02261 VRVal = VRVal/3600.0;
02262 }
02263 }
02264
02265 if (LWValue.tabf != NULL){
02266 if (LWValue.tabf->asize>0){
02267 LWVal = LWValue.getval (time);
02268 }
02269 }
02270
02271 if (SLWValue.tabf != NULL){
02272 if (SLWValue.tabf->asize>0){
02273 SLWVal = SLWValue.getval (time);
02274 }
02275 }
02276
02277 if (NorRainValue.tabf != NULL){
02278 if (NorRainValue.tabf->asize>0){
02279 NorRainVal = NorRainValue.getval (time);
02280 }
02281 }
02282
02283 if (GasFluxValue.tabf != NULL){
02284 if (GasFluxValue.tabf->asize>0){
02285 GasFluxVal = GasFluxValue.getval (time);
02286 }
02287 }
02288
02289 if (HeatFluxValue.tabf != NULL){
02290 if (HeatFluxValue.tabf->asize>0){
02291 HeatFluxVal = HeatFluxValue.getval (time);
02292 }
02293 }
02294
02295 if (ShWRadValue.tabf != NULL){
02296 if (ShWRadValue.tabf->asize>0){
02297 ShWRadVal = ShWRadValue.getval (time);
02298 }
02299 }
02300
02301 if (LoWRadValue.tabf != NULL){
02302 if (LoWRadValue.tabf->asize>0){
02303 LoWRadVal = LoWRadValue.getval (time);
02304 }
02305 }
02306
02307 if (RadValue.tabf != NULL){
02308 if (RadValue.tabf->asize>0){
02309 RadVal = RadValue.getval (time);
02310 }
02311 }
02312
02313 if (TempbValue.tabf != NULL){
02314 if (TempbValue.tabf->asize>0){
02315 TempbVal = TempbValue.getval (time);
02316 }
02317 }
02318
02319 }
02320
02321
02322
02323
02324
02325
02326
02327
02328 void climatcond::prepare_values (double time)
02329 {
02330
02331 while (time >= 31536000){
02332 time = time - 31536000;
02333 }
02334
02335 if (fabs(time-last_time) < zero){
02336
02337
02338 }else{
02339 compute_values (time);
02340 }
02341
02342 }
02343
02344
02345
02346
02347
02348
02349
02350
02351 double climatcond::condit_vapour_diffusion_moisture_flux_salt (long VDKind,double w, double cf, double temp, long nid)
02352 {
02353 double VapPres, WindFactVD,VapPresF;
02354 double rh,betaWind, psat, betaExch, pv, smc;
02355
02356
02357 switch (mattyp){
02358 case 203:{
02359 rh = Tt->nodes[nid].eqother[1];
02360 smc = Tt->nodes[nid].eqother[3];
02361 cfmax = Tt->nodes[nid].eqother[4];
02362 break;
02363 }
02364 case 202:
02365 case 201:
02366 case 200:{
02367
02368 temp = 283.15;
02369 switch (MatCharC[0]){
02370 case 30:{
02371
02372
02373 break;}
02374 case 2:{
02375
02376 break;
02377 }
02378 }
02379
02380 break;
02381 }
02382 default:{
02383 print_err("unknown definition of Condit_Vapour_Diffusion_moisture_flux_salt is required",__FILE__,__LINE__,__func__);
02384 }
02385 }
02386
02387
02388 switch (VDKind){
02389 case 0:{
02390
02391 betaExch = ExchCoeffVD ;
02392 break;
02393 }
02394 case 61:{
02395
02396
02397
02398
02399
02400 if (fabs (basicdirect-WDVal) <= M_PI){
02401 betaWind = fabs (basicdirect-WDVal);
02402 }else{
02403 betaWind = 2*M_PI - fabs(basicdirect-WDVal);
02404 }
02405
02406 if (betaWind < (M_PI/2)){
02407 betaExch = ExchCoeffVD +WindFactVD * sqrt(WSVal);
02408 }else{
02409 betaExch = ExchCoeffVD ;
02410 }
02411 break;
02412 }
02413 case 62:{
02414
02415 betaExch = ExchCoeffVD*w + alfaK;
02416 break;
02417 }
02418 case 63:{
02419
02420 betaExch = ExchCoeffVD*exp(alfaK*w);
02421 break;
02422 }
02423 case 64:{
02424
02425 betaExch = Ax*w*w + Bx*w + Cx;
02426 break;
02427 }
02428 case 65:{
02429
02430 betaExch = exchangecoeff.getval (w);
02431 break;
02432 }
02433 case 67:{
02434
02435
02436
02437 betaExch = ExchCoeffVD;
02438
02439 TepVal = rhtemp;
02440
02441 RHVal=rhint/100;
02442 break;
02443 }
02444 default:{
02445 print_err("unknown definition of Condit_Vapour_Diffusion_moisture_flux_salt-VDKind is required",__FILE__,__LINE__,__func__);
02446 }
02447 }
02448
02449
02450
02451
02452 psat = exp(23.5771 - 4042.9/(TepVal - 37.58));
02453
02454 if (RHorVP == 0){
02455 VapPres = psat * RHVal;
02456 }else{
02457 VapPres = VapPresF;
02458 }
02459
02460 psat = exp(23.5771 - 4042.9/(temp - 37.58));
02461 pv = psat * rh ;
02462
02463
02464
02465
02466 double tok;
02467 if (cf >= cfmax){
02468
02469 if ((VapPres - pv) <0){
02470 if(Tt->nodes[nid].eqother[11] == 1){
02471 tok = betaExch *(pv - VapPres)*cf/1000;
02472 Tt->nodes[nid].eqother[11]= Tt->nodes[nid].eqother[11] + 1;
02473 Tt->nodes[nid].eqother[12]= Tt->nodes[nid].eqother[12] + tok;
02474 Tt->nodes[nid].eqother[13]= tok;
02475 return(tok);
02476 }else{
02477 Tt->nodes[nid].eqother[11]= Tt->nodes[nid].eqother[11]+1;
02478 if(Tt->nodes[nid].eqother[11] ==5)
02479 Tt->nodes[nid].eqother[11] = 1;
02480 return (Tt->nodes[nid].eqother[13]);
02481 }
02482 }else{
02483 if(Tt->nodes[nid].eqother[11] == 1){
02484 tok = ((smc -w)/smc)*betaExch *(pv -VapPres)*cf/1000;
02485 Tt->nodes[nid].eqother[11]= Tt->nodes[nid].eqother[11]+1;
02486 Tt->nodes[nid].eqother[12]= Tt->nodes[nid].eqother[12]+ tok;
02487 Tt->nodes[nid].eqother[13]= tok;
02488 return(tok);
02489 }else{
02490 Tt->nodes[nid].eqother[11]= Tt->nodes[nid].eqother[11]+1;
02491 if(Tt->nodes[nid].eqother[11] ==5)
02492 Tt->nodes[nid].eqother[11] = 1;
02493 return(Tt->nodes[nid].eqother[13]);
02494 }
02495 }
02496 }else{
02497 if( Tt->nodes[nid].eqother[12] > 1.0e-20){
02498 if ((VapPres - pv) <0){
02499 if(Tt->nodes[nid].eqother[11] == 1){
02500 tok =betaExch *(pv - VapPres)*cf/1000;
02501 if (fabs(tok) > Tt->nodes[nid].eqother[12])
02502 tok = Tt->nodes[nid].eqother[12];
02503 Tt->nodes[nid].eqother[11]= Tt->nodes[nid].eqother[11]+1;
02504 Tt->nodes[nid].eqother[12]= Tt->nodes[nid].eqother[12]-tok;
02505 Tt->nodes[nid].eqother[13]= tok;
02506 return(tok);
02507 }else{
02508 Tt->nodes[nid].eqother[11]= Tt->nodes[nid].eqother[11]+1;
02509 if(Tt->nodes[nid].eqother[11] ==5)
02510 Tt->nodes[nid].eqother[11] = 1;
02511 return(Tt->nodes[nid].eqother[13]);
02512 }
02513 }else{
02514 if(Tt->nodes[nid].eqother[11] == 1){
02515 tok = ((smc -w)/smc)*betaExch *(pv -VapPres)*cf/1000;
02516 if (fabs(tok) > Tt->nodes[nid].eqother[12])
02517 tok = Tt->nodes[nid].eqother[12];
02518 Tt->nodes[nid].eqother[11]= Tt->nodes[nid].eqother[11]+1;
02519 Tt->nodes[nid].eqother[12]= Tt->nodes[nid].eqother[12]-tok;
02520 Tt->nodes[nid].eqother[13]= tok;
02521 return(tok);
02522 }else{
02523 Tt->nodes[nid].eqother[11]= Tt->nodes[nid].eqother[11]+1;
02524 if(Tt->nodes[nid].eqother[11] ==5)
02525 Tt->nodes[nid].eqother[11] = 1;
02526 return(Tt->nodes[nid].eqother[13]);
02527 }
02528 }
02529 }else{
02530 Tt->nodes[nid].eqother[11] = 1;
02531 Tt->nodes[nid].eqother[12] = 0.0;
02532 Tt->nodes[nid].eqother[13] = 0.0;
02533 return(0);
02534 }
02535 }
02536 }
02537
02538