Parameter packs for (de)serialization

Developers
1

I wanted to understand

and figured the best way to do that would be to take parameter packs out for a test spin.
Encouraged by the positive feedback from previous tutorial-like posts, I’ll going to drop a little suggestion for a Stan feature (not so much as a suggestion as a request for comments and that someone volunteer to do it—this is a fun feature!).

Serialization

The application I have in mind is converting a sequence of one or more Stan numerical data types (e.g., real, int, matrix, array of vectors, etc.) into a single array of real numbers and then converting back again to structured data without the user having to do any index fiddling. This should be particularly useful in dealing with arguments and system functions in our higher-order functionals like ODE solvers, where generic parameters are passed as an array.

It’s called “serialization” because structured objects are converted to a serial order.

Parameter packs

Parameter packs are the new C++ way of dealing with variadic arguments. In C parlance, “variadic” means accepting any number of arguments. We’ll be able to use that to write a function that accepts an arbitrary sequence of Stan types.

Serialization specification

We want a user-facing Stan function that serializes an arbitrary sequence of numerical data types x1, …, xN to a single array.

real[] serialize(T1 x1, ..., TN xN);

For example, it can be invoked with arguments serialize(real, real[], matrix) to produce the following array result.

real u = 17.2;
real v[2] = { 1.3, 9 };
matrix[2, 2] = [[1, 0], [0, 1]];

serialize(17.2, { 1.3, 9 }, [[1, 0], [0, 1]])
  == { 17.2, 1.3, 9, 1, 0, 0, 1 }

I’m then going to make an unusual proposal for deserialization. We actually pass by reference the things to be filled. I just want to call out that so far, I’ve managed to keep things simple so that as of version 2.18,

All arguments to Stan functions are passed by constant reference.

That makes it feel to users like everything is pass-by-value. The proposal I’m about to make will change that in this one particular case (unless we don’t treat it like a function; see below).

real a;
real b[2];
matrix[2, 2] c;

deserialize({ 17.2, 1.3, 9, 1, 0, 0, 1 }, a, b, c);
    a == 17.2
    b == { 1.3, 9 }
    c == [[1, 0], 
          [0, 1]]

Here, the arguments a, b, and c are passed by non-constant reference and are assigned by the deserialize function.

If “(de)serialize” is too comp-sci jargony, we could use “(un)pack”.

Special syntax?

We could use a special assignment statement to call this one out:

pack a, b, c into x;
unpack a, b, c from x;

Or we could use a specialized symbol

x <== a, b, c;

and one of

x => a, b, c;

or

a, b, c <= x;

We can get away with the overload of the less-than-or-equal operator <= because the type of x is an array. We probably want at least one element in the a, b, c sequence, though zero makes sense as a boundary condition.

Serialization and parameter packs

The key to understanding parameter packs is that they are written like ML or Prolog head/tail matching code. They consist of a base case for zero arguments and then are coded typically to peel off arguments one by one in a series of recursive cases.

The base case is easy:

template <typename T>
inline void serialize(std::vector<T>& x) { }

If there are no arguments, nothing happens. It works for any type T of container element, but we only care about int and real for Stan, as those are the argument types to our functionals.

Next, we can handle the scalar case,

template <typename T, typename... Txs>
inline void serialize(std::vector<T>& s,
                      double x, const Txs&... xs) {
  push_back(x);
  serialize(s, xs...);
}

The template pack parameter Txs and its variable xs require three syntactic conventions:

  • typename... in the template declaration
  • const Txs&... in the argument declaration
  • xs... as an argument

The key is that const Txs&... xs will absorb any number of distinct arguments, including zero. This function follows the usual functional programming pattern of tail recursion, which:

  • does something with the first (head) variadic argument
  • recursively calling the remaining (tail) variadic arguments for processing

With what we have so far, we’ll be able to serialize arbitrary sequences of scalars.

To keep things simple, here’s a function that’ll also allow 1D arrays to be part of the mix:

