Following 2009 Trangenstein "Numeric Simulations of Hyperbolic Conservation Laws", section 7.4:
Notice that his book coincides with anholonomic orthonormal basis.
So that's what I'm using here. $u_1 ... u_n$ are my coordinate chart parameters. Indexes are in anholonomic non-coordinates.
$\bar{j} =$ the coordinate basis index, ${e^\bar{j}}_j =$ change of basis from coordinate (non-Cartesian) to non-coordinate basis.

Conservation of our field F is a $\pmatrix{ p+1 \\ q}$ in space and time:
$\nabla \cdot F = 0$
with indexes:
$\nabla_\mu \cdot {F^{\mu \bar{i}_1 ... \bar{i}_p}}_{\bar{k}_1 ... \bar{k}_q} = 0$
separating time:
$ \nabla_t \cdot {F^{t \bar{i}_1 ... \bar{i}_p}}_{\bar{k}_1 ... \bar{k}_q} + \nabla_\bar{j} \cdot {F^{\bar{j} \bar{i}_1 ... \bar{i}_p}}_{\bar{k}_1 ... \bar{k}_q} = 0$
assume all ${\Gamma^t}_{\alpha\beta} = 0$ and all ${\Gamma^\alpha}_{t\mu} = 0$
$ \partial_t \cdot {F^{t \bar{i}_1 ... \bar{i}_p}}_{\bar{k}_1 ... \bar{k}_q} + \nabla_\bar{j} \cdot {F^{\bar{j} \bar{i}_1 ... \bar{i}_p}}_{\bar{k}_1 ... \bar{k}_q} = 0$
Let our state ${U^{\bar{i}_1 ... \bar{i}_p}}_{\bar{k}_1 ... \bar{k}_q} = {F^{t \bar{i}_1 ... \bar{i}_p}}_{\bar{k}_1 ... \bar{k}_q}$ be whatever is travelling through time:
$ \partial_t \cdot {U^{\bar{i}_1 ... \bar{i}_p}}_{\bar{k}_1 ... \bar{k}_q} + \nabla_{\bar{j}} \cdot {F^{\bar{j} \bar{i}_1 ... \bar{i}_p}}_{\bar{k}_1 ... \bar{k}_q} = 0$

We will only consider $U$ to be a single $\left( \begin{matrix} p \\ q \end{matrix} \right)$-form. You can repeat for collections of tensors.

The indexes $\bar{i}_1 ... \bar{i}_p, \bar{j}_1 ... \bar{j}_q$ span our chart coordinates $u_1 ... u_n$.

Include the change-of-coordinates:
$ \partial_t \cdot ( {U^{i_1 ... i_p}}_{k_1 ... k_q} \cdot {e^{\bar{i}_1}}_{i_1} \cdot ... \cdot {e^{\bar{i}_p}}_{i_p} \cdot {e_{\bar{k}_1}}^{k_1} \cdot ... \cdot {e_{\bar{k}_q}}^{k_p} ) + \nabla_{\bar{j}} \cdot ( {F^{\bar{j} i_1 ... i_p}}_{k_1 ... k_q} \cdot {e^{\bar{i}_1}}_{i_1} \cdot ... \cdot {e^{\bar{i}_p}}_{i_p} \cdot {e_{\bar{k}_1}}^{k_1} \cdot ... \cdot {e_{\bar{k}_q}}^{k_p} ) = 0$
Assume the grid is constant, so $\partial_t {e_\bar{i}}^i = 0$:
$ \partial_t \cdot ( {U^{i_1 ... i_p}}_{k_1 ... k_q} ) \cdot {e^{\bar{i}_1}}_{i_1} \cdot ... \cdot {e^{\bar{i}_p}}_{i_p} \cdot {e_{\bar{k}_1}}^{k_1} \cdot ... \cdot {e_{\bar{k}_q}}^{k_p} + \nabla_{\bar{j}} \cdot ( {F^{\bar{j} i_1 ... i_p}}_{k_1 ... k_q} \cdot {e^{\bar{i}_1}}_{i_1} \cdot ... \cdot {e^{\bar{i}_p}}_{i_p} \cdot {e_{\bar{k}_1}}^{k_1} \cdot ... \cdot {e_{\bar{k}_q}}^{k_p} ) = 0$
Note that the connection of this $\nabla_{\bar{j}}$ is associated with the coordinate metric and therefore is symmetric.
You could equivalently use the connection $\nabla_j$ associated with the anholonomic components, which is antisymmetric, and then you would not need to worry about rescaling the components (because that information would already be stored in the connection).
That would look like this:
$( \partial_t \cdot ( {U^{i_1 ... i_p}}_{k_1 ... k_q} ) + \nabla_j \cdot {F^{j i_1 ... i_p}}_{k_1 ... k_q} ) \cdot {e^{\bar{i}_1}}_{i_1} \cdot ... \cdot {e^{\bar{i}_p}}_{i_p} \cdot {e_{\bar{k}_1}}^{k_1} \cdot ... \cdot {e_{\bar{k}_q}}^{k_p} = 0$

Is it safe to integrate the normalized-coordinate components across the volume, or should the components be converted to coordinate components inside the integral and then converted back to non-coordinates outside the integral?
$ \partial_t \cdot ( {U^{i_1 ... i_p}}_{k_1 ... k_q} ) + \nabla_j \cdot {F^{j i_1 ... i_p}}_{k_1 ... k_q} = 0$
...?



For $\pmatrix{0 \\ 0}$ case:
For scalar conservation laws you don't have to worry about renormalizing indexes.
Also in this case the divergence turns out to be $\nabla F = \frac{1}{V} \partial (V F)$, where $V = det({e_i}^a)$, the determinant of the derivative of the coordinate basis with the Cartesian basis (assuming the non-coordinate basis is normalized).

$\frac{\partial}{\partial t} U + \frac{1}{V} \partial_j (V F^j) = 0$

Next integrate the divergence over the volume.

$0 = \int_{t_1}^{t_2} \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \frac{\partial}{\partial t} U + \frac{1}{V} \partial_j (V F^j) \right) V du_n ... du_1 dt $

$0 = \int_{t_1}^{t_2} \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( V \frac{\partial}{\partial t} U + \partial_j (V F^j) \right) du_n ... du_1 dt $

$0 = \int_{t_1}^{t_2} \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( V \frac{\partial}{\partial t} U \right) du_n ... du_1 dt + \int_{t_1}^{t_2} \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \partial_j (V F^j) du_n ... du_1 dt $

Assume the grid volume and the flux doesn't change over time. Assume U is constant throughout the region

$0 = \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} V du_n ... du_1 \cdot \left( U|_{t=t_2} - U|_{t=t_1} \right) + (t_2 - t_1) \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \partial_j (V F^j) du_n ... du_1 $

