{'fraction_non_alphanumeric': 0.0439608081158415, 'fraction_numerical': 0.014263418017861788, 'mean_word_length': 3.7988210818307904, 'pattern_counts': {'":': 0, '<': 12, '<?xml version=': 0, '>': 15, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 17, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Undecidability of the fate of relaxation in one-dimensional quantum systems

We investigate the relaxation dynamics in an isolated quantum many-body system. The stationary value of an observable after relaxation is a topic of researches in the field of quantum thermalization, since thermalization is a relaxation phenomena where this stationary value coincides with the equilibrium value. Therefore, computing the stationary value in quantum many-body systems is regarded as an important problem. We, however, prove that the stationary value in quantum many-body systems is incomputable. More precisely, we show that whether the stationary value is in the vicinity of a given value or not is an undecidable problem. Our undecidable result is still satisfied when we restrict our system to a one-dimensional shift-invariant system with nearest-neighbor interaction, our initial state to a product state of a state on a single site, and our observable to a shift-sum of a one-body observable. This result clearly shows that there is no general theorem or procedure to decide the presence or absence of thermalization in a given quantum many-body system.

In this paper, we consider this problem by a completely different approach from conventional ones. We employ the ideas from theoretical computer science and quantum information and prove that the above problem of thermalization is undecidable. More precisely, we prove that the expectation value of an observable A after relaxation with a given Hamiltonian H of a one-dimensional system is undecidable. Our result of undecidability is still valid if we restrict the class of the observable A as a shift sum of a given one-body observable, the Hamiltonian H as a shift-invariant nearest-neighbor interaction, and the initial state as a product state of identical states on a single site (except the first site). This result directly implies that there is no general procedure to decide the presence or absence of thermalization, obviously no general theorem.
Undecidability sometimes appears in physics. Examples are dynamical systems [56], repeated quantum measurement [57], and the spectral gap in quantum many-body systems [58], and our finding tells that quantum thermalization stands in this line. The following two tools inspire our proof: the reduction to the halting problem of Turing machine [59,60] and Feynman-Kitaev type quantum emulation of classical machines [61,62]. By adopting several technical ideas, we successfully construct a quantum many-body system where the halting problem of Turing machine determines the destination after relaxation.
This paper is organized as follows: In Sec. 2 and Sec. 3, we provide a pedagogical review of quantum thermalization and theoretical computer science. In Sec. 4, we state three main theorems and a technical lemma. The latter lemma is the most important result from a theoretical aspect. In Sec. 5, we derive two theorems, undecidability of relaxation, from the technical lemma. Most of the remainder of this paper is devoted to proving this lemma. In Sec. 6, we briefly sketch the proof strategy. In Sec. 7, we introduce a classical universal reversible Turing machine, which is responsible for the halting problem. In Sec. 8, we construct the Hamiltonian of the quantum system emulating the classical dynamics of the Turing machine. Since the details of dynamics of our quantum system are a little complicated, in Sec. 9 we introduce an analogous setting, which is easier to treat, and solve the dynamics. In Sec. 10, we go back to the original setting and construct the initial state with which the expectation value after relaxation is indeed undecidable. In Sec. 11, we prove the remaining theorem, undecidability of thermalization, by extending the proof techniques for the previous lemma. In Sec. 12, we briefly comment on what happens if we numerically simulate the constructed system. This paper also serves as the supplementary note for the short letter [63]. Two main theorems and a key lemma in Ref. [63] correspond to Theorem 1b, Theorem 2, and Lemma 1 in this paper, respectively.

Background of quantum thermalization
Before going to our main result, we first summarize the problem of thermalization and theoretical computer science in this and next sections for readers who are not familiar with these topics. If a reader is familiar with them, one can skip this and the next sections.
In the research field of quantum thermalization, we mainly consider whether a quantum many-body system with a Hamiltonian H at the initial state |ψ thermalizes (with respect to an observable A) or not. A state |φ is called thermal with respect to an observable A if its expectation value of A is close to its equilibrium value: φ|A|φ where ρ MC is a microcanonical state with energy φ|H|φ . The symbol means that both-hand sides coincide in the thermodynamic limit (i.e., = 0). We call that an initial state |ψ under the Hamiltonian H thermalizes with respect to A = O(1) 1 if |ψ(t) := e −iHt |ψ is thermal with respect to A for almost all t, that is, is satisfied for any ε > 0. Here χ{·} is the indicator function which takes 1 (resp. 0) if the statement inside the bracket is true (resp. false). In the above case, χ{·} takes one if ψ(t)|A|ψ(t) − Tr[Aρ MC ] < ε is satisfied at the time t, and takes zero otherwise. Note that due to the quantum recurrence theorem [64], for any T there exists τ > T such that the state at time τ , |ψ(τ ) , and the initial state |ψ(0) is arbitrarily close. Our definition of thermalization allows recurrence, while recurrence time should become extremely long in a large system.
LetĀ := lim T →∞ 1 T T 0 dt ψ(t)|A|ψ(t) be the long-time average of A. Then, an initial state thermalizes if the following two conditions are satisfied: • (Relaxation): The time-series fluctuation around the long-time averageĀ converges to zero: The former condition, relaxation, is proven for initial states under some condition (a diverging effective dimension) [36][37][38] and this condition is shown to be fulfilled analytically in physically plausible initial states [39]. Thus, the remaining hard task is to handle the latter condition, convergence to the equilibrium value. This is why various concepts and arguments raised in the field of quantum thermalization concern the long-time average and the equilibrium value.
The main goal of this paper is to prove the incomputability of the value of the long-time average A. In other words, we have no general procedure to detect the behavior of quantum many-body systems after relaxation.

Background of theoretical computer science
We here briefly review some basic notions of theoretical computer science; the Turing machines, the Church-Turing thesis, decision problems, and the halting problem of Turing machines. Detailed explanations can be seen in e.g., the textbook of Moore and Mertens [60].  We employ the description with Turing machines (TM) as computation. The TM is a very simple computation system which consists of a one-dimensional tape (an infinite line of cells) filled with symbols and a finite control with a head and its own internal state (see Supplementary Figure 1). The head reads a single cell on the tape, and can rewrite the symbol of the cell, and can move left or right one cell. At the same time, the internal state of the finite control may change. The transition function specify the rule of dynamics in the following manner as an example: If the internal state of the finite control is q i and the symbol in the cell read by the head is g m , then the machine changes the internal state to q j , rewrites the cell to g n , and moves the head left.
A TM is defined as a set of tape alphabets, a set of internal states, and a transition function 2 . One may feel that TMs are very primitive and too simple. However, there exist universal Turing machines which can implement all possible TMs, and these universal Turing machines can implement almost all computational tasks in our world (e.g., computation by C++ and Python). Thus, it is declared that computational functions are those computable by TMs, which is known as the Church-Turing thesis. By accepting the Church-Turing thesis, our computational ability is equivalent to that of a TM, and therefore we can safely restrict our attention to a universal Turing machine.
We next explain decision problems. In decision problems, a set of infinite inputs with a Yes/No assignment is given, and our task is to answer Yes/No correctly for all inputs. In the case of the primary test, for example, inputs are natural numbers, and we need to answer Yes if and only if the input is a prime number. If there exists a TM which outputs 1 with any input assigned to Yes and outputs 0 with any input assigned to No, this decision problem is called decidable. If a decision problem is not decidable (i.e., no TM answers this problem correctly), this decision problem is called undecidable. We note that the speed of computation does not matter to decidability. In fact, combinatorial optimization problems, which are regarded as hard tasks in practice, are categorized as decidable problems because we can solve them by brute force methods. We also note that undecidable problems must accompany infinite possible inputs because there always exists a TM whose output coincides with the correct answer for finite inputs accidentally. This is one of the reasons why standard setups of theoretical computer science usually accompany infinite inputs.
A famous undecidable problem is the halting problem of TMs. In this decision problem, inputs are input codes for a fixed universal Turing machine 3 , and we need to answer whether this TM halts at some time or moves forever. It is proved that this decision problem is undecidable. In our proof of the undecidability of thermalization, we reduce our problem to the halting problem of TMs.
We sometimes treat decision problems with promise. In the case with promise, the inputs are restricted to a subclass of possible inputs where this promise is satisfied. In other words, we allow not to answer or to answer incorrectly if the input does not satisfy this promise. We again take the primary test as an example. Suppose a promise that the input is a prime or a product of two primes. Then, if the input is a product of three primes, we may incorrectly answer Yes (a TM may output 1) for this input.

Setup and main results
We first state the undecidability of relaxation and then state the undecidability of thermalization in this section. All the essential ideas for the proof of the undecidability of thermalization have already appeared in that of relaxation. Therefore, in this paper we first treat that of relaxation, which spends most of this paper, and then briefly discuss some additional ideas to extend the result of relaxation to thermalization.
Consider a one-dimensional chain of d-level quantum systems with the periodic boundary condition, whose underlying local Hilbert space is denoted by H. Suppose that the Hamiltonian of this chain H is shift-invariant and contains only 1-body terms and nearestneighbor 2-body interactions. We set the length of the system as L. Let A be a non-negative observable of H, and A L be the normalized spatial average of A's: where A i is the operator A acting on the site i. Our interest is a long-time average of A L : where ρ L is the initial state of the system 4 with length 5 L. Whenever without confusion, the dependency of A on H and ρ is dropped. We argue, roughly, that A is incomputable 6 . Our statement still holds even in the case 3 Recall that a single universal Turing machine emulates any TM with any input. 4 More precisely, we first introduce a state ρ with infinite length, and define ρ L as the restriction of the state ρ to the first L consecutive sites. 5 We note that how the system size increases (e.g., 1 ≤ i ≤ L or −L/2 + 1 ≤ i ≤ L/2) does not matter to our result of undecidability because this system has the periodic boundary condition and almost uniform initial states in the form of |φ 0 |φ 1 ⊗L−1 suffices to show the undecidability. 6 To make such an assertion appropriately, some vocabulary from computer science is necessary: We say the operator A is computable if a Turing machine can compute any component (with respect to some standard basis) of A with any specified accuracy (So the inputs of a Turing machine are the indices of the components and the integer specifying the accuracy). Also, the Hamiltonian and the density operator over the one-dimensional chain are computable if their restriction to the finite system size (L) is computable. In our case, the Hamiltonian is computable if and only if its 1-and 2-body terms are computable.
Whenever we say an operator A is an input to the problem, we mean that the bit string describing A is given to a Turing machine as the input. Here, the map of the bit string to an arbitrary component should be computable up to any given accuracy. When the input is an observable or a state over H ⊗∞ , the bit string description is given again. The computability of the map of the bit string to the observable or the that the Hamiltonian is shift-invariant and nearest-neighbor interaction, and the initial state is an almost uniform product state in the following form: To state our claim in a rigorous manner, we define two decision problems with a promise; STA (state time average) and HTA (Hamiltonian time average).
[STA]: The dimension of the local Hilbert space d = dim H, a Hamiltonian H, and a one-body observable A are fixed. These two, H and A, are parameters of the problem 7 .
Input: a density operator ρ.

Promise: Either
Decision problem: Decide which of the above two holds.
[HTA]: The dimension of the local Hilbert space d = dim H, a one-body observable A, and an initial state ρ are fixed. These two, A and ρ, are parameters of the problem 8 .
Input: a Hamiltonian H.

Promise: Either
Decision problem: Decide which of the above two holds.
We now state our main theorems, which claim that these decision problems are undecidable. Below, the necessary dimension d 0 is a fixed number, whose rough estimation is presented in Sec. 10.5 as d 0 120.
We first state the undecidability of relaxation in the form of STA: Theorem 1a Fix the local dimension d ≥ d 0 . Then, there exists a family of shift-invariant Hamiltonians H with 1-body terms and 2-body nearest-neighbor interaction terms having the following properties: Let A be an arbitrary observable on a single site and ε > 0 be an arbitrary error margin. We put no requirement on the observable. Then, for any M > 1 there exist an operator A with A − A ≤ ε, a value A * , and proper errors ε 1 and ε 2 with ε 2 ≤ ε and ε 2 /ε 1 ≥ M such that the STA with the aforementioned A , H, c = A * , ε 2 and ε 1 is undecidable. This remains valid even if the initial state ρ (inputs for STA) is restricted to a product state of a pure state of a single site in the form of Eq. (7), where |φ 0 and |φ 1 are orthogonal to each other.
We here clarify the dependency of quantities: The system Hamiltonian H is fixed independent of the observable A and a small parameter ε. In contrast, the target value A * and state is defined considering restriction to the finite-size system. 7 We also require that H and A are computable. 8 We also require that A and ρ are computable.
two parameters for errors ε 1 and ε 2 depend on H, A, and the parameters ε and M . The task of STA is to determine A ∈ [c + ε 1 , c − ε 1 ] or A / ∈ [c + ε 2 , c − ε 2 ] with c = A * for all inputs ρ with the form (7), and thus to show undecidability we pick up a family of bad inputs depending on A, ε 1 , and ε 2 .
We slightly modified the observable from A to A in order to avoid several unwanted cases that the observable A is an identity operator and that the basis determined by the Hamiltonian H and the eigenbasis of A are accidentally aligned in an unwanted direction.
We next state the undecidability of relaxation in the form of HTA. The next theorem on HTA, Theorem 1b, is stronger than Theorem 1a in that two of three main quantities, the initial state, the observable, and the Hamiltonian, are arbitrary.

Theorem 1b
Fix the local dimension d ≥ d 0 . Let A be an arbitrary observable on a single site and ρ be an arbitrary pure state in the form of Eq. (7). We require that there exists a state |φ 2 orthogonal to |φ 0 and |φ 1 such that φ 2 |A|φ 2 = φ 1 |A|φ 1 . Then, for any large M > 1, there exist proper ε 2 , ε 1 with ε 2 /ε 1 ≥ M and A * such that the HTA for A and ρ with c = A * , ε 2 and ε 1 , is undecidable. This remains valid even if the Hamiltonian H (inputs of HTA) is restricted to a shift-invariant Hamiltonian H with 1-body terms and 2-body nearest-neighbor interaction terms.
Here, we put the assumption on |φ 2 in order to exclude the case that A is close to an identity operator.
We again clarify the dependency of quantities: Both the observable A and the initial state ρ are arbitrary (as far as the condition for |φ 2 is satisfied). Two parameters for the error, ε 1 and ε 2 , depend on the choice of A and ρ. In addition, the target value A * also depends on the choice of A and ρ. The task of HTA is to determine A ∈ [c + ε 1 , c − ε 1 ] or A / ∈ [c + ε 2 , c − ε 2 ] with c = A * for all inputs H in the form of shift-invariant and nearest-neighbor interaction, and thus to show undecidability we pick up a family of bad inputs depending on A, ρ, ε 1 , and ε 2 .
These two theorems are, in fact, different rewritings of the following technical lemma. This lemma claims that there exists a shift-invariant one-dimensional Hamiltonian which decodes the input code from the Hamiltonian itself and emulates the dynamics of a universal reversible Turing machine (URTM) properly. Note that in this lemma, u takes all possible binary bit strings: u ∈ ∪ ∞ n=1 {0, 1} n .

Lemma 1
Fix the dimension of the local Hilbert space at d ≥ d 0 . A complete orthonormal system (CONS) of the local Hilbert space {|e i } d−1 i=0 and an observable A over H with e 1 |A|e 1 = 0 and e 2 |A|e 2 > 0 are given arbitrarily. A universal reversible Turing machine (URTM) on a single tape is given arbitrarily. Then, for any η > 0 there exist a shift-invariant Hamiltonian H (which depends on the CONS) and a set of computable unitary operators {V u } over H (which depends on the CONS, A, and η > 0) with the following properties: • H consists of 1-body terms and 2-body nearest-neighbor interaction terms.
• For any u, V u |e 0 = |e 0 is satisfied.
• By setting the initial state as is satisfied if and only if the URTM halts with the input u, and is satisfied if and only if the URTM does not halt with the input u.
Here, we bound both A and VAV † in order to treat STA and HTA simultaneously. By this lemma, we can prove that the halting problem, a well-known undecidable problem, is not easier than STA nor HTA: So the latter is also undecidable.
We have discussed the undecidability of relaxation, where we consider whether the longtime average is close to a given value A * . In the case of thermalization, we consider whether the long-time average is close to the equilibrium value Tr[Aρ MC ]. Modifying slightly the proof techniques for Lemma 1, we can obtain HTA-type undecidability of thermalization. Here, we set d 1 = d 0 + 5.
The conditions of |φ 2 and |φ 3 mean that the initial state |φ 1 and |φ 0 is not at the edge of the spectrum of A. Except for such anomalous cases, almost all setups are expected to satisfy these conditions.
We clarify the dependency of quantities, which is similar to that of Theorem 1b: Both the observable A and the initial state ρ are arbitrary (as far as the condition for |φ 2 and |φ 3 is satisfied). Two parameters for the error, ε 1 and ε 2 , depend on the choice of A and ρ.

The task of HTA is to determine
for all inputs H in the form of shift-invariant and nearest-neighbor interaction, and thus to show undecidability we pick up a family of bad inputs depending on A, ρ, ε 1 , and ε 2 .

Proof of Theorem 1a and Theorem 1b (undecidability of relaxation) from Lemma 1
Below, without loss of generality, we suppose e 1 |A|e 1 = 0, which can be met by subtracting a constant factor. We also set the target value c in STA and HTA to zero.
Proof of Theorem 1a. We fix the CONS {|e i } arbitrarily. For any given ε > 0 and A, we set the operator A as This choice guarantees e 1 |A |e 1 = 0 and | e 2 |A |e 2 | ≥ ε. We set the system Hamiltonian H as what is given in Lemma 1 with η = ε/8M . With these H and η, we claim that the STA for the following state family for u with ε 1 = η and ε 2 = (1/4 − η)ε is undecidable. The ratio of two errors is bounded as where we used η ≤ 1/8. The undecidability follows from Lemma 1, which suggests that holds if and only if the URTM with the input u halts, and holds if and only if the URTM with the input u does not halt. Since the halting problem of the URTM is undecidable, this STA is also undecidable.
Proof of Theorem 1b. We set a CONS {|e i } such that |e 0 , |e 1 , and |e 2 are the given states |φ 0 , |φ 1 , and |φ 2 , respectively, which ensures e 2 |A|e 2 = 0. Without loss of generality, we can suppose e 2 |A|e 2 > 0. (In case of e 2 |A|e 2 < 0, we obtain an analogous conclusion by setting A → −A.) Let H be the Hamiltonian whose existence is guaranteed by Lemma 1 with the above CONS, the above observable A, and η = e 2 |A|e 2 /8M . Using this H, we construct an A-dependent family of Hamiltonians for u given as We now set η = min( e 2 |A|e 2 /8M, 1/8) and ε 1 = η and holds if and only if the URTM with the input u halts, and holds if and only if the URTM with the input u does not halt. Since the halting problem of the URTM is undecidable, the HTA with the initial state |ψ I = |φ 0 ⊗ |φ 1 ⊗L−1 , the observable A, ε 1 and ε 2 is also undecidable.
6 Strategy for the proof of Lemma 1 Most of the remainder of this paper, from Sec. 7 to Sec. 10, is devoted to proving Lemma 1. Before going to the details of the proof, we here sketch the strategy of the proof. We first introduce a classical reversible Turing machine which our quantum system emulates (Sec. 7). As explained in Sec. 7.2, in order to implement it by a one-dimensional quantum system we generalize conventional Turing machines such that the finite control sits in the tape. We construct a quantum system which emulates the above classical reversible Turing machine (Secs. [8][9][10]. The Hamiltonian is constructed in a similar manner to the Feynman-Kitaev Hamiltonian without a clock: One-step time evolution of the classical Turing machine corresponds to an application of a quantum partial isometry as a quantum walk (Sec. 8). The idea of the Feynman-Kitaev Hamiltonian is first demonstrated with a simple toy example (Sec. 8.1), and then we construct the quantum isometry for our classical Turing machine (Sec. 8.2).
Since the dynamics of our system are a little complicated, we first introduce an analogous setting to our original one but easier to analyze, where the initial state is a single computational basis state representing the input code (Sec. 9). In this analogous setting, we demonstrate how the value of A changes depending on the halting/non-halting of the URTM. In this case, the dynamics are fully solvable since the effective Hamiltonian for this initial state can be expressed as a tridiagonal matrix, which can be exactly diagonalized.
We then proceed to the original setting, a shift-invariant initial state (except site 1), where the initial state is a superposition of computational basis states (Sec. 10). The first two subsections (Sec. 10.1 and Sec. 10.2) are devoted to describing the decoding process of the input code for the URTM from an almost uniform initial state. With the help of the law of large numbers, we demonstrate that the input code is successfully decoded with arbitrarily high probability amplitude by choosing proper parameters. In the subsequent two subsections, we describe how the value of A changes in the case of halting and nonhalting. If the URTM halts, it starts flipping spins in A-cells to increase the value of A. The stable state of the corresponding quantum walk is uniformly distributed over all the states traversed by this discrete-time dynamics, so the expectation value of A becomes large if the URTM halts (Sec. 10.3). On the other hand, since the URTM is simulated on relatively small numbers of sites (M-cells), which are sparsely located in the 1D-chain, if the URTM does not halt, an overwhelming part of the system (A-cells) remains unchanged (Sec. 10.4).

Proof of Lemma 1: (1) Classical Turing machine
In this section, we construct a classical machine which is emulated by quantum many-body systems. Our classical machine consists of three TMs, TM1, TM2, and TM3, working on two types of cells, M-cells and A-cells. TM2 simulates a given URTM, and TM1 decodes its input from bit sequences. TM3 is responsible for changing the value of A when TM2 halts.
An M-cell has a three-layered structure: The first layer is a working space of TM1 and TM2, and the second and third layers, which take 0 or 1, contain the information of the input code for TM2. An A-cell has two distinct states; a 1 and a 2 .

Universal reversible Turing machine
We fix a URTM characterized by the set (Q, Γ, q 0 , q f , s 0 , δ). Here, Q is a set of states of the finite control, Γ is a set of tape alphabet symbols, q 0 and q f are the initial and the final state of the finite control respectively, s 0 is the blank symbol, and δ is the transition function. The set Γ contains the symbol , which stands for the leftmost cell. Without loss of generality, we suppose the followings: 1. We essentially use the quadruple form, that is, split each step into the two sub-steps: The first sub-step is rewriting the internal state and the symbol in the cell, and the second sub-step is the move of the tape head. To realize such moves of the machine, we set Q := Q m × Q u with Q m := {m 0 , m 1 }, where q m ∈ Q m is used to represent the two distinct sub-steps. It equals m 0 when rewriting and m 1 when moving the tape  head. At every step, the value of q m is flipped so that rewriting (other than the Q m part) and moving on the tape occurs alternatingly. This is for the sake of clarity and at the same time to realize the movement of the tape head to the left by a two-body interaction.
2. We suppose that TM2 has unique direction property [65]: Q u is split into three disjoint subsets, Q + , Q − , and Q 0 (i.e., if and only if the tape head has moved to the right (resp. left), and in Q 0 if and only if the tape head has not moved. We promise that q 0 ∈ Q m ⊗ Q + and q f ∈ Q m ⊗ Q − . This property confirms the reversibility of this TM.
3. The tape head of the initial configuration is at the leftmost cell, and the Q m part equals m 1 . Hence, at the first step, it moves to the right.

Generalized URTM
To implement classical TMs by one-dimensional quantum many-body systems, we need some modifications to the above URTM. First, in conventional TMs the finite control sits outside the tape, while we will set the finite control in the tape. Second, we add two (non-universal) TMs, TM1 and TM3, which work before and after TM2 runs, respectively. We first generalize the TM such that the finite control sits in the tape. Any RTM (reversible Turing machine) can be simulated by the following generalized URTM, as long as the length of the tape suffices. The generalized TM is almost the same as the URTM in the previous subsection, except for the following respects: • The finite control sits between the two cells of the tape, and can read only one of the two adjacent cells in a single step. If it reads the cell right (left) to it, it moves to the right (left) or does not move. Each element q of the extended state space Q (we employ the same symbol for this extended space for brevity) tells which cell will interact with the finite control in this step.
• (Unique pairing property) For the sake of reversibility, we assume the following: Each element q of the state space Q tells the previous position of the finite control and the cell it has interacted with. In analogy with the unique direction property of an RTM, we call this property a unique pairing property.
• There are several configurations with no successor.
We next add TM1 and TM3. Since the detailed rule of TM1 will be presented in Sec. 10.1, we here briefly explain the dynamics of TM1 and TM3.
As explained, the second and third layers consist of 01-bit sequences. Let β be a real number whose binary expansion is equal to the input code u, and set the relative frequency of 1 in the second layer to β. We also set the third layer such that almost all cells are filled with 0, and the state 1 rarely appears. TM1 counts the relative frequency of 1 in the second layer until we encounter the first 1 in the third layer.
If and only if TM2 halts, TM3 starts working. TM3 moves right, and if TM3 encounters an A-cell with the state a 1 , it flips the state to a 2 . If TM3 encounters an A-cell with the state a 2 , which means that TM3 has moved around the system and goes back to the start point, then TM3 stops.
We call the composite classical machine of these TMs as generalized URTM and denote it by M G .

The states, symbols, and tape of M G
We define the set of states and symbols of M G as  where we set Here, Q m × Q u represents the states of the finite control of TM1 and TM2, and Q r = {r} is the sign for TM3, which flips the states in A-cells. The alphabets Γ, Γ in and Γ A appear in the first layer of M-cells, the second and third layers of M-cells, and A-cells, respectively. The symbol s 1 and ψ 0 are signs that this cell is an A-cell. Note that Γ contains the blank symbol s 0 . We compile these symbols in Supplementary Figure 3. The rate α will be taken sufficiently small to observe the gap between the halting case and the non-halting case (i.e., almost all cells are set as A-cells, though the number of M-cells is also sufficiently large in L → ∞ limit).
Initially, M G 's state is set to (m 0 , q 0 ), and the Γ A -part of the cell is set to a 1 . Cells in the left of the finite control are set to the blank cell s 0 . Below, for brevity we leave off the symbols · and · and also drop −1, as L is large enough (i.e., we write αL and (1 − α)L instead of α(L − 1) and (1 − α)(L − 1) ).

The move of M G : the first step
When M G is in the initial state (m 0 , q 0 ), this RTM not only mimics the move of TM1 but also prepares the cell with . The finite control reads the site on its right and rewrites the cell as s 0 → and s 1 → (There is no other possibility in legal configurations). We employ two different symbols and only for satisfying the reversibility, and and plays completely the same role: telling the left/right end of the tape. Both TM1 and TM2 stop if they hit these cells, while TM3 just passes this cell transparently. Therefore, in the remainder of this paper, we do not distinguish from . At the first step, the state and the position of the finite control will be changed in the same manner as the TM1's first step. Hence, first, it evolves to (m 1 , q), and then it moves to the right and evolves to (m 0 , q). • If q m = m 0 , it reads the cell in its right. If it is an M-cell, the machine updates its state of Q u part, q u say, and the symbol in the cell according to the transition function δ. At the same time, Q m part of the finite control, q m say, evolves to m 1 . If it is an A-cell, q u does not change and q m simply evolves to m 1 .
• If q m = m 1 and q u ∈ Q + , the finite control is swapped with the cell in its right, and q m evolves to m 0 .
• If q m = m 1 and q u ∈ Q − , the finite control is swapped with the cell in its left, and q m evolves to m 0 .
• If q m = m 1 and q u ∈ Q 0 , the finite control does not move, and q m simply evolves to m 0 .
If q m = m 1 , the previous partner of the finite control is the cell on its right, and its previous position is the same as the present position. If q m = m 0 and q u ∈ Q + (resp. q u ∈ Q − ), the previous partner of the finite control is the cell on the left (resp. right) of it, and the previous position is at the left (resp. right) of the previous partner.
Eventually, the machine may run out of the tape: In case of the periodic boundary condition, the finite control reads the symbol with q u ∈ Q + \{q 0 }. Since these states are illegal configurations, no successor is defined to them.
If the machine runs into the halting state q f , then the finite control evolves into the state r ∈ Q r . The dynamics after halting are given in the next subsection.

The move of M G : simulating TM3
If TM2 halts and the state of the finite control becomes in q ∈ Q r (i.e., q = r), then TM3 starts flipping the state in A-cells from a 1 to a 2 (Supplementary Figure 2.(d)). If the cell on the right of the finite control is an A-cell with its state a 1 , then the head flips the state to a 2 and moves right. If the cell on the right of the finite control is an M-cell, then the finite control just moves right. If the cell on the right of the finite control is an A-cell with its state a 2 , which implies that all of the A-cells have already been flipped to a 2 and the finite control has gone around the periodic system, then TM3 stops.
8 Proof of Lemma 1: (2) Quantum partial isometry corresponding to M G

Toy example of Feynman-Kitaev type Hamiltonian
Before constructing the quantum partial isometry for M G , for readers who are not familiar with the emulation of classical machines by quantum systems, we here demonstrate a toy example of a quantum partial isometry emulating a very simple classical machine. If a reader is familiar with this topic, one can skip this subsection and start reading the subsection entitled Construction of quantum partial isometry for M G .
We consider a simple classical machine, flipping bits with moving from left to right on a single cell line. The finite control also settles in this line. The finite control takes three possible states, q 0 , q 1 and q 2 , and the transition rule of this machine is the following (see also Supplementary Figure 4): 1. If the state of the finite control is q 0 , then the state of the finite control evolves to q 1 .
At the same time, the cell on the left of the finite control is flipped to (if the cell is 0) or˜ (if the cell is 1).
2. If the state of the finite control is q 1 and the cell on the right of the finite control is neither nor˜ , then the machine flips the bit (i.e., 0 → 1 and 1 → 0) in the cell on the right of the finite control. At the same time, the state of the finite control evolves to q 2 .
3. If the state of the finite control is q 2 , then the finite control moves to the right. At the same time, the state of the finite control evolves to q 1 .
4. If the state of the finite control is q 1 and the cell on the right of the finite control is or˜ , then the machine stops.
We emulate the above classical machine by a one-dimensional quantum system with nearest-neighbor interaction. The local Hilbert space is spanned by the following seven states: {|0 , |1 , |q 0 , |q 1 , |q 2 , | , |˜ }.
The local quantum isometry on the sites i and i + 1 is given as with which correspond to the transition rules 1,2,3, respectively. The quantum isometry of the total system is given by The initial state is restricted to the form that only a single site takes |q 0 and other sites are in the local Hilbert subspace spanned by {|0 , |1 }, which ensures the emulation of legal states of the classical machine. By denoting by |ψ n the state representing the state of the classical machine at the n-th step, then is satisfied, which means that the quantum isometry U induces one step evolution of the quantum state.

Construction of quantum partial isometry for M G
Now we consider our original Turing machine. By defining then the sequence of the elements x = (x 1 , x 2 , . . . , x L ) (x i ∈ X) represents an instantaneous configuration of M G , where only one of x i 's is a member of Q ∪ Q r . We here dropped the second and third layers of M-cells (i.e., Γ in from {Γ × Γ A }) for brevity, which are explicitly treated in Sec. 10 A quantum state |x ∈ (H X ) ⊗L represents x, a configuration of M G . Since the states in H X are only the label (except for Γ A ), we assign the CONS given in Lemma 1 to the CONS {|x |x ∈ X} as With this assignment, the value of A indeed varies between i = 1 and i = 2 with the state s 1 : s 1 , a 2 |A|s 1 , a 2 = e 2 |A|e 2 > 0.
Each step of the move (time evolution) of M G can be described by a partial isometry acting on at most two sites simultaneously. Thus, the sum of the partial isometries acting on at most two sites properly emulates the time evolution of M G . Suppose that in M G the tape head is at the i-th cell, or equivalently x i ∈ Q ∪ Q r . We set the partial isometry U i acting on the i − 1, i, and i + 1-th sites as where each summand acts on at most two sites. If the state of the finite control is in (m 0 , q), the move of the machine is emulated 9 by U 0 i . If the i-th and i + 1-th cells in M G evolve from (x i , x i+1 ) to (x i , x i+1 ), the corresponding isometry on the sites i and i + 1 is given by |x i , x i+1 x i , x i+1 |. If the state of the finite control is in (m 1 , q), the move is emulated by one of U 1+ i , U 1− i , or U 10 i , depending on whether q ∈ Q + , Q − , or Q 0 . They swap the sites and change m 1 to m 0 . U r i mimics the move of TM3 if the state is r ∈ Q r and the finite control is at i-th site.
We note that if the i-th site does not correspond to the finite control (i.e., x i / ∈ Q ∪ Q r ), then U i |x = 0, which ensures the fact that only the vicinity of the finite control can change. Owing to this, the isometry on (H X ) ⊗L given by properly emulates the dynamics of the generalized RTM M G .

Proof of Lemma 1: (3) Evaluating A for computational basis state
Since the computation of A from the almost uniform initial state |φ 0 ⊗ |φ 1 ⊗L−1 is a little complicated, in this section we first consider an easier setting with a computational basis initial state. With this initial state, We treat almost uniform initial states in the next section.

General expression of long-time average
We first derive a general expression of the long-time average A. Let H be a Hamiltonian of the system, and E i and |E i be eigenenergy and corresponding energy eigenstate. We expand the initial state |ψ with the energy eigenbasis as |ψ = i c i |E i . Then, the long-time average of an observable A is calculated as where χ(E j = E i ) takes 1 if E j = E i and takes zero otherwise. We set the Planck constant to unity. In the third line, we used the fact that lim T →∞ In particular, if the Hamiltonian H has no degeneracy, we have a simple expression A = i |c i | 2 E i |A|E i . In contract, if there exist degeneracy (E i = E j with i = j), we need to handle off-diagonal elements E 1 |A|E j .

Hamiltonian emulating M G and effective Hamiltonian
We now construct the Hamiltonian of the quantum system emulating M G . We define the Hamiltonian as where U i is introduced in Eq. (37). Let x j be the configuration of M G at the j-th step, and J be the number of time steps until the machine stops. Then, its corresponding quantum state |x j satisfies and our Hamiltonian restricted to the Hilbert space spanned by {x i } J i=1 can be expressed as Note that x J does not have a successor, and hence the evolution of M G stops at this point.
Since H eff is a tridiagonal matrix with the basis {|x i } J i=1 , we can fully solve its eigenenergies and eigenstates. The k-th eigenenergy is written as E k = 2 cos kπ J + 1 (43) (k = 1, 2, . . . , J) with the corresponding energy eigenstate Note that 0 < kπ J+1 < π guarantees the absence of degeneracy in this Hamiltonian.

Computing long-time average of A
The initial state |x 1 is expanded by the energy eigenstates as and thus A is calculated as Here, in the first line we used Eq. (39), and in the last line we used the following relation: Now, we shall show that the last term x j |A L |x j is sufficiently small. Since x j = x j and A L is a sum of one-body observables; A L = 1 L i A i , we find that x j |A L |x j can take a nonzero value only if x j and x j differs only in a single site. By denoting this site by i * , we bound the last term x j |A L |x j as Thus, the last term of Eq. (46) is bounded above as whose right-hand side vanishes in the L → ∞ limit. Therefore, in the following, we dropped this term for brevity. We finally evaluate the first two terms of Eq. (46). Let N A (x j ) be the number of A-cells in x j filled with a 2 . Then, we have where the difference between x j |A L |x j and N A (x j ) L e 2 |A|e 2 comes from the presence of M-cells. In the case of halting, by taking L sufficiently large, we can make M G halt before J/2 steps. In this condition, (51) By taking α sufficiently small, we arrive at the relation A ≥ 1 4 − η e 2 |A|e 2 in case of halting.
In contrast, in case of non-halting, since N A (x j ) = 0 for any j, we have By taking α sufficiently small, we arrive at the relation A ≤ η in case of non-halting.
10 Proof of Lemma 1: (4) Evaluating A for superposition of computational basis states 10.1 Setting of the initial state and decoding of the input We now consider our original setting where the initial state is shift-invariant except the first site, which takes the form of Here, |e 0 represents the state corresponding to the initial state of the finite control, and others represent cells. As announced in Sec. 8.2, we first elongate the local Hilbert subspace to H Γ ⊗ H Γ A ⊗ H in in order to treat the alphabets in Γ in in M G , which represent the second and third layers of M-cells. The Hilbert space H in is a 5-dimensional space which is a sum of a 2 qubit space and a single state space {|ψ 0 }: The input code for TM2 is encoded into |ψ in ∈ {|0 , |1 } ⊗2 ⊂ H in , which sits in all M-cells.
If the cell is an A-cell, this part is blank (constant independent of the input) denoted by |ψ 0 ∈ H in . In the previous section, we set αL cells to M-cells and (1 − α)L cells to A-cells deterministically. In this section, instead of this, we set V |e 1 as a superposition of these two types of cells. We first assign |e 1 and |e 2 as and then apply an operator V which slightly rotates |e k (k = 1, 2) as Here, two states, |s 0 , a 1 |ψ in and |e 1 = |s 1 , a 1 |ψ 0 , correspond to the initial states of M-cells and A-cells, respectively (see also Supplementary Figure 3). The symbol s 0 serves as a blank cell of M-cells in the first layer, and |ψ in ∈ H in stores the input u for TM2 in the second and third layers of M-cells. In the state |s 0 , a k |ψ in , the symbol a k plays no role.
In the state |s 1 , a k |ψ 0 , the symbol s 1 is a sign to be A-cells, and |ψ 0 plays no role. The symbols a 1 and a 2 distinguish two states of A-cells, which change the value of A. We set the operator V as an identity operator on the orthogonal component of In particular, V stabilizes |e 0 (i.e., V |e 0 = |e 0 ). The state |ψ in takes the form of where the binary expansion of β is set to be equal to the input code u in the form of a binary bit string. The amount of β is guessed by the relative frequency of 1's in the second layer (see also Supplementary Figure 2.(a)). The second layer is a superposition of computational basis states, and thus TM1 runs in each computational basis state as a quantum superposition. Consider m copies of √ β |1 + √ 1 − β |0 , which is expanded as Here, N 1 (w) is the number of 1's in the binary sequence w. The probability amplitude for a computational basis state |w is . Due to the law of large numbers, the probability amplitude for states where the relative frequency of 1's is close to β converges to 1 in the large m limit: lim m→∞ w: where the precise meaning of the symbol N1(w) m β is clarified shortly (in Eq. (61)). Hence, if m is sufficiently large compared to the length of the input code, TM1 guesses β correctly from the frequency of 1's.
The length of qubit m is determined by another bit sequence √ γ |1 + √ 1 − γ |0 in the third layer. For any given accuracy 0 < ξ < 1, the information of m is encoded to γ as satisfying As γ set to extremely close to 0, almost all qubits are |0 in this sequence, and the state |1 rarely appears. In particular, |1 appears only after m-th digit with probability larger than 1 − ξ. Owing to this, if |1 appears at the m -th digit for the first time, this is taken as the sign of m ≤ m . Based on the observed value m , the length of the output by TM1 (i.e., the presumed length of the digit of β) is determined as n = 1 4 log 2 m , which ensures With this choice of the output length n , guessing m larger than the true value does not affect the correctness of the estimation of β.

Case when TM2 halts 1: decoding and expression of states
Suppose that TM2 halts on the encoded input. We claim that in this case, with overwhelming probability amplitude the dynamics by U before TM2 halts are described by the configuration of the first L 0 of M-cells, where L 0 is sufficiently large in the following sense: • It is not smaller than the space used by TM1 and TM2.
• It should be large enough to encounter at least a single |1 in the third layer in the first L 0 of M-cells with high probability. In other words, for a given 0 < ξ, ξ < 1 we set γ and L 0 such that where the first line is the same as Eq. (60). By setting ξ, ξ 1, the first |1 in the second layer appears between the first m of M-cells and the first L 0 of M-cells with high probability.
We emphasize that L 0 is independent of the system size L.
Let C be the cluster of the first |C| sites with Here, the o(L 0 ) term is chosen so that C contains at least L 0 of M-cells with high probability. The cells outside of the cluster C may interact with the finite control only when q ∈ Q r , which means that only Γ A -parts may touch. This fact motivates the following representation of the initial state: Let y be the first |C| components of x, that is, the configuration of sites in the cluster C, and denote by Y the set of all possible y's. Then, the initial state of the total system can be expressed as We call y "good" initial configurations if (i) TM1 correctly decode the input code u for TM2, and (ii) in the dynamics induced by the isometry U the finite control does not come out from the cluster C before TM2 reaches the halting state. We write as Y * the set of good initial configurations in the above sense, and define the following unnormalized state: By construction, states outside Y * have negligibly small probability weight, and thus for any δ > 0 we can prepare ψ V,L |ψ * Since δ can be arbitrarily small, in the following we regard |ψ * as the initial state itself. For each good initial configuration y, we define the j-th state as |j, y := U j−1 |y ⊗ (V |e 1 ) ⊗L−|C| .
Recalling that U updates the Γ-part of the M-cell in simulating TM1 and TM2, and Γ A part after the halting, we find that U either changes y into another configuration y , or interchanges V |e k 's (k = 1, 2). The above observation implies that |j, y is written in the form of where k i ∈ {1, 2}. This confirms that two states at different steps with the same initial state are orthogonal to each other: j, y|j , y = 0 (69) for any j = j .
where E k,y := 2 cos kπ J y + 1 , |E k,y := 2 J y + 1 Jy j=1 sin jkπ J y + 1 |j, y are the k-th energy eigenvalue and corresponding energy eigenstate of the effective Hamiltonian H eff = j |j + 1, y j, y| + c.c., and J y is the total number of steps for the termination of M G starting from the configuration in the cluster C as y . Using Eq. (39), the long-time average of A from the initial state |ψ * V,L given in Eq. (65) is calculated as where χ(E k,y = E k ,y ) takes 1 if E k,y = E k ,y and takes zero otherwise. The contribution from the case of y = y (diagonal elements) has already been calculated and shown to be a finite amount in Sec. 9. We shall prove that the contribution to A from the case of y = y (off-diagonal elements) is sufficiently small.
Due to the form of Eq. (71), the condition E k,y = E k ,y implies Let G be the greatest common divisor of J y + 1 and J y + 1, and define k 0 := Jy+1 G and k 0 := J y +1 G . Then, k and k with Eq. (74) are expressed as k = lk 0 and k = lk 0 with l = 1, 2, . . . , G, and thus is satisfied. Hence, for any y = y we have In the sixth line, we used the fact that G l=1 cos (j+j )lπ G is equal to G if j + j is a multiple of 2G and is equal to zero otherwise. In the seventh line, we used the following relation similar to Eq. (48) for any y = y . As is the case of Eq. (48), this relation follows from the fact that the configurations (j, y) and (j , y ) differs at least a single site. Substituting Eq. (76) into Eq. (73) in case of y = y , we arrive at the upper bound for the off-diagonal sum: k,k ,y,y y =y Notably, |Y * | is independent of the total system size L, and hence, by taking the thermodynamic limit L → ∞ the right-hand side becomes arbitrarily small. Thus, using the result in Sec. 9, we have and By taking α sufficiently small, we arrive at the desired result for any η > 0.

Case when TM2 does not halt
In this subsection, we consider the case that TM2 does not halt with the input u. To bound the off-diagonal elements from above, we need a completely different treatment from the case of halting, because the size of the cluster C now becomes the entire system (The cluster C should contain the working space of TM2, which is unlimited in case of non-halting).
To treat the non-halting case, we focus on the fact that most of the cells in the initial state are A-cells. By construction, the number of A-cells is invariant under the time evolution, and if the decoding of input succeeds, all of A-cells are filled with the symbol a 1 at all times.
We first restrict the state space to "good" configurations of x. In this subsection, we employ the word "good" initial configurations with a slightly different definition from the previous subsection. We say that x is a "good" initial configuration if (i) the input for TM2 is correctly decoded by TM1, and (ii) the fraction of A-cells is larger than 1 − α where α is set properly as slightly larger than but close to α. We do not require the size of the working space for TM2. We write as X * the set of good initial configurations in the above sense, and define the following unnormalized state: For a similar reason to the case of Y * in Sec. 10.2, the dynamics starting from |ψ * * L denoted by |ψ * * L (t) is arbitrarily close to the actual dynamics |ψ L (t) . Therefore, in the following, we consider the state |ψ * * L (t) instead of the actual one. Let P := |e 1 e 1 | be the projector onto the A-cell filled with the symbol a 1 , and P i be the aforementioned projector acting on the i-th site. If |x represents a single configuration of M G , Tr[P i |x x|] equals 1 (resp, 0) iff the i-th cell is (resp. is not) an A-cell with the symbol a 1 . Hence, we have for any t. Let ρ i (t) denote the reduced density operator of |ψ * * L (t) ψ * * L (t)| onto the i-th site, and define the averaged density operator over all sites as Then, using a relation we arrive at the desired result: Since we can set α < 2α, we arrive at the desired relation A < 2 √ 2α A . To show Lemma 1, we also need to evaluate the state with V rotation, which is evaluated as The first term has readily obtained through Eq. (86) with replacing A to V AV † , which bounds the first term from above by 2 √ 2α A . The second term is calculated as where in the last line we used Eqs. (56) and (55). Combining them, we find whose right-hand side can become arbitrarily small by taking α sufficiently small. In particular, max{A, V AV † } ≤ η is fulfilled, which completes the whole proof of Lemma 1.

Note: dimension of the local Hilbert space
We have not discussed how large the sufficient dimension of the local Hilbert space is. We here present a rough estimation of it instead of rigorous calculation.
To reduce the dimension, it is useful to enlarge the tape alphabet. More symbols the tape alphabet has, fewer states the finite control needs to have. It has been established that there are URTMs with 10-state and 8-symbol, 15-state and 6-symbol, 24-state and 4symbol, 32-state and 3-symbol, and 138-state and 2-symbol [65]. By employing the 24-state 4-symbol URTM, 24 states in the finite control and (4 + 1) × 2 × 2 2 = 40 symbols in the tape alphabet suffice to simulate TM2 (see also Supplementary Figure 3). Although we do not evaluate the sufficient number of states and symbols for TM1, since the tasks of TM1 are all elementary; counting the number of 1's and 0's. computing logarithm, and division, we roughly estimate that 50 states and 4 symbols suffice for TM1. If the above estimation is correct, the sufficient dimension of the local Hilbert space is (24+50+1)+(40+2+2) = 119.
11 Proof of Theorem 2 (undecidability of thermalization) Our main idea to prove Theorem 2 is tuning the equilibrium value Tr[Aρ MC ] to the desired target value A * by changing the Hamiltonian. To this end, we need some refinements and modifications of the proof of Lemma 1. In Sec.11.1, we first specify the long-time average of A when the URTM halts. In Sec.11.2, we modify the dynamics of TM3, which enables us to push the value of A in case of halting away from that in case of non-halting. In Sec.11.3, we construct a family of Hamiltonians whose microcanonical average of A is tuned to the target value A * .

Value of A in thermodynamic limit when TM2 halts
In the previous section, we only derive an inequality for the value of A when TM2 halts. We here briefly calculate the value of A in the thermodynamic limit when TM2 halts. We start from Eq. (50) and the equality part of Eq. (51) for non-uniform (computational basis) initial states. By taking L to infinity, the number of steps before TM2 halts is negligibly small compared to the total number of steps J. Thus the right-hand side in the equality of Eq. (51) is calculated as lim L→∞ 3 2(J + 1) for computational basis initial states, where we assumed that M-cells and A-cells are uniformly distributed. By combining Eq. (50), this relation leads to the bound for A for computational basis initial states: A − 1−α 2 e 2 |A|e 2 ≤ α A . In the case of uniform initial states (excepts for the first site), we have already shown that off-diagonal terms and various other correction terms vanish in the thermodynamic limit. Thus, by using Eqs. (73) and (78), we readily have for uniform initial states. This means that both A and V AV † converge close to 1−α 2 e 2 |A|e 2 in the thermodynamic limit with errors α A and 2(α + α(1 − α)) A . Since α can be taken arbitrarily small, these relations roughly state that both A and V AV † converge around

Pushing the value of A away from zero
In the setup in the proof of Lemma 1, the value of A when TM2 halts is shown to approach 1 2 e 2 |A|e 2 . This value, more precisely 1 2 ( e 1 |A|e 1 + e 2 |A|e 2 ), comes from the following calculation: In the thermodynamic limit, the time evolution of A(t) in the classical TM is roughly described by a linear function in time, and the dynamics stops when A(t) reaches e 2 |A|e 2 (see the middle graph in Supplementary Figure 5). Thus, the time average of A(t) is 1 2 e 2 |A|e 2 .
In order to inflate A, we modify the rule of TM3 such that the machine spends most of the time with A(t) at e 2 |A|e 2 . The modified TM3, named TM3+, has two additional states in A-cells, and three additional internal states in the finite control. We first recall the states and the rule of TM3 10 : In TM3, A-cells take two different states, a 1 and a 2 , and the finite control has a single internal state r. The rule of the dynamics is as follows: • If the head reads a 1 with the internal state r, then the state of the A-cell is changed to a 2 and the finite control moves right.
• If the head reads a 2 with the internal state r, then the finite control stops. Now, we describe the states and the rule of TM3+. In TM3+, A-cells take four different states, a 1 , a 2 , b l and b r , and the finite control has four different internal states, r, r l , r r , and r c . The rule of the dynamics is as follows: • If the head reads a 1 with the internal state r, then the state of the A-cell is changed to a 2 , and the finite control moves right.
• If the head reads a 2 with the internal state r, then the state of the A-cell is changed to b l , the internal state is changed to r c , and the finite control moves right.
• If the head reads a 2 with the internal state r c , then the state of the A-cell is changed to b r , the internal state is changed to r l , and the finite control does not move (state 1-2 to 1-3 in Supplementary Figure 6).  • If the head reads a 2 with the internal state r l , then the state of the A-cell and the internal state are kept unchanged, and the finite control moves left (state 1-3 to1-4).
• If the head reads b l with the internal state r l , then the state of the A-cell is kept unchanged, the internal state is changed to r r , and the finite control moves left (state 2-1 to 2-2).
• If the head reads a 2 or b l with the internal state r r , then the state of the A-cell and the internal state is kept unchanged, and the finite control moves right (state 2-2 to 2-3 and 2-3 to 3-1).
• If the head reads b r with the internal state r r , then the state of the A-cell is changed to a 2 , the internal state is changed to r c , and the finite control moves right (state 1-1 to 1-2).
• If the head reads b l with the internal state r c , then the finite control stops. This rule yields the following dynamics: After flipping all A-cells to a 2 , then TM3+ shuttles between b l and b r with pushing b r . AT the beginning of this shuttle, the cell b r sits the right of b l . During a single shuttle between b l and b r , b l does not move on average and b r moves right one cell (see Supplementary Figure 6). This shuttle finishes when b r hits b l from the left. The flipping from a 1 to a 2 takes O(L) steps, while the shuttle process after flipping takes O(L 2 ) steps. Namely, at almost all times, all A-cells take the state a 2 . Thus, following the argument in the previous subsection, we find that in the corresponding quantum system the value of A approaches e 2 |A|e 2 when TM2 halts.

Tuning the microcanonical average
We now demonstrate how to tune the long-time average A in case of halting to the microcanonical average Tr[Aρ MC ]. As seen in the previous subsections, by employing TM3+ the long-time average A in case of halting (resp. non-halting) approaches e 2 |A|e 2 (resp. e 1 |A|e 1 = 0) in the thermodynamic limit. We shall show the existence of a proper orthonormal set 11 {|e i } where Tr[Aρ MC ] = e 2 |A|e 2 holds.
We introduce the Hilbert subspace with states orthogonal to |φ 0 , |φ 1 , A |φ 0 and A |φ 1 denoted by S ⊂ H. Let P S be the projection operator onto S. We consider an operator P S AP S on S, and express its eigenstates as {|σ i }. Since P S |σ i = |σ i and P S |σ j = |σ j , we have We set the indices of {|σ i } in decreasing order of expectation value of A; where we denoted the dimension of S by d . Thanks to the assumption φ 2 |A|φ 2 > max |ψ ∈span{|φ0 ,|φ1 } ψ|A|ψ and φ 3 |A|φ 3 < min |ψ ∈span{|φ0 ,|φ1 } ψ|A|ψ , we find 12 Let ρ MC 1 be the reduced density matrix of a microcanonical state to a single site with energy ψ 0 |H|ψ 0 (|ψ 0 := |φ 0 ⊗ |φ 1 ⊗ · · · ⊗ |φ 1 ). Due to the translation invariance of H, Tr[Aρ MC ] = Tr[Aρ MC 1 ] is satisfied. Note that the Hamiltonian H depends on the input code u and several parameters, and thus ρ MC 1 also depends on them. We formally expand it with the set {|e i } as ρ MC 1 = i,j |e i e j | e i |ρ MC 1 |e j =: i,j |e i e j | ρ MC ij . We have already set |e 0 = |φ 0 and |e 1 = |φ 1 .
We first consider the case that σ 2 |A|σ 2 = σ d −1 |A|σ d −1 . We set the orthonormal set 11 This set is not CONS because this set does not span the whole Hilbert space. However, since the dynamics are closed in this subspace spanned by this set, this point does not matter to our argument. (We assign extremely large energy to states outside this subspace in order to keep the contribution from these states to microcanonical average negligible). 12 Since |φ 2 ∈ S, we can expand |φ 2 as |φ 2 = i c i |σ i . Due to the diagonal property of A with the orthonormal set {|σ i }, φ 2 |A|φ 2 = i |c i | 2 σ i |A|σ i < σ 1 |A|σ 1 . A similar argument holds for |φ 3 .
where |e , |e ∈ T := span{|φ 0 , |φ 1 } are proper states diagonalizing P T AP T in this subspace T . By construction, for any p and q, is satisfied. We now change the parameter p from p = 0 to p = 1 continuously with a fixed q. Accordingly the long-time average e 2 |A|e 2 changes from σ d |A|σ d to σ 1 |A|σ 1 continuously. Hence, there exists a proper 0 ≤ p * ≤ 1 which fulfills e 2 |A|e 2 = Tr[Aρ MC 1 ] (see Supplementary Figure 7). Note that by tuning q if needed, we can safely avoid the undesired situation that the above coincidence happens at Tr[Aρ MC 1 ] = e 2 |A|e 2 = 0. In summary, by employing the above p * , these orthonormal set {|e i } realizes the Hamiltonian with which the long-time average A is equal to Tr[Aρ MC 1 ] = 0 if and only if the URTM with the input u halts.
We next consider the case that Figure 8: A possible situation in numerical exact diagonalization. Although the long-time average A L appears to converge to some value, it may suddenly deviate from this value and converge to another value as (a). We cannot distinguish this long metastable case (a) from a true convergence case (b).
σ 1 |A|σ 1 13 . In this case, we set with i ≥ 4, where we again supposed ρ MC 33 = 0 without loss of generality. If q is close to 1, we have e 2 |A|e 2 > Tr[Aρ MC 1 ] by setting p = 1. Thus, by changing p from p = 0 to p = 1, we again find that there exists a proper p * such that e 2 |A|e 2 = Tr[Aρ MC 1 ] is satisfied. We can avoid the unwanted situation Tr[Aρ MC 1 ] = e 2 |A|e 2 = 0 by tuning q. In summary, by employing the above p * , these orthonormal set {|e i } realizes the Hamiltonian with which the long-time average A is equal to Tr[Aρ MC 1 ] = 0 if and only if the URTM with the input u halts. This completes the proof of Theorem 2.
12 Remarks on the infinite limit 12.1 What happens if we numerically simulate this system? Our result shows that the long-time average of A is uncomputable even with the unlimited computational resource. Therefore, one may wonder about the relation between the uncomputability and the fact that we can numerically simulate any system with finite size. We here comment on this point.
Although there are several possible scenarios, we here describe the most intuitive one. We consider the case that TM2 halts after very long steps. If the system size is insufficient, the finite control passes the periodic boundary before halting and stops its movement. With such a small system size, A L takes a value close to zero. On the other hand, if the system size becomes sufficiently large for TM2 to halt, then A L start approaching around 1 2 e 2 |A|e 2 . The problem lies in the fact that we never know what system size is sufficiently large. The undecidability of Turing machines tells that we cannot exclude the possibility that a Turing machine runs for very long steps and suddenly halts. In terms of numerical simulations, any large numerical simulation might be fooled due to insufficient system size, and thus even with the unlimited computational resource we cannot distinguish the case with a long metastable state (Supplementary Figure 8 (a)) and the case that the system has already relaxed ( Supplementary Figure 8 (b)). This shows clear contrast to the numerical simulation of conventional many-body systems, where by taking sufficiently large system size we can make the amount of error from the true value in the thermodynamic limit arbitrarily small.

Difference from the behavior of near-integrable systems
The above explanation may convey the impression that our constructed system is essentially the same as near-integrable systems; H = H int +εV , with an integrable Hamiltonian H int and a small perturbation εV . In near-integrable systems, the small parameter ε characterizes the necessary system size to reach the true long-time average, which can become extremely large by taking ε close to zero. By preparing a sufficiently large system, where its sufficiency is characterized by ε, we succeed in reaching the correct long-time average within an arbitrarily small error.
We claim that our system does not have such a small parameter going to an arbitrarily small number, and as its consequence sufficiently large system size does not exist. One may feel that the length of an input code u, which induces the change in the input Hamiltonian H, serves as a small parameter. If this guess is true, by preparing a large system whose size is determined by the input length, we safely observe the presence or absence of thermalization in this large system. However, unfortunately, this is not the case.
To explain why, we introduce the busy beaver function BB(n), which is the maximum number of steps taken by a halting Turing machine (TM) with n different internal states and an empty input code. We remark that a TM with n internal states and an empty input can be implemented by a fixed URTM with a small number of internal states and an input code with length at most l(n), and vice versa. Thus, the busy beaver function can be regarded as the maximum number of steps the URTM can take with a halting input code with a given length. This fact directly implies that the estimation of the sufficient size for thermalization with a given input code u is essentially the same as the computation of the busy beaver function.
However, the busy beaver function is proven not to be a computable function. More surprisingly, it is also proven that we cannot compute the busy beaver function BB(n) with n ≥ 748 as far as the Zermelo-Fraenkel set theory with the axiom of choice (ZFC), which is roughly equivalent to our mathematics, is consistent [66,67]. This striking result comes from the fact that there exists a 748-state TM which halts if and only if the ZFC is inconsistent. The Gödel's incompleteness theorem reveals that the ZFC cannot prove the consistency of the ZFC itself. Hence, if BB(748) is computable, then by running the above TM for ZFC consistency until the step determined by BB(748) we can prove the consistency of the ZFC, which is a contradiction.
The incomputability of BB(748) leads to the absence of small parameters in our system of thermalization. We collect all Hamiltonians which correspond to the TMs with 748 internal states and empty input codes. Since there are only a finite number of Hamiltonians, any two Hamiltonians in this set have a finite difference between them (i.e., there is no small parameter going to zero). However, as we have seen, we do not have sufficiently large system size in spite of no small parameter in these Hamiltonians, as long as our mathematics is consistent.

Concluding remarks
We have shown that in one-dimensional shift-invariant quantum many-body systems, whether the long-time average A is in the vicinity of c or not is undecidable. The main idea of our proof is a map from a classical URTM to the Hamiltonian of quantum many-body systems. To make the interaction as nearest-neighbor and shift-invariant, we employ the Feynman-Kitaev type construction of a Hamiltonian.
Our method to implement the halting problem to dynamics of quantum many-body systems will apply to various other problems related to thermalization. However, to encode the halting problem, two cases, the halting case and the non-halting case, should be mapped onto two distinct situations (e.g., thermalizing in one case and not thermalizing in the other case). Hence, the presence of both two phenomena should be proved in some concrete systems. This prevents the extension of our method to the eigenstate thermalization hypothesis because no concrete model is proven to satisfy the eigenstate thermalization hypothesis. In contrast, a recent study proves the non-integrability of some concrete models in the sense that the system has no local conserved quantity [54]. Therefore, the decision problem on (non-)integrability for a given Hamiltonian may serve as a stage to extend our method.
We remark that our result does not exclude the possibility that a particular model is proven to thermalize or not to thermalize. Our result excludes the possibility to obtain a general and ultimate criterion to judge the presence or absence of thermalization, but not a model-specific result. In fact, integrable systems are known not to thermalize. We also remark that our model is highly artificial, which is a limitation of our result. We thus have two interesting open problems: One is to find a more natural model exhibiting undecidability. The other is to construct a set of a more restricted class of Hamiltonians whose fate of thermalization is decidable. Both problems merit further research.