Stochastic quantization associated with the $\exp(\Phi)_2$-quantum field model driven by space-time white noise on the torus in the full $L^1$-regime

The present paper is a continuation of our previous work on the stochastic quantization of the $\exp(\Phi)_2$-quantum field model on the two-dimensional torus. Making use of key properties of Gaussian multiplicative chaos and refining the method for singular SPDEs introduced in the previous work, we construct a unique time-global solution to the corresponding parabolic stochastic quantization equation in the full"$L^{1}$-regime"$\vert\alpha\vert<\sqrt{8\pi}$ of the charge parameter $\alpha$. We also identify the solution with an infinite-dimensional diffusion process constructed by the Dirichlet form approach.


Background
In the present paper, we study stochastic quantization associated with space-time quantum fields with interactions of exponential type, called the exp(Φ) 2 -quantum field model in the Euclidean quantum field theory, in finite volume. The exp(Φ) 2 -quantum field (or the exp(Φ) 2 -measure) µ (α) is a probability measure on D ′ (Λ), the space of distributions on the two-dimensional torus Λ = T 2 = (R/2πZ) 2 , which is given by where the massive Gaussian free field µ 0 is the Gaussian measure on D ′ (Λ) with zero mean and the covariance operator (1 − △) −1 , △ being the Laplacian in L 2 (Λ) with the periodic boundary conditions, α(∈ R) is called the charge parameter, the Wick exponential exp ⋄ (αφ)(x) is formally introduced by the expression and is the normalizing constant. We remark that the diverging term E µ 0 [φ(x) 2 ] plays a role of the Wick renormalization. Since this quantum field model was first introduced by Høegh-Krohn [Høe71] in the "L 2 -regime" it is also called the Høegh-Krohn model. For a physical background and related early works of this model, see e.g., [AH73,AH74,Sim74] and references therein. Kahane [Kah85] constructed a random measure ν called the Gaussian mulptiplicative chaos, in the "L 1 -regime" |α| < √ 8π.
It implies the existence of the exp(Φ) 2 -measure µ (α) under |α| < √ 8π, which gives a generalization of the early works mentioned above. After that, the relevance of both the Gaussian multiplicative chaos and the exp(Φ) 2 -quantum field model has been received much attention by many people in connection with topics like the Liouville conformal field theory and the stochastic Ricci flow. See e.g., [Kah85, DS11, RV14, Ber17, JS17, DS19, Bis20, BP21] and references therein. We should also mention that Kusuoka [Kus92] independently studied the exp(Φ) 2 -quantum field model under |α| < √ 8π.
By heuristic calculations, we observe that the exp(Φ) 2 -measure µ (α) is an invariant measure of the following two-dimensional parabolic stochastic partial differential equation (SPDE in short) involving exponential nonlinearity: where (Ẇ t ) t≥0 is an R-valued Gaussian space-time white noise, that is, the time derivative of a standard L 2 (Λ)-cylindrical Brownian motion {W t = (W t (x)) x∈Λ } t≥0 . We call (1.1) the exp(Φ) 2stochastic quantization equation associated with µ (α) . Due to the singularity of the nonlinear drift term, the interpretation and construction of a solution to this singular-SPDE have been a challenging problem over the past years. For a concise overview on stochastic quantization equations, we refer to [AMR15, AK20, ADVG21, AKMR20] and references therein. Albeverio and Röckner [AR91] first solved (1.1) (in the case where Λ is replaced by R 2 ) weakly under |α| < √ 4π by using the Dirichlet form theory. Inspired by recent quick developments of singular SPDEs based on Hairer's groundbreaking work on regularity structures [Hai14] and the related work, called paracontrolled calculus, due to Gubinelli, Imkeller and Perkowski [GIP15], Garban [Gar20] constructed a unique strong solution to (1.1) (for the case where (△ − 1) is replaced by △, i.e., massless case) in a more restrictive condition than |α| < √ 4π. In our previous paper [HKK21], we constructed the time-global and pathwise-unique solution to the SPDE (1.1) under |α| < √ 4π by splitting the original equation (1.1) into the Ornstein-Uhlenbeck process ∂ t X t (x) = 1 2 (△ − 1)X t (x) +Ẇ t (x), and the shifted equation This split is based on the idea of Da Prato and Debussche [DPD03], which is now called the Da Prato-Debussche trick. By the uniqueness of the solution, we also obtained the identification with the limit of the solutions to the stochastic quantization equations generated by the approximating measures to the exp(Φ) 2 -measure µ (α) , and with the process obtained by the Dirichlet form approach. Our construction of the solution to the shifted equation (1.2) is different from the standard fixed-point argument applied in [DPD03,Gar20].
To be more precise, we proved convergence of solutions to approximating equations of (1.1) by using compact embedding, and then identified the limit as the solution. We should mention that, after [HKK21], Oh, Robert and Wang [ORW21] independently constructed the time-global unique solution to (1.1) in the same regime in [HKK21]. Later in [ORTW20], together with Tzvetkov, they studied the massless case on two-dimensional compact Riemannian manifolds in the L 2 -regime. Besides, elliptic SPDEs, which also realize the exp(Φ) 2 -quantum field model have been studied in e.g., [ADVG21].
The main purpose of the present paper is to construct the time-global and pathwise-unique solution to the parabolic SPDE (1.1) in the full "L 1 -regime" |α| < √ 8π. Although the present paper builds on our previous work [HKK21], we significantly improve the arguments of [HKK21] in several ways. To apply the Da Prato-Debussche trick, we need to construct the Wick exponential of the Ornstein-Uhlenbeck process {exp ⋄ (αX t )} t≥0 as a driving noise of the shifted equation (1.2). Since the Gaussian free field µ 0 is the stationary measure of the Ornstein-Uhlenbeck process {X t } t≥0 , this problem is reduced to the construction of the Wick exponential exp ⋄ (αφ). In [HKK21, Theorem 2.2], we constructed it under |α| < √ 4π by combining the Wick calculus of the Gaussian free field µ 0 with the standard Fourier expansion on a negative order L 2 -Sobolev space H s (Λ) (s < 0). However, this kind of argument does not work beyond the L 2 -regime. Refining existing results on the convergence of the Gaussian multiplicative chaos ν RV14,Ber17], we construct the Wick exponential exp ⋄ (αφ) on a suitable Besov space under |α| < √ 8π (see Theorem 2.1). This is one of the important contributions of the present paper. On the other hand, in this case, since the Wick exponential exp ⋄ (αφ) does not have L 2integrability with respect to µ 0 unlike the case of |α| < √ 4π, we need to modify our arguments mentioned above into L p -setting for the construction of the time-global and pathwise-unique solution to (1.1). Besides, due to the lack of the L 2 -integrability, we cannot follow the argument as in [AR90, AR91, HKK21, AKMR20] to show the closability of the associated Dirichlet form. To overcome this difficulty, in Corollary 2.4, we prove that the Wick exponential exp ⋄ (αφ) has the L 2 -integrability with respect not to µ 0 , but to µ (α) . This key property plays a significant role not only for the closability of the Dirichlet form, but also for the identication of the diffusion process obtained by the Dirichlet form approach with the solution to the SPDE (1.1).
We should mention here that our model is closely connected with the sine-Gordon model (or cos(Φ) 2 -quantum field model), which has been studied for a long period by many authors. See e.g., [Frö76,FS76,FP77,AH79] for the early works. Since the sine-Gordon model is formally obtained by replacing the nonlinearity e αφ by e √ −1αφ , it has some similarities with the exp(Φ) 2model. Indeed, it can be constructed rigorously in the same way as the exp(Φ) 2 -model in the case |α| < √ 4π. On the other hand, for large values of |α| up to √ 8π, further renormalization by counter-terms is required (see [BGN82,DH93] for details). To make a rigorous meaning to stochastic quantization equations associated with both the Φ 4 3 -model and the sine-Gordon model, we require further renormalization procedures beyond the Wick renormalization, and recent developments of regularity structure and paracontrolled calculus enable us to study such singular SPDEs rigorously. In [HS16,CHS18], Hairer, Shen and Chandra proved local well-posedness of (the massless version of) the sine-Gordon stochastic quantization equation by applying regularity structure. Hence, at first sight, one might guess that regularity structure or paracontrolled calculus is applicable to the exp(Φ) 2 -stochastic quantization equation (1.1) beyond the L 2regime. To apply such general theories, we usually assume that the inputs of the solution map to the shifted equation of a given singular SPDE take values in a Besov space B s ∞,∞ (Λ) with some s < 0. (We mention here that the reconstruction theorem in B s p,q (R d ) was also studied by Hairer and Labbé [HL17], but they considered only the models with B s ∞,∞ -type bounds.) In contrast, the Wick exponential of the Ornstein-Uhlenbeck process {exp ⋄ (αX t )} t≥0 , which plays a role of an input in our case, belongs to B s p,p (Λ) for some p ∈ [1, 2), but does not to B s ∞,∞ (Λ) (see Theorem 3.2). Moreover, since the nonlinear term of the SPDE (1.1) has exponential growth, it is out of results by these general theories. Alternatively, by making use of the nonnegativity of exp ⋄ (αX t ), we may define a product between two rough objects exp(αY t ) and exp ⋄ (αX t ) on the right-hand side of the shifted equation (1.2) (see Theorem 4.3). This is the most crucial point in our argument. We remark that the nonnegativity of the Wick exponential is a remarkable and useful property, and is also applied in proofs of previous results (see e.g [AH74, Gar20, ORW21, ORTW20, HKK21]).
The organization of the rest of the present paper is as follows: In Section 1.2, we present the framework and state the main theorems (Theorems 1.1, 1.5 and 1.7). In Section 1.3, we fix some notations and summarize several basic properties on Besov spaces. In Section 2, we introduce an approximation of the Wick exponentials of the Gaussian free fields and show its almost-sure convergence in an appropriate Besov space (see Theorem 2.1). For later use, we modify the argument of [Ber17] to obtain a stronger estimate than existing results. Moreover, we also prove that the exp(Φ) 2 -measure µ (α) is well-defined and the Wick exponential exp ⋄ (αφ) has the L 2 -integrability with respect to µ (α) (see Corollaries 2.3 and 2.4). In Section 3, we prove the almost-sure convergence of the Wick exponential of the infinite-dimensional Ornstein-Uhlenbeck process (see Theorem 3.2). In Section 4, we prove Theorem 1.1 using the result of Section 3. In Sections 5 and 6, we prove Theorems 1.5 and 1.7, respectively. Since some parts of Sections 4, 5 and 6 go in very similar ways to the arguments of the previous paper [HKK21], we sometimes omit the details in the present paper. Finally, in Appendix, we give several estimates on the approximation of the Green function of (1 − △), which is used in Section 2.

