# Reduction on two dimensional arrays

**URL:** <https://pyfr.discourse.group/t/reduction-on-two-dimensional-arrays/1228>\
**Category:** Development\
**Created:** [4 March 2025 18:17 UTC](https://pyfr.discourse.group/t/reduction-on-two-dimensional-arrays/1228 "2025-03-04T18:17:55Z")\
**Posts on this page:** 5\
**Page:** 1

<div class="post-metadata">

**Author:** ![Frankx9](https://avatars.discourse-cdn.com/v4/letter/f/3e96dc/32.png) [@Frankx9](https://pyfr.discourse.group/u/Frankx9)\
**Post date:** [4 March 2025 18:17 UTC](https://pyfr.discourse.group/t/reduction-on-two-dimensional-arrays/1228/1 "2025-03-04T18:17:55Z")

</div>

Hi,

I would like to know if the reduction kernel could be adapted to work on arrays whose ioshape is not three-dimensional (e.g. nupts, nvars, neles) but rather bidimensional (e.g. nvars, neles) and still perform a reduction over neles to get a nvars shape as a resulting array.  
If so could you point out to the major changes that should be set for indexing the matrix appropriately?

> <https://github.com/PyFR/PyFR/blob/11785be6921777b36f332ddee471737abcc27206/pyfr/backends/cuda/blasext.py#L58C1-L66C1>

Best

---

<div class="post-metadata">

**Author:** ![fdw](https://avatars.discourse-cdn.com/v4/letter/f/a587f6/32.png) [@fdw](https://pyfr.discourse.group/u/fdw)\
**Post date:** [4 March 2025 22:12 UTC](https://pyfr.discourse.group/t/reduction-on-two-dimensional-arrays/1228/2 "2025-03-04T22:12:04Z")

</div>

The overall logic should be simple to adapt (basically just removing the inner loop which does the accumulation over the first dimension). I would suggest starting with the OpenMP backend as its reduction code is the simplest.

Regards, Freddie.

---

<div class="post-metadata">

**Author:** ![Frankx9](https://avatars.discourse-cdn.com/v4/letter/f/3e96dc/32.png) [@Frankx9](https://pyfr.discourse.group/u/Frankx9)\
**Post date:** [5 March 2025 10:23 UTC](https://pyfr.discourse.group/t/reduction-on-two-dimensional-arrays/1228/3 "2025-03-05T10:23:00Z")

</div>

Hi @fdw,

thanks for the hint, I changed the inner loop as you suggested to:

```
if (i < ncolb)
{
    ixdtype_t idx = SOA_IX(i, blockIdx.y, gridDim.y);
    r = data[idx];
    acc += r;
}

```

but when I test it with a sum of all elements I get this as results

```
    test = np.ones((5, 49055)).astype(np.float32)
    test= self._be.matrix(test, initval= test, tags={'align'})
    red = self._be.kernel('sum', test)
    self._be.run_kernels([red], wait=True)

```

with a result which is

`[49054. 49054. 49054. 49054. 49055.]`

Is there anything else that should be changed?

Best

---

<div class="post-metadata">

**Author:** ![fdw](https://avatars.discourse-cdn.com/v4/letter/f/a587f6/32.png) [@fdw](https://pyfr.discourse.group/u/fdw)\
**Post date:** [5 March 2025 12:49 UTC](https://pyfr.discourse.group/t/reduction-on-two-dimensional-arrays/1228/4 "2025-03-05T12:49:37Z")

</div>

As per my suggestion, start with the OpenMP backend first for testing.

Regards, Freddie.

---

<div class="post-metadata">

**Author:** ![Frankx9](https://avatars.discourse-cdn.com/v4/letter/f/3e96dc/32.png) [@Frankx9](https://pyfr.discourse.group/u/Frankx9)\
**Post date:** [6 March 2025 14:58 UTC](https://pyfr.discourse.group/t/reduction-on-two-dimensional-arrays/1228/5 "2025-03-06T14:58:56Z")

</div>

> (basically just removing the inner loop which does the accumulation over the first dimension)

In the openMP kernel is it

```auto
                #pragma omp simd
                for (ixdtype_t _xj = 0; _xj < SOA_SZ; _xj++)
                {
                 ...
                }

```

?

Best