What Reap/Sow tricks do you use?

Question

I admit that I have been reluctant to use Reap and Sow for years and that I used it only when I had no choice with EvaluationMonitor.

That is partly because the documentation was difficult for me and partly because I could not think of any scenario where that would be useful. My impression is that I am not the only one that avoids Reap/Sow constructs.

However, I am noticing that it can be useful. Here are some examples:

Sub expressions

Obtaining a list of sub-expressions :

Reap[expression /. 
  x : pattern :> Sow[x]]

Depending on the pattern that might differ from

Cases[expression,pattern,All]

Example:

Cases[{x, x^2, m[x]*m[y]*e}, m[a_]*m[b_], All]

{}

Reap[{x, x^2, m[x]*m[y]*e} /. x : m[a_]*m[b_] :> Sow[x]][[2, 1]]

{m[x] m[y]}

Substitution for AppendTo.

Compare :

l = {};
Do[AppendTo[l, RandomReal[]], 10^5] // AbsoluteTiming

10.8 seconds

with

l = {};
Reap[Do[Sow[RandomReal[]], 10^5]][[2, 1]]; // AbsoluteTiming

0.06 seconds

AppendTo is slow because it copies the list each time. See the following links for further details on AppendTo:

https://mathematica.stackexchange.com/a/72625/86543

and point 3 and 5 in this answer:

https://mathematica.stackexchange.com/a/29351/86543

See also the following for a way to reproduce and enhance Reap/Sow:

https://mathematica.stackexchange.com/a/244033/86543

EvaluationMonitor or StepMonitor

Nothing much to say other than check out the documentation but I am guilty of not using this to check why FindRoot or Plot is not converging or is slow.

Debugging

I use Echo usually but if one needs to extract a lot of intermediate parts in a code and check how they relate to each other using another code, then copy pasting Echo outputs would be overly complicated.

Also, Reap and Sow could be helpful to check what Map, MapAt, MapAll, SubsetMap is working on especially if one included a level specification. That said, I prefer taking a small example and using Framed

Finding specific terms using functions that take a function as an argument

For example one may have wrote a code that acts only on specific terms using a function argument. A quick way to extract those specific terms is to use Reap and Sow.

Arbitrary depth and position categorization

One can use Reap/Sow with tags in combination with MapAt, MapAll,MapBlock, MapIndexed, ArrayReduce etc to generalize Select and Gather to different depths and positions.

Example 1:

Gather sub-expressions of an expression by LeafCount:

Reap[MapAll[Sow[#, LeafCount[#]] &, 
    Integrate[1/(1 + Tan[x]^20), x]]][[2]] // Map[DeleteDuplicates]

Example 2:

Gather and select columns of a matrix depending on whether the sum of elements of the column is prime:

(
Reap[
    ArrayReduce[
                Sow[{#,PrimeQ@Total@#}
                    ,
                    PrimeQ[Total@#]
                   ]&
        ,
        EchoFunction[MatrixForm]@
        Array[RandomInteger[40]&,{3,4}]
        ,
        1
    ]
][[2]]
//Map[MatrixForm,#,{2,Infinity}]&
)

The matrix:

$$ \left( \begin{array}{cccc} 13 & 13 & 25 & 28 \\ 21 & 26 & 7 & 28 \\ 9 & 35 & 21 & 4 \\ \end{array} \right) $$

The gathered columns:

$$\left( \begin{array}{cc} \{\{13,21,9\},\text{True}\} & \{\{25,7,21\},\text{True}\} \\ \{\{13,26,35\},\text{False}\} & \{\{28,28,4\},\text{False}\} \\ \end{array} \right)$$

If we want to select only the column where the sum of elements is prime than we can add an extra argument True at the end in the Reap command.

I should also mention SelectEquivalents from this answer:

https://stackoverflow.com/a/6245166

but I am not sure how to use it in practice.

---

What else ?

Answer 1

Ok, I'll bite. Here is one example where I found Reap/Sow useful. It's also an excuse to document my Backtrack function that I find useful but hasn't seen many example of Wolfram Language code implementing it in a most generic form.

---

Backtracking is a very useful technique for solving combinatorial problems, and Donald Knuth dedicated a substantial part of TAOCP 4A to it. Let's define a simple backtrack solver:

Backtrack[start_, next_, visitQ_, visit_] := Module[{
   level = 1, (* current search level *)
   stack = {{start}}, (* list of states to try at each level *)
   current = {1}, (* 
   list of indexes of states that need to be tried at each level*)
   state (* current state *)},
  While[level > 0,
   If[current[[level]] > Length[stack[[level]]], (* 
    We must backtrack *)
    level--,
    state = stack[[level, current[[level]]]];
    current[[level]]++;
    If[visitQ@state, visit@state];
    With[{nextStates = next@state},
     If[
      Length@nextStates > 0, (* We must go deeper one level *)
      level++;
      If[Length@stack < level, stack = Append[stack, nextStates], 
       stack[[level]] = nextStates];
      If[Length@current < level, current = Append[current, 1], 
       current[[level]] = 1]
      ]]]]]

The function Backtrack uses the following parameters:

• start - the initial state
• next - a function that, given the current state, produces a list of next states
• visitQ - a function that checks whether the current state should be visited
• visit - a function that visits the current state

It's the most generic version, and we can define more specialized versions with more convenient interface by modifying the visit function. Notice the use of Reap/Sow in BacktrackCollect - I find it idiomatic to use it for data collection, esp. given very non-linear nature of backtracking search.

BacktrackCollect[start_, next_, visitQ_, n_ : -1] := 
 Module[{count = 0, r},
  r = Reap@
    Catch@Backtrack[start, next, 
      visitQ, (Sow[#]; count++; If[count == n, Throw["Enough"]]) &];
  If[Length@r[[2]] > 0, r[[2, 1]], {}]]
BacktrackCount[start_, next_, visitQ_] := 
 Module[{count = 0}, Backtrack[start, next, visitQ, count++ &]; count]

---

Now let's use our freshly minted BacktrackCollect to solve n-queens problem:

queenProblem[n_] := Module[{next, visitQ},
  next[queens_] := Map[
    Append[queens, #] &,
    Select[
     Complement[Range[n], queens],
     q |-> Module[{m = Length@queens, i},
       AllTrue[Range[m], Abs[queens[[#]] - q] != m + 1 - # &]
       ]]];
  visitQ[queens_] := Length@queens == n;
  {next, visitQ}]
{next, visitQ} = queenProblem[5];
BacktrackCount[{}, next, visitQ] // Timing (* {0.003299, 10} *)
displayQueens[queens_] := Module[{n = Length[queens], board, i},
  board = ConstantArray[" ", {n, n}];
  For[i = 1, i <= n, i++,
   board[[queens[[i]], i]] = [BlackQueen]];
  Grid[board, Frame -> All, 
   Background -> {Automatic, Automatic, 
     Flatten[Table[{i, j} -> 
        If[EvenQ[i + j], Darker[White], White], {i, n}, {j, n}]]}]]
Multicolumn[displayQueens /@ BacktrackCollect[{}, next, visitQ], 5]