[guided]The delicate point is that $\mathcal H$ may be unbounded, so differentiability cannot be proved by treating $\mathcal H$ as an ordinary bounded matrix. We use the Spectral Theorem for Self-Adjoint Operators to reduce the difference quotient to scalar multiplication. Since $\psi_0\in D(\mathcal H)$, the spectral domain characterization obtained from the unbounded functional calculus for $\lambda\mapsto\lambda$ gives
\begin{align*}
\int_{\mathbb R}\lambda^2\,d\mu_{\psi_0}(\lambda)<\infty.
\end{align*}
For $h\ne0$, the spectral isometry identity gives
\begin{align*}
\left\|\frac{U(h)\psi_0-\psi_0}{h}+i\mathcal H\psi_0\right\|_H^2=\int_{\mathbb R}\left|\frac{e^{-ih\lambda}-1}{h}+i\lambda\right|^2\,d\mu_{\psi_0}(\lambda).
\end{align*}
For each fixed $\lambda\in\mathbb R$, the scalar derivative of $h\mapsto e^{-ih\lambda}$ at $0$ is $-i\lambda$, so
\begin{align*}
\frac{e^{-ih\lambda}-1}{h}\to -i\lambda
\end{align*}
as $h\to0$. To justify passing the limit through the spectral integral, we need a dominating function independent of $h$. The mean value estimate for $r\mapsto e^{-ir\lambda}$ gives
\begin{align*}
\left|\frac{e^{-ih\lambda}-1}{h}\right|\le |\lambda|.
\end{align*}
Hence
\begin{align*}
\left|\frac{e^{-ih\lambda}-1}{h}+i\lambda\right|^2\le4\lambda^2.
\end{align*}
The function $\lambda\mapsto4\lambda^2$ is integrable with respect to $\mu_{\psi_0}$ by the domain characterization, so the [Dominated Convergence Theorem](/theorems/4) gives
\begin{align*}
\lim_{h\to0}\left\|\frac{U(h)\psi_0-\psi_0}{h}+i\mathcal H\psi_0\right\|_H=0.
\end{align*}
Thus the orbit is differentiable at $0$ with derivative $-i\mathcal H\psi_0$.
We also need to know that the orbit remains inside $D(\mathcal H)$. This is where the modulus-one property of the multiplier is used a second time. For every Borel set $A\subset\mathbb R$, the projection $E(A)$ is the functional-calculus operator associated to the indicator function $\mathbb{1}_A:\mathbb R\to\{0,1\}$. Since the functional calculus is multiplicative, $E(A)U(t)$ is the functional-calculus operator associated to $\mathbb{1}_A(\lambda)e^{-it\lambda}$. Applying the spectral isometry identity to this bounded Borel multiplier gives
\begin{align*}
\mu_{U(t)\psi_0}(A)=\|E(A)U(t)\psi_0\|_H^2=\int_{\mathbb R}\mathbb{1}_A(\lambda)|e^{-it\lambda}|^2\,d\mu_{\psi_0}(\lambda).
\end{align*}
Because $|e^{-it\lambda}|=1$ for every $\lambda\in\mathbb R$, this becomes
\begin{align*}
\mu_{U(t)\psi_0}(A)=\int_A 1\,d\mu_{\psi_0}(\lambda)=\mu_{\psi_0}(A).
\end{align*}
Therefore
\begin{align*}
\int_{\mathbb R}\lambda^2\,d\mu_{U(t)\psi_0}(\lambda)=\int_{\mathbb R}\lambda^2\,d\mu_{\psi_0}(\lambda)<\infty,
\end{align*}
so $U(t)\psi_0\in D(\mathcal H)$. On this domain, multiplication by $\lambda$ and by $e^{-it\lambda}$ commute in the spectral representation, giving
\begin{align*}
\mathcal H U(t)\psi_0=U(t)\mathcal H\psi_0.
\end{align*}
Finally, the group law gives
\begin{align*}
\frac{U(t+h)\psi_0-U(t)\psi_0}{h}=U(t)\frac{U(h)\psi_0-\psi_0}{h}.
\end{align*}
Since $U(t)$ is bounded, it may be passed through the norm limit. Hence
\begin{align*}
\psi'(t)=U(t)(-i\mathcal H\psi_0)=-i\mathcal H U(t)\psi_0.
\end{align*}
Equivalently, $i\psi'(t)=\mathcal H\psi(t)$ for every $t\in\mathbb R$, and $\psi(0)=\psi_0$ follows from $U(0)=I_H$.[/guided]