# Multivariate Function with Known Gradients - RStan

**URL:** <https://discourse.mc-stan.org/t/multivariate-function-with-known-gradients-rstan/2704>\
**Category:** Developers\
**Tags:** rstan\
**Created:** [December 6, 2017, 4:13pm UTC](https://discourse.mc-stan.org/t/multivariate-function-with-known-gradients-rstan/2704 "2017-12-06T16:13:11Z")\
**Posts on this page:** 4\
**Page:** 1

<div class="post-metadata">

**Author:** ![Jonas\_K](https://avatars.discourse-cdn.com/v4/letter/j/bbe5ce/32.png) [@Jonas\_K](https://discourse.mc-stan.org/u/Jonas_K)\
**Post date:** [December 6, 2017, 4:13pm UTC](https://discourse.mc-stan.org/t/multivariate-function-with-known-gradients-rstan/2704/1 "2017-12-06T16:13:11Z")

</div>

Hi,

As in the issue [Custom Function with Known Gradients](http://discourse.mc-stan.org/t/custom-function-with-known-gradients/2682), I got a function to work which maps `double -> double`. Eventually I want to use a function that maps `vector -> vector`, but for now I have difficulties to implement a function that maps `vector -> double` and is callable from a Stan file.

The function is declared as

```
namespace tests {

    double vecToDouble(const std::vector<double>& x)
    {
        double res = 0;
        for (int i=0; i<x.size(); ++i) res += x[i] * x[i];
        return res;
    }

    std::vector<double> grad_vecToDouble(const std::vector<double>& x)
    {
        std::vector<double> res;
        for (int i=0; i<x.size(); ++i) res.push_back(2 * x[i]);
        return res;
    }
}

```

and for the Stan interface I currently have

```
namespace stan {
    namespace math {

        double vecToDouble(const std::vector<double>& x) {
            return tests::vecToDouble(x);
        }

        var vecToDouble(const Eigen::Matrix<var, Eigen::Dynamic, 1>& x) {
            std::vector<var> b(x.data(), x.data() + x.rows() * x.cols());

            std::vector<double> a = value_of(b);
            double fa = tests::vecToDouble(a);
            std::vector<double> grad_fa = tests::grad_vecToDouble(a);

            return precomputed_gradients(fa, b, grad_fa);
        }

    }
}

```

It works in pure C++ as in

```
int main() {

    int N = 3;

    Eigen::Matrix<stan::math::var, Eigen::Dynamic, 1> x(N);
    std::vector<double> xx;
    for (int i=0; i<N; ++i)
    {
        stan::math::var xi = (double) i;
        x(i, 0) = xi;
        xx.push_back((double) i);
    }

    printf("function value: %.2f\n", stan::math::vecToDouble(xx));

    stan::math::var y = stan::math::vecToDouble(x);
    printf("vecToDouble(x) = %.2f\n", y.val());

    std::vector<double> g;
    std::vector<stan::math::var> in(x.data(), x.data() + x.rows() * x.cols());
    y.grad(in, g);
    for (int i=0; i<g.size(); ++i)
    {
        printf("g[%d] = %.2f\n", i, g[i]);
    }

    return 0;
}

```

However, once I want to call it from Stan it leads to this error

_…/Library/R/3.4/library/RcppEigen/include/Eigen/src/Core/functors/AssignmentFunctors.h:24:104: error: assigning to ‘double’ from incompatible type 'const stan::math::var’_  
_EIGEN\_DEVICE\_FUNC EIGEN\_STRONG\_INLINE void assignCoeff(DstScalar& a, const SrcScalar& b) const { a = b; }_

and a couple of errors follow after that which seem to be related to this error.

As a side note, I first declared the Stan interface as described in [https://github.com/stan-dev/math/wiki/Adding-a-new-function-with-known-gradients](https://github.com/stan-dev/math/wiki/Adding-a-new-function-with-known-gradients) like:

```
namespace stan {
    namespace math {

        double vecToDouble(const std::vector<double>& x);

        var vecToDouble(const std::vector<var>& x);

    }
}

```

But this lead to an error saying that the function `vecToDouble` wasn’t declared. I could see in the generated Cpp file that the variable passed into the function is declared as `Eigen::Matrix<stan::math::var, Eigen::Dynamic, 1>`, so I changed the types and the function can apparently be found.

Any help will be greatly appreciated!

Regards,  
Jonas

Operating System: OS X 10.12  
Interface Version: RStan 2.16.2  
Compiler/Toolkit: clang-900.0.38

---

<div class="post-metadata">

**Author:** ![Jonas\_K](https://avatars.discourse-cdn.com/v4/letter/j/bbe5ce/32.png) [@Jonas\_K](https://discourse.mc-stan.org/u/Jonas_K)\
**Post date:** [December 12, 2017, 12:14pm UTC](https://discourse.mc-stan.org/t/multivariate-function-with-known-gradients-rstan/2704/2 "2017-12-12T12:14:02Z")

</div>

The solution to this issue is quite simple as I got know in the thread [Barriers to New Developers](http://discourse.mc-stan.org/t/barriers-to-new-developers/163/4). Basically, one should not expose a vector-valued function based on `std::vector`, but always based on `Eigen::Matrix`. Additionally, there always has to be an overloaded function definition for both `double` and `var`.

```
namespace stan {
    namespace math {

        inline double
        vecToDoubleTwo(const Eigen::Matrix<double, Eigen::Dynamic, 1>& m) {

            double val = 0;
            for (int i=0; i<m.size(); ++i) val += m(i) * m(i);

            return val;
        }

        inline var
        vecToDoubleTwo(const Eigen::Matrix<var, Eigen::Dynamic, 1>& m) {
            using Eigen::Matrix;

            Matrix<double, Eigen::Dynamic, 1> m_d(m.rows());
            for (int i = 0; i < m.size(); ++i)
                m_d(i) = m(i).val();

            double val = stan::math::vecToDouble(m_d);

            vari** varis
            = ChainableStack::memalloc_.alloc_array<vari*>(m.size());
            for (int i = 0; i < m.size(); ++i)
                varis[i] = m(i).vi_;

            double* gradients
            = ChainableStack::memalloc_.alloc_array<double>(m.size());
            for (int i = 0; i < m.size(); ++i)
                gradients[i] = 9.9 * m_d(i);

            return var(new precomputed_gradients_vari(val, m.size(),
                                                      varis, gradients));
        }

    }
}
```

---

<div class="post-metadata">

**Author:** ![wds15](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/wds15/32/908_2.png) [@wds15](https://discourse.mc-stan.org/u/wds15)\
**Post date:** [December 12, 2017, 12:22pm UTC](https://discourse.mc-stan.org/t/multivariate-function-with-known-gradients-rstan/2704/3 "2017-12-12T12:22:09Z")

</div>

The alternative to providing overloads for var and double is to write the function templated in a way that it “does the right” thing automagically depending on the types passed into the function. For these types of functions, the stan-math library has invented the operands\_and\_partials thing. It is quite heavy in the use of template stuff, but I figured it out by looking at how it is used in the stan-math library. For example the definition of the normal lpdf density can serve as an example. The code gets a lot cleaner when using this.

---

<div class="post-metadata">

**Author:** ![Jonas\_K](https://avatars.discourse-cdn.com/v4/letter/j/bbe5ce/32.png) [@Jonas\_K](https://discourse.mc-stan.org/u/Jonas_K)\
**Post date:** [December 12, 2017, 12:40pm UTC](https://discourse.mc-stan.org/t/multivariate-function-with-known-gradients-rstan/2704/4 "2017-12-12T12:40:36Z")

</div>

I’ll have a look at that. Thanks!