A Crank-Nicolson type minimization scheme
for a hyperbolic free boundary problem

In this paper we treat a variational problem related to the following hyperbolic free boundary problem:

We observe that $u$ satisfies the wave equation where $u>0$, and that $u$ is harmonic when $u<0.$ It can be formally shown that, in the energy-preserving regime, solutions fulfill the following free boundary condition:

This kind of problem is a natural prototype for explaining phenomena involving oscillations in the presence of an obstacle, e.g., an elastic string hitting a desk or soap bubbles moving atop water. However, due to the lack of rigorous mathematical results for hyperbolic free boundary problems, numerical studies of this problem are very limited. On the other hand, interesting results have been already obtained for first-order hyperbolic free boundary problems, see, e.g., [3] for the analysis of a model of tumor growth, or [6] for the analysis of wave-structure interactions.

In this paper, we propose a numerical scheme with good energy conservation properties and provide its theoretical background at least in the case of space dimension 1. We remark that similar types of problems have been treated in [8], [13], [7] [5], and that the recent paper [11] has established a precise mathematical formulation. These papers revealed that the discrete Morse flow (also known as minimizing movements), a variational method based on time-discretized functionals, is an effective tool not only for problems of elliptic and parabolic type but also in the hyperbolic setting.

We now briefly review those previous results. Kikuchi and Omata [8] studied the problem in the one-dimensional domain $\Omega=(0,\infty)$. They showed the existence and uniqueness of a strong solution $u\in C^{2}(\Omega\times(0,\infty)\cap\{u>0\})$, the regularity of its free boundary $\partial\{u>0\}$ and the well-posedness of the problem under suitable compatibility conditions. Yoshiuchi et al. [13] addressed a similar problem to Problem 1 including a damping term $\alpha u_{t}$. Using the discrete Morse flow, they derived an energy estimate for approximate solutions, and provided numerical results. The following problem, stated here without initial and boundary conditions, has been treated by Kikuchi in [7]:

He constructed a weak solution to this problem using a minimizing method in the spirit of the discrete Morse flow. Moreover, two of the present authors, Ginder and Svadlenka, dealt with a hyperbolic free boundary problem under a volume constraint in [5]. They constructed a weak solution in the one dimensional setting, again using the discrete Morse flow. The analysis of hyperbolic obstacle problems based on this variational framework of discrete Morse flow was significantly extended to nonlocal problems (with fractional Laplacian) and to semilinear problems in a vector-valued setting by the group of Bonafini, Novaga and Orlandi in the papers [1, 2].

The technique used in this paper is similar to the discrete Morse flow used in the above papers, but it has some distinct differences. We now explain our approach, together with the organization of the paper. The first step is to consider the minimization, within the set

of the following time-discretized, Crank-Nicolson type functional:

where $\mathcal{S}_{m}(u)$ denotes the set $\{u>0\}\cup\{u_{m-1}>0\}\cup\{u_{m-2}>0\}$. In the above, $u_{m}$ represents an approximation of the solution to the original problem at a fixed time $t_{m}=mh$, where $m\in{\mathbb{N}}$, and $h>0$ denotes the time step, obtained by dividing the time interval $(0,T)$ into $M$ equal parts. The difference with the previous approaches (e.g., [13], [5]) appears in the gradient term. As will be shown, our gradient term yields the energy conservation property even in the time discretized setting when there is no free boundary (i.e., when the restriction to the set $\mathcal{S}_{m}(u)$ is omitted). This is a significant improvement to the previous approaches, where energy conservation does not hold and leads to inaccurate numerical solutions. We believe that combination of the new energy-preserving scheme with the more advanced analytical techniques developed in [1, 2] will lead to a powerful framework. The details of the new variational scheme are presented in Section 2.

The sequence $\{u_{m}\}$ is constructed using the initial conditions $u_{0}\in\mathcal{K}$ and $v_{0}\in H^{1}_{0}(\Omega)$. In particular, a forward differencing is used to define $u_{1}:=u_{0}+hv_{0}\in{\mathcal{K}}$. Then, for any integer $m\geq 2$, a minimizer $\widetilde{u}_{m}$ of $J_{m}$ exists and has the subsolution property (Theorem 3.1, Proposition 3.2 in Section 3), which can be shown by a standard argument, such as in [13]. We then define the cut-off minimizer $u_{m}:=\max\{\widetilde{u}_{m},0\}$ and find that the minimizers $u_{m}$ satisfy an energy estimate (Theorem 3.3). In particular, for any integer $k\geq 1$, we have

where $C$ is independent of $k$ and $h$. The regularity of $u_{m}$ can then be obtained and this leads to the first variation formula (Proposition 3.7). The next step is to define approximate weak solutions of the free boundary problem, as well as a weak solution of Problem 1. This is done in Section 4. In the case of a 1-dimensional domain, the energy estimate and the embedding theorem allow us to pass the time step size $h$ to zero in the equation of approximate weak solutions to obtain a weak solution of Problem 1.

In Section 4, we also formally discuss the role of free boundary condition on a more general setting for this problem, namely a hyperbolic free boundary problem with an adhesion term. When energy is not conserved, we cannot calculate the first variation in the usual way. In this case, we adopt a formal approach and derive the equation from a measure theoretic point of view. We conclude by presenting numerical results and their analysis in Section 5, emphasizing the wide applicability of the proposed numerical method.

In this section, we derive the energy preserving property of the functional $J_{m}$ in (3), which has not been achieved in previous research. To this end, let us consider the following modified functional:

