> I was able to get the idea of the problem here, although not the solution. It's a little harder to test and understand
> fwd diff as it's not yet used in Stan.

We actually want the forward and reverse-mode auto-diff to
be usable outside Stan.  So they all have to be tested
independently of their use in Stan.

Peter has reasonable tests in place for everything as of
now, so that can be an example of how to test.  The easiest
way to understand how it all works is to look at those tests.
But it seems like you have the basic idea given the implementation.


In generla, every fvar op is a pair consisting of the value and
the derivative with respect to some specified variable (or
more than one in the case of directional derivatives).
Derivatives propagate like this through any operation (just do
the obvious generalization or restriction for 3 args or 1 arg):

    fvar<T> OP(fvar<T> x, fvar<T> y) {
       return fvar<T>(OP(x.val_, y.val_),
                      d_ * d/dx OP(x,y) + y.d_ * d/dy OP(x,y) );
    }

> Basically you something like this: template <typename T> inline fvar<typename
> stan::return_type<T>::type> log_diff_exp(const fvar<T>& x1, const fvar<T>& x2) { using stan::math::log_diff_exp; using
> stan::math::NOT_A_NUMBER; using std::exp; if (x1.val_ <= x2.val_) return fvar<T>(NOT_A_NUMBER, NOT_A_NUMBER); return
> fvar<T>(log_diff_exp(x1.val_, x2.val_), x1.d_ / (1 - exp(x2.val_ - x1.val_) ) + x2.d_ / (1 - exp(x1.val_ - x2.val_) ) );
> } ??? (not tested)

The text got smushed here, but yes, that's basically right assuming those
are the right partials.

> Great. Could you coordinate with Peter Li before jumping in on this? I just suggested some fixes for Peter Li (on the
> dev list) for the matrix ops, which really need to be improved for efficiency and numerical stability.

> I was going to open a discussion at stan-dev before proceeding, I just already done these changes so I wouldn't forget
> them later.

Thanks.

- Bob