Let $\mathcal{V} = \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} V du_n ... du_1 $
Let $\Delta t = t_2 - t_1$

$U|_{t=t_2} = U|_{t=t_1} - \frac{\Delta t}{\mathcal{V}} \int_{u_{1,L}}^{u_{1,R}} \overset{\hat{j}}{...} \int_{u_{n,L}}^{u_{n,R}} \left( (V F^j)|_{u_j=u_{j,R}} - (V F^j)|_{u_j=u_{j,L}} \right) du_n \overset{\hat{j}}{...} du_1 $

Assume the flux is constant across the surface:

$U|_{t=t_2} = U|_{t=t_1} - \frac{\Delta t}{\mathcal{V}} \int_{u_{1,L}}^{u_{1,R}} \overset{\hat{j}}{...} \int_{u_{n,L}}^{u_{n,R}} \left( (V F^j)|_{u_j=u_{j,R}} ) - (V F^j)|_{u_j=u_{j,L}} ) \right) du_n \overset{\hat{j}}{...} du_1 $

Assume $F^j$ is constant throughout the surfaces $u = u_{j,L}$ and $u = u_{j,R}$:

$U|_{t=t_2} = U|_{t=t_1} - \frac{\Delta t}{\mathcal{V}} \left( F^j|_{u_j=u_{j,R}} \int_{u_{1,L}}^{u_{1,R}} \overset{\hat{j}}{...} \int_{u_{n,L}}^{u_{n,R}} V|_{u_j=u_{j,R}} du_n \overset{\hat{j}}{...} du_1 - F^j|_{u_j=u_{j,L}} \int_{u_{1,L}}^{u_{1,R}} \overset{\hat{j}}{...} \int_{u_{n,L}}^{u_{n,R}} V|_{u_j=u_{j,L}} du_n \overset{\hat{j}}{...} du_1 \right) $

Let $\Sigma_j = \int_{u_{1,L}}^{u_{1,R}} \overset{\hat{j}}{...} \int_{u_{n,L}}^{u_{n,R}} V|_{u_j=u_{j,R}} du_n \overset{\hat{j}}{...} du_1 $
$U|_{t=t_2} = U|_{t=t_1} - \frac{\Delta t}{\mathcal{V}} \left( (F^j \Sigma_j)|_{u_j=u_{j,R}} - (F^j \Sigma_j)|_{u_j=u_{j,L}} \right) $

Now for specific coordinate system examples of $\pmatrix{0 \\ 0}$:

Cartesian:
$g_{ij} = \delta_{ij}$
${e_{\bar{i}}}^i = \delta_\bar{i}^i$
$V = 1$
$\mathcal{V} = \Pi_{k=1}^n \Delta u_k$
$\Sigma_j = \Pi_{k=1,k\ne j}^n \Delta u_k$

$U|_{t=t_2} = U|_{t=t_1} - \frac{\Delta t}{ \Pi_{k=1}^n (\Delta u_k) } \left( ( F^j \Pi_{k=1,k\ne j}^n (\Delta u_k) )|_{u_j=u_{j,R}} - ( F^j \Pi_{k=1,k\ne j}^n (\Delta u_k) )|_{u_j=u_{j,L}} \right) $

$U|_{t=t_2} = U|_{t=t_1} - \frac{\Delta t}{\Delta u_j} \left( F^j|_{u_j=u_{j,R}} - F^j|_{u_j=u_{j,L}} \right)$

Polar, $u = \{ r, \phi \}$
$g_{\bar{i}\bar{j}} = diag(1, r^2)$
${e_{\bar{i}}}^i = diag(1, r)$
$V = r$
$\mathcal{V} = \frac{1}{2} (r_R^2 - r_L^2) (\phi_R - \phi_L) = \frac{1}{2} \Delta (r^2) \Delta \phi$
$\Sigma_r = \int_{\phi=\phi_L}^{\phi=\phi_R} r d\phi = r (\phi_R - \phi_L) = r \Delta \phi$
$\Sigma_\phi = \int_{r=r_L}^{r=r_R} r dr = \frac{1}{2} (r_R^2 - r_L^2) = \frac{1}{2} \Delta (r^2)$

$U|_{t=t_2} = U|_{t=t_1} - \frac{\Delta t}{ \frac{1}{2} \Delta (r^2) \Delta \phi } \left( F^r|_{u_r=u_{r,R}} \cdot r_R \Delta \phi - F^r|_{u_r=u_{r,L}} \cdot r_L \Delta \phi + F^\phi|_{u_\phi=u_{\phi,R}} \cdot \frac{1}{2} \Delta (r^2) - F^\phi|_{u_\phi=u_{\phi,L}} \cdot \frac{1}{2} \Delta (r^2) \right) $

$U|_{t=t_2} = U|_{t=t_1} - \Delta t \left( \frac{2}{ \Delta (r^2) } \left( F^r|_{u_r=u_{r,R}} \cdot r_R - F^r|_{u_r=u_{r,L}} \cdot r_L \right) + \frac{1}{ \Delta \phi } \left( F^\phi|_{u_\phi=u_{\phi,R}} - F^\phi|_{u_\phi=u_{\phi,L}} \right) \right) $

Spherical, $u = \{r, \theta, \phi\}$
$g_{\bar{i}\bar{j}} = diag(1, r^2, r^2 sin(\theta)^2)$
${e_{\bar{i}}}^i = diag(1, r, r sin(\theta))$
$V = r^2 sin(\theta)$
$\mathcal{V} = \frac{1}{3} (r_R^3 - r_L^3) (\phi_R - \phi_L) (cos(\theta_R) - cos(\theta_L)) = \frac{1}{3} \Delta (r^3) \Delta cos(\theta)$
$\Sigma_r = \int_{\theta=\theta_L}^{\theta=\theta_R} \int_{\phi=\phi_L}^{\phi=\phi_R} r^2 sin(\theta) d\phi d\theta = r^2 (cos(\theta_R) - cos(\theta_L)) (\phi_R - \phi_L) = r^2 \Delta cos(\theta) \Delta \phi$
$\Sigma_\theta = \int_{r=r_L}^{r=r_R} \int_{\phi=\phi_L}^{\phi=\phi_R} r^2 sin(\theta) d\phi dr = \frac{1}{3} (r_R^3 - r_L^3) sin(\theta) (\phi_R - \phi_L) = \frac{1}{3} \Delta (r^3) sin(\theta) \Delta \phi$
$\Sigma_\phi = \int_{r=r_L}^{r=r_R} \int_{\theta=\theta_L}^{\theta=\theta_R} r^2 sin(\theta) d\theta dr = \frac{1}{3} (r_R^3 - r_L^3) (cos(\theta_R) - cos(\theta_L)) = \frac{1}{3} \Delta (r^3) \Delta cos(\theta)$