Statement of the main theorems
We begin with introducing some notations and objects. Let Λ = T 2 = (R/2πZ) 2 be the twodimensional torus equipped with the Lebesgue measure dx. Let L p (Λ) (p ∈ [1, ∞]) be the usual real-valued Lebesgue space. In particular, L 2 (Λ) is a Hilbert space equipped with the usual inner product Let C ∞ (Λ) be the space of real-valued smooth functions on Λ equipped with the topology given by the convergence f n → f in C ∞ (Λ): , the L 2 -inner product ·, · is naturally extended to the pairing of C ∞ (Λ) and its dual space D ′ (Λ). For k = (k 1 , k 2 ) ∈ Z 2 and x = (x 1 , x 2 ) ∈ Λ, we write |k| = (k 2 1 + k 2 2 ) 1/2 and k · x = k 1 x 1 + k 2 x 2 . Although we work in the framework of real-valued functions, it is sometimes easier to do computations by using a system of complex-valued functions {e k } k∈Z 2 defined by For s ∈ R, we define the real L 2 -Sobolev space of order s with periodic boundary condition by This space is a Hilbert space equipped with the inner product Note that H 0 (Λ) coincides with L 2 (Λ) and we regard H −s (Λ) as the dual space of H s (Λ) through the standard chain H s (Λ) ⊂ L 2 (Λ) ⊂ H −s (Λ) for s ≥ 0. We now define the massive Gaussian free field measure µ 0 by the centered Gaussian measure on D ′ (Λ) with covariance (1 − △) −1 , that is, determined by the formula where △ is the Laplacian acting on L 2 (Λ) with periodic boundary condition. This formula implies for any ε > 0, and thus the Gaussian free field measure µ 0 has a full support on H −ε (Λ). For a charge parameter α ∈ (− √ 8π, √ 8π), we then define the exp(Φ) 2 -quantum field (or the exp(Φ) 2 -measure) µ (α) on D ′ (Λ) by where Z (α) > 0 is the normalizing constant and exp ⋄ (αφ) is the Wick exponential which will be rigorously constructed in Section 2. We prove in Theorem 2.1 that the function φ → Λ exp ⋄ (αφ)(x) dx is a positive measurable function for all |α| < √ 8π. Hence, we may also regard µ (α) as a probability measure on D ′ (Λ) (see Corollary 2.3).
In the present paper, we consider the stochastic quantization equation (1.1) associated with exp(Φ) 2 -measure µ (α) , that is a parabolic SPDE where W = {W t (x); t ≥ 0, x ∈ Λ} is an L 2 (Λ)-cylindrical Brownian motion defined on a filtered probability space (Ω, F, (F t ) t≥0 , P) and (Ẇ t ) t≥0 is its time derivative in weak sense. This driving noise has a convenient Fourier series representation where {e k } k∈Z 2 is a real-valued complete orthonormal system (CONS) of L 2 (Λ) defined by The exponential term of the SPDE (1.1) is difficult to treat as it is, because the solution Φ t is expected to take values in D ′ (Λ) \ C(Λ). For this reason, we need to give a rigorous meaning of this SPDE by the renormalization. We assume some properties for the multiplier function.
Hypothesis 1. ψ : R 2 → [0, 1] is a function satisfying the following properties: Note that ψ does not need to be continuous except the origin. For a function ψ satisfying Hypothesis 1, we define the Fourier cut-off operator P N on D ′ (Λ) by From Hypothesis 1, we have the following.
• lim N →∞ P N f − f H s = 0 for any s ∈ R and f ∈ H s (Λ).
By introducing approximating equations driven by the regularized white noise (P NẆt ) t≥0 , we obtain the following theorem in the full L 1 -regime of the charge parameter α. See Section 1.3 below for the definition of the Besov space B −ε p,p (Λ).
In this paper we call this Φ the strong solution of the SPDE (1.1), because in view of Theorem 1.1 we have the mapping from the initial value φ and the driving noiseẆ t to the process Φ.
Remark 1.2. The key ingredient of the proof is Theorem 2.1 below is the almost-sure convergence of Gaussian multiplicative chaos (GMC in short). The law of GMC was first constructed by Kahane [Kah85], and Robert and Vargas [RV10] extended it for general convolution approximations of the covariance kernel. Although these results give only convergence in law, some stronger convergence results were also obtained: almost-sure convergence for the circle average and standard Fourier projection [DS11] and the convergence in probability for general convolution approximations [Ber17]. See [RV14,BP21] for the reviews of these theories. Our proof of Theorem 2.1 is a modification of [Ber17]. We remark that Hypothesis 1 is prepared for the main theorems on singular SPDEs (for e.g. Theorem 1.1), and the circle average approximation contained in [Ber17] does not satisfy Hypothesis 1. However, our construction of Wick exponentials of the Gaussian free field in Section 2 includes the case of the approximations by averaging treated in [Ber17], in particular the circle average approximation, because the estimates (2.5) and (2.6) below hold also for the approximations by averaging (see Section A.3). See Section 2 for our construction of Wick exponentials (GMC). Remark 1.3. Since the exp(Φ) 2 -measure µ (α) is absolutely continuous with respect to µ 0 (see Corollary 2.3), "µ 0 -almost every φ" can be replaced by "µ (α) -almost every φ".
Remark 1.4. We can refine the state space of the strong solution obtained in Theorem 1.1. Precisely, the strong solution is in C([0, T ]; H −ε (Λ)) almost surely (see Corollary 1.6 for detail).
To introduce another approach to the SPDE (1.1), we define the regularized exp(Φ) 2 -measure where Z  Hypothesis 2. The operators P N defined by (1.5) satisfy the following properties.
(ii) For any p ∈ (1, 2), s ∈ R, there exists a constant C > 0 such that for any f ∈ B s p,p (Λ). If ψ is a Schwartz function and the inverse Fourier transform of ψ is a nonnegative function, then Hypothesis 2 holds. See e.g., [BCD11, Proposition 2.78].
Theorem 1.5. Assume that ψ satisfies Hypotheses 1 and 2. Let |α| < √ 8π and ε > 0. For any Let ξ N be a random variable with the law µ is a stationary process and converges in law as N → ∞ to the strong solution Φ stat of (1.1) with an initial law µ (α) , on the space C([0, T ]; H −ε (Λ)) for any T > 0. Moreover, the law of the random variable Φ stat t is µ (α) for any t ≥ 0.
Finally, we discuss a connection between the SPDE (1.1) and the Dirichlet form theory. Let s ∈ (0, 1) be an exponent fixed later (see Corollary 2.4) and set H = L 2 (Λ) and E = H −s (Λ). Recall that {e k } k∈Z 2 is a real-valued CONS of H defined by (1.4). We then denote by FC ∞ b the space of all smooth cylinder functions F : E → R having the form with n ∈ N, f ∈ C ∞ b (R n ; R) and l 1 , . . . , l n ∈ Span{e k ; k ∈ Z 2 }. Since supp(µ (α) ) = E, two different functions in FC ∞ b are also different in L p (µ (α) )-sense. Moreover, We then consider a pre-Dirichlet form (E, FC ∞ b ) which is given by where (·, ·) H is the inner product of H. Applying the integration by parts formula for µ (α) , we obtain that (E, FC ∞ b ) is closable on L 2 (µ (α) ) (see Proposition 6.1 below for detail), so we can define D(E) as the completion of FC ∞ b with respect to E 1/2 1 -norm. Thus, by directly applying the general methods in the theory of Dirichlet forms (cf. [MR92,CF12]), we can prove quasi-regularity of (E, D(E)) and the existence of a diffusion process M = (Θ, G, (G t ) t≥0 , (Ψ t ) t≥0 , (Q φ ) φ∈E ) properly associated with (E, D(E)).
The following theorem says that the diffusion process Ψ = (Ψ t ) t≥0 coincides with the strong solution Φ obtained in Theorem 1.1.

