A Monotone Discretization for the Fractional Obstacle Problem

We introduce a novel monotone discretization method for addressing obstacle problems involving the integral fractional Laplacian with homogeneous Dirichlet boundary conditions over bounded Lipschitz domains. This problem is prevalent in mathematical finance, particle systems, and elastic theory. By leveraging insights from the successful monotone discretization of the fractional Laplacian, we establish uniform boundedness, solution existence, and uniqueness for the numerical solutions of the fractional obstacle problem. We employ a policy iteration approach for efficient solution of discrete nonlinear problems and prove its finite convergence. Our improved policy iteration, adapted to solution regularity, demonstrates superior performance by modifying discretization across different regions. Numerical examples underscore the method's efficacy.


Introduction
In this work, we consider the obstacle problem associated with the integral form of the fractional Laplacian, referred to as the fractional obstacle problem.
For n ≥ 1, let Ω ⊂ R n be a bounded Lipschitz domain satisfying the exterior ball condition, with boundary denoted as ∂Ω.Additionally, Ω c represents the complement of Ω.The specific form of the fractional obstacle problem is as follows: Given two prescribed functions f : Ω → R and the obstacle function ψ : R n → R with the nondegeneracy condition ψ| Ω c < 0, and s ∈ (0, 1), we seek a function u : R n → R that satisfies the following nonlinear equation with homogeneous Dirichlet boundary conditions: (1) Here, the integral fractional Laplacian of order s ∈ (0, 1), defined by (2) (−∆) s u(x) := C n,s P.V.
The normalization constant is given by C n,s := The fractional obstacle problem finds applications in various fields, such as mathematical finance theory, systems of particles and elasticity problems [25,9,22,23,13].For a comprehensive overview of obstacle problems, including traditional obstacle problems and those related to integral-differential operators (of which the fractional obstacle problem is a particular case), we refer to the survey article [19].
For fractional operators, it has been shown that the maximum principle holds [18].In terms of numerical stability, it is desirable to preserve this property at the discrete level, known as the monotonicity of the scheme.Constructing monotone schemes is a natural advantage of finite difference methods.For the Laplace operator (−∆), the monotone difference scheme has been paid attention to by researchers for a long time, [5,26].In the context of the finite element method, achieving monotonicity in the discretization process depends on specific grid conditions [27].
Furthermore, considering the regularity of the solution is essential for the convergence of numerical schemes.The fractional obstacle problem has received significant attention in mathematics, particularly in the field of PDE regularity of solutions and the free boundary.In this regard, a series of works [24,7,8] have established the Hölder regularity of C 1,s (R n ) for solutions to the s-order obstacle problem.For the case of a bounded domain, it is expected that the interior regularity exhibits the same C 1,s Hölder regularity.However, lower regularity may arise near the boundary of the domain.For the fractional Laplace equation, a well-known result is that when the problem domain satisfies the exterior ball condition and the righthand side belongs to L ∞ , the solution exhibits singularity near the boundary ∂Ω in terms of dist(x, ∂Ω) s , resulting in global C s Hölder regularity.This characteristic is observed widely in a large class of nonlocal elliptic equations, as discussed in the survey article [18].Therefore, the global regularity (boundary regularity) of the fractional obstacle problem on bounded domains still requires careful investigation.
In recent years, numerous numerical algorithms have been developed for solving fractional-order partial differential equations (PDEs).For the fractional Laplace equation, the finite element method can be employed, utilizing the variational formulation and accounting for the low regularity near the boundary [1].Similarly, the fractional obstacle problem can be formulated as a variational inequality problem, where the solution u is sought to minimize a functional subject to the constraint u ≥ ψ [24].To address this variational inequality, [4] utilizes the finite element method, establishing interior and boundary regularity results for bounded fractional obstacle problems, and providing error estimates based on these regularity results.
However, for fractional operators, due to their non-local nature, it is challenging to establish general conditions that guarantee the monotonicity of finite element discretization.This challenge is further amplified when dealing with the finite element discretization of nonlinear problems.
On the other hand, several finite difference discretizations that satisfy monotonicity have been proposed [15,10,11] for the fractional operator (−∆) s .However, most of these works focus on problems solved on structured grids and make overly strong assumptions about the Hölder regularity of the true solution.In our recent work [14], we propose a monotone difference scheme based on quadrature formulas, which can be applied to solve the fractional Laplace equation on general bounded Lipschitz domains.This work provides a comprehensive analysis of the scheme's consistency, taking into account the actual regularity of the problem.Furthermore, we obtain rigorous pointwise error estimates by utilizing discrete barrier functions.
The objective of this study is to develop a monotone scheme for solving the problem defined by equation ( 1) on a bounded Lipschitz domain and to devise an efficient solution solver.Building upon our previous work [14], we extend this discretization to the fractional obstacle problem.A crucial aspect of our work is the introduction of the enhanced discrete comparison principle, as demonstrated in Lemma 3.4, for the discrete fractional Laplacian operator.By utilizing this result, we establish that the discretization scheme maintains monotonicity and satisfies the discrete maximum principle, guaranteeing the uniqueness of the discrete solution.Additionally, we employ the discrete Perron method to establish the existence and uniform boundedness of the discrete solution.
In addition, this study explores the policy iteration method for solving nonlinear discrete problems.Leveraging the enhanced discrete comparison principle, we establish the finite convergence of the standard policy iteration method.Furthermore, the paper provides an in-depth discussion on the relationship between the discretization scheme for the fractional Laplace equation and the regularity results of the problem.Based on this discussion, we propose an improved policy iteration method that considers the low regularity near the contact set of the problem and selects different discrete scales for different regions.Numerical results demonstrate the superior performance of the improved method compared to the standard policy iteration method.
The remaining part of the paper is organized as follows.In Section 2, we introduce some necessary preliminary results.In Section 3, we present the monotone discretization for the fractional obstacle problem and discuss the properties of the scheme.In Section 4, we apply the policy iteration method for solving nonlinear discrete problems and prove the convergence of the iteration.Additionally, taking into account the specific regularity of the problem, we propose an improved policy iteration method that exhibits better numerical performance in practical computations.Numerical examples are presented in Section 5 to support our theoretical results.Furthermore, for the boundary Hölder regularity of the problem, a simple proof is provided in Appendix A.

