Diffusive Limit with Geometric Correction of Unsteady Neutron Transport Equation

We consider the diffusive limit of an unsteady neutron transport equation in a two-dimensional plate with one-speed velocity. We show the solution can be approximated by the sum of interior solution, initial layer, and boundary layer with geometric correction. Also, we construct a counterexample to the classical theory in \cite{Bensoussan.Lions.Papanicolaou1979} which states the behavior of solution near boundary can be described by the Knudsen layer derived from the Milne problem.

The study of neutron transport equation dates back to 1950s. The main methods include the explicit formula and spectral analysis of the transport operators(see [5], [4], [6], [7], [8], [9], [10], [11], [12]). In the classical paper [1], a systematic construction of boundary layer was provided via Milne problem. However, this construction was proved to be problematic for steady equation in [13] and a new boundary layer construction based on ǫ-Milne problem with geometric correction was presented. In this paper, we extend this result to unsteady equation and consider a more complicated case with initial layer involved. Then we can show the diffusive limit of the equation (1.1). Theorem 1.2. Assume g(t, x 0 , w) ∈ C 4 ([0, ∞) × Γ − ) and h( x, w) ∈ C 4 (Ω × S 1 ). Then for the unsteady neutron transport equation (1.1), the unique solution u ǫ (t, x, w) ∈ L ∞ ([0, ∞) × Ω × S 1 ) satisfies (1), (1.8) where the interior solution U ǫ 0 is defined in (3.56), the initial layer U ǫ I,0 is defined in (3.55), and the boundary layer U ǫ B,0 is defined in (3.54). Moreover, if g(t, θ, φ) = t 2 e −t cos φ and h( x, w) = 0, then there exists a C > 0 such that when ǫ is sufficiently small, where the interior solution U 0 is defined in (6.19), the initial layer U I,0 is defined in (6.18), and the boundary layer U ǫ B,0 is defined in (6.17). Remark 1.3. θ and φ are defined in (3.34) and (3.39).
It is easy to see, by a similar argument, the results in Theorem 1.1 and Theorem 1.2 also hold for the one-dimensional unsteady neutron transport equation, where the temporal domain is [0, ∞), spacial domain is [0, L] for fixed L > 0, and velocity domain is [−1/2, 1/2].
1.3. Notation and Structure of This Paper. Throughout this paper, C > 0 denotes a constant that only depends on the parameter Ω, but does not depend on the data. It is referred as universal and can change from one inequality to another. When we write C(z), it means a certain positive constant depending on the quantity z. We write a b to denote a ≤ Cb.
Our paper is organized as follows: in Section 2, we establish the L ∞ well-posedness of the equation (1.1) and prove Theorem 1.1; in Section 3, we present the asymptotic analysis of the equation (1.1); in Section 4, we give the main results of the ǫ-Milne problem with geometric correction; in Section 5, we prove the first part of Theorem 1.2; finally, in Section 6, we prove the second part of Theorem 1.2.

