57 if( strcmp(
rt.chEscFunSubord,
"SIMP") == 0 )
60 escinc_v = 1./(1. + tau);
73 escinc_v = 1./(1. + 1.6*tau);
81 b = 1.6 + (3.*pow(2.*a,-0.12))/(1. + atau);
85 b = 1.6 + (3.*pow(2.*a,-0.12))/(1. + 1./sqrt(atau));
89 escinc_v = 1./(1. + b*tau);
114 esccom_v += 0.333*sqrt(a/(
SQRTPI*tau));
203 *dest = (dstin + dstout)/2.f;
213 ASSERT( escla_v >=0. && *dest>=0. && *esin>=0. );
277 if( tau_out <0 || tau_in < 0. )
280 tt =
MIN2( tau_out , tau_in );
285 tt = tau_out - tau_in;
330 if( tau_out <0 || tau_in < 0. )
333 tt =
MIN2( tau_out , tau_in );
338 tt = tau_out - tau_in;
377 static double a[5]={1.00,-0.1117897,-0.1249099917,-9.136358767e-3,
379 static double b[6]={1.00,0.1566124168,9.013261660e-3,1.908481163e-4,
380 -1.547417750e-7,-6.657439727e-9};
381 static double c[5]={1.000,19.15049608,100.7986843,129.5307533,-31.43372468};
382 static double d[6]={1.00,19.68910391,110.2576321,169.4911399,-16.69969409,
400 else if( tau < 0.01 )
402 esca0k2_v = 1. - 2.*tau;
405 else if( tau <= 11. )
407 suma = a[0] + tau*(a[1] + tau*(a[2] + tau*(a[3] + a[4]*tau)));
408 sumb = b[0] + tau*(b[1] + tau*(b[2] + tau*(b[3] + tau*(b[4] +
410 esca0k2_v = tau/2.5066283*log(tau/SQRTPI) + suma/sumb;
416 arg = 1./log(tau/SQRTPI);
417 sumc = c[0] + arg*(c[1] + arg*(c[2] + arg*(c[3] + c[4]*arg)));
418 sumd = d[0] + arg*(d[1] + arg*(d[2] + arg*(d[3] + arg*(d[4] +
420 esca0k2_v = (sumc/sumd)/(2.*tau*sqrt(log(tau/SQRTPI)));
464 if(
HFLines[i].Emis->TauIn < -1. )
492 escmase_v = 1. - tau*(0.5 + tau/6.);
494 else if( tau > -30. )
496 escmase_v = (1. - exp(-tau))/tau;
500 fprintf(
ioQQQ,
" DISASTER escmase called with 2big tau%10.2e\n",
502 fprintf(
ioQQQ,
" This is zone number%4ld\n",
nzone );
508 ASSERT( escmase_v >= 1. );
551 remain = (1.998069357 + t*(66.4037741 + t*107.2041376))/(1. +
552 t*(37.4009646 + t*105.0388805));
557 remain = (1.823707708 + t*2.395042899)/(1. + t*(2.488885899 -
566 cone2_v = remain*tln/(2. + t);
599 else if( hnukt > 1. && tau > 100. )
614 escpcn_v = sumesc/sumrec;
649 esc0 = 1./((0.6451 + tau)*(0.47 + 1.08/(1. + 7.3e-6*tau)));
651 esc0 =
MAX2(0.,esc0);
652 esc0 =
MIN2(1.,esc0);
656 taulog = log10(
MIN2(1e8,tau));
670 taucon =
MIN2(2.,beta*tau);
674 fac1 = -1.25 + 0.475*taulog;
675 fac2 = -0.485 + 0.1615*taulog;
676 fac = -fac1*taucon + fac2*
POW2(taucon);
687 *dest = (
realnum)(beta/(0.30972 -
MIN2(.28972,0.03541667*taulog)));
696 *dest =
MIN2(*dest,1.f-*esc);
697 *dest =
MAX2(0.f,*dest);
709 fprintf(
ioQQQ,
"scatdest tau %.2e beta%.2e 1-al%.2e al%.2e dest%.2e \n",
777 if( abund <= 0. || conopc <= 0. )
785 beta = conopc/(abund*
SQRTPI*crsec/widl + conopc);
798 eovrlp_v =
MIN2(1e-3,8.5*beta);
804 eovrlp_v =
MIN2(1e-3,8.5*beta);
811 eovrlp_v =
MIN2(1e-3,8.5*beta);
815 fprintf(
ioQQQ,
" chCore of %i not understood by RT_DestProb.\n",
821 eovrlp_v /= 1. + eovrlp_v;
824 eovrlp_v *= 1. - escp;
833 enum {DEBUG_LOC=
false};
840 fprintf(
ioQQQ,
"%li RT_DestProb\t%g\n",
868 double RT_LineWidth_v,
930 else if( tau <= 20. )
935 aa = 4.8 + 5.2*tau + (4.*tau - 1.)*atau;
937 RT_LineWidth_v = vth*0.8862*aa/b*(1. -
MIN2( 1.,
940 if( Pesc < 10. * FLT_EPSILON )
946 atau = log(
MAX2(0.0001,tau));
947 aa = 1. + 2.*atau/pow(1. + 0.3*t->
Emis->
damp*tau,0.6667) + pow(6.5*
949 b = 1.6 + 1.5/(1. + 0.20*t->
Emis->
damp*tau);
950 RT_LineWidth_v = vth*0.8862*aa/b*(1. -
MIN2( 1. ,
955 RT_LineWidth_v *= 2.;
964 RT_LineWidth_v = vth*sqrt(log(
MAX2(1.,tau))*.2821);
969 if( r*vth <= RT_LineWidth_v )
971 RT_LineWidth_v = vth*r*log(RT_LineWidth_v/(r*vth));
976 ASSERT( RT_LineWidth_v >= 0. );
977 return( RT_LineWidth_v );
1006 fhummr_v = 3.8363 - 0.56329*x;
1010 fhummr_v = 2.79153 - 0.75325*x;
1014 fhummr_v = 1.8446 - 1.0238*x;
1018 fhummr_v = 0.72500 - 1.5836*x;