Preliminary results
In this section, some preliminary knowledge and results will be introduced, including notation, definitions of function spaces, and regularity results.Given an open set U ⊂ R n with ∂U ̸ = ∅, the C β (U ) Hölder seminorm, with β > 0, is denoted by | • | C β (U ) .More precisely, for β = k + t with k integer and t ∈ (0, 1], define As usual, the nonessential constant denoted by C may vary from line to line.X ≲ Y means that there exists a constant C > 0 such that X ≤ CY , and X ≂ Y means that X ≲ Y and Y ≲ X.
2.1.Regularity: the linear problem.Here are some well-known conclusions regarding the regularity of the integral fractional Laplacian: The first lemma pertains to the interior regularity of fractional harmonic functions.This result establishes that the solution of (1) is smooth within the region unaffected by ψ, or in other words, where there is no contact with ψ.
For the global regularity result, however, it can be proven that the solution belongs to C s (R n ), and this result is sharp [20].This is due to the limited regularity at the boundary and suggests the need for a more comprehensive discussion using weighted Hölder spaces, as discussed in [20,1].

2.2.
Regularity: the fractional obstacle problem.Next, some regularity results for the fractional obstacle problem are presented.It is worth noting that, when discussing regularity, the assumption that the forcing term f is zero can be made without loss of generality.This is because the general solution can be decomposed as follows: u = u 1 + w f , where u 1 solves (1) with zero forcing (f = 0) and obstacle ψ 1 := ψ − w f .Here, w f is the solution to the linear problem (3).
The contact set is defined as follows: The first result on the interior Hölder regularity of the fractional obstacle problem was established by Silvester in R n [24].When considering the fractional obstacle problem in bounded domains, Nochetto, Borthagaray, and Salgado in [4] derived a corresponding result under the assumption: The main technique employed in their work is the use of truncation functions.The result is presented as follows.