template <typename T, typename U, typename... Txs>
inline void serialize(std::vector<T>& s,
                      const std::vector<U>& x, const Txs&... xs) {
  s.insert(s.end(), x.begin(), x.end());
  serialize(s, xs...);
}

The usage within insert requires that values of type U are assignable to type T. This is only going to work for scalars. And then I could only get it to work with pure declarations to enable recursive calling; this has to go before the above so that everything’s declared before it’s ever instantiated.

template <typename T>
inline void serialize(std::vector<T>& x);

template <typename T, typename... Txs>
inline void serialize(std::vector<T>& s,
                      double x, const Txs&... xs);

template <typename T, typename U, typename... Txs>
inline void serialize(std::vector<T>& s,
                      const std::vector<U>& x, const Txs&... xs);

Deserialization

Easy peasy. We just reverse everything, including the const declarations.

template <typename T>
inline void deserialize(const std::vector<T>& s, std::size_t n) { }

template <typename T, typename... Txs>
inline void deserialize(const std::vector<T>& s, std::size_t n,
                        double& x, Txs&... xs) {
  x = s[n];
  deserialize(s, n + 1, xs...);
}

template <typename T, typename U, typename... Txs>
inline void deserialize(const std::vector<T>& s, std::size_t n,
                        std::vector<U>& x, Txs&... xs) {
  for (std::size_t i = 0; i < x.size(); ++i) x[i] = s[n++];
  deserialize(s, n, xs...);
}

These also need to be declared before being used recursively like this. The only tricky bit is that we need to provide an argument n to indicate the current position in the vector being deserialized.

Test drive

Here’s an example of how it runs:

int main() {
  std::vector<double> x = { 2, 3.5 };
  double y = 150;
  std::vector<double> z = { 1.9, 15.2, -1000000000000 };
  std::vector<double> s;
  serialize(s, x, y, z);
  print_vec("s", s);

  std::vector<double> x2(2);
  double y2;
  std::vector<double> z2(3);
  deserialize(s, 0u, x2, y2, z2);
  print_vec("x2", x2);
  std::cout << "y2 = " << y2 << std::endl;
  print_vec("z2", z2);
}

And here’s the output:

$ clang++ -std=c++1y serialize.cpp
$ ./a.out
s = (2, 3.5, 150, 1.9, 15.2, -1e+12)
x2 = (2, 3.5)
y2 = 150
z2 = (1.9, 15.2, -1e+12)

Exercises

  1. Generalize the third definition for serialize to allow arbitrary dimensional arrays. Hint: call serialize recursively on the elements of x.

  2. Generalize to matrix, vector, and row vector types.

The whole enchilada

Here’s everything if you want to run it:

#include <iostream>
#include <vector>


template <typename T>
inline void serialize(std::vector<T>& x);

template <typename T, typename... Txs>
inline void serialize(std::vector<T>& s,
                      double x, const Txs&... xs);

template <typename T, typename U, typename... Txs>
inline void serialize(std::vector<T>& s,
                      const std::vector<U>& x, const Txs&... xs);


template <typename T>
inline void serialize(std::vector<T>& x) { }

template <typename T, typename... Txs>
inline void serialize(std::vector<T>& s,
                      double x, const Txs&... xs) {
  s.push_back(x);
  serialize(s, xs...);
}

template <typename T, typename U, typename... Txs>
inline void serialize(std::vector<T>& s,
                      const std::vector<U>& x, const Txs&... xs) {
  s.insert(s.end(), x.begin(), x.end());
  serialize(s, xs...);
}



template <typename T>
inline void deserialize(const std::vector<T>& s, std::size_t n);

template <typename T, typename... Txs>
inline void deserialize(const std::vector<T>& s, std::size_t n,
                        double& x, Txs&... xs);

template <typename T, typename U, typename... Txs>
inline void deserialize(const std::vector<T>& s, std::size_t n,
                        std::vector<U>& x, Txs&... xs);

template <typename T>
inline void deserialize(const std::vector<T>& s, std::size_t n) { }

template <typename T, typename... Txs>
inline void deserialize(const std::vector<T>& s, std::size_t n,
                        double& x, Txs&... xs) {
  x = s[n];
  deserialize(s, n + 1, xs...);
}