$U|_{t=t_2} = U|_{t=t_1} - \frac{\Delta t}{\mathcal{V}} \left( (F^r \Sigma_r)|_{u_r=u_{r,R}} - (F^r \Sigma_r)|_{u_r=u_{r,L}} + (F^\theta \Sigma_\theta)|_{u_\theta=u_{\theta,R}} - (F^\theta \Sigma_\theta)|_{u_\theta=u_{\theta,L}} + (F^\phi \Sigma_\phi)|_{u_\phi=u_{\phi,R}} - (F^\phi \Sigma_\phi)|_{u_\phi=u_{\phi,L}} \right) $




TODO fix this: $\frac{1}{V} \partial_j V$ is not correct for all $\pmatrix{p \\ q}$ tensors.
$0 = \int_{t_1}^{t_2} \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \frac{\partial}{\partial t}({U^{i_1 ... i_p}}_{k_1 ... k_q}) + \frac{1}{V} ( \partial_j ({F^{j i_1 ... i_p}}_{k_1 ... k_q} V) ) \right) V du_n ... du_1 dt $

assuming the grid is not dependent on time ($\frac{dV}{dt} = 0$)
$0 = \int_{t_1}^{t_2} \frac{\partial}{\partial t} \left( \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( {U^{i_1 ... i_p}}_{k_1 ... k_q} V \right) du_n ... du_1 \right) dt + \int_{t_1}^{t_2} \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \partial_j ({F^{j i_1 ... i_p}}_{k_1 ... k_q} V) \right) du_n ... du_1 dt $
$0 = \left( \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( {U^{i_1 ... i_p}}_{k_1 ... k_q} V \right) du_n ... du_1 \right)|^{t_2}_{t_1} + \int_{t_1}^{t_2} \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \partial_j ({F^{j i_1 ... i_p}}_{k_1 ... k_q} V) \right) du_n ... du_1 dt $
assuming the flux is not dependent on time:
$0 = \left( \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( {U^{i_1 ... i_p}}_{k_1 ... k_q} V \right) du_n ... du_1 \right)|^{t_2}_{t_1} + \int_{t_1}^{t_2} dt \cdot \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \partial_j ({F^{j i_1 ... i_p}}_{k_1 ... k_q} V) \right) du_n ... du_1 $
$0 = \left( \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( {U^{i_1 ... i_p}}_{k_1 ... k_q} V \right) du_n ... du_1 \right)|^{t_2}_{t_1} + (t_2 - t_1) \cdot \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \partial_j ({F^{j i_1 ... i_p}}_{k_1 ... k_q} V) \right) du_n ... du_1 $
assuming the state is not dependent on space -- that it is constant throughout the integral region
$ {U^{i_1 ... i_p}}_{k_1 ... k_q}|^{t_2}_{t_1} \cdot \left( \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} V du_n ... du_1 \right) = - (t_2 - t_1) \cdot \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \partial_j ({F^{j i_1 ... i_p}}_{k_1 ... k_q} V) \right) du_n ... du_1 $
$ {U^{i_1 ... i_p}}_{k_1 ... k_q}|^{t_2}_{t_1} = - \frac{ (t_2 - t_1) }{ \left( \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} V du_n ... du_1 \right) } \cdot \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \partial_j ({F^{j i_1 ... i_p}}_{k_1 ... k_q} V) \right) du_n ... du_1 $
(But why would you integrate the volume over the volume? Shouldn't you omit the $V$?)
For coordinate basis, ${e_i}^a = \delta_i^a$ and V = 1...
$ {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_2) - {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_1) = - \frac{ (t_2 - t_1) }{ \left( \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} du_n ... du_1 \right) } \cdot \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \partial_j ({F^{j i_1 ... i_p}}_{k_1 ... k_q}) \right) du_n ... du_1 $
For a Cartesian basis, where $\nabla = \partial$ and $g_{ab} =\delta_{ab}$
$ {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_2) - {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_1) = - \frac{ (t_2 - t_1) }{ \left( \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} du_n ... du_1 \right) } \cdot \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \partial_j ({F^{j i_1 ... i_p}}_{k_1 ... k_q}) \right) du_n ... du_1 $
$ {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_2) - {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_1) = - \frac{ (t_2 - t_1) }{ \left( \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} du_n ... du_1 \right) } \cdot \Sigma_{j=1}^n \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{j-1,L}}^{u_{j-1,R}} \int_{u_{j+1,L}}^{u_{j+1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \int_{u_{j,L}}^{u_{j,R}} \left( \partial_j ({F^{j i_1 ... i_p}}_{k_1 ... k_q}) \right) du_j du_n ... du_{j+1} du_{j-1} .. du_1 $
...using the FTC ...
$ {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_2) - {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_1) = - \frac{ (t_2 - t_1) }{ \Pi_{j=1}^n \left( u_{j,R} - u_{j,L} \right) } \cdot \Sigma_{j=1}^n \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{j-1,L}}^{u_{j-1,R}} \int_{u_{j+1,L}}^{u_{j+1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( {F^{j i_1 ... i_p}}_{k_1 ... k_q}(u_j = u_{j,R}) - {F^{j i_1 ... i_p}}_{k_1 ... k_q}(u_j = u_{j,L}) \right) du_n ... du_{j+1} du_{j-1} .. du_1 $
Assuming the flux is not dependent on space:
$ {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_2) - {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_1) = - \frac{ (t_2 - t_1) }{ \Pi_{j=1}^n \left( u_{j,R} - u_{j,L} \right) } \cdot \Sigma_{j=1}^n \left( {F^{j i_1 ... i_p}}_{k_1 ... k_q}(u_j = u_{j,R}) - {F^{j i_1 ... i_p}}_{k_1 ... k_q}(u_j = u_{j,L}) \right) \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{j-1,L}}^{u_{j-1,R}} \int_{u_{j+1,L}}^{u_{j+1,R}} ... \int_{u_{n,L}}^{u_{n,R}} du_n ... du_{j+1} du_{j-1} .. du_1 $
$ {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_2) - {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_1) = - \frac{ (t_2 - t_1) }{ \Pi_{j=1}^n \left( u_{j,R} - u_{j,L} \right) } \cdot \Sigma_{j=1}^n \left( \left( {F^{j i_1 ... i_p}}_{k_1 ... k_q}(u_j = u_{j,R}) - {F^{j i_1 ... i_p}}_{k_1 ... k_q}(u_j = u_{j,L}) \right) \cdot \Pi_{k=1,k\ne j}^{k=n} \left( u_{k,R} - u_{k,L} \right) \right) $
$ {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_2) - {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_1) = - (t_2 - t_1) \cdot \Sigma_{j=1}^n \left( \frac{ 1 }{ u_{j,R} - u_{j,L} } \cdot \left( {F^{j i_1 ... i_p}}_{k_1 ... k_q}(u_j = u_{j,R}) - {F^{j i_1 ... i_p}}_{k_1 ... k_q}(u_j = u_{j,L}) \right) \right) $

Go back before the assumption of $\nabla = \partial$
... $ {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_2) - {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_1) = - \frac{ (t_2 - t_1) }{ \left( \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} du_n ... du_1 \right) } \cdot \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \nabla_j ({F^{j i_1 ... i_p}}_{k_1 ... k_q}) \right) du_n ... du_1 $
$ {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_2) - {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_1) = - \frac{ (t_2 - t_1) }{ \Pi_{j=1}^n \left( u_{j,R} - u_{j,L} \right) } \cdot \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \nabla_j ({F^{j i_1 ... i_p}}_{k_1 ... k_q}) \right) du_n ... du_1 $
...expand covariant derivative into connections...
$ {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_2) - {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_1) = - \frac{ (t_2 - t_1) }{ \Pi_{j=1}^n \left( u_{j,R} - u_{j,L} \right) } \cdot \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \partial_j \left( {F^{j i_1 ... i_p}}_{k_1 ... k_q} \right) + \Sigma_l \left( {F^{j i_1 ... i_{l-1} m i_{l+1} ... i_p}}_{k_1 ... k_q} {\Gamma^{i_l}}_{mj} \right) - \Sigma_l \left( {F^{j i_1 ... i_p}}_{k_1 ... k_{l-1} m k_{l+1} ... k_q} {\Gamma^m}_{k_l j} \right) \right) du_n ... du_1 $
$ {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_2) - {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_1) = - \frac{ (t_2 - t_1) }{ \Pi_{j=1}^n \left( u_{j,R} - u_{j,L} \right) } \cdot \left( \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \partial_j \left( {F^{j i_1 ... i_p}}_{k_1 ... k_q} \right) \right) du_n ... du_1 + \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \Sigma_l \left( {F^{j i_1 ... i_{l-1} m i_{l+1} ... i_p}}_{k_1 ... k_q} {\Gamma^{i_l}}_{mj} \right) \right) du_n ... du_1 - \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \Sigma_l \left( {F^{j i_1 ... i_p}}_{k_1 ... k_{l-1} m k_{l+1} ... k_q} {\Gamma^m}_{k_l j} \right) \right) du_n ... du_1 \right) $
...using the FTC and assuming flux is not dependent on space...
$ {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_2) - {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_1) = - \frac{ (t_2 - t_1) }{ \Pi_{j=1}^n \left( u_{j,R} - u_{j,L} \right) } \cdot \left( \Sigma_{j=1}^n \left( \left( {F^{j i_1 ... i_p}}_{k_1 ... k_q}(u_j = u_{j,R}) - {F^{j i_1 ... i_p}}_{k_1 ... k_q}(u_j = u_{j,L}) \right) \cdot \Pi_{k=1,k\ne j}^{k=n} \left( u_{k,R} - u_{k,L} \right) \right) + \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \Sigma_l \left( {F^{j i_1 ... i_{l-1} m i_{l+1} ... i_p}}_{k_1 ... k_q} {\Gamma^{i_l}}_{mj} \right) \right) du_n ... du_1 - \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \Sigma_l \left( {F^{j i_1 ... i_p}}_{k_1 ... k_{l-1} m k_{l+1} ... k_q} {\Gamma^m}_{k_l j} \right) \right) du_n ... du_1 \right) $
$ {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_2) - {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_1) = - (t_2 - t_1) \Sigma_{j=1}^n \left( \frac{ 1 }{ u_{j,R} - u_{j,L} } \cdot \left( {F^{j i_1 ... i_p}}_{k_1 ... k_q}(u_j = u_{j,R}) - {F^{j i_1 ... i_p}}_{k_1 ... k_q}(u_j = u_{j,L}) \right) + \frac{ 1 }{ \Pi_{j=1}^n \left( u_{j,R} - u_{j,L} \right) } \cdot \left( \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \Sigma_l \left( {F^{j i_1 ... i_{l-1} m i_{l+1} ... i_p}}_{k_1 ... k_q} {\Gamma^{i_l}}_{mj} \right) \right) du_n ... du_1 - \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \Sigma_l \left( {F^{j i_1 ... i_p}}_{k_1 ... k_{l-1} m k_{l+1} ... k_q} {\Gamma^m}_{k_l j} \right) \right) du_n ... du_1 \right) \right) $

