# Faster isapprox()

**URL:** <https://discourse.julialang.org/t/faster-isapprox/101202>\
**Category:** General Usage\
**Tags:** question, performance, linearalgebra, loopvectorization, preallocation\
**Created:** [July 5, 2023, 10:12am UTC](https://discourse.julialang.org/t/faster-isapprox/101202 "2023-07-05T10:12:04Z")\
**Posts on this page:** 15\
**Page:** 1

<div class="post-metadata">

**Author:** ![user\_231578](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/user_231578/32/24170_2.png) [@user\_231578](https://discourse.julialang.org/u/user_231578)\
**Post date:** [July 5, 2023, 10:12am UTC](https://discourse.julialang.org/t/faster-isapprox/101202/1 "2023-07-05T10:12:04Z")

</div>

I have to run (possibly) millions of times the function isapprox(A, A\_, rtol=1e-4), where A and A\_ are two matrices. To improve efficiency, I defined a custom isapprox() function using a more efficient norm:

```julia
using LoopVectorization

function is_approx(A::Array{Float64,2}, A_::Array{Float64,2}, rtol::Float64)
	# similar to built-in isapprox() but faster
	return norm(A - A_) < rtol * max(norm(A), norm(A_))
end

function norm(A::Array{Float64,2})
	# fast L(2,2) norm
	# discussion: https://discourse.julialang.org/t/performance-of-norm-function/14709/13
	x = zero(eltype(A))
	@avx for i in eachindex(A)
		@fastmath x += A[i] * A[i]
	end
	@fastmath sqrt(x)
end

```

If you want, you can find more context about what I use this for [here](https://github.com/massimilianofurlan/ql_agent/blob/main/ql_agent.jl) .

Despite this being considerably more efficient than isapprox(), is\_approx() is still the main performance bottleneck of my code. Do you have any ideas on how to make is\_approx() more efficient?

If it is relevant, profiling the function,

```julia
using StatProfilerHTML 
A, A_ = rand(10,10), rand(10,10)
@profilehtml for i in 1:10^6 is_approx(A,A_,1e-4) end

```

gives the following report:

 ![Screenshot 2023-07-05 at 10.58.57 AM](https://global.discourse-cdn.com/julialang/original/3X/a/9/a9e29b24f28fb9bf55debbb83840e21cb8ea1d44.png)

---

<div class="post-metadata">

**Author:** ![Dan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dan/32/42581_2.png) [@Dan](https://discourse.julialang.org/u/Dan)\
**Post date:** [July 5, 2023, 10:20am UTC](https://discourse.julialang.org/t/faster-isapprox/101202/2 "2023-07-05T10:20:08Z")

</div>

`is_approx` is calling `norm` 3 times (for `A - A_` , `A` , and `A_`).  
Fusing `norm` with `is_approx` can allow calculation of the norms together (keeping three accumulators for the relevant squared numbers), and will reduce memory access 3x leading to significant improvement to performance.

---

<div class="post-metadata">

**Author:** ![user\_231578](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/user_231578/32/24170_2.png) [@user\_231578](https://discourse.julialang.org/u/user_231578)\
**Post date:** [July 5, 2023, 10:25am UTC](https://discourse.julialang.org/t/faster-isapprox/101202/3 "2023-07-05T10:25:08Z")

</div>

What does fusing `norm` with `is_approx` mean? Would you be able to propose an example for your solution? Thanks for your help!

---

<div class="post-metadata">

**Author:** ![jondea](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jondea/32/39086_2.png) [@jondea](https://discourse.julialang.org/u/jondea)\
**Post date:** [July 5, 2023, 10:25am UTC](https://discourse.julialang.org/t/faster-isapprox/101202/4 "2023-07-05T10:25:09Z")

</div>

Given that you are comparing norms, can you skip all the sqrts and square your rtol instead?

---

<div class="post-metadata">

**Author:** ![Dan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dan/32/42581_2.png) [@Dan](https://discourse.julialang.org/u/Dan)\
**Post date:** [July 5, 2023, 10:38am UTC](https://discourse.julialang.org/t/faster-isapprox/101202/5 "2023-07-05T10:38:34Z")

</div>

Especially, you should skip calculating `A - A_` which is also allocating memory (slow). This calculation should not allocate any memory.

---

<div class="post-metadata">

**Author:** ![user\_231578](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/user_231578/32/24170_2.png) [@user\_231578](https://discourse.julialang.org/u/user_231578)\
**Post date:** [July 5, 2023, 10:39am UTC](https://discourse.julialang.org/t/faster-isapprox/101202/6 "2023-07-05T10:39:11Z")

</div>

Thanks for helping. This seems to make it just slightly faster.

```julia

function is_approx2(A::Array{Float64,2}, A_::Array{Float64,2}, rtol::Float64)
	# similar to built-in isapprox() but faster
	return norm2(A - A_) < abs2(rtol) * max(norm2(A), norm2(A_))
end

function norm2(A::Array{Float64,2})
	# fast L(2,2) norm
	# discussion: https://discourse.julialang.org/t/performance-of-norm-function/14709/13
	x = zero(eltype(A))
	@avx for i in eachindex(A)
		@fastmath x += A[i] * A[i]
	end
	return x
end

@btime is_approx($A,$A_,$1e-4)
  96.901 ns (1 allocation: 896 bytes)

@btime is_approx2($A,$A_,$1e-4)
  93.640 ns (1 allocation: 896 bytes)

```

---

<div class="post-metadata">

**Author:** ![gdalle](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/gdalle/32/27854_2.png) [@gdalle](https://discourse.julialang.org/u/gdalle)\
**Post date:** [July 5, 2023, 10:45am UTC](https://discourse.julialang.org/t/faster-isapprox/101202/7 "2023-07-05T10:45:13Z")

</div>

The main bottleneck of your code is allocation when you create the new matrix `A - A_`.

---

<div class="post-metadata">

**Author:** ![Dan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dan/32/42581_2.png) [@Dan](https://discourse.julialang.org/u/Dan)\
**Post date:** [July 5, 2023, 10:46am UTC](https://discourse.julialang.org/t/faster-isapprox/101202/8 "2023-07-05T10:46:27Z")

</div>

Something efficient might look like:

```julia
function is_approx3(A::Array{Float64,2}, A_::Array{Float64,2}, rtol::Float64)
    # make sure sizes of A and A_ are the same here (?)
    s1, s2, s3 = 0.0, 0.0, 0.0
    for i in eachindex(A)
        a, a_ = A[i], A_[i]
        s1 += (a - a_)^2
        s2 += a^2
        s3 += a_^2
    end
    return s1 < rtol^2*max(s2,s3)
end

```

---

<div class="post-metadata">

**Author:** ![gdalle](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/gdalle/32/27854_2.png) [@gdalle](https://discourse.julialang.org/u/gdalle)\
**Post date:** [July 5, 2023, 10:52am UTC](https://discourse.julialang.org/t/faster-isapprox/101202/9 "2023-07-05T10:52:50Z")

</div>

With LoopVectorization:

```julia
julia> function is_approx3(A::Matrix{T1}, B::Matrix{T2}, rtol) where {T1,T2}
           normA = zero(T1)
           normB = zero(T2)
           normdiff = zero(promote_type(T1, T2))
           @avx for i in eachindex(A, B)
               a, b = A[i], B[i]
               normA += abs2(a)
               normA += abs2(b)
               normdiff += abs2(a - b)
           end
           return normdiff < abs2(rtol) * max(normA, normB)
       end
is_approx3 (generic function with 1 method)

julia> A, B = rand(10, 10), rand(10, 10);

julia> @btime is_approx($A, $B, $1e-4);
  147.392 ns (1 allocation: 896 bytes)

julia> @btime is_approx2($A, $B, $1e-4);
  144.150 ns (1 allocation: 896 bytes)

julia> @btime is_approx3($A, $B, $1e-4);
  17.524 ns (0 allocations: 0 bytes)

```

Note that `eachindex(A, B)` makes sure that the sizes are compatible

---

<div class="post-metadata">

**Author:** ![user\_231578](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/user_231578/32/24170_2.png) [@user\_231578](https://discourse.julialang.org/u/user_231578)\
**Post date:** [July 5, 2023, 10:54am UTC](https://discourse.julialang.org/t/faster-isapprox/101202/10 "2023-07-05T10:54:17Z")

</div>

your solution, together with `@avx`, makes it x5 faster!

```julia
function is_approx3(A::Array{Float64,2}, A_::Array{Float64,2}, rtol::Float64)
    s1, s2, s3 = 0.0, 0.0, 0.0
    @avx for i in eachindex(A)
        a, a_ = A[i], A_[i]
        s1 += (a - a_)^2
        s2 += a^2
        s3 += a_^2
    end
    return s1 < rtol^2*max(s2,s3)
end

@btime is_approx3($A,$A_,$1e-4)
  18.829 ns (0 allocations: 0 bytes)

```

---

<div class="post-metadata">

**Author:** ![DNF](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dnf/32/10191_2.png) [@DNF](https://discourse.julialang.org/u/DNF)\
**Post date:** [July 5, 2023, 2:07pm UTC](https://discourse.julialang.org/t/faster-isapprox/101202/11 "2023-07-05T14:07:43Z")

</div>

`LoopVectorization.@avx` is the old name for `LoopVectorization.@turbo`. I thought it changed several years ago. If `@avx` is still around, I suggest you use the newer, official, name, `@turbo`.

---

<div class="post-metadata">

**Author:** ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)\
**Post date:** [July 5, 2023, 2:29pm UTC](https://discourse.julialang.org/t/faster-isapprox/101202/12 "2023-07-05T14:29:11Z")

</div>

> [@user\_231578](#):
>
> I have to run (possibly) millions of times the function isapprox(A, A\_, rtol=1e-4), where A and A\_ are two matrices.

I’m curious if you can share information about your application? It’s pretty rare for `isapprox` to be a performance-critical function. Even if you’re using it as a convergence test in an iterative solver of some kind, usually there are other parts of your iteration that will dominate.

---

<div class="post-metadata">

**Author:** ![user\_231578](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/user_231578/32/24170_2.png) [@user\_231578](https://discourse.julialang.org/u/user_231578)\
**Post date:** [July 5, 2023, 2:39pm UTC](https://discourse.julialang.org/t/faster-isapprox/101202/13 "2023-07-05T14:39:11Z")

</div>

Thanks for letting me know!

---

<div class="post-metadata">

**Author:** ![user\_231578](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/user_231578/32/24170_2.png) [@user\_231578](https://discourse.julialang.org/u/user_231578)\
**Post date:** [July 5, 2023, 2:41pm UTC](https://discourse.julialang.org/t/faster-isapprox/101202/14 "2023-07-05T14:41:30Z")

</div>

> [@user\_231578](#):
>
> If you want, you can find more context about what I use this for [here](https://github.com/massimilianofurlan/ql_agent/blob/main/ql_agent.jl) .

Sure! I wrote the code example in the link specifically for this question.

---

<div class="post-metadata">

**Author:** ![Elrod](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/elrod/32/22461_2.png) [@Elrod](https://discourse.julialang.org/u/Elrod)\
**Post date:** [July 5, 2023, 2:51pm UTC](https://discourse.julialang.org/t/faster-isapprox/101202/15 "2023-07-05T14:51:14Z")

</div>

> [@DNF](#):
>
> `LoopVectorization.@avx` is the old name for `LoopVectorization.@turbo`. I thought it changed several years ago. If `@avx` is still around, I suggest you use the newer, official, name, `@turbo`.

Yes, use `@turbo`. I haven’t made a breaking release since then, which is why `@avx` still works.  
`@turbo` replaced `@avx` early in the 0.12 series. We’re now at 0.12.162.