# Improve performance for multiple For-loops

**URL:** <https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289>\
**Category:** Performance\
**Tags:** geo, io\
**Created:** [January 13, 2021, 4:17pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289 "2021-01-13T16:17:40Z")\
**Posts on this page:** 20\
**Page:** 1

<div class="post-metadata">

**Author:** ![EuRoXy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/euroxy/32/13970_2.png) [@EuRoXy](https://discourse.julialang.org/u/EuRoXy)\
**Post date:** [January 13, 2021, 4:17pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/1 "2021-01-13T16:17:40Z")

</div>

I am trying to optimize a calculation of correlation on a geo-grid of size 16x16:

```julia
using NetCDF
function corr1vall(s, lon1, lat1) # scenario, longitude 1, latitude 1
    pathPv = myPath * "pv/" * ssps[s] * "/day/EU/regr/"
    rmcp(pathPv)
    
    corr = zeros(16,16,28)
    for m in 1:28 # model
        fnam = pathPv * readdir(pathPv)[m]
        pv = ncread(fnam, "pv") # size (16,16,7300)

        @views pv1 = pv[lon1,lat1,:] # 20-yr time series of grid 1

        for lat in 1:16
            for lon in 1:16
                @views pv2 = pv[lon,lat,:] # 20-yr time series of grid 2
                corr[lon,lat,m] = cor(pv1/areaLat[1], pv2/areaLat[lat]) # normalize
            end
        end
    end
    return mean(corr, dims=3)[:,:,1] # multi-model mean
end

```

One grid takes currently 6.6 s:

```julia
@time corr1vall(1,1,1)
> 6.625234 seconds (89.95 k allocations: 1.171 GiB, 0.90% gc time)

```

Further I have 256 grids and 3 scenarios:

```julia
corrAllGrids = zeros(16,16,16*16)
g = 1 
for lat in 1:16
    for lon in 1:16
        corrAllGrids[:,:,g] = corr1vall(1, lon, lat)
        g += 1
    end
end

```

I hope to reduce the time for the first function before trying parallelizing, is there evident modifications I could make?

---

<div class="post-metadata">

**Author:** ![Oscar\_Smith](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/oscar_smith/32/25343_2.png) [@Oscar\_Smith](https://discourse.julialang.org/u/Oscar_Smith)\
**Post date:** [January 13, 2021, 4:23pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/2 "2021-01-13T16:23:31Z")

</div>

Replacing ` cor(pv1/areaLat[1], pv2/areaLat[lat])` with ` cor(pv1*(1/areaLat[1]), pv2*(1/areaLat[lat]))` should be a major improvement. Division is about 6x slower than multiplication.

---

<div class="post-metadata">

**Author:** ![EuRoXy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/euroxy/32/13970_2.png) [@EuRoXy](https://discourse.julialang.org/u/EuRoXy)\
**Post date:** [January 13, 2021, 4:37pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/3 "2021-01-13T16:37:54Z")

</div>

Thanks for the hint. However, I still get 6.6-6.7 s after this.

---

<div class="post-metadata">

**Author:** ![cchderrick](https://avatars.discourse-cdn.com/v4/letter/c/ecd19e/32.png) [@cchderrick](https://discourse.julialang.org/u/cchderrick)\
**Post date:** [January 13, 2021, 4:38pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/4 "2021-01-13T16:38:40Z")

</div>

Are you sure it is the loop has performance issue, not an IO issue?

> [@EuRoXy](#):
>
> ` pv = ncread(fnam, "pv") # size (16,16,7300)`

Does this allocate a 1,868,800 element 28 times? I think there is an `read!` so you could allocate only once.  
([High-level interface · NetCDF.jl](https://juliageo.org/NetCDF.jl/dev/highlevel/#NetCDF.ncread))!

---

<div class="post-metadata">

**Author:** ![Skoffer](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/skoffer/32/378_2.png) [@Skoffer](https://discourse.julialang.org/u/Skoffer)\
**Post date:** [January 13, 2021, 4:53pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/5 "2021-01-13T16:53:39Z")

</div>

You should check for type stability. It looks like `areaLat` and `myPath` are defined outside of the function (global variables?)

---

<div class="post-metadata">

**Author:** ![EuRoXy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/euroxy/32/13970_2.png) [@EuRoXy](https://discourse.julialang.org/u/EuRoXy)\
**Post date:** [January 13, 2021, 4:56pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/6 "2021-01-13T16:56:30Z")

</div>

One problem is that not every model has the same time span: some 7300 days, some 7305, etc. Is it still feasible to allocate only once?

---

<div class="post-metadata">

**Author:** ![EuRoXy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/euroxy/32/13970_2.png) [@EuRoXy](https://discourse.julialang.org/u/EuRoXy)\
**Post date:** [January 13, 2021, 4:58pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/7 "2021-01-13T16:58:25Z")

</div>

Right, I declared area weights `areaLat` and `myPath` globally.

---

<div class="post-metadata">

**Author:** ![Skoffer](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/skoffer/32/378_2.png) [@Skoffer](https://discourse.julialang.org/u/Skoffer)\
**Post date:** [January 13, 2021, 4:59pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/8 "2021-01-13T16:59:25Z")

</div>

Also, I barely remember statistics, but I guess `cor(a*c1, b*c2)` can be rewritten somehow, maybe `c1*c2*cor(a, b)` or something like that. If it is the case, constants should be factored out, since each `pv1/areaLat[1]` allocate new array.

---

<div class="post-metadata">

**Author:** ![Skoffer](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/skoffer/32/378_2.png) [@Skoffer](https://discourse.julialang.org/u/Skoffer)\
**Post date:** [January 13, 2021, 5:03pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/9 "2021-01-13T17:03:12Z")

</div>

You can avoid by declaring

```julia
function corr1vall(s, lon1, lat1, pathPv, areaLat)
...
end

```

It increases the number of arguments but makes everything type stable. If a large number of parameters is not what you want, you can wrap extra parameters in named tuple and use `UnPack.jl` to extract them in the function body

```julia
opts = (; pathPv = pathPv, areaLat = areaLat)
corr1vall(s, lon1, lat1, opts)
   @unpack pathPv, areaLat = opts
....
end

```

---

<div class="post-metadata">

**Author:** ![baggepinnen](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/baggepinnen/32/693_2.png) [@baggepinnen](https://discourse.julialang.org/u/baggepinnen)\
**Post date:** [January 13, 2021, 5:20pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/10 "2021-01-13T17:20:25Z")

</div>

> [@Skoffer](#):
>
> c1_c2_cor(a, b)

This is not quite correct, I think you just remove the multipliers, they don’t change the correlation except for the sign?

---

<div class="post-metadata">

**Author:** ![EuRoXy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/euroxy/32/13970_2.png) [@EuRoXy](https://discourse.julialang.org/u/EuRoXy)\
**Post date:** [January 13, 2021, 5:31pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/11 "2021-01-13T17:31:12Z")

</div>

Based on the suggestions so far, I modified to

```julia
function corr1vall(s, lon1, lat1, pathPv, areaLat) # scenario, longitude 1, latitude 1
    pathPv = myPath * "pv/" * ssps[s] * "/day/EU/regr/"
    rmcp(pathPv)
    
    corr = zeros(16,16,28)
    pv = zeros()
    for m in 1:28 # model
        fnam = pathPv * readdir(pathPv)[m]
        pv = ncread(fnam, "pv") # size (16,16,7300 or 7305 or 7160)

        @views pv1 = pv[lon1,lat1,:] * (1/areaLat[lat1]) # normalize

        for lat in 1:16
            for lon in 1:16
                @views pv2 = pv[lon,lat,:] * (1/areaLat[lat])
                corr[lon,lat,m] = cor(pv1, pv2)
            end
        end
    end
    return mean(corr, dims=3)[:,:,1]
end

```

Current runtime for one grid:

```julia
@time corr1vall(1,1,1, pathPv, areaLat)
> 6.357588 seconds (143.15 k allocations: 805.057 MiB, 0.60% gc time)

```

But I am still not sure about how to deal with I/O issue from `ncread`, since the time dimension of each model is not identical.

---

<div class="post-metadata">

**Author:** ![baggepinnen](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/baggepinnen/32/693_2.png) [@baggepinnen](https://discourse.julialang.org/u/baggepinnen)\
**Post date:** [January 13, 2021, 5:35pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/12 "2021-01-13T17:35:16Z")

</div>

How long time does each read take, the calculations don’t seem that expensive, while reading data can take quite some time.

---

<div class="post-metadata">

**Author:** ![EuRoXy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/euroxy/32/13970_2.png) [@EuRoXy](https://discourse.julialang.org/u/EuRoXy)\
**Post date:** [January 13, 2021, 5:38pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/14 "2021-01-13T17:38:30Z")

</div>

Reading each m takes around 0.20-0.22 s

---

<div class="post-metadata">

**Author:** ![cchderrick](https://avatars.discourse-cdn.com/v4/letter/c/ecd19e/32.png) [@cchderrick](https://discourse.julialang.org/u/cchderrick)\
**Post date:** [January 13, 2021, 5:42pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/15 "2021-01-13T17:42:43Z")

</div>

0.22 sec \* 28 models = 6.16 sec total for `corr1vall`, looks like a IO time to me.

---

<div class="post-metadata">

**Author:** ![Skoffer](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/skoffer/32/378_2.png) [@Skoffer](https://discourse.julialang.org/u/Skoffer)\
**Post date:** [January 13, 2021, 5:54pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/16 "2021-01-13T17:54:57Z")

</div>

You are right, instead of

```julia
cor(pv1/areaLat[1], pv2/areaLat[lat])

```

it should be just

```julia
sign(areaLat[1]) * sign(areaLat[lat]) * cor(pv1, pv2)

```

---

<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:** [January 13, 2021, 5:58pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/17 "2021-01-13T17:58:43Z")

</div>

I was going to say, don’t use global variables, move this

> [@EuRoXy](#):
>
> `readdir(pathPv)[m]`

outside the loop, use `dropdims` instead of this

> [@EuRoXy](#):
>
> `mean(corr, dims=3)[:,:,1] `

which makes an unnecessary copy.

But clearly, it’s all about file io.

(But also, making paths by string concatenation is not robust, use `joinpath`)

---

<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:** [January 13, 2021, 6:08pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/18 "2021-01-13T18:08:30Z")

</div>

> [@EuRoXy](#):
>
> `pv[lon,lat,:]`

You seem to be working along the last dimension, which is not advantageous. And that makes me wonder: how are data actually stored on file?

Is there some chance that those data are contiguous on file, but are reorganized during reading, perhaps to mimic how they work in python or C? That could cause the file reading itself to become slow, and later also hampers the correlation calculations.

---

<div class="post-metadata">

**Author:** ![EuRoXy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/euroxy/32/13970_2.png) [@EuRoXy](https://discourse.julialang.org/u/EuRoXy)\
**Post date:** [January 13, 2021, 6:16pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/19 "2021-01-13T18:16:30Z")

</div>

@fabiangans

---

<div class="post-metadata">

**Author:** ![fabiangans](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fabiangans/32/2624_2.png) [@fabiangans](https://discourse.julialang.org/u/fabiangans)\
**Post date:** [January 14, 2021, 6:36am UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/21 "2021-01-14T06:36:14Z")

</div>

FIrst of all, did you profile this to find out where the time is spent? In case all the time is spent in `ncread`, there is not a lot you can do currently, however, 6s read time still is unusually high so it would be worth filing an issue and share an example file if you can confirm that time is spent in IO (just time reading a single file).

Another thing to remember is that `ncread` is and can not be type-stable because you never know what the data types inside a netcdf file are. So you have to make sure that the type of pv is known in some way to make an efficient inner loop. This can be done either by using a function barrier as `NetCDF.open(function, filename, varname)` provides or by using the mutating form of `ncread!` and pre-allocating the array as suggested already, for example like this:

```julia
corr = zeros(16,16,28)
pv = zeros(16,16,7300)

    for m in 1:28 # model
        fnam = pathPv * readdir(pathPv)[m]
        ncread!(fnam, "pv", pv) # size (16,16,7300)

```

This would reduce allocations and make it easier for the compiler to predict types in your second loop. However, as mentioned this only makes sense to change if not all time is spent inside `ncread`.

---

<div class="post-metadata">

**Author:** ![EuRoXy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/euroxy/32/13970_2.png) [@EuRoXy](https://discourse.julialang.org/u/EuRoXy)\
**Post date:** [January 14, 2021, 4:09pm UTC](https://discourse.julialang.org/t/improve-performance-for-multiple-for-loops/53289/23 "2021-01-14T16:09:40Z")

</div>

Thanks a lot for all your helpful insights! After moving reading .nc file to the outermost loop, I get an acceptable computation time.