34 #define RC_INI(rs) (rs[_r].rc>1 ? DEC_RC_(rs) : (rs[_r].rc==1 ? INC_NDX_(rs) : rs[_r].ini ))
35 #define DEC_RC_(rs) (rs[_r].rc--,rs[_r].ini)
36 #define INC_NDX_(rs) (_r++,rs[_r-1].ini)
62 for( i=0; i <
LIMELM; i++ )
80 for( i=0; i < iup; i++ )
89 fprintf(
ioQQQ,
" 3BOD rate:" );
90 for( i=1; i <= 14; i++ )
94 fprintf(
ioQQQ,
"\n" );
100 fprintf(
punch.
ioRecom,
" 3-body rec coef vs charge \n" );
101 for( i=0; i < iup; i++ )
143 double alogte , alogne;
178 alogt = alogte - 2.*zlog;
179 alogn = alogne - 7.*zlog;
183 alogt = alogte -
zlog2[nz-1];
184 alogn = alogne -
zlog7[nz-1];
204 if( alogt <= 2.1760913 )
206 if( alogn < (3.5*alogt - 8.) )
210 else if( alogn > (3.5*alogt - 2.) )
213 alogn = 3.5*alogt - 2.;
217 else if( alogt <= 2.4771213 )
226 else if( alogt <= 5.1139434 )
244 nt = (long)(9.9657843*alogt + 1.);
248 nt = (long)(19.931568*alogt - 19.);
254 if( fabs(alogt-
tz[nt-1]) >= fabs(alogt-
tz[
MIN2(83,nt+1)-1]) )
258 else if( fabs(alogt-
tz[nt-1]) > fabs(alogt-
tz[
MAX2(1,nt-1)-1]) )
264 if( fabs(alogt-
tz[nt-1]) < 0.00001 )
283 xnc =
xinvrs(alogn,c,d,u,v,&jfail);
284 if( xnc <= 0. || jfail != 0 )
307 if( (nt <= 21) && (z == 1.0) )
321 if( alogt <= 2.1760913 )
332 if( alogt <= 2.4771213 )
353 yyb[0] = 6.93410742e03;
359 yyb[0] = 1.36653821e03;
371 c =
xmap(xx,yya,alogt);
374 d =
xmap(xx,yyb,alogt);
377 u =
xmap(xx,yyx,alogt);
386 yyb[0] = 2.25774918e01;
392 yyb[0] = 6.70408707e02;
397 yya[0] =
a2[nt0-(21)];
398 yyb[0] =
b2[nt0-(21)];
399 yyx[0] =
x2[nt0-(21)];
402 yya[1] =
a2[nt-(21)];
403 yya[2] =
a2[nt1-(21)];
404 c =
xmap(xx,yya,alogt);
405 yyb[1] =
b2[nt-(21)];
406 yyb[2] =
b2[nt1-(21)];
407 d =
xmap(xx,yyb,alogt);
408 yyx[1] =
x2[nt-(21)];
409 yyx[2] =
x2[nt1-(21)];
410 u =
xmap(xx,yyx,alogt);
414 xnc =
xinvrs(alogn,c,d,u,v,&jfail);
415 if( xnc <= 0. || jfail != 0 )
426 yya[0] = -8.04963875;
427 yyb[0] = 1.07205127e03;
432 yya[0] = -8.54721069;
433 yyb[0] = 4.70450195e02;
445 a =
xmap(xx,yya,alogt);
448 b =
xmap(xx,yyb,alogt);
452 y =
xmap(xx,yyy,alogt);
455 expp = a - y*alognc + b*pow(xnc,x);
458 da_v = z*pow(10.,expp);
488 {
static struct{
long rc;
double ini; } _rs0[] = {
519 for(_i=_r=0L; _i < 28; _i++)
521 {
static struct{
long rc;
double ini; } _rs1[] = {
552 for(_i=_r=0L; _i < 28; _i++)
554 {
static struct{
long rc;
double ini; } _rs2[] = {
640 for(_i=_r=0L; _i < 83; _i++)
642 {
static struct{
long rc;
double ini; } _rs3[] = {
728 for(_i=_r=0L; _i < 83; _i++)
730 {
static struct{
long rc;
double ini; } _rs4[] = {
812 for(_i=_r=0L; _i < 79; _i++)
814 {
static struct{
long rc;
double ini; } _rs5[] = {
821 for(_i=_r=0L; _i < 4; _i++)
823 {
static struct{
long rc;
double ini; } _rs6[] = {
909 for(_i=_r=0L; _i < 83; _i++)
912 {
static struct{
long rc;
double ini; } _rs7[] = {
998 for(_i=_r=0L; _i < 83; _i++)
1000 {
static struct{
long rc;
double ini; } _rs8[] = {
1082 for(_i=_r=0L; _i < 79; _i++)
1084 {
static struct{
long rc;
double ini; } _rs9[] = {
1091 for(_i=_r=0L; _i < 4; _i++)
1093 {
static struct{
long rc;
double ini; } _rs10[] = {
1179 for(_i=_r=0L; _i < 83; _i++)
1181 {
static struct{
long rc;
double ini; } _rs11[] = {
1247 for(_i=_r=0L; _i < 63; _i++)
1249 {
static struct{
long rc;
double ini; } _rs12[] = {
1315 for(_i=_r=0L; _i < 63; _i++)
1317 {
static struct{
long rc;
double ini; } _rs13[] = {
1383 for(_i=_r=0L; _i < 63; _i++)
1418 x12m = xmapx1 - xmapx2;
1419 y13m = yinit - y[2];
1420 x3 = (xmapx1 + x3)*x13m;
1421 xmapx2 = (xmapx1 + xmapx2)*x12m;
1422 b = ((yinit - y[1])*x3 - y13m*xmapx2)/(x12m*x3 - x13m*xmapx2);
1423 a = (y13m - x13m*b)/x3;
1424 c = yinit - a*xmapx1*xmapx1 - b*xmapx1;
1426 xmap_v = a*xmapx0*xmapx0 + b*xmapx0 + c;
1449 static long itmax = 10;
1462 for( i=0; i < itmax; i++ )
1465 fx = y - a - bxu + v*xlog;
1466 dfx = v*.4342945 - bxu*u;
1470 fxdfx = fabs(fx/dfx);
1471 fxdfx =
MIN2(0.2,fxdfx);
1472 xx = x*(1. -
sign(fxdfx,fx/dfx));
1478 xx = x*(1. -
sign(0.2,fx));
1481 if( (fabs(xx-x)/x) < 1.e-4 )