
(FPCore (w0 M D h l d) :precision binary64 (* w0 (sqrt (- 1.0 (* (pow (/ (* M D) (* 2.0 d)) 2.0) (/ h l))))))
double code(double w0, double M, double D, double h, double l, double d) {
return w0 * sqrt((1.0 - (pow(((M * D) / (2.0 * d)), 2.0) * (h / l))));
}
module fmin_fmax_functions
implicit none
private
public fmax
public fmin
interface fmax
module procedure fmax88
module procedure fmax44
module procedure fmax84
module procedure fmax48
end interface
interface fmin
module procedure fmin88
module procedure fmin44
module procedure fmin84
module procedure fmin48
end interface
contains
real(8) function fmax88(x, y) result (res)
real(8), intent (in) :: x
real(8), intent (in) :: y
res = merge(y, merge(x, max(x, y), y /= y), x /= x)
end function
real(4) function fmax44(x, y) result (res)
real(4), intent (in) :: x
real(4), intent (in) :: y
res = merge(y, merge(x, max(x, y), y /= y), x /= x)
end function
real(8) function fmax84(x, y) result(res)
real(8), intent (in) :: x
real(4), intent (in) :: y
res = merge(dble(y), merge(x, max(x, dble(y)), y /= y), x /= x)
end function
real(8) function fmax48(x, y) result(res)
real(4), intent (in) :: x
real(8), intent (in) :: y
res = merge(y, merge(dble(x), max(dble(x), y), y /= y), x /= x)
end function
real(8) function fmin88(x, y) result (res)
real(8), intent (in) :: x
real(8), intent (in) :: y
res = merge(y, merge(x, min(x, y), y /= y), x /= x)
end function
real(4) function fmin44(x, y) result (res)
real(4), intent (in) :: x
real(4), intent (in) :: y
res = merge(y, merge(x, min(x, y), y /= y), x /= x)
end function
real(8) function fmin84(x, y) result(res)
real(8), intent (in) :: x
real(4), intent (in) :: y
res = merge(dble(y), merge(x, min(x, dble(y)), y /= y), x /= x)
end function
real(8) function fmin48(x, y) result(res)
real(4), intent (in) :: x
real(8), intent (in) :: y
res = merge(y, merge(dble(x), min(dble(x), y), y /= y), x /= x)
end function
end module
real(8) function code(w0, m, d, h, l, d_1)
use fmin_fmax_functions
real(8), intent (in) :: w0
real(8), intent (in) :: m
real(8), intent (in) :: d
real(8), intent (in) :: h
real(8), intent (in) :: l
real(8), intent (in) :: d_1
code = w0 * sqrt((1.0d0 - ((((m * d) / (2.0d0 * d_1)) ** 2.0d0) * (h / l))))
end function
public static double code(double w0, double M, double D, double h, double l, double d) {
return w0 * Math.sqrt((1.0 - (Math.pow(((M * D) / (2.0 * d)), 2.0) * (h / l))));
}
def code(w0, M, D, h, l, d): return w0 * math.sqrt((1.0 - (math.pow(((M * D) / (2.0 * d)), 2.0) * (h / l))))
function code(w0, M, D, h, l, d) return Float64(w0 * sqrt(Float64(1.0 - Float64((Float64(Float64(M * D) / Float64(2.0 * d)) ^ 2.0) * Float64(h / l))))) end
function tmp = code(w0, M, D, h, l, d) tmp = w0 * sqrt((1.0 - ((((M * D) / (2.0 * d)) ^ 2.0) * (h / l)))); end
code[w0_, M_, D_, h_, l_, d_] := N[(w0 * N[Sqrt[N[(1.0 - N[(N[Power[N[(N[(M * D), $MachinePrecision] / N[(2.0 * d), $MachinePrecision]), $MachinePrecision], 2.0], $MachinePrecision] * N[(h / l), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
w0 \cdot \sqrt{1 - {\left(\frac{M \cdot D}{2 \cdot d}\right)}^{2} \cdot \frac{h}{\ell}}
\end{array}
Sampling outcomes in binary64 precision:
Herbie found 4 alternatives:
| Alternative | Accuracy | Speedup |
|---|
(FPCore (w0 M D h l d) :precision binary64 (* w0 (sqrt (- 1.0 (* (pow (/ (* M D) (* 2.0 d)) 2.0) (/ h l))))))
double code(double w0, double M, double D, double h, double l, double d) {
return w0 * sqrt((1.0 - (pow(((M * D) / (2.0 * d)), 2.0) * (h / l))));
}
module fmin_fmax_functions
implicit none
private
public fmax
public fmin
interface fmax
module procedure fmax88
module procedure fmax44
module procedure fmax84
module procedure fmax48
end interface
interface fmin
module procedure fmin88
module procedure fmin44
module procedure fmin84
module procedure fmin48
end interface
contains
real(8) function fmax88(x, y) result (res)
real(8), intent (in) :: x
real(8), intent (in) :: y
res = merge(y, merge(x, max(x, y), y /= y), x /= x)
end function
real(4) function fmax44(x, y) result (res)
real(4), intent (in) :: x
real(4), intent (in) :: y
res = merge(y, merge(x, max(x, y), y /= y), x /= x)
end function
real(8) function fmax84(x, y) result(res)
real(8), intent (in) :: x
real(4), intent (in) :: y
res = merge(dble(y), merge(x, max(x, dble(y)), y /= y), x /= x)
end function
real(8) function fmax48(x, y) result(res)
real(4), intent (in) :: x
real(8), intent (in) :: y
res = merge(y, merge(dble(x), max(dble(x), y), y /= y), x /= x)
end function
real(8) function fmin88(x, y) result (res)
real(8), intent (in) :: x
real(8), intent (in) :: y
res = merge(y, merge(x, min(x, y), y /= y), x /= x)
end function
real(4) function fmin44(x, y) result (res)
real(4), intent (in) :: x
real(4), intent (in) :: y
res = merge(y, merge(x, min(x, y), y /= y), x /= x)
end function
real(8) function fmin84(x, y) result(res)
real(8), intent (in) :: x
real(4), intent (in) :: y
res = merge(dble(y), merge(x, min(x, dble(y)), y /= y), x /= x)
end function
real(8) function fmin48(x, y) result(res)
real(4), intent (in) :: x
real(8), intent (in) :: y
res = merge(y, merge(dble(x), min(dble(x), y), y /= y), x /= x)
end function
end module
real(8) function code(w0, m, d, h, l, d_1)
use fmin_fmax_functions
real(8), intent (in) :: w0
real(8), intent (in) :: m
real(8), intent (in) :: d
real(8), intent (in) :: h
real(8), intent (in) :: l
real(8), intent (in) :: d_1
code = w0 * sqrt((1.0d0 - ((((m * d) / (2.0d0 * d_1)) ** 2.0d0) * (h / l))))
end function
public static double code(double w0, double M, double D, double h, double l, double d) {
return w0 * Math.sqrt((1.0 - (Math.pow(((M * D) / (2.0 * d)), 2.0) * (h / l))));
}
def code(w0, M, D, h, l, d): return w0 * math.sqrt((1.0 - (math.pow(((M * D) / (2.0 * d)), 2.0) * (h / l))))
function code(w0, M, D, h, l, d) return Float64(w0 * sqrt(Float64(1.0 - Float64((Float64(Float64(M * D) / Float64(2.0 * d)) ^ 2.0) * Float64(h / l))))) end
function tmp = code(w0, M, D, h, l, d) tmp = w0 * sqrt((1.0 - ((((M * D) / (2.0 * d)) ^ 2.0) * (h / l)))); end
code[w0_, M_, D_, h_, l_, d_] := N[(w0 * N[Sqrt[N[(1.0 - N[(N[Power[N[(N[(M * D), $MachinePrecision] / N[(2.0 * d), $MachinePrecision]), $MachinePrecision], 2.0], $MachinePrecision] * N[(h / l), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
w0 \cdot \sqrt{1 - {\left(\frac{M \cdot D}{2 \cdot d}\right)}^{2} \cdot \frac{h}{\ell}}
\end{array}
w0\_m = (fabs.f64 w0)
w0\_s = (copysign.f64 #s(literal 1 binary64) w0)
(FPCore (w0_s w0_m M D h l d)
:precision binary64
(let* ((t_0 (* (/ D d) (/ M 2.0)))
(t_1
(* w0_m (sqrt (- 1.0 (* (pow (/ (* M D) (* 2.0 d)) 2.0) (/ h l)))))))
(*
w0_s
(if (<= t_1 5e+186) t_1 (* w0_m (sqrt (- 1.0 (/ (* t_0 (* t_0 h)) l))))))))w0\_m = fabs(w0);
w0\_s = copysign(1.0, w0);
double code(double w0_s, double w0_m, double M, double D, double h, double l, double d) {
double t_0 = (D / d) * (M / 2.0);
double t_1 = w0_m * sqrt((1.0 - (pow(((M * D) / (2.0 * d)), 2.0) * (h / l))));
double tmp;
if (t_1 <= 5e+186) {
tmp = t_1;
} else {
tmp = w0_m * sqrt((1.0 - ((t_0 * (t_0 * h)) / l)));
}
return w0_s * tmp;
}
w0\_m = private
w0\_s = private
module fmin_fmax_functions
implicit none
private
public fmax
public fmin
interface fmax
module procedure fmax88
module procedure fmax44
module procedure fmax84
module procedure fmax48
end interface
interface fmin
module procedure fmin88
module procedure fmin44
module procedure fmin84
module procedure fmin48
end interface
contains
real(8) function fmax88(x, y) result (res)
real(8), intent (in) :: x
real(8), intent (in) :: y
res = merge(y, merge(x, max(x, y), y /= y), x /= x)
end function
real(4) function fmax44(x, y) result (res)
real(4), intent (in) :: x
real(4), intent (in) :: y
res = merge(y, merge(x, max(x, y), y /= y), x /= x)
end function
real(8) function fmax84(x, y) result(res)
real(8), intent (in) :: x
real(4), intent (in) :: y
res = merge(dble(y), merge(x, max(x, dble(y)), y /= y), x /= x)
end function
real(8) function fmax48(x, y) result(res)
real(4), intent (in) :: x
real(8), intent (in) :: y
res = merge(y, merge(dble(x), max(dble(x), y), y /= y), x /= x)
end function
real(8) function fmin88(x, y) result (res)
real(8), intent (in) :: x
real(8), intent (in) :: y
res = merge(y, merge(x, min(x, y), y /= y), x /= x)
end function
real(4) function fmin44(x, y) result (res)
real(4), intent (in) :: x
real(4), intent (in) :: y
res = merge(y, merge(x, min(x, y), y /= y), x /= x)
end function
real(8) function fmin84(x, y) result(res)
real(8), intent (in) :: x
real(4), intent (in) :: y
res = merge(dble(y), merge(x, min(x, dble(y)), y /= y), x /= x)
end function
real(8) function fmin48(x, y) result(res)
real(4), intent (in) :: x
real(8), intent (in) :: y
res = merge(y, merge(dble(x), min(dble(x), y), y /= y), x /= x)
end function
end module
real(8) function code(w0_s, w0_m, m, d, h, l, d_1)
use fmin_fmax_functions
real(8), intent (in) :: w0_s
real(8), intent (in) :: w0_m
real(8), intent (in) :: m
real(8), intent (in) :: d
real(8), intent (in) :: h
real(8), intent (in) :: l
real(8), intent (in) :: d_1
real(8) :: t_0
real(8) :: t_1
real(8) :: tmp
t_0 = (d / d_1) * (m / 2.0d0)
t_1 = w0_m * sqrt((1.0d0 - ((((m * d) / (2.0d0 * d_1)) ** 2.0d0) * (h / l))))
if (t_1 <= 5d+186) then
tmp = t_1
else
tmp = w0_m * sqrt((1.0d0 - ((t_0 * (t_0 * h)) / l)))
end if
code = w0_s * tmp
end function
w0\_m = Math.abs(w0);
w0\_s = Math.copySign(1.0, w0);
public static double code(double w0_s, double w0_m, double M, double D, double h, double l, double d) {
double t_0 = (D / d) * (M / 2.0);
double t_1 = w0_m * Math.sqrt((1.0 - (Math.pow(((M * D) / (2.0 * d)), 2.0) * (h / l))));
double tmp;
if (t_1 <= 5e+186) {
tmp = t_1;
} else {
tmp = w0_m * Math.sqrt((1.0 - ((t_0 * (t_0 * h)) / l)));
}
return w0_s * tmp;
}
w0\_m = math.fabs(w0) w0\_s = math.copysign(1.0, w0) def code(w0_s, w0_m, M, D, h, l, d): t_0 = (D / d) * (M / 2.0) t_1 = w0_m * math.sqrt((1.0 - (math.pow(((M * D) / (2.0 * d)), 2.0) * (h / l)))) tmp = 0 if t_1 <= 5e+186: tmp = t_1 else: tmp = w0_m * math.sqrt((1.0 - ((t_0 * (t_0 * h)) / l))) return w0_s * tmp
w0\_m = abs(w0) w0\_s = copysign(1.0, w0) function code(w0_s, w0_m, M, D, h, l, d) t_0 = Float64(Float64(D / d) * Float64(M / 2.0)) t_1 = Float64(w0_m * sqrt(Float64(1.0 - Float64((Float64(Float64(M * D) / Float64(2.0 * d)) ^ 2.0) * Float64(h / l))))) tmp = 0.0 if (t_1 <= 5e+186) tmp = t_1; else tmp = Float64(w0_m * sqrt(Float64(1.0 - Float64(Float64(t_0 * Float64(t_0 * h)) / l)))); end return Float64(w0_s * tmp) end
w0\_m = abs(w0); w0\_s = sign(w0) * abs(1.0); function tmp_2 = code(w0_s, w0_m, M, D, h, l, d) t_0 = (D / d) * (M / 2.0); t_1 = w0_m * sqrt((1.0 - ((((M * D) / (2.0 * d)) ^ 2.0) * (h / l)))); tmp = 0.0; if (t_1 <= 5e+186) tmp = t_1; else tmp = w0_m * sqrt((1.0 - ((t_0 * (t_0 * h)) / l))); end tmp_2 = w0_s * tmp; end
w0\_m = N[Abs[w0], $MachinePrecision]
w0\_s = N[With[{TMP1 = Abs[1.0], TMP2 = Sign[w0]}, TMP1 * If[TMP2 == 0, 1, TMP2]], $MachinePrecision]
code[w0$95$s_, w0$95$m_, M_, D_, h_, l_, d_] := Block[{t$95$0 = N[(N[(D / d), $MachinePrecision] * N[(M / 2.0), $MachinePrecision]), $MachinePrecision]}, Block[{t$95$1 = N[(w0$95$m * N[Sqrt[N[(1.0 - N[(N[Power[N[(N[(M * D), $MachinePrecision] / N[(2.0 * d), $MachinePrecision]), $MachinePrecision], 2.0], $MachinePrecision] * N[(h / l), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]}, N[(w0$95$s * If[LessEqual[t$95$1, 5e+186], t$95$1, N[(w0$95$m * N[Sqrt[N[(1.0 - N[(N[(t$95$0 * N[(t$95$0 * h), $MachinePrecision]), $MachinePrecision] / l), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]]), $MachinePrecision]]]
\begin{array}{l}
w0\_m = \left|w0\right|
\\
w0\_s = \mathsf{copysign}\left(1, w0\right)
\\
\begin{array}{l}
t_0 := \frac{D}{d} \cdot \frac{M}{2}\\
t_1 := w0\_m \cdot \sqrt{1 - {\left(\frac{M \cdot D}{2 \cdot d}\right)}^{2} \cdot \frac{h}{\ell}}\\
w0\_s \cdot \begin{array}{l}
\mathbf{if}\;t\_1 \leq 5 \cdot 10^{+186}:\\
\;\;\;\;t\_1\\
\mathbf{else}:\\
\;\;\;\;w0\_m \cdot \sqrt{1 - \frac{t\_0 \cdot \left(t\_0 \cdot h\right)}{\ell}}\\
\end{array}
\end{array}
\end{array}
if (*.f64 w0 (sqrt.f64 (-.f64 #s(literal 1 binary64) (*.f64 (pow.f64 (/.f64 (*.f64 M D) (*.f64 #s(literal 2 binary64) d)) #s(literal 2 binary64)) (/.f64 h l))))) < 4.99999999999999954e186Initial program 91.7%
if 4.99999999999999954e186 < (*.f64 w0 (sqrt.f64 (-.f64 #s(literal 1 binary64) (*.f64 (pow.f64 (/.f64 (*.f64 M D) (*.f64 #s(literal 2 binary64) d)) #s(literal 2 binary64)) (/.f64 h l))))) Initial program 45.0%
lift-*.f64N/A
lift-pow.f64N/A
lift-*.f64N/A
lift-*.f64N/A
lift-/.f64N/A
lift-/.f64N/A
associate-*r/N/A
lower-/.f64N/A
lower-*.f64N/A
unpow2N/A
lower-*.f64N/A
times-fracN/A
lower-*.f64N/A
lower-/.f64N/A
lower-/.f64N/A
times-fracN/A
lower-*.f64N/A
lower-/.f64N/A
lower-/.f6471.6
Applied rewrites71.6%
lift-*.f64N/A
lift-*.f64N/A
associate-*l*N/A
lower-*.f64N/A
lift-*.f64N/A
lift-/.f64N/A
*-commutativeN/A
lower-*.f64N/A
lift-/.f64N/A
lower-*.f6478.1
lift-*.f64N/A
lift-/.f64N/A
*-commutativeN/A
lower-*.f64N/A
lift-/.f6478.1
Applied rewrites78.1%
w0\_m = (fabs.f64 w0)
w0\_s = (copysign.f64 #s(literal 1 binary64) w0)
(FPCore (w0_s w0_m M D h l d)
:precision binary64
(let* ((t_0
(* w0_m (sqrt (- 1.0 (* (pow (/ (* M D) (* 2.0 d)) 2.0) (/ h l)))))))
(* w0_s (if (<= t_0 INFINITY) t_0 (* w0_m (sqrt 1.0))))))w0\_m = fabs(w0);
w0\_s = copysign(1.0, w0);
double code(double w0_s, double w0_m, double M, double D, double h, double l, double d) {
double t_0 = w0_m * sqrt((1.0 - (pow(((M * D) / (2.0 * d)), 2.0) * (h / l))));
double tmp;
if (t_0 <= ((double) INFINITY)) {
tmp = t_0;
} else {
tmp = w0_m * sqrt(1.0);
}
return w0_s * tmp;
}
w0\_m = Math.abs(w0);
w0\_s = Math.copySign(1.0, w0);
public static double code(double w0_s, double w0_m, double M, double D, double h, double l, double d) {
double t_0 = w0_m * Math.sqrt((1.0 - (Math.pow(((M * D) / (2.0 * d)), 2.0) * (h / l))));
double tmp;
if (t_0 <= Double.POSITIVE_INFINITY) {
tmp = t_0;
} else {
tmp = w0_m * Math.sqrt(1.0);
}
return w0_s * tmp;
}
w0\_m = math.fabs(w0) w0\_s = math.copysign(1.0, w0) def code(w0_s, w0_m, M, D, h, l, d): t_0 = w0_m * math.sqrt((1.0 - (math.pow(((M * D) / (2.0 * d)), 2.0) * (h / l)))) tmp = 0 if t_0 <= math.inf: tmp = t_0 else: tmp = w0_m * math.sqrt(1.0) return w0_s * tmp
w0\_m = abs(w0) w0\_s = copysign(1.0, w0) function code(w0_s, w0_m, M, D, h, l, d) t_0 = Float64(w0_m * sqrt(Float64(1.0 - Float64((Float64(Float64(M * D) / Float64(2.0 * d)) ^ 2.0) * Float64(h / l))))) tmp = 0.0 if (t_0 <= Inf) tmp = t_0; else tmp = Float64(w0_m * sqrt(1.0)); end return Float64(w0_s * tmp) end
w0\_m = abs(w0); w0\_s = sign(w0) * abs(1.0); function tmp_2 = code(w0_s, w0_m, M, D, h, l, d) t_0 = w0_m * sqrt((1.0 - ((((M * D) / (2.0 * d)) ^ 2.0) * (h / l)))); tmp = 0.0; if (t_0 <= Inf) tmp = t_0; else tmp = w0_m * sqrt(1.0); end tmp_2 = w0_s * tmp; end
w0\_m = N[Abs[w0], $MachinePrecision]
w0\_s = N[With[{TMP1 = Abs[1.0], TMP2 = Sign[w0]}, TMP1 * If[TMP2 == 0, 1, TMP2]], $MachinePrecision]
code[w0$95$s_, w0$95$m_, M_, D_, h_, l_, d_] := Block[{t$95$0 = N[(w0$95$m * N[Sqrt[N[(1.0 - N[(N[Power[N[(N[(M * D), $MachinePrecision] / N[(2.0 * d), $MachinePrecision]), $MachinePrecision], 2.0], $MachinePrecision] * N[(h / l), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]}, N[(w0$95$s * If[LessEqual[t$95$0, Infinity], t$95$0, N[(w0$95$m * N[Sqrt[1.0], $MachinePrecision]), $MachinePrecision]]), $MachinePrecision]]
\begin{array}{l}
w0\_m = \left|w0\right|
\\
w0\_s = \mathsf{copysign}\left(1, w0\right)
\\
\begin{array}{l}
t_0 := w0\_m \cdot \sqrt{1 - {\left(\frac{M \cdot D}{2 \cdot d}\right)}^{2} \cdot \frac{h}{\ell}}\\
w0\_s \cdot \begin{array}{l}
\mathbf{if}\;t\_0 \leq \infty:\\
\;\;\;\;t\_0\\
\mathbf{else}:\\
\;\;\;\;w0\_m \cdot \sqrt{1}\\
\end{array}
\end{array}
\end{array}
if (*.f64 w0 (sqrt.f64 (-.f64 #s(literal 1 binary64) (*.f64 (pow.f64 (/.f64 (*.f64 M D) (*.f64 #s(literal 2 binary64) d)) #s(literal 2 binary64)) (/.f64 h l))))) < +inf.0Initial program 88.9%
if +inf.0 < (*.f64 w0 (sqrt.f64 (-.f64 #s(literal 1 binary64) (*.f64 (pow.f64 (/.f64 (*.f64 M D) (*.f64 #s(literal 2 binary64) d)) #s(literal 2 binary64)) (/.f64 h l))))) Initial program 0.0%
Taylor expanded in M around 0
Applied rewrites77.5%
w0\_m = (fabs.f64 w0)
w0\_s = (copysign.f64 #s(literal 1 binary64) w0)
(FPCore (w0_s w0_m M D h l d)
:precision binary64
(*
w0_s
(if (<= (* (pow (/ (* M D) (* 2.0 d)) 2.0) (/ h l)) -1e+116)
(fma (/ (* (pow (* M D) 2.0) (* h w0_m)) (* d (* d l))) -0.125 w0_m)
(* w0_m (sqrt 1.0)))))w0\_m = fabs(w0);
w0\_s = copysign(1.0, w0);
double code(double w0_s, double w0_m, double M, double D, double h, double l, double d) {
double tmp;
if ((pow(((M * D) / (2.0 * d)), 2.0) * (h / l)) <= -1e+116) {
tmp = fma(((pow((M * D), 2.0) * (h * w0_m)) / (d * (d * l))), -0.125, w0_m);
} else {
tmp = w0_m * sqrt(1.0);
}
return w0_s * tmp;
}
w0\_m = abs(w0) w0\_s = copysign(1.0, w0) function code(w0_s, w0_m, M, D, h, l, d) tmp = 0.0 if (Float64((Float64(Float64(M * D) / Float64(2.0 * d)) ^ 2.0) * Float64(h / l)) <= -1e+116) tmp = fma(Float64(Float64((Float64(M * D) ^ 2.0) * Float64(h * w0_m)) / Float64(d * Float64(d * l))), -0.125, w0_m); else tmp = Float64(w0_m * sqrt(1.0)); end return Float64(w0_s * tmp) end
w0\_m = N[Abs[w0], $MachinePrecision]
w0\_s = N[With[{TMP1 = Abs[1.0], TMP2 = Sign[w0]}, TMP1 * If[TMP2 == 0, 1, TMP2]], $MachinePrecision]
code[w0$95$s_, w0$95$m_, M_, D_, h_, l_, d_] := N[(w0$95$s * If[LessEqual[N[(N[Power[N[(N[(M * D), $MachinePrecision] / N[(2.0 * d), $MachinePrecision]), $MachinePrecision], 2.0], $MachinePrecision] * N[(h / l), $MachinePrecision]), $MachinePrecision], -1e+116], N[(N[(N[(N[Power[N[(M * D), $MachinePrecision], 2.0], $MachinePrecision] * N[(h * w0$95$m), $MachinePrecision]), $MachinePrecision] / N[(d * N[(d * l), $MachinePrecision]), $MachinePrecision]), $MachinePrecision] * -0.125 + w0$95$m), $MachinePrecision], N[(w0$95$m * N[Sqrt[1.0], $MachinePrecision]), $MachinePrecision]]), $MachinePrecision]
\begin{array}{l}
w0\_m = \left|w0\right|
\\
w0\_s = \mathsf{copysign}\left(1, w0\right)
\\
w0\_s \cdot \begin{array}{l}
\mathbf{if}\;{\left(\frac{M \cdot D}{2 \cdot d}\right)}^{2} \cdot \frac{h}{\ell} \leq -1 \cdot 10^{+116}:\\
\;\;\;\;\mathsf{fma}\left(\frac{{\left(M \cdot D\right)}^{2} \cdot \left(h \cdot w0\_m\right)}{d \cdot \left(d \cdot \ell\right)}, -0.125, w0\_m\right)\\
\mathbf{else}:\\
\;\;\;\;w0\_m \cdot \sqrt{1}\\
\end{array}
\end{array}
if (*.f64 (pow.f64 (/.f64 (*.f64 M D) (*.f64 #s(literal 2 binary64) d)) #s(literal 2 binary64)) (/.f64 h l)) < -1.00000000000000002e116Initial program 63.3%
lift-*.f64N/A
lift-pow.f64N/A
lift-*.f64N/A
lift-*.f64N/A
lift-/.f64N/A
lift-/.f64N/A
associate-*r/N/A
lower-/.f64N/A
lower-*.f64N/A
unpow2N/A
lower-*.f64N/A
times-fracN/A
lower-*.f64N/A
lower-/.f64N/A
lower-/.f64N/A
times-fracN/A
lower-*.f64N/A
lower-/.f64N/A
lower-/.f6461.9
Applied rewrites61.9%
lift-*.f64N/A
lift-*.f64N/A
associate-*l*N/A
lower-*.f64N/A
lift-*.f64N/A
lift-/.f64N/A
*-commutativeN/A
lower-*.f64N/A
lift-/.f64N/A
lower-*.f6463.2
lift-*.f64N/A
lift-/.f64N/A
*-commutativeN/A
lower-*.f64N/A
lift-/.f6463.2
Applied rewrites63.2%
Taylor expanded in M around 0
+-commutativeN/A
*-commutativeN/A
lower-fma.f64N/A
lower-/.f64N/A
associate-*r*N/A
lower-*.f64N/A
pow-prod-downN/A
lower-pow.f64N/A
*-commutativeN/A
lower-*.f64N/A
lower-*.f64N/A
lower-*.f64N/A
unpow2N/A
lower-*.f6444.8
Applied rewrites44.8%
lift-*.f64N/A
lift-*.f64N/A
associate-*l*N/A
lower-*.f64N/A
lower-*.f6447.7
Applied rewrites47.7%
if -1.00000000000000002e116 < (*.f64 (pow.f64 (/.f64 (*.f64 M D) (*.f64 #s(literal 2 binary64) d)) #s(literal 2 binary64)) (/.f64 h l)) Initial program 87.1%
Taylor expanded in M around 0
Applied rewrites94.5%
w0\_m = (fabs.f64 w0) w0\_s = (copysign.f64 #s(literal 1 binary64) w0) (FPCore (w0_s w0_m M D h l d) :precision binary64 (* w0_s (fma (/ (* (pow (* M D) 2.0) (* h w0_m)) (* d (* d l))) -0.125 w0_m)))
w0\_m = fabs(w0);
w0\_s = copysign(1.0, w0);
double code(double w0_s, double w0_m, double M, double D, double h, double l, double d) {
return w0_s * fma(((pow((M * D), 2.0) * (h * w0_m)) / (d * (d * l))), -0.125, w0_m);
}
w0\_m = abs(w0) w0\_s = copysign(1.0, w0) function code(w0_s, w0_m, M, D, h, l, d) return Float64(w0_s * fma(Float64(Float64((Float64(M * D) ^ 2.0) * Float64(h * w0_m)) / Float64(d * Float64(d * l))), -0.125, w0_m)) end
w0\_m = N[Abs[w0], $MachinePrecision]
w0\_s = N[With[{TMP1 = Abs[1.0], TMP2 = Sign[w0]}, TMP1 * If[TMP2 == 0, 1, TMP2]], $MachinePrecision]
code[w0$95$s_, w0$95$m_, M_, D_, h_, l_, d_] := N[(w0$95$s * N[(N[(N[(N[Power[N[(M * D), $MachinePrecision], 2.0], $MachinePrecision] * N[(h * w0$95$m), $MachinePrecision]), $MachinePrecision] / N[(d * N[(d * l), $MachinePrecision]), $MachinePrecision]), $MachinePrecision] * -0.125 + w0$95$m), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
w0\_m = \left|w0\right|
\\
w0\_s = \mathsf{copysign}\left(1, w0\right)
\\
w0\_s \cdot \mathsf{fma}\left(\frac{{\left(M \cdot D\right)}^{2} \cdot \left(h \cdot w0\_m\right)}{d \cdot \left(d \cdot \ell\right)}, -0.125, w0\_m\right)
\end{array}
Initial program 80.6%
lift-*.f64N/A
lift-pow.f64N/A
lift-*.f64N/A
lift-*.f64N/A
lift-/.f64N/A
lift-/.f64N/A
associate-*r/N/A
lower-/.f64N/A
lower-*.f64N/A
unpow2N/A
lower-*.f64N/A
times-fracN/A
lower-*.f64N/A
lower-/.f64N/A
lower-/.f64N/A
times-fracN/A
lower-*.f64N/A
lower-/.f64N/A
lower-/.f6486.2
Applied rewrites86.2%
lift-*.f64N/A
lift-*.f64N/A
associate-*l*N/A
lower-*.f64N/A
lift-*.f64N/A
lift-/.f64N/A
*-commutativeN/A
lower-*.f64N/A
lift-/.f64N/A
lower-*.f6488.1
lift-*.f64N/A
lift-/.f64N/A
*-commutativeN/A
lower-*.f64N/A
lift-/.f6488.1
Applied rewrites88.1%
Taylor expanded in M around 0
+-commutativeN/A
*-commutativeN/A
lower-fma.f64N/A
lower-/.f64N/A
associate-*r*N/A
lower-*.f64N/A
pow-prod-downN/A
lower-pow.f64N/A
*-commutativeN/A
lower-*.f64N/A
lower-*.f64N/A
lower-*.f64N/A
unpow2N/A
lower-*.f6458.5
Applied rewrites58.5%
lift-*.f64N/A
lift-*.f64N/A
associate-*l*N/A
lower-*.f64N/A
lower-*.f6462.3
Applied rewrites62.3%
herbie shell --seed 2025066
(FPCore (w0 M D h l d)
:name "Henrywood and Agarwal, Equation (9a)"
:precision binary64
(* w0 (sqrt (- 1.0 (* (pow (/ (* M D) (* 2.0 d)) 2.0) (/ h l))))))