Well-posedness of Unsteady Neutron Transport Equation
In this section, we consider the well-posedness of the unsteady neutron transport equation    ǫ 2 ∂ t u + ǫ w · ∇ x u + u −ū = f (t, x, w) in [0, ∞) × Ω, u(0, x, w) = h( x, w) in Ω u(t, x 0 , w) = g(t, x 0 , w) for w · n < 0 and x 0 ∈ ∂Ω, where the initial and boundary data satisfy the compatibility condition h( x 0 , w) = g(0, x 0 , w) for w · n < 0 and x 0 ∈ ∂Ω. (2.2) We define the L 2 and L ∞ norms in Ω × S 1 as usual: |f ( x, w)| . (2.4) Define the L 2 and L ∞ norms on the boundary as follows: 2.1. Preliminaries. In order to show the L 2 and L ∞ well-posedness of the equation (2.1), we start with some preparations of the penalized neutron transport equation.
Then for the penalized neutron transport equation x 0 , w) for x 0 ∈ ∂Ω and w · n < 0.
for k ≥ 1. Therefore, u k λ converges strongly in L ∞ to a limit solution u λ satisfying Since u 1 λ can be rewritten along the characteristics as based on Lemma 2.1, we can directly estimate where dγ = ( w · n)ds on the boundary.
Proof. We divide the proof into several steps: Step 1: Applying Lemma 2.3 to the solution of the equation (2.16). Then for any φ ∈ L ∞ ([0, ∞) × Ω × S 1 ) satisfying Our goal is to choose a particular test function φ. We first construct an auxiliary function ζ(t). Since u λ (t) ∈ L ∞ (Ω × S 1 ), it naturally impliesū λ (t) ∈ L ∞ (Ω) which further leads toū λ (t) ∈ L 2 (Ω). We define ζ(t, x) on Ω satisfying In the bounded domain Ω, based on the standard elliptic estimate, we have Step 2: Without loss of generality, we only prove the case with s = 0. We plug the test function into the weak formulation (2.30) and estimate each term there. Naturally, we have We estimate the two term on the right-hand side of (2.35) separately. By (2.31) and (2.33), we have In the second equality, above cross terms vanish due to the symmetry of the integral over S 1 . On the other hand, for the second term in (2.35), Hölder's inequality and the elliptic estimate imply Based on (2.32), (2.34), the boundary condition of the penalized neutron transport equation (2.16), the trace theorem, Hölder's inequality and the elliptic estimate, we have Note that we will take Then the only remaining term is Now we have to tackle ∂ t ∇ζ L 2 ([0,t]×Ω×S 1 ) .
Step 3: For test function φ( x, w) which is independent of time t, in time interval [t − δ, t] the weak formulation in (2.30) can be simplified as Taking difference quotient as δ → 0, we know Then (2.44) can be simplified into Then the left-hand side of (2.46) is actually . By a similar argument as in Step 2 and the Poincaré inequality, the right-hand side of (2.46) can be bounded as Therefore, we have For all t, we can further integrate over [0, t] to obtain Step When 0 ≤ λ < 1 and 0 < ǫ < 1, we get the desired uniform estimate with respect to λ. Theorem 2.5. Assume e λ0t f (t, x, w) ∈ L 2 ([0, ∞) × Ω × S 1 ), h( x, w) ∈ L 2 (Ω × S 1 ) and e λ0t g(t, x 0 , w) ∈ L 2 ([0, ∞) × Γ − ) for some λ 0 > 0. Then for the unsteady neutron transport equation (2.1), there exists λ * 0 satisfying 0 < λ * 0 ≤ λ 0 and a unique solution u(t, x, w) ∈ L 2 ([0, ∞) × Ω × S 1 ) satisfying 1 for any 0 ≤ λ ≤ λ * 0 . When λ 0 = 0, we have λ * 0 = 0. Proof. We divide the proof into several steps: Step 1: Weak formulation. In the weak formulation (2.30), we may take the test function φ = u λ to get the energy estimate Hence, this naturally implies On the other hand, we can square on both sides of (2.28) to obtain Multiplying a sufficiently small constant on both sides of (2.56) and adding it to (2.55) to absorb u λ Hence, we have A simple application of Cauchy's inequality leads to Taking C sufficiently small, we can divide (2.58) by ǫ 2 to obtain Step 2: Convergence. Since above estimate does not depend on λ, it gives a uniform estimate for the penalized neutron transport equation (2.16). Thus, we can extract a weakly convergent subsequence u λ → u as λ → 0. The weak lower semi-continuity of norms · L 2 ([0,t]×Ω×S 1 ) and · L 2 ([0,t]×Γ + ) implies u also satisfies the estimate (2.60). Hence, in the weak formulation (2.30), we can take λ → 0 to deduce that u satisfies equation (2.1). Also u λ − u satisfies the equation By a similar argument as above, we can achieve When λ → 0, the right-hand side approaches zero, which implies the convergence is actually in the strong sense. The uniqueness easily follows from the energy estimates.
Step 3: x 0 , w) for x 0 ∈ ∂Ω and w · n < 0. (2.63) Similar to the argument in Step 1, we can obtain Then when λ is sufficiently small, we have which implies exponential decay of u.
Then for the unsteady neutron transport equation for any 0 ≤ λ ≤ λ * 0 . When λ 0 = 0, we have λ * 0 = 0. Proof. We divide the proof into several steps to bootstrap an L 2 solution to an L ∞ solution: Step 1: Double Duhamel iterations. The characteristics of the equation (2.1) is given by (2.11). Hence, we can rewrite the equation (2.1) along the characteristics as where the backward exit time t b is defined as (2.14). For the convenience of analysis, we transform it into a simpler form Here a ∧ b denotes min{a, b}. Note we have replacedū by the integral of u over the dummy velocity variable w t . For the last term in this formulation, we apply the Duhamel's principle again to u(ǫ 2 s, x−ǫ(t/ǫ 2 −s) w, w t ) and obtain where we introduce another dummy velocity variable w s and Step 2: Estimates of all but the last term in (2.69). We can directly estimate as follows: Step 3: Estimates of the last term in (2.69). We can first transform the last term I in (2.69) into by substitution s → s * = (s − t/ǫ 2 + t b ) and r → r * = (r − s + s b ). Now we decompose the right-hand side in (2.77) as for some δ > 0. We can estimate I 1 directly as Then we can bound I 2 as (2.80) By the definition of t b and s b , we always have ǫ 2 (r * +s * +t/ǫ Hence, we may interchange the order of integration and apply Hölder's inequality to obtain (2.81) Note w t ∈ S 1 , which is essentially a one-dimensional variable. Thus, we may write it in a new variable ψ as w t = (cos ψ, sin ψ). Then we define the change of (2.82) Therefore, for s b − r * ≥ δ, we can directly compute the Jacobian Thus, in order to guarantee the Jacobian is strictly positive, we may further decompose I 2 into I 2,1 + I 2,2 where in I 2,1 , we have w 1 cos ψ + w 2 sin ψ > 1 − δ and in I 2,2 , we have w 1 cos ψ + w 2 sin ψ ≤ 1 − δ. Since w = (w 1 , w 2 ) ∈ S 1 , based on trigonometric identity, we obtain (2.84) Hence, we may simplify (2.81) as Then we may further utilize L 2 estimate of u in Theorem 2.5 to obtain Step 4: L ∞ estimate. In summary, collecting (2.71), (2.72), (2.73), (2.74), (2.75), (2.76), (2.79) and (2.86), for fixed 0 < δ < 1, we have Then we may take 0 < δ ≤ 1/2 to obtain Taking supremum of u over all (t, x, w), we have Step 5: By a similar argument as in Step 3 and Step 4, combined with the L 2 decay, we can finally show the desired estimate.

