Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

I need to make a list of some elements from a 4D array with dimensions 4x4x4x4. I need to select the elements based on their position in the following way: name an element e(x,y,z,w); I need all the elements where x≥y and z≥w. This basically excludes all the upper triangular part of each submatrix and all the submatrices in the upper triangular part of the tensor. The array in question is this one:

b={{{{12, 0, 0, 0}, {-18 Sqrt[2], 0, 0, 0}, {12 Sqrt[3], 0, 0, 0}, {-6, 
0, 0, 0}}, {{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 
0}}, {{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 
0}}, {{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 
0}}}, {{{-18 Sqrt[2], 0, 0, 0}, {84, 0, 0, 0}, {-48 Sqrt[6], 0, 0,
 0}, {54 Sqrt[2], 0, 0, 0}}, {{0, 0, 0, 0}, {0, 60, 0, 
0}, {0, -90, 0, 0}, {0, 18 Sqrt[10], 0, 0}}, {{0, 0, 0, 0}, {0, 0,
 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}}, {{0, 0, 0, 0}, {0, 0, 0, 
0}, {0, 0, 0, 0}, {0, 0, 0, 0}}}, {{{12 Sqrt[3], 0, 0, 
0}, {-48 Sqrt[6], 0, 0, 0}, {276, 0, 0, 0}, {-186 Sqrt[3], 0, 0, 
0}}, {{0, 0, 0, 0}, {0, -90, 0, 0}, {0, 240, 0, 
0}, {0, -90 Sqrt[10], 0, 0}}, {{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 
168, 0}, {0, 0, -84 Sqrt[6], 0}}, {{0, 0, 0, 0}, {0, 0, 0, 0}, {0,
 0, 0, 0}, {0, 0, 0, 0}}}, {{{-6, 0, 0, 0}, {54 Sqrt[2], 0, 0, 
0}, {-186 Sqrt[3], 0, 0, 0}, {648, 0, 0, 0}}, {{0, 0, 0, 0}, {0, 
18 Sqrt[10], 0, 0}, {0, -90 Sqrt[10], 0, 0}, {0, 600, 0, 0}}, {{0,
 0, 0, 0}, {0, 0, 0, 0}, {0, 0, -84 Sqrt[6], 0}, {0, 0, 504, 
0}}, {{0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 0}, {0, 0, 0, 360}}}}

I tried in many ways but I can't figure it out. Thanks in advance.

share|improve this question
2  
You can try MapIndexed for level {-1} in combination with Reap/Sow. –  Spawn1701D May 16 '13 at 23:06

4 Answers 4

up vote 7 down vote accepted

From your definition:

This basically excludes all the upper triangular part of each submatrix and all the submatrices in the upper triangular part of the tensor.

With[{lt = LowerTriangularize, ca = ConstantArray},
  Pick[b, # * ca[#, {4, 4}] & @ lt @ ca[1, {4, 4}], 1]
] ~Flatten~ 3
share|improve this answer
    
I think this is what Im going to use, because very soon I will have to do it for a 20x20x20x20 array and this seems pretty quick. –  user24273 May 17 '13 at 10:56

Do these give what you want?

rv1[b_] := 
  Module[{tag}, 
   SparseArray[{x_, y_, z_, w_} /; x >= y && z >= w :> 
      b[[x, y, z, w]], Dimensions@b, tag]["NonzeroValues"]];

rv2[b_] := 
  Extract[b, 
   Cases[Range@4~Tuples~{4}, {x_, y_, z_, w_} /; x >= y && z >= w]];

rv3[b_] := 
  Extract[b, 
   Table[{x, y, z, w}, {y, 4}, {x, y, 4}, {w, 4}, {z, w, 4}]~Flatten~
    3];

rv4[b_] := 
  With[{extract = 
     With[{parts = Table[{x, y}, {y, 4}, {x, y, 4}]~Flatten~1}, 
      Extract[#, parts] &]}, Flatten[extract /@ extract@b, 1]];

(*This is @Spawn1701D's suggestion*)
spawn[b_] := 
 Reap[MapIndexed[Sow[#1, #2[[1]] >= #2[[2]] && #2[[3]] >= #2[[4]]] &, 
    b, {4}], True][[-1, 1]]

(* from @MrWizard *)
wiz1[b_] := 
 With[{lt = LowerTriangularize, ca = ConstantArray}, 
   Pick[b, #*ca[#, {4, 4}] &@lt@ca[1, {4, 4}], 1]]~Flatten~3

Let's check the results

contenders = {rv1, rv2, rv3, rv4, spawn, wiz1};

Equal @@ Sort /@ Through@contenders[b]
(* True *)

Some timings

Through@(Composition @@@ 
    Tuples@{{First}, {AbsoluteTiming}, contenders})[b]

(* {0.019003, 0., 0.00100, 0., 0.002001, 0.} *)

Using this function TimingAverage that I took from somewhere sometime and never looked into

TimingAverage[expr_, time_ : 1.`] := 
 Module[{t = Timing[expr;][[1]], tries = 1}, 
     While[t < time, tries *= 2;
        t = Timing[Do[expr, {tries}];][[1]];];
    t/tries]

You get

Through@(Composition @@@ 
    Tuples@{{TimingAverage}, contenders})[b]

(* {0.0185251, 0.000716020, 0.000154249, 0.0000656987, 0.00237658, 0.000148536} *)

So, from top to bottom, rv4, wiz1, rv3, rv2, spawn, rv1

share|improve this answer
    
Thanks. I'm curious, do you personally use TimingAverage over the timeAvg function (originally from Timo) I've posted many times willfully, or was this an ad hoc function? –  Mr.Wizard May 16 '13 at 23:40
    
@Mr.Wizard because I never even looked into any of them to see pros and cons and I have it packaged up from a long time ago. If you think the other one is better, edit it in and I'll trust you and change my packaging ;) –  Rojo May 16 '13 at 23:41
    
No need; it does the same thing I believe, only I chose a multiplier of five. –  Mr.Wizard May 16 '13 at 23:48
    
Could you add timings for the method I just posted, please? –  Mr.Wizard May 16 '13 at 23:50
    
@Mr.Wizard already done. +1 btw –  Rojo May 16 '13 at 23:53

A flexible solution :

Reap[ReplacePart[
      b,
      ({x_, y_, z_, w_} /; x >= y && z >= w) :> Sow[Extract[b, {x, y, z, w}]]
    ];][[2, 1]]   

So far I know, ReplacePart[] is the only instruction that accepts a pattern of indexes as argument. That's the reason why I use ReplacePart[]. In fact we don't care of the result of the replacements.

This solution is very flexible but not efficient in terms of quickness and memory occupation.

share|improve this answer
    
thank you very much –  user24273 May 17 '13 at 10:55

Two Pick's; first picks out matrices, second elements:

With[{s = ConstantArray[True, {3, 3}] // UpperTriangularize},
  Map[
   Pick[#, s] &,
   Pick[m, s], {-3}]] // Flatten

Pick also works with SparseArray objects, so:

Pick[b, SparseArray[{x_, y_, z_, w_} /; x >= y && z >= w -> True,
   Dimensions@b]] // Flatten
share|improve this answer
    
I like this. +1 –  Mr.Wizard May 21 '13 at 13:39

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.