\section{Differentiable Light Transport in 1D}
\label{sec:rte1d}
%
In this section, we focus on 1D light transport in homogeneous non-emissive media.
For such a material with extinction coefficient~$\sigma_t$, scattering coefficient~$\sigma_s$, and the isotropic phase function, the radiative transfer equation (RTE) reduces to
%
\begin{equation}
	\label{eqn:rte1d}
	\omega \left(\frac{\partial}{\partial x}L_{\omega}\right)(x) + \sigma_t L_{\omega}(x) = \frac{\sigma_s}{2} (L_{\omega}(x) + L_{-\omega}(x)),
\end{equation}
%
where $\omega \in \{1, -1\}$, and $L_{\omega}$ denotes the radiance along $\omega$.
%
It can be shown that Eq.~\eqref{eqn:rte1d} has analytical solutions of the form
%
\begin{equation}
	\label{eqn:rte1d_sol}
	L_{\omega}(x) = a_{\omega} \,\E^{-\lambda x}+ b_{\omega} \,\E^{\lambda x},
	%\label{eqn:rte1d_sol_diff}
	%\left(\frac{\partial}{\partial x} L_{\omega}\right)(x) = - a_{\omega} \lambda \,\E^{-\lambda x} + b_{\omega} \lambda \,\E^{\lambda x},
\end{equation}
%
for some $a_{\omega}, b_{\omega} \in \Real$ and $\lambda := \sqrt{\sigma_t (\sigma_t - \sigma_s)}$.
Assuming boundary conditions
%
\begin{equation}
	\label{eqn:rte1d_boundary}
	L_1(x_{-1}) = v_{-1}, \qquad L_{-1}(x_1) = v_1,
\end{equation}
%
for some $x_{-1} < x_1$, then the coefficients $a_{\omega}$ and $b_{\omega}$ in Eq.~\eqref{eqn:rte1d_sol} can be obtained by solving a linear system that includes Eq.~\eqref{eqn:rte1d_boundary} and
%
\begin{align*}
	%L_1(x_{-1}) &= v_{-1},\\
	-\left( \frac{\partial}{\partial x} L_{-1} \right)(x_{-1}) + \sigma_t L_{-1}(x_{-1}) &= \frac{\sigma_s}{2} (L_{-1}(x_{-1}) + L_1(x_{-1})),\\
	\left( \frac{\partial}{\partial x} L_1 \right)(x_1) + \sigma_t L_1(x_1) &= \frac{\sigma_s}{2} (L_{-1}(x_1) + L_1(x_1)).%\\
	%L_{-1}(x_1) &= v_1,
\end{align*}
%
\iffalse
which yields
%
\begin{equation}
	\label{eqn:rte1d_sol_coeffs}
	a_{\omega} = \frac{\E^{\lambda (x_{-1} + x_1)} \left(-\omega\,\E^{\lambda x_{-\omega}} \, v_{\omega} \, \sigma_s + \E^{\lambda x_{\omega}} \, v_{-\omega} \, c_{-\omega}\right)}{d},
	\quad
	b_{\omega} = \frac{\omega \, \E^{\lambda x_{\omega}} \, v_{\omega} \, \sigma_s + \E^{\lambda x_{-\omega}} \, v_{-\omega} \, c_{\omega}}{d},
\end{equation}
%
\begin{comment}
\begin{equation}
	\begin{aligned}
		a_0 &= \frac{\E^{\lambda (x_0 + x_1)} \left(-\E^{\lambda x_0} \, v_1 \, \sigma_s + \E^{\lambda x_1} \, v_0 \, c_0\right)}{d}, &
		b_0 &= \frac{\E^{\lambda x_1} \, v_1 \, \sigma_s + \E^{\lambda x_0} \, v_0 \, c_1}{d};\\
		a_1 &= \frac{\E^{\lambda (x_0 + x_1)} \left(\E^{\lambda x_1} \, v_0 \, \sigma_s + \E^{\lambda x_0} \, v_1 \, c_1 \right)}{d}, &
		b_1 &= \frac{-\E^{\lambda x_0} \, v_0 \, \sigma_s + \E^{\lambda x_1} \, v_1 \, c_0}{d};
	\end{aligned}
\end{equation}
\end{comment}
%
where
%
\begin{equation}
	\label{eqn:rte1d_hvars}
	c_{\omega} := 2\lambda + \omega (\sigma_s - 2\sigma_t), \quad
	d := \E^{2 \lambda x_{-1}} \, c_1 - \E^{2 \lambda x_1} \, c_{-1}.
\end{equation}
\fi

\paragraph{Coupled system}
%
We now consider a more complex problem involving a couple system comprised of two homogeneous media with material parameters $\sigma_t^{(i)}, \sigma_s^{(i)}$ and radiance fields $L_{\omega}^{(i)}$ for $i = 1, 2$.
Geometrically, the first medium is contained in $[0, u]$ while the second in $[u, 1]$ for some $0 < u < 1$.
Additionally, at the interface of $x = u$, we assume the interface conditions
%
\begin{equation}
	\label{eqn:rte1d_coupled_ic}
	\begin{aligned}
		L_{-1}^{(1)}(u) &= F \, L_1^{(1)}(u) + (1 - F) L_{-1}^{(2)}(u) \frac{1}{\eta},\\
		L_1^{(2)}(u) &= F \, L_{-1}^{(2)}(u) + (1 - F) L_1^{(1)}(u) \, \eta,
	\end{aligned}