Now go back to when we are using a non-coordinate system...
$ {U^{i_1 ... i_p}}_{k_1 ... k_q}|^{t_2}_{t_1} = - (t_2 - t_1) \frac{ 1 }{ \left( \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} V du_n ... du_1 \right) } \cdot \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( \nabla_j ({F^{j i_1 ... i_p}}_{k_1 ... k_q} V) \right) du_n ... du_1 $
...expand covariant derivative into connections...
$ {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_2) - {U^{i_1 ... i_p}}_{k_1 ... k_q}(t=t_1) = - \frac{ (t_2 - t_1) }{ \left( \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} V du_n ... du_1 \right) } \cdot \int_{u_{1,L}}^{u_{1,R}} ... \int_{u_{n,L}}^{u_{n,R}} \left( {e^\bar{j}}_j \partial_\bar{j} \left( {F^{j i_1 ... i_p}}_{k_1 ... k_q} V \right) + \Sigma_l \left( {F^{j i_1 ... i_{l-1} m i_{l+1} ... i_p}}_{k_1 ... k_q} {\Gamma^{i_l}}_{mj} V \right) - \Sigma_l \left( {F^{j i_1 ... i_p}}_{k_1 ... k_{l-1} m k_{l+1} ... k_q} {\Gamma^m}_{k_l j} V \right) \right) du_n ... du_1 $

Specific examples:

Cylindrical holonomic (grid metric)
connections of the holonomic basis:
connection coefficients / 2nd kind Christoffel: ${{{{ \Gamma}^{\phi}}_{\phi}}_r} = {\frac{1}{r}}$ ; ${{{{ \Gamma}^{\phi}}_r}_{\phi}} = {\frac{1}{r}}$ ; ${{{{ \Gamma}^r}_{\phi}}_{\phi}} = {-{r}}$