on the set $\mathcal{K}:=\{u\in H^{1}(\Omega)\,;\,u=u_{0}\,\hbox{on}\,\partial\Omega\}$. This functional can be regarded as the no-free boundary version of $J_{m}$. We note that a unique minimizer exists for each $I_{m}$ whenever $I_{m}(u_{0})<\infty$ since the functional is convex and lower semicontinuous.

For any integer $m\geq 2$, we introduce the following functional:

where $\mathcal{S}_{m}(u):=\{u>0\}\cup\{u_{m-1}>0\}\cup\{u_{m-2}>0\}$.

We determine a sequence of functions $\{u_{m}\}$ iteratively by taking $u_{0}\in{\mathcal{K}}$ and $u_{1}=u_{0}+hv_{0}\in{\mathcal{K}}$, defining $\widetilde{u}_{m}$ as a minimizer of $J_{m}$ in $\mathcal{K}$, and setting $u_{m}:=\max\{\widetilde{u}_{m},0\}$.

We now study the existence and regularity of minimizers which guarantees the possibility of applying the first variation formula to $J_{m}$.

Proof. Given $u_{m-1},u_{m-2}$, we show the existence of $\widetilde{u}_{m}$. Since the infimum of $J_{m}$ is non-negative, we have only to show the lower semi-continuity of $J_{m}$. Take any minimizing sequence $\{u^{j}\}\subset\mathcal{K}$ such that $J_{m}(u^{j})\to\inf_{u\in\mathcal{K}}J_{m}(u)$ as $j\to\infty$. Since the sequence $\{u^{j}-u_{0}\}\subset H^{1}_{0}(\Omega)$ is bounded in $H^{1}(\Omega)$, there exist $\widetilde{u}\in H^{1}_{0}(\Omega)$ and $\gamma\in L^{p}(\Omega)$, $p\in[1,\infty)$ such that, up to extracting a subsequence,

Since $0\leq\gamma\leq 1$ a.e. on $\Omega$, and $\gamma=1$ a.e. on $\{u>0\}$, where $u:=\widetilde{u}+u_{0}\in\mathcal{K}$,

where the second inequality follows from (3).

The minimizers of $J_{m}$ have the following subsolution property.

Proof. Fixing $\zeta\in C_{0}^{\infty}(\Omega)$ with $\zeta\geq 0$, and $\varepsilon>0$, we have

Noting that,

we continue the estimate as

Dividing by $\varepsilon$, letting $\varepsilon$ decrease to zero from above, and applying a density argument concludes the proof.

We now derive an energy estimate satisfied by the minimizers of $J_{m}$.

The following theorem is obtained by a standard argument from elliptic regularity theory. For the sake of completeness, we shall briefly demonstrate it.

To prove this, we prepare two lemmas.

Proof. We use mathematical induction for $m\geq 2$. For $m=2$, setting $\psi_{\delta}(u):=u-\delta(u+u_{0}-k)^{+}\in\mathcal{K}$, where $u:=\widetilde{u}_{2}$, $(u+u_{0}-k)^{+}:=\max\{u+u_{0}-k,0\}$, $\delta>0$, $k\geq\max\{2\max_{\partial\Omega}u_{0},1\}$, we calculate the quantity $J_{2}(\psi_{\delta}(u))-J_{2}(u)$, which is non-negative by the minimality of $u$. Noting that $S_{m}(\psi_{\delta}(u))\subset S_{m}(u)$, we have

Dividing by $\delta$, letting $\delta\to 0+$, and setting $A_{k}:=\{u+u_{0}>k\}$, we get

where we have used Young’s inequality at the last line. Since $k\geq 1$, we get

Therefore, by [9, Theorem 2.5.1], we find that $u+u_{0}\in L^{\infty}(\Omega)$, and hence $u=\widetilde{u}_{2}\in L^{\infty}(\Omega)$.

Next, we assume that $\widetilde{u}_{k}\in L^{\infty}(\Omega)$ for all $k=2,...,m-1$. Since $u_{k}=\max\{\widetilde{u}_{k},0\}$ belongs to $L^{\infty}(\Omega)$ for all $k=2,...,m-1$, by repeating the above argument with $\widetilde{u}_{2},u_{1},u_{0}$ replaced by $\widetilde{u}_{m},u_{m-1},u_{m-2}$, respectively, we get $\widetilde{u}_{m}+u_{m-2}\in L^{\infty}(\Omega)$. Therefore, $\widetilde{u}_{m}\in L^{\infty}(\Omega)$.

The previous result implies that there exists $\mu>0$, which depends only on $\Omega,u_{0},u_{1},h$ but not on $m$, such that $\sup_{\Omega}|\widetilde{u}_{m}|\leq\mu$.

Proof. For fixed $m\geq 2$, first we show the statement for $U=\widetilde{u}_{m}+u_{m-2}$. We set $\zeta=\eta^{2}\max\{u+u_{m-2}-k,0\}$ in Proposition 3.2, where $u:=\widetilde{u}_{m}$, $k$ is a real number with $k\geq\max_{B_{r}}(u+u_{m-2})-d$, $\eta$ is smooth function with $\hbox{spt}\,\eta\subset B_{r}$, $0\leq\eta\leq 1$, $\eta\equiv 1$ on $B_{s}$, $|\nabla\eta|\leq 2/(r-s)$ in $B_{r}\setminus B_{s}$, and $s=r-\sigma r\in(0,r)$, $\sigma\in(0,1)$. Then, using the boundedness of $u,u_{m-1},u_{m-2}$, we get