1

I am box filtering data from a DNS simulation. When filtering the original data read from a file, this is fast. When filtering a function of these data after applying Map, it is 10 times slower. If I write the data after Map operation to a file and read in again, the filtering is fast again.

What is the difference of data after Map operation to simple numbers? How can I convert the result to numbers without writing / reading to a file?

Here is a test example:


    idiff = 100; jdiff = 100; kdiff = 400;
cfct3Dr = Array[0.5 &, {idiff + 1, jdiff + 1, kdiff + 1}];

omcc[c_, m_] = -c^(1 + m) (-1 + c^m) (1 + m) // Simplify

m0 = 49/11

omm0[c_] = omcc[c, m0]


omDNScalc = ParallelMap[omm0[#] &, cfct3Dr, {3}];

Export["omDNSdata.h5", omDNScalc];
omDNS = Import["omDNSdata.h5", "Dataset1"];
Delta = 32;shift = 15;niLES =4;njLES=4;nkLES=24;

For[omDNSLES = Array[0.*# &, {niLES, njLES, nkLES}]; i = 1, 
  i <= niLES, i++, 
  For[j = 1, j <= njLES, j++, 
   For[k = 1, k <= nkLES, k++, 
    omDNSLES[[i]][[j]][[k]] = 
     Total[Total[
        Total[omDNS[[shift*(i - 1) + 1 ;; shift*(i - 1) + Delta, 
          shift*(j - 1) + 1 ;; shift*(j - 1) + Delta, 
          shift*(k - 1) + 1 ;; 
           shift*(k - 1) + Delta]]]]]/(Delta^3)]]] // AbsoluteTiming
{0.0406014, Null}

For[omDNSLEScalc = Array[0.*# &, {niLES, njLES, nkLES}]; i = 1, 
  i <= niLES, i++, 
  For[j = 1, j <= njLES, j++, 
   For[k = 1, k <= nkLES, k++, 
    omDNSLEScalc[[i]][[j]][[k]] = 
     Total[Total[
        Total[omDNScalc[[shift*(i - 1) + 1 ;; shift*(i - 1) + Delta, 
           shift*(j - 1) + 1 ;; shift*(j - 1) + Delta, 
           shift*(k - 1) + 1 ;; 
            shift*(k - 1) + Delta]]]]]/(Delta^3)]]] // AbsoluteTiming
{1.41565, Null}

In this example, niLES=4,njLES=4,nkLES=27,shift=12,Delta=32, cfct3Dr and cfomDNS are 3D (100x100x340) arrays, cfct3Dr is read in from a .h5 file. The second loop need 10 times more CPU time, I could not find a difference when inspecting omDNS and omDNScalc.

Also: Is there a more elegant way to box filter data from a bigger 3D array, maybe using Map instead of the nested For (or Do) lists?

Michael E2
  • 235,386
  • 17
  • 334
  • 747
  • 1
    Please make your code self-contained. There are several undefined variables. – Roman Dec 19 '20 at 16:20
  • 2
    omm0[c] is a syntax error. It should read: omm0[c_]. Therefore omDNScalc will not work as expected. – Daniel Huber Dec 19 '20 at 16:39
  • Your code does not run for me: https://i.stack.imgur.com/CgegL.png -- It looks like you're just constructing Table[] so that seems better that For[]. Your three nested Total[] seem like the same as Total[..., 3]. These are minor compared to Daniel's remark. – Michael E2 Dec 19 '20 at 18:48

1 Answers1

3

The reason for the difference in speed is that omDNScalc is not a packed array. This is a faster way to compute it that constructs a packed array:

omDNScalc = omm0@cfct3Dr;

With that change, the times are comparable and fast.


Alternative: I had to experiment a little to figure out you don't sum over the whole array omDNS (that is, idiff etc. are not evenly covered by shift and Delta). That is the reason we don’t take all of omDNS below.

See Developer`PartitionMap for a useful function and Partition of which it is a generalization. The code is no faster (except for omm0@cfct3Dr), but it's shorter.

idiff = 100; jdiff = 100; kdiff = 400;
cfct3Dr = ConstantArray[0.5, {idiff + 1, jdiff + 1, kdiff + 1}];
omcc[c_, m_] = -c^(1 + m) (-1 + c^m) (1 + m) // Simplify;
m0 = 49/11;
omm0[c_] = omcc[c, m0];
omDNScalc = omm0@cfct3Dr;
omDNS = omDNScalc;
Delta = 32; shift = 15; niLES = 4; njLES = 4; nkLES = 24;

omDNSLES = Developer`PartitionMap[ Total[#, Infinity] &, omDNS[[;; -shift, ;; -shift, ;; -shift]], {Delta, Delta, Delta}, {shift, shift, shift}] / Delta^3;

Michael E2
  • 235,386
  • 17
  • 334
  • 747