When a Dirichlet problem on a smooth domain is approximated by a polygon, an error occurs that is suboptimal for quadratic approximation [1, 10, 11]. However, this can be corrected by a modification of the variational form . Here we review this approach and suggest a new one.
Let be a smooth, bounded, two-dimensional domain. Consider the Poisson equation with Dirichlet boundary conditions:
We assume that and are sufficiently smooth that can be extended to be in , where contains a neighborhood of the closure of .
One way to discretize (1) is to approximate the domain by polygons , where the edge lengths of are of order in size. Then conventional finite elements can be employed, with the Dirichlet boundary conditions being approximated by the assumption that on , with
appropriately defined. For example, let us suppose for the moment thatand we take as well. In particular, we assume that is triangulated with a quasi-uniform mesh of maximum triangle size , and the boundary vertices of are in . We define where
Then the standard finite element approximation finds satisfying
where . Here we assume that is extended smoothly outside of .
This approach for
(piecewise linear approximation) leads to the error estimate
However, when this approach is applied with piecewise quadratic polynomials (), the best possible error estimate is
which is less than optimal order by a factor of . The reason of course is that we have made only a piecewise linear approximation of . Table 1 summarizes some computational experiments for the test problem in Section 2.1. We see a significant improvement for quadratics over linears, but there is almost no improvement with cubics. Moreover, we will see that a significant improvement using quadratics can be obtained using simple approaches that modify the variational form.
There have been many techniques introduced to circumvent the loss of accuracy with quadratics (and higher-order piecewise polynomials) [11, 6]. However, all of them require some modification of the quadrature for the elements at the boundary.
Here we review an approach by Bramble et al.  that solves directly on , but with a modified variational form based on the method of Nitsche . The method  has been modified and applied in many ways . However, the method in  leads to a non symmetric bilinear form. Given this shortcoming we define a new method that is symmetric and solves the problem on that has similar convergence results. As we will see in the next section, one main idea in  is that one uses a Taylor series of the solution near the boundary to define appropriate boundary conditions on . We should mention that this idea has been used recently (see for example [5, 8]).
|L2 err||rate||H1 err||rate||seg||hmax|
2. The Bramble-Dupont-Thomée approach
The method  of Bramble-Dupont-Thomée (BDT) achieves high-order accuracy by modifying Nitsche’s method  applied on . We assume that and we do not necessarily assume that the boundary vertices of belong to . The bilinear form used in  is
Here, denotes the outward-directed normal to and
Contrast the definition of to the closely related function defined by
For simplicity the assume that . Then the BDT method will find such that
If were 0, this would be Nitsche’s method on .
Corrections of arbitrary order, involving terms for are studied in , but for simplicity we restrict attention to the first-order correction to Nitsche’s method given in (4). The error estimates obtained in  are as follows
Thus using the variational form (4) leads to an approximation that is optimal-order with quadratics and cubics and is only suboptimal for quartics by a factor of .
2.1. An example of a circle
We consider a numerical example. Consider the case where is a disc of radius centered at the origin, in which case we have . However, it is more difficult to evaluate . We have for , where denotes the outward normal to . We can write , and . Since , we have
Note that for , and . Since , we must pick the plus sign, so
It is not hard to see that in this case.
This problem is simple to implement using the FEniCS Project code dolfin . We take , , and in the computational experiments described subsequently. Computational results for this example are given in Table 2 where we see optimal order approximation for , improvement for over (suboptimal by a factor ), and no improvement for quintics. These errors are depicted in Figure 2.
|hmax||L2 error||rate||H1 error||rate|
3. A new method based on a Robin-type approach
One issue with the BDT method is that the resulting linear system is not symmetric, although it is possible to symmetrize the method as we discuss in Section 8. Here we develop a technique that leads to a symmetric system. Moreover, this method does not require the parameter(s) from Nitsche’s method. For Nitsche’s method to succeed, must be chosen appropriately .
We first separate to its piecewise linear part and its curvilinear part. We will assume that where is a piecewise linear segment and are and no where linear. We let the end points of to be .
For the method in this section we assume that the vertices of belong to and hence might not be a subdomain of . Hence, we need to define in this case. We assume that for every that is there is a unique smallest number in absolute value such that
We assume that the approximate domain boundary can be decomposed into three parts, as follows. Let be the edges of .
where denotes the interior of . Let . We assume the following.
We assume that all the vertices of belong to and that each (for ) is a vertex of . Finally, we assume that
Our method is based on a Robin type of boundary condition on . In fact, our method will be based on the closely related problem:
Here we define for and not a vertex of . The key here is that, using that vanishes on , for ( not a vertex of ) we have
where lies in the line segment with end points and .
Now we can write the method. We start by defining the finite element space we will use
where is a suitable approximation of and is a piecewise polynomial of degree at most on .
The bilinear form is given by
Then the method solves:
Find such that
is the linear interpolant onto. Note that we can define only knowing on . Alternatively, if we have an analytic representation of we can define as a smooth extension of outside of .
4. Error Analysis
4.1. Stability Analysis
Unfortunately, the bilinear form is not positive definite. However, we will be able to prove stability of method. In order to do so, we need to decompose the space into its boundary contribution and interior contribution. More precisely, we can write
where . We will define a norm on :
and a semi-norm
Note that is in fact a norm on .
The following crucial lemma will allow us to prove stability.
There exists a constant such that
Let be the collection of edges that are a subset of and let be triangles such that has an edge in . Then, if and using inverse estimates we have
The result is complete after we use that for . ∎
We note that may not be well defined for all . Therefore, we need to make an assumption on such that this is not the case.
We assume that is such that
For example, if has a lower bound as follows, then (9) will hold. Suppose that the end points of are and . Then we assume that there exists a constant and a such that
where is independent of . Under these conditions, Assumption 2 holds.
We can now prove the stability result.
We know we can write where and . Define by
Note that . Now we can estimate .
Hence, we have
Here we used (8) twice. In particular, we used . Assuming we have
We therefore have
Hence, we obtain using (10)
Thus we arrive at
From this and (10) we get
We can now prove error estimates after we make an assumption more stringent than Assumption 2.
Suppose that the end points of are and . Then we assume that there exists a constant such that
where is independent of .
We let . Then we see that
where , and .
Note that using integration by parts we have
First consider then we have
However, we have
Therefore, we get
Now consider .
Here we used (6).
Now lets consider . If we let then
Now let we then have
Let , with end points and . Then, we have . Hence, using Assumption 3 we get
We then obtain the following estimate, after summing over all edges ,
We get the following inequality after using approximation properties of the Lagrange interpolant:
Therefore, we have
Combining the above results we get
The result now follows from Theorem 1.
One feature of Nitsche’s method, that is preserved with BDT, is that one uses the full space of piecewise polynomials without restriction at the boundary. The modification of