# Optimizing sums of products (dot products)

**URL:** <https://discourse.julialang.org/t/optimizing-sums-of-products-dot-products/119515>\
**Category:** Performance\
**Tags:** simd, linearalgebra, sum\
**Created:** [September 17, 2024, 4:40pm UTC](https://discourse.julialang.org/t/optimizing-sums-of-products-dot-products/119515 "2024-09-17T16:40:15Z")\
**Posts on this page:** 1\
**Showing post:** 18

<div class="post-metadata">

**Author:** ![mbauman](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mbauman/32/31082_2.png) [@mbauman](https://discourse.julialang.org/u/mbauman)\
**Post date:** [September 24, 2024, 7:40pm UTC](https://discourse.julialang.org/t/optimizing-sums-of-products-dot-products/119515/18 "2024-09-24T19:40:11Z")

</div>

> [@mbauman](#):
>
> I’d love to build a realistic example that demonstrates different values in an `if` than the `sqrt`

OK, this was just too juicy to resist. Here’s a barely-contrived example (Julia v1.10.5, Apple M1):

```julia-repl
julia> discriminant(a,b,c) = @fastmath b^2 - 4*a*c
discriminant (generic function with 1 method)

julia> function roots(a,b,c)
           if discriminant(a,b,c) < 0
               "no real roots"
           else
               roots = [(-b ± sqrt(discriminant(a,b,c)))/(2*a) for (±) in (+, -)]
               "root at $(unique(roots))"
           end
       end
roots (generic function with 1 method)

julia> roots(.2, .1, .0125)
"no real roots"

```

That works just fine. Sure, you might miss an exact root here or there, but those are the edge of numerical stability anyhow. What you probably wouldn’t expect is how easy it is for the above answer to change (and error!) depending on what you happen to do **_outside_** that function!

```julia-repl
julia> f1(a,b,c) = roots(a,b,c), b^2.0
f1 (generic function with 1 method)

julia> f1(.2, .1, .0125)
("no real roots", 0.010000000000000002)

julia> f2(a,b,c) = roots(a,b,c), b*b
f2 (generic function with 1 method)

julia> f2(.2, .1, .0125)
ERROR: DomainError with -8.326672684688674e-19:
sqrt was called with a negative real argument but will only return a complex result if called with a complex argument. Try sqrt(Complex(x)).

```

Is that scary enough? That said, clang and gcc do this sort of thing by default and you don’t see many folks complaining of issues like this.

---

_[View the full topic](https://discourse.julialang.org/t/optimizing-sums-of-products-dot-products/119515)._