Maximum Principle.
Theorem 2.7. When f = 0, the solution u(t, x, w) to the unsteady neutron transport equation (2.1) satisfies the maximum principle, i.e.
Proof. We claim that it suffices to show u(t, x, w) ≤ 0 whenever g(t, x 0 , w) ≤ 0 and h( x, w) ≤ 0. Suppose the claim is justified. Then define x 0 , w) for w · n < 0 and x 0 ∈ ∂Ω, Hence, m − h ≤ 0 and m − g ≤ 0 implies u 2 ≤ 0, which is actually u ≥ m. Therefore, the maximum principle is established.
x 0 , w) ≤ 0 and h( x, w) ≤ 0. In the proof of Lemma 2.2, we have shown u k λ → u λ in L ∞ as k → ∞. Therefore, we have u λ (t, x, w) ≤ 0. Based on the proof of Lemma 2.5, we know u λ → u in L 2 as λ → 0, where u is the solution of the equation (2.1). Then we naturally obtain u ≤ 0. Also, this is the unique solution to the equation (2.1). This justifies the claim and completes the proof.
Then for the unsteady neutron transport equation (2.1), there exists a unique solution Proof. Since the equation (2.1) is a linear equation, then we can utilize the superposition property, i.e. we can separate the solution u = u 1 + u 2 where u 1 satisfies the equation x 0 , w) for w · n < 0 and x 0 ∈ ∂Ω, (2.103) and u 2 satisfies the equation x 0 , w) = 0 for w · n < 0 and x 0 ∈ ∂Ω, Note that the data in (2.103) and (2.104) satisfy the compatibility condition (1.3). Therefore, we can apply the previous results in this section. Corollary 2.8 yields Also, Theorem 2.6 leads to Combining (2.105) and (2.106), we have the desired result.
Finally, we can apply Theorem 2.9 to the equation (1.1) and obtain Theorem 1.1.
. Then for the unsteady neutron transport equation