Cylindrical, anholonomic, normalized basis
${e_\bar{j}}^j = \left[ \begin{matrix} 1 & 0 & 0 \\ 0 & r & 0 \\ 0 & 0 & 1 \end{matrix} \right]$
${e^\bar{j}}_j = \left[ \begin{matrix} 1 & 0 & 0 \\ 0 & \frac{1}{r} & 0 \\ 0 & 0 & 1 \end{matrix} \right]$
$V = det({e_\bar{j}}^j) = r$
connections of the anholonomic normalized basis:
commutation coefficients: ${{{{ c}_{\hat{\theta}}}_r}^{\hat{\theta}}} = {\frac{1}{r}}$ ; ${{{{ c}_r}_{\hat{\theta}}}^{\hat{\theta}}} = {-{\frac{1.0}{r}}}$
connection coefficients / 2nd kind Christoffel: ${{{{ \Gamma}^{\hat{\theta}}}_{\hat{\theta}}}_r} = {\frac{1}{r}}$ ; ${{{{ \Gamma}^r}_{\hat{\theta}}}_{\hat{\theta}}} = {-{\frac{1.0}{r}}}$


Spherical, anholonomic, normalized
${e_\bar{j}}^j = \left[ \begin{matrix} 1 & 0 & 0 \\ 0 & r & 0 \\ 0 & 0 & r sin (\theta) \end{matrix} \right]$
${e^\bar{j}}_j = \left[ \begin{matrix} 1 & 0 & 0 \\ 0 & \frac{1}{r} & 0 \\ 0 & 0 & \frac{1}{r sin(\theta)} \end{matrix} \right]$
$V = det({e_\bar{j}}^j) = r^2 sin(\theta)$
connections of the anholonomic normalized basis:
${\Gamma^r}_{\phi \phi} = -\frac{1}{r}$
${\Gamma^r}_{\theta \theta} = -\frac{1}{r}$
${\Gamma^\theta}_{\theta r} = \frac{1}{r}$
${\Gamma^\theta}_{\phi \phi} = -\frac{1}{r} \frac{cos(\theta)}{sin(\theta)}$
${\Gamma^\phi}_{\phi r} = \frac{1}{r}$
${\Gamma^\phi}_{\phi \theta} = \frac{1}{r} \frac{cos(\theta)}{sin(\theta)}$
commutation coefficients: ${{{{ c}_{\hat{\theta}}}_r}^{\hat{\theta}}} = {\frac{1}{r}}$ ; ${{{{ c}_r}_{\hat{\theta}}}^{\hat{\theta}}} = {-{\frac{1}{r}}}$ ; ${{{{ c}_{\hat{\phi}}}_r}^{\hat{\phi}}} = {\frac{1}{r}}$ ; ${{{{ c}_{\hat{\phi}}}_{\hat{\theta}}}^{\hat{\phi}}} = {\frac{cos\left( \theta\right)}{{{r}} {{sin\left( \theta\right)}}}}$ ; ${{{{ c}_r}_{\hat{\phi}}}^{\hat{\phi}}} = {-{\frac{1.0}{r}}}$ ; ${{{{ c}_{\hat{\theta}}}_{\hat{\phi}}}^{\hat{\phi}}} = {\frac{-{cos\left( \theta\right)}}{{{r}} {{sin\left( \theta\right)}}}}$
connection coefficients / 2nd kind Christoffel: ${{{{ \Gamma}^{\hat{\theta}}}_{\hat{\theta}}}_r} = {\frac{1}{r}}$ ; ${{{{ \Gamma}^{\hat{\phi}}}_{\hat{\phi}}}_r} = {\frac{1}{r}}$ ; ${{{{ \Gamma}^r}_{\hat{\theta}}}_{\hat{\theta}}} = {-{\frac{1.0}{r}}}$ ; ${{{{ \Gamma}^{\hat{\phi}}}_{\hat{\phi}}}_{\hat{\theta}}} = {\frac{cos\left( \theta\right)}{{{r}} {{sin\left( \theta\right)}}}}$ ; ${{{{ \Gamma}^r}_{\hat{\phi}}}_{\hat{\phi}}} = {-{\frac{1.0}{r}}}$ ; ${{{{ \Gamma}^{\hat{\theta}}}_{\hat{\phi}}}_{\hat{\phi}}} = {\frac{-{cos\left( \theta\right)}}{{{r}} {{sin\left( \theta\right)}}}}$

for the scalar case...
$0 = \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( \frac{\partial}{\partial t}(U) + \frac{1}{ r^2 sin(\theta) } ( \nabla_j ( r^2 sin(\theta) F^j ) ) \right) r^2 sin(\theta) d\phi d\theta dr dt $
$0 = \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( r^2 sin(\theta) \frac{\partial}{\partial t}(U) + \nabla_r ( r^2 sin(\theta) F^r ) + \nabla_\theta ( r^2 sin(\theta) F^\theta ) + \nabla_\phi ( r^2 sin(\theta) F^\phi ) \right) d\phi d\theta dr dt $
$0 = \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( r^2 sin(\theta) \frac{\partial}{\partial t}(U) + sin(\theta) \partial_r ( r^2 F^r ) + r^2 \partial_\theta ( sin(\theta) F^\theta ) + r^2 sin(\theta) \partial_\phi ( F^\phi ) \right) d\phi d\theta dr dt $
$0 = \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( r^2 sin(\theta) \frac{\partial}{\partial t}(U) \right) d\phi d\theta dr dt + \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( sin(\theta) \partial_r ( r^2 F^r ) \right) d\phi d\theta dr dt + \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( r^2 \partial_\theta ( sin(\theta) F^\theta ) \right) d\phi d\theta dr dt + \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( r^2 sin(\theta) \partial_\phi ( F^\phi ) \right) d\phi d\theta dr dt $
$0 = \int_{\phi_L}^{\phi_R} d\phi \cdot \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{t_1}^{t_2} \frac{\partial}{\partial t}(U) dt r^2 sin(\theta) d\theta dr + \int_{t_1}^{t_2} dt \cdot \int_{\phi_L}^{\phi_R} d\phi \cdot \int_{\theta_L}^{\theta_R} \int_{r_L}^{r_R} \partial_r ( r^2 F^r ) dr sin(\theta) d\theta + \int_{t_1}^{t_2} dt \cdot \int_{\phi_L}^{\phi_R} d\phi \cdot \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \partial_\theta ( sin(\theta) F^\theta ) d\theta r^2 dr + \int_{t_1}^{t_2} dt \cdot \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \partial_\phi ( F^\phi ) d\phi sin(\theta) d\theta r^2 dr $
$0 = (U|_{t=t_2} - U|_{t=t_1}) \cdot \frac{1}{3} ( (r_R)^3 - (r_L)^3 ) \cdot (-cos(\theta_R) + cos(\theta_L)) \cdot (\phi_R - \phi_L) + (t_2 - t_1) \cdot ( (r_R)^2 F^r|_{r=r_R} - (r_L)^2 F^r|_{r=r_L} ) \cdot (-cos(\theta_R) + cos(\theta_L)) \cdot (\phi_R - \phi_L) + (t_2 - t_1) \cdot \frac{1}{3} ( (r_R)^3 - (r_L)^3 ) \cdot ( sin(\theta_R) F^\theta|_{\theta=\theta_R} - sin(\theta_L) F^\theta|_{\theta=\theta_L} ) \cdot (\phi_R - \phi_L) + (t_2 - t_1) \cdot \frac{1}{3} ( (r_R)^3 - (r_L)^3 ) \cdot (-cos(\theta_R) + cos(\theta_L)) \cdot ( F^\phi|_{\phi=\phi_R} - F^\phi|_{\phi=\phi_L} ) $
$U|_{t=t_2} = U|_{t=t_1} - (t_2 - t_1) \left( \frac{ (r_R)^2 F^r|_{r=r_R} - (r_L)^2 F^r|_{r=r_L} }{ \frac{1}{3} ( (r_R)^3 - (r_L)^3 ) } + \frac{ sin(\theta_R) F^\theta|_{\theta=\theta_R} - sin(\theta_L) F^\theta|_{\theta=\theta_L} }{ cos(\theta_L) - cos(\theta_R) } + \frac{ F^\phi|_{\phi=\phi_R} - F^\phi|_{\phi=\phi_L} }{ \phi_R - \phi_L } \right) $

