vector/matrix ctors for var

[Bob Carpenter]

It's a bit different for vars. With these three functions:

|
double fun1(const Matrix<var,Dynamic,1>& x) {
Matrix<var,Dynamic,1> y(x);
return y(2).val();
}
double fun2(const Matrix<var,Dynamic,1>& x) {
Matrix<var,Dynamic,1> y(x.size());
std::memcpy(&y(0),&x(0),x.size() * sizeof(double));
return y(2).val();
}
double fun3(const Matrix<var,Dynamic,1>& x) {
std::vector<var> y(&x(0),&x(0) + x.size());
return y[2].val();
}
|

the proportion of time taken is:

fun1: 46%
fun2: 30%
fun3: 24%

Oddly, memcpy is only taking about 20% of the time
spent in fun2. So I can't figure out why fun3 isn't
a whole lot faster.

This is all using clang++ -O3, by the way, with the
Stan libs also compiled with the same flags.

And you can see it's using the sse42 move
functionality.

I'm having a bit of trouble understanding the profiler output.
I'm only guessing at fun3 time because it doesn't show up at
all --- it's just that fun1 + fun2 take 76% of the time.
The other time seems to be missing. I'm attaching a screen
shot to show what I mean.

- Bob

[Bob Carpenter]

OK, I think I sorted out what's going on with
these by selecting the "invert call tree" checkbox
under "Call Tree" on the left and then sorting
by Self time. It's still confusing that it lists
some of the time under main() and some under each
of the functions.

So I just went and instrumented the code to make sure
I was interpreting the profile correctly. I get:

fun1: 27s Eigen vector copy ctor
fun2: 17s Eigen vector size ctor + memcpy
fun3: 13s std::vector iterator ctor

for the following:

> double fun1(const Matrix<var,Dynamic,1>& x) {
> Matrix<var,Dynamic,1> y(x);
> return y(2).val();
> }
> double fun2(const Matrix<var,Dynamic,1>& x) {
> Matrix<var,Dynamic,1> y(x.size());
> std::memcpy(&y(0),&x(0),x.size() * sizeof(double));
> return y(2).val();
> }
> double fun3(const Matrix<var,Dynamic,1>& x) {
> std::vector<var> y(&x(0),&x(0) + x.size());
> return y[2].val();
> }
> |

I also checked and the size ctor in fun2() only allocates
the memory, it doesn't actually fill in values. So we
don't get extraneous values on the autodiff stack using var.

So if copy ctors for Eigen vectors are taking up too much time,
I'd just rewrite them as in fun2. It's still not as fast
as std::vector, but there's no Eigen Matrix ctor that takes in
either the memory for the vector elements or an iterator.

- Bob

[Ben Goodrich]

On Wednesday, May 29, 2013 5:22:10 PM UTC-4, Bob Carpenter wrote:

So if copy ctors for Eigen vectors are taking up too much time,
I'd just rewrite them as in fun2.  It's still not as fast
as std::vector, but there's no Eigen Matrix ctor that takes in
either the memory for the vector elements or an iterator.

This whole thing strikes me as odd, so maybe it would be good to ask on the Eigen list? But since we are utilizing stan/math/matrix/EigenDenseBaseAddons.hpp anyway, if it is more convenient for us to write additional Eigen::Matrix constructors, then we could go ahead and do that.

Ben

[Michael Betancourt]

Wow, memcpy absolutely smokes it.  Each run constructed the given object 50000 times.

Eigen::VectorXd temp(copy.size());

std::memcpy(&temp(0),&copy(0), copy.size() * sizeof(double));

     2.7e-05 (ms)

Eigen::VectorXd temp(copy);

     0.004721 (ms)

std::vector<double> temp(cont_params);

     3.9e-05 (ms)

I'll revamp the code as soon as I can.

--
You received this message because you are subscribed to the Google Groups "stan development mailing list" group.
To unsubscribe from this group and stop receiving emails from it, send an email to stan-dev+u...@googlegroups.com.
For more options, visit https://groups.google.com/groups/opt_out.
 
 

[Bob Carpenter]

I followed Ben's advice and wrote to the Eigen mailing
list. Here's what I got back, with attachment. I'm going
to do some more profiling in a second, but thought I should
send this now as is.

- Bob

==============from Eigen mailing list=======================
Hi,