(boundary Hölder regularity).
Let Ω be a bounded Lipschitz domain satisfying the exterior ball condition, and ψ ∈ C( Ω) ∩ C 2,1 (Ω) with nondegeneracy condition ψ| Ω c < 0. Let u solves (1), then In fact, in the boundary Hölder regularity estimate, the smoothness assumption on ψ can be relaxed to ψ ∈ L ∞ (Ω).The main idea is to revert back to the fundamental theoretical framework of boundary regularity [20] and explore its extension to the fractional obstacle problem.Additionally, compared to [4], the proof is more direct.The detailed proof of this new approach is included in the Appendix A.
Before concluding this section, let us provide an intuitive summary of the implications of the regularity results for the design of numerical schemes: (1) Lemma 2.1 (balayage) indicates that the solutions are smooth within both Λ and Ω \ Λ if ψ is smooth.Therefore, the standard discretization methods would be suitable for these regions.(2) By combining the interior regularity discussed in Proposition 2.2, it suggests that the smoothness across the contact boundary decreases to C 1,s .As a result, discretization methods that rely on higher-order derivatives would be inappropriate for this region.(3) The boundary regularity result, Proposition 2.3, is similar to that of linear problems.Thus the techniques used for linear problems, such as the regularity under weighted norms, are still applicable.From a numerical perspective, this phenomenon motivates the use of graded grids to enhance the convergence order, as explored in the work by [1,14].This issue will be revisited in Section 5.

A monotone discretization
For simplicity, it is assumed that Ω is a polygon in 2D and a polyhedron in 3D.Let T h be a triangulation of the computation domain Ω, i.e. ∪ T ∈T h T = Ω.For T ∈ T h , let h T denote the diameter of element T , ρ T denote the radius of the largest inscribed ball contained in T .The triangulation is referred to as local quasi-uniform if there exists a constant λ 1 such that Meanwhile, the triangulation T h is called shape regular, if there is λ 2 > 0 such that h T ≤ λ 2 ρ T .For n ≥ 2, the shape regular property implies the local quasi-uniform properties, see [6].Let N h denote the node set of T h and N b h := {x i ∈ N h : x i ∈ ∂Ω} be the collection of boundary nodes, and interior node set is denoted by , where P 1 (T ) denotes the linear function collection defined on element T .In [14], a monotone discretization for the integral fractional Laplacian was proposed and a corresponding pointwise error estimate was conducted under realistic Hölder regularity assumption.In this work, the monotone discretization for the integral fractional Laplacian is applied, denoted by (−∆) s h : ( Here, Ω i ⊂ Ω is a star-shaped domain centered at x i satisfying some symmetrical condition [14, Equ.(3.5) and (3.6)] with a typical scaling (6) , where h i is the mesh size around x i , and δ i := dist(x i , ∂Ω).A concrete example that satisfies the condition is the n-cubic domain which is also utilized in our numerical experiments.The parameter α i is carefully selected to address the regularity concern of integral fractional Laplacian [20,14].The singular integral within Ω i is approximated using the finite difference method, denoted as ∆ , where e j represents the unit vector of the jth coordinate.The coefficient κ n,s,i is a known positive constant, provided in [14,Equ. (3.10)].
3.1.Properties of (−∆) s h .In this subsection, several fundamental properties of the (−∆) s h operator are revisited and established.These properties will be utilized in subsequent discussions.One of the key features is its monotonicity, as proven in [14, Lemma 3.4], which is closely related to the matrix discretization being an M -matrix.An even stronger property of (−∆) s h will be explored: the discretization matrix exhibits a strong diagonal dominance.Additionally, the consistency of (−∆) s h will be discussed.To begin with, let us revisit the discrete barrier function and monotone property given in [14, Lemma 3.4 & Lemma 3.5].
, where the constant C depends only on s and Ω.
The above monotonicity result is insufficient for the obstacle problem, as the fractional order operator equation in the obstacle problem holds only in certain regions of the domain.Let us denote the resulting discrete matrix of (−∆) s h as L ∈ R N ×N , where N is the number of interior nodes.Below, it will be shown that L is a strongly diagonal dominant M -matrix.
Proof.The property (8a) follows directly from the definition of (−∆) s h .In fact, for the equation ( 5) at x i , it is straightforward to observe that n+2s) dy.Furthermore, all the other contributions to the off-diagonal terms are non-positive.
Utilizing Lemma 3.1, the discrete barrier function can be expressed as the vector B ∈ R N , where all entries are equal to one.Therefore, (7) implies that which leads to (8b).□ Now, the enhanced discrete comparison principle applicable to the obstacle problem will be presented.This principle plays a crucial role in the subsequent analysis.