now for the vector case...
$0 = \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( r^2 sin(\theta) \frac{\partial}{\partial t}(U^i) + \nabla_r ( r^2 sin(\theta) F^{r i} ) + \nabla_\theta ( r^2 sin(\theta) F^{\theta i} ) + \nabla_\phi ( r^2 sin(\theta) F^{\phi i} ) \right) d\phi d\theta dr dt $
$0 = \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( r^2 sin(\theta) \frac{\partial}{\partial t}(U^i) + \partial_r ( r^2 sin(\theta) F^{r i} ) + r^2 sin(\theta) F^{r k} {\Gamma^i}_{r k} + \partial_\theta ( r^2 sin(\theta) F^{\theta i} ) + r^2 sin(\theta) F^{\theta k} {\Gamma^i}_{\theta k} + \partial_\phi ( r^2 sin(\theta) F^{\phi i} ) + r^2 sin(\theta) F^{\phi k} {\Gamma^i}_{\phi k} \right) d\phi d\theta dr dt $
$0 = \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( r^2 sin(\theta) \frac{\partial}{\partial t}(U^i) + \partial_r ( r^2 sin(\theta) F^{r i} ) + r^2 sin(\theta) (F^{r r} {\Gamma^i}_{r r} + F^{r \theta} {\Gamma^i}_{r \theta} + F^{r \phi} {\Gamma^i}_{r \phi} ) + \partial_\theta ( r^2 sin(\theta) F^{\theta i} ) + r^2 sin(\theta) (F^{\theta r} {\Gamma^i}_{\theta r} + F^{\theta \theta} {\Gamma^i}_{\theta \theta} + F^{\theta \phi} {\Gamma^i}_{\theta \phi}) + \partial_\phi ( r^2 sin(\theta) F^{\phi i} ) + r^2 sin(\theta) (F^{\phi r} {\Gamma^i}_{\phi r} + F^{\phi \theta} {\Gamma^i}_{\phi \theta} + F^{\phi \phi} {\Gamma^i}_{\phi \phi}) \right) d\phi d\theta dr dt $
$\left\{ \begin{matrix} 0 = \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( r^2 sin(\theta) \frac{\partial}{\partial t}(U^r) + \partial_r ( r^2 sin(\theta) F^{r r} ) + r^2 sin(\theta) (F^{r r} {\Gamma^r}_{r r} + F^{r \theta} {\Gamma^r}_{r \theta} + F^{r \phi} {\Gamma^r}_{r \phi} ) + \partial_\theta ( r^2 sin(\theta) F^{\theta r} ) + r^2 sin(\theta) (F^{\theta r} {\Gamma^r}_{\theta r} + F^{\theta \theta} {\Gamma^r}_{\theta \theta} + F^{\theta \phi} {\Gamma^r}_{\theta \phi}) + \partial_\phi ( r^2 sin(\theta) F^{\phi r} ) + r^2 sin(\theta) (F^{\phi r} {\Gamma^r}_{\phi r} + F^{\phi \theta} {\Gamma^r}_{\phi \theta} + F^{\phi \phi} {\Gamma^r}_{\phi \phi}) \right) d\phi d\theta dr dt \\ 0 = \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( r^2 sin(\theta) \frac{\partial}{\partial t}(U^\theta) + \partial_r ( r^2 sin(\theta) F^{r \theta} ) + r^2 sin(\theta) (F^{r r} {\Gamma^\theta}_{r r} + F^{r \theta} {\Gamma^\theta}_{r \theta} + F^{r \phi} {\Gamma^\theta}_{r \phi} ) + \partial_\theta ( r^2 sin(\theta) F^{\theta \theta} ) + r^2 sin(\theta) (F^{\theta r} {\Gamma^\theta}_{\theta r} + F^{\theta \theta} {\Gamma^\theta}_{\theta \theta} + F^{\theta \phi} {\Gamma^\theta}_{\theta \phi}) + \partial_\phi ( r^2 sin(\theta) F^{\phi \theta} ) + r^2 sin(\theta) (F^{\phi r} {\Gamma^\theta}_{\phi r} + F^{\phi \theta} {\Gamma^\theta}_{\phi \theta} + F^{\phi \phi} {\Gamma^\theta}_{\phi \phi}) \right) d\phi d\theta dr dt \\ 0 = \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( r^2 sin(\theta) \frac{\partial}{\partial t}(U^\phi) + \partial_r ( r^2 sin(\theta) F^{r \phi} ) + r^2 sin(\theta) (F^{r r} {\Gamma^\phi}_{r r} + F^{r \theta} {\Gamma^\phi}_{r \theta} + F^{r \phi} {\Gamma^\phi}_{r \phi} ) + \partial_\theta ( r^2 sin(\theta) F^{\theta \phi} ) + r^2 sin(\theta) (F^{\theta r} {\Gamma^\phi}_{\theta r} + F^{\theta \theta} {\Gamma^\phi}_{\theta \theta} + F^{\theta \phi} {\Gamma^\phi}_{\theta \phi}) + \partial_\phi ( r^2 sin(\theta) F^{\phi \phi} ) + r^2 sin(\theta) (F^{\phi r} {\Gamma^\phi}_{\phi r} + F^{\phi \theta} {\Gamma^\phi}_{\phi \theta} + F^{\phi \phi} {\Gamma^\phi}_{\phi \phi}) \right) d\phi d\theta dr dt \\ \end{matrix} \right\}$
$\left\{ \begin{matrix} 0 = \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( r^2 sin(\theta) \frac{\partial}{\partial t}(U^r) + \partial_r ( r^2 sin(\theta) F^{r r} ) + \partial_\theta ( r^2 sin(\theta) F^{\theta r} ) + \partial_\phi ( r^2 sin(\theta) F^{\phi r} ) + r^2 sin(\theta) (-F^{\theta \theta} \frac{1}{r}) + r^2 sin(\theta) (-F^{\phi \phi} \frac{1}{r}) \right) d\phi d\theta dr dt \\ 0 = \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( r^2 sin(\theta) \frac{\partial}{\partial t}(U^\theta) + \partial_r ( r^2 sin(\theta) F^{r \theta} ) + \partial_\theta ( r^2 sin(\theta) F^{\theta \theta} ) + \partial_\phi ( r^2 sin(\theta) F^{\phi \theta} ) + r^2 sin(\theta) (F^{\theta r} \frac{1}{r}) + r^2 sin(\theta) (-F^{\phi \phi} \frac{cos(\theta)}{r sin(\theta)}) \right) d\phi d\theta dr dt \\ 0 = \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( r^2 sin(\theta) \frac{\partial}{\partial t}(U^\phi) + \partial_r ( r^2 sin(\theta) F^{r \phi} ) + \partial_\theta ( r^2 sin(\theta) F^{\theta \phi} ) + \partial_\phi ( r^2 sin(\theta) F^{\phi \phi} ) + r^2 sin(\theta) (F^{\phi r} \frac{1}{r} + F^{\phi \theta} \frac{cos(\theta)}{r sin(\theta)}) \right) d\phi d\theta dr dt \\ \end{matrix} \right\}$
$\left\{ \begin{matrix} 0 = \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( r^2 sin(\theta) \frac{\partial}{\partial t}(U^r) + \partial_r ( r^2 sin(\theta) F^{r r} ) + \partial_\theta ( r^2 sin(\theta) F^{\theta r} ) + \partial_\phi ( r^2 sin(\theta) F^{\phi r} ) - r sin(\theta) F^{\theta \theta} - r sin(\theta) F^{\phi \phi} \right) d\phi d\theta dr dt \\ 0 = \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( r^2 sin(\theta) \frac{\partial}{\partial t}(U^\theta) + \partial_r ( r^2 sin(\theta) F^{r \theta} ) + \partial_\theta ( r^2 sin(\theta) F^{\theta \theta} ) + \partial_\phi ( r^2 sin(\theta) F^{\phi \theta} ) + r sin(\theta) F^{\theta r} - r cos(\theta) F^{\phi \phi} \right) d\phi d\theta dr dt \\ 0 = \int_{t_1}^{t_2} \int_{r_L}^{r_R} \int_{\theta_L}^{\theta_R} \int_{\phi_L}^{\phi_R} \left( r^2 sin(\theta) \frac{\partial}{\partial t}(U^\phi) + \partial_r ( r^2 sin(\theta) F^{r \phi} ) + \partial_\theta ( r^2 sin(\theta) F^{\theta \phi} ) + \partial_\phi ( r^2 sin(\theta) F^{\phi \phi} ) + r sin(\theta) F^{\phi r} + r cos(\theta) F^{\phi \theta} \right) d\phi d\theta dr dt \\ \end{matrix} \right\}$
$\left\{ \begin{matrix} 0 = (U^r|_{t=t_2} - U^r|_{t=t_1}) \cdot (\frac{1}{3} ( (r_R)^3 - (r_L)^3 )) \cdot (-cos(\theta_R) + cos(\theta_L)) \cdot (\phi_R - \phi_L) + (t_2 - t_1) \cdot ( (r_R)^2 F^{rr}|_{r=r_R} - (r_L)^2 F^{rr}|_{r=r_L} ) \cdot (-cos(\theta_R) + cos(\theta_L)) \cdot (\phi_R - \phi_L) + (t_2 - t_1) \cdot (\frac{1}{3} ( (r_R)^3 - (r_L)^3 )) \cdot (sin(\theta_R) F^{\theta r}|_{\theta=\theta_R} - sin(\theta_L) F^{\theta r}|_{\theta=\theta_L}) \cdot (\phi_R - \phi_L) + (t_2 - t_1) \cdot (\frac{1}{3} ( (r_R)^3 - (r_L)^3 )) \cdot (-cos(\theta_R) + cos(\theta_L)) \cdot (F^{\phi r}|_{\phi=\phi_R} - F^{\phi r}|_{\phi=\phi_L}) - (t_2 - t_1) \cdot (\frac{1}{2} ( (r_R)^2 - (r_L)^2 )) \cdot (-cos(\theta_R) + cos(\theta_L)) \cdot (\phi_R - \phi_L) \cdot F^{\theta \theta} - (t_2 - t_1) \cdot (\frac{1}{2} ( (r_R)^2 - (r_L)^2 )) \cdot (-cos(\theta_R) + cos(\theta_L)) \cdot (\phi_R - \phi_L) \cdot F^{\phi \phi} \\ 0 = (U^\theta|_{t=t_2} - U^\theta|_{t=t_1}) \cdot (\frac{1}{3} ( (r_R)^3 - (r_L)^3 )) \cdot (-cos(\theta_R) + cos(\theta_L)) \cdot (\phi_R - \phi_L) + (t_2 - t_1) \cdot ( (r_R)^2 F^{r \theta}|_{r=r_R} - (r_L)^2 F^{r \theta}|_{r=r_L} ) \cdot (-cos(\theta_R) + cos(\theta_L)) \cdot (\phi_R - \phi_L) + (t_2 - t_1) \cdot (\frac{1}{3} ( (r_R)^3 - (r_L)^3 )) \cdot (sin(\theta_R) F^{\theta \theta}|_{\theta=\theta_R} - sin(\theta_L) F^{\theta \theta}|_{\theta=\theta_L}) \cdot (\phi_R - \phi_L) + (t_2 - t_1) \cdot (\frac{1}{3} ( (r_R)^3 - (r_L)^3 )) \cdot (-cos(\theta_R) + cos(\theta_L)) \cdot (F^{\phi \theta}|_{\phi=\phi_R} - F^{\phi \theta}|_{\phi=\phi_L}) + (t_2 - t_1) \cdot (\frac{1}{2} ( (r_R)^2 - (r_L)^2 )) \cdot (-cos(\theta_R) + cos(\theta_L)) \cdot (\phi_R - \phi_L) \cdot F^{\theta r} - (t_2 - t_1) \cdot (\frac{1}{2} ( (r_R)^2 - (r_L)^2 )) \cdot (sin(\theta_R) - sin(\theta_L)) \cdot (\phi_R - \phi_L) \cdot F^{\phi \phi} \\ 0 = (U^\phi|_{t=t_2} - U^\phi|_{t=t_1}) \cdot (\frac{1}{3} ( (r_R)^3 - (r_L)^3 )) \cdot (-cos(\theta_R) + cos(\theta_L)) \cdot (\phi_R - \phi_L) + (t_2 - t_1) \cdot ( (r_R)^2 F^{r \phi}|_{r=r_R} - (r_L)^2 F^{r \phi}|_{r=r_L} ) \cdot (-cos(\theta_R) + cos(\theta_L)) \cdot (\phi_R - \phi_L) + (t_2 - t_1) \cdot (\frac{1}{3} ( (r_R)^3 - (r_L)^3 )) \cdot (sin(\theta_R) F^{\theta \phi}|_{\theta=\theta_R} - sin(\theta_L) F^{\theta \phi}|_{\theta=\theta_L}) \cdot (\phi_R - \phi_L) + (t_2 - t_1) \cdot (\frac{1}{3} ( (r_R)^3 - (r_L)^3 )) \cdot (-cos(\theta_R) + cos(\theta_L)) \cdot (F^{\phi \phi}|_{\phi=\phi_R} - F^{\phi \phi}|_{\phi=\phi_L}) + (t_2 - t_1) \cdot (\frac{1}{2} ( (r_R)^2 - (r_L)^2 )) \cdot (-cos(\theta_R) + cos(\theta_L)) \cdot (\phi_R - \phi_L) \cdot F^{\phi r} + (t_2 - t_1) \cdot (\frac{1}{2} ( (r_R)^2 - (r_L)^2 )) \cdot (sin(\theta_R) - sin(\theta_L)) \cdot (\phi_R - \phi_L) \cdot F^{\phi \theta} \end{matrix} \right\}$
$\left\{ \begin{matrix} U^r|_{t=t_2} = U^r|_{t=t_1} - (t_2 - t_1) \cdot \left( \frac{ (r_R)^2 F^{rr}|_{r=r_R} - (r_L)^2 F^{rr}|_{r=r_L} }{ \frac{1}{3} ( (r_R)^3 - (r_L)^3 ) } + \frac{ sin(\theta_R) F^{\theta r}|_{\theta=\theta_R} - sin(\theta_L) F^{\theta r}|_{\theta=\theta_L} }{ -cos(\theta_R) + cos(\theta_L) } + \frac{ F^{\phi r}|_{\phi=\phi_R} - F^{\phi r}|_{\phi=\phi_L} }{ \phi_R - \phi_L } - F^{\theta \theta} - F^{\phi \phi} \right) \\ U^\theta|_{t=t_2} = U^\theta|_{t=t_1} - (t_2 - t_1) \cdot \left( \frac{ (r_R)^2 F^{r \theta}|_{r=r_R} - (r_L)^2 F^{r \theta}|_{r=r_L} }{ \frac{1}{3} ( (r_R)^3 - (r_L)^3 ) } + \frac{ sin(\theta_R) F^{\theta \theta}|_{\theta=\theta_R} - sin(\theta_L) F^{\theta \theta}|_{\theta=\theta_L} }{ -cos(\theta_R) + cos(\theta_L) } + \frac{ F^{\phi \theta}|_{\phi=\phi_R} - F^{\phi \theta}|_{\phi=\phi_L} }{ \phi_R - \phi_L } + F^{\theta r} - \frac{ sin(\theta_R) - sin(\theta_L) }{ -cos(\theta_R) + cos(\theta_L) } \cdot F^{\phi \phi} \right) \\ U^\phi|_{t=t_2} = U^\phi|_{t=t_1} - (t_2 - t_1) \cdot \left( \frac{ (r_R)^2 F^{r \phi}|_{r=r_R} - (r_L)^2 F^{r \phi}|_{r=r_L} }{ \frac{1}{3} ( (r_R)^3 - (r_L)^3 ) } + \frac{ (sin(\theta_R) F^{\theta \phi}|_{\theta=\theta_R} - sin(\theta_L) F^{\theta \phi}|_{\theta=\theta_L}) }{ (-cos(\theta_R) + cos(\theta_L)) } + \frac{ (F^{\phi \phi}|_{\phi=\phi_R} - F^{\phi \phi}|_{\phi=\phi_L}) }{ \phi_R - \phi_L } + F^{\phi r} + \frac{ (sin(\theta_R) - sin(\theta_L)) }{ (-cos(\theta_R) + cos(\theta_L)) } \cdot F^{\phi \theta} \right) \end{matrix} \right\}$
$\left[ \begin{matrix} U^r \\ U^\theta \\ U^\phi \end{matrix} \right]_{t=t_2} = \left[ \begin{matrix} U^r \\ U^\theta \\ U^\phi \end{matrix} \right]_{t=t_1} - (t_2 - t_1) \cdot \left( \frac{1}{ \frac{1}{3} ( (r_R)^3 - (r_L)^3 ) } \cdot \left[ \begin{matrix} (r_R)^2 F^{rr}|_{r=r_R} - (r_L)^2 F^{rr}|_{r=r_L} \\ (r_R)^2 F^{r \theta}|_{r=r_R} - (r_L)^2 F^{r \theta}|_{r=r_L} \\ (r_R)^2 F^{r \phi}|_{r=r_R} - (r_L)^2 F^{r \phi}|_{r=r_L} \end{matrix} \right] + \frac{1}{ -cos(\theta_R) + cos(\theta_L) } \cdot \left[ \begin{matrix} sin(\theta_R) F^{\theta r}|_{\theta=\theta_R} - sin(\theta_L) F^{\theta r}|_{\theta=\theta_L} \\ sin(\theta_R) F^{\theta \theta}|_{\theta=\theta_R} - sin(\theta_L) F^{\theta \theta}|_{\theta=\theta_L} \\ sin(\theta_R) F^{\theta \phi}|_{\theta=\theta_R} - sin(\theta_L) F^{\theta \phi}|_{\theta=\theta_L} \end{matrix} \right] + \frac{1}{ \phi_R - \phi_L } \cdot \left[ \begin{matrix} F^{\phi r}|_{\phi=\phi_R} - F^{\phi r}|_{\phi=\phi_L} \\ F^{\phi \theta}|_{\phi=\phi_R} - F^{\phi \theta}|_{\phi=\phi_L} \\ F^{\phi \phi}|_{\phi=\phi_R} - F^{\phi \phi}|_{\phi=\phi_L} \end{matrix} \right] + \left[ \begin{matrix} - F^{\theta \theta} - F^{\phi \phi} \\ F^{\theta r} - \frac{ sin(\theta_R) - sin(\theta_L) }{ -cos(\theta_R) + cos(\theta_L) } \cdot F^{\phi \phi} \\ F^{\phi r} + \frac{ (sin(\theta_R) - sin(\theta_L)) }{ (-cos(\theta_R) + cos(\theta_L)) } \cdot F^{\phi \theta} \end{matrix} \right] \right) $

TODO check this out too:
https://www.groundai.com/project/fifth-order-finite-volume-weno-in-general-orthogonally-curvilinear-coordinates4982/1