Asymptotic Analysis
In this section, we construct the asymptotic expansion of the equation (1.1).
3.1. Discussion of Compatibility Condition. The initial and boundary data satisfy the compatibility condition h( x 0 , w) = g(0, x 0 , w) for w · n < 0. (3.1) Then in the half-space w · n < 0 at (0, x 0 , w), the equation is valid, which implies In order to show the diffusive limit, the condition (3.3) holds for arbitrary ǫ. Since g and h are all independent of ǫ, we must have for w · n < 0, The relation (3.6) implies the improved compatibility condition h( x 0 , w) = g(0, x 0 , w) = C 0 for w · n < 0, (3.7) for some constant C 0 . This fact is of great importance in the following analysis.

Interior Expansion.
We define the interior expansion as follows: where U ǫ k can be defined by comparing the order of ǫ via plugging (3.8) into the equation (1.1). Thus, we have The following analysis reveals the equation satisfied by U ǫ k : Plugging (3.9) into (3.10), we obtain Integrating (3.14) over w ∈ S 1 , we achieve the final form Note that in order to determine U ǫ k , we need to define the initial data and boundary data.

(3.21)
We define the initial layer expansion as follows: where U ǫ I,k can be determined by comparing the order of ǫ via plugging (3.22) into the equation (3.21). Thus, we have . . .
The following analysis reveals the equation satisfied by U ǫ I,k : Integrate (3.23) over w ∈ S 1 , we have  Similarly, we can derive U ǫ I,k (τ, x, w) for k ≥ 1 satisfies (3.31)

Boundary Layer Expansion with Geometric Correction.
In order to determine the boundary condition for U ǫ k , we need to define the boundary layer expansion. Hence, we need several substitutions:
3.5. Initial-Boundary Layer Expansion. Above construction of initial layer and boundary layer yields an interesting fact that at the corner point (t, x) = (0, x 0 ) for x 0 ∈ ∂Ω, the initial layer starting from this point has a contribution on the boundary data, and the boundary layer starting from this point has a contribution on the initial data. Therefore, we have to find some additional functions to compensate these noises. The classical theory of asymptotic analysis requires the so-called initial-boundary layer, where both the temporal scaling and spacial scaling should be used simultaneously. Fortunately, based on our analysis, the improved compatibility condition (3.7) implies the value at this corner point is a constant for w · n < 0. Then These contribution must be zero at the zeroth order, i.e.
Therefore, the zeroth order initial-boundary layer is absent.
3.6. Construction of Asymptotic Expansion. The bridge between the interior solution, the initial layer, and the boundary layer is the initial and boundary condition of (1.1). To avoid the introduction of higher order initial-boundary layer, we only require the zeroth order expansion of initial and boundary data be satisfied, i.e. we have The construction of U ǫ k . U ǫ I,k and U ǫ B,k are as follows: Assume the cut-off function ψ and ψ 0 are defined as and define the force as Step 1: Construction of zeroth order terms. The zeroth order boundary layer solution is defined as f ǫ 0 (t, 0, θ, φ) = g(t, θ, φ) for sin φ > 0, lim η→∞ f ǫ 0 (t, η, θ, φ) = f ǫ 0 (t, ∞, θ).
The zeroth order initial layer is defined as
Step 2: Construction of first order terms. Define the first order boundary layer solution as (3.57) Define the first order initial layer as (3.58) Define the first order interior solution as (3.59) Step 3: Construction of U ǫ 2 and U ǫ 2 . Define the second order boundary layer solution as Define the second order initial layer as

