00001
00002
00003
00004
00005
00006
00007
00008
00009
00010 #include <stdio.h>
00011 #include <stdlib.h>
00012 #include <math.h>
00013 #include "constrel.h"
00014 #include "C60baroghel.h"
00015 #include "globalt.h"
00016
00017 C60barmat::C60barmat()
00018 {
00019 mw = 18.01528;
00020 ma = 28.9645;
00021 gasr = 8314.41;
00022
00023 t0 = 273.15;
00024 p0 = 101325.0;
00025 tcr = 647.3;
00026
00027 w1 = 163.0;
00028 c1 = 450.0;
00029
00030
00031 phi0 = 0.0825468;
00032 aphi = 0.00017417;
00033
00034
00035 k0 = 2.0e-18;
00036 ak = 6.29e-3;
00037
00038
00039 scr = 1.0;
00040 ag = 1.0;
00041
00042
00043 sir = 0.2;
00044 aw = 2.0;
00045 bw = 6.0;
00046
00047
00048 rhos = 2564.0;
00049
00050 fs = 1.0;
00051
00052
00053 lambdas0 = 1.9215152;
00054 alam = -0.0012525253;
00055
00056
00057 ac = -0.22626263;
00058 cps0 = 855.25758;
00059
00060
00061 hydren = 0.5e+6;
00062
00063 finv = 0.65;
00064
00065 fste = 0.36;
00066
00067
00068 dld = 8.0E-2;
00069
00070
00071 emod0 = 3.452222e+10;
00072
00073
00074 vcoeff = 0.18;
00075
00076
00077 betas = 3.0*10.2e-6;
00078
00079
00080 alpha = 0.5;
00081
00082
00083
00084
00085 at = 1.0;
00086 bt = 15000.0;
00087
00088 acc = 1.876;
00089
00090 ddbw0 = 1.0e-20;
00091 }
00092 C60barmat::~C60barmat()
00093 {}
00094
00095
00096
00097
00098
00099
00100
00101
00102 double C60barmat::sat(double pc,double t)
00103 {
00104 double sw,n,a,b,tt,q0,q2,q3,z,e,e0,g;
00105
00106 b = 2.27;
00107 q3 = 18.62e6 ;
00108 q2 = 7.0e6;
00109 n = 1.2;
00110 z = 0.5;
00111
00112 if (t <= 373.15)
00113 a = q3;
00114 else{
00115 tt=(t-373.15)/(647.15-373.15);
00116 q0=(q3-q2)*((2.*(tt*tt*tt))-(3.0*(tt*tt))+1.0);
00117 q2 = 25.0e6;
00118 a = q0 + q2;
00119 }
00120
00121 if (t <= tcr)
00122 e = pow(((tcr - t0)/(tcr - t)),n);
00123 else{
00124 e0=pow(((tcr-293.15)/z),n);
00125 e = n/z*e0*t + e0 - n/z*e0*(tcr - z);
00126 }
00127
00128 g = pow((e/a*pc),(b/(b-1.0)));
00129 sw = pow((g + 1.0),(-1.0/b));
00130
00131
00132
00133 return(sw);
00134 }
00135
00136
00137
00138
00139
00140
00141
00142
00143 double C60barmat::dsat_dpc(double pc,double t)
00144 {
00145 double dsw_dpc;
00146 double n,a,b,dg_dpc,tt,q0,q2,q3,z,e,e0,g;
00147
00148 b = 2.27;
00149 q3 = 18.62e6 ;
00150 q2 = 7.0e6;
00151 n = 1.2;
00152 z = 0.5;
00153
00154 if (t <= 373.15)
00155 a = q3;
00156 else{
00157 tt=(t-373.15)/(647.15-373.15);
00158 q0=(q3-q2)*((2.*(tt*tt*tt))-(3.0*(tt*tt))+1.0);
00159 q2 = 25.0e6;
00160 a = q0 + q2;
00161 }
00162
00163 if (t <= tcr){
00164 e = pow(((tcr - t0)/(tcr - t)),n);
00165 }
00166 else{
00167 e0=pow(((tcr-293.15)/z),n);
00168 e = n/z*e0*t + e0 - n/z*e0*(tcr - z);
00169 }
00170 g = pow((e/a*pc),(b/(b-1.0)));
00171 dg_dpc = (b/(b-1.0))*e/a*pow((e/a*pc),(b/(b-1.0)-1.0));
00172 dsw_dpc = (-1.0/b)*pow((g + 1.0),(-1.0-1.0/b))*dg_dpc;
00173
00174
00175
00176 return(dsw_dpc);
00177 }
00178
00179
00180
00181
00182
00183
00184
00185
00186 double C60barmat::dsat_dt(double pc,double t)
00187 {
00188 double dsw_dt,e,b,n,de_dt,g,a,da_dt,dea_dt,dg_dt,tt,dtt_dt,q0,dq0_dt,q2,q3,e0,z;
00189
00190 b = 2.27;
00191 q3 = 18.62e6 ;
00192 q2 = 7.0e6;
00193 n = 1.2;
00194 z = 0.5;
00195
00196 if (t <= 373.15){
00197 a = q3;
00198 da_dt = 0.0;
00199 }
00200 else{
00201 tt=(t-373.15)/(647.15-373.15);
00202 q0=(q3-q2)*((2.*(tt*tt*tt))-(3.0*(tt*tt))+1.0);
00203 dtt_dt = 1.0/(647.15-373.15);
00204 dq0_dt = (q3-q2)*(6.0*tt*tt*dtt_dt - 6.0*tt*dtt_dt);
00205 q2 = 25.0e6;
00206 a = q0 + q2;
00207 da_dt = dq0_dt;
00208 }
00209
00210 if (t <= tcr){
00211 e = pow(((tcr - t0)/(tcr - t)),n);
00212 de_dt = n*pow(((tcr - t0)/(tcr - t)),(n-1.0))*(tcr - t0)/(tcr - t)/(tcr - t);
00213 }
00214 else{
00215 e0=pow(((tcr-293.15)/z),n);
00216 e = n/z*e0*t + e0 - n/z*e0*(tcr - z);
00217 de_dt = n/z*e0;
00218 }
00219
00220 g = pow((e/a*pc),(b/(b-1.0)));
00221 dea_dt = (de_dt*a - da_dt*e)/a/a;
00222 dg_dt = (b/(b-1.0))*pc*pow((e/a*pc),(b/(b-1.0)-1.0))*dea_dt;
00223 dsw_dt = (-1.0/b)*pow((g + 1.0),(-1.0-1.0/b))*dg_dt;
00224
00225
00226
00227 return(dsw_dt);
00228 }
00229
00230
00231
00232
00233
00234
00235 double C60barmat::ssp()
00236 {
00237 return(0.55);
00238 }
00239
00240
00241
00242
00243
00244
00245
00246
00247 double C60barmat::C60bar_phi(double t)
00248 {
00249 double phi;
00250
00251 phi = phi0 + aphi*(t - t0);
00252
00253 return(phi);
00254 }
00255
00256
00257
00258
00259
00260
00261
00262
00263
00264
00265 double C60barmat::C60bar_kintr(double pc,double pg,double t,double dam)
00266 {
00267 double kintr,bk;
00268 state_eq tt;
00269
00270 bk = log(5.0/3.0)/log(4.0);
00271
00272 kintr = k0*exp(log(10.0)*ak*(t-t0))*pow((pg/p0),bk);
00273
00274
00275
00276
00277
00278
00279
00280
00281
00282
00283
00284
00285 return(kintr);
00286 }
00287
00288
00289
00290
00291
00292
00293
00294
00295 double C60barmat::C60bar_krg(double pc,double t)
00296 {
00297 double krg,s;
00298
00299 s = sat(pc,t);
00300
00301 krg = 1.0 - pow((s/scr),ag);
00302
00303 return(krg);
00304 }
00305
00306
00307
00308
00309
00310
00311
00312
00313
00314
00315 double C60barmat::C60bar_krw(double pc,double t,double rh)
00316 {
00317 double krw,help,s;
00318
00319 s = sat(pc,t);
00320
00321
00322 if (rh < 0.75)
00323 krw = pow(((s-sir)/(1.0-sir)),aw);
00324 else{
00325 help = 1.0 + pow(((1.0-rh)/0.25),bw);
00326 krw = pow(help,-1.0)*pow(s,aw);
00327 }
00328
00329 return(krw);
00330 }
00331
00332
00333
00334
00335
00336
00337
00338
00339 double C60barmat::C60bar_dd(double pc,double t)
00340 {
00341 double tau,dd;
00342
00343 tau = C60bar_tau(pc,t);
00344
00345 dd = tau*mw*ma/gasr;
00346
00347
00348 return(dd);
00349 }
00350
00351
00352
00353
00354
00355
00356
00357
00358
00359 double C60barmat::C60bar_deff(double pc,double pg,double t)
00360 {
00361
00362 double deff,dd,phi,s,cdiff;
00363 state_eq tt;
00364
00365 dd = C60bar_dd(pc,t);
00366 phi = C60bar_phi(t);
00367 s = sat(pc,t);
00368 cdiff = tt.get_cdiff(pc,pg,t);
00369
00370 deff = dd*phi*(1.0 - s)*fs*cdiff;
00371
00372 return(deff);
00373 }
00374
00375
00376
00377
00378
00379
00380
00381 double C60barmat::C60bar_cps(double t)
00382 {
00383 double cps;
00384
00385 if ((t-t0) < 600.0)
00386 cps = cps0 + ac*(t - t0);
00387 else
00388 cps = cps0 + ac*(600.0);
00389
00390 return(cps);
00391 }
00392
00393
00394
00395
00396
00397
00398
00399
00400
00401 double C60barmat::C60bar_rhocp(double pc,double pg,double t)
00402 {
00403 double s,phi,rhocp,rhow,rhog,rhogw,cps,cpw,cpga,cpgw,rhos;
00404 state_eq tt;
00405
00406 s = sat(pc,t);
00407 phi = C60bar_phi(t);
00408 rhow = tt.get_rhow(t);
00409 rhog = tt.get_rhog(pc,pg,t);
00410 rhogw = tt.get_rhogw(pc,t);
00411 cps = C60bar_cps(t);
00412 cpw = tt.get_cpw();
00413 cpga = tt.get_cpga();
00414 cpgw = tt.get_cpgw();
00415 rhos = C60bar_rhos();
00416
00417 rhocp = (1.0-phi)*rhos*cps + phi*(s*rhow*cpw + (1.0-s)*(rhog*cpga + rhogw*(cpgw-cpga)));
00418
00419 return(rhocp);
00420 }
00421
00422
00423
00424
00425
00426
00427
00428
00429 double C60barmat::C60bar_tau(double pc,double t)
00430 {
00431 double phi,s,tau;
00432
00433 s = sat(pc,t);
00434 phi = C60bar_phi(t);
00435
00436 tau = pow(phi,(1.0/3.0))*pow((1.0-s),(7.0/3.0));
00437
00438 return(tau);
00439 }
00440
00441
00442
00443
00444
00445
00446
00447
00448
00449
00450 double C60barmat::C60bar_lambdas(double t)
00451 {
00452 double lambdas;
00453
00454 if ((t-t0) < 600.0)
00455 lambdas = lambdas0 + alam*(t-t0);
00456 else
00457 lambdas = lambdas0 + alam*(600.0);
00458
00459 return(lambdas);
00460 }
00461
00462
00463
00464
00465
00466
00467
00468
00469 double C60barmat::C60bar_lambdaeff(double pc,double pg,double t)
00470 {
00471 double lambdaeff,lambdas,s,phi,rhow,rhos;
00472 state_eq tt;
00473
00474 s = sat(pc,t);
00475 phi = C60bar_phi(t);
00476 rhow = tt.get_rhow(t);
00477 lambdas = C60bar_lambdas(t);
00478 rhos = C60bar_rhos();
00479
00480 lambdaeff = lambdas*(1.0 + 4.0*phi*rhow*s/(1.0-phi)/rhos);
00481
00482 return(lambdaeff);
00483 }
00484
00485
00486
00487
00488
00489
00490 double C60barmat::C60bar_rhos()
00491 {
00492 return(rhos);
00493 }
00494
00495
00496
00497
00498
00499
00500 double C60barmat::C60bar_betas()
00501 {
00502 return(betas);
00503 }
00504
00505
00506
00507
00508
00509
00510
00511
00512
00513
00514
00515
00516
00517 double C60barmat::C60bar_hydw(double pc,double pg,double t)
00518 {
00519 double fhy,hydw;
00520
00521 if((t-t0)< 105.0){
00522 fhy = 0.0;
00523 }
00524 else{
00525 fhy = (1.0+sin(3.1416/2.0*(1.0-2.0*exp(-0.004*((t-t0)-105.0)))))/2.0;
00526 }
00527
00528 hydw = fste*finv*c1*fhy;
00529
00530 return(hydw);
00531 }
00532
00533
00534
00535
00536
00537
00538
00539
00540
00541 double C60barmat::C60bar_dehydw_dt(double pc,double pg,double t)
00542 {
00543 double dfhyt,dehydw_dt;
00544
00545 if((t-t0)< 105.0){
00546 dfhyt = 0.0;
00547 }
00548 else{
00549 dfhyt = (3.1416*0.004/2.0)*cos(3.1416/2.0*(1.0-2.0*exp(-0.004*((t-t0)-105.0))))*exp(-0.004*((t-t0)-105.0));
00550 }
00551
00552 dehydw_dt = fste*finv*c1*dfhyt;
00553
00554 return(dehydw_dt);
00555 }
00556
00557
00558
00559
00560
00561
00562
00563
00564
00565
00566 double C60barmat::C60bar_hydren(double pc,double pg,double t)
00567 {
00568 return(hydren);
00569 }
00570
00571
00572
00573
00574
00575
00576
00577
00578
00579 double C60barmat::C60bar_fste(double pc,double pg,double t)
00580 {
00581 return(fste);
00582 }
00583
00584
00585
00586
00587
00588
00589
00590
00591
00592 double C60barmat::C60bar_ddbw(double pc,double pg,double t)
00593 {
00594 double ddbw;
00595
00596 ddbw = ddbw0*exp(-t/(273.15+23));
00597 if (t > tcr)
00598 ddbw = ddbw0*exp(-tcr/(273.15+23));
00599
00600 ddbw=0.0;
00601
00602 return(ddbw);
00603 }
00604
00605
00606
00607
00608
00609
00610
00611
00612
00613
00614
00615
00616 double C60barmat::C60bar_emod(double pc,double pg,double t)
00617 {
00618 double emod,ttc;
00619
00620 ttc = t - 273.15;
00621
00622 if (ttc <= 500.0)
00623 emod = 3.5604e10 - 54.089*1.0e6*ttc;
00624 else
00625 emod = 0.8560e10*exp((-54.089/8560.0)*(ttc-500.0));
00626
00627 return(emod);
00628 }
00629
00630
00631
00632
00633
00634
00635
00636
00637
00638
00639 double C60barmat::C60bar_fct(double pc,double pg,double t)
00640 {
00641 double fct,ttc;
00642
00643 ttc = t - 273.15;
00644
00645 if (ttc <= 600.0)
00646 fct = (6.0 - 8.56e-3*ttc)*1.0e6;
00647 else
00648 fct = 0.864*1.0e6;
00649
00650 return(fct);
00651 }
00652
00653
00654
00655
00656
00657
00658
00659
00660
00661
00662
00663 double C60barmat::C60bar_xk0(double pc,double pg,double t)
00664 {
00665 double fct,emod,xk0;
00666
00667 fct = C60bar_fct(pc,pg,t);
00668 emod = C60bar_emod(pc,pg,t);
00669
00670 xk0 = fct/emod;
00671
00672 return(xk0);
00673 }
00674
00675
00676
00677
00678
00679
00680
00681
00682
00683
00684 double C60barmat::C60bar_bcc(double pc,double pg,double t)
00685 {
00686 double bcc,ttc;
00687
00688 ttc = t - 273.15;
00689
00690 bcc = 1550.0*exp(-0.0015154*ttc);
00691
00692 return(bcc);
00693 }
00694
00695
00696
00697
00698
00699
00700 double C60barmat::C60bar_alpha()
00701 {
00702 return(alpha);
00703 }
00704
00705
00706
00707
00708
00709
00710
00711 double C60barmat::C60bar_nu()
00712 {
00713 return(vcoeff);
00714 }
00715
00716
00717
00718
00719
00720
00721 void C60barmat::read(XFILE *in)
00722 {}
00723
00724
00725
00726
00727
00728
00729 void C60barmat::print(FILE *out)
00730 {}
00731
00732
00733
00734
00735
00736
00737
00738
00739
00740
00741
00742
00743
00744
00745 void C60barmat::give_reqntq(long *antq)
00746 {
00747
00748 antq[scal_iso_damage-1] = 1;
00749 }