User: rhaas
Date: 2014/03/12 10:01 PM
Modified:
/trunk/src/nuc_eos_cxx/
helpers.hh, nuc_eos_dpdrhoe_dpderho.cc, nuc_eos_full.cc, nuc_eos_press.cc, nuc_eos_press_cs2.cc, nuc_eos_short.cc, readtable.cc
Log:
EOS_Omni: low-level optimizations
* use log-rules to transform base 10 logs to natural ones
faster and hopefully more accurate
* optimize away a division when computing cs2
* optimize bracketing test to (a-a1)*(a-a2)<0
* replace divisions by multiplication by inverse
* use CCTK_BUILTIN_EXPECT to indicate likely outcome
this should help branch prediction since the compiler can put the the
less likely branch into the unfavorable location.
* remove superfluous if statements
whenever we leave a itmax loop regularly we iterated for too long
* explicitly replace x/y by x*1./y
From: Roland Haas <rhaas(a)tapir.caltech.edu>
File Changes:
Directory: /trunk/src/nuc_eos_cxx/
==================================
File [modified]: helpers.hh
Delta lines: +49 -52
===================================================================
--- trunk/src/nuc_eos_cxx/helpers.hh 2014-03-13 03:01:52 UTC (rev 100)
+++ trunk/src/nuc_eos_cxx/helpers.hh 2014-03-13 03:01:55 UTC (rev 101)
@@ -18,24 +18,24 @@
// 105 -- rho too high
// 106 -- rho too low
- if(xrho > eos_rhomax) {
+ if(CCTK_BUILTIN_EXPECT(xrho > eos_rhomax, false)) {
return 105;
}
- if(xrho < eos_rhomin) {
+ if(CCTK_BUILTIN_EXPECT(xrho < eos_rhomin, false)) {
return 106;
}
- if(xye > eos_yemax) {
+ if(CCTK_BUILTIN_EXPECT(xye > eos_yemax, false)) {
return 101;
}
- if(xye < eos_yemin) {
+ if(CCTK_BUILTIN_EXPECT(xye < eos_yemin, false)) {
// this is probably not pure and should be removed
fprintf(stderr,"xye: %15.6E eos_yemin: %15.6E\n",xye,eos_yemin);
return 102;
}
- if(xtemp > eos_tempmax) {
+ if(CCTK_BUILTIN_EXPECT(xtemp > eos_tempmax, false)) {
return 103;
}
- if(xtemp < eos_tempmin) {
+ if(CCTK_BUILTIN_EXPECT(xtemp < eos_tempmin, false)) {
return 104;
}
return 0;
@@ -53,16 +53,16 @@
// 105 -- rho too high
// 106 -- rho too low
- if(xrho > eos_rhomax) {
+ if(CCTK_BUILTIN_EXPECT(xrho > eos_rhomax, false)) {
return 105;
}
- if(xrho < eos_rhomin) {
+ if(CCTK_BUILTIN_EXPECT(xrho < eos_rhomin, false)) {
return 106;
}
- if(xye > eos_yemax) {
+ if(CCTK_BUILTIN_EXPECT(xye > eos_yemax, false)) {
return 101;
}
- if(xye < eos_yemin) {
+ if(CCTK_BUILTIN_EXPECT(xye < eos_yemin, false)) {
return 102;
}
return 0;
@@ -331,10 +331,12 @@
// then interpolate in y
// assume rectangular grid
- double t1 = (fs[1]-fs[0])/(xs[1]-xs[0]) * (x - xs[0]) + fs[0];
- double t2 = (fs[3]-fs[2])/(xs[1]-xs[0]) * (x - xs[0]) + fs[2];
+ double dxi = 1./(xs[1]-xs[0]);
+ double dyi = 1./(ys[1]-ys[0]); // x*1./y uses faster instructions than x/y
+ double t1 = (fs[1]-fs[0])*dxi * (x - xs[0]) + fs[0];
+ double t2 = (fs[3]-fs[2])*dxi * (x - xs[0]) + fs[2];
- return (t2 - t1)/(ys[1]-ys[0]) * (y-ys[0]) + t1;
+ return (t2 - t1)*dyi * (y-ys[0]) + t1;
}
static inline __attribute__((always_inline))
@@ -359,6 +361,8 @@
const double dltp = log(1.2);
const double dltm = log(0.8);
+ const double leps0_prec = leps0*prec;
+
// temporary local vars
double lt, lt1, lt2;
double ltmin = logtemp[0];
@@ -404,7 +408,7 @@
#endif
bcount++;
- if(bcount >= maxbcount) {
+ if(CCTK_BUILTIN_EXPECT(bcount >= maxbcount, false)) {
*keyerrt = 667;
return;
}
@@ -432,22 +436,20 @@
nuc_eos_C_linterp_one(idx,delx,dely,delz,&f2a,iv);
fmid=f2a-leps0;
- if(fmid <= 0.0) lt=ltmid;
+ if(CCTK_BUILTIN_EXPECT(fmid <= 0.0, false)) lt=ltmid;
#if DEBUG
fprintf(stderr,"bisection step 2 it %d, fmid: %15.6E f2a: %15.6E lt: %15.6E ltmid: %15.6E dlt: %15.6E\n",
it,fmid,f2a,lt,ltmid,dlt);
#endif
- if(fabs(1.0-f2a/leps0) <= prec) {
+ if(CCTK_BUILTIN_EXPECT(fabs(leps0-f2a) <= leps0_prec, false)) {
*ltout = ltmid;
return;
}
} // for it = 0
- if(it >= itmax-1) {
- *keyerrt = 667;
- return;
- }
+ *keyerrt = 667;
+ return;
} // bisection
@@ -465,7 +467,7 @@
// local variables
const int itmax = 200; // use at most 10 iterations, then go to bisection
- double dlepsdlt; // derivative dlogeps/dlogT
+ double dlepsdlti; // 1 / derivative dlogeps/dlogT
double ldt;
double leps,leps0; // temp vars for eps
double ltn, lt; // temp vars for temperature
@@ -489,6 +491,7 @@
fprintf(stderr,"it: %d t: %15.6E leps: %15.6E eps0: %15.6E del: %15.6E\n",
it,lt, leps,leps0,fabs(leps-leps0)/(fabs(leps0)));
#endif
+ // TODO: profile this to see which outcome is more likely
if(fabs(leps-leps0) < prec*fabs(leps0)) {
*ltout = lt0;
return;
@@ -548,8 +551,7 @@
// Check if we are already bracketing the input internal
// energy. If so, interpolate for new T.
- if( (leps0 >= epst1 && leps0 <= epst2)
- || (leps0 <= epst1 && leps0 >=epst2)) {
+ if(CCTK_BUILTIN_EXPECT((leps0 - epst1) * (leps0 - epst2) <= 0., false)) {
*ltout = (logtemp[itemp]-logtemp[itemp-1]) / (epst2 - epst1) *
(leps0 - epst1) + logtemp[itemp-1];
@@ -559,8 +561,8 @@
// well, then do a Newton-Raphson step
// first, guess the derivative
- dlepsdlt = (epst2-epst1)/(logtemp[itemp]-logtemp[itemp-1]);
- ldt = -(leps - leps0) / dlepsdlt * fac;
+ dlepsdlti = (logtemp[itemp]-logtemp[itemp-1])/(epst2-epst1);
+ ldt = -(leps - leps0) * dlepsdlti * fac;
ltn = MIN(MAX(lt + ldt,ltmin),ltmax);
lt = ltn;
@@ -578,7 +580,7 @@
if(oerr < err) fac *= 0.9;
oerr = err;
- if(err < prec*fabs(leps0)) {
+ if(CCTK_BUILTIN_EXPECT(err < prec*fabs(leps0), false)) {
*ltout = lt;
return;
}
@@ -587,20 +589,18 @@
} // while(it < itmax)
- if(it >= itmax - 1) {
- // try bisection
- // bisection(lr, lt0, ye, leps0, ltout, 1, prec, keyerrt);
+ // try bisection
+ // bisection(lr, lt0, ye, leps0, ltout, 1, prec, keyerrt);
#if DEBUG
- fprintf(stderr, "Failed to converge. This is bad. Trying bisection!\n");
+ fprintf(stderr, "Failed to converge. This is bad. Trying bisection!\n");
#endif
- bisection(lr,lt0,ye,leps0,prec,ltout,1,keyerrt);
+ bisection(lr,lt0,ye,leps0,prec,ltout,1,keyerrt);
#if DEBUG
- if(*keyerrt==667) {
- fprintf(stderr,"This is worse. Bisection failed!\n");
- abort();
- }
+ if(*keyerrt==667) {
+ fprintf(stderr,"This is worse. Bisection failed!\n");
+ abort();
+ }
#endif
- }
return;
@@ -620,7 +620,7 @@
// local variables
const int itmax = 200; // use at most 10 iterations, then go to bisection
- double dentdlt; // derivative dentropy/dlogT
+ double dentdlti; // 1 / derivative dentropy/dlogT
double ldt;
double ent,ent0; // temp vars for eps
double ltn, lt; // temp vars for temperature
@@ -644,7 +644,7 @@
fprintf(stderr,"it: %d t: %15.6E leps: %15.6E eps0: %15.6E del: %15.6E\n",
it,lt,ent,ent0,fabs(ent-ent0)/(fabs(ent0)));
#endif
- if(fabs(ent-ent0) < prec*fabs(ent0)) {
+ if(CCTK_BUILTIN_EXPECT(fabs(ent-ent0) < prec*fabs(ent0), false)) {
*ltout = lt0;
return;
}
@@ -703,8 +703,7 @@
// Check if we are already bracketing the input internal
// energy. If so, interpolate for new T.
- if( (ent0 >= ent1 && ent0 <= ent2)
- || (ent0 <= ent1 && ent0 >= ent2)) {
+ if(CCTK_BUILTIN_EXPECT((ent0 - ent1) * (ent0 - ent2) <= 0., false)) {
*ltout = (logtemp[itemp]-logtemp[itemp-1]) / (ent2 - ent1) *
(ent0 - ent1) + logtemp[itemp-1];
@@ -714,8 +713,8 @@
// well, then do a Newton-Raphson step
// first, guess the derivative
- dentdlt = (ent2-ent1)/(logtemp[itemp]-logtemp[itemp-1]);
- ldt = -(ent - ent0) / dentdlt * fac;
+ dentdlti = (logtemp[itemp]-logtemp[itemp-1])/(ent2-ent1);
+ ldt = -(ent - ent0) * dentdlti * fac;
ltn = MIN(MAX(lt + ldt,ltmin),ltmax);
lt = ltn;
@@ -733,7 +732,7 @@
if(oerr < err) fac *= 0.9;
oerr = err;
- if(err < prec*fabs(ent0)) {
+ if(CCTK_BUILTIN_EXPECT(err < prec*fabs(ent0), false)) {
*ltout = lt;
return;
}
@@ -742,19 +741,17 @@
} // while(it < itmax)
- if(it >= itmax - 1) {
- // try bisection
+ // try bisection
#if DEBUG
- fprintf(stderr, "Failed to converge. This is bad. Trying bisection!\n");
+ fprintf(stderr, "Failed to converge. This is bad. Trying bisection!\n");
#endif
- bisection(lr,lt0,ye,ent0,prec,ltout,2,keyerrt);
+ bisection(lr,lt0,ye,ent0,prec,ltout,2,keyerrt);
#if DEBUG
- if(*keyerrt==667) {
- fprintf(stderr,"This is worse. Bisection failed!\n");
- abort();
- }
+ if(*keyerrt==667) {
+ fprintf(stderr,"This is worse. Bisection failed!\n");
+ abort();
+ }
#endif
- }
return;
File [modified]: nuc_eos_dpdrhoe_dpderho.cc
Delta lines: +3 -3
===================================================================
--- trunk/src/nuc_eos_cxx/nuc_eos_dpdrhoe_dpderho.cc 2014-03-13 03:01:52 UTC (rev 100)
+++ trunk/src/nuc_eos_cxx/nuc_eos_dpdrhoe_dpderho.cc 2014-03-13 03:01:55 UTC (rev 101)
@@ -30,7 +30,7 @@
// Note that this code now requires that the
// temperature guess be within the table bounds
keyerr[i] = checkbounds_kt0_noTcheck(rho[i], ye[i]);
- if(keyerr[i] != 0) {
+ if(CCTK_BUILTIN_EXPECT(keyerr[i] != 0, false)) {
*anyerr = 1;
}
}
@@ -45,7 +45,7 @@
const double lt = log(MIN(MAX(temp[i],eos_tempmin),eos_tempmax));
double ltout;
const double epstot = eps[i]+energy_shift;
- if(epstot>0.0e0) {
+ if(CCTK_BUILTIN_EXPECT(epstot>0.0e0, true)) {
// this is the standard scenario; eps is larger than zero
// and we can operate with logarithmic tables
const double lxeps = log(epstot);
@@ -61,7 +61,7 @@
keyerr[i] = 667;
} // epstot > 0.0
- if(keyerr[i]) {
+ if(CCTK_BUILTIN_EXPECT(keyerr[i] != 0, false)) {
// now try negative temperature treatment
double eps0, eps1;
int idx[8];
File [modified]: nuc_eos_full.cc
Delta lines: +6 -6
===================================================================
--- trunk/src/nuc_eos_cxx/nuc_eos_full.cc 2014-03-13 03:01:52 UTC (rev 100)
+++ trunk/src/nuc_eos_cxx/nuc_eos_full.cc 2014-03-13 03:01:55 UTC (rev 101)
@@ -39,7 +39,7 @@
for(int i=0;i<n;i++) {
// check if we are fine
keyerr[i] = checkbounds(rho[i], temp[i], ye[i]);
- if(keyerr[i] != 0) {
+ if(CCTK_BUILTIN_EXPECT(keyerr[i] != 0, false)) {
*anyerr = 1;
}
}
@@ -148,7 +148,7 @@
prs[i] = exp(prs[i]);
eps[i] = exp(eps[i]) - energy_shift;
#if HAVEGR
- cs2[i] = cs2[i] / (1.0 + eps[i] + prs[i]/rho[i]);
+ cs2[i] = rho[i] * cs2[i] / (rho[i] + rho[i] * eps[i] + prs[i]);
#endif
}
@@ -194,7 +194,7 @@
// Note that this code now requires that the
// temperature guess be within the table bounds
keyerr[i] = checkbounds_kt0_noTcheck(rho[i], ye[i]);
- if(keyerr[i] != 0) {
+ if(CCTK_BUILTIN_EXPECT(keyerr[i] != 0, false)) {
*anyerr = 1;
}
}
@@ -209,7 +209,7 @@
const double lt = log(MIN(MAX(temp[i],eos_tempmin),eos_tempmax));
double ltout;
const double epstot = eps[i]+energy_shift;
- if(epstot>0.0e0) {
+ if(CCTK_BUILTIN_EXPECT(epstot>0.0e0, true)) {
// this is the standard scenario; eps is larger than zero
// and we can operate with logarithmic tables
const double lxeps = log(epstot);
@@ -225,7 +225,7 @@
} // epstot > 0.0
- if(keyerr[i]) {
+ if(CCTK_BUILTIN_EXPECT(keyerr[i] != 0, false)) {
*anyerr=1;
} else {
temp[i] = exp(ltout);
@@ -318,7 +318,7 @@
for(int i=0;i<n;i++) {
prs[i] = exp(prs[i]);
#if HAVEGR
- cs2[i] = cs2[i] / (1.0 + eps[i] + prs[i]/rho[i]);
+ cs2[i] = rho[i]*cs2[i] / (rho[i] + rho[i]*eps[i] + prs[i]);
#endif
}
File [modified]: nuc_eos_press.cc
Delta lines: +4 -4
===================================================================
--- trunk/src/nuc_eos_cxx/nuc_eos_press.cc 2014-03-13 03:01:52 UTC (rev 100)
+++ trunk/src/nuc_eos_cxx/nuc_eos_press.cc 2014-03-13 03:01:55 UTC (rev 101)
@@ -25,7 +25,7 @@
for(int i=0;i<n;i++) {
// check if we are fine
keyerr[i] = checkbounds(rho[i], temp[i], ye[i]);
- if(keyerr[i] != 0) {
+ if(CCTK_BUILTIN_EXPECT(keyerr[i] != 0, false)) {
*anyerr = 1;
}
}
@@ -87,7 +87,7 @@
// Note that this code now requires that the
// temperature guess be within the table bounds
keyerr[i] = checkbounds_kt0_noTcheck(rho[i], ye[i]);
- if(keyerr[i] != 0) {
+ if(CCTK_BUILTIN_EXPECT(keyerr[i] != 0, false)) {
*anyerr = 1;
}
}
@@ -103,7 +103,7 @@
const double lt = log(MIN(MAX(temp[i],eos_tempmin),eos_tempmax));
double ltout;
const double epstot = eps[i]+energy_shift;
- if(epstot>0.0e0) {
+ if(CCTK_BUILTIN_EXPECT(epstot>0.0e0, true)) {
// this is the standard scenario; eps is larger than zero
// and we can operate with logarithmic tables
const double lxeps = log(epstot);
@@ -119,7 +119,7 @@
keyerr[i] = 667;
} // epstot > 0.0
- if(keyerr[i]) {
+ if(CCTK_BUILTIN_EXPECT(keyerr[i] != 0, false)) {
// now try negative temperature treatment
double eps0, eps1;
int idx[8];
File [modified]: nuc_eos_press_cs2.cc
Delta lines: +6 -6
===================================================================
--- trunk/src/nuc_eos_cxx/nuc_eos_press_cs2.cc 2014-03-13 03:01:52 UTC (rev 100)
+++ trunk/src/nuc_eos_cxx/nuc_eos_press_cs2.cc 2014-03-13 03:01:55 UTC (rev 101)
@@ -26,7 +26,7 @@
for(int i=0;i<n;i++) {
// check if we are fine
keyerr[i] = checkbounds(rho[i], temp[i], ye[i]);
- if(keyerr[i] != 0) {
+ if(CCTK_BUILTIN_EXPECT(keyerr[i] != 0, false)) {
*anyerr = 1;
}
}
@@ -63,7 +63,7 @@
prs[i] = exp(prs[i]);
eps[i] = exp(eps[i]) - energy_shift;
#if HAVEGR
- cs2[i] = cs2[i] / (1.0 + eps[i] + prs[i]/rho[i]);
+ cs2[i] = rho[i] * cs2[i] / (rho[i] + rho[i] * eps[i] + prs[i]);
#endif
}
@@ -96,7 +96,7 @@
// Note that this code now requires that the
// temperature guess be within the table bounds
keyerr[i] = checkbounds_kt0_noTcheck(rho[i], ye[i]);
- if(keyerr[i] != 0) {
+ if(CCTK_BUILTIN_EXPECT(keyerr[i] != 0, false)) {
*anyerr = 1;
}
}
@@ -112,7 +112,7 @@
const double lt = log(MIN(MAX(temp[i],eos_tempmin),eos_tempmax));
double ltout;
const double epstot = eps[i]+energy_shift;
- if(epstot>0.0e0) {
+ if(CCTK_BUILTIN_EXPECT(epstot>0.0e0, true)) {
// this is the standard scenario; eps is larger than zero
// and we can operate with logarithmic tables
const double lxeps = log(epstot);
@@ -124,7 +124,7 @@
keyerr[i] = 667;
} // epstot > 0.0
- if(keyerr[i]) {
+ if(CCTK_BUILTIN_EXPECT(keyerr[i] != 0, false)) {
// now try negative temperature treatment
double eps0, eps1;
int idx[8];
@@ -175,7 +175,7 @@
for(int i=0;i<n;i++) {
prs[i] = exp(prs[i]);
#if HAVEGR
- cs2[i] = cs2[i] / (1.0 + eps[i] + prs[i]/rho[i]);
+ cs2[i] = rho[i] * cs2[i] / (rho[i] + rho[i]*eps[i] + prs[i]);
#endif
}
File [modified]: nuc_eos_short.cc
Delta lines: +9 -9
===================================================================
--- trunk/src/nuc_eos_cxx/nuc_eos_short.cc 2014-03-13 03:01:52 UTC (rev 100)
+++ trunk/src/nuc_eos_cxx/nuc_eos_short.cc 2014-03-13 03:01:55 UTC (rev 101)
@@ -31,7 +31,7 @@
for(int i=0;i<n;i++) {
// check if we are fine
keyerr[i] = checkbounds(rho[i], temp[i], ye[i]);
- if(keyerr[i] != 0) {
+ if(CCTK_BUILTIN_EXPECT(keyerr[i] != 0, false)) {
*anyerr = 1;
}
}
@@ -94,7 +94,7 @@
prs[i] = exp(prs[i]);
eps[i] = exp(eps[i]) - energy_shift;
#if HAVEGR
- cs2[i] = cs2[i] / (1.0 + eps[i] + prs[i]/rho[i]);
+ cs2[i] = rho[i] * cs2[i] / (rho[i] + rho[i]*eps[i] + prs[i]);
#endif
}
@@ -131,7 +131,7 @@
// Note that this code now requires that the
// temperature guess be within the table bounds
keyerr[i] = checkbounds_kt0_noTcheck(rho[i], ye[i]);
- if(keyerr[i] != 0) {
+ if(CCTK_BUILTIN_EXPECT(keyerr[i] != 0, false)) {
*anyerr = 1;
}
}
@@ -146,7 +146,7 @@
const double lt = log(MIN(MAX(temp[i],eos_tempmin),eos_tempmax));
double ltout;
const double epstot = eps[i]+energy_shift;
- if(epstot>0.0e0) {
+ if(CCTK_BUILTIN_EXPECT(epstot>0.0e0, true)) {
// this is the standard scenario; eps is larger than zero
// and we can operate with logarithmic tables
const double lxeps = log(epstot);
@@ -163,7 +163,7 @@
} // epstot > 0.0
- if(keyerr[i]) {
+ if(CCTK_BUILTIN_EXPECT(keyerr[i] != 0, false)) {
*anyerr = 1;
} else {
temp[i] = exp(ltout);
@@ -211,7 +211,7 @@
for(int i=0;i<n;i++) {
prs[i] = exp(prs[i]);
#if HAVEGR
- cs2[i] = cs2[i] / (1.0 + eps[i] + prs[i]/rho[i]);
+ cs2[i] = rho[i] * cs2[i] / (rho[i] + rho[i] * eps[i] + prs[i]);
#endif
}
@@ -247,7 +247,7 @@
// Note that this code now requires that the
// temperature guess be within the table bounds
keyerr[i] = checkbounds_kt0_noTcheck(rho[i], ye[i]);
- if(keyerr[i] != 0) {
+ if(CCTK_BUILTIN_EXPECT(keyerr[i] != 0, false)) {
*anyerr = 1;
}
}
@@ -264,7 +264,7 @@
nuc_eos_findtemp_entropy(lr,lt,ye[i],ent[i],*prec,
(double *restrict)(<out),&keyerr[i]);
- if(keyerr[i]) {
+ if(CCTK_BUILTIN_EXPECT(keyerr[i] != 0, false)) {
*anyerr = 1;
} else {
temp[i] = exp(ltout);
@@ -313,7 +313,7 @@
prs[i] = exp(prs[i]);
eps[i] = exp(eps[i]) - energy_shift;
#if HAVEGR
- cs2[i] = cs2[i] / (1.0 + eps[i] + prs[i]/rho[i]);
+ cs2[i] = rho[i] * cs2[i] / (rho[i] + rho[i] * eps[i] + prs[i]);
#endif
}
File [modified]: readtable.cc
Delta lines: +8 -4
===================================================================
--- trunk/src/nuc_eos_cxx/readtable.cc 2014-03-13 03:01:52 UTC (rev 100)
+++ trunk/src/nuc_eos_cxx/readtable.cc 2014-03-13 03:01:55 UTC (rev 101)
@@ -194,11 +194,15 @@
// pressure
energy_shift = energy_shift * EPSGF;
for(int i=0;i<nrho;i++) {
- logrho[i] = log(pow(10.0,logrho[i]) * RHOGF);
+ // rewrite:
+ //logrho[i] = log(pow(10.0,logrho[i]) * RHOGF);
+ // by using log(a^b*c) = b*log(a)+log(c)
+ logrho[i] = logrho[i] * log(10.) + log(RHOGF);
}
for(int i=0;i<ntemp;i++) {
- logtemp[i] = log(pow(10.0,logtemp[i]));
+ //logtemp[i] = log(pow(10.0,logtemp[i]));
+ logtemp[i] = logtemp[i]*log(10.0);
}
// allocate epstable; a linear-scale eps table
@@ -214,12 +218,12 @@
{ // pressure
int idx = 0 + NTABLES*i;
- alltables[idx] = log(pow(10.0,alltables[idx])*PRESSGF);
+ alltables[idx] = alltables[idx] * log(10.0) + log(PRESSGF);
}
{ // eps
int idx = 1 + NTABLES*i;
- alltables[idx] = log(pow(10.0,alltables[idx])*EPSGF);
+ alltables[idx] = alltables[idx] * log(10.0) + log(EPSGF);
epstable[i] = exp(alltables[idx]);
}