Your version #3 is not correct, you have to use y.insert(y.begin(),
x.data(), x.data()+x.size()); With gcc I get the exact same same
result for each version. See attached benchmark file. I also observed
that in some cases clang is able to optimize away memcpy when the copy
is actually not necessary. This is the case in the attached file.

gael

On Thu, May 30, 2013 at 4:58 AM, Bob Carpenter <ca...@alias-i.com> wrote:
> We refactored our application to use Eigen column vectors
> instead of std::vector and observed a fairly significant
> slowdown. Profiling pointed to the copy constructor for
> Eigen vectors, so I'm wondering if there's something we're
> missing. Here's what I tried.
>
> To copy a column vector
>
> Matrix<double,Dynamic,1> x;
>
> there's the obvious copy constructor:
>
> 1. Matrix<double,Dynamic,1> y(x);
>
> or the direct sizing + memcpy approach:
>
> 2. Matrix<double,Dynamic,1> y(x.size());
> std::memcpy(&y(0),&x(0),x.size() * sizeof(double));
>
> By comparison, we can use an iterator based on the
> memory used for x's values to construct a standard
> vector:
>
> 3. std::vector<double> y(&x(0),&x(0) + x.size());
>
> I timed each version by repeatedly calling it, and found
> the times varied pretty widely:
>
> 1. 26s
> 2. 17s
> 3. 14s
>
> For Matrix<double,Dynamic,Dynamic> on the other hand, the
> copy constructor and the memcpy versions took about the same
> time.
>
> Is there a faster way to copy vectors? I'm guessing we'll have
> the same issue if I were to pass in references, resize them and
> then assign (like copy ctor) or memcpy.
>
> This was using clang++ -O3 on a Macbook Pro with i7 CPU and clang
> version:
>
> Apple clang version 4.1 (tags/Apple/clang-421.11.66) (based on LLVM
> 3.1svn)
> Target: x86_64-apple-darwin12.3.0
>
> - Bob
>
>

ctor.cpp

#include <Eigen/Dense>
#include <vector>
#include <iostream>
#include <bench/BenchTimer.h>

using namespace Eigen;

EIGEN_DONT_INLINE
double* cp1(const VectorXd& A, int i = 0)
{
VectorXd B(A);
return &B(i);
}

EIGEN_DONT_INLINE
double* cp2(const VectorXd& A, int i = 0)
{
VectorXd B(A.size());
std::memcpy(B.data(), A.data(), A.size()*sizeof(double));
return &B(i);
}

EIGEN_DONT_INLINE
double* cp3(const VectorXd& A, int i = 0)
{
std::vector<double> B;
B.insert(B.begin(), A.data(), A.data()+A.size());
return &B[i];
}

int main(int argc, char** argv)
{
using namespace std;

BenchTimer t1,t2,t3;
int tries = 10;
int rep = 1000;
int n = 100000;
VectorXd A(n);
A.setRandom();
BENCH(t1,tries,rep,cp1(A));
BENCH(t2,tries,rep,cp2(A));
BENCH(t3,tries,rep,cp3(A));
cout << t1.best() << " ; " << t2.best() << " ; " << t3.best() << endl;
}

[Michael Betancourt]

I think there's some confusion.  The problem is not the speed of copying an Eigen::VectorXd into another Eigen::VectorXd verses a std::vector, but rather an Eigen::VectorXd into an Eigen::VectorXd verses a std::vector into a std::vector.  The pure copy constructor of the Eigen::VectorXd is definitely ass-slow.

}

--
You received this message because you are subscribed to the Google Groups "stan development mailing list" group.

To unsubscribe from this group and stop receiving emails from it, send an email to stan-dev+unsubscribe@googlegroups.com.

[Michael Betancourt]

I updated the constructors and ran through the first volume of BUGS examples again only to see....no change.  The copy constructor does seem to be using less space, but now a huge amount of time is spent freeing objects.

v1.3.0 Time    New Time     Name

3 %                11.3 %         szone_free_definite_size

6 %                9.4 %           grad_log_prob

6 %                9.4 %           chain()

6 %                7.5 %           free

3 %                7.5 %           tiny_malloc_from_free_list

9%                 5.6%            check_not_nan

NA                 3.7%            ps_point(ps_point const&)

I do use slightly more temporary objects in the new code, but grad_log_prob and chain() are noticeably slower, too.  Has anyone benchmarked the autodiff code by itself?