template <typename T, typename U, typename... Txs>
inline void deserialize(const std::vector<T>& s, std::size_t n,
                        std::vector<U>& x, Txs&... xs) {
  for (std::size_t i = 0; i < x.size(); ++i) x[i] = s[n++];
  deserialize(s, n, xs...);
}

template <typename T>
void print_vec(const std::string& var, const std::vector<T>& y) {
  std::cout << var << " = (";
  for (int i = 0; i < y.size(); ++i) {
    if (i > 0) std::cout << ", ";
    std::cout << y[i];
  }
  std::cout << ")" << std::endl;
}

int main() {
  std::vector<double> x = { 2, 3.5 };
  double y = 150;
  std::vector<double> z = { 1.9, 15.2, -1000000000000 };
  std::vector<double> s;
  serialize(s, x, y, z);
  print_vec("s", s);

  std::vector<double> x2(2);
  double y2;
  std::vector<double> z2(3);
  deserialize(s, 0u, x2, y2, z2);
  print_vec("x2", x2);
  std::cout << "y2 = " << y2 << std::endl;
  print_vec("z2", z2);
}
5 Likes
2

This is cool, but won’t tuples replace this?

The line of thought that got me there was, “how is this going to check dimensions match?”, “maybe there should be something in the language”, “isn’t that tuples”?

3

I think this is how we could implement something like tuples in just the math library (well, the first bit before the language extensions), except that they also get serialized into a single byte array instead of a C++ type. This matters when sending data elsewhere like with MPI and GPU. In that form I don’t think it checks types at compile time on deserialization.

If we wanted to go the tuple route here, we’d first get those in, then add some serialization methods that do a similar thing using parameter packs for tuples specifically, then figure out how we want to hook up the tuples to the serialization methods (providing additional specializations on all relevant math library MPI/GPU/etc functions? or having the Stan compiler emit them when appropriate? or some other C++ magic?).

4

I would like to say that even for a complete coding disaster such as myself, parameter packs and variadic templates are quite straightforward to use.

There’s not much left to do beyond this - I think you’ve done the hard bits…

The one potential concern would be packs that mix int and double. If you cast from an int to a double to an int are you guaranteed to get the same int?

My case study now as always is sparse matrices. I’m about half way through a lme4 sparse matrix implementation and the ability to pack the sparsity structrue with the values would be convenient from an interface point of view. (Right now the interface is fine but a touch ugly)

2 Likes
5

Oh yeah, that is true. This is probly what the C++ implementation looks like regardless.

If this comes into the language before tuples, is there any reason to do anything other than placeholder functions?

I guess, which is easier to deprecate?

And another couple probably not so well thought out ideas,

  1. But what about functions that have access to the data, transformed data, parameters, and transformed parameters blocks? Maybe they’d be functor types or something? So that methods that accepted them could call things like get_flat_list_of_parameters(), get_flat_list_of_data() or whatever?

  2. What about parameterpacks in Stan?

But yeah those sound way harder than what you’ve already designed here in like 40 lines of code hahaha. But for the fuuuuuuture.

As a side note for #1, a problem came up recently where a user accidentally packed in parameters they weren’t using and it caused the algebra solver to break: Algebra solver has a side effect on log probability

Might be avoidable if these autodiff links were handled automagically?

6

Can you say a little more about this? Are you saying that you could end up putting in parameters that aren’t used, which then get weird autodiff stuff happening to them that shouldn’t, and there might be a way to notice that automatically?

Are you sort of suggesting that we should just allow essentially closures in the functions passed to map_rect and GPU stuff and handle it automatically? I think we could do that in the language but the implementation might be made easier by these serialization and deserialization functions.

7

Now that I read this again and think about it, I don’t think that guy should have had the problem he had. If the variable wasn’t used, it shouldn’t have been affected. I’m going to make an issue. This should be investigated more.

Yeah. That wouldn’t really solve everything, cause there’d still be stuff you’d want to pass to these functions manually, but I think it’d solve a lot of the annoying packing.