Lemma 3.4 (enhanced discrete comparison principle for (−∆)
where Λ h ⊂ N 0 h is an arbitrary subset.Then, v h ≥ w h in Ω. Proof.Since v h , w h are piecewise linear functions, it suffices to prove v h (x i ) ≥ w h (x i ) for all x i ∈ N 0 h .Let z h := v h − w h , and Z ∈ R N denotes its coefficient under the basis function.Our objective is to demonstrate that Z ≥ 0. Note that the resulting discrete matrix of (−∆) s h , denoted as L, exhibits strong diagonal dominance.Then, (9) has the matrix representation Here, the indices of all interior points in Λ h are combined into the first group.The remaining indices form the second group.This implies that L11Z 1 + L12Z 2 ≥ 0 and Z 2 ≥ 0. By above lemma (strongly diagonal dominant M -matrix), the entries in the offdiagonal block L 12 are non-positive, which yields Since L is a strongly diagonal dominant M -matrix, so do L 11 , that is (L −1 11 ) ij ≥ 0. Therefore, it can be deduced that Z 1 ≥ 0 and Z ≥ 0, which completes the proof.□ 3.2.Numerical scheme.Now, the numerical scheme for solving the fractional obstacle problem is proposed: It is noted that, although not explicitly stated, there is a parameter α in the discretization of the integral fractional Laplacian that is influenced by the regularity of the solution u.The discrete comparison principle for G h will be demonstrated, from which the uniqueness of (10) follows directly.

Lemma 3.5 (discrete comparison principle for
h .The interior points will be classified into two sets based on the values of G h [w h ]: Furthermore, for any x i ∈ N 0 h \ Λ w h , the following inequalities are true: ).This leads to the desired result by taking Λ h = Λ w h in Lemma 3.4 (enhanced discrete comparison principle for (−∆) s h ).□ Corollary 3.6 (uniqueness).For the problem (10), the discrete solution is unique.

Stability and existence.
The existence and stability of (10) will be shown in this section.Theorem 3.7 (existence and stability).There exists a unique u h ∈ V h that solves (10).The solution u h is stable in the sense that ∥u h ∥ L ∞ (Ω) ≤ C, where the constant C is independent of h.

Proof.
To establish the existence, a monotone sequence of discrete functions {{u k h } ∞ k=0 } is constructed, starting with the initial guess u 0 h ∈ V h that satisfies the following condition: where the discrete barrier function b h is defined in (7), and the constant E > 0 will be specified later.By the Lemma 3.1 (discrete barrier function), the following relation holds: .
Step 2 (Perron construction).Induction is employed.Suppose that there already exists a discrete function h and satisfaction of ( 12) is as follows: All interior nodes are considered sequentially, and auxiliary functions u k,i−1 h ∈ V h are constructed using the first i − 1 nodes, starting with h ](x i ) = 0.This is achievable since the discrete matrix of (−∆) s h satisfies L ii > 0, as stated in (8a).Moreover, at other nodes This process is repeated with the remaining nodes x j for i < j ≤ N , and u k+1 h := u k,N h is set as the last intermediate function.By construction, the following is obtained: k=1 is monotonically decreasing and clearly clearly bounded from below by ψ(x i ).Hence, the sequence converges, and the limit

