\documentclass[acmtog]{acmart}
\acmSubmissionID{126}
\usepackage{booktabs} % For formal tables
\usepackage{amsmath}
\usepackage{nicematrix}
\usepackage{tikz}
\usetikzlibrary {arrows.meta}
\usepackage[normalem]{ulem}
\usepackage{wrapfig}
\usepackage{algpseudocodex}
\usepackage{multirow}
\usetikzlibrary{decorations.pathreplacing}
\usetikzlibrary{spy} % LATEX and plain TEX
% TOG prefers author-name bib system with square brackets
\citestyle{acmauthoryear}
%\setcitestyle{nosort,square} % nosort to allow for manual chronological ordering
\usepackage{overpic}
\usepackage[ruled]{algorithm2e} % For algorithms
\renewcommand{\algorithmcfname}{ALGORITHM}
\SetAlFnt{\small}
\SetAlCapFnt{\small}
\SetAlCapNameFnt{\small}
\SetAlCapHSkip{0pt}
\DeclareMathSymbol{:}{\mathord}{operators}{"3A}
\DeclareMathOperator{\diag}{diag}
\DeclareMathOperator{\Id}{Id}
% Metadata Information
\setcopyright{cc}
\setcctype{by}
\acmJournal{TOG}
\acmYear{2025} \acmVolume{44} \acmNumber{4} \acmArticle{} \acmMonth{8} \acmPrice{}\acmDOI{10.1145/3730821}
\usepackage[many]{tcolorbox}
\tcbuselibrary{hooks, breakable}
\newtcolorbox{prop}[2][]{boxsep=0mm,bottom=1pt,top=0pt,toptitle=0mm, bottomtitle=0mm,detach title,after upper={\par\texttt{\tcbtitle}\hfill},arc=0mm,outer arc=0mm, boxrule=1pt, toprule=1pt,boxsep=0pt,colbacktitle=white,colback=white,coltitle=black,
title={#2},#1}
%% https://tex.stackexchange.com/questions/568782/allowing-paragraph-continuation-after-a-tcolorbox
\makeatletter
\tcbset{
after app={%
\ifx\tcb@drawcolorbox\tcb@drawcolorbox@breakable
\else
% add only when not breakabel
\@endparenv
\fi
}
}
% for breakable
\appto\tcb@use@after@lastbox{\@endparenv\@doendpe}
\makeatother
% Document starts
\begin{document}
% Title portion
\title{Sobol' Sequences with Guaranteed-Quality 2D Projections}
% DO NOT ENTER AUTHOR INFORMATION FOR ANONYMOUS TECHNICAL PAPER SUBMISSIONS TO SIGGRAPH 2019!
\author{Nicolas Bonneel}
\orcid{0000-0001-5243-4810}
\affiliation{
\institution{CNRS, Université Claude Bernard Lyon 1, INSA Lyon}
\country{France}
}
\email{nicolas.bonneel@liris.cnrs.fr}
\author{David Coeurjolly}
\orcid{0000-0003-3164-8697}
\affiliation{
\institution{CNRS, Université Claude Bernard Lyon 1, INSA Lyon}
\country{France}
}
\email{david.coeurjolly@cnrs.fr}
\author{Jean-Claude Iehl}
\orcid{0000-0001-6877-2398}
\affiliation{
\institution{Université Claude Bernard Lyon 1, CNRS, INSA Lyon}
\country{France}
}
\email{jean-claude.iehl@liris.cnrs.fr}
\author{Victor Ostromoukhov}
\orcid{0009-0004-3123-9388}
\affiliation{
\institution{Université Claude Bernard Lyon 1, CNRS, INSA Lyon}
\country{France}
}
\email{victor.ostromoukhov@liris.cnrs.fr}
%\renewcommand\shortauthors{Bonnneel, N. et al}
\setlength{\fboxsep}{.2pt}
\newcommand{\be}{{{\small [e]}}}
\newcommand{\btwoe}{{{\small [2e]}}}
\newcommand{\bfoure}{{{\small [4e]}}}
\newcommand{\bet}{{{\small [\tilde{e}]}}}
\newcommand{\btwoet}{{{\small [2\tilde{e}]}}}
\newcommand{\bfouret}{{{\small [4\tilde{e}]}}}
\definecolor{OliveGreen}{RGB}{63, 126, 49}
\definecolor{BurntOrange}{RGB}{255, 142, 39}
\newcommand{\dummyfig}[3]{
\begin{tikzpicture}
\draw[color=blue,fill=blue!20] (0,0) rectangle (#2,#3) node[pos=.5] {#1};
\end{tikzpicture}}
\begin{teaserfigure}
\centering
\includegraphics[width=12cm]{figures_siggraph25/teaser_e5_10.pdf}
\begin{overpic}[width=17cm]{figures_siggraph25/JKvsOneTwoSeq.png}
\put(-2,8.5) {\rotatebox{90}{Sobol'}}
\put(-2,2) {\rotatebox{90}{Ours}}
\put(-2.2,-3){\rotatebox{90}{dims.}}
\put(1.5,-2){$(0,1)$}
\put(8.61538,-2){$(2,3)$}
\put(15.7308,-2){$(4,5)$}
\put(22.8462,-2){$(6,7)$}
\put(29.9615,-2){$(8,9)$}
\put(37.0769,-2){$(10,11)$}
\put(44.1923,-2){$(12,13)$}
\put(51.3077,-2){$(14,15)$}
\put(58.4231,-2){$(16,17)$}
\put(65.5385,-2){$(18,19)$}
\put(72.6538,-2){$(20,21)$}
\put(79.7692,-2){$(22,23)$}
\put(86.8846,-2){$(24,25)$}
\put(94,-2){$(26,27)$}
\end{overpic}
\vspace{0.5cm}
\caption{We demonstrate that 2D Sobol' sequences constructed with polynomials $p$ and $p^2+p+1$ have a characteristic matrix $K = M_{p^2+p+1} M_{p}^{-1}$ that can be obtained with a simple recursive algorithm.
This is illustrated using polynomials of degrees $e=5$ and $2e=10$, where each colored block has dimensions $e \times e$.
They produce high-quality $(1, 2)$-sequences (with quality factor $t=1$) under mild conditions on $K$'s blocks and $p$ (bottom). The quality factor $t$ of each point set is indicated in the upper-right corner of each patch (only 256 points are shown).
We use these $(1, 2)$-sequences to construct higher-dimensional low-discrepancy sequences with high-quality 2D and 4D projections. }
\label{fig:teaser}
\end{teaserfigure}
\begin{abstract}
Low-discrepancy sequences, and more particularly Sobol' sequences are gold standard for drawing highly uniform samples for quasi-Monte Carlo applications.
They produce so-called $(t,s)$-sequences, that is, sequences of $s$-dimensional samples whose uniformity is controlled by a non-negative integer quality factor $t$.
The Monte Carlo integral estimator has a convergence rate that improves as $t$ decreases. Sobol' construction in base 2 also provides extremely fast sampling point generation using efficient xor-based arithmetic.
Computer graphics applications, such as rendering, often require high uniformity in consecutive 2D projections and in higher-dimensional projections at the same time.
However, it can be shown that, in the classical Sobol' construction, only a single 2D sequence of points (up to scrambling), constructed using irreducible polynomials $x$ and $x+1$, achieves the ideal $t=0$ property.
Reusing this sequence in projections necessarily loses high dimensional uniformity.
We prove the existence and construct many 2D Sobol' sequences having $t=1$ using irreducible polynomials $p$ and $p^2+p+1$.
They can be readily combined to produce higher-dimensional low discrepancy sequences with a high-quality $t=1$, guaranteed in consecutive pairs of dimensions.
We provide the initialization table that can be directly used with any existing Sobol' implementation, along with the corresponding generator matrices, for an optimized 692-dimensional Sobol' construction.
In addition to guaranteeing the $(1,2)$-sequence property for all consecutive pairs, we ensure that $t \leq 4$ for consecutive 4D projections up to $2^{15}$ points.
\end{abstract}
%
% The code below should be generated by the tool at
% http://dl.acm.org/ccs.cfm
% Please copy and paste the code instead of the example below.
%
\begin{CCSXML}
10002950.10003648.10003670.10003687
Mathematics of computing~Random number generation
500
10010147.10010371.10010372
Computing methodologies~Rendering
300
10002950.10003714.10003715.10003718
Mathematics of computing~Computations in finite fields
100
10002950.10003714.10003740
Mathematics of computing~Quadrature
300
\end{CCSXML}
\ccsdesc[500]{Mathematics of computing~Random number generation}
\ccsdesc[300]{Computing methodologies~Rendering}
\ccsdesc[100]{Mathematics of computing~Computations in finite fields}
\ccsdesc[300]{Mathematics of computing~Quadrature}
%
% End generated code
%
\keywords{Quasi-Monte Carlo, Sobol', low discrepancy sequences, irreducible polynomials}
\maketitle
\section{Introduction}
Integrating a function using quasi-Monte Carlo consists in evaluating it at well-chosen uniformly distributed sample points, and averaging these values.
Sobol' sequences arise as the cornerstone of quasi-Monte Carlo, by producing extremely well-distributed sampling points whose uniformity drastically improves the convergence rate, compared to classical Monte Carlo methods.
These points are constructed using a fast and compact recursive algorithm involving polynomials and matrices.
Sobol' sequences have thus been widely adopted in computer graphics, notably for rendering where efficiency is paramount. Their mathematical beauty, their connections to Pascal matrices (or Sierpinsky triangle) and Galois theory also make them appealing to the mathematician, but their difficulty comprehending may discourage others. This paper gives new fundamental mathematical insights on Sobol' sequences, exploring particular pairs of polynomials and new Sobol' matrix constructions, with practical and provable benefits in terms of quasi-Monte Carlo integration error.
Sobol' sequences produce a sequence of samples in arbitrary dimension by multiplying Sobol' matrices with a base-$b$ representation of the sample index, where one Sobol' matrix is used per dimension. To compute a Sobol' matrix, a polynomial $p$ of degree $e$ and an initialization matrix of size $e \times e$ (often called \textit{direction vectors} in prior work) are required. The (triangular) Sobol' matrix is recursively constructed column by column, where the next column is computed as a linear combination of the previous $e$ columns
weighted by the polynomial coefficients (plus a shifted column), with all operations performed modulo $b$ (or over Galois Field $GF(b)$).
The uniformity of the produced sample points is determined by a non-negative integer parameter $t$, where $t=0$ corresponds to the best achievable quality, and the quasi-Monte Carlo integral estimator converges with a rate roughly proportional to $b^t$~(see \cite{Lemieux2009MCQMCSampling} page 157, and \cite{niederreiter1988low}).
An $s$-dimensional set with $b^m$ samples with quality $t$ is called a $(t,m,s)$-net (see Fig.~\ref{fig:t1}). If such a point set is a $(t,m,s)$-net for all $m$, then it is a $(t,s)$-sequence.
We are particularly interested in the case $s=2$, not only for producing 2D points for 2-dimensional problems, but most importantly to control 2-dimensional projections of higher-dimensional problems~\cite{joe2008constructing} arising for example in computer graphics such as rendering.
It has been demonstrated that, in base $b=2$, and when considering 2D points, matrices that generate $(0,2)$-sequences are inherently related by a Pascal matrix~\cite{HoferSuzuki2019,ahmed2023analysis}.
Consequently, the only pair of polynomials that produce $t=0$ using Sobol' construction are $x$ and $ x+1$.
The space of $t=0$ sequences is thus extremely limited, when $b=2$.
One possible solution is to increase $b$, as suggested in~\cite{Quad2024}, which allows for additional polynomials that generate $t=0$ sequences.
However, the effect on the integration error for $t\neq 0$ becomes more significant due to the $b^t$ factor in the convergence rate.
In addition, base $b=2$ allows for extremely fast implementations using vectorized xor-based arithmetic, which is not the case when $b>2$. For practical reasons, Sobol' polynomials and initialization matrices have thus been optimized mostly in base $b=2$~\cite{joe2008constructing}, though consecutive dimensions typically produce increasing $t$ values as dimension increases which limits their use for rendering~\cite{christensen2018renderman}.
Focusing on the case $b=2$, in this paper we show that there exist many $(1,2)$-sequences that can be constructed from pairs of irreducible polynomials $p$ and $p^2+p+1$.
These 2D sequences can be combined to produce higher dimensional $(t,s)$-sequences of high-quality $t=1$ in consecutive 2D projections, which is the best quality achievable for $b=2$ aside from the single $t=0$ pair mentioned above.
We also observed that, in practical integration problems, the quasi-Monte Carlo convergence rate did not differ significantly between $t=0$ and $t=1$ in base $2$ (see Fig.~\ref{fig:correlation}), making $t=1$ a compelling compromise.
Our use of a standard $b=2$ Sobol' framework makes our sequences readily usable in production renderers already using Sobol' sequences, by simply replacing existing polynomials and initializations with ours.
In our quest to prove the quality of our polynomials, we discovered a new recursive construction of Sobol' matrices, derived by iteratively squaring polynomials.
This, in turn, led to the identification of interesting patterns that characterize these sequences.
This paper explores the depths of this new construction and demonstrates how our sequences can be applied to quasi-Monte Carlo rendering. Code and data are available at \href{https://github.com/liris-origami/OneTwoSobolSequences}{https://github.com/liris-origami/OneTwoSobolSequences}.
\begin{figure}[!h]
\begin{overpic}[width=1.5cm]{figures_siggraph25/tms/t0/0_0_4.pdf}
\put(-25,0){\rotatebox{90}{$(0,4,2)-$net}}
\put(-25,-130){\rotatebox{90}{$(1,4,2)-$net}}
\put(-25,-255){\rotatebox{90}{$(2,4,2)-$net}}
\end{overpic}
\includegraphics[width=1.5cm]{figures_siggraph25/tms/t0/0_1_3.pdf}
\includegraphics[width=1.5cm]{figures_siggraph25/tms/t0/0_2_2.pdf}
\includegraphics[width=1.5cm]{figures_siggraph25/tms/t0/0_3_1.pdf}
\includegraphics[width=1.5cm]{figures_siggraph25/tms/t0/0_4_0.pdf}\\
\vspace{.2cm}
\hrule
\vspace{.2cm}
\includegraphics[width=1.5cm]{figures_siggraph25/tms/t1/1_0_3.pdf}
\includegraphics[width=1.5cm]{figures_siggraph25/tms/t1/1_1_2.pdf}
\includegraphics[width=1.5cm]{figures_siggraph25/tms/t1/1_2_1.pdf}
\includegraphics[width=1.5cm]{figures_siggraph25/tms/t1/1_3_0.pdf}
\vspace{.2cm}
\hrule
\vspace{.2cm}
\includegraphics[width=1.5cm]{figures_siggraph25/tms/t2/2_0_2.pdf}
\includegraphics[width=1.5cm]{figures_siggraph25/tms/t2/2_1_1.pdf}
\includegraphics[width=1.5cm]{figures_siggraph25/tms/t2/2_0_2.pdf}
\caption{The definition of a $(t,m,2)$-net of $n=2^m$ samples is that all dyadic intervals of area $2^t/2^m$ contain $2^t$ points.
As such, a $(0,4,2)$-net (top) has exactly one (=$2^0$) point in each interval of volume $1/{2^4}$, a $(1,4,2)$-net has exactly two points (=$2^1$) in each interval of volume $1/{2^3}$ (middle), and a $(2,4,2)$-net has four points on intervals of volume $1/{2^2}$ (bottom).
Increasing $t$ thus reduces uniformity constraints and produces larger gaps and clusters in the distribution.
\label{fig:t1}}
\end{figure}
\begin{figure*}[!ht]
\includegraphics[width=\textwidth]{scripts/t-correlation.pdf}
\caption{Experimental validation in 2D of the impact of the $t$ value of a Sobol' sequence on various metrics (from left to right, the generalized $L_2$ discrepancy \cite{hickernell1998generalized} and Monte Carlo integration errors on random Gaussians and Heaviside functions). We randomly select a pair of Sobol' polynomials from the first thousand entries of Joe and Kuo \shortcite{joe2008constructing}, we evaluate the metrics for each sample count and plot error distribution (box plots) per values $t$ (box plot colors) observed for that sample count (64 realizations for each $t$). We observe a strong correlation between the observed $t$ value of the point set and its quality for Monte Carlo integration, while $t=1$ appears as a good compromise in quality compared to $t=0$. \label{fig:correlation}}
\end{figure*}
\section{Related work}
We summarize fundamentals about Sobol' construction~\cite{Sobol1967distribution} and refer the reader to reference books~\cite{Niederreiter1992random,Lemieux2009MCQMCSampling,Dick_Pillichshammer_2010} for in-depth discussions.
\paragraph{Sobol' construction}
To construct $(t,s)$-sequences, Sobol' proposed a solution based on primitive polynomials. Given $s$ primitive polynomials $p_0, \dots, p_{s-1}$ in the Galois Field of prime base $b$ called $GF(b)$ (think of it as the set of integers modulo $b$), Sobol' constructs $s$ upper triangular matrices $M_{p_0}, \dots, M_{p_{s-1}}$, one per dimension, which are used to obtain each coordinate of the $i^{th}$ sample point $\mathbf{x}_i$ in the sequence. Specifically, the $k^{th}$ coordinate of the $i^{th}$ sample point, $\mathbf{x}_i^k$, is obtained using a matrix-vector product between matrix $M_{p_k}$ and the base $b$ decomposition of the sample point index $i$ interpreted as a column vector, denoted $\bar{i}$ hereafter, $\bar{q} = M_{p_k} \bar{i}$ with $i = \sum_j \bar{i}_{j} b^j$. The sample point coordinate is now
$$\mathbf{x}_i^k = \sum_j \bar{q}_j b^{-j}\,.$$
The construction of the upper triangular invertible matrix $M_{p_k}$ using the primitive polynomial $p_k$ uses a recursive formula to obtain a column given its $e$ previous columns, where $e$ is the degree of the polynomial $p_k$.
The initialization of this recursion, an $e \times e$ upper triangular matrix at the top-left corner of matrix $M_{p_k}$, provides additional degrees of freedom in addition to the chosen polynomial.
We base our construction on matrix blocks instead of the more common column-wise recursion, as proposed by Faure and Lemieux~\shortcite{faure2016irreducible}, that we briefly describe in Sec.~\ref{sec:blocklemieux}.
Joe and Kuo numerically optimized top-left $e\times e$ blocks, resulting in improved Sobol' sequences on consecutive projections~\cite{joe2008constructing}.
Matrices can be directly obtained without Sobol' recurrence using an integer linear program solver, but this limits their use to only moderately large problem~\cite{matbuilder22}.
Faure and Lemieux showed that the larger set of irreducible polynomials can be used instead of primitive polynomials~\cite{Sloane2001,faure2016irreducible}.
Irreducible polynomials are similar to prime numbers, meaning they cannot be factored into products of other non-constant polynomials.
Faure and Lemieux showed that the parameter $t$ of the resulting $(t,s)$-sequence is bounded by the sum of the polynomial degrees minus one.
A simple way to obtain $(0,b)$-sequences in base $b$ consists of using the first $b$ irreducible polynomials $p_0(x)=x$, $p_1(x)=x+1$, $\dots$, $p_{b-1}(x)=x+b-1$, each of degree $1$.
The theorem by Faure and Lemieux then shows that $0 \leq t \leq \sum_i \left(deg(p_i)-1\right) = 0$.
However, this produces a unique sequence (up to sample permutations), related to those produced by Faure~\shortcite{Faure1982}, which limits its use in more general settings that require sample diversity.
Ostromoukhov et al.~\shortcite{Quad2024} used a construction with quadruplets of irreducible polynomials in base $b=3$ to achieve progressive point sets of excellent consecutive projections.
\paragraph{LDS and projective LDS in Computer Graphics.} In rendering applications, low-discrepancy sequences can have a significant impact on path-tracing performance~\cite{keller2004myths,keller2013quasi,christensen2018renderman,jarosz2019orthogonal}.
When the sampling pattern defined on the canonical domain $[0,1)^s$ is mapped to a pixel (or a group of pixels), decorrelating the pattern across different pixels typically requires a scrambling procedure.
Owen's scrambling is usually considered, as it preserves the $t$ value of the point set \cite{owen1995randomly}.
Due to the nature of the rendering equation, several authors have explored projective strategies aimed at achieving highly uniform consecutive 2D projections.
Achieving high-quality in 2D projections often comes at the cost of degrading uniformity in higher dimensions~\cite{kollig2002efficient,perrier2018sequences,paulin2021cascaded,ahmed2020screen}. Notably, the ZeroTwo sequence uses the first two Sobol' dimensions repeatedly with random permutations~\cite{pharr2023physically}, while padded 4D Sobol' repeats and shuffles samples of the first four dimensions~\cite{burley2020practical}. These provide ideal behavior in consecutive 2 or 4D projections, but behave similarly to white noise in higher dimensions.
Some methods are dedicated to generating point sets rather than sequences~\cite{matbuilder22,Quad2024}, or are not low discrepancy \cite{reinert2016projective,paulin2020sliced}.
Our new construction enables the definition of complete $(t,s)-$sequences with guaranteed high-quality (\emph{i.e.} $(1,2)$-sequences) 2D projections.
\section{Overview of our construction}
We first focus on 2-dimensional Sobol' sequences. Our goal is thus to obtain Sobol' matrices $M_{p_0}$ and $M_{p_1}$ for two polynomials $p_0$ and $p_1$. In our framework, we may use standard Sobol' construction to generate these matrices using respectively polynomials $p_0$ and $p_1$ and initialization matrices $D_{p_0}$ and $D_{p_1}$: our contribution is to provide simple conditions for these initialization matrices to yield high-quality parameter $t=1$ Sobol' sequences (see Fig.~\ref{fig:mainIteration}).
This section describes an overview of this process, while Sec.~\ref{sec:mathconstruction} details proofs. In the following, we assume modulo 2 arithmetic, and all values involved are binary. % (we then have $x+x=0$ and $-x=x$).
We consider a matrix $K = M_{p_1} M_{p_0}^{-1}$ that uniquely represents a 2D Sobol' sequence up to a permutation of points, called characteristic matrix (denoted $C$ by Ahmed et al.~\shortcite{ahmed2023analysis}).
We show that when $p_1 = p_0^2+p_0+1$, where $p_0$ is a degree $e$ polynomial and $p_1$ thus has degree $2e$, $K$ has a peculiar form, and can be obtained recursively from a decomposition into square blocks $A$, $B$ and $C$: % i starts at 1 here
$$
K^{(i)} =
\left[\begin{array}{cc}
A & B\\
0 & C
\end{array}\right]
\rightarrow K^{(i+1)} =
\left[\begin{array}{cc|cc}
A & B & A+B & A\\
0& C & C & 0 \\ \hline
0&0& A & A+B\\
0& 0&0 & C
\end{array}\right]\,.
$$
where each block $A$, $B$ or $C$ has size $2^{i-1}e$. This produces a matrix $K$ of arbitrary size, doubling its size at each iteration. Our paper introduces conditions on $K^{(1)}$ and $K^{(2)}$, that lead to conditions on the Sobol' initialization matrices $D_{p_0}$ and $D_{p_1}$, for our Sobol' sequences to be $(1,2)$-sequences.
We note that the initialization $K^{(1)}$, of size $2e \times 2e$, is $ K^{(1)}= D_{p_1} \widetilde{D}_{p_0}^{-1}$, where we denote $\widetilde{D}_{p_0}$ the $2e \times 2e$ matrix obtained by extending the $e \times e$ matrix $D_{p_0}$ to $2e \times 2e$ using standard Sobol' iterations.
To obtain Sobol' initialization matrices, we thus first generate a random invertible triangular $e \times e$ matrix $D_{p_0}$ which we extend to $2e \times 2e$ using Sobol' iterations. We then use a matrix $K^{(1)}$ satisfying our conditions (see next), and easily obtain $2e \times 2e$ initialization matrix $D_{p_1} = K^{(1)} \widetilde{D}_{p_0}$. Matrices $M_{p_0}$ and $M_{p_1}$, and 2D sample coordinates are then obtained using standard Sobol' procedures, from $p_0$, $p_1=p_0^2+p_0+1$, $D_{p_0}$ and $D_{p_1}$.
\begin{figure}
\definecolor{dx}{RGB}{255, 0, 255}
\definecolor{ddx}{RGB}{255, 128, 128}
\definecolor{dy}{RGB}{128, 0, 128}
\begin{tikzpicture}
\node (itera) at (0,0) {\includegraphics[width=\linewidth]{figures_siggraph25/mainIteration.pdf}};
\node[text=dx] (dp) at (-4, -1) {$D_{p}$};
\node[text=ddx] (dpp) at (-2.7, -1) {$D_{p^2}$};
\node (dppp) at (-1.7, -2) {${\color{dy}D_{p^2+p+1}} = K^{(1)}D_{p^2}$};
\node (k1) at (-3.5,2.2) {$K^{(1)}$};
\node (k2) at (-1.6,2.2) {$K^{(2)}$};
\node (k3) at (1.8,2.2) {$K^{(3)}$};
\draw [draw=gray,decorate,decoration={brace,amplitude=5pt,raise=4ex}]
(-0.2,2) -- (4.5,2) node[color = gray,midway,yshift=3em]{Only used for proofs (Theorem \ref{theorem:constructionRanks})};
\draw [draw=gray,very thick] (-4.4,2.5) rectangle (-.4,-.3);
\node[text=gray] at (-.85,2.3) {Alg.~\ref{algo:CharacteristicMatrices}};
\draw [draw=gray,very thick] (-4.4,-.6) rectangle (-.4,-2.5);
\node[text=gray] at (-.85,-.8) {Alg.~\ref{algo:ConstructingDxDy}};
\draw [-{Stealth[length=2mm]}](-3.6,.8) to [bend right=15] (dp);
%\draw [-{Stealth[length=2mm]}] (-3.4,.8) to [bend left=15] (dppp);
\draw [-{Stealth[length=2mm]}] (dp) -- (dpp);
\draw [-{Stealth[length=2mm]}] (dpp) -- (dppp);
\draw [draw=gray,very thick] (-4.4,-2.7) rectangle (4,-7);
\node[text=gray] at (2,-2.9) {Classical Sobol' construction};
\node (mp) at (-3.8, -3.5) {$M_{p}$};
\node (m1) at (-2,-5.2) {\includegraphics[width=3cm]{figures_siggraph25/matrixM1.pdf}};
\node (mppp) at (0.2, -3.5) {$M_{p^2+p+1}$};
\node (m2) at (2,-5.2) {\includegraphics[width=3cm]{figures_siggraph25/matrixM2.pdf}};
\draw [-{Stealth[length=2mm]}] (dp) to [bend right=15] (mp);
\draw [-{Stealth[length=2mm]}] (dppp) -- (mppp);
\node (pts) at (0,-9) {\includegraphics[width=3cm]{figures_siggraph25/pts56.pdf}};
\draw [very thick] (1.55,-10.45) -- (1.55,-7.55);
\draw [-{Stealth[length=2mm]}] (m2) to [bend left=25] (1.65,-9);
\draw [very thick] (-1.5,-7.45) -- (1.5,-7.45);
\draw [-{Stealth[length=2mm]}] (m1) to [bend left=15] (0,-7.4);
\end{tikzpicture}
\caption{Overview. We introduce a recursive construction of the characteristic matrix associated with a pair of polynomials ($p,p^2+p+1$). We use it in proofs to obtain conditions for generating $(1,2)$-sequences based on the first iteration alone of this recursion.
From characteristic matrices meeting these conditions, we
derive Sobol' initialization matrices $D_{p}$ and $D_{p^2+p+1}$, which in turn allows to construct the corresponding Sobol' matrices $M_{p}$ and $M_{p^2+p+1}$ generating $(1,2)-$sequences in base 2.}\label{fig:mainIteration}
\end{figure}
To generate higher-dimensional Sobol' sequences, we combine pairs of dimensions but further require that $p_0$ and $p_1 = p_0^2+p_0+1$ be irreducible polynomials so as to guarantee that the resulting $s$-dimensional sequence is a $(t,s)$-sequence~\cite{faure2016irreducible}.
We claim the following contributions.
\begin{theorem}\label{theorem:construction}
The sequence of iterations
\begin{equation}
K^{(i)} =
\left[\begin{array}{cc}
A & B\\
0 & C
\end{array}\right]
\rightarrow K^{(i+1)} =
\left[\begin{array}{cc|cc}
A & B & A+B & A\\
0& C & C & 0 \\ \hline
0&0& A & A+B\\
0& 0&0 & C
\end{array}\right],
\label{eq:KiKiplusone}
\end{equation}
where each matrix block $A$, $B$, $C$, is of size $2^{i-1}e \times 2^{i-1}e$, produces the characteristic matrix of a 2D Sobol' sequence given by a degree $e$ polynomial $p$, and degree $2e$ polynomial $p^2+p+1$.
\end{theorem}
\begin{corollary}\label{corollary:coeffs}
A 2D Sobol' sequence given by polynomials $p$ and $p^2+p+1$ only depends on the degree of the polynomial and initialization matrices, and does not depend on the coefficients of $p$ themselves, up to a permutation of samples.
\end{corollary}
This corollary is readily justified since the recursive construction of theorem~\ref{theorem:construction} does not involve polynomial coefficients, but merely polynomial degrees.
As such, a pair $p_0$ of degree $e$ and $p_1 = p_0^2+p_0+1$ with given initialization matrices would result in the same sequence as $p_0'$ of degree $e$ and $p_1' = p_0'^2+p_0'+1$ for \textit{another pair} of initialization matrices.
%$p_0(x)=x^e$ and $p_1(x)=x^{2e}+x^e+1$ produces the same Sobol' sequence of samples as any other pair of polynomials of the same degree satisfying $p_1=p_0^2+p_0+1$ if initialization matrices are the same. \nb{PAS BON! A REVOIR PROPREMENT}
\begin{theorem}\label{theorem:constructionRanks}
Iterations $K^{(i)} \rightarrow K^{(i+1)}$ characterize Sobol' $(1,2)$-sequences if and only if both conditions are met:
\begin{itemize}
\item $(\mathcal{P})$: $\emph{\text{corank}}(T) \leq 1$ for any rectangular $(w-1) \times w$ submatrix $T$ of $K^{(2)}$ anchored at its first row
\item $(\mathcal{Q})$: $\emph{\text{corank}}(C') \leq 1$ for any square submatrix $C'$ of block $C$ in $K^{(1)}$ obtained by removing any consecutive set of $k$ columns and the last $k$ rows of $C$.
\end{itemize}
\end{theorem}
These results allow us to pre-compute many matrices $K^{(1)}$ satisfying the conditions of theorem~\ref{theorem:construction} for a given polynomial degree $e$ and infer pairs of irreducible polynomials and initialization matrices that produce Sobol' $(1, 2)$-sequences.
In the process of proving these theorems, we discovered new insights on Sobol' constructions for more general polynomials:
\begin{itemize}
\item Polynomials $p$ and $p^2$ will produce the same Sobol' matrices, given proper initialization matrices (see lemma~\ref{Lemma1:doubling}).
\item This property allows building general Sobol' matrices in a recursive way by doubling their size at each iteration via polynomial squaring.
\item The characteristic matrix can also be constructed recursively.
\end{itemize}
The goal of Sec.~\ref{sec:mathconstruction} is to provide mathematical proofs for the claims we summarized in this overview section. The practitioner may thus skip Sec.~\ref{sec:mathconstruction} at first read and jump to the description of our algorithms in Sec.~\ref{sec:experimental}.
Specifically, sections ~\ref{sec:blocklemieux} and \ref{sec:squaring} prove lemma~\ref{Lemma1:doubling} related to polynomials $p$ and $p^2$ yielding the same matrices. This in turns helps proving in Sec.~\ref{sec:blockchar} that characteristic matrices can be constructed recursively. When applied to polynomials $p$ and $p^2+p+1$, Sec.~\ref{sec:pp2p1} proves that the characteristic matrix has the recursive construction of Eq.~\ref{eq:KiKiplusone} and thus proves Theorem~\ref{theorem:construction}. Finally, Sec.~\ref{sec:ranks} and our supplementary materials prove Theorem~\ref{theorem:constructionRanks} that explicits conditions under which the characteristic matrices generate $(1, 2)$-sequences.
\section{Construction of $(1,2)$-sequences in base $2$}
\label{sec:mathconstruction}
We first recall a construction of Sobol' matrices based on matrix blocks by Faure and Lemieux~\shortcite{faure2016irreducible} in Sec.~\ref{sec:blocklemieux}, which will serve as a basis for the next sections describing our proof. Our proof first consists of introducing a new recursive construction of Sobol' matrices by squaring polynomials in Sec.~\ref{sec:squaring}. We then show that a similar squaring procedure can be obtained for characteristic matrices in Sec.~\ref{sec:blockchar}. We then show, using this construction, that the characteristic matrix for polynomials $p$ and $p^2+p+1$ has a specific form exhibiting a self-similar pattern in Sec.~\ref{sec:pp2p1}. We finally show that ranks of characteristic matrices with this self-similar pattern are necessarily such that the produced 2D sequences are $(1,2)$-sequences in Sec.~\ref{sec:ranks}.
\subsection{Block-based recursive construction}
\label{sec:blocklemieux}
We first briefly describe a block-based 1-D Sobol' construction as described by Faure and Lemieux~\shortcite{faure2016irreducible}.
For a given irreducible polynomial $p(x) = x^e + \sum_{i=0}^{e-1} a_i x^i$ and upper $e\times e$ invertible triangular initialization matrix $D_p$, Faure and Lemieux~\shortcite{faure2016irreducible} rewrite Sobol' iterations in terms of block matrices:
$$
M_p =
\begin{pmatrix}
B_{1,1} & B_{1,2} & B_{1,3} & \dots \\
0 & B_{2,2} & B_{2,3} & \dots \\
0 & 0 & B_{3,3} & \dots \\
\vdots & \vdots & \vdots & \ddots
\end{pmatrix},
$$
where the blocks $B_{i,j}$, of size $e\times e$, are defined according to the following recursive procedure:
$$
B_{1,1} = D_p ; \hspace{10pt} B_{i,i} = D_p F_p^{i-1} \\
$$
\begin{equation}
B_{i,j} = \left\{
\begin{array}{lll}
B_{i,j-1} Q_p F_p && \text{when} \hspace{5pt} i=1 \\
B_{i,j-1} Q_p F_p + B_{i-1,j-1} F_p && \text{elsewhere}.
\end{array}
\right.
\label{eq:Blr}
\end{equation}
Here, $Q_p$ is an $e\times e$ lower triangular Toeplitz matrix involving the coefficients $(a_i)_{i}$ of polynomial $p$:
\begin{equation}
Q_p =
\begin{pmatrix}
a_{0} & 0 & 0 & \dots & 0 \\
a_{1} & a_{0} & 0 & \dots & 0 \\
a_{2} & a_{1} &a_{0} & \dots & 0 \\
\vdots & \vdots & \vdots & \ddots & \vdots \\
a_{e-1} & a_{e-2} & a_{e-3} & \dots & a_{0}
\end{pmatrix}.
\label{eq:Qp}
\end{equation}
Matrix $F_p$ of size $e\times e$ is defined as
\begin{equation}
F_p = ({Id}_{\small e} + R_{p,2} ) ({Id}_{\small e} + R_{p,3} ) \dots ({Id}_{\small e} + R_{p,e} ),
\label{eq:F}
\end{equation}
where $Id_{\small e}$ is an identity matrix of size $e\times e$, and $R_{p, k}$ are matrices of size $e\times e$
with zeros everywhere except in the first $k-1$ entries of the $k$-th column, given by the coefficients $(a_{e-(k-1)},\dots ,a_{e-1})$ of polynomial $p$.
We introduce a new recursion to build Sobol' matrices inspired by the construction of Faure and Lemieux~\shortcite{faure2016irreducible}.
\subsection{Sobol' construction by squaring polynomials}
\label{sec:squaring}
In the following, we introduce a squared superscript notation to clarify matrix sizes when appropriate, e.g., $M_{p}^\btwoe$ denoting the Sobol' matrix of polynomial $p$ restricted to the first $2e$ rows and $2e$ columns.
We observe that the Sobol' matrix $M_{p^2}$ of a squared polynomial (although not irreducible) is identical to the Sobol' matrix $M_p$ of the original polynomial, provided that the initialization matrix $D_{p^2}$ coincides with the top-left corner of $M_p$. We formalize this:
\begin{lemma}\label{Lemma1:doubling}
The Sobol' matrix $M_p$ generated by a polynomial $p$ of degree $e$ and initialization matrix $D_p$ is identical to the Sobol' matrix $M_{p^2}$ of polynomial $p^2$ and initialization matrix
\begin{equation}
D_{p^2} = M^\btwoe_{p} = \begin{pmatrix}
D_p & 0 \\
0 & D_p
\end{pmatrix} \begin{pmatrix}
Id_{\small e} & Q_p F_p \\
0 & F_p
\end{pmatrix},\label{eq:lemma41}
\end{equation}
corresponding to the top-left $2e \times 2e$ submatrix of $M_p$.
Matrices $Q_{p^2}$ and $F_{p^2}$ of Faure and Lemieux~\shortcite{faure2016irreducible} can be obtained by applying a Kronecker product (or tensor product) with the matrix ${Id}_2 = {\begin{pmatrix} 1 & 0 \\ 0 & 1 \end{pmatrix}}$ to matrices $Q_p$ and $F_p$:
\begin{equation}
Q_{p^2} = Q_{p} \otimes {Id}_2, \hspace{20pt} F_{p^2} = F_{p} \otimes {Id}_2.
\label{eq:FQ-prime}
\end{equation}
where $Q_{p^2}$ and $F_{p^2}$ are of size $2e \times 2e$.
\end{lemma}
\begin{figure*}
\begin{tikzpicture}
\draw[very thick] (2,-1.5) rectangle (3.5,0) node[pos=.5] {0};
\draw[very thick] (3.5,-1.5) rectangle (5,0) node[pos=.5] {${\color{teal}M_p^\be} Y^{(1)}$};
\draw [very thick](3.5,0) rectangle (5,1.5) node[pos=.5] {${\color{teal}M_p^\be} X^{(1)}$};
\draw[fill=teal,color=teal,very thick] (2,0) rectangle (3.5,1.5) node[pos=.5] {\color{white}\small$M_p^\be=D_p$};
\draw[very thick] (5,-1.5) rectangle (8,1.5) node[pos=.5] {${\color{magenta}M_p^\btwoe} X^{(2)}$};
\draw[very thick] (5,-4.5) rectangle (8,-1.5) node[pos=.5] {${\color{magenta}M_p^\btwoe} Y^{(2)}$};
\draw[very thick] (2,-4.5) rectangle (5,-1.5) node[pos=.5] {0};
\node at (3.5,2) {${\color{magenta}M_p^\btwoe}$};
\draw[color=magenta,very thick] (2,-1.5) rectangle (5,1.5);
\node at (8.5,1.5) {\ldots\ldots};
\node at (2.0,-4.9) {\rotatebox{90}{\ldots\ldots}};
\node at (3.5,-5) {$X^{(1)} = Q_pF_p$};
\node at (6.8,-5) {$X^{(i)} = X^{(i-1)}\otimes Id_2$};
\node at (3.3,-5.5) {$Y^{(1)} = F_p$};
\node at (6.8,-5.5) {$Y^{(i)} = Y^{(i-1)}\otimes Id_2$};
\draw [color=gray,decorate,decoration={brace,amplitude=5pt,raise=.5ex}]
(2,0) -- (2,1.5) node[midway,xshift=-1em]{$e$};
\draw [color=gray,decorate,decoration={brace,amplitude=5pt,raise=1ex}]
(1.8,-1.5) -- (1.8,1.5) node[midway,xshift=-1.5em]{$2e$};
\draw [color=gray,decorate,decoration={brace,amplitude=5pt,raise=1ex}]
(1.4,-4.5) -- (1.4,1.5) node[midway,xshift=-1.5em]{$4e$};
\end{tikzpicture}\hspace{.5cm}
\begin{tikzpicture}
\draw[very thick] (2,-1.5) rectangle (3.5,0) node[pos=.5] {0};
\draw[very thick] (3.5,-1.5) rectangle (5,0) node[pos=.5] {${\color{teal}K_{q,r}^\bet} Y'^{(1)}$};
\draw[very thick] (3.5,0) rectangle (5,1.5) node[pos=.5] {${\color{teal}K_{q,r}^\bet} X'^{(1)}$};
\draw[fill=teal,color=teal,very thick] (2,0) rectangle (3.5,1.5) node[pos=.5] {\tiny\color{white}$K_{q,r}^\bet=D_rD_q^{-1}$};
\draw[very thick] (5,-1.5) rectangle (8,1.5) node[pos=.5] {${\color{magenta}K_{q,r}^\btwoet} X'^{(2)}$};
\draw[very thick] (5,-4.5) rectangle (8,-1.5) node[pos=.5] {${\color{magenta}K_{q,r}^\btwoet} Y'^{(2)}$};
\draw[very thick] (2,-4.5) rectangle (5,-1.5) node[pos=.5] {0};
\node at (3.5,2) {${\color{magenta}K_{q,r}^\btwoet}$};
\draw[very thick,color=magenta,very thick] (2,-1.5) rectangle (5,1.5);
\node at (8.5,1.5) {\ldots\ldots};
\node at (2.0,-4.9) {\rotatebox{90}{\ldots\ldots}};
\node at (4.2,-5) {$X'^{(1)} = D_{q} (Q_{q} + Q_{r}F_{r}F_{q}^{-1} ) D_{q}^{-1}$};
\node at (7.8,-5) {$X'^{(i)} = X'^{(i-1)}\otimes Id_2$};
\node at (3.7,-5.5) {$Y'^{(1)} = D_{q} ( F_{r}F_{q}^{-1} ) D_{q}^{-1}$};
\node at (7.8,-5.5) {$Y'^{(i)} = Y'^{(i-1)}\otimes Id_2$};
\draw [color=gray,decorate,decoration={brace,amplitude=5pt,raise=.5ex}]
(2,0) -- (2,1.5) node[midway,xshift=-1em]{$e$};
\draw [color=gray,decorate,decoration={brace,amplitude=5pt,raise=1ex}]
(1.8,-1.5) -- (1.8,1.5) node[midway,xshift=-1.5em]{$2e$};
\draw [color=gray,decorate,decoration={brace,amplitude=5pt,raise=1ex}]
(1.4,-4.5) -- (1.4,1.5) node[midway,xshift=-1.5em]{$4e$};
\end{tikzpicture}
\caption{In contrast to the column-by-column Sobol’ approach~\cite{Sobol1967distribution}, and the block formulation of Faure and Lemieux~\shortcite{faure2016irreducible} (Sec.~\ref{sec:blocklemieux}), we present a novel recursive construction for Sobol’ matrices (left) and characteristic matrices (right) by squaring polynomials.
The construction of matrix $M_p$ follows from lemma~\ref{Lemma1:doubling}, Eq.~(\ref{eq:lemma41}), while the construction of $K^\btwoet_{q,r}$ is obtained from Eq.~(\ref{eq:C}). Matrices $X^{(i)}$, $Y^{(i)}$, $X'^{(i)}$ and $Y'^{(i)}$ provide compact representations of expressions involving $D_p$, $F_p$, and $Q_p$, using the distributivity of the Kronecker product over matrix multiplication.}
\label{fig:doublingMK}
\end{figure*}
This lemma brings a new recursive construction of Sobol' matrices, doubling their sizes at each iteration by squaring polynomials, illustrated in Fig.~\ref{fig:doublingMK}.
%To double the size of an $m \times m$ Sobol' matrix, where $m=2^k e$, we consider this matrix as an initialization matrix $D^{\small [m]}_{p^{m/e}}$ for a polynomial of degree $m$, and use $Q^{\small [m]}_{p^{m/e}} = Q^\be_{p} \otimes {Id}_{m/e}$ and $F^{\small [m]}_{p^{m/e}} = F^\be_{p} \otimes {Id}_{m/e}$ in a single iteration of the block-based formulation of Faure and Lemieux~\shortcite{faure2016irreducible} to obtain the corresponding three $m \times m$ blocks required to double its size.
In our development, we first note that, while $F_p$ can have a complicated form, its inverse can be expressed much more easily, as an upper triangular Toeplitz matrix:
\begin{equation}
F_p^{-1} =
\begin{pmatrix}
1 & a_{e-1} & a_{e-2} & \dots & a_{1} \\
0 & 1 & a_{e-1} & \dots & a_{2} \\
0 & 0 &1 & \dots & a_{3} \\
\vdots & \vdots & \vdots & \ddots & \vdots \\
0 & 0 & 0 & \dots & 1
\end{pmatrix},
\label{eq:Q}
\end{equation}
where we ignore signs in modulo 2 arithmetic. This is obtained by observing that matrices ${Id}_e + R_{p,k}$ are their own inverses, and that inverting $F_p$ then merely reverses the order of the multiplication: $F_p^{-1} = ({Id}_{\small e} + R_{p,e} ) \dots ({Id}_{\small e} + R_{p,2} )$.
\begin{proof}[Proof of Lemma \ref{Lemma1:doubling}]
The relation between $D_{p^2}$ and $D_{p}$ is simply obtained by running Faure and Lemieux' iterations for one block of columns.
We may now assume that the top-left $2e \times 2e$ submatrices of $M_p$ and $M_{p^2}$ coincide.
We also note that for a polynomial $p(x) = x^e + \sum_{i=0}^{e-1} a_i x^i$, given modulo-2 arithmetic cancelling odd degrees, $p^2(x) = x^{2e} + \sum_{i=0}^{e-1} a_i x^{2i}$.
This, in turn, leads to $Q_{p^2} = Q_{p} \otimes Id_2$, and similarly, to $F^{-1}_{p^2} = F^{-1}_{p} \otimes {Id}_2$ and thus $F_{p^2} = F_{p} \otimes {Id}_2$ (the inverse of the Kronecker product being the Kronecker product of the inverse matrices), and finally, $Q_{p^2} F_{p^2} = (Q_{p} F_{p})\otimes {Id}_2$.
The recursive construction of Eq.~(\ref{eq:Blr}) thus produces the same matrices whether using $Q_{p^2}$ and $F_{p^2}$ or $Q_{p}$ and $F_{p}$.
\end{proof}
\subsection{Block-based characteristic matrices}
\label{sec:blockchar}
For a pair of dimensions with Sobol' matrices $M_{q}$ and $M_{r}$, we base our analysis on the characteristic matrix $K = M_{r} M_{q}^{-1}$ as defined by Ahmed et al.~\shortcite{ahmed2023analysis}
, which uniquely characterizes Sobol' 2D sequences up to a permutation of samples.
When the polynomials $q$ and $r$ are of the same degree $\tilde{e}$, we may build a characteristic matrix of size $2\tilde{e}$ by applying one iteration of Faure and Lemieux' block construction:
\begin{align}
K^{\btwoet}_{q, r} &= M^\btwoe_{r} \left(M^{\btwoe}_{q}\right)^{-1}\nonumber \\
&= {\begin{pmatrix} D_{r} & 0 \\ 0 & D_{r}\end{pmatrix}} {\begin{pmatrix} Id_{\small \tilde{e}} & Q_{r} F_{r} \\ 0 & F_{r}\end{pmatrix}} \left( {\begin{pmatrix} D_{q} & 0 \\ 0 & D_{q}\end{pmatrix}} {\begin{pmatrix} Id_{\small \tilde{e}} & Q_{q} F_{q} \\ 0 & F_{q}\end{pmatrix}} \right)^{-1} \nonumber \\
&= {\begin{pmatrix} D_{r} & 0 \\ 0 & D_{r}\end{pmatrix}} {\begin{pmatrix} Id_{\small \tilde{e}} & Q_{q} + Q_{r} F_{r} F_{q}^{-1} \\ 0 & F_{r} F_{q}^{-1}\end{pmatrix}} {\begin{pmatrix} D_{q} & 0 \\ 0 & D_{q}\end{pmatrix}}^{-1} \nonumber \\
%= {\begin{pmatrix} D_{r} & 0 \\ 0 & D_{r}\end{pmatrix}} {\begin{pmatrix} R_{11} & R_{12} \\ 0 & R_{22} \end{pmatrix}} {\begin{pmatrix} D_{q}^{-1} & 0 \\ 0 & D_{q}^{-1}\end{pmatrix}} \\
&= {\begin{pmatrix} K^{\bet}_{q, r} & 0 \\ 0 & K^{\bet}_{q, r}\end{pmatrix}} {\begin{pmatrix} Id_{\small \tilde{e}} & D_{q} (Q_{q} + Q_{r}F_{r}F_{q}^{-1} ) D_{q}^{-1} \\ 0 & D_{q} ( F_{r}F_{q}^{-1} ) D_{q}^{-1}\end{pmatrix}.}
\label{eq:C}
\end{align}
where we used the identity $
\begin{array}{lll}
\begin{pmatrix} A & B \\ 0 & C\end{pmatrix}^{-1} =
\begin{pmatrix} A^{-1} & -A^{-1} B C^{-1} \\
0 & C^{-1} \end{pmatrix}
\end{array}
$ for invertible $A$ and $C$, and where $K^{\bet}_{q, r} = D_{r} D_{q}^{-1}$ is the $\tilde{e} \times \tilde{e}$ initialization block for this characteristic matrix, obtained from the initializations of $M_{q}$ and $M_{r}$.
This construction considers polynomials $q$ and $r$ of the same degree, but we can use lemma~\ref{Lemma1:doubling} to square our first polynomial and then consider $q = p^2$ and $r = p^2+p+1$ of the same degree. We thus consider this construction using $\tilde{e}=2e$.
Also, this construction merely produces a matrix of size $2\tilde{e}$ given polynomials of degree $\tilde{e}$. However, it can be used recursively by considering lemma~\ref{Lemma1:doubling}, doubling the size of the matrix at each iteration by squaring polynomials. By lemma~\ref{Lemma1:doubling}, the effect of squaring polynomials on all matrices involved is merely a tensor product with $Id_2$. We also illustrate this process in Fig.~\ref{fig:doublingMK}.
\subsection{The special case $(p, p^2+p+1)$}
\label{sec:pp2p1}
The construction for characteristic matrices in Sec.~\ref{sec:blockchar} is general, and applies to any pair of polynomials. In this section, we show the special case when $q = p^2$, and $r = p^2+p+1$.
First, we note that applying lemma~\ref{Lemma1:doubling} for using ${q} = p^2$ in place of ${q} = p$ leads to an extended initialization matrix $D_{p^2}$ which inverse can be expressed:
\begin{equation}
D_{p^2}^{-1} = \left(M^{\btwoe}_p\right)^{-1} = \begin{pmatrix}
Id_{\small e} & Q_p \\
0 & F_p^{-1}
\end{pmatrix}\begin{pmatrix}
D_p^{-1} & 0 \\
0 & D_p^{-1}
\end{pmatrix}. \label{eq:invDp2}
\end{equation}
We thus now consider that our polynomials have degree $2e$ and we seek to apply our characteristic matrix construction to obtain a matrix $K^{\bfoure}_{q, r}$ of size $4e \times 4e$.
It can then be shown (see Appendix \ref{app:proofsR}) that $F_{p^2+p+1} F_{p^2}^{-1}$ has a particular form:
\begin{equation}
F_{p^2+p+1} F_{p^2}^{-1} = \begin{pmatrix}
Id_{\small e} & F_p \\
0 & Id_{\small e}
\end{pmatrix}.\label{eq:FrinvFp}
\end{equation}
The resulting matrix is also its own inverse: $F_{p^2+p+1} F_{p^2}^{-1} = F_{p^2} F_{p^2+p+1}^{-1}$.
Combining Eq.~(\ref{eq:invDp2}) and (\ref{eq:FrinvFp}) (see Appendix \ref{app:proofs1}), we have:
\begin{align}
F_{p^2+p+1} F_{p^2}^{-1} D_{p^2}^{-1} &= D_{p^2}^{-1} \begin{pmatrix}
Id_{\small e} & Id_{\small e} \\
0 & Id_{\small e}
\end{pmatrix}.\label{eq:FrinvFpDp2}
\end{align}
Similarly, we have
\begin{equation}
(Q_{p^2} + Q_{p^2+p+1} F_{p^2+p+1} F_{p^2}^{-1}) D_{p^2}^{-1} = D_{p^2}^{-1} \begin{pmatrix}
Id_{\small e} & Id_{\small e} \\
Id_{\small e} & 0
\end{pmatrix},\label{eq:Qp2plus}
\end{equation}
where we further use the properties of the product of our Toeplitz matrices, see Appendix \ref{app:proofs2} for details.
Putting all together, and considering the initially given matrix $K^{\btwoe}_{p^2, p^2+p+1}$ has size $2e$ since we considered $p^2$ and $p^2+p+1$ of degree $2e$, we see that the recursive construction of $K$ becomes
\begin{equation}
\begin{array}{lll}
K^{\bfoure}_{p^2, p^2+p+1}
= {\begin{pmatrix} K^{\btwoe}_{p^2, p^2+p+1} & 0 \\ 0 & K^{\btwoe}_{p^2, p^2+p+1}\end{pmatrix}} {\begin{pmatrix} Id_{\small 2e} &\begin{pmatrix}
Id_{\small e} & Id_{\small e} \\
Id_{\small e} & 0
\end{pmatrix} \\ 0 & \begin{pmatrix}
Id_{\small e} & Id_{\small e} \\
0 & Id_{\small e}
\end{pmatrix} \end{pmatrix}}\,.
\end{array}
\label{eq:C2}
\end{equation}
Lemma~\ref{Lemma1:doubling} indicates that this procedure becomes recursive by squaring polynomials, allowing to double the size of $K$ at each iteration. We rewrite these iterations using the following notation involving blocks $A$, $B$, $C$ of size $2^ie \times 2^ie$, doubling their size at each iteration:
$$
K^{(i)} =
\left[\begin{array}{cc}
A & B\\
0 & C
\end{array}\right]
\rightarrow K^{(i+1)} =
\left[\begin{array}{cc|cc}
A & B & A+B & A\\
& C & C & 0 \\ \hline
& & A & A+B\\
& & & C
\end{array}\right]
$$
with the $2e \times 2e$ initial matrix
\begin{align*}
K^{(1)} &= K^{\btwoe}_{p^2, p^2+p+1} = D_{p^2+p+1} D_{p^2}^{-1} \\
&=
D_{p^2+p+1}
\left(\begin{pmatrix}
D_{p} & 0 \\
0 & D_{p}
\end{pmatrix} \begin{pmatrix}
Id_{\small e} & Q_p F_p \\
0 & F_p
\end{pmatrix}\right)^{-1}, \end{align*}
for any given $e \times e$ upper triangular invertible matrix $D_p$ and $2e \times 2e$ upper triangular invertible matrix $D_{p^2+p+1}$ giving the degrees of freedom for the generated sequences. $D_{p^2}$ is here obtained by running standard Sobol' iterations to extend the $e \times e$ matrix $D_{p}$ to obtain the next $e$ rows and columns.
\subsection{Rank of submatrices}
\label{sec:ranks}
Let us denote by $\bar{T}^{j,w}$ a square $w \times w$ submatrix of $K$
starting at column $1\leq j2$ within our framework to discover $t=0$ sequences in higher dimensions, which remains a gold standard for numerical integration.
\begin{acks}
This work was partially funded by ANR-20-CE45-0025 (MoCaMed), by ANR-22-CE46-000 (StableProxies), and donations from Adobe Inc.
\end{acks}
% Bibliography
\bibliographystyle{ACM-Reference-Format}
\bibliography{references}
% Appendix
\appendix
\section{Additional derivations}
\label{app:details}
\subsection{Proofs of Eq.~(\ref{eq:FrinvFp})}
\label{app:proofsR}
Starting from eq.~\ref{eq:F} and $F_p^{-1} = ({Id}_{\small e} + R_{p,e} ) \dots ({Id}_{\small e} + R_{p,2} )$, we have:
\begin{align*}
&F_{p^2+p+1} F_{p^2}^{-1}=\\
&({Id}_{\small 2e} + R_{p^2+p+1,2} ) \dots ({Id}_{\small 2e} + R_{p^2+p+1,2e} )({Id}_{\small 2e} + R_{p^2,2e} ) \dots ({Id}_{\small 2e} + R_{p^2,2} )
\end{align*}
To simplify the notations, we denotes $R'_{k}=({Id}_{\small 2e} + R_{p^2+p+1,k})$ and $R''_{k}=({Id}_{\small 2e} + R_{p^2,k})$. Hence, we have:
\begin{align}
&F_{p^2+p+1}F_{p^2}^{-1}= \nonumber\\
&\quad \underbrace{R'_{2} \dots R'_{e-1}}_{(i)}\underbrace{R'_{e} \dots R'_{2e-1}R'_{2e} R''_{2e} R''_{2e-1}\dots R''_{e}}_{(ii)}\underbrace{R''_{e-1}\dots R''_{2}}_{(iii)}\label{eq:appcases}
\end{align}
In the following, we will use this illustration for $R'_{k}$ (the column of index $k$ contains the $(k-1)$ highest degree coefficients of $p^2+p+1$):
\begin{center}
\begin{overpic}[width=2cm]{figures_siggraph25/blockmatbis.pdf}
\put(30,60){$Id$}
\put(80,10){$Id$}
\put(82,60){0}
\put(61,49){\rotatebox{90}{$p^2+p$}}
\put(30,12){0}
\put(65,12){0}
\put(65,27.5){1}
\end{overpic}
\end{center}
Note that by definition of $R'_k$ and $R''_k$ matrices, we can only consider polynomials $p^2+p$ and $p^2$ respectively, as the constant factor is dropped by construction.
Let us first consider the first innermost product in part $(ii)$ of Eq.~(\ref{eq:appcases}) involving the $p^2+p$ and $p^2$ polynomials $R'_{2e} R''_{2e}$:
\begin{center}
\begin{overpic}[width=2cm]{figures_siggraph25/blockmat.pdf}
\put(40,60){$Id$}
\put(84,40){\rotatebox{90}{$p^2+p$}}
\put(90,5){1}
\put(40,5){0}
\end{overpic}
\begin{overpic}[width=2cm]{figures_siggraph25/blockmat.pdf}
\put(40,60){$Id$}
\put(84,50){\rotatebox{90}{$p^2$}}
\put(90,5){1}
\put(40,5){0}
\end{overpic}
\raisebox{1cm}{=} \begin{overpic}[width=2cm]{figures_siggraph25/blockmatpplusone.pdf}
\put(40,60){$Id$}
\put(87,68){\rotatebox{90}{$p$}}
\put(89,30){0}
\put(90,5){1}
\put(40,5){0}
\put(103,68){$e$}
\end{overpic}
\end{center}
as $p^2+p^2$ coefficients cancel out for rows greater or equal to $e$. We denote by $U_1$ the resulting matrix. Let us now consider the product $U_2 = R'_{2e-1}\,U_1\,R''_{2e-1}$:
\begin{center}
\begin{overpic}[width=2cm]{figures_siggraph25/blockmater.pdf}
\put(40,60){$Id$}
\put(74,38){\rotatebox{90}{$p^2+p$}}
\put(90,2){1}
\put(79,13){1}
\put(79,2){0}
\put(90,13){0}
\put(40,2){0}
\put(40,2){0}
\put(90,50){0}
\end{overpic}
\begin{overpic}[width=2cm]{figures_siggraph25/blockmatpplusone.pdf}
\put(40,60){$Id$}
\put(87,72){\rotatebox{90}{$p$}}
\put(89,30){0}
\put(90,5){1}
\put(40,5){0}
\end{overpic}
\begin{overpic}[width=2cm]{figures_siggraph25/blockmater.pdf}
\put(40,60){$Id$}
\put(74,50){\rotatebox{90}{$p^2$}}
\put(40,13){0}
\put(40,2){0}
\put(90,2){1}
\put(79,13){1}
\put(79,2){0}
\put(90,13){0}
\put(90,50){0}
\end{overpic}
\end{center}
which simplifies to
\begin{center}
\begin{overpic}[width=2cm]{figures_siggraph25/blockmater2.pdf}
\put(40,60){$Id$}
\put(87,64){\rotatebox{90}{$p$}}
\put(72,72){\rotatebox{90}{$p$}}
\put(90,2){1}
\put(75,18){1}
\put(90,18){0}
\put(90,38){0}
\put(75,38){0}
\put(75,2){0}
\put(40,2){0}
\put(40,18){0}
\end{overpic}
\end{center}
If we repeat this process for all triplets of matrices $R'_{k}\, U_{2e-k}\,R''_k$ for the $k$ indices of $(ii)$, we end up with the matrix $U_e$:
\begin{center}
\begin{overpic}[width=2cm]{figures_siggraph25/blockmatU.pdf}
\put(20,70){$Id$}
\put(65,22){$Id$}
\put(65,70){$F^{-1}_p$}
\put(22,22){0}
\end{overpic}
\end{center}
Indeed, for each product $R'_{k}\, U_{2e-k}\,R''_k$, all $p^2$ coefficents vanish, leading to a triangular upper-right block with shifted $p$ coefficients as in Eq.~(\ref{eq:Q}).
Let us now consider the product between $(ii)$ and the $(i)$ and $(iii)$ parts in Eq.~(\ref{eq:appcases}).
First, we observe that
\begin{center}
\raisebox{.9cm}{$R''_{e-1}R''_{e-2}=$}~~~~ \begin{overpic}[width=2cm]{figures_siggraph25/blockmati.pdf}
\put(20,70){$Id$}
\put(65,22){$Id$}
\put(75,70){0}
\put(22,22){0}
\put(40,68){\rotatebox{90}{\small $p^2$}}
\put(35,85){\rotatebox{90}{\small $p^2$}}
\put(48,51){\tiny 1}
\put(40,51){\tiny 0}
\put(40,57){\tiny 1}
\end{overpic}
\end{center}
By doing such products for all matrices of $(iii)$, we obtain an upper-left block which corresponds to the first $e\times e$ entries of $ F^{-1}_{p^2+p}$.
$$
R''_{e-1}\ldots R''_2 = \begin{pmatrix}
F^{-1}_{p^2} & 0 \\
0 & Id_e
\end{pmatrix}\,.
$$
For products in $(i)$, we use the fact that
$$
\left(R'_{2}\ldots {R'}_{e-1}\right)^{-1} = {R'}_{e-1}^{-1}\ldots {R'}_{2}^{-1} = R'_{e-1}\ldots R'_{2} =\begin{pmatrix}
F^{-1}_{p^2+p} & 0 \\
0 & Id_e
\end{pmatrix}\,,
$$
as $R'_{k}$ is its own inverse and using a similar construction as for $(iii)$. Thus, using the inverse of a block matrix, we obtain
$$
R'_{2}\ldots {R'}_{e-1} = \begin{pmatrix}
F_{p^2+p} & 0 \\
0 & Id_e
\end{pmatrix}\,,$$
which equals to $\begin{pmatrix}
F_{p^2} & 0 \\
0 & Id_e
\end{pmatrix}$ as no coefficients for the $p$ term are present in the upper-left block.
We finally have
\begin{align*}
F_{p^2+p+1}F_{p^2}^{-1}&= \begin{pmatrix}
F_{p^2} & 0 \\
0 & Id_e
\end{pmatrix}
\begin{pmatrix}
Id_e & F^{-1}_p \\
0 & Id_e
\end{pmatrix}
\begin{pmatrix}
F^{-1}_{p^2} & 0 \\
0 & Id_e
\end{pmatrix}\\
&=\begin{pmatrix}
F_{p^2} & 0 \\
0 & Id_e
\end{pmatrix}
\begin{pmatrix}
F^{-1}_{p^2} & F^{-1}_p \\
0 & Id_e
\end{pmatrix}\\
&= \begin{pmatrix}
Id_{\small e} & F_p \\
0 & Id_{\small e}
\end{pmatrix}.
\end{align*}
which concludes Eq.~(\ref{eq:FrinvFp}). The last step uses the observation that, for the upper-right block, $F_{p^2}F_p^{-1} = \left(\left(F_{p^2} F^{-1}_{p}\right)^{-1}\right)^{-1} = \left(F_p F_{p^2}^{-1}\right)^{-1} = \left(F_p F_{p}^{-1}F_{p}^{-1}\right)^{-1} = F_p$.
\subsection{Proofs of Eq.~(\ref{eq:FrinvFpDp2})}
\label{app:proofs1}
First, combining Eq.~(\ref{eq:invDp2}) and (\ref{eq:FrinvFp}), we have:
\begin{align}
F_{p^2+p+1} F_{p^2}^{-1} D_{p^2}^{-1} &=\begin{pmatrix}
Id_{\small e} & F_p \\
0 & Id_{\small e}
\end{pmatrix}\begin{pmatrix}
Id_{\small e} & Q_p \\
0 & F_p^{-1}
\end{pmatrix}\begin{pmatrix}
D_p^{-1} & 0 \\
0 & D_p^{-1}
\end{pmatrix} \\
&= \begin{pmatrix}
Id_{\small e} & Q_p + Id_{\small e} \\
0 & F_p^{-1}
\end{pmatrix} \begin{pmatrix}
D_p^{-1} & 0 \\
0 & D_p^{-1}
\end{pmatrix} \label{eq:Gx} \\
&= \begin{pmatrix}
Id_{\small e} & Q_p \\
0 & F_p^{-1}
\end{pmatrix} \begin{pmatrix}
D_p^{-1} & D_p^{-1} \\
0 & D_p^{-1}
\end{pmatrix} \\
&=\begin{pmatrix}
Id_{\small e} & Q_p \\
0 & Id_{\small e}
\end{pmatrix}\begin{pmatrix}
D_p^{-1} & 0 \\
0 & D_p^{-1}
\end{pmatrix}\begin{pmatrix}
Id_{\small e} & Id_{\small e} \\
0 & Id_{\small e}
\end{pmatrix} \\
&= D_{p^2}^{-1} \begin{pmatrix}
Id_{\small e} & Id_{\small e} \\
0 & Id_{\small e}
\end{pmatrix},
\end{align}
leading to Eq.~(\ref{eq:FrinvFpDp2}). Starting from Eq.~(\ref{eq:Gx}), we also have
\begin{align}
F_{p^2+p+1} F_{p^2}^{-1} D_{p^2}^{-1} &= \begin{pmatrix}
Id_{\small e} & Q_p + Id_{\small e} \\
0 & F_p^{-1}
\end{pmatrix} \begin{pmatrix}
D_p^{-1} & 0 \\
0 & D_p^{-1}
\end{pmatrix} \\
&= \begin{pmatrix}
Id_{\small e} & Q_{p+1} \\
0 & F_p^{-1}
\end{pmatrix} \begin{pmatrix}
D_p^{-1} & 0 \\
0 & D_p^{-1}
\end{pmatrix} \label{eq:GxBis}\,,
\end{align}
that will be used later.
\subsection{Proofs of Eq. (\ref{eq:Qp2plus})}
\label{app:proofs2}
Let us now prove the following statement:
$$
(Q_{p^2} + Q_{p^2+p+1} F_{p^2+p+1} F_{p^2}^{-1}) D_{p^2}^{-1} = D_{p^2}^{-1} \begin{pmatrix}
Id_{\small e} & Id_{\small e} \\
Id_{\small e} & 0
\end{pmatrix}\,.
$$
From now on, we make explicit the size of the matrices using $[e]$ or $[2e]$ superscripts. First, by definition, $Q^\btwoe_{p^2}$ is
\begin{equation}
Q^\btwoe_{p^2} =
\begin{pmatrix}
a_{0} & 0 & 0 & \dots & 0 \\
0 & a_{0} & 0 & \dots & 0 \\
a_{1} & 0 &a_{0} & \dots & 0 \\
0 & a_1 & 0 & \dots & 0 \\
\vdots & \vdots & \vdots & \ddots & \vdots\\
a_{e-1} & 0 & a_{e-2} & \dots & a_0
\end{pmatrix}\,.
\end{equation}
We decompose $Q^\btwoe_{p^2}$ into $e\times e$ blocks:
\begin{equation}
Q^\btwoe_{p^2} = \begin{pmatrix}
Q^\be_{p^2} & 0 \\
\widetilde{Q}^\be_{p^2} & Q^\be_{p^2}
\end{pmatrix}.\label{eq:decQtwoe}
\end{equation}
Similar to $Q^\btwoe_{p^2}$, $\widetilde{Q}^\be_{p^2}$ is also Toeplitz.
Furthermore, we have
$$
Q^\be_{p+q} = Q^\be_p + Q^\be_q\quad\text{and}\quad Q^\be_{pq} = Q^\be_p Q^\be_q\,,
$$
for any polynomial $p$ and $q$ of degree $e$. The same holds for $\widetilde{Q}^\be_{p+q}$ and $\widetilde{Q}^\be_{pq}$ matrices.
Now,
%
%
\begin{align}
(Q_{p^2} &+ Q_{p^2+p+1} F_{p^2+p+1} F^{ -1}_{p^2}) D_{p^2}^{-1} =\nonumber\\
&=Q_{p^2} D_{p^2}^{-1} + Q_{p^2+p+1} F_{p^2+p+1} F^{ -1}_{p^2}D_{p^2}^{-1}\nonumber\\
\text{using Eq.~(\ref{eq:invDp2})}& \text{ and \ref{eq:GxBis} with $T=\begin{pmatrix}
D_p^{-1} & 0 \\
0 & D_p^{-1}
\end{pmatrix}$ }\nonumber\\
&=\left(\begin{pmatrix}
Q^\be_{p^2} & 0 \\
\widetilde{Q}^\be_{p^2} & Q^\be_{p^2}
\end{pmatrix}\begin{pmatrix}
Id_{\small e} & Q^{\be}_{p}
\\
0 & \widetilde{Q}^{\be}_{p}
\end{pmatrix} \right.\nonumber\\
&\left.\quad + \begin{pmatrix}
Q^\be_{p^2+p+1} & 0 \\
\widetilde{Q}^\be_{p^2+p+1} & Q^\be_{p^2+p+1}
\end{pmatrix}\begin{pmatrix}
Id_{\small e} & Q^{\be}_{p+1} \\
0 & \widetilde{Q}^{\be}_{p+1}
\end{pmatrix}\right)T\,. \label{eq:Tfactor}
\end{align}
First, we observe that $\widetilde{Q}^{\be}_{p} = \widetilde{Q}^{\be}_{p+1}$. The first factor can be rewritten
\begin{align*}
&\begin{pmatrix}
Q^\be_{p^2} & Q^\be_{p^3} \\
\widetilde{Q}^\be_{p^2} & \widetilde{Q}^\be_{p^2}Q^\be_{p}+Q^\be_{p^2}\widetilde{Q}^\be_{p}
\end{pmatrix}\\
&\quad\quad+
\begin{pmatrix}
Q^\be_{p^2+p+1} & Q^\be_{p^3+1} \\
\widetilde{Q}^\be_{p^2+p+1} & \widetilde{Q}^\be_{p^2+p+1}Q^\be_{p+1}+{Q}^\be_{p^2+p+1}\widetilde{Q}^{\be}_{p+1}
\end{pmatrix}\,,
\end{align*}
since $Q^\be_{p^2}Q^\be_{p}=Q^\be_{p^3}$ and $Q^\be_{p^2+p+1}Q^\be_{p}=Q^\be_{p^3+p^2+p}$. Furthermore, for any polynomial $p$ and $q$ of degree $e$, we have
\begin{align*}
Q^\btwoe_{pq} &= Q^\btwoe_{p}Q^\btwoe_{q}\\
&= \begin{pmatrix}
Q^\be_{p} & 0 \\
\widetilde{Q}^\be_{p} & Q^\be_{p}
\end{pmatrix}\begin{pmatrix}
Q^\be_{q} & 0 \\
\widetilde{Q}^\be_{q} & Q^\be_{q}
\end{pmatrix}\\
&= \begin{pmatrix}
{Q}^\be_{pq} &0 \\
\widetilde{Q}^\be_{p}{Q}^\be_{q} + {Q}^\be_{p}\widetilde{Q}^\be_{q}&{Q}^\be_{pq}
\end{pmatrix}\,.
\end{align*}
Hence, the first factor of Eq.~(\ref{eq:Tfactor}) is
\begin{align*}
\begin{pmatrix}
Q^\be_{p^2} & Q^\be_{p^3} \\
\widetilde{Q}^\be_{p^2} & \widetilde{Q}^\be_{p^3}
\end{pmatrix}
+
\begin{pmatrix}
Q^\be_{p^2+p+1} & Q^\be_{p^3+1} \\
\widetilde{Q}^\be_{p^2+p+1} & \widetilde{Q}^\be_{p^3+1}
\end{pmatrix}
&= \begin{pmatrix}
Q^\be_{p+1} & Q^\be_{1} \\
\widetilde{Q}^\be_{p+1} & \widetilde{Q}^\be_{p^3}+\widetilde{Q}^\be_{p^3+1}
\end{pmatrix}\\
&= \begin{pmatrix}
Q^\be_{p+1} & Id_e \\
F_p^{-1}& 0
\end{pmatrix}\,,
\end{align*}
using the fact that $\widetilde{Q}^{\be}_{p}= F_p^{-1}$ from the construction of both matrices.
Finally,
\begin{align*}
(Q_{p^2} + Q_{p^2+p+1}&F_{p^2+p+1} F^{ -1}_{p^2}) D_{p^2}^{-1} =\begin{pmatrix}
Q^\be_{p+1} & Id_e \\
F_p^{-1} & 0
\end{pmatrix} T\\
&=\begin{pmatrix}
Q^\be_{p+1} & Id_e \\
F_p^{-1} & 0
\end{pmatrix} \begin{pmatrix}
D_p^{-1} & 0 \\
0 & D_p^{-1}
\end{pmatrix}\\
&=\begin{pmatrix}
Q^\be_{p} + Id_e & Id_e \\
F_p^{-1} & 0
\end{pmatrix} \begin{pmatrix}
D_p^{-1} & 0 \\
0 & D_p^{-1}
\end{pmatrix}\\
&=\begin{pmatrix}
Id_e&Q^\be_{p} + Id_e \\
0& F_p^{-1}
\end{pmatrix}\begin{pmatrix}
0 &1 \\
1&0
\end{pmatrix} \begin{pmatrix}
D_p^{-1} & 0 \\
0 & D_p^{-1}
\end{pmatrix}\\
&=D_{p^2}^{-1} \begin{pmatrix}
Id_{\small e} & Id_{\small e} \\
Id_{\small e} &0
\end{pmatrix} \,,
\end{align*}
which concludes the proof of Eq.~(\ref{eq:Qp2plus}).
\subsection{$(1,2)-$sequences and corank 1 submatrices of $K$}
\label{app:corank}
Let us consider two Sobol' matrices $M_p$ and $M_q$ of size $m\times m$ forming a $(t,m,2)$-net. We denote $K=M_qM_p^{-1}$. First we remind that the pairs of matrices $(M_p,M_q)$ and $(Id_m, K)$ generate the same point set (up to indices permutation). Let $\mathcal{K}_k^t$ the $(m-t)\times m$ matrix consisting of the first $k$ rows of $Id_m$ and the first $m-k-t$ rows of $K$:
\vspace{.5cm}
\begin{center}
\raisebox{1cm}{$\mathcal{K}_k^t$ =}
\begin{overpic}[width=3cm]{figures_siggraph25/tmat.pdf}
\put(10,50) {$Id_k$}
\put(10,15) {$K'$}
\put(60,15) {$K''$}
\put(60,50) {0}
\put(10,70) {$k$}
\put(105,50) {$k$}
\put(50,70) {$m-k$}
\put(105,15) {$m-k-t$}
\end{overpic}
\end{center}
\begin{lemma}[Niederreiter \shortcite{Niederreiter1992random} (p. 73) and Paulin et al \shortcite{paulin2022generator}]
$M_p$ and $M_q$ is a $(t,m,2)-$net if and only if for all $k\in \{1,\dots,m\}$, $\mathcal{K}_k^t$ has corank $t$.
\end{lemma}
From block-wise rank computation (the corank of a block triangular matrix with one full rank diagonal block is the corank of the other diagonal block \cite{meyer1973generalized}), we have
$$
\text{corank}(\mathcal{K}_k^t) = \text{corank}(K'')\,.
$$
Focusing on $(1,2)$-sequences, matrices $K''$ for all $m$ and all $k$ of size $(m-k-1)\times(m-k)$ are exactly the $T^{j,w}$ matrices involved in the property $\mathcal{P}$ (see Sect.~\ref{sec:ranks}). As a consequence, if each such matrices $T$ has corank 1, we can conclude that $K$ characterizes a $(1,2)-$sequence.
\end{document}