\end{equation}
%
where $0 \leq F \leq 1$ is the Fresnel reflection term, and $\eta > 0$ captures the difference in refractive indices.
Then, the coefficients $a_{\omega}^{(i)}, b_{\omega}^{(i)}$ can be determined uniquely by solving the linear system given by Eq.~\eqref{eqn:rte1d_coupled_ic} and
%
\begin{align*}
	L_1^{(1)}(0) &= v_{-1},\\
	-\left( \frac{\partial}{\partial x} L_{-1}^{(1)} \right)(0) + \sigma_t^{(1)} L_{-1}^{(1)}(0) &= \frac{\sigma_s^{(1)}}{2} \left(L_{-1}^{(1)}(0) + L_1^{(1)}(0)\right),\\
	\left( \frac{\partial}{\partial x} L_1^{(1)} \right)(u) + \sigma_t^{(1)} L_1^{(1)}(u) &= \frac{\sigma_s^{(1)}}{2} \left(L_{-1}^{(1)}(u) + L_1^{(1)}(u)\right),\\
	-\left( \frac{\partial}{\partial x} L_{-1}^{(2)} \right)(u) + \sigma_t^{(2)} L_{-1}^{(2)}(u) &= \frac{\sigma_s^{(2)}}{2} \left(L_{-1}^{(2)}(u) + L_1^{(2)}(u)\right),\\
	\left( \frac{\partial}{\partial x} L_1^{(2)} \right)(1) + \sigma_t^{(2)} L_1^{(2)}(1) &= \frac{\sigma_s^{(2)}}{2} \left(L_{-1}^{(2)}(1) + L_1^{(2)}(1)\right),\\
	L_{-1}^{(2)}(1) &= v_1.
\end{align*}
%
Figure~\ref{fig:rte1d_sol} demonstrates the solution radiance fields of an example problem.

\begin{figure}[t]
	\centering
	\includegraphics[height=1.8in]{images/rte1d_sol.pdf}
	\caption{\label{fig:rte1d_sol}
		Solution radiance field for $\sigma_s^{(1)} = 0.8$, $\sigma_t^{(1)} = 1.0$, $\sigma_s^{(2)} = 0.5$, $\sigma_t^{(2)} = 1.0$, $v_{-1} = 5.0$, $v_1 = 1.0$, $F = 0.7$, $\eta = 1.5$, and $u = 0.5$.
		The orange and blue curves illustrate $L_1^{(i)}$ and $L_{-1}^{(i)}$ functions, respectively.
	}
\end{figure}
%
\begin{figure}[b]
	\centering
	\includegraphics[width=\columnwidth]{images/rte1d_sol_diff.pdf}
	\caption{\label{fig:rte1d_sol_diff}
		Comparing $\partial L_{-1}^{(1)}(x_0)/\partial u$ calculated directly (blue) and via Eq.~\protect\eqref{eqn:rte1d_sol_intd} (orange) for $x_0 = 0.2$  under the same configuration as Figure~\protect\ref{fig:rte1d_sol}.
		The two agree perfectly.
	}
\end{figure}

\paragraph{Differentiating $L_{\omega}^{(i)}$}
Provided $x_0 \in (0, 1)$, the solution radiance can be differentiated with respect to $u$, the geometric location of the interface.
Without loss of generality, we focus on calculating $\partial L_{-1}^{(1)}(x_0)/\partial u$.

On one hand, since $L_{-1}^{(1)}(x_0)$ can be expressed analytically as a function of $u$, it an be differentiated symbolically.
On the other hand, the IRTE gives that
%
\begin{equation}
	\label{eqn:rte1d_sol_int}
	\begin{split}
		L_{-1}^{(1)}(x_0) =\;& \int_{x_0}^u \E^{-\sigma_t^{(1)}(x - x_0)} \frac{\sigma_s^{(1)}}{2} (L_{-1}^{(1)}(x) + L_1^{(1)}(x)) \,\D x\\
		& + \E^{-\sigma_t^{(1)} (u - x_0)} L_{-1}^{(1)}(u),
	\end{split}
\end{equation}
%
which can be differentiated using \thmref{reynolds}, leading to
%
\begin{equation}
	\label{eqn:rte1d_sol_intd}
	\begin{split}
		\frac{\partial}{\partial u} L_{-1}^{(1)}(x_0) =\;& \int_{x_0}^u \E^{-\sigma_t^{(1)}(x - x_0)} \frac{\sigma_s^{(1)}}{2} \frac{\partial}{\partial u} \left( L_{-1}^{(1)}(x) + L_1^{(1)}(x) \right) \,\D x\\
		& + \E^{-\sigma_t^{(1)}(u - x_0)} \frac{\sigma_s^{(1)}}{2} \left( L_{-1}^{(1)}(u) + L_1^{(1)}(u) \right)\\
		& + \frac{\partial}{\partial u} \left( \E^{-\sigma_t^{(1)} (u - x_0)} L_{-1}^{(1)}(u) \right).
	\end{split}
\end{equation}
%
In the RHS of Eq.~\eqref{eqn:rte1d_sol_intd}, radiances from the same region, namely $L_{-1}^{(1)}(u)$ and $L_1^{(1)}(u)$, are used when evaluating the integrand at $x = u$.
Figure~\ref{fig:rte1d_sol_diff} shows illustrations of derivative functions calculated with both methods. 