Also, why is check_not_nan so slow?  In the sampler code I use the (var != var) trick which seems much faster.

[Marcus Brubaker]

There have been some changes in the autodiff code for multivariate model that may have made chain() a little slower, but it should balance out with reduced computation elsewhere.  Plus, most of the BUGS models don't use the multivariate distributions, so it shouldn't have such a consistent across the board hit.

I would benchmark grad_log_prob() for a single model to see if it's there.

And yeah, between calls to check_range, check_not_nan, etc, there is a big speedup to be had if we can ever reduce the need for repeatedly checking variables like this.

Cheers,

Marcus

To unsubscribe from this group and stop receiving emails from it, send an email to stan-dev+u...@googlegroups.com.

[Marcus Brubaker]

Also, is it possible to use an Eigen::Map instead of copying the vectors around?  You may need to rewrite some functions using the techniques from here: http://eigen.tuxfamily.org/dox/TopicFunctionTakingEigenTypes.html but it could avoid all these temporary allocations.

Cheers,

Marcus

To unsubscribe from this group and stop receiving emails from it, send an email to stan-dev+u...@googlegroups.com.

[Bob Carpenter]

On 5/31/13 9:19 AM, Michael Betancourt wrote:
> I think there's some confusion. The problem is not the speed of copying an Eigen::VectorXd into another Eigen::VectorXd
> verses a std::vector, but rather an Eigen::VectorXd into an Eigen::VectorXd verses a std::vector into a std::vector.
> The pure copy constructor of the Eigen::VectorXd is definitely ass-slow.

That's why I was suggesting allocating a VectorXd by size,
then copying into its memory. Now I'm not sure it's legal.

I'm about to go follow up all of the stuff from the Eigen
list yesterday.

- Bob

> On Thu, May 30, 2013 at 9:57 PM, Bob Carpenter <ca...@alias-i.com <mailto:ca...@alias-i.com>> wrote:
>
> I followed Ben's advice and wrote to the Eigen mailing
> list. Here's what I got back, with attachment. I'm going
> to do some more profiling in a second, but thought I should
> send this now as is.
>
> - Bob
>
> ==============from Eigen mailing list=======================
> Hi,
>
> Your version #3 is not correct, you have to use y.insert(y.begin(),
> x.data(), x.data()+x.size()); With gcc I get the exact same same
> result for each version. See attached benchmark file. I also observed
> that in some cases clang is able to optimize away memcpy when the copy
> is actually not necessary. This is the case in the attached file.
>
> gael
>

> On Thu, May 30, 2013 at 4:58 AM, Bob Carpenter <ca...@alias-i.com <mailto:ca...@alias-i.com>> wrote:
> > We refactored our application to use Eigen column vectors
> > instead of std::vector and observed a fairly significant
> > slowdown. Profiling pointed to the copy constructor for
> > Eigen vectors, so I'm wondering if there's something we're
> > missing. Here's what I tried.
> >
> > To copy a column vector
> >
> > Matrix<double,Dynamic,1> x;
> >
> > there's the obvious copy constructor:
> >
> > 1. Matrix<double,Dynamic,1> y(x);
> >
> > or the direct sizing + memcpy approach:
> >
> > 2. Matrix<double,Dynamic,1> y(x.size());

> > std::memcpy(&y(0),&x(0),x.__size() * sizeof(double));

> stan-dev+unsubscribe@__googlegroups.com <mailto:stan-dev%2Bunsu...@googlegroups.com>.
> For more options, visit https://groups.google.com/__groups/opt_out <https://groups.google.com/groups/opt_out>.

>
>
>
> --
> You received this message because you are subscribed to the Google Groups "stan development mailing list" group.

> To unsubscribe from this group and stop receiving emails from it, send an email to stan-dev+u...@googlegroups.com.

[Michael Betancourt]

The temporaries all need to be in memory at the same time -- in particular you need something like \sum_{i = 1}^{n} i = n (n  + 1) / 2 objects in memory when building a tree of depth n.

[Michael Betancourt]

That's why I was suggesting allocating a VectorXd by size,

then copying into its memory.  Now I'm not sure it's legal.

It looks like it works to me, and achieves the speed of the std::vector copy constructor.  I was just noting that Gael was benchmarking Eigen::Vector -> { Eigen::Vector, std::vector } which is something else entirely. 

