I've adapted the KyberSlash countermeasure to also work for the compression instances with d = 10 and d = 11. Not sure it's necessary but doesn't hurt either. I'd definitely be happy about some other sharp pair of eyes.

Polynomials store their coefficients as a std::vector<>. Some of those store secret values and should therefore use secure_vector. Also, it might be beneficial to add first-level support for Polynomials (and friends) as Strong types. (Most likely this would be another pull request, though.)

I feel like this is a top-level concept we might want to establish deeper in the library. Perhaps as an (internal) replacement for the DataSource class.

Namely: A thing that can provide a concrete number of bytes where I don't really know how many bytes will be needed. The source may be infinite or not.

Things that could fall in this abstraction may be:

• the RandomNumberGenerator
• a XOF
• read from some buffer (think BufferSlicer)
• from some I/O source (e.g. the network)

In Kyber/Dithithium this is used to either sample some random polynomial from a XOF or to unmarshal it from a bitpacked encoding.

Other use cases may be:

• parsing in general (TLS, X.509 and friends)
• data loading, e.g. like:
  load_be<uint64_t>(rng) to get a random 64 bit integer. (See also Internal: load_be/le should accept a BufferSlicer #4068)
• among other things

From a technical perspective, one might even implement it like shown here: As a C++20 concept that tries to build a thin wrapper around "the source". That way, new data sources could just hook themselves into this "helper" and decide for themselves to be a "data source" or not.

Cool idea 👍. I'm not sure if we want to commit to a pure concept construction, though. In some cases, it may be more convenient to have a class like DataSource. We could still extend the DataSource class to fulfill the respective concept if one wants to juggle with templates.

This was different in the reference implementation. @FAlbertDev we did look into this a while ago, I feel we should take this up again and get to the bottom of it.

I've just investigated the usage of Montgomery in Kyber. First, a short recap on Montgomery:
The main observation is that every Montgomery multiplication $c=a \cdot b$ adds a factor of $R^{-1}$ to the result c where $R$ is the Montgomery factor (i.e., we compute $a \cdot b \cdot R^{-1}$) At some point we need to multiply $R$ to neutralize this $R^{-1}$ factor. Therefore, we could multiply $R$ to either $a$ or $b$. However, we can also multiply $R$ to the result $c$.
In our Kyber and Dilithium code, the following operations add or remove the Montgomery factor.

• the function montgomery adds a factor $R$
• the function inv_ntt adds a factor $R$. Note that ntt doesn't add a Montgomery factor.
• As already said, a Montgomery multiplication adds a factor $R^{-1}$ to each poly, vector, or matrix element. This holds for all multiplication types, i.e., polynomial-polynomial, vector-vector, matrix-vector multiplication, etc.

For encryption and decryption, this setup works well since there is always one multiplication (adding factor $R^{-1}$) followed by inv_ntt (adding factor $R$). This does not work here since we store the value t in NTT domain and, therefore, do not call inv_ntt. We have a factor of $R^{-1}$ to eliminate. Currently, we do that by adding the factor $R$ to the matrix beforehand using montgomery. Instead, I propose to add this factor to t after the multiplication. The fewer elements multiplied by $R$, the faster it should be.

   auto A = kyber_sample_matrix(rho, false /* not transposed */, mode);
   auto s = ntt(ps.sample_polynomial_vector_cbd_eta1());
   const auto e = ntt(ps.sample_polynomial_vector_cbd_eta1());

   auto t = montgomery(A * s);
   t += e;
   t.reduce();

   // End Algorithm 12 ---------------------------

Nice! Thanks for double-checking and the detailed write-up! It would probably be really nice to track the montgomery parameter in the type system as well, but it does feel like overkill.

Also: by applying montgomery() to the t vector, we can remove the overload for PolynomialMatrix altogether. :)

As mentioned before: it would be good if these could be the basis for strong types to disambiguate them in interfaces E.g. to replace the multiple poly_pack_t0, ..._t1, ... methods in dilithium_algos.cpp with an overloaded poly_pack() that just "does the right thing" ™ based on the strong type.

Dilithium has many integer constants that change depending on the instantiation of the algorithm. I found it worthwhile to restrict them with explicit enums. I'm okay with this concept, also given the restricted (algorithm-internal) scope. But I feel there must be a better way for this. :)

I also like this idea. While the implied switch-cases are annoying, they guarantee we do not forget a value.

Those became ubiquitous and at this point I do start wondering whether we should just put that in the destructors of Stuffer and Slicer (perhaps explicitly disengagable). What do y'all think?

Sounds good to me.

I'll do that after 3.5.0 is cut, to avoid unwanted surprises.

I've just investigated the usage of Montgomery in Kyber. First, a short recap on Montgomery:
The main observation is that every Montgomery multiplication $c=a \cdot b$ adds a factor of $R^{-1}$ to the result c where $R$ is the Montgomery factor (i.e., we compute $a \cdot b \cdot R^{-1}$) At some point we need to multiply $R$ to neutralize this $R^{-1}$ factor. Therefore, we could multiply $R$ to either $a$ or $b$. However, we can also multiply $R$ to the result $c$.
In our Kyber and Dilithium code, the following operations add or remove the Montgomery factor.

• the function montgomery adds a factor $R$
• the function inv_ntt adds a factor $R$. Note that ntt doesn't add a Montgomery factor.
• As already said, a Montgomery multiplication adds a factor $R^{-1}$ to each poly, vector, or matrix element. This holds for all multiplication types, i.e., polynomial-polynomial, vector-vector, matrix-vector multiplication, etc.

For encryption and decryption, this setup works well since there is always one multiplication (adding factor $R^{-1}$) followed by inv_ntt (adding factor $R$). This does not work here since we store the value t in NTT domain and, therefore, do not call inv_ntt. We have a factor of $R^{-1}$ to eliminate. Currently, we do that by adding the factor $R$ to the matrix beforehand using montgomery. Instead, I propose to add this factor to t after the multiplication. The fewer elements multiplied by $R$, the faster it should be.

   auto A = kyber_sample_matrix(rho, false /* not transposed */, mode);
   auto s = ntt(ps.sample_polynomial_vector_cbd_eta1());
   const auto e = ntt(ps.sample_polynomial_vector_cbd_eta1());

   auto t = montgomery(A * s);
   t += e;
   t.reduce();

   // End Algorithm 12 ---------------------------

Could the gist of this be useful beyond the CRYSTALS scope? E.g. FrodoKEM has a pack/unpack routine that -- on first glance at least -- looks like it could use it as well. @randombit perhaps there are more examples?