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 tutoriallike 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 higherorder 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 userfacing 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 passbyvalue. 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 nonconstant reference and are assigned by the deserialize
function.
If “(de)serialize” is too compsci 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 lessthanorequal 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

Generalize the third definition for
serialize
to allow arbitrary dimensional arrays. Hint: callserialize
recursively on the elements ofx
. 
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);
}