[Bob Carpenter]

Great. That's all we care about. We don't need
this to work for arbitrary types.

The arithmetic types (double, int) and our scalar
autodiff types (var, fvar) are all safe to copy
bitwise.

- Bob

[Bob Carpenter]

On 5/31/13 12:41 PM, Michael Betancourt wrote:
> The temporaries all need to be in memory at the same time -- in particular you need something like \sum_{i = 1}^{n} i =
> n (n + 1) / 2 objects in memory when building a tree of depth n.

Could you remind me why we need more than log n?

- Bob

> On Fri, May 31, 2013 at 5:23 PM, Marcus Brubaker <marcus.a...@gmail.com <mailto:marcus.a...@gmail.com>> wrote:
>
> Also, is it possible to use an Eigen::Map instead of copying the vectors around? You may need to rewrite some
> functions using the techniques from here: http://eigen.tuxfamily.org/dox/TopicFunctionTakingEigenTypes.html but it
> could avoid all these temporary allocations.
>
> Cheers,
> Marcus
>
>
> On Fri, May 31, 2013 at 12:03 PM, Michael Betancourt <betan...@gmail.com <mailto:betan...@gmail.com>> wrote:
>
> I updated the constructors and ran through the first volume of BUGS examples again only to see....no change.
> The copy constructor does seem to be using less space, but now a huge amount of time is spent freeing objects.

I wonder if we could allocate a pool of these objects and
just re-use them? They're all double vectors of the same size.
Boost has a pool package:

http://www.boost.org/doc/libs/1_53_0/libs/pool/doc/html/index.html

I'd rather not introduce a static pool, but we could attach
it to one of the sampler objects.

> v1.3.0 Time New Time Name
> 3 % 11.3 % szone_free_definite_size
> 6 % 9.4 % grad_log_prob
> 6 % 9.4 % chain()
> 6 % 7.5 % free
> 3 % 7.5 % tiny_malloc_from_free_list
> 9% 5.6% check_not_nan
> NA 3.7% ps_point(ps_point const&)
>
> I do use slightly more temporary objects in the new code, but grad_log_prob and chain() are noticeably slower,
> too. Has anyone benchmarked the autodiff code by itself?

I haven't, but we should. The only change is that there's
a new virtual delete() function, but it shouldn't be called.
And there's a bit of cleanup for the new stack Marcus introduced,
but that should also be trivial.

The "tiny_malloc_from_free_list" also seems to indicate
that htere are more smaller objects getting allocated and
de-allocated. That's a problem, but it won't be coming
from auto-diff.

> Also, why is check_not_nan so slow? In the sampler code I use the (var != var) trick which seems much faster.

No idea --- it's a Boost implementation.

We need to fix it for the next release. I just created an issue.

- Bob

[Michael Betancourt]

You need to store rho and the state (at the very least momentum) at the boundary of each subtree -- maybe that scales as log(n)? I'm pretty sure the code is as tight as it can be without appealing to pools or the like.

[Bob Carpenter]

On 5/31/13 4:15 PM, Michael Betancourt wrote:
> You need to store rho and the state (at the very least momentum) at the boundary of each subtree -- maybe that scales as log(n)? I'm pretty sure the code is as tight as it can be without appealing to pools or the like.

Imisread n as the number of leapfrog steps,
not the depth.

Let's switch to the NUTS arXiv paper notation, where
j is tree depth and 2^j the number of leaf nodes (leapfrog
steps are in one-one correspondence with leaf nodes).

Section 3.1.2 of the paper provides a O(j) memory algorithm. It's a pretty
simple induction that assumes you only store the sample that
would be produced (if any) from each subtree and its (position,momentum)
at the endpoints. Because NUTS only checks that the ends of
each subtree don't make a U-turn, once you've completed a subtree
and chosen a sample from it, you never need to look inside it again.

Of course, it's still O(2^j) time because we do have to do each
leapfrog step.

If you're using O(j^2) for some reason, that could account
for the increased number of objects you're finding in the profiles.

- Bob

[Michael Betancourt]

It's a little more complicated than that: because you're building the tree recursively you have many incomplete trees in memory at any given time. It's easier to draw than to describe, but in any case it might end up being O(j) and I just made a bad order calculation. Again, I'm pretty sure there's no way to simplify the memory at all but I'll go over it again -- other eyes on the NUTS algorithm would be welcome.