and if the last inequality were strict,
Step 2 could be applied to further improve u h .This demonstrates the existence of a discrete solution to (10).The above proof also implies which is the uniform bound, as asserted.□ Remark 3.8 (discussion on the consistency and convergence).An important concept in numerical methods is its consistency, which involves investigating the prop- | tending to zero as the mesh size decreases.In fact, for the integral fractional Laplacian operator, a detailed analysis in [14] reveals that the discretization achieves consistency at interior points that are a constant number of mesh sizes away from the boundary.The inconsistency near boundary points can be addressed by incorporating a discrete barrier function (7).This analysis technique is commonly used in the convergence analysis of semi-Lagrangian or twoscale method [12,17] and can be viewed as an extension to the Barles-Souganidis [2] analysis framework for nonlinear problems.Due to space limitations, elaboration on the corresponding results will not be provided here.

Solver
In this section, the policy iteration (Howard's algorithm) is employed to solve the discrete fractional obstacle problem (10).It is worth noting that policy iteration is a well-established and extensively studied technique in dynamic control problems.For the problem (10), the convergence of the policy iteration is established by leveraging the monotonicity of the discrete operator.Furthermore, an improved policy iteration will be proposed by incorporating prior knowledge, such as the regularity of the solution to the fractional obstacle problem discussed in Section 2.

Policy iteration.
Policy iteration is utilized to solve the discrete problem (10).The algorithm is outlined as follows.
Algorithm 1: Policy iteration for the fractional obstacle problem (10) Input: right hand side f , obstacle function ψ.
Output: solution u h .
Initialize u h , then stop; Otherwise, continue to Step 1.
There exists a convergence result for policy iteration when solving the nonlinear problem (10), which was initially proven in [3].Here, a brief overview of the proof using our notation is provided for clarity.It is worth noting that the monotonicity property plays a crucial role in the convergence analysis.

Theorem 4.1 (convergence of policy iteration). The sequence {u (k)
h } ∞ k=0 given by above algorithm satisfies: (2) The algorithm converges in at most N iterations.
Proof.With u h can be defined as N 0 h to ensure that (13) holds for the initialization.
Step 1 (increasing sequence).(13) implies that for k ≥ 0 Then, the definition of Next, using the update condition ( 13), the following relation holds: .
By combining Lemma 3.4 (enhanced discrete comparison principle for (−∆) s h ), the increasing sequence is shown as Step 2 (strictly decreasing discrete contact set).First, the increasing sequence {u h , the update condition (13) and the above property implies (−∆) s h [u for the definition of C (k+1) h in Step 1.This means h , then the iteration clearly stops at (k + 1)-th step, due to the update condition (13).That is, the discrete contact set is strictly decreasing unless the iteration stops, which makes at most N iterations.□ Remark 4.2 (convergence profile of Algorithm 1).The N -step iteration is proven to be sharp in the general case, as demonstrated by a specific example presented in [21].However, in numerical experiments, convergence is often observed within just a few iterations.

4.2.
An improved policy iteration.Next, improvements will be introduced to the discrete operator in the policy iteration algorithm.It is observed that the scaling H i of the singular integral part of (−∆) s h over the region Ω i can be adjusted, as shown in expression (6).In the policy iteration, modifications are made to the selection of H i , primarily for the following reasons: (1) During the update of u h in Step 3, if the discrete contact set C h has been updated, equation ( 13) can be seen as the discrete fractional Laplacian on the non-contact points N 0 h \ C h .Instead, the values at the contact points C h can be viewed as the Dirichlet "boundary" data.
(2) In the discretization of the fractional Laplacian, the singular part is approximated by a scaled Laplace operator [14].However, it is important to note that the regularity estimate for the obstacle problem, Proposition 2.2, indicates that the true solution has only C 1,s continuity over the boundary of the contact set.This limitation restricts the accuracy of such approximation.Based on the above, when discretizing the fractional Laplacian on N 0 h \ C h , it is necessary to restrict the size of H i to ensure that Ω i does not intersect with the region relating to the discrete contact set C h .Therefore, the following improved strategy is proposed: ( 14) , where θ is a constant that depends on the shape regularity of the mesh.In the numerical experiments, θ is set to 1/4.
Since the discrete contact set in policy iteration is continuously updated, the improvement in the algorithm primarily lies in updating the discretization of the fractional Laplacian through ( 14) after updating the contact set.The improved algorithm is outlined as follows.
Algorithm 2: Improved policy iteration for the fractional obstacle problem (10) Input: right hand side f , obstacle function ψ.
by choosing h , then stop; Otherwise, continue to Step 1.
Remark 4.3 (convergence profile of Algorithm 2).In the improved algorithm, it is important to note that the solutions in the consecutive steps correspond to different discretized fractional Laplacian operators, as seen in (15).Therefore, it is not possible to directly apply Theorem 4.1 to provide a rigorous convergence analysis.However, it is worth noting that although the discretized fractional Laplacian operators vary across steps, they are all appropriate discretizations of the continuous operator and maintain monotonicity.From this perspective, the convergence behavior of the improved policy iteration algorithm should be similar to the original version.Indeed, this observation has been confirmed in numerical experiments, where the improved algorithm exhibits better performance, particularly near the contact interface.

Numerical Experiments
In this section, several numerical experiments are conducted to test the stability and effectiveness of the methods, as well as the convergence of the improved policy iteration.
According to the discussion in Section 2.2, the fractional obstacle problem exhibits low regularity at the domain boundary, which is consistent with the fractional linear problem.To address this issue, various approaches have been proposed, including the use of graded grids.In this context, as discussed in [1,4,14], the concept of graded grids is introduced with a parameter h.Let µ ≥ 1 be a constant such that for any T ∈ T h , ( 16) It is worth noting that the quasi-uniform grid corresponds to the case where µ = 1.Additionally, the choice of µ has an impact on both the grid quality and the overall computational complexity.The value of µ will be specified in each experiment.
In the discretization of the fractional Laplacian operator, the scale H i of the singular approximation part not only depends on the grid but also relies on α i , as shown in Equation ( 6).The optimal choice of α depends on the local smoothness of the solution or the smoothness of the forcing term f , as discussed in [14].In the tests conducted in this section, a smooth f is assumed, and based on the findings in [14, Theorem 6.1], the optimal choice of α is used, i.e., α = 1 2 .5.1.Convergence order test.In this 1D test, the convergence rate is tested on both uniform and graded grids.The domain considered is Ω = (−1, 1), and an explicit solution for problem (1) is constructed as follows, A direct calculation shows that Next, the following forcing and obstacle terms are considered: 0) and (−∆) s u − f = 0 outside of this set.Consequently, u represents the solution of problem (1), and the contact set is B 1 2 (0). Figure 1 illustrates the convergence rates for both uniform and graded grids.According to [14, Section 7.1], the optimal choice for the pointwise estimate of the linear fractional Laplacian is µ = 2−s s , which will also be used in this test.The convergence rate over uniform grids is observed to be approximately s, which is consistent with the behavior observed in the quasi-uniform grids for the linear case.Moreover, the graded grid yields a significantly improved convergence rate.It is worth noting that for the linear case, the convergence order under the optimal choice is 2 − s, which is also observed for some cases (e.g., s = 0.6).However, for other cases, the strategy of choosing H i in (14) based on the local distance to the  In Figure 2, the numerical solutions u h obtained for different values of s are presented below.A clear qualitative difference is observed between the solutions for different choices of s.When s = 0.9, the discrete solution resembles the expected solution of the classical obstacle problem, where the operator (−∆) s is effectively replaced by −∆.As s decreases, the solution approaches ψ + .Throughout the experiments, a consistent observation is made that as s increases, the contact set decreases and always contains the point x 0 .
The difference between the improved method and the original method is remarkable.Observations reveal that when using the original method with larger values of s, oscillations occur near the free boundary due to the handling of the singular part across it.On the other hand, the improved method yields a numerically smoother solution without oscillations.As s decreases, the solutions obtained by the two methods become closer to each other, and the differences diminish.As mentioned in Remark 4.3, it is highly likely that the improved algorithm possesses the property of an increasing solution sequence (or decreasing discrete contact set) as stated in Theorem 4.1 for the standard algorithm.In this experiment, as demonstrated in Figure 3 1.Number of iterations for s = 0.3, 0.6, 0.9 using graded grids (µ = 2−s s ), with different number of degrees of freedom (DOFs).
Furthermore, the actual number of iterations obtained in the experiments is presented in Table 1.It can be observed that the actual number of iterations required is significantly smaller than the theoretical upper bound N .Additionally, as the grid is refined, the algorithm's performance remains stable.This demonstrates the efficiency of our algorithm.5.4.2D test.Finally, problem (1) is considered in the domain B 1 (0) ⊂ R 2 with f = 0, and the obstacle function is given by For s = 0.1, the numerical solution requires 13,395 degrees of freedom, and the number of improved policy iterations is 4. Conversely, for s = 0.9, the numerical solution involves 4,253 degrees of freedom, and the number of improved policy iterations is 10.The number of iterations increases as s grows, but significantly less than N , as also observed in the 1D case.As expected, the 2D numerical results demonstrate similar qualitative characteristics as observed in 1D.

Conclusion
In this paper, a discrete scheme for the fractional obstacle problem is proposed, which is based on the monotone discretization for the fractional Laplace operator.Utilizing the distinctive structure of the problem (10) and conducting a comprehensive study of the discrete operator, the study reveals that the nonlinear discrete operator G h upholds the discrete comparison principle.Based on this property, the existence, uniqueness, and uniform boundedness of numerical solutions are established.
Due to the limited regularity of the true solution near the domain boundary, a graded grid has been introduced to capture this behavior.For solving nonlinear problems, the policy iteration method is used.Benefiting from the monotonicity of the discrete scheme of the fractional Laplace operator, the policy iterative process can converge in finite iterations.Moreover, considering the reduced regularity of the actual solution near the contact set boundary, the discretization of the fractional Laplacian has been adaptively refined by iteratively updating contact nodes.This refinement has led to an improved policy iteration approach.In contrast to the Appendix A. proof of Proposition 2.3 under a weak assumption of ψ Proposition 2.3 is proven, assuming ψ ∈ L ∞ (Ω).The proof follows closely the global Hölder estimate for the linear problem [20], so only the main thread will be sketched.First, some elementary tools are recalled.
The proof of Lemma A.2 follows from the comparison principle.As a direct consequence, the L ∞ estimate of the fractional Laplace equation implies the following inequality for f ∈ L ∞ (Ω) and ψ ∈ L ∞ (Ω): Lemma A.3 (supersolution, Lemma 2.6 in [20]).There exist C 1 > 0 and a radial continuous function φ 1 satisfying By the Lemma A.2 and Lemma A.3, one could construct an upper barrier for |u| by scaling and translating the supersolution.The proof of Lemma A.4 is similar to the Lemma 2.7 in [20].The only difference is for each point x 0 ∈ ∂Ω, when constructing the upper barrier on the ball touched x 0 from outside, the radius of the ball not only depends on the exterior ball condition of Ω, but also depends on r 0 which appeared in the assumption (4) for the contact set.
Furthermore, by employing the inequality |u(x in Ω (see Lemma A.4), the following result can be deduced: Further, (18) and Lemma A.4 also imply that Now, one can use the Proposition A.1, which the C β seminorm of ũ can be bounded by (17), (18) and (19), and obtain

Figure 2 .
Figure 2. Numerical solutions to the fractional obstacle problem for s = 0.3, 0.6, 0.9, computed using the standard policy iteration (Algorithm 1) and its improved version (Algorithm 2), respectively.

5. 3 .
Convergence history of improved policy iteration.Next, the convergence history of the improved policy iteration (Algorithm 2) is examined.The experiment is conducted with f = 1, and the obstacle function is defined as follows:ψ(x) = 3 − 6|x − 1/4|.

Figure 3 .
Figure 3. Convergence history of the improved policy iteration (Algorithm 2) with s = 0.6 and 953 free vertices, including residual decay versus iterations (last picture).The red points indicate the contact set, and the blue points represent the non-contact set.

Figure 4 .
Figure 4. Discrete solutions to the fractional obstacle problem for s = 0.1 (left) and s = 0.9 (right).Top: lateral view.Bottom: top view, with the discrete contact set highlighted.
, the convergence profile is observed through the values of |C