Pythagorean Triple finder
Here are some comments on your revised code, many of which apply to your original code as well.

for (auto&amp; T: U) transform(T).find(max_perim, proc);


You're using the names T and U to refer to things that are not template type parameters. In fact, U is a reference to a static data member of the enclosing class, whose only declaration appears many lines below this use. This is extremely surprising code. I might instead write

for (auto&amp;&amp; formula : this-&gt;child_formulas) {
    this-&gt;transform(formula).find(max_perim, proc);
}


Notice the use of braces around the loop body; the use of auto&amp;&amp;, Which Always Works; the use of this-&gt; to clarify for the human reader what scope we're expecting to find these names in; and especially the descriptive name child_formulas for U. (child_formulas says what U is; this-&gt; says whose it is. Both of those things were in question, in the old code.)



class triple : public std::array&lt;long, 3&gt; {


Don't inherit from standard types. It's legal C++, but it's bad practice.


  If you (or your project, or your company) didn’t write class Foo, then class Foo should not be granted control over the API of your own class. And that’s what you’re doing when you inherit from a base class: you’re granting that class control over your API.


I like your original x, y, z data members much better.



template &lt;typename F&gt;
void find(long max_perim, F&amp;&amp; proc) {


You're passing proc by perfect-forwarding, but then when you use it, you're using it as a plain old lvalue. This is fine (in the wild-west anything-goes world of templates), but it doesn't express the meaning of proc quite as well as I'd like. proc is supposed to be a callback that gets called for each triple, right? Calling a callback doesn't modify the callback. So, we can and should require that proc be const-callable, and then we just write

template&lt;class F&gt;
void find(int max_perim, const F&amp; proc) {


(Drive-by irrelevant style adjustments to save typing. long varies in size and is the same size as int on many platforms, so let's just use int until we're ready to step all the way up to a well-defined int64_t or __int128_t.)

If you really want to support non-const and/or rvalue-callable Fs, you can perfect-forward proc all over the place, but trust me, it's not worth it. (The std::forward&lt;F&gt;(proc)s will clutter your code, and all you're doing is enabling your users to write confusing and unintuitive client code.)

Vice versa, re your question about concepts: You can constrain this template to SFINAE away in C++17 (or C++11) like this:

template&lt;class F, class = std::enable_if_t&lt;std::is_invocable_v&lt;const F&amp;, triple&amp;&gt;&gt;&gt;
void find(int max_perim, const F&amp; proc) {


This is not much more boilerplate than the C++2a concepts version:

template&lt;class F&gt; requires std::is_invocable&lt;const F&amp;, triple&amp;&gt;
void find(int max_perim, const F&amp; proc) {


Notice that in C++2a you will also be able to write "simply"

void find(int max_perim, const std::invocable&lt;triple&amp;&gt; auto&amp; proc) {


but (A) this doesn't express exactly the same SFINAE-constraint as the other two, and (B) IMHO it looks more confusing and scary.

But should you constrain this function to SFINAE away when it's given a "bad" lambda type? IMHO, no, you shouldn't. There's no reason to SFINAE here. SFINAE is for when you need this version of find to drop out of the overload set so that some other more general version can pick up the call. That's not the situation that you have, here.

Compare the error messages you get from (A) (or the C++2a Concepts version (A2))

template&lt;class F, class = std::enable_if_t&lt;std::is_invocable_v&lt;const F&amp;, triple&amp;&gt;&gt;&gt;
void find(int max_perim, const F&amp; proc) {
    proc(*this);
}

this-&gt;find(42, 7);


versus (B)

template&lt;class F&gt;
void find(int max_perim, const F&amp; proc) {
    static_assert(std::is_invocable_v&lt;const F&amp;, triple&amp;&gt;);
    proc(*this);
}

this-&gt;find(42, 7);


versus (C)

template&lt;class F&gt;
void find(int max_perim, const F&amp; proc) {
    proc(*this);
}

this-&gt;find(42, 7);


and see which version feels most "user-friendly" for the client programmer. (Remember that someone using find correctly will never see any of these error messages, so you should optimize to help the guy who doesn't know how to use it.)



  cnt_next_depth += U.size(); // always 3
  if (--cnt_this_depth == 0) {
    if (++depth &gt; max_depth) return;
    cnt_this_depth = cnt_next_depth;
    cnt_next_depth = 0;
  }
  for (auto&amp; T: U) q.push(t.transform(T));


Code comments are usually a red flag, at least on CodeReview. ;) It seems that this += corresponds to the three calls to q.push below; i.e., your algorithm requires that cnt_next_depth exactly track q.size(). But you can't use q.size() because you are reusing q to store both the elements at the current level and the elements at the next level.

It would make more sense to use two different queues:

std::vector&lt;triple&gt; this_level = ...;
std::vector&lt;triple&gt; next_level;
for (int i=0; i &lt; max_depth; ++i) {
    for (auto&amp;&amp; t : this_level) {
        proc(t);
        for (auto&amp;&amp; formula : child_formulas) {
            next_level.push_back(t.transform(formula));
        }
    }
    this_level = std::move(next_level);
    next_level.clear();
}


As a bonus, it turns out that we don't even need them to be queues anymore; they can be plain old vectors, and we save a bunch of heap traffic.



friend std::ostream&amp; operator&lt;&lt;(std::ostream&amp; stream, const triple&amp; t) noexcept {
  // Frustrating boiler plate. Any terse alternatives, that do it quickly and correctly?
  char comma[] = {'\0', ' ', '\0'}; // NOLINT
  for (auto d: t) {
    stream &lt;&lt; comma &lt;&lt; d;
    comma[0] = ',';
  }
  return stream;
}


No good answer. The idiomatic way would be

friend std::ostream&amp; operator&lt;&lt;(std::ostream&amp; stream, const triple&amp; t) {
    bool first = true;
    for (auto&amp;&amp; d : t) {
        if (!first) stream &lt;&lt; ", ";
        stream &lt;&lt; d;
        first = false;
    }
    return stream;
}


Notice that this function is not noexcept, because any of these stream operations might throw. (For example, if you're writing to a stringstream and it throws bad_alloc.)

Also, in real life I'd write for (int d : t) to indicate exactly what type of "elements" I expect to be getting out of a triple. I'm using auto&amp;&amp; here only because I saw your comment that you want to stick with "Almost Always Auto" style.



  triple transform(const trans&amp; T) const noexcept {
    auto res = triple{};
    std::transform(T.begin(), T.end(), res.begin(), [this](vec3 V) {
      return std::inner_product(V.begin(), V.end(), this-&gt;begin(), 0L);
    });
    return res;
  }


This is interesting use of STL algorithms, but algorithms are really meant for operating on big anonymous ranges of data elements, not for constant-size situations like we have here. In order to use STL algorithms, you've been forced to anonymize your formulas into faceless data ranges:

  static constexpr auto U = std::array&lt;trans, 3&gt;{{
      // https://en.wikipedia.org/wiki/Pythagorean_triple#Parent.2Fchild_relationships
      {{{{{1, -2, 2}},    // vec3 U[0][0]
         {{2, -1, 2}},    // vec3 U[0][1]
         {{2, -2, 3}}}}}, // vec3 U[0][1]

      {{{{{1, 2, 2}},    // vec3 U[1][0]
         {{2, 1, 2}},    // vec3 U[1][1]
         {{2, 2, 3}}}}}, // vec3 U[1][2]

      {{{{{-1, 2, 2}},    // vec3 U[2][0]
         {{-2, 1, 2}},    // vec3 U[2][1]
         {{-2, 2, 3}}}}}, // vec3 U[2][2]
  }};


Compare that code to a more "code-driven" version:

triple transform(const Formula&amp; f) const noexcept {
    return f(x, y, z);
}

auto formulas[] = {
    +[](int x, int y, int z){ return triple{ x - 2*y + 2*z, 2*x - y + 2*z, 2*x - 2*y + 3*z}; },
    +[](int x, int y, int z){ return triple{ x + 2*y + 2*z, 2*x + y + 2*z, 2*x + 2*y + 3*z}; },
    +[](int x, int y, int z){ return triple{ -x + 2*y + 2*z, -2*x + y + 2*z, -2*x + 2*y + 3*z}; },
};


In this version, trans is gone already, and t.transform(f) is hardly pulling its weight.



if (perimeter() &gt; max_perim) return;


Consider what you'd do if you wanted to find all triples up to some maximum side length, or up to some maximum hypotenuse length. Does this approach generalize to

triple.find(callback_on_each_triple, stopping_condition);


?