Yeah I think they’re the right things too – just trying to poke holes.

8

This really sounds like we should use straight std::tuple in the math library and maybe borrow some code from the (header -only) c++ cereal lib to handle the few things we want to serialize. I don’t think a home baked approach is the way to go here(?)

9

I like cereal’s packaging! We actually only need to serialize very few things, and we have specific desires about how we do it. I could see it going either way - importing and configuring cereal for our specific use case might be just as much work as the short examples Bob has above. If everything else is close enough to equal, I’d go for less dependencies. This came up with the spdlog part of the I/O proposal and that depencency weighs heavily on me, but it seems like that library offers quite a lot of what we need and has cheap enough config hooks… Could still be wrong there.

10

Right, which means documenting those desires before looking for a solution might be the place to start! For example, are we talking about serializing var-types or not? (can’t tell from the posts above).

1 Like
11

I agree this is important, even with header-only libs. Another option with something like cereal or arrow (which was meant for shared-memory computations) is to borrow the parts we need only. The advantage there is that when we deal with cerealization on the interface end (e.g.-writing binary output) we could use the same lib and kill two birds with one stone.

12

For what? We don’t have tuples in the Stan language yet.

We need the Stan types:

  • int, real, vector, row_vector, matrix (default Eigen column major)
  • arrays of same (all row major)

I also strongly prefer fewer dependencies.

Synchronized logging with config is a whole lot more complex than serialization of Stan types, so I think the dependency on spdlog makes sense there.

13

This makes no sense—the Stan language implementation doesn’t have to control how the math lib implements something. The Stan language is just generating some C++ code—generate code that uses std::tuple. When thinking about the stated use-case you don’t even need to serialize:

Why not just make the C++-level functional variadic to begin with so that you don’t have to generate packing code in the Stan language at all?

14

To me this looks great what you are proposing here. There is just one caveat in terms of efficiency: When you deserialize the data, then we are going to assign this to locally delcared variables. Thus, the serialized data inside an ODE functor will be cast into a var which is rather costly to deal with. I know that this is an orthogonal problem… but can we do something about it?

Other than that I like this suggestion a lot, because:

  • it will fit into our existing code base nicely
  • tuples are probably superior, but those take time until they are there and even if they have landed we would still need to augment all our ODE/algebra_solver/… stuff
  • as I understand the only limitation is that double and int things need to be separate, but I think that’s fine as those are anyway always separate arguments to our functions.

Going one step further… how do you suggest to deal with ragged arrays and serialization? So what would be nice to also cover is map_rect style arrays. For that you need to pack J units and for each you create a 1D array. Would you recommend that the user then does J packings himself in a loop like

real x_r[J,big number];
for(j in 1:J) {
    x_r[j] = serialize(...stuff...);
}

How much effort do you think this is to do… and does someone has resources for it?

15

… also: Would this land in stan-math or in stan? It got nothing to do with autodiff and is a language convenience feature.

16

It’s essentially a C++ convenience tool for the MPI and ODE style functions in Math. Could be exposed to Math API users and Stan language users.

17

Now that you mention it, that would be way better than serializing.

The main thing we need to do is separate the parameters from the data, which I think we can do with traits. Then we need to serialize the parameters into a vector to use with the underlying solvers.

Then the calls look like this:

vector[N, K] y
 = ode_integrate(foo_sys, vector y0, real t0, vector ts, T1 x1, ..., TN xN);

where the Tn are any Stan type. Then the function foo_sys would have to look like:

vector foo_sys(vector y, real t, T1 x1, ..., TN xN);

Let’s do this.

This is an orthogonal optimization. We should be able to tighten up our variable declarations. One approach is to just use auto pretty much everywhere. I’m pinging @mitzimorris as she’s the one most familiar with the new generator code.

1 Like
18

Yes yes yes!

19

Off topic, but while we’re at it, what are our chances of optional/named parameters?

1 Like
20

At some point you write something very similar to the variadic serialization code you have above, no? But agreed - don’t make the user do it! Good point @sakrejda.

1 Like