Notations and preliminaries
Throughout this paper, we use the notation A B for two functions A = A(λ) and B = B(λ) of a variable λ, if there exists a constant c > 0 independent of λ such that A(λ) ≤ cB(λ) for any λ. We write A ≍ B if A B and B A. We write A µ B if we want to emphasize that the constant c depends on another variable µ. We collect several basic facts on function spaces used through this paper. Below we usually denote L p (Λ), H s (Λ) and B s p,q (Λ) by L p , H s and B s p,q , respectively, for the sake of brevity. Denote by S(R 2 ) for the space of real-valued Schwartz functions on R 2 and denote its dual by S ′ (R 2 ), which is the space of tempered distributions. The Fourier transform F is defined by and so the inverse Fourier transform is given by , the usual generalization of the Fourier transform is considered. Let (χ, ρ) be a dyadic partition of unity, that is, they satisfy the following: • χ, ρ : R 2 → [0, 1] are smooth radial functions, where B(x, r) stands for the open ball in R 2 centered at x and with radius r. We then set ρ −1 := χ and ρ j := ρ(2 −j ·) for j ≥ 0. We define the Littlewood-Paley blocks (or the Littlewood- We then define the inhomogeneous Besov norm · B s p,q and the Besov space We recall mainly from [BCD11] some basic properties of Besov spaces. We remark that the setting in [BCD11] is not on a torus but on the Euclidean spaces. However, it is known that most results in [BCD11] also follow in the case of function spaces on a torus, and are proved by a parallel argument or by extending functions on a torus to those on the Euclidean spaces periodically (see e.g. [GIP15, Appendix A]). In view of this fact we refer associate results in [BCD11] below, though there is a difference between a torus and the Euclidean spaces. The following embeddings are immediate consequences of the definition.
It is important to note that B s 2,2 coincides with the Sobolev space H s for any s ∈ R, and B s ∞,∞ coincides with the Hölder space C s (Λ) for any s ∈ R \ N with the equivalent norms ([BCD11, Page 99]). The second and third properties above implies that H s ⊂ B s−ε p,p for any p ∈ [1, 2] and ε > 0.
The following is an immediate consequence of the interpolations of L p -spaces and of ℓ p -spaces.
The following equivalence of norms plays an important role in Corollary 2.4.
p,q be the set of all nonnegative elements in B s p,q . Thanks to the following theorem, a nonnegative distribution is regarded as a nonnegative Borel measure. This fact plays a crucial role in Section 4.
Consequently, the domain of ξ can be extended to whole C(Λ).
In the present paper, we consider the approximation with general Fourier multiplier operators as in (1.5). Since we need a stronger convergence for our purpose, we give a self-contained proof of the construction in this section. As mentioned in Remark 2.6 below, our arguments work more general approximations than previous results.

GFFs and Wick exponentials
Recall that µ 0 is the centered Gaussian measure on D ′ (Λ) with covariance (1 − △) −1 . On the probability space (Ω, F, P), a D ′ (Λ)-valued random variable X with the law µ 0 is called a (massive) Gaussian free field. Recalling (1.3), we have the covariance formula of the random field X: The aim of this section is to define the formal exponential exp(αX) for any GFF X and any α with |α| < √ 8π. Since X is D ′ (Λ)-valued, we need a renormalization procedure to give a rigorous meaning to it. Recall that ψ satisfies Hypothesis 1, and the Fourier cut-off operator P N on D ′ (Λ) is defined by (1.5): where ψ N := ψ(2 −N ·). Since P N maps H −1−ε to C(Λ) for small ε > 0 as mentioned before (after Hypothesis 1), the approximation X N := P N X is a continuous function, so the exponential exp(αX N ) is well-defined. However, to take a limit as N → ∞, we need an approximation with renormalization The following is the main theorem of this section.
Theorem 2.1. Assume that ψ satisfies Hypothesis 1. Let |α| < √ 8π and choose parameters p, β such that Then the sequence {exp ⋄ N (αX)} N ∈N converges in the space B −β p,p , P-almost surely and in L p (P). Moreover, by regarding exp ⋄ N (αX) as the random nonnegative Borel measure exp ⋄ N (αX)(x) dx on Λ for N ∈ N, one has the weak convergence of {exp ⋄ N (αX)} N ∈N almost surely. The limits obtained by different ψ's coincide with each other almost surely.
Remark 2.2. The conclusion of Theorem 2.1 holds under the estimates (2.5) and (2.6) in Proposition 2.5 below, even without Hypothesis 1. See Remark 2.6 below for details. In most references, approximations with continuous parameter are used for the convergence in probability and in L p (P). It is associated to adopt ψ ε := ψ(ε·) instead of ψ N for the approximation. For almost-sure convergence we need discretization of the approximation parameter and sufficiently high speed of the approximation with respect to the parameter in order to control the P-null sets. This is the reason why we choose approximation with discrete parameter as appeared in the definition of ψ N in Theorem 2.1. Here, we remark that for the convergence in L p (P) (in particular the convergence in probability), we do not need to discretize the approximation parameter. Furthermore we remark that we choose the exponential speed 2 −N for the definition ψ N because of the simplicity of the proof, and N −r with sufficiently large r > 0 instead of 2 −N is also sufficient for the almost-sure convergence. See the proof of Theorem 2.1 in the last part of Section 2.4.
in short, and regard them as the corresponding random nonnegative Borel measures on Λ, according to Theorem 1.11.
Although the proof of (i) is completely the same as [HKK21, Corollary 2.3], we note the fact on the uniform positivity of the normalizing constants which is used in the next corollary. By Jensen's inequality, Here we used the fact that D ′ (Λ) M (α) φ,N (x)µ 0 (dφ) = 1 for any x ∈ Λ, which follows from the definition.
Next we show (ii). Let p and β be as in Theorem 2.1. For any n ∈ N, we have N and the dominated convergence theorem, this implies dµ (α) dµ 0 is bounded and strictly positive µ 0 -almost everywhere.
Even though Theorem 2.1 and Corollary 2.3 imply that the random variable φ → exp ⋄ (αφ) belongs to L p (µ (α) ; B −β p,p ), the state space can be chosen smaller. The following fact plays a crucial role in Sections 5 and 6.
By the interpolation between Besov spaces (Proposi- φ,N (Λ), t ∈ (0, 1] by the bound of the heat kernel in spacial component. By Proposition 1.10 we have φ , e k for any k ∈ Z 2 almost everywhere, by using Fatou's lemma we have Thus we complete the proof. Below, we give a self-contained proof of Theorem 2.1. For the proof we prepare a lot of technical results, and in the end of Section 2, Theorem 2.1 is proved.

Approximation of the Green function
By definition, the random field X N = P N X has the covariance function N (x, x). The function G M,N approximates the Green function G Λ defined by (2.1). In the following proposition, we summarize the properties of the function G M,N used in the proof of Theorem 2.1. We regard G M,N as a periodic function on R 2 × R 2 , rather than a function on Λ × Λ.
Proposition 2.5. Assume that ψ satisfies Hypothesis 1. Then for any x, y ∈ R 2 with |x−y| < 1 and any M, N ∈ N, where the remainder term R M,N (x, y) is uniformly bounded over x, y, M, N . Moreover, there exist constants C > 0 and θ > 0 such that, for any M, N ∈ N, Since the proof of Proposition 2.5 is long and technical, we provide it in Appendix A. We remark that (2.6) can be improved by L p -estimate for all p ∈ [1, ∞) (see Proposition A.5).
Remark 2.6. Theorem 2.1 holds true for any multiplier ψ such that the function G M,N defined from ψ satisfies the estimates (2.5) and (2.6). Indeed, in the proof of Theorem 2.1 after Proposition 2.5, we use only (2.5) and (2.6), but do not use Hypothesis 1. The class of approximations satisfying (2.5) and (2.6) is quite large, and includes the approximations by averaging, treated in [Ber17], in particular the circle average approximation (see Section A.3). Moreover, our proofs would go similarly even if we replace the torus Λ with the Lebesgue measure dx and the Gaussian field X generated by free field measure, by a two-dimensional compact Riemannian manifold M with its volume measure σ and a Gaussian random field X M on M with covariance function G M satisfying (2.5) and (2.6) with replacement of |x − y| by the metric d(x, y) in M, respectively. However, in the case of M and 2) for renormalization, which is a constant in the case of the torus with the Lebesgue measure dx, will depend on x ∈ M generally. We are also able to extend it to compact Riemannian manifold with other dimensions. In the case the range of the charge constant α should be changed according to the dimension. Even though we have such extensions, for simplicity, we discuss our problem only on the torus Λ with the Lebesgue measure dx in the present paper.

Uniform integrability
Using the first property (2.5) of Proposition 2.5, we first prove the uniform bound of in short. At the beginning, we present Kahane's convexity inequality (cf. [Kah85]), which plays a significant role in the proof.
Lemma 2.7 (see [Bis20, Proposition 5.6]). Let D be an open and bounded subset of R 2 . Let ϕ 1 , ϕ 2 be continuous Gaussian random fields on D with mean zero and with covariance functions for any x, y ∈ D, then one has The following estimate is useful to determine the regularity of M (α) N . The estimate is called a multifractal property and is proved also in previous results (see e.g. [BP21, Theorem 3.23], [Gar20, Proposition 3.9] and [RV14, Theorem 2.14]). As mentioned in Remark 2.6, our arguments work in the case of more general approximations than those treated in the previous results.
The following lemmas are useful to show the uniform integrability of Λ M (α) Lemma 2.9. For α ∈ R and p ∈ [1, 2] there exists a constant C > 0 such that, for any N ∈ N and δ ∈ (0, 1/4], Proof. For any δ ∈ (0, 1/4] we can choose {x i ; i = 1, 2, . . . , n δ } such that where c is an absolute constant. Since we have by the elementary inequality (a + b) p/2 ≤ a p/2 + b p/2 for a, b ≥ 0 and the shift invariance of the law of M Hence Proposition 2.8 yields the conclusion.
Lemma 2.10. For any α ∈ R there exists a constant C > 0 such that, for any N ∈ N and δ ∈ (0, 1/4], Proof. By the estimate (2.5), Proof. Choosing finite points {x i } such that Λ = [−π, π) 2 ⊂ i B(x i , 1/2) and using the shift invariance of the law of M for an absolute constant C > 0. Let δ ∈ (0, 1/4] and we decompose In the second inequality, we use p ≤ 2 and the nonnegativity of M (α) N . Applying Lemmas 2.9 and 2.10, we have where the constant C ′ is independent of N and δ. Since α 2 p < 8π, by choosing sufficiently small δ, we complete the proof. Proof. By definition of the Besov norm, By the shift invariance of the law of M . Hence by Proposition 2.8,

Almost-sure convergence
In this subsection, we show the almost-sure weak convergence of M (α) N as N → ∞ in the space of positive Borel measures on Λ, and we finally complete the proof of Theorem 2.1. We apply the following proposition several times, which follows from direct computation.
Proposition 2.13. Let X be an n-dimensional centered Gaussian random vector with a covariance matrix V . Then, for a ∈ R n and a Borel function f on R n , E e a·X f (X) = e a·(V a)/2 E [f (X + V a)] .
The following theorem plays a crucial role to prove Theorem 2.1.
Theorem 2.14. Let |α| < √ 8π. Then, there exist positive constants c and C such that Proof. Our proof is based on the same spirit as [Ber17, Sections 3 and 4]. It is well-known that the limiting measure M (α) must be supported on the points x such that lim N →∞ called α-thick points. An essential point of [Ber17] is to decompose M for some fixed n 0 ∈ N and δ > 0. Then, L 1 contribution of M > N can be eliminated, while M < N has a good control in L 2 depending on the choice of n 0 and δ. However, we need the following modifications to obtain the stronger estimate (2.7).
• Let n 0 = δ 3 N be a variable depending on δ and N .
• Replace the indicator function 1 with some Lipschitz function.
Now we start the proof of this theorem. Denote byB(x, r) the open ball in Λ centered at x and with radius r under the canonical metric of Λ. It is sufficient to show (2.7) for f ∈ C(Λ) with suppf ⊂B(0, 1/2). Indeed, we obtain the assertion for general f ∈ C(Λ), once we apply the finite decomposition f = k f k with f k supported in some ballB(x k , 1/2) and the shift invariance of the law of M  N (x, x) for M, N ∈ N. By the estimate (2.5), for any x, y ∈ R 2 with |x| ∨ |y| < 1/2 and any M, N ∈ N with M ≤ N , we have These yield the following: for any sufficiently small δ > 0, there exists an integer N ′ δ depending on δ such that, for any N ′ δ ≤ M ≤ N and |x| ∨ |y| < 1/2 The parameter δ is to be chosen later, as a sufficiently small number compared with 1 − α 2 /8π and the exponent θ in the estimate (2.6). Furthermore, let χ δ be a function on R such that Then we define for each N, i ∈ N such that N ≤ i (actually we will let i = N or N + 1), Let N δ be an integer such that N δ ≥ N ′ δ /δ 3 . From (2.8) we have that, if N ≥ N δ , then for any integers m, n with δ 3 N ≤ m ≤ n and |x| ∨ |y| < 1/2, (1) The terms M > N,N +1 and M > N,N . For any fixed i ∈ {N, N + 1} and x ∈ Λ, we apply Proposition 2.13 to the (i − [δ 3 N ] + 1)-dimensional random vector X = (X n (x)) δ 3 N ≤n≤i and a fixed vector a = (0, . . . , 0, α). Then, since V a = (αC n,i ) δ 3 N ≤n≤i , we have where we used the elementary inequality (1 − a n ), a 1 , . . . , a K ∈ [0, 1].
Since X n (x) has a variance C n,n and (2.9) implies C n,n = (1 + o(δ))C n , we have by the tail estimate of the normal distribution, for some positive constants C δ and C ′ δ depending on δ. Therefore, we obtain the exponential decay  . For any fixed x, y ∈ Λ, we apply Proposition 2.13 to the multidimensional Gaussian random variable X = (X n (x)) δ 3 N ≤n≤i , (X m (y)) δ 3 N ≤m≤j and a fixed vector a ∈ R (i−[δ 3 N ]+1)+(j−[δ 3 N ]+1) such that a · X = α(X i (x) + X j (y)). Since the covariance matrix V of X is given by
(2-2) The integral over |x − y| ≥ 2 −δ 3 N . We have to consider combinations of I terms. We consider only I N +1,N − I N,N , since the other difference I N +1,N +1 − I N,N +1 is estimated by a similar way. For simplicity, we write

Now we decompose
In the region |x − y| ≥ 2 −δ 3 N , we have no choice but to do However, we can use the estimate (2.6). Indeed, Hence by the estimate (2.6) we have Since θ > 0, this decays exponentially if δ is chosen sufficiently small. Finally we consider J 2 . The estimate (2.5) implies that for 2 −δ 3 N ≤ |x − y| ≤ 1,

Hence we have
Since c δ is positive for sufficiently small δ, this completes the proof.  SinceP n X N converges as N → ∞ toP n X uniformly in x ∈ Λ almost surely for each n, we have Letting n → ∞, we have f,M Then, for ω ∈ Ω \ N , it is easy to see thatM (α) ∞ (ω) can be extended to a linear operator on the space linearly spanned by D. Moreover, since D includes the constant function 1, and hence, for f ∈ D In view of these facts, for ω ∈ Ω \ N ,M Proof of Theorem 2.1. Since convergence of the corresponding measures follows from Corollary 2.16, we prove convergence in the Besov space and independence of the limit in ψ.

Wick exponentials of Ornstein-Uhlenbeck processes
For φ ∈ D ′ (Λ) and an L 2 (Λ)-cylindrical Brownian motion W , let X = X(φ) be the unique solution of the initial value problem (3.1) In this section, we consider the Wick exponential of the infinite-dimensional Ornstein-Uhlenbeck (OU in short) process X. First we recall the basic estimate of X in [HKK21].
It is known that the GFF measure µ 0 is the invariant measure of the process X (see e.g., [DPZ96, Theorem 6.2.1]). Therefore, the random variable is also a GFF under the probability measure P ⊗ µ 0 for any t > 0. Thus the existence of the Wick exponential of X is an immediate consequence of Theorem 2.1.
Proof. Using the invariance of µ 0 with respect to X t and using Theorems 2.1 and 2.14, we have the exponential decay for some positive constants c and C, where X is a GFF under the probability P. Then the assertion is obtained by a similar way to the proof of Theorem 2.1.
Denote by X ∞ := lim N →∞ X N the P ⊗ µ 0 -almost-sure limit. The following result is an immediate consequence of the P ⊗ µ 0 -almost-sure convergence in Theorem 3.2.

Global well-posedness of the strong solution
In this section, we consider the approximating equation (1.6): Φ N 0 = P N φ, and prove Theorem 1.1. The proof goes in a similar way to [HKK21, Section 3] with a slight modification. Similarly to the previous paper, we use the Da Prato-Debussche trick, that is, we decompose the solution of (1.6) by Note that X N = P N X(φ), where X(φ) is the solution of (3.1) with the initial value φ. Hence the renormalized exponential of X N in the latter equation (4.2) is equal to Since X N converges to X ∞ in L p ([0, T ]; B −β p,p ) as stated in Corollary 3.3, we consider the solution map of the deterministic equation for any generic nonnegative X ∈ L p ([0, T ]; B −β p,p ).

Products of continuous functions and nonnegative distributions
Since any nonnegative distribution is regarded as a nonnegative Borel measure by Theorem 1.11, the product of a function f ∈ C(Λ) and a nonnegative distribution ξ ∈ D ′ (Λ) is well-defined as a Borel measure.
Definition 4.1. For any f ∈ C(Λ) and any nonnegative ξ ∈ D ′ (Λ), we define the signed Borel measure where µ ξ (dx) is the Borel measure associated with ξ, as in Theorem 1.11.
We recall some properties of the product map M from [HKK21, Section 3.1]. Recall that B s,+ p,q (Λ) denotes the set of nonnegative elements in B s p,q (Λ).
Then the correspondence (Y, X ) → M(f (Y ), X ) is well-defined as a map for any r ′ ∈ [1, r]. Moreover, if r ′ < r, this map is continuous.

Global well-posedness of Υ
In this part, we can consider more general parameters than those in (2.3). We fix such parameters p, β and the time interval [0, T ]. We consider the initial value problem for any δ ∈ (0, 2 p (p − 1) − β), and the mapping is continuous.
Remark 4.6. We remark that, if △ − 1 is replaced by △, then is the right t-uniform estimate ([KOT03, Lemma 2.2]). The constant 1 comes from the bound of e t△/2 ∆ −1 u. In the above proposition, we can omit this constant by using the factor e −t .
We first show the uniqueness of the solution, by following [HKK21, Lemma 3.8]. Since the function x → |x| p is not twice differentiable if p < 2, we need to modify the previous argument. Proof. Let Υ, Υ ′ ∈ Y T be two solutions of (4.4) with the same X and υ. Then Z = Υ − Υ ′ solves the equation where D ∈ L p ([0, T ]; B −β p,p ), because of definition of Y T and Theorem 4.2. Let ε > 0 and define Z ε = e ε△ Z. Then Z ε solves the equation By the regularizing effect of the heat semigroup (Proposition 4.5), e ε△ D belongs to L p ([0, T ]; C ∞ (Λ)). Then by the Schauder estimate (Proposition 4.7), we have that Z ε belongs to C([0, T ]; C ∞ (Λ)).
Hence for any f ∈ C 2 (R), we have where the first equality is justified as a Riemann-Stieltjes integral, because Next we show the existence of the solution, by following [HKK21, Lemma 3.10]. Since the only difference is that we use Besov spaces instead of Sobolev spaces, we omit some details in this part. The following embedding theorem is frequently used below.
. This yields the bound (4.7) for Υ, thus in particular Υ ∈ Y T .

Stationary solution
In this section, we prove Theorem 1.5 and Corollary 1.6 by assuming that ψ satisfies Hypotheses 1 and 2. Recall that Φ N = Φ N (φ) is a unique solution of the SPDE (1.8): and Φ = Φ(φ) is the strong solution obtained by Theorem 1.1. Since the nonlinear term of (1.8) is given by the log-derivative of the approximating measure µ for any t ≥ 0. This means that Φ stat t has a law µ (α) for any t > 0. Corollary 1.6 is obtained as follows. Since for µ (α) -almost every φ ∈ D ′ (Λ). Since µ (α) and µ 0 are absolutely continuous with respect to each other (Corollary 2.3), "µ (α) -almost every φ" can be replaced by "µ 0 -almost every φ". We now turn to proofs of (i) and (ii). The proofs go in very similar ways to [HKK21, Section 4].
Proof of (i). By the definition (3.1) of the OU process X, we can decompose Φ N, , for sufficiently small δ, ε > 0, by the a priori estimate of the OU process (Proposition 3.1) and the uniform bound of Radon-Nikodym derivatives Then by a similar argument to [HKK21, Theorem 4.2], we have the tightness of Proof of (ii). By a similar argument to the proof of [HKK21, Theorem 1.3], we can assume that ξ N converges to ξ in H −ε almost surely. Then we can complete the proof of (ii) by showing that To do this, we decompose Φ N,stat = X(ξ N ) + Y N , as in the proof of (i), and decompose Φ stat = X(ξ) + Y, where Y = S(0, X ∞ (ξ)). Since  ). This is obtained by a similar way to Lemma 4.10. Indeed, the a priori estimate (4.7) holds for Υ N uniformly over N , since {P N } are nonnegative and uniformly bounded as operators on B −β p,p , in view of Hypothesis 2. If {Υ N k } k∈N is a convergent subsequence, then the limit Υ solves the equation (4.5) as a consequence of the continuity of P N as N → ∞, which is assumed by Hypothesis 2.

Relation with Dirichlet form theory
In this section, we prove Theorem 1.7. Although the proof goes in a very similar way to one in [HKK21], we provide a sketch of the proof for readers' convenience.
We fix the parameter s ∈ (0, 1) appearing in Corollary 2.4 and set D = Span{e k ; k ∈ Z 2 }, H = L 2 and E = H −s . In what follows, ·, · stands for the pairing of E and its dual space E * = H s . By Corollary 2.4, the map φ → exp ⋄ (αφ) can be regarded as a B(E)/B(E)-measurable map. Let (E, FC ∞ b ) be the pre-Dirichlet form defined by (1.9), that is, Then we obtain the following: Proposition 6.1. It holds that where LF ∈ L 2 (µ (α) ) is given by Proof. Let ψ = 1 [−1,1] 2 , which satisfies Hypothesis 1. Applying the Gaussian integration by parts formula with respect to µ 0 (see [GJ86,page 207 for all F ∈ FC ∞ b , h ∈ D and N ∈ N. Now we recall Theorem 2.1, Corollary 2.4 and lim N →∞ Z (α) N = Z (α) > 0. Taking the limit N → ∞ on both sides of the above equality, we obtain and this leads us to the desired integration by parts formula (6.1). Besides, applying Corollary 2.4 again, we obtain LF ∈ L 2 (µ (α) ). This completes the proof.

A Green functions and their approximation
In this appendix, we provide some properties of Green functions and their approximation on the whole space and the torus. In the end, we prove a proposition, which yields Proposition 2.5.

A.1 Green function on the whole plane
Recall that ψ is a function satisfying Hypothesis 1, ψ N = ψ(2 −N ·), and G M,N (x, y) = 1 2π We regard G M,N as a periodic function on R 2 × R 2 , rather than a function on Λ × Λ. Then by the Poisson summation formula, we can write it as an infinite sum of decreasing functions Hence we need to observe the behavior of K M,N for our purpose. Setting ρ M,N = 1 2π F −1 (ψ M ψ N ), we can write K M,N as a convolution where △ R 2 is the Laplacian on R 2 , and K is the Green function of 1 − △ R 2 .
Proposition A.1. The function K : R 2 \ {0} → R is positive and has the estimates Proof. By the relation between the heat semigroup and the resolvent kernel, we have dt t for x = 0. Since the integral over (0, |x|/2) and (|x|/2, ∞) are equal in view of the change of variables by s = |x| 2 /4t, we have Hence we observe the behavior of the function Since the integrand is bounded by e −rs/2 on s ≥ 1, we have g(r) e −r/2 for r ≥ 1, so the latter part of (A.1) follows. To consider the estimate on r < 1, we decompose g(r) = Thus we have the former part of (A.1).
Next we consider the convolution of K and a function with sufficient decay. for some C > 0 and γ > 0. Set ρ N = 2 2N ρ(2 N ·) for N ∈ N. Then for any |x| < 1 and N ∈ N, where R is a bounded function with rapid decay as |x| → ∞. Since R * ρ N is bounded, it is sufficient to show that (ρ N * log(| · | ∧ 1)) (x) = log |x| ∨ 2 −N + O(1).

A.2 Green function on the torus
We return to the proof of Proposition 2.5.
Lemma A.4. Let ψ be a function satisfying Hypothesis 1. Then there exists a smooth function ψ with the following properties: •ψ satisfies Hypothesis 1.
• For any k ∈ N 2 there exists a constant C k such that for any x ∈ R 2 , where κ is a constant as in Hypothesis 1(ii) and |k| 1 := |k 1 | + |k 2 | for each k = (k 1 , k 2 ) ∈ N 2 .
Now we prove the following proposition, which yields Proposition 2.5. The estimate (A.8) in the following proposition is better than (2.6), because (A.8) is L p -estimate for all p ∈ [1, ∞).

A.3 Approximations by averaging
We introduce a class of approximations of the Gaussian free field, which contains the circle average (see e.g. [Ber17,BP21,DS11]), and show that the associated kernels also satisfy (2.5) and (2.6) in Proposition 2.5. This implies that our construction of Wick exponentials of the Gaussian free field in Section 2 includes the circle averaging approximation. Let X be the Gaussian free field on Λ = T 2 as defined in Section 2.1, and extend X on R 2 periodically. Let ν be a probability measure on R 2 supported in the unit ball B(0, 1) such that (A.11) sup |x|≤2 R 2 | log(x − y)|ν(dy) < ∞.
For N ∈ N denote by ν N the measure given by ν N (A) = ν(2 N A) for a Borel set A. Define the approximation X N of X by Then the random field X N has the covariance function From these inequalities it follows that Fν is bounded and ζ-Hölder continuous for any ζ ∈ (0, 1]. Hence, the estimate (2.6) is obtained in the same way as the proof of (A.8) in Proposition A.5.