[Bob Carpenter]

I'm looking at it NUTS now. I'm still trying to figure out
if we can get away from static casting ps points.

To see that memory is O(j) for a depth j tree, just use
a simple induction.

Base step: You can build a depth 0 tree (single node) in
constant time --- there's nothing to build.

Inductive step: Suppose you can build a depth n tree in
memory size O(n). To build a depth (n + 1) tree, you need to
first build a depth n tree for the first half, which will
take O(n) memory. Then you trim it down to its edges, so you're
now using a constant amount of memory and have maxed out at O(n).
Now build the depth n tree making up the second half of the
depth (n + 1) tree. That will take O(n) memory. Adding the
constant memory required to be stored for the first half, we
get O(n + 1) total.

- Bob

[Michael Betancourt]

> I'm looking at it NUTS now. I'm still trying to figure out
> if we can get away from static casting ps points.

Keep an eye on the fast_constructor branch, I've been cleaning up base_nuts.hpp here, too. I think the right approach to avoiding static casts is to introduce copy constructors and operator= for the base class and then just call them explicitly when manipulating a derived class.

> To see that memory is O(j) for a depth j tree, just use
> a simple induction.
>
> Base step: You can build a depth 0 tree (single node) in
> constant time --- there's nothing to build.
>
> Inductive step: Suppose you can build a depth n tree in
> memory size O(n). To build a depth (n + 1) tree, you need to
> first build a depth n tree for the first half, which will
> take O(n) memory. Then you trim it down to its edges, so you're
> now using a constant amount of memory and have maxed out at O(n).
> Now build the depth n tree making up the second half of the
> depth (n + 1) tree. That will take O(n) memory. Adding the
> constant memory required to be stored for the first half, we
> get O(n + 1) total.

Right. That's what I meant? :-p

[Bob Carpenter]

On 6/1/13 2:36 PM, Michael Betancourt wrote:
>> I'm looking at it NUTS now. I'm still trying to figure out
>> if we can get away from static casting ps points.
>
> Keep an eye on the fast_constructor branch, I've been cleaning up base_nuts.hpp here, too.

That's what I'm looking at.

> I think the right approach to avoiding static casts is to introduce copy constructors
> and operator= for the base class and then just call them explicitly when manipulating a derived class.

Why not just hold references to the base class. Then
all the virtual functions do what you want them to do rather
than being sliced as the static cast does. And it gets away
from all the ugly pointers in the code, too. For example, consider

ps_point* z = 0;
Eigen::VectorXd* rho = 0;

if (this->_rand_uniform() > 0.5)
{
z = &z_plus;
rho = &rho_plus;
util.sign = 1;
}
else
{
z = &z_minus;
rho = &rho_minus;
util.sign = -1;
}

this->_z.copy_base(*z);

int n_valid_subtree = build_tree(_depth, *rho, 0, z_propose, util);

*z = static_cast<ps_point>(this->_z);


You don't need the init to 0 here, as they'll get set whichever
way you go through the branch.

The above could be recoded as follows using references (and less
vertical space) as:

ps_point& z = &z_plus;
Eigen::VectorXd& rho = &rho_plus
if (this->_rand_uniform() > 0.5) {
util.sign = 1;
} else {
z = &z_minus;
rho = &rho_minus;
util.sign = -1;
}
this->_z = z;
int n_valid_subtree = build_tree(_depth, rho, 0, z_propose, util);
z = &(this->z_);

And there's still something fishy about having to use "this->" that
I don't understand.

If util.sign is used to multiply doubles, it should probably be
a double to avoid having to do the promotion on the fly.

- Bob

[Michael Betancourt]

> Why not just hold references to the base class. Then
> all the virtual functions do what you want them to do rather
> than being sliced as the static cast does. And it gets away
> from all the ugly pointers in the code, too. For example, consider
>
> ps_point* z = 0;
> Eigen::VectorXd* rho = 0;
>
> if (this->_rand_uniform() > 0.5)
> {
> z = &z_plus;
> rho = &rho_plus;
> util.sign = 1;
> }
> else
> {
> z = &z_minus;
> rho = &rho_minus;
> util.sign = -1;
> }
>
> this->_z.copy_base(*z);
>
> int n_valid_subtree = build_tree(_depth, *rho, 0, z_propose, util);
>
> *z = static_cast<ps_point>(this->_z);
>
> You don't need the init to 0 here, as they'll get set whichever
> way you go through the branch.