(3.61)
Define the first order interior solution as (3.62) Step 4: Generalization to arbitrary k. Similar to above procedure, we can define the k th order boundary layer solution as (3.64) Define the k th order interior solution as (3.65) When g and h are sufficiently smooth, then all the functions defined above are well-posed. The key point here is in the boundary layer, the source term including ∂ θ U ǫ B,k is in L ∞ due to the substitution (3.39).

Diffusive Limit
In this section, we prove the first part of Theorem 1.2.
Theorem 5.1. Assume g(t, x 0 , w) ∈ C 4 ([0, ∞) × Γ − ) and h( x, w) ∈ C 4 (Ω × S 1 ). Then for the unsteady neutron transport equation (1.1), the unique solution u ǫ (t, x, w) ∈ L ∞ ([0, ∞) × Ω × S 1 ) satisfies where the interior solution U ǫ 0 is defined in (3.56), the initial layer U ǫ I,0 is defined in (3.55), and the boundary layer U ǫ B,0 is defined in (3.54). Proof. We divide the proof into several steps: Step 1: Remainder definitions. We may rewrite the asymptotic expansion as follows: The remainder can be defined as Noting the equation is equivalent to the equations (3.21) and (3.40), we write L to denote the neutron transport operator as follows: Step 2: Estimates of LQ N . The interior contribution can be estimated as Similarly, for higher order term, we can estimate Step 3: Estimates of LQ I,N . The initial layer contribution can be estimated as Similarly, we have Step 4: Estimates of LQ B,N . The boundary layer solution is . Notice ψ 0 ψ = ψ 0 , so the boundary layer contribution can be estimated as (5.21) It is easy to see Since ψ 0 = 1 when η ≤ 1/(4ǫ), the effective region of ∂ η ψ 0 is η ≥ 1/(4ǫ) which is further and further from the origin as ǫ → 0. By Theorem 4.2, the first term in (5.21) can be controlled as For the second term in (5.21), we have Similarly, for higher order term, we can estimate It is obvious that Away from the origin, the first term in (5.26) can be controlled as For the second term in (5.26), we have Consider the asymptotic expansion to N = 4, then the remainder R 4 satisfies the equation x 0 , w) for w · n < 0 and x 0 ∈ ∂Ω. Note that the initial data and boundary data are nonzero due the contribution of initial layer and boundary data at the point (t, x) = (0, x 0 ). By Theorem 2.9, we have Hence, we have Since it is easy to see

Counterexample for Classical Approach
In this section, we present the classical approach in [1] to construct asymptotic expansion, especially the boundary layer expansion, and give a counterexample to show this method is problematic in unsteady equation.
U k (t, η, θ, ξ)dξ. (6.14) 6.3. Classical Approach to Construct Asymptotic Expansion. Similarly, we require the zeroth order expansion of initial and boundary data be satisfied, i.e. we have The construction of U k , U I,k , and U B,k by the idea in [1] can be summarized as follows: Assume the cut-off function ψ and ψ 0 are defined as (3.51) and (3.52).

(6.17)
The zeroth order initial layer is defined as Then we can define the zeroth order interior solution as where (t, x, w) is the same point as (τ, η, θ, ξ).
Step 2: Construction of first order terms.
Define the first order boundary layer solution as Define the first order initial layer as Define the first order interior solution as Step 3: Construction of second order terms. Define the second order boundary layer solution as x 0 , w) for sin(θ + ξ) > 0, lim η→∞ f 2 (t, η, θ, ξ) = f 2 (t, ∞, θ).

(6.23)
Define the second order initial layer as
Similar to above procedure, we can define the k th order boundary layer solution as Define the k th order initial layer as (6.28) By the idea in [1], we should be able to prove the following result: