optimization

Calculate u**2 in advance. Manual common sub expression.

source/f90/fint1.f90

@@ -21,7 +21,7 @@ double precision function fint1(u, eps, d, layers, nper, eps_size)
   logical :: rational, user
   double precision :: acoef, bcoef, ccoef, cospsi, den, dif, dlimf, e, elleps
   double precision :: pi, rom, rop, ru, sinpsi, sum, um, usurlo, t
-  double precision :: tanpsi, wn
+  double precision :: tanpsi, wn, u2
 
   common / param / acoef, bcoef, ccoef, elleps, cospsi, sinpsi, tanpsi,  &
       ru, um, dlimf, wn, user, rational
@@ -37,14 +37,15 @@ double precision function fint1(u, eps, d, layers, nper, eps_size)
     return
   endif
   e = tanpsi * u
-  rom = (1.0d0 - e)**2 + u**2
-  rop = (1.0d0 + e)**2 + u**2
+  u2 = u**2
+  rom = (1.0d0 - e)**2 + u2
+  rop = (1.0d0 + e)**2 + u2
   sum = rop + rom
   rom = sqrt(rom)
   rop = sqrt(rop)
   dif = rop - rom
   den = sqrt((2.0d0 - dif) * (2.0d0 + dif)) * rop * rom
-  fint1 = pi * u**2 * (4.0d0 * sum - dif**2 * (sum - rop * rom)) / den**3
+  fint1 = pi * u2 * (4.0d0 * sum - dif**2 * (sum - rop * rom)) / den**3
   if (rational) then
     return
   endif

source/f90/o1.f90

@@ -3,8 +3,10 @@ double precision function o1(u)
   implicit none
   
   double precision, intent(in) :: u
+  double precision :: u2
 
-  o1 = -sinh(u) * log(u**2) + u * ((0.03114d0 * u**2 + 0.41666d0) * u**2 + 0.84557d0)
+  u2 = u**2
+  o1 = -sinh(u) * log(u2) + u * ((0.03114d0 * u2 + 0.41666d0) * u2 + 0.84557d0)
 
   return
 end function o1

source/f90/o2.f90

@@ -3,9 +3,11 @@ double precision function o2(u)
   implicit none
   
   double precision, intent(in) :: u
+  double precision :: u2
 
-  o2 = (((202.91d0 / u**2 + 932.21d0) / u**2 + 41.740d0) / u**2 + 2.0d0) /  &
-       (((540.88d0 / u**2 + 345.67d0) / u**2 + 18.961d0) / u**2 + 1.0d0) / u
+  u2 = 1 / u**2
+  o2 = (((202.91d0 * u2 + 932.21d0) * u2 + 41.740d0) * u2 + 2.0d0) /  &
+       (((540.88d0 * u2 + 345.67d0) * u2 + 18.961d0) * u2 + 1.0d0) / u
 
   return
 end function o2