(Right, but ps_point is also used in many other files. Rather safe than sorry.)

> The above could be recoded as follows using references (and less
> vertical space) as:
>
> ps_point& z = &z_plus;
> Eigen::VectorXd& rho = &rho_plus
> if (this->_rand_uniform() > 0.5) {
> util.sign = 1;
> } else {
> z = &z_minus;
> rho = &rho_minus;
> util.sign = -1;
> }
> this->_z = z;
> int n_valid_subtree = build_tree(_depth, rho, 0, z_propose, util);
> z = &(this->z_);

Does the second call work? i.e. once z has been created as a reference, does z = &z_minus change the reference or call operator equals on z_plus? Setting z_plus = z_minus would obviously break the algorithm.

I tried to get everything to work with references but failed, hence the pointers.

> And there's still something fishy about having to use "this->" that
> I don't understand.

My best understanding is that there's some collision with the templates and the extension of the base class so you have to use this-> to reference any protected members of the base class.

> If util.sign is used to multiply doubles, it should probably be
> a double to avoid having to do the promotion on the fly.

Good call.

[Bob Carpenter]

On 6/1/13 3:59 PM, Michael Betancourt wrote:
>> Why not just hold references to the base class. Then
>> all the virtual functions do what you want them to do rather
>> than being sliced as the static cast does. And it gets away
>> from all the ugly pointers in the code, too. For example, consider
>>
>> ps_point* z = 0;
>> Eigen::VectorXd* rho = 0;
>>
>> if (this->_rand_uniform() > 0.5)
>> {
>> z = &z_plus;
>> rho = &rho_plus;
>> util.sign = 1;
>> }
>> else
>> {
>> z = &z_minus;
>> rho = &rho_minus;
>> util.sign = -1;
>> }
>>
>> this->_z.copy_base(*z);
>>
>> int n_valid_subtree = build_tree(_depth, *rho, 0, z_propose, util);
>>
>> *z = static_cast<ps_point>(this->_z);
>>
>> You don't need the init to 0 here, as they'll get set whichever
>> way you go through the branch.
>
> (Right, but ps_point is also used in many other files. Rather safe than sorry.)

The point here is that you're safe not setting it at all, because
it's set in one of the branches of the if/then. You have to be able to rely
on if/then doing the right thing!

My larger point was that it's NOT safe to just assume a "reasonable" default.
Quite the opposite --- it might cover up some potential problem.
It's much safer to init with NaN, so you know if you fail to override
the init you'll have a problem.

>> The above could be recoded as follows using references (and less
>> vertical space) as:
>>
>> ps_point& z = &z_plus;
>> Eigen::VectorXd& rho = &rho_plus
>> if (this->_rand_uniform() > 0.5) {
>> util.sign = 1;
>> } else {
>> z = &z_minus;
>> rho = &rho_minus;
>> util.sign = -1;
>> }
>> this->_z = z;
>> int n_valid_subtree = build_tree(_depth, rho, 0, z_propose, util);
>> z = &(this->z_);
>
> Does the second call work? i.e. once z has been created as a reference, does z = &z_minus change the reference

Yes, it changes the reference.

> or call operator equals on z_plus? Setting z_plus = z_minus would obviously break the algorithm.

Just setting z_plus = z_minus would employ the copy constructor
to set whatever z_plus refers to to hold copy of z_minus.

> I tried to get everything to work with references but failed, hence the pointers.
>
>> And there's still something fishy about having to use "this->" that
>> I don't understand.
>
> My best understanding is that there's some collision with the templates and the extension of the base class so you have to use this-> to reference any protected members of the base class.

I'll try to sort this out. All the this's are very
distracting to me.

- Bob

[Michael Betancourt]

> Yes, it changes the reference.

Sweet, that makes me very happy.

>> or call operator equals on z_plus? Setting z_plus = z_minus would obviously break the algorithm.
>
> Just setting z_plus = z_minus would employ the copy constructor
> to set whatever z_plus refers to to hold copy of z_minus.

Right -- equating the values kind of destroys the tree you were trying to build by collapsing the boundaries.

> I'll try to sort this out. All the this's are very
> distracting to me.

Agreed.