[proofplan]
We rewrite the [wave equation](/page/Wave%20Equation) as a first-order evolution equation on the Dirichlet energy space. For smooth compatible data and smooth forcing, differentiating the proposed formula shows that it solves the first-order inhomogeneous equation with the correct initial value, so uniqueness gives the formula. For general energy data and $L^1((0,T);L^2(U))$ forcing, we approximate by smooth data and pass to the limit using the energy estimate, strong continuity of $S(t)$, and convergence of Bochner integrals.
[/proofplan]
[step:Write the wave equation as a first-order evolution on the energy space]
Define the $c$-energy norm on $\mathcal{H}(U)$ by
\begin{align*}
\|(v,w)\|_{\mathcal{H}_c}:=\left(c^2\|\nabla v\|_{L^2(U)}^2+\|w\|_{L^2(U)}^2\right)^{1/2}
\end{align*}
for $(v,w)\in H^1_0(U)\times L^2(U)$, where $c>0$ is the wave speed from the theorem statement. Since $c>0$, this norm is equivalent to the product norm on $\mathcal{H}(U)$; explicitly,
\begin{align*}
\min\{c,1\}\|(v,w)\|_{H^1_0(U)\times L^2(U)}\leq \|(v,w)\|_{\mathcal{H}_c}\leq \max\{c,1\}\|(v,w)\|_{H^1_0(U)\times L^2(U)}
\end{align*}
when the $H^1_0(U)$ norm is taken to be $\|\nabla v\|_{L^2(U)}$. Define the operator
\begin{align*}
A:D(A)\subset \mathcal{H}(U)\to \mathcal{H}(U)
\end{align*}
by
\begin{align*}
D(A):=\bigl(H^2(U)\cap H^1_0(U)\bigr)\times H^1_0(U),
\end{align*}
and
\begin{align*}
A(v,w):=(w,c^2\Delta v).
\end{align*}
For a sufficiently regular solution
\begin{align*}
u:[0,T]\to H^1_0(U),
\end{align*}
define its phase-space trajectory
\begin{align*}
Y:[0,T]\to \mathcal{H}(U), \qquad t\mapsto (u(t),\partial_t u(t)).
\end{align*}
Also define the forcing curve
\begin{align*}
F:(0,T)\to \mathcal{H}(U), \qquad s\mapsto (0,f(s)).
\end{align*}
Then the equation $\partial_t^2u-c^2\Delta u=f$ is equivalent to
\begin{align*}
Y'(t)=AY(t)+F(t)
\end{align*}
for those times $t$ at which the classical derivatives exist.
Equip $\mathcal{H}(U)$ with the [inner product](/page/Inner%20Product) associated to $\|\cdot\|_{\mathcal{H}_c}$. We verify the generator statement directly. For $(v,w),(\phi,\psi)\in D(A)$, [integration by parts](/theorems/210) on $U$ is valid because $v,\phi\in H^2(U)\cap H^1_0(U)$ and the traces of $v$ and $\phi$ vanish on $\partial U$. It gives
\begin{align*}
(A(v,w),(\phi,\psi))_{\mathcal{H}_c}= -((v,w),A(\phi,\psi))_{\mathcal{H}_c}.
\end{align*}
Thus $A$ is skew-symmetric. To verify maximality, fix $\lambda\in\mathbb{R}\setminus\{0\}$ and fix an arbitrary pair $(g,h)\in H^1_0(U)\times L^2(U)=\mathcal{H}(U)$. The resolvent equation $(\lambda I-A)(v,w)=(g,h)$ for an unknown pair $(v,w)\in D(A)$ is equivalent to $w=\lambda v-g$ and
\begin{align*}
(\lambda^2-c^2\Delta)v=\lambda g+h.
\end{align*}
The [bilinear form](/page/Bilinear%20Form)
\begin{align*}
B_\lambda[v,\phi]:=\lambda^2\int_U v\phi\,d\mathcal{L}^n(x)+c^2\int_U \nabla v\cdot\nabla \phi\,d\mathcal{L}^n(x)
\end{align*}
on $H^1_0(U)$ is continuous and coercive because $\lambda^2>0$ and $c^2>0$. The [Lax-Milgram theorem](/theorems/91) gives a unique weak solution $v\in H^1_0(U)$, and elliptic regularity for the smooth bounded Dirichlet domain upgrades $v$ to $H^2(U)\cap H^1_0(U)$. Then $w=\lambda v-g\in H^1_0(U)$, so $\lambda I-A$ is onto $\mathcal{H}(U)$. The same argument for $-\lambda I-A$ gives maximality on both signs. Hence $A$ is skew-adjoint on the energy space by the standard maximal skew-symmetric operator criterion, so Stone's theorem for skew-adjoint generators generates a strongly continuous unitary group, which is precisely the homogeneous wave solution group $S(t)$. Therefore, for each $Z\in D(A)$,
\begin{align*}
\frac{d}{dt}S(t)Z=AS(t)Z=S(t)AZ.
\end{align*}
Unitarity gives preservation of the energy norm:
\begin{align*}
\|S(t)Z\|_{\mathcal{H}_c}=\|Z\|_{\mathcal{H}_c}
\end{align*}
for every $Z\in\mathcal{H}(U)$ and every $t\in\mathbb{R}$.
[/step]
[step:Verify the formula for smooth compatible data]
Assume first that
\begin{align*}
(u_0,u_1)\in D(A)
\end{align*}
and that the forcing curve
\begin{align*}
F:[0,T]\to \mathcal{H}(U), \qquad s\mapsto (0,f(s)),
\end{align*}
belongs to $C^1([0,T];\mathcal{H}(U))$ and satisfies $F(s)\in D(A)$ for every $s\in[0,T]$, with
\begin{align*}
AF:[0,T]\to\mathcal{H}(U), \qquad s\mapsto A(F(s)),
\end{align*}
continuous. These hypotheses imply that $s\mapsto S(t-s)F(s)$ is Bochner integrable in $D(A)$ with the graph norm of $A$. We use the Banach-space Leibniz rule for parameter-dependent Bochner integrals: if the integrand and its parameter derivative are continuous as Banach-space-valued maps, then differentiation under the integral sign is valid and the moving upper endpoint contributes the boundary value. Thus differentiating $S(t-s)F(s)$ in $t$ and moving $A$ through the Bochner integral are justified. Define
\begin{align*}
Z:[0,T]\to \mathcal{H}(U)
\end{align*}
by
\begin{align*}
Z(t):=S(t)(u_0,u_1)+\int_0^t S(t-s)F(s)\,d\mathcal{L}^1(s).
\end{align*}
At $t=0$, the Bochner integral is over the empty interval, hence
\begin{align*}
Z(0)=(u_0,u_1).
\end{align*}
Differentiating the first term gives
\begin{align*}
\frac{d}{dt}S(t)(u_0,u_1)=AS(t)(u_0,u_1).
\end{align*}
For the integral term, the Leibniz rule for parameter-dependent Bochner integrals gives
\begin{align*}
\frac{d}{dt}\int_0^t S(t-s)F(s)\,d\mathcal{L}^1(s)
=
F(t)+\int_0^t AS(t-s)F(s)\,d\mathcal{L}^1(s).
\end{align*}
Since $A$ commutes with $S(t-s)$ on the regular data under consideration, this equals
\begin{align*}
F(t)+A\int_0^t S(t-s)F(s)\,d\mathcal{L}^1(s).
\end{align*}
Therefore
\begin{align*}
Z'(t)=AZ(t)+F(t).
\end{align*}
Writing $Z(t)=(z(t),\partial_t z(t))$, this first-order equation is exactly
\begin{align*}
\partial_t^2z-c^2\Delta z=f,
\end{align*}
with homogeneous Dirichlet boundary condition and initial data
\begin{align*}
(z(0),\partial_tz(0))=(u_0,u_1).
\end{align*}
By the standard uniqueness theorem for the smooth homogeneous Dirichlet wave problem, applied to the difference of two smooth solutions with zero initial data, zero forcing, and homogeneous Dirichlet trace, $Z(t)=(u(t),\partial_tu(t))$ for every $t\in[0,T]$.
[guided]
The goal is to check that the right-hand side of the proposed formula is not merely plausible, but actually satisfies the same differential equation and initial condition as the forced wave solution.
Assume the data are smooth enough that all differentiations are classical in the following precise sense: $(u_0,u_1)\in D(A)$, the curve $F:[0,T]\to\mathcal{H}(U)$ belongs to $C^1([0,T];\mathcal{H}(U))$, $F(s)\in D(A)$ for every $s\in[0,T]$, and $s\mapsto A(F(s))$ is continuous as a map into $\mathcal{H}(U)$. We first recall why the operator-theoretic input is available. With respect to the $c$-energy inner product, [integration by parts](/theorems/2098) on $H^2(U)\cap H^1_0(U)$ and the homogeneous Dirichlet trace show that $A$ is skew-symmetric. Maximality follows by solving $(\lambda I-A)(v,w)=(g,h)$: the equation reduces to $w=\lambda v-g$ and the coercive Dirichlet elliptic problem $(\lambda^2-c^2\Delta)v=\lambda g+h$ on $H^1_0(U)$, whose weak solution is in $H^2(U)\cap H^1_0(U)$ by elliptic regularity on the smooth bounded domain. Thus $A$ is skew-adjoint, and Stone's theorem gives the strongly continuous unitary group $S(t)$ generated by $A$. Hence, for every $Z\in D(A)$, the map $t\mapsto S(t)Z$ is differentiable in $\mathcal{H}(U)$ and satisfies $\frac{d}{dt}S(t)Z=AS(t)Z=S(t)AZ$. Since $F(s)\in D(A)$ and $s\mapsto AF(s)$ is continuous, these assumptions ensure that $S(t-s)F(s)$ is differentiable in $t$ with derivative $AS(t-s)F(s)$ and that $A$ can be passed through the Bochner integral. The precise analytic tool is the Banach-space Leibniz rule for parameter-dependent Bochner integrals: continuity of the integrand and of its parameter derivative permits differentiating under the integral sign, and the variable upper limit contributes the value of the integrand at the endpoint. Define
\begin{align*}
F:[0,T]\to \mathcal{H}(U), \qquad s\mapsto (0,f(s)),
\end{align*}
and define
\begin{align*}
Z:[0,T]\to \mathcal{H}(U)
\end{align*}
by
\begin{align*}
Z(t):=S(t)(u_0,u_1)+\int_0^t S(t-s)F(s)\,d\mathcal{L}^1(s).
\end{align*}
The initial condition is immediate from the definition of the Bochner integral over a degenerate interval:
\begin{align*}
Z(0)=S(0)(u_0,u_1)=(u_0,u_1).
\end{align*}
Now differentiate. The homogeneous group satisfies
\begin{align*}
\frac{d}{dt}S(t)(u_0,u_1)=AS(t)(u_0,u_1).
\end{align*}
For the inhomogeneous term, two effects occur when $t$ changes: the upper limit changes, producing the boundary contribution $F(t)$, and the kernel $S(t-s)$ changes, producing the generator $A$. Thus
\begin{align*}
\frac{d}{dt}\int_0^t S(t-s)F(s)\,d\mathcal{L}^1(s)
=
F(t)+\int_0^t AS(t-s)F(s)\,d\mathcal{L}^1(s).
\end{align*}
On the present smooth class, $A$ may be passed outside the integral, so
\begin{align*}
\frac{d}{dt}\int_0^t S(t-s)F(s)\,d\mathcal{L}^1(s)
=
F(t)+A\int_0^t S(t-s)F(s)\,d\mathcal{L}^1(s).
\end{align*}
Combining this with the derivative of the homogeneous term gives
\begin{align*}
Z'(t)=AZ(t)+F(t).
\end{align*}
Finally, if $Z(t)=(z(t),q(t))$, then the first component of $Z'=AZ+F$ says $q(t)=\partial_tz(t)$, while the second component says
\begin{align*}
\partial_t q(t)=c^2\Delta z(t)+f(t).
\end{align*}
Hence
\begin{align*}
\partial_t^2z-c^2\Delta z=f.
\end{align*}
The boundary condition is homogeneous Dirichlet in the trace sense because the first component $z(t)$ lies in $H^1_0(U)$ for each $t\in[0,T]$, and for the present smooth compatible data this trace statement agrees with the classical boundary condition $z(t)|_{\partial U}=0$. Since $Z(0)=(u_0,u_1)$, the standard uniqueness theorem for the smooth homogeneous Dirichlet wave problem, applied to the difference of $Z$ and the smooth forced solution with zero initial data, zero forcing, and zero boundary trace, identifies $Z$ with that solution.
[/guided]
[/step]
[step:Approximate the energy data and the forcing]
Let $C_c^\infty(U)$ denote the [vector space](/page/Vector%20Space) of smooth functions $\varphi:U\to\mathbb{R}$ whose support is compact and contained in $U$. Because $U$ is a bounded smooth domain, $C_c^\infty(U)$ is dense in $H^1_0(U)$ and $H^1_0(U)$ is dense in $L^2(U)$. Hence $D(A)=\bigl(H^2(U)\cap H^1_0(U)\bigr)\times H^1_0(U)$ is dense in $\mathcal{H}(U)$ for the energy topology. Choose a sequence
\begin{align*}
(u_{0,k},u_{1,k})\in D(A), \qquad k\in\mathbb{N},
\end{align*}
such that
\begin{align*}
(u_{0,k},u_{1,k})\to (u_0,u_1)
\end{align*}
in $\mathcal{H}(U)$. Since $L^2(U)$ is a separable [Hilbert space](/page/Hilbert%20Space) and $H^1_0(U)$ is dense in $L^2(U)$, we construct the forcing approximation as follows. First approximate $f$ in $L^1((0,T);L^2(U))$ by a [simple function](/page/Simple%20Function) whose values lie in $H^1_0(U)$. Extend that simple function by zero outside $(0,T)$, convolve in the time variable with a one-dimensional smooth compactly supported mollifier, and then restrict back to $[0,T]$. A final polynomial-in-time approximation, if needed, gives a sequence
\begin{align*}
f_k:[0,T]\to L^2(U), \qquad k\in\mathbb{N},
\end{align*}
with $f_k\in C^1([0,T];H^1_0(U))$ such that
\begin{align*}
f_k\to f
\end{align*}
in $L^1((0,T);L^2(U))$. Define
\begin{align*}
F_k:(0,T)\to \mathcal{H}(U), \qquad s\mapsto (0,f_k(s)).
\end{align*}
Then $F_k(s)\in D(A)$ for every $s\in[0,T]$, because $0\in H^2(U)\cap H^1_0(U)$ and $f_k(s)\in H^1_0(U)$. Moreover
\begin{align*}
AF_k:[0,T]\to\mathcal{H}(U), \qquad s\mapsto (f_k(s),0),
\end{align*}
is continuous because $f_k\in C^1([0,T];H^1_0(U))\subset C([0,T];L^2(U))$.
For each $k\in\mathbb{N}$, classical well-posedness for the smooth homogeneous Dirichlet wave problem applies because $(u_{0,k},u_{1,k})\in D(A)$ and $f_k\in C^1([0,T];H^1_0(U))$ satisfy the compatibility needed for the operator-domain formulation. This is the standard smooth well-posedness theorem for the Dirichlet wave equation in the generator domain, obtained from the strongly continuous group generated by $A$ and classical differentiation in $D(A)$. Let
\begin{align*}
u_k:[0,T]\to H^1_0(U)
\end{align*}
be this corresponding smooth forced solution, with
\begin{align*}
\partial_tu_k:[0,T]\to L^2(U).
\end{align*}
By the smooth case already proved,
\begin{align*}
(u_k(t),\partial_tu_k(t))
=
S(t)(u_{0,k},u_{1,k})
+
\int_0^t S(t-s)F_k(s)\,d\mathcal{L}^1(s)
\end{align*}
for every $t\in[0,T]$.
[/step]
[step:Pass to the limit in the homogeneous and inhomogeneous terms]
Fix $t\in[0,T]$. Since $S(t):\mathcal{H}(U)\to\mathcal{H}(U)$ is bounded and strongly continuous,
\begin{align*}
S(t)(u_{0,k},u_{1,k})\to S(t)(u_0,u_1)
\end{align*}
in $\mathcal{H}(U)$.
For the limiting forcing term, the map $s\mapsto F(s)=(0,f(s))$ is strongly measurable as a map into $\mathcal{H}(U)$ because $f$ is strongly measurable as an $L^2(U)$-valued function. Since $r\mapsto S(r)$ is strongly continuous and $S(r)$ is an isometry in $\|\cdot\|_{\mathcal{H}_c}$, the map $s\mapsto S(t-s)F(s)$ is strongly measurable and satisfies
\begin{align*}
\|S(t-s)F(s)\|_{\mathcal{H}_c}=\|f(s)\|_{L^2(U)}
\end{align*}
for $\mathcal{L}^1$-a.e. $s\in(0,t)$. Thus $s\mapsto S(t-s)F(s)$ is Bochner integrable on $(0,t)$ because $f\in L^1((0,T);L^2(U))$. For convergence of the forcing terms, use the same energy isometry of $S(t-s)$. For each $s\in(0,t)$,
\begin{align*}
\|S(t-s)(F_k(s)-F(s))\|_{\mathcal{H}_c}
=
\|F_k(s)-F(s)\|_{\mathcal{H}_c}
=
\|f_k(s)-f(s)\|_{L^2(U)}.
\end{align*}
Therefore
\begin{align*}
\left\|
\int_0^t S(t-s)(F_k(s)-F(s))\,d\mathcal{L}^1(s)
\right\|_{\mathcal{H}_c}
\le
\int_0^t \|f_k(s)-f(s)\|_{L^2(U)}\,d\mathcal{L}^1(s).
\end{align*}
Since $f_k\to f$ in $L^1((0,T);L^2(U))$, the right-hand side tends to $0$. Hence
\begin{align*}
\int_0^t S(t-s)F_k(s)\,d\mathcal{L}^1(s)
\to
\int_0^t S(t-s)F(s)\,d\mathcal{L}^1(s)
\end{align*}
in $\mathcal{H}(U)$.
[/step]
[step:Identify the limit with the energy solution]
Let
\begin{align*}
Y_k:[0,T]\to\mathcal{H}(U), \qquad t\mapsto (u_k(t),\partial_tu_k(t)),
\end{align*}
and
\begin{align*}
Y:[0,T]\to\mathcal{H}(U), \qquad t\mapsto (u(t),\partial_tu(t)).
\end{align*}
The difference $Y_k-Y$ is the energy solution of the forced Dirichlet wave equation with initial datum $(u_{0,k},u_{1,k})-(u_0,u_1)$ and forcing $f_k-f\in L^1((0,T);L^2(U))$. The standard Dirichlet wave energy estimate for energy solutions applies because both solutions have homogeneous Dirichlet boundary condition, the initial difference lies in $\mathcal{H}(U)$, and the forcing difference lies in $L^1((0,T);L^2(U))$. Here an energy solution means a curve in $C([0,T];H^1_0(U))$ with time derivative in $C([0,T];L^2(U))$ satisfying the wave equation in the usual weak distributional sense. The estimate is obtained by differentiating the Dirichlet energy for smooth differences, using the boundary condition to remove the boundary term, applying Cauchy-Schwarz in $L^2(U)$ to the forcing term, and then passing to energy solutions by density through the energy-solution well-posedness theorem. It gives, for every $t\in[0,T]$,
\begin{align*}
\|Y_k(t)-Y(t)\|_{\mathcal{H}_c}
\le
\|(u_{0,k},u_{1,k})-(u_0,u_1)\|_{\mathcal{H}_c}
+
\int_0^t \|f_k(s)-f(s)\|_{L^2(U)}\,d\mathcal{L}^1(s).
\end{align*}
Both terms on the right tend to $0$, so
\begin{align*}
Y_k(t)\to Y(t)
\end{align*}
in $\mathcal{H}(U)$.
Taking limits in the identity for $Y_k(t)$ from the smooth case and using the convergence of the homogeneous and inhomogeneous terms proved above yields
\begin{align*}
Y(t)
=
S(t)(u_0,u_1)
+
\int_0^t S(t-s)F(s)\,d\mathcal{L}^1(s).
\end{align*}
Since $F(s)=(0,f(s))$, this is precisely
\begin{align*}
(u(t),\partial_tu(t))
=
S(t)(u_0,u_1)
+
\int_0^t S(t-s)(0,f(s))\,d\mathcal{L}^1(s).
\end{align*}
The identity holds for the fixed $t\in[0,T]$, and since $t$ was arbitrary, it holds for every $t\in[0,T]$.
[guided]
The smooth calculation proves the formula only for regular data, so we still have to justify that the identity survives the passage to energy data and $L^1$ forcing. Define
\begin{align*}
Y_k:[0,T]\to\mathcal{H}(U), \qquad t\mapsto (u_k(t),\partial_tu_k(t)),
\end{align*}
and define
\begin{align*}
Y:[0,T]\to\mathcal{H}(U), \qquad t\mapsto (u(t),\partial_tu(t)).
\end{align*}
The difference $Y_k-Y$ solves the Dirichlet wave equation with initial datum $(u_{0,k},u_{1,k})-(u_0,u_1)$ and forcing $f_k-f$. The Dirichlet energy estimate applies because the boundary condition is homogeneous for both solutions, the initial difference belongs to $\mathcal{H}(U)$, and $f_k-f\in L^1((0,T);L^2(U))$. In this passage, an energy solution is a curve in $C([0,T];H^1_0(U))$ with time derivative in $C([0,T];L^2(U))$ satisfying the wave equation weakly. For smooth differences, the estimate comes from differentiating the Dirichlet energy, cancelling the boundary term by the zero trace condition, and applying Cauchy-Schwarz in $L^2(U)$ to the forcing contribution; the energy-solution version follows by density. Therefore, for every $t\in[0,T]$,
\begin{align*}
\|Y_k(t)-Y(t)\|_{\mathcal{H}_c}
\le
\|(u_{0,k},u_{1,k})-(u_0,u_1)\|_{\mathcal{H}_c}
+
\int_0^t \|f_k(s)-f(s)\|_{L^2(U)}\,d\mathcal{L}^1(s).
\end{align*}
The first term tends to $0$ by the chosen convergence of the initial data, and the integral term tends to $0$ because $f_k\to f$ in $L^1((0,T);L^2(U))$. Hence
\begin{align*}
Y_k(t)\to Y(t)
\end{align*}
in $\mathcal{H}(U)$.
Now we pass to the limit in the formula for $Y_k(t)$. The homogeneous term converges because $S(t):\mathcal{H}(U)\to\mathcal{H}(U)$ is bounded:
\begin{align*}
S(t)(u_{0,k},u_{1,k})\to S(t)(u_0,u_1).
\end{align*}
For the Duhamel term, the energy isometry gives
\begin{align*}
\left\|
\int_0^t S(t-s)(F_k(s)-F(s))\,d\mathcal{L}^1(s)
\right\|_{\mathcal{H}_c}
\le
\int_0^t \|f_k(s)-f(s)\|_{L^2(U)}\,d\mathcal{L}^1(s),
\end{align*}
which tends to $0$. Taking limits in the smooth identity gives
\begin{align*}
Y(t)
=
S(t)(u_0,u_1)
+
\int_0^t S(t-s)F(s)\,d\mathcal{L}^1(s).
\end{align*}
Since $F(s)=(0,f(s))$, this is exactly
\begin{align*}
(u(t),\partial_tu(t))
=
S(t)(u_0,u_1)
+
\int_0^t S(t-s)(0,f(s))\,d\mathcal{L}^1(s).
\end{align*}
Because the argument was made for an arbitrary fixed $t\in[0,T]$, the formula holds for every $t\in[0,T]$.
[/guided]
[/step]