The Tamed MHD Equations

We study a regularised version of the magnetohydrodynamics (MHD) equations, the tamed MHD (TMHD) equations. They are a model for the flow of electrically conducting fluids through porous media. We prove existence and uniqueness of TMHD on the whole space $\mathbb{R}^{3}$, that smooth data give rise to smooth solutions, and show that solutions to TMHD converge to a suitable weak solution of the MHD equations as the taming parameter $N$ tends to infinity. Furthermore, we adapt a regularity result for the Navier-Stokes equations to the MHD case.

1. Introduction 1.1. Magnetohydrodynamics. The magnetohydrodynamics (MHD) equations describe the dynamic motion of electrically conducting fluids. They combine the equations of motion for fluids (Navier-Stokes equations) with the field equations of electromagnetic fields (Maxwell's equations), coupled via Ohm's law. In plasma physics, the equations are a macroscopic model for plasmas in that they deal with averaged quantities and assume the fluid to be a continuum with frequent collisions. Both approximations are not met in hot plasmas. Nonetheless, the MHD equations provide a good description of the low-frequency, long-wavelength dynamics of real plasmas. In this thesis, we consider the incompressible, viscous, resistive equations with homogeneous mass density, and regularised variants of it. In dimensionless formulation, the MHD equations are of the following form: (1.1) Here, v = v(x, t), B = B(x, t) denote the velocity and magnetic fields, p = p(x, t) is the pressure, Re > 0, Rm > 0 are the Reynolds number and the magnetic Reynolds number and S > 0 denotes the Lundquist number (all of which are dimensionless constants). The two last equations concerning the divergence-freeness of the velocity and magnetic field are the incompressibility of the flow and Maxwell's second equation. For simplicity, in the remainder of the paper, we set S = Rm = Re = 1. Mathematical treatment of the deterministic MHD equations reaches back to the works of G. Duvaut and J.-L. Lions [12] and M. Sermange and R. Temam [41]. Since then, a large amount of papers have been devoted to the subject. We only mention several interesting regularity criteria [8,21,22,27] and the more recent work on nonresistive MHD equations (Rm = ∞) by C.L. Fefferman, D.S. McCormick J.C. Robinson and J.L. Rodrigo on local existence via higher-order commutator estimates [15,16].
In this paper, we want to study a regularised version of the MHD equations on the whole space R 3 , which we call the tamed MHD equations (TMHD), following M. Röckner and X.C. Zhang [39]. They arise from (1.1) by adding two extra terms (the taming terms) that act as restoring forces: The taming terms are discussed in more detail in Section 1.2, and we discuss the results of this work in Section 1.3. The extra term ∇π, which we call the magnetic pressure, will be explained in Section 1.2.3. However, before we study the tamed equations, we want to give an overview of regularisation schemes for the Navier-Stokes and the MHD equations to put our model into the broader context of the mathematical literature.

Damped Navier-Stokes Equations (or BFeD Models).
A related model to the tamed Navier-Stokes equations are the so-called (nonlinearly) damped Navier-Stokes equations: with α > 0 and β ≥ 1. The damping term −α|v| β−1 v models the resistence to the motion of the flow resulting from physical effects like porous media flow, drag or friction or other dissipative mechanisms (cf. [7] and Section 1.2.1). It represents a restoring force, which for β = 1 assumes the form of classical, linear damping, whereas β > 1 means a restoring force that grows superlinearly with the velocity (or magnetic field). X.J. Cai and Q.S. Jiu [7] first proved existence and uniqueness of a global strong solution for 7 2 ≤ β ≤ 5. This range was lowered down to β ∈ (3, 5] by Z.J. Zhang, X.L. Wu and M. Lu in [53]. Furthermore, they considered the case β = 3 to be critical [53,Remark 3.1]. Y. Zhou in [55] proved the existence of a global solution for all β ∈ [3,5]. For the case β ∈ [1, 3), he established regularity criteria that ensure smoothness. Uniqueness holds for any β ≥ 1 in the class of weak solutions. Existence, decay rates and qualitative properties of weak solutions were also investigated by S.N. Antontsev and H.B. de Oliveira [1].
The Brinkman-Forchheimer-extended Darcy model (BFeD model, cf. Section 1.2.1) is a related model for flow of fluids through porous media and uses the damping terms for β ∈ {1, 2, 3} The first problems studied were continuous dependence of the solutions on their parameters, e.g. in F. Franchi, B. Straughan [20]. V.K. Kalantarov and S. Zelik [26] and P.A. Markowich, E.S. Titi and S. Trabelsi [34] proved existence and uniqueness of a weak solution for Dirichlet and periodic boundary conditions, respectively. Long-time behaviour and existence of global attractors have been studied by several authors [35,45,46,49]. An anisotropic version of the equations was studied by H. Bessaih, S. Trabelsi and H. Zorgati [2]. The flow of electrically conducting fluids through porous media, modelled by MHD equations with damping, was studied first by Z. Ye in [48]. He considered the system with nonlinear damping in the equations for both the velocity field (with nonlinear damping parameter α) and the magnetic field (with paramter β) and he proved existence and uniqueness of global strong solutions in the full space case for several ranges of parameters, most interestingly for our purposes for α, β ≥ 4. Z.J. Zhang and X. Yang [54] tried to improve this to α, β > 3, but apparently made a mistake in their proof ( [52, Remark after Equation (9), p. 2]). Z.J. Zhang, C.P. Wu, Z.A. Yao [52] then improved the range to α ∈ [3, 27 8 ], β ≥ 4. The present paper, in a way, deals with the "critical" case α = β = 3, see the discussion of the results below. Furthermore, E.S. Titi and S. Trabelsi [44] proved global well-posedness for an MHD model with nonlinear damping only in the velocity field. They thus avoid the magnetic pressure problem outlined in Section 1.2.3, as opposed to the above papers which seem to have overlooked this issue.
1.2. The Tamed Equations. We first motivate the tamed equations from a physical point of view by pointing out situations where similar models arise naturally in applications. The tamed Navier-Stokes equations are in a sense a variant of the Navier-Stokes equations with damping in the critical case β = 3, combined with a cutoff.

Physical Motivation -Flows Through Porous
Media. Since the tamed equations are closely related to the damped equations of Section 1.1.1, which are much more well-studied, we focus on the occurence of these in the physics literature.
A system with possibly nonlinear damping is considered as a model for the flow of a fluid through porous media. To make this clear, consider the following 1-D compressible Euler equations with damping: The interpretation that this equation models the flow through porous media is in line with the result that as t → ∞, the density ρ converges to the solution of the porous medium equation (cf. e.g. F.M. Huang, R.H. Pan [25]). The momentum, on the other hand, is described in the limit by Darcy's law : which represents a simple linear relationship between the flow rate and the pressure drop in a porous medium.
Here, k is the permeability of the porous medium and µ is the dynamic viscosity. The velocity v is called Darcy's seepage velocity.
In the interface region between a porous medium and a fluid layer, C.T. Hsu and P. Cheng [24,Equation (31), p. 1591] proposed the following equation 1 : where v is the so-called volume-averaged Darcy seepage velocity and p is the volume-averaged pressure. This equation is motivated by a quadratic correction of P. Forchheimer to Darcy's law, called Forchheimer's law or Darcy-Forchheimer law (cf. for example P.A. Markowich, E.S. Titi, S. Trabelsi [34]): with the Forchheimer coefficient γ > 0, the Forchheimer velocity v F as well as the density ρ F . Furthermore, this correction becomes necessary at higher flow rates through porous media, see below for a more detailed discussion. The question arises whether there are cases where a nonlinear correction of yet higher degree is necessary, i.e., where the flow obeys a cubic Forchheimer's law : Indeed, this seems to be the case. P. Forchheimer [18] himself suggested several corrections to Darcy's law at higher flow velocities, one of them being the cubic law (1.3). M. Firdaouss, J.-L. Guermond and P. Le Quéré [17] revisited several historic data sets, amongst them the ones used by Darcy and Forchheimer (who did not correct for Reynolds numbers) and found that the data are actually better described by a linear and cubic Darcy-Forchheimer law (i.e., where γ = 0), in the regime of low to moderate Reynolds numbers, which, as they note [17, p. 333], includes most practical cases: At higher Reynolds numbers, the correct behaviour seems to be quadratic, i.e., Forchheimer's law, in accordance with numerical simulations, e.g. in the work of M. Fourar et al. [19]. The point at which this behaviour changes seems to be dimension-dependent: it occurs much earlier in the numerical simulations of [19,Figure 7] in the 3D case than in the 2D case. Another instance where a cubic Forchheimer law is observed is the high-rate flow in a radial fracture with corrugated walls, cf. M. Buès The tamed Navier-Stokes equations model the behaviour of the flow through porous media in the regime of relatively low to moderate Reynolds numbers, assuming that the higher-order behaviour is much more significant than the linear Darcy behaviour. For a more physically accurate model, one should also include the linear damping term, but we want to focus on the nonlinear effects here and thus for simplicity have omitted this term. The fact that the onset of nonlinear behaviour occurs at higher flow rates is modelled by the cutoff function g N which is nonzero only for sufficiently high velocity. (1.5) The "taming function" g N (defined below) allowed them to obtain stronger estimates than for the untamed Navier-Stokes equations, and hence regularity results that are out of reach for the Navier-Stokes equations. Furthermore, they could show that bounded solutions to the Navier-Stokes equations, if they exist, coincide with the solutions to the tamed Navier-Stokes equations, as shown in [39]. This is a feature that most regularisations of the Navier-Stokes equations do not share. The deterministic case was further studied by X.C. Zhang on uniform C 2 -domains in [51]. In a series of subsequent papers, various properties of the stochastic version of the equations were studied: existence and uniqueness to the stochastic equation as well as ergodicity in [39], Freidlin-Wentzell type large deviations in [38] as well as the case of existence, uniqueness and small time large deviation principles for the Dirichlet problem in bounded domains [37] (both with T.S. Zhang). More recently, there has been resparked interest in the subject, with contributions by Z. Dong and R.R. Zhang [11] (existence and uniqueness for multiplicative Lévy noise) as well as Z. Brzeźniak and G. Dhariwal [5] (existence, uniqueness and existence of invariant measures in the full space R 3 for a slightly simplified system and by different methods).
The taming function was subsequently simplified by changing the expression of g N as well as replacing the argument of the function g N by the square of the spatial L ∞ norm of the velocity, i.e., g N v 2 L ∞ , see W. Liu and M. Röckner [33, pp. 170 ff]. This leads to simpler assumptions on g N as well as easier proofs, especially when spatial derivatives are concerned (which then act only on the remaining factor v). However, this only seems to work within the framework of locally monotone operators, which cannot be applied here due to the crucial assumption of compact embeddings. Thus we do not use this simplification in this work.

The Magnetic Pressure Problem.
From the form of the MHD equations, it would seem like there should also be a "pressure" term ∇π in the equation for the magnetic field. That this is not the case is due to the structure of the nonlinear term in the equation, as was noted already in the work of M. Sermange and R. Temam [41, p. 644].
To make this precise, note that i.e., the nonlinear terms in the magnetic field equation combine to an expression that is manifestly divergence-free. If there existed a magnetic pressure π such that taking the divergence of this equation, observing that div B = 0, would give where ∇π(t, x) ∈ L 2 loc (R + ; L 2 (R 3 )), which implies ∇π = 0. Thus, a careful balancing in the two nonlinear terms leads to the "magnetic pressure" being zero. Now, if we introduce further nonlinearities into the equation for the magnetic field, we might offset this cancellation and thus we will get an artificial "magnetic pressure" in our tamed equations. We can show that this pressure converges to zero as N → ∞, but for the tamed equations, it is undeniably present. We will informally name this phenomenon the magnetic pressure problem: Introducing extra terms N (v, B) that are not divergence-free into the equation for the magnetic field B in the MHD equations will lead to the appearance of an artificial, possibly unphysical "magnetic pressure" π, i.e., it will be of the form This term, being of gradient type, does not manifest itself in the weak formulation of the problem, which is most often studied. Our system is no exception here, so when talking about the pointwise form of the equation, we have to include the magnetic pressure term π, as above. This fact is easily overlooked when introducing regularising terms into the equation for the magnetic field. To give an example, in the work of Z.J. Zhang, C.P. Wu and Z.A.
Yao [52], the authors introduce a damping term |B| β−1 B into the magnetic field equation, but forgot to include a "magnetic pressure" in the strong form of this equation. Note that in other regularisations of the MHD equations, such as the Leray-α model, this problem is avoided by only introducing terms that preserve the structure of the nonlinearities (1.6).
Ideally, one should thus introduce taming terms for the velocity field only. For mathematical reasons, however, at this point we have to content ourselves with taming terms in both components, for otherwise, in the crucial H 1 -estimate (2.10) we could not cancel all four nonlinearities.
1.2.4. The Magnetic Field: To Regularise or Not to Regularise? There seems to be no clear answer, even for schemes which do not introduce magnetic pressure, to the question of whether in the MHD equations the magnetic field should be regularised as well, or whether one should restrict oneself to only regularising the velocity field. A mathematical criticism formulated in J.S. Linshiz and E.S. Titi [32, p. 3] is that regularising the magnetic part as well might add an unnecessary amount of dissipativity to the system. However, for the mathematical reasons discussed in the previous section, we add a taming term to the magnetic field equation as well.
1.2.5. The Tamed MHD Equations. We investigate the case of the deterministic version of tamed magnetohydrodynamics (TMHD) equations in this paper. They can be understood as a model of an electrically conducting fluid in a porous medium at low to moderate Reynolds numbers (cf. P.A. Markowich, E.S. Titi and S. Trabelsi [34]). Following the approach of M. Röckner and X.C. Zhang, we study the following equations: (1.7) If we write y := (v, B), the equations differ from the "untamed" MHD equations by the taming term which is a direct generalisation of the term in (1.5). The norm is defined in equation (2.2) below. The taming function g N : R + → R + is defined by Here, the constant C taming is defined by C taming := 2 max{Re, Rm} = 2.
For the Navier-Stokes case, M. Röckner and X.C. Zhang in [39] set C = 1 ν ∝ Re, so the fact that C taming ∝ Re is not surprising. The factor 2 arises from the fact that we need to tame more terms here. The dependency on Rm seems natural as well.
The idea of the taming procedure remains very clear: try to counteract the nonlinear terms of which there are four in the case of the MHD equations. To pinpoint the exact place where the power of the taming function unfolds, see the discussion after Lemma 2.2.

1.3.
Results and Structure of This Paper. We follow the ideas of [39]. However, the proof of the regularity of the solution requires an MHD adaptation of a result from E.B. Fabes, B.F. Jones and N.M. Rivière [14], which the author could not find in the literature. See Appendix A for a discussion and a proof of this result.
Our main results can be summarised as follows: . For any N > 0, there exists a unique weak solution y to the TMHD equation in the sense of Definition 2.4, depending continuously on the initial data, such that (i) For all t ≥ 0, (1.10) (ii) The solution satisfies y ∈ C(R + ; H 1 ) ∩ L 2 loc (R + ; H 2 ), ∂ t y ∈ L 2 loc (R + ; H 0 ) and for all t ≥ 0, (1.11) (iii) There exist real-valued functions p(t, x) and π(t, x), satisfying the conditions ∇p ∈ L 2 loc (R + ; L 2 (R 3 ; R 3 )), ∇π ∈ L 2 loc (R + ; L 2 (R 3 ; R 3 )), such that for almost all t ≥ 0, in L 2 (R 3 ; In the case of smooth data, we can prove smoothness of the solutions to the TMHD equations: a unique smooth magnetic field and smooth pressure functions which are defined up to a time-dependent constant. Furthermore, the quadruplet (v N , B N , p N , π N ) solves the tamed MHD equations (1.7). Moreover, we have the following estimates: for any T, N > 0 (1.12) sup (1.14) sup .
To be precise, the constant C T,y0,f depends on y 0 H 1 and´T 0 f H 0 ds and goes to zero as both these quantities tend towards zero. The constant C ′ T,y0,f depends on T , y 0 H 2 and sup t∈[0,T ] f (t) H 0 as well as´T 0 ∂ s f 2 H 0 ds. Finally, we have the following convergence result for vanishing taming terms, i.e., in the limit N → ∞.
the unique solutions to the tamed equations (1.7) with initial value y N 0 given by Theorem 1.2.
Then there is a subsequence (N k ) k∈N such that y N k converges to a y in L 2 ([0, T ]; L 2 loc ) and p N k converges weakly to some p in L 9/8 ([0, T ]; L 9/5 (R 3 )). The magnetic pressure π N k converges to zero, weakly in L 9/8 ([0, T ]; L 9/5 (R 3 )). Furthermore, (y, p) is a weak solution to (1.1) such that the following generalised energy inequality holds: Remark 1.5. Compared to the Navier-Stokes case, we get a different type of term in the inequality, namely the last one on the right-hand side of the inequality. Note also that the "magnetic pressure" disappears as N k → ∞.
We have been able to extend all the results of [39] to the case of tamed MHD equations. This posed several technical obstacles: we had to extend the regularity result of [14] to the MHD case, which the author could not find in the literature. Moreover, we describe the magnetic pressure problem in regularised MHD equations. Furthermore, our work basically provides the critical case α = β = 3 of the model considered in [48,52,54].
The paper is organised as follows: we start in Section 2.1 by introducing the functional framework of the problem. Then we state and prove a number of elementary lemmas regarding estimates as well as (local) convergence results for the operators appearing in the tamed MHD equations. Existence and uniqueness of a weak solution is shown in Section 2.2 via a Faedo-Galerkin approximation procedure. Employing the results of Appendix A, we then show in Section 2.3 that for smooth data the solution to the tamed MHD equations remains smooth. Finally, in Section 2.4 we show that as N → ∞, the solution to the tamed MHD equations converges to a suitable weak solution of the (untamed) MHD equations.
The results of this paper are part of the dissertation [40] of the author.
1.4. Notation. Let G ⊂ R 3 be a domain and denote the divergence operator by div. We use the following notational hierarchy for L p and Sobolev spaces: (1) For the spaces L p (G, R) of real-valued integrable (equivalence classes of) functions we use the notation L p (G) or L p if no confusion can arise. These are the spaces of the components v i , B i of the velocity and magnetic field vector fields. (2) We sometimes use the notation L p (G) := L p (G; R 3 ) to denote 3-D vector-valued integrable quantities, especially the velocity vector field and magnetic vector field v and B. For Sobolev spaces, we use the same notational conventions, so for example Finally, if the domain of the functions is not in R 3 , in particular if it is a real interval (for the time variable), then we use the unchanged L p notation. For brevity, we use the following terminology when discussing the terms on the right-hand side of the tamed MHD equations: the terms involving the Laplace operator are called the linear terms, the terms involving the taming function g N are called taming terms and the other terms are called the nonlinear terms. Furthermore, we refer to the initial data y 0 = (v 0 , B 0 ) and the force f = (f v , f B ) collectively as the data of the problem.

Auxiliary Results.
We define the following spaces: the closure being with respect to the norm u m,p := This norm is equivalent to the Sobolev norm given by where ∇ j u denotes the j-th total weak derivative of u of order j. We define the solenoidal spaces by where the divergence is taken in the sense of Schwartz distributions.
To handle the velocity and the magnetic field of the MHD equations at the same time, we will need to define a norm on the space H m := H m × H m . We will define the scalar products in the usual way (see [41], p. 7): for the vector field y = (v, B) define and similarly, for y ∈ H m × H m , we set and accordingly for the norms. They behave just like an ℓ 2 -type product norm.
In a similar fashion we define Lebesgue norms by We will often employ the following Gagliardo-Nirenberg-Sobolev-type interpolation inequality: Let p, q, r ≥ 1 and 0 ≤ j < m. Assume the following three conditions: Then for any u ∈ W m,p ∩ L q (R 3 ; R 3 ), we have the following estimate: Applying it to each component of the norm for y = (v, B), the same estimate carries over to yield (2.6) ∇ j y L r ≤ C m,j,p,q,r y α m,p y 1−α L q . Define the space of (solenoidal) test functions by . We then have the following lemma.
Let P : L 2 (R 3 ; R 3 ) → H 0 be the Leray-Helmholtz projection. Then P commutes with derivative operators ( [36, Lemma 2.9, p. 52]) and can be restricted to a bounded linear operator Furthermore, consider the tensorised projection Then P : L 2 → H 0 is a bounded linear operator. We define the following operator for the terms on the right-hand side of the TMHD equations, projected on the space of divergence free functions: For y := (v, B) and a test functionỹ := (ṽ,B) ∈ H 1 , consider (using the self-adjointness of the projection P) Let us give names to the linear, nonlinear and taming terms of (2.8): The following lemma provides elementary estimates on the terms defined above.
To prove (i), using also the Sobolev embedding theorem we find which implies the assertion. Equality (2.9) follows from the divergence-freeness. To prove (2.10), we start with the inequality The nonlinear terms can be estimated by Young's inequality: The taming terms are estimated using g N (r) ≥ C(r − N ): Since C taming = 2, we get (2.10) by combining the above three estimates.

Existence and Uniqueness of Weak Solutions.
In this section we will study the well-posedness of the weak formulation of the TMHD equations. We start by stating our notion of weak solution. We proceed to show uniqueness first and then existence of weak solutions via a Faedo-Galerkin approximation scheme.
This definition deals with purely spatial test functions, but it can be extended to the case of test functions that depend also on time, as the next proposition demonstrates.
be a weak solution and let T > 0. Then for allỹ ∈ C 1 ([0, T ]; H 1 ) such that Moreover, the following energy equality holds: Proof. To prove that the right-hand side of (2.12) is well-defined, we proceed as in [  Proof. We fix a T > 0 and set z(t) := y 1 (t) − y 2 (t) =: v B . Then we find We first investigate the linear term and find, using integration by parts and the definition of the norms · H 1 that The nonlinear term I N L (t) consists of four terms of the same structure which are estimated by Young's inequality: Thus we need to estimate the latter term which we do by applying the Sobolev embedding theorem as well as Young's inequality: We collect the four terms and use the previous estimates to find (again using the definition of the Sobolev norms) where, by our definition of weak solutions Concerning the taming term, I T (t), we have, using the mean-value theorem of calculus (and the fact that |g ′ | ≤ C), the inequality g N (r) ≤ Cr, Sobolev's embedding theorem and Young's inequality (for (p, q) = (4, 4/3))

This implies that
Hence, altogether we have the inequality Choosing ε = 1 10 , we find that Our next step is to establish existence of weak solutions, i.e., the existence part of Theorem 1.2.
Proof of existence of a weak solution. We use a Faedo-Galerkin approximation on [0, T ]. Take an orthonormal basis {e k | k ∈ N} ⊂ V of H 1 such that span{e k } is dense in H 3 . Fix an n ∈ N. For z = z 1 , . . . , z n ∈ R n and e = e 1 , . . . , e n ∈ V n set where the ρ n are a family of mollifiers such that Now we consider the ordinary differential equation Noting that z · e ∈ H 3 ⊂ H 2 , the estimate (2.10) on A(y), y H 1 from Lemma 2.2 then yields where the constant C n,N contains the norms of all the e i for i = 1, . . . , n (as all the terms on the right-hand side of (2.10) are quadratic in y). Moreover, since is a polynomial in the components of z in each component, it is a smooth map. Hence, the differential equation (2.14) has a unique solution z n (t) such that z n (t) = z n (0) +ˆt 0 b n (z n (s))ds +ˆt 0 f n (s)ds, t ≥ 0.

Now we set
Then the function y n satisfies the differential equation Using the definition of H 1 , the self-adjointness of (I − ∆) and Young's inequality for the last term as well as (2.10) for the second term (dropping the nonlinear terms, all of which have negative signs), we find that This implies Dropping the second term on the left-hand side and using Gronwall's lemma, we find that Again invoking (2.16), we see that sup n∈N sup t∈[0,T ] y n (t) 2 H 1 ≤ C. Since closed balls of H 1 are weakly compact by the Banach-Alaoglu theorem, we find that for almost all t ∈ [0, T ] we have the H 1 -weak convergence y n (t) ⇀ y(t) as n → ∞. To conclude that this holds true for all t ∈ [0, T ] we note that G (k) , as the uniform limit of continuous functions, is continuous, and that on the other hand by the weak convergence just mentioned, n (t) = y n (t), e k H 1 → y(t), e k H 1 . Hence, t → y(t), e k H 1 is continuous for all k ∈ N and by the density of the {e k } k∈N , we find that t → y(t),ỹ H 1 is continuous for allỹ ∈ H 1 . We can thus conclude that t → y(t) is weakly continuous in H 1 , and that for allỹ ∈ H 1 : We next invoke the Helmholtz-Weyl decomposition L 2 = H 0 ⊕ H 0 ⊥ . Since y n (t) − y(t) ∈ H 1 ⊂ H 0 , this allows us to conclude that To be more precise, we can use the Helmholtz-Weyl decomposition for L 2 (R 3 ; R 3 ) for both the velocity component and the magnetic field component of y n − y, and putting these two together yields (2.17). Note that this equation implies the component-wise convergence This can be seen by takingz of the formz = (0, . . . , 0, z, 0, . . . , 0), where z ∈ L 2 (R 3 ; R).
From this equation as well as (2.16) and Fatou's lemma, we get Next we want to show that y is indeed a solution of the tamed MHD equations (1.5).
As this equation is linear in e k , it holds for linear combinations and since span{e k } forms a dense subset in H 3 , we conclude We are left to prove (i) -(iii). We will start with (iii). In equation (2.11) we setỹ = (ṽ, 0) and we find for almost all t ≥ 0 that and infer from [43, Proposition 1.1.] the existence of a functionp with ∇p ∈ L 2 loc (R + ; L 2 (R 3 )) such that (for almost all t ≥ 0) Now we define the pressure p byp and we observe that since ∇ |B| 2 2 = (B · ∇)B ∈ L 2 loc (R + ; L 2 (R 3 ; R 3 )) due to (1.11), and by the regularity ofp, p also satisfies the right regularity condition ∇p ∈ L 2 loc (R + ; L 2 (R 3 )). In the same way (cf. Section 1.2.3), by testing againstỹ = (0,B) ∈ V to get we can find a π such that ∇π ∈ L 2 loc (R + ; L 2 (R 3 )) (note that this function in general cannot be expected to be equal to zero for our equation, as we argued in Section 1.2.3) such that for almost all t ≥ 0 These two statements imply (iii). Next we want to prove (i). We take the scalar product in H 0 of the system (2.21), (2.22) with y(t) as well as (2.9) to find which yields (as ∂ t y(t), y(t) Integrating this inequality gives The right-hand side is continuous in t and thus together with the weak continuity of t → y(t) ∈ H 1 by [50, Proposition 21.23 (d), p. 258] we get that t → y(t) ∈ H 1 is strongly continuous. We then apply (2.10), (i) and Young's inequality to find Hence we can conclude that which implies (1.11).

Existence, Uniqueness and Regularity of a Strong Solution.
In this section, we will show that for smooth initital data, the TMHD equations admit a smooth solution, i.e., Theorem 1.3. To prove this, we have to prove their regularity, which is done via the regularity result of Appendix A. We use the notation from Appendix A. We denote the space-time L p norms by y L p (ST ) . Let (F t ) t≥0 be the Gaussian heat semigroup on R 3 . We define its action on a function by the space-convolution We can thus rewrite the operatorB as Then by Appendix A, Theorem A.4, the weak solution y constructed in Theorem 1.2 satisfies the integral equation The Riesz projection term vanishes here because the Helmholtz-Leray projection P ensures that the divergence of the taming term is zero, and the forcing term has zero divergence by assumption.
Proof of Theorem 1.3. We first prove the regularity statement. To this end, we show the following for all k ∈ N: Hence we find that g N (|y| 2 )y ∈ L 10/3 (S T ). Now, as y 0 and f are smooth, by Lemma A.2, we find that D α x D j t f N ∈ L 10/3 (S T ), |α| + 2j ≤ 2. An application of Theorem A.5 then yields D α x D j t y ∈ L 10/3 (S T ), |α| + 2j ≤ 2.
k → k + 1. Assume (2.26). We want to apply the Sobolev embedding theorem, which is justified as Therefore, and by another application of Lemma A.2, this yields and hence, by Theorem A.5, We have thus shown that The next step of the proof consists of another induction on the number of derivatives. Namely we want to show that We have shown the base case m = 2 already. So we are left to show the induction step m → m + 1.
There are two cases to consider: (a) There is at least one spatial derivative, i.e., we have Applying linearity and the triangle inequality, we see that the term containing the initial condition y 0 is bounded. For the other term we get the upper bound The first term is bounded by the induction hypothesis, since The second term is bounded by Young's convolution inequality and the fact that D x k F t ∈ L 1 (S T ) just like in the proof of Lemma A.1. (b) There are only derivatives with respect to time, i.e., The term containing the initial condition is again not a problem. In a similar way as before we find where we have used that F t solves the heat equation. Now by using integration by parts, in the last term, we transfer one spatial derivative from the Laplacian to the second factor and since 2j − 1 = m, we conclude the boundedness as before by Young's convolution inequality and the L 1 -boundedness of D x k F t−s .
By the Sobolev embedding theorem we now find that y is smooth. Thus we get for every t ≥ 0 that We take this equation and apply 3 different inner products to both sides of the equations: (a) ·, ∂ t y(t) H 0 , which will lead to an estimate for´T 0 ∂ t y 2 H 0 dt (b) First apply ∂ t , then apply ·, ∂ t y(t) H 0 . This will lead to an estimate for ∂ t y 2 H 0 . (c) ·, y(t) H 1 , which gives an estimate for y(t) 2 H 2 .
(a) We find by using Young's inequality Here, we denote by G N a primitive function of g N . Since g N (r) ≤ 2r, we find 0 ≤ G N (r) :=´r 0 g N (s)ds ≤ 2 r 2 2 = r 2 . Integrating from 0 to T yields -estimating the nonpositive terms by zero -the following: T,N,y0,f .
(2.28) (b) We first differentiate Equation (2.27) with respect to t and then take the inner product with ∂ t y in H 0 . Note that for θ, φ, ψ ∈ {v, B}, we get By the (anti-)symmetry of the nonlinear terms, if φ = ψ, the second term vanishes, which accounts for the (v · ∇)v and (B · ∇)B terms. The other two nonlinear terms cancel each other, so we are left with four variations of the first term of the right-hand side of Equation (2.29), which can be simplified further using the divergence-freeness to yield Taking this into account and applying Young's inequality, we find 2.4. Convergence to the Untamed MHD Equations. In this section we stress the dependence of the solution to the tamed equation on N by writing y N . We will prove that as N → ∞, the solutions to the tamed equations converge to weak solutions of the untamed equations, i.e., Theorem 1.4.
Proof. The proof follows along the same lines as that of Theorem 1.2 in [39]. Let y N 0 ∈ H 1 with y N 0 → y 0 in H 0 and (y N , p N ) be the associated unique strong solution given by Theorem 1.2. Combining (2.24) with (2.25) yields (2.31) sup For q ∈ [2, ∞), r ∈ (2, 6] such that 3 r by an application of the Sobolev embedding (2.6) and (2.31) we find Employing the Arzelà-Ascoli theorem and the Helmholtz-Weyl decomposition in the same way as in the proof of Theorem 1.2, we find a subsequence y N k (again denoted by y N ) and a y = v, In fact, we can even prove that for every bounded open set G ⊂ R 3 , we have To this end, let G ⊂Ḡ ⊂ Q for a cuboid Q, and ρ be a smooth, non-negative cutoff function with ρ ≡ 1 on G, ρ ≡ 0 on R 3 \Q, and by Friedrichs' inequality (2.19) The terms I 1 (N, ε), I 2 (N, ε) vanish for N → ∞ as using (2.33) we get lim N →∞ and an analogous computation yields lim N →∞ I 2 (N, ε) = 0. The other two terms can be bounded by and the arbitrariness of ε > 0 implies the claim.
This can be seen as follows: where we have used (2.32) and Chebychev's inequality.
As in the proof of Theorem 1.2, there exist pressure functions p N , π N ∈ L 2 ([0, T ]; L 2 loc (R 3 )) such that ∇p N , ∇π N ∈ L 2 ([0, T ]; L 2 (R 3 ; R 3 )), and we have for almost all t ≥ 0 that To derive the generalised energy inequality, we take a non-negative φ ∈ C ∞ 0 ((0, T ) × R 3 ) and then take the inner products with 2y N φ in H 0 of this equation. Let us use the abbreviation˜=´T 0´R 3 dxdt. Then we geẗ For the second line (2.38), we proceed along similar lines, this time with respect to the space variable. To avoid confusion, we will write the equation in components: We will focus on the velocity terms, the magnetic field works in exactly the same way.
The last term on the right-hand side equals, after another application of integration by parts, and thus Hence, (2.38) can be rewritten as The third line, (2.39) will be dealt with term by term. By the incompressibility condition The second term, in a similar fashion, becomes (again using the divergence-freeness) The first term of the last line here cancels with the second term of (2.40). Thus we only have to deal with the first term of (2.40): and therefore The last terms that we have to treat are the pressure terms of (2.42). For the first term, we find, again by integration by parts and the incompressibility and similarly for the second term Thus, altgether we find that (2.43) Since 0 ≤ φ ∈ C ∞ c , it acts as a density, and thus, by [10, Theorem 1.2.1], the map y →´T 0´R 3 |∇y| 2 φ dxds is lower semi-continuous in L 2 ([0, T ]; H 0 ). Thus, the limit of the left-hand side of (2.43) as N → ∞ is greater than or equal to lim inf On the other hand, the limit of the right-hand side as N → ∞ consists of four terms, which we treat individually. We denote G := supp φ.
For the first term, by Cauchy-Schwarz-Buniakowski, (2.31) and (2.34), we find The second term can be treated in a similar fashion: For the termˆT (2.34) implies convergence in measure for the finite measure µ := φdx ⊗ ds. One can also prove uniform integrability with respect to this measure. Then by the generalised Lebesgue dominated convergence theorem we get Moving on with the energy inequality, the last term of (2.43) can be treated in the same way as just discussed. We are left with the pressure term. As in [39], pp. 547 f., we take the divergence of (2.35) to find (2.44) Similarly, we take the divergence of (2.36) and obtain 3 We note that for N sufficiently large (g N (r)) 9/8 · r 9/16 ≤ Cg N (r) · r. This is obviously true on the set {r | g N (r) = 0}. If r > 0 is such that g N (r) > 0 (which implies r ≥ 1), we have (g N (r)) 9/8 · r 9/16 ≤ g N (r)(g N (r)) 1/8 · r 9/16 ≤ g N (r)2 1/8 r 1/8 · r 9/16 ≤ g N (r)2 1/8 · r 11/16 ≤ 2 1/8 g N (r) · r, where the factor of 2 appears due to the definition of the taming function. Using this inequality and (2.31), we find For the first three nonlinear terms on the right-hand side of (2.44), we note that by Hölder's inequality (first for the product measure dx ⊗ dt and with p = 16/7, q = 16/9, then for dt with p = 14/6, q = 14/8) and the Sobolev embedding (2.6) we havê Recall that ∆p N = div R N , where R N is defined by (2.44). Then we have Here we used the L p theory for singular integrals, e.g. [ , there is a subsequence (p N k ) k and a function p ∈ L 9/8 ([0, T ]; L 9/5 (R 3 ; R 3 )) such that for k → ∞ (2.49) p N k → p weakly in L 9/8 ([0, T ]; L 9/5 (R 3 ; R 3 )).
Thus, in the same way as above, we can employ the generalised Lebesgue dominated convergence theorem to conclude that for φ ∈ C ∞ 0 ((0, In exactly the same way we find a subsequence (N k ) k∈N such that π N k → π weakly in L 9/8 ([0, T ]; L 9/5 (R 3 ; R 3 )).
The limit π satisfies the equation ∆π = 0, which, combined with the integrability property of π yields π ≡ 0, thus eliminating the "magnetic pressure" from the resulting weak equation as well as the generalised energy inequality. Hence we have shown that the solutions to the tamed MHD equations converge to suitable weak solutions to the MHD equations.
They used this formulation to prove regularity estimates in mixed space-time L p spaces. Their results play an important role in showing smoothness for smooth initial data for the weak solution of the tamed Navier-Stokes equations in [39] and in this section we attempt to derive analogous results for the MHD equations.
The main idea of [14] is threefold: 1. Find a divergence-free solution to the heat equation with the initial data of the Navier-Stokes problem via Fourier analysis. 2. Use this solution as a test function in the weak formulation to derive the integral equation.
3. Prove regularity of the solution to the integral equation (which amounts to estimating the nonlinear term B(u, u).
We follow their steps with the necessary modifications of the MHD case. Fabes, Jones and Riviere consider mixed space-time L p,q -norms on S T := R d × [0, T ] for p, q ≥ 2, defined by The space of functions that have finite L p,q -norm is denoted by L p,q (S T ). As we are only interested in the case p = q, we will occasionally assume this for simplicity in the following. All the results that follow, however, are true also in the more general case.
A.1. A Divergence-Free Solution to the Heat Equation on the Whole Space. The first step consists in constructing a symmetric d × d matrix-valued function (t, x) → (E ij (x, t)) d i,j=1 with the following properties: , g ∈ L p (R d ) and ∇ · g = 0 in the sense of distributions), the following convergence holds:ˆR The function is given by where Γ(x, t) := 1 (4πt) n/2 e −|x| 2 /4t denotes the Weierstraß kernel and R j denotes the Riesz transformation, Then one can show that E ij (x, t) ∈ C ∞ (R d × (0, ∞)), and for 1 ≤ p < ∞ and g ∈ L p (R d ), We now want to define the nonlinear operator B. Recall that E i denotes the i-th row of E ij . We denote by u(y, s), , and define the operatorB(u, w) by Note that even though E ij (·, 1) / ∈ L 1 (R d ), since its Fourier transform is not continuous at the origin, and L 1 functions have uniformly continuous Fourier transform, cf. [47, Satz V.2.2, p. 212], we still have D x k E ij ∈ L 1 (S T ). This implies the following: Lemma A.1. Let u, w ∈ L p (S T ), p ≥ 2. ThenB(u, w) ∈ L p/2 (S T ). 4 Note that E is a d × d-matrix, so taking its gradient we get a tensor of rank 3. By multiplying with the vector u, we obtain a tensor of rank 2, i.e., a matrix.
Proof. As this statement is not entirely obvious and we will need it below, we prove it here for the reader's convenience. We want to estimate To simplify notations, we let We denote by f ⊛ g the convolution in space and time, i.e., Then, using that s ∈ [0, t] if and only if t − s ∈ [0, t], the inner integral of the expression we want to estimate can be written and estimated as Thus by Young's convolution inequality and the Cauchy-Schwarz-Buniakowski inequality where we used Hölder's inequality for p = 1, p ′ = ∞ for the sum j,k and the embedding ℓ ∞ (d) ⊂ ℓ 1 (d). The proof is complete.
Given this operator, we define the operator B(y 1 , y 2 ) by The operator B occurs naturally when we (formally) use the function E i as a test function in the weak formulation. This will be the subject of the next section. Before moving on with the theory, let us give a useful Sobolev version of the classical Schauder estimates (for lack of a better name) for the heat equation. To this end, we introduce another short-hand notation for the convolution appearing in the definition ofB, Equation (A.4): Lemma A.2. Let γ ∈ R, p, q ∈ (1, ∞), and f ∈ L q ((0, T ); W γ,p ) be a function. Then the convolution of f with the heat kernel Γ lies in the space L q ((0, T ); W γ+2,p ). More precisely, These can also be rewritten into estimates with respect to temporal derivatives. By S(R d+1 ) we denote the Schwartz space of rapidly decreasing functions, and by S ′ (R d+1 ) we denote its dual space, the space of tempered distributions. Our space of test functions for the weak formulation of the MHD equations is The weak solution satisfies a scalar equation. We cast this scalar equation, by choosing suitable test functions, into an equivalent vector-valued equation. More precisely, we have the following result, corresponding to Theorem (2.1) in [14] for f = 0 and Theorem (4.4) for f = 0.

Theorem A.4 (Integral Equation)
. Let y 0 ∈ L r (R d ), 1 ≤ r < ∞, p, q ≥ 2, p < ∞. If y ∈ L p,q (S T ) is a weak solution of the MHD equations with initial value y 0 , then y solves the integral equation Proof. This can be proven in exactly the same way as in [14].
A.3. Regularity of Solutions to the Integral Equation.
Theorem A.5 (Regularity). Let y be a solution to the equation y + B(y, y) = f and y ∈ L p,q (S T ), 2 q + d p ≤ 1. Let k be a positive integer such that k + 1 < p, q < ∞. If D α x D j t f ∈ L p/(|α|+2j+1),q/(|α|+2j+1) (S T ) whenever |α| + 2j ≤ k, then also D α x D j t y ∈ L p/(|α|+2j+1),q/(|α|+2j+1) (S T ) f or |α| + 2j ≤ k. Proof. The proof proceeds along the same lines as that of Theorem 3.4 in [14], with the necessary modifications to the MHD case. Since y = f − B(y, y), we only have to show that D α x D j t B(y, y) ∈ L p/(|α|+2j+1),q/(|α|+2j+1) (S T ). For k = 1, our assumptions read f ∈ L p,q (S T ) and D xi f ∈ L p/2,q/2 (S T ) for all i. We split the argument into two parts: that for the v-part of the equation, i.e., for the equation v + B 1 (y, y) = f 1 and that for the B-part of the equation, i.e., B + B 2 (y, y) = f 2 . Since the termsB i (u, v)(x, t) are of the form (with summation convention) and as our assumption on the coefficients p, q implies 1 < p 2 , q 2 < ∞, we can apply Lemma A.2 as well as the L p,q -boundedness of the Riesz transform (cf. J.E. Lewis [31,Theorem 4,p. 226]) to find that and therefore D xi B 1 (y, y) ∈ L p/2,q/2 (S T ). This in turn implies D xi v ∈ L p/2,q/2 (S T ). The same argument yields D xi B ∈ L p/2,q/2 (S T ).
This implies that , which is finite by the induction hypothesis.
Case 2: j > 0. In this case, we place all the spatial derivatives on the functions u 1 , u 2 , so we can write Since by definition and integration by parts we havē B(u, w) =ˆt by applying D t we get two kinds of terms from the product rule: (i) If the derivative hits the indicator function, we (formally) get terms of the form i.e., in the first term both integrals shrink to a point due to the delta functions and we are left with one spatial derivative of u and w. In the second term, we are left with an unproblematic spatial integral. (ii) If the derivative operator hits the function E ij , we use the definition of E ij to find D t E ij = ∆E ij . So each temporal derivative is transformed into two spatial derivatives. We can then use integration by parts again to transfer all but one of these (spatial) derivatives to u and w, so this term becomes proportional to Proceeding inductively, we see that if we apply D t for j > 1 times, we have two cases: (i) The time derivative hits the indicator function at least once. In this case we get a term where s ≤ r and |ν| + |η| + 2r = 2(j − 1) + 1 = 2j − 1. Here we get the factor j − 1 because we "lose" one time derivative to the δ-distribution, but we get one more spatial derivative (with the scaling factor 1) due to the derivative from ∇E. (ii) All the derivatives hit ∇E. In this case, by continuing as in case 2 above, transferring all but one derivative onto u and w, we get a term where |β ′ | + |γ ′ | = 2j − 1.
With regard to (A.8), we replace u by D β x u and w by D γ x w, |β| + |γ| = |α| and apply Lemma A.2 to find D j t D α xB (u, w) .