figurec pnassupportinginfo Corresponding author Predrag Popović.
E-mail<EMAIL_ADDRESS>
# Sediment load determines the shape of rivers

Predrag Popović  Olivier Devauchelle  Anaïs Abramian  and Eric Lajeunesse

##  S1 Experimental parameters and uncertainties

In this section we summarize the parameters used in the experiments of
Abramian et al. (1), show the results for each of the experimental runs, and
describe how we estimate the river properties and their uncertainties.

###  S.1.1 Details of the experimental runs

In Table 1, we show the experimental parameters and their uncertainties
(uncertainties are estimated in Abramian et al. (1)). In Fig. S1, we show the
depth and sediment flux profiles for each of the experiments of Abramian et
al. (1) and compare them to our model predictions. In Table 2 and Fig. S2, we
show the properties of these rivers.

We can see that our model predicts rivers that are slightly wider than the
experiments (by about 20%), while the depth, $D_{\text{max}}$, is not biased
in an obvious way. The width of the sediment flux profiles is captured quite
accurately, but the experimental profiles seem to be somewhat higher (by about
30% for the maximum sediment flux).

In our model, the bed slope at the river bank approximately equals the
friction coefficient, $D^{\prime}|_{D=0}\approx\mu_{t}$, since the stress,
$\tau$, vanishes when $D=0$ and the sediment flux,
$q_{\mu}\text{exp}[-\xi/\lambda]$, is negligible (Eq. 13 of the main text).
However, some of the experiments show a high bank angle, unusual for granular
material with irregular grains. This could be due to capillary forces acting
near the bank, to rivers being not fully in equilibrium, to surface armoring
by removing the loose grains (2), or to the fact that the fluid stress does
not necessarily vanish near the bank. In fact, using the value $\mu_{t}=1.2$,
combined with a length-scale, $L_{s}$, smaller by about 20% (which is within
the experimental uncertainty), makes our model agree with the experiment well
in all metrics. However, to avoid treating $\mu_{t}$ as a tuning parameter, we
decided to use $\mu_{t}=0.9$, which is the largest value estimated in an
independent experiment (3).

###  S.1.2 Estimating the river properties and their uncertainties

We estimated all of the river properties based on the cross-sections of Fig.
S1. In particular, we estimated the total sediment discharge, $Q_{s}$, as the
integral of the sediment flux profile found by grain tracking. We identify the
uncertainty in the sediment discharge from the fluctuations of $Q_{s}$ about
the mean, once the equilibrium is reached (Fig. 2 in Abramian et al. (1)).
Next, we estimate the width, $W$, as the distance from one bank to the other,
where we identify the banks by the sudden change of slope. In doing so, we
introduce an error that is of the order of the grain size. However, we can
also estimate the natural variability of the river width along its path from
the overhead images of the experiment. We find it to be about $5\%$ of the
mean — an error larger than that introduced by the cross-section measurement.
For this reason, we identify the uncertainty for the width with this natural
variability. Next, we estimate the depth, $D_{\text{max}}$, by taking the
minimum of a parabola we fit to the the river bottom, in order to minimize the
error due to bed roughness. Abramian et al. (1) estimate the measurement error
of depth due to laser inaccuracy to be about $0.5\text{ mm}$. We do not have
access to the natural variability of the depth, although we can say that it is
at least of the order of the grain size. Again, this natural variability is
greater than the measurement error estimated by Abramian et al. (1), so we
take the error of $D_{\text{max}}$ to be $d_{s}=0.83\text{ mm}$ for all
experiments. We follow a similar procedure to find the maximum sediment flux,
$q_{s,\text{max}}$ — we fit a parabola around the center of the sediment flux
profile, and take its maximum to be $q_{s,\text{max}}$. Again, we do not have
access to the natural variability of $q_{s,\text{max}}$ along the river’s
path, so we take the error of $q_{s,\text{max}}$ to be the difference between
the observed maximum sediment flux profile and the height of the fitted
parabola.

As we noted in section Dependence on Water and Sediment Discharge of the main
text, the downstream slope cannot be measured directly. However, we can
estimate it indirectly. The sediment flux at the river center,
$q_{s,\text{max}}$, is driven only by the fluid stress, which is about
$\tau\approx\rho_{f}gD_{\text{max}}S$ in the shallow-water approximation.
Therefore, dropping the contributions of momentum diffusion and gravity from
Eq. 13 of the main text, we find
$D_{\text{max}}S/L_{s}\approx\mu_{t}+q_{s,\text{max}}/q_{\mu}$, and we can
estimate the slope in the experiment as

$S\approx\frac{L_{s}}{D_{\text{max}}}\left(\mu_{t}+\frac{q_{s,\text{max}}}{q_{\mu}}\right)\text{
. }$ (S1)

This is the estimate we show in Fig. S2. This estimate is not independent from
our model. We can, nevertheless, use it to check the consistency between the
parameters in our model and the experiments. We estimate the error for the
slope calculated in this way as a combination of the errors for the quantities
that enter Eq. S1.

The quantities we described above are all point-measurements and, therefore,
are sensitive to the roughness of the bed and measurement error. On the other
hand, integral quantities such as the mean sediment flux, $\langle
q_{s}\rangle$, and transport width, $W_{T}$, are robust against such errors.
In this case, we expect that the uncertainty of $\langle q_{s}\rangle$ and
$W_{T}$ is mostly due to their natural variability in time and space, rather
than by any measurement error. We, however, do not have access to measurements
that would allow us to estimate the variability, so we assume that the natural
variability of $\langle q_{s}\rangle$ is of the same order as the variability
of the sediment flux profile across the channel. Therefore, we take the
uncertainty of $\langle q_{s}\rangle$ to be equal to the uncertainty of
$q_{s,\text{max}}$ — a conservative estimate. Finally, since
$W_{T}=Q_{s}/\langle q_{s}\rangle$, we estimate its uncertainty as a
combination of the uncertainty of the sediment discharge and that of the mean
sediment flux.

Table 1: Table of parameters used in the experiments and the model. Definition |  Notation |  Value |  Unit
---|---|---|---
Grain diameter |  $d_{s}$ |  $0.83\pm 0.2$ |  mm
Fluid viscosity |  $\nu$ |  $10^{-5}$ |  m2s-1
Fluid density |  $\rho_{f}$ |  $1160\pm 5$ |  kg m-3
Sediment density |  $\rho_{s}$ |  $1490$ |  kg m-3
Fluid discharge |  $Q_{w}$ |  $0.97\pm 0.05$ |  l min-1
Sediment discharge |  $Q_{s}$ |  $(0,60)$ |  grains $\text{s}^{-1}$
Threshold Shields parameter |  $\theta_{t}$ |  $0.167\pm 0.003$ |  None
Friction coefficient |  $\mu_{t}$ |  $0.9\pm 0.2$ |  None
Sediment diffusion length |  $\lambda$ |  $0.10\pm 0.03$ |  mm
Sediment flux scale |  $q_{\mu}$ |  $107\pm 30$ |  grains cm-1s-1
Figure S1: Depth and sediment profiles for all experiments of Abramian et al.
(1). The left column are the measured (brown lines) and modeled (blue dashed
lines) depth profiles. The right column are the measured (red lines) and
modeled (blue lines) sediment flux profiles. Table 2: River properties for
each of the experiments performed by Abramian et al. (1).

Experiment label |  Sediment discharge |  Unit |
---|---|---|---
$1$ | $0$ | $\text{grains s}^{-1}$ |
$2$ | $12.6\pm 1.2$ |
$3$ | $24\pm 2.4$ |
$4$ | $44.9\pm 4.5$ |
$5$ | $59.7\pm 6.0$ |
Quantity |  Experiment label |  Value |  Unit
Width, $W$ | $1$ | $3.20\pm 0.16$ | cm
$2$ | $3.36\pm 0.17$
$3$ | $3.53\pm 0.18$
$4$ | $4.83\pm 0.24$
$5$ | $4.56\pm 0.23$
Depth, $D_{\text{max}}$ | $1$ | $0.744\pm 0.083$ | cm
$2$ | $0.748\pm 0.083$
$3$ | $0.530\pm 0.083$
$4$ | $0.474\pm 0.083$
$5$ | $0.539\pm 0.083$
Maximum sediment flux, $q_{s,\text{max}}$ | $1$ | $0$ | $\text{grains s}^{-1}\text{cm}^{-1}$
$2$ | $15.96\pm 0.09$
$3$ | $24.65\pm 1.44$
$4$ | $26.97\pm 3.07$
$5$ | $27.42\pm 2.44$
Transport width, $W_{T}$ | $1$ | n/a | cm
$2$ | $1.01\pm 0.11$
$3$ | $1.24\pm 0.22$
$4$ | $2.16\pm 0.53$
$5$ | $2.46\pm 0.49$
Downstream slope, $S$ | $1$ | $0.005\pm 0.002$ | None
$2$ | $0.006\pm 0.002$
$3$ | $0.009\pm 0.004$
$4$ | $0.011\pm 0.005$
$5$ | $0.009\pm 0.004$

Figure S2: Comparison of various river properties in dimensional units between
our model and experiments of Abramian et al. (1). Slope in panel (c) is
estimated using Eq. S1. We only show the weak transport regime for the
properties related to the sediment flux.

##  S2 Fluid induced stress, $\tau$

In this section we will first derive Eq. 3 of the main text for the fluid
stress, $\tau$. Then, we will show that this stress can be written as a series
expansion that assumes the aspect ratio of the channel is large. We will show
that our model, Eq. 13 of the main text, is the first term in this expansion.
Finally, we will test our approximation on several tractable examples.

###  S.2.1 Deriving the fluid stress equation

In this section, following Devauchelle et al. (4), we derive equation 3 of the
main text for the fluid stress, $\tau$. To that end, we first integrate the
Stokes’ flow equation, Eq. 2 of the main text, over the vertical to get

$\eta_{v}\int_{-D}^{0}\frac{\partial^{2}u}{\partial
y^{2}}\text{d}z-\tau_{z}=-\rho_{f}gSD\text{ , }$ (S2)

where the $z$-component of the stress is
$\tau_{z}\equiv-\eta_{v}\left.\frac{\partial u}{\partial z}\right|_{z=-D}$,
and we have used the boundary condition that the stress vanishes at the
surface, $\tau_{z}(z=0)=0$. To get Eq. 3 of the main text, we need to pull the
$y$-derivatives in Eq. S2 outside the integral. Doing this, we get

$\int_{-D}^{0}\frac{\partial^{2}u}{\partial
y^{2}}\text{d}z=\frac{\partial^{2}}{\partial
y^{2}}\int_{-D}^{0}u\text{d}z+\frac{\text{d}D}{\text{d}y}\left.\frac{\partial
u}{\partial y}\right|_{z=-D}\text{ , }$ (S3)

where we have used the boundary condition that the velocity vanishes at the
boundary, $u(z=-D)=0$. We can use this same condition to relate the term
$\left.\frac{\partial u}{\partial y}\right|_{z=-D}$ to the shear stress,
$\tau_{z}$. Namely, differentiating the boundary condition, $u(z=-D)=0$,
yields

$\left.\frac{\partial u}{\partial
y}\right|_{z=-D}=\frac{\text{d}D}{\text{d}y}\left.\frac{\partial u}{\partial
z}\right|_{z=-D}\text{ . }$ (S4)

Substituting relations Eqs. S3 and S4 into Eq. S2, we find

$\eta_{v}\frac{\text{d}^{2}}{\text{d}y^{2}}(D\bar{u})-\tau_{z}\left(1+\left(\frac{\text{d}D}{\text{d}y}\right)^{2}\right)=-\rho_{f}gSD\text{
, }$ (S5)

where we have introduced $\bar{u}(y)\equiv\frac{1}{D}\int_{-D}^{0}u\text{d}z$.
Equation S4 relates the $y$ and $z$ components of the stress as
$\tau_{y}=D^{\prime}\tau_{z}$. This means that the total stress, $\tau$, is
related to $\tau_{z}$ as

$\tau=\left(\tau_{z}^{2}+\tau_{y}^{2}\right)^{1/2}=\tau_{z}\left(1+D^{\prime
2}\right)^{1/2}{,}$ (S6)

where prime stands for $\text{d}/\text{d}y$. Substituting Eq. S6 into Eq. S5,
we can express $\tau$ as

$\tau=\frac{\eta_{v}(D\bar{u})^{\prime\prime}+\rho_{f}gSD}{\left(1+D^{\prime
2}\right)^{1/2}}\text{ , }$ (S7)

Identifying $\left(1+D^{\prime 2}\right)^{-1/2}$ with $\cos{\phi}$, we finally
retrieve Eq. 3 of the main text

$\tau=\left(\eta_{v}(D\bar{u})^{\prime\prime}+\rho_{f}gSD\right)\cos{\phi}\text{
, }$ (S8)

So far, we made no approximation to get here from the original Stokes equation
(Eq. 2 of the main text).

###  S.2.2 The stress approximation as a series expansion

In the main text, we used the shallow-water velocity to approximate the fluid
stress (Eqs. 12 and 13 of the main text). In this section, we will show that
this approximation is the first term in a series expansion that assumes the
aspect ratio of the river is large. In this way, this approximation may be
systematically improved, assuming that the series converges.

We begin with the Stokes equation

$\frac{\partial^{2}u}{\partial y^{2}}+\frac{\partial^{2}u}{\partial
z^{2}}=\frac{-gS}{\nu}$ (S9)

If the channel depth, $D(y)$, varies smoothly, then the second $y$-derivative
of velocity scales as $U/W^{2}$ while the second $z$-derivative scales as
$U/D_{\text{max}}^{2}$, where $U$ is the velocity scale. Therefore, if the
aspect ratio, $W/D_{\text{max}}$, is large, the $z$-derivative is much larger
than the $y$-derivative. We can make this obvious by rescaling the variables
in the Stokes equation as

$\hat{{y}}=\frac{y}{W}\quad\text{ ,
}\quad\hat{{z}}=\frac{z}{D_{\text{max}}}\quad\text{ , }\quad\hat{{u}}\equiv
u\frac{\nu}{gSD_{\text{max}}^{2}}$ (S10)

In these new coordinates, the channel has unit width and depth, and the fluid
is driven with a unit forcing. The Stokes equation then becomes

$\frac{D_{\text{max}}^{2}}{W^{2}}\frac{\partial^{2}\hat{{u}}}{\partial\hat{{y}}^{2}}+\frac{\partial^{2}\hat{{u}}}{\partial\hat{{z}}^{2}}=-1$
(S11)

We can then expand the velocity as

$\hat{{u}}=\hat{{u}}^{(0)}+\hat{{u}}^{(1)}+\hat{{u}}^{(2)}+\text{... ,}$ (S12)

where $\hat{{u}}^{(0)}$ is the term of the order unity, $\hat{{u}}^{(1)}$ is
the term of the order $D_{\text{max}}^{2}/W^{2}$, $\hat{{u}}^{(2)}$ is the
term of the order $D_{\text{max}}^{4}/W^{4}$, etc. Equating orders in the
Stokes equation, we get (back in the physical coordinates)

$\displaystyle\frac{\partial^{2}u^{(0)}}{\partial z^{2}}=\frac{-gS}{\nu}\text{
, }$ (S13) $\displaystyle\frac{\partial^{2}u^{(0)}}{\partial
y^{2}}+\frac{\partial^{2}u^{(1)}}{\partial z^{2}}=0\text{ , ...}$ (S14)

The boundary conditions for each term are the same — velocity at each order
must vanish on the channel bed and its $z$-derivative must vanish on the fluid
surface. From Eq. S13, we can see that the zeroth order term is in fact equal
to the shallow-water velocity, $u^{(0)}\equiv u_{\text{sw}}$. All of the terms
in the expansion are polynomial in $z$, and we can easily find each term of
order $(n+1)$ based on the previous term of order $(n)$ by performing
integration over $z$. In this way, we find for the first two terms:

$\displaystyle u^{(0)}=\frac{gS}{2\nu}\left(D^{2}-z^{2}\right)\text{ , }$
(S15) $\displaystyle
u^{(1)}=\frac{gS}{4\nu}\left(D^{2}-z^{2}\right)(D^{2})^{\prime\prime}\text{ .
}$ (S16)

The term $u^{(1)}$ is inversely proportional to the square of the aspect ratio
since the correcting factor, $(D^{2})^{\prime\prime}$, is of order
$D_{\text{max}}^{2}/W^{2}$. From here, we can find the vertically averaged
velocity as an expansion:

$\displaystyle\bar{u}=\bar{u}^{(0)}+\bar{u}^{(1)}+\text{...}$ (S17)
$\displaystyle\bar{u}^{(0)}=\frac{gSD^{2}}{3\nu}\text{ , }$ (S18)
$\displaystyle\bar{u}^{(1)}=\frac{gSD^{2}}{6\nu}(D^{2})^{\prime\prime}\text{ .
}$ (S19)

Equation S18 for the zeroth order velocity is the same as Eq. 12 of the main
text for the shallow-water velocity. With this, we can find the fluid stress.
Recall the exact equation for the stress, Eq. 3 of the main text,

$\tau=\left(\rho_{f}gSD+\rho_{f}\nu(\bar{u}D)^{\prime\prime}\right)\cos\phi\text{
. }$ (S20)

We can also write the stress as an expansion

$\tau=\tau^{(0)}+\tau^{(1)}+\text{... ,}$ (S21)

where $\tau^{(0)}$ is the leading order term, and $\tau^{(1)}$ is proportional
to $D_{\text{max}}^{2}/W^{2}$. Using the expansion for $\bar{u}$, we see that

$\displaystyle\tau^{(0)}=\rho_{f}gSD\cos\phi$ (S22)
$\displaystyle\tau^{(1)}=\rho_{f}\nu(\bar{u}^{(0)}D)^{\prime\prime}\cos\phi=\left[\rho_{f}gS\frac{1}{3}\left(D^{3}\right)^{\prime\prime}\right]\cos\phi\text{
. }$ (S23)

Since $(D^{3})^{\prime\prime}\sim D_{\text{max}}^{3}/W^{2}$, the first-order
term is about $D_{\text{max}}^{2}/W^{2}$ times smaller than the leading order
term: $\tau^{(1)}\sim\tau^{(0)}D_{\text{max}}^{2}/W^{2}$. The cosine term in
the above equations is $\cos\phi=\left(1+D^{\prime 2}\right)^{-1/2}$, and,
since $D^{\prime 2}\sim D_{\text{max}}^{2}/W^{2}$, this term could also be
expanded in a series but, for convenience, we keep the entire term here.

The zeroth order term, $\tau^{(0)}$, with $\cos\phi\approx 1$, is the shallow-
water stress that was used previously to estimate the shape of inert rivers
(5). Our approximation (Eq. 13 of the main text) amounts to truncating the
$\tau$-series at first order, $\tau\approx\tau^{(0)}+\tau^{(1)}$. In the main
text, we showed that the river is formed by the combined action of fluid
stress and gravity. The gravity term in the force ratio, $(D^{\prime})^{2}$,
is of the order of $D_{\text{max}}^{2}/W^{2}$, and is, therefore, of the same
order of magnitude as $\tau^{(1)}$. For this reason, keeping only the shallow-
water term in the force ratio, $\mu$, is inconsistent, since it neglects a
stress term that is comparable to gravity.

Although the expansion above may be corrected to arbitrary order, it is always
local — the stress is always a function of the depth and its derivatives. To
get non-local effects, whereby the stress depends on the entire bed shape, we
would have to add infinitely many terms in the expansion.

###  S.2.3 Testing the approximation

In the previous section we showed that the fluid velocity and stress can be
found using a series expansion when the flow is shallow enough. In this
section, we will explicitly test this expansion on several examples for which
exact solutions exist.

As a first example, we look at the flow over an inclined surface with an
opening angle $\phi$, that we assume to be small (Fig. S3a). In this case, the
depth profile is given by

$D(y)=y\tan\phi\text{ . }$ (S24)

This problem is not well-posed since we do not specify the boundary condition
on the open edge of the wedge. Nevertheless, if the conditions far away do not
matter, the solution of Stokes equation and the corresponding stress should be

$\displaystyle u=\frac{gS}{2\nu}\frac{D^{2}-z^{2}}{1-(\tan\phi)^{2}}\text{ ,
}$ (S25) $\displaystyle\tau_{z}=\frac{gS\rho_{f}D}{1-(\tan\phi)^{2}}\text{ ,
}$ (S26)

where $\tau_{z}$ is the vertical component of the stress. We can see that the
exact solution has the same form as the shallow-water approximation with the
correction factor $1-(\tan\phi)^{2}$. Our approximation to first order yields

$\displaystyle
u=\frac{gS}{2\nu}\left(D^{2}-z^{2}\right)\left(1+(\tan\phi)^{2}\right)\text{ ,
}$ (S27) $\displaystyle\tau_{z}=gS\rho_{f}D\left(1+(\tan\phi)^{2}\right)\text{
, }$ (S28)

where we have expanded the cosine terms in Eq. S22 and S23 to first order in
$\phi$. Equations S27 and S28 show that our approximation is correct to first
order in $\phi$. We compare the stress in our approximation to the exact
solution, and to the shallow-water stress in Fig. S3d. For $\tan\phi\geq 1$,
our approximation breaks down. In fact, the exact solution given by Eqs. S25
and S28 also breaks down at $\tan\phi=1$. For large angles, the exact solution
we presented above is not viable because, in that case, the boundary condition
at the far end cannot be ignored. This is an example of a non-local effect we
mentioned in the previous section — for $\tan\phi\geq 1$, our approximation
breaks down because it can only relate the flow to the local bed shape.

As our second example, we consider the flow in an elliptic channel with a
ratio of semi-major to semi-minor axes equal to $R=W/(2D_{\text{max}})$, that
we assume to be large (Fig. S3b). The depth profile is given by

$D(y)=\frac{1}{R}\sqrt{\left(\frac{W}{2}\right)^{2}-y^{2}}\text{ , }$ (S29)

The exact solution of the Stokes equation in such a channel is

$\displaystyle
u=\frac{gS}{2\nu}\frac{\left(D^{2}-z^{2}\right)R^{2}}{(1+R^{2})}$ (S30)
$\displaystyle\tau_{z}=\frac{gS\rho_{f}DR^{2}}{1+R^{2}}\text{ , }$ (S31)

Our approximation to first order gives

$\displaystyle
u=\frac{gS}{2\nu}\left(D^{2}-z^{2}\right)\left(1-\frac{1}{R^{2}}\right)$ (S32)
$\displaystyle\tau_{z}=gS\rho_{f}D\left(1-\frac{1}{R^{2}}\right)$ (S33)

Our approximation provides the correct first order when $R$ is large, but
fails when $R\geq 1$. The reason is the same as before — in a narrow channel,
rather than being controlled by the local depth configuration, the fluid
velocity becomes dominated by the side walls. We compare our approximation
with the exact solution and the shallow-water approximation in Fig. S3e.

Finally, we consider the flow the flow over a sinusoidally perturbed bed (Fig.
S3c). The depth profile in this case is given by

$D(y)=D_{0}+\delta\sin(ky)\text{ , }$ (S34)

where $\delta$ is the amplitude and $k$ is the wavenumber of the perturbation.
If $\delta$ is small, the full Stokes equation can be linearized to find the
stress. Abramian et al. (6) derived the expression for the stress in this
case:

$\displaystyle\tau_{z}=gS\rho_{f}D_{0}+gS\rho_{f}\delta\left(1-kD_{0}\text{tanh}(kD_{0})\right)\sin(ky)\text{
. }$ (S35)

We expect our approximation to work when the wavelength of the perturbation is
large compared with the flow depth (i.e. when $kD_{0}$ is small), but should
not be limited by $\delta$. Therefore, we can compare our approximation to the
result of Abramian et al. (6) when both $\delta/D_{0}$ and $kD_{0}$ are small.
In this case, our approximation yields:

$\displaystyle\tau_{z}=gS\rho_{f}D_{0}+gS\rho_{f}\delta\left(1-(kD_{0})^{2}\right)\sin(ky)$
(S36)

Again, this is the correct first order expansion of Eq. S35 for small
$kD_{0}$. In the classical shallow-water theory, the stress,
$\tau_{z}=gS\rho_{f}D$, is always in phase with the depth perturbation.
Conversely, in the linearization of Abramian et al. (6), Eq. S35, the phase of
the stress can reverse when the wavelength of the perturbation is small
enough, so that the stress maximum is where the flow is shallowest. This
feature of the full linearized equation is reproduced in our approximation.
This shows that our approximation can capture qualitative effects of the
cross-stream diffusion of momentum. We compare our first order approximation,
linearization of Abramian et al. (6) (Eq. S35), and the shallow-water
approximation for a sinusoidal perturbation in Fig. S3f.

It is clear that our approximation cannot always work. For example, in a
rectangular channel, our approximation predicts uniform stress and velocity
above the entire flat bottom, to all orders of the approximation. This is
clearly not the case in reality, since the fluid must slow down near the
vertical walls. The effect of the walls is, again, an example of a non-
locality — the velocity of the fluid is not only determined by the local depth
and its derivatives, but is also affected by the far away walls. Therefore, in
all cases we considered above, the failure of our approximation was related to
the non-locality of the velocity field.

Figure S3: Comparison of stress between our first order approximation (blue
lines), shallow-water approximation (black dashed lines), and a analytical
reference solution (red lines) in examples for which an analytical solution is
available. The upper row shows the channel shape, along with the velocity
field of the analytical solution (lighter blue stands for faster moving
fluid). The lower row shows the normalized stress profile,
$\tau(y/L)/\rho_{f}gSL$, where $L$ is a reference length scale that is
different for each example. The spatial coordinates, $y$ and $z$, in upper and
lower panels are also normalized by $L$. (a) and (d) Flow over an inclined
surface of an opening angle $\phi=25^{\circ}$ (Eqs. S24 to S28). The reference
length scale in this case can be arbitrarily chosen. (b) and (e) Flow in an
elliptic channel with a ratio of the semi-major to semi-minor axis $R=1.8$
(Eqs. S29 to S33). The reference length scale is the channel width, $L=W$. (c)
and (f) Flow over a surface with a small-amplitude, large-wavelength
perturbation, $\delta/D_{0}=0.2$ and $kD_{0}=0.9$ (Eqs. S34 to S36). The
reference length scale is the unperturbed channel depth, $L=D_{0}$.

##  S3 River profile as a function of water and sediment discharges

As we mentioned in the main text, our model is completely determined by
setting five parameters — $\mu_{t}$, $\lambda$, $L_{s}$, $S$, and $\xi$. Of
those, $\mu_{t}$, $\lambda$, and $L_{s}$ are directly measurable in the
experiments, while $S$ and $\xi$ depend implicitly on discharges $Q_{w}$ and
$Q_{s}$. In this section we describe how we numerically find this relation to
get the river profiles as a function of the discharges of water and sediment.

###  S.3.1 Solving the boundary value problem

To get the river as a function of the discharges, we first need to be able to
solve our model, Eq. 13 of the main text, for given values of the parameters
$\mu_{t}$, $\lambda$, $L_{s}$, $S$, and $\xi$. This presents a slight
challenge. We described this model as a second order boundary value problem
with boundary conditions $D=0$ when $y=0$, and $D^{\prime}=0$ when $y=-W/2$.
There are two problems here — first, the second derivative,
$D^{\prime\prime}$, given by Eq. 13 of the main text, diverges when
$D\rightarrow 0$, and, second, the width of the river, $W$, is unknown a
priori, so the second boundary condition is ill-defined. To deal with the
first issue, we expand the solution near $D=0$, so that we only consider
depths larger than some small value, $\epsilon$. Thus, we change the first
boundary condition to $D=\epsilon$ when $y=0$ (we can arbitrarily shift the
solution along $y$ due to translational invariance). To deal with the second
issue, we replace the second boundary condition with a new condition for
$D^{\prime}$ that, like the first boundary condition, starts at $y=0$. In this
way, instead of dealing with a boundary value problem, we can solve a simpler
initial value problem. To do this, we need to approximate the slope near the
river bank, $D^{\prime}$, when the depth, $D=\epsilon$, is small. We use the
fact that $D^{\prime\prime}$, which is a known function of $D$ and
$D^{\prime}$, diverges when $D\rightarrow 0$. In particular, on the left river
bank (when $D=\epsilon$ and $D^{\prime}>0$),
$D^{\prime\prime}\rightarrow\infty$ as $\epsilon\rightarrow 0$ for all
solutions with $D^{\prime}$ smaller than that of the river, while
$D^{\prime\prime}\rightarrow-\infty$ for all solutions with $D^{\prime}$
greater than that of the river. Therefore, for small $D$, the river
approximately lies on the curve $D^{\prime\prime}(D,D^{\prime})=0$. We can
thus invert the relation $D^{\prime\prime}(\epsilon,D^{\prime})=0$ to find
$D^{\prime}$ when $y=0$, therefore defining the second initial condition. In
summary, we solve the initial value problem, $D=\epsilon$ and
$D^{\prime\prime}(\epsilon,D^{\prime})=0$ when $y=0$, by forward stepping and
we stop when we reach the river center, $D^{\prime}=0$. Once we find the left
bank, the right bank follows by symmetry.

###  S.3.2 Interpolating the model

Once we can find a river for given model parameters, we fix $\mu_{t}$,
$\lambda$, and $L_{s}$ to their experimental values and then solve our model
(Eq. 13 of the main text) for multiple values of $S$ and $\xi$ to create a
grid of solutions for a range of model parameters (Fig. S4a). For each
numerical solution in this grid, we can find $Q_{w}$ and $Q_{s}$ using Eq. 15
of the main text. Interpolating over this grid, we then find $S$ and $\xi$ (or
any other property of the river rivers such as the width or depth) as
functions of $Q_{w}$ and $Q_{s}$.

One problem is how to choose the range of $S$ and $\xi$ for the solution grid.
Namely, increases in the sediment discharge are controlled by minute
differences of $\tilde{\xi}$ (on the order of $10^{-9}$) from
$\tilde{\xi}_{c}$, where $\tilde{\xi}\equiv\xi S/L_{s}$ (see section Inert,
active, and limiting river of the main text). To probe a significant range of
sediment discharge, we need to explore $\tilde{\xi}$ in a very narrow range
around $\tilde{\xi}_{c}$. We do not know the value of $\tilde{\xi}_{c}$ a
priori, and we have to find it numerically. Moreover, the value of
$\tilde{\xi}_{c}$ changes with $S$, so for each $S$, we need to independently
estimate $\tilde{\xi}_{c}$ to high precision. Since there are no river-like
solutions below $\tilde{\xi}_{c}$, we can easily identify whether a solution
with given parameters $S$ and $\xi$ has $\tilde{\xi}>\tilde{\xi}_{c}$ or
$\tilde{\xi}<\tilde{\xi}_{c}$. Therefore, we can find $\tilde{\xi}_{c}$ by a
shooting method. Moreover, in the process of approaching $\tilde{\xi}_{c}$ by
interval halving, we solve our model for multiple values of $\xi$, so we can
immediately use the solutions with $\tilde{\xi}<\tilde{\xi}_{c}$ in our
interpolation grid. This also has the useful property that exponentially
approaching $\tilde{\xi}_{c}$ leads to approximately linearly increasing
$Q_{s}$. Therefore, by finding $\tilde{\xi}_{c}$, we probe the space of
sediment discharge more or less uniformly. Once we find $\tilde{\xi}_{c}$ for
a given $S$, we add the solutions to the interpolation grid, we choose
randomly a different $S$, and we repeat the procedure.

To probe the relevant range of $S$ that corresponds roughly to the fluid
discharge in the experiments, $Q_{w}$, we look at inert rivers ($Q_{s}=0$).
Using Eq. 15 of the main text, we can write the fluid discharge as

$Q_{w}=\frac{gL_{s}^{4}}{\nu
S^{3}}\int_{-\tilde{W}/2}^{\tilde{W}/2}\frac{\tilde{D}^{3}}{3}\text{d}\tilde{y}\text{
, }$ (S37)

where tildes stand for quantities made non-dimensional with the length scale,
$L_{s}/S$, as defined in Eq. 16 of the main text. The term
$\tilde{Q}_{w}=\int\tilde{D}^{3}/3\text{d}\tilde{y}$ is the non-dimensional
fluid discharge that depends on the shape of the river, but not its size. As
we have discussed in section Inert, active, and limiting river of the main
text, the shape of non-dimensional inert rivers only depends on $\mu_{t}$. So,
for an inert river, $\tilde{Q}_{w,0}(\mu_{t})$ is independent of the fluid
discharge. Therefore, we find that

$S_{0}=\left(\frac{gL_{s}^{4}\tilde{Q}_{w,0}(\mu_{t})}{\nu
Q_{w}}\right)^{1/3}\quad\text{ when }\quad Q_{s}=0\text{ . }$ (S38)

Based on this expression, we estimate the order of magnitude of the
experimental slope. Then, taking values between, for example, $S_{0}/5$ and
$5S_{0}$, ensures that we cover the range $S$ that is relevant for our
experiment, even when the sediment discharge is finite.

###  S.3.3 Dependence on water and sediment discharges

In Fig. S5, we show how several river properties depend on the discharges of
fluid and sediment in our model. Namely, the aspect ratio increases with the
sediment discharge, and only depends on the fluid discharge for large values
of the sediment discharge.The shape of an inert river ($Q_{s}=0$) is
independent of the fluid discharge. The maximum sediment flux,
$q_{s,\text{max}}$, increases with the sediment discharge and is largely
independent of the fluid discharge. In fact, for large fluid discharge,
$Q_{w}$, the maximum sediment flux saturates at
$q_{\mu}(\tilde{D}_{\text{max},0}-\mu_{t})$, as predicted by the Parker
regime. Figure S5f shows that the downstream slope, $S$, scales approximately
as $Q_{w}^{-1/3}$ (for an inert river, this scaling is exact, Eq. S38).
Therefore, the scale of the river, given by $L_{s}/S$ roughly increases with
fluid discharge as $Q_{w}^{1/3}$. The sediment discharge affects the slope,
$S$, and the size of the river only slightly (Fig. S5c).

Figure S4: An interpolation grid of solutions used to find model dependence on
discharges of fluid and sediment. Each grey point is a solution to our model
with different values of $S$ and $\xi$ (and fixed $\mu_{t}$, $\lambda$,
$L_{s}$). For visual clarity, we only show every fifth point. (a)
Interpolation grid in the space of model parameters, $\xi$ and $S$. The
horizontal axis shows the logarithm of the difference between $\tilde{\xi}$
and the limiting value $\tilde{\xi}_{c}$ which is approximately proportional
to the sediment discharge. Solid lines correspond to constant sediment
discharge (red lines; numbers denote the value in $\text{grains s}^{-1}$) and
fluid discharge (blue lines; numbers denote the value in $\text{l min}^{-1}$)
found after interpolation. The non-smooth appearance of some of these lines
for small and large values of $S$ is due to the sparseness of the grid. (b)
Interpolation grid in the space of fluid and sediment discharge, $Q_{w}$ and
$Q_{s}$. Each point on this grid has a corresponding point in panel (a). Solid
lines correspond to constant aspect ratio (red lines) and width (blue lines;
numbers denote the value in cm).  Figure S5: (a-c) River properties in our
model as a function of sediment discharge, $Q_{s}$, for various values of
fluid discharge, $Q_{w}$. (a) Aspect ratio, $D_{\text{max}}/W$. (b) Maximum
sediment flux, $q_{s,\text{max}}/q_{\mu}$. (c) Downstream slope, $S$. (d-f)
River properties in our model as a function of fluid discharge, $Q_{w}$, for
various values of sediment discharge, $Q_{s}$. (d) Aspect ratio,
$D_{\text{max}}/W$. (e) Maximum sediment flux, $q_{s,\text{max}}/q_{\mu}$. (f)
Downstream slope, $S$. The black dashed line corresponds to a curve
$S=(Q_{w}^{*}/Q_{w})^{1/3}$, where $Q_{w}^{*}\equiv
gL_{s}^{4}\tilde{Q}_{w,0}/\nu\approx 2.6\times 10^{-7}\text{ l}\text{
min}^{-1}$ is the characteristic fluid discharge.

##  S4 Non-dimensional river model

In section Inert, active, and limiting river of the main text, we discussed
our model in non-dimensional form (Eq. 17 of the main text). In this section,
we discuss the mathematical properties of Eq. 17, such as its fixed points and
its phase portrait. These properties clarify why an infinite, limiting river
exists in our model. In addition, we support claims made in sections Inert,
active, and limiting river and The Parker regime of the main text about inert
and infinite rivers.

###  S.4.1 Fixed points

Since Eq. 17 of the main text is of second order, the second derivative of the
depth, $\tilde{D}^{\prime\prime}$, is a function of $\tilde{D}$ and
$\tilde{D}^{\prime}$. Thus, we can view it as a dynamical system that can be
integrated by forward-stepping, starting from a given value of $\tilde{D}$ and
$\tilde{D}^{\prime}$.

Fixed points of this equation are defined as points in the
$\tilde{D}$-$\tilde{D}^{\prime}$ space (the phase space), that give constant
solutions under integration. This is satisfied when $\tilde{D}^{\prime}=0$ and
$\tilde{D}^{\prime\prime}=0$. Inserting these conditions into Eq. 17 of the
main text yields

$\tilde{D}=\mu_{t}+e^{(\tilde{D}-\tilde{\xi})/\tilde{\lambda}}\text{ . }$
(S39)

Solving this transcendental equation for $\tilde{D}$ yields the depths of the
fixed points. For given values of the parameters $\mu_{t}$, $\tilde{\lambda}$,
and $\tilde{\xi}$, Eq. S39 has zero, one, or two solutions. We find that when
$\tilde{\xi}=\tilde{\xi}_{\text{bif}}$, Eq. S39 has only one solution at depth
$\tilde{D}_{\text{bif}}$ (the subscript “bif” stands for “bifurcation”), where

$\displaystyle\tilde{\xi}_{\text{bif}}=\mu_{t}+\tilde{\lambda}(1-\ln{\tilde{\lambda}})\text{
, }$ (S40) $\displaystyle\tilde{D}_{\text{bif}}=\mu_{t}+\tilde{\lambda}\text{
. }$ (S41)

For $\tilde{\xi}<\tilde{\xi}_{\text{bif}}$, there is no fixed point, while for
$\tilde{\xi}>\tilde{\xi}_{\text{bif}}$ there are two (at depths
$\tilde{D}_{1}$ and $\tilde{D}_{2}$; brown circle and star in Figs. S6 and
S7). These fixed points represent flat solutions that extend to infinity in
the $\tilde{y}$-direction (Fig. S6c).

###  S.4.2 Phase portrait

We now explore how the shape of the river depends on the model parameters
$\mu_{t}$, $\tilde{\lambda}$, and $\tilde{\xi}$ by looking at the phase
portrait, i.e. the trajectories described by Eq. 17 in the
$\tilde{D}$-$\tilde{D}^{\prime}$ space (7). In addition, we use the phase
portrait to unambiguously show that there exists an infinite, limiting river.

In Fig. S6a, we show an example of a phase portrait for given parameters
$\mu_{t}$, $\tilde{\lambda}$, and $\tilde{\xi}$. Each trajectory in the phase
portrait corresponds to a bed profile that satisfies Eq. 17 under different
boundary conditions (black lines in Fig. S6a). Depending on the boundary
conditions, Eq. 17 has multiple, qualitatively different solutions (Fig.
S6d-g). However, for a given set of parameters, $\mu_{t}$, $\tilde{\lambda}$,
and $\tilde{\xi}$, there exists only one channel-like solution of Eq. 17 with
two banks where the depth vanishes (blue line in Fig. S6a and b). Such a
solution is the only one that can represent a river transporting finite
amounts of water and sediment.

In Fig. S7, we show how the phase portrait changes as we change $\tilde{\xi}$,
but keep $\mu_{t}$ and $\tilde{\lambda}$ fixed. For
$\tilde{\xi}\rightarrow\infty$ (Fig. S7a), the river is inert and reaches a
maximum depth $\tilde{D}_{\mathrm{max},0}>\mu_{t}$, while there exists only
one fixed point at $\tilde{D}_{1}\rightarrow\mu_{t}$. Decreasing $\tilde{\xi}$
to finite values (Fig. S7b), the second fixed point appears, and the river
depth lies between the two points,
$\tilde{D}_{1}<\tilde{D}_{\mathrm{max}}<\tilde{D}_{2}$. For a particular value
$\tilde{\xi}=\tilde{\xi}_{c}$ (Fig. S7c), the river solution passes through
the second fixed point, $\tilde{D}_{\mathrm{max},c}=\tilde{D}_{2}$. Therefore,
the river becomes infinite in the $\tilde{y}-\tilde{z}$ space. The value of
$\xi_{c}$ depends on $\mu_{t}$ and $\tilde{\lambda}$, but it exists for any
value of these parameters.

Reducing $\tilde{\xi}$ below this critical value, $\tilde{\xi}_{c}$, changes
the phase portrait such that no river solution can exist anymore (Fig. S7d) —
the river solution does not exist since the two banks (solutions starting at
$\tilde{D}=0$) do not join at the center. Reducing $\tilde{\xi}$ even further,
below $\tilde{\xi}_{\text{bif}}$ (Fig. S7e), the two fixed points merge and
disappear so that there are neither river nor fixed point solutions.

We show the positions of the fixed points, $\tilde{D}_{1}$ and
$\tilde{D}_{2}$, and the river depth, $\tilde{D}_{\text{max}}$, as a function
of $\tilde{\xi}$ in the bifurcation diagram, Fig. S7f. There, we can see that,
as $\tilde{\xi}$ decreases from $\tilde{\xi}\rightarrow\infty$ to finite
values, the river depth increases until it meets with $\tilde{D}_{2}$ at
$\tilde{\xi}=\tilde{\xi}_{c}$.

The existence of $\tilde{\xi}_{c}$ implies a limiting flux $q_{s,c}$ of the
river (section Inert, active, and limiting river of the main text). By the
same reasoning, one would be tempted to conclude that the existence of
$\tilde{\xi}_{\text{bif}}$ indicates a maximum sediment flux for a flat bed
(that could be realized, for example, in a flume experiment). This is,
however, not true — when $\tilde{\xi}\rightarrow\infty$, the fixed points at
$\tilde{D}_{1}=\mu_{t}$ and $\tilde{D}_{2}\rightarrow\infty$ have sediment
fluxes $q_{s,1}=0$ and $q_{s,2}\rightarrow\infty$. Decreasing $\tilde{\xi}$ to
finite values, $q_{s,1}$ increases while $q_{s,2}$ decreases until they meet,
thereby covering the entire range of possible sediment flux values from 0 to
$\infty$. A flat bed in a flume can, therefore, carry any sediment flux, at
least in principle. In practice, however, some of these solutions may become
unstable to perturbations (6).

###  S.4.3 Inert river

In the Parker and weak transport regimes, the shape of the river is determined
by that of an inert river. We cannot find the properties of an inert river in
our model analytically; instead, to find its profile,
$\tilde{D}_{0}(\tilde{y})$, we have to solve Eq. 17 of the main text with
$\tilde{\xi}\rightarrow\infty$:

$\sqrt{\left(\tilde{D}_{0}+\tfrac{1}{3}(\tilde{D}_{0}^{3})^{\prime\prime}\right)^{2}+\tilde{D}_{0}^{\prime
2}}-\mu_{t}=0\text{ . }$ (S42)

The only parameter in this equation is $\mu_{t}$, so the inert river shape
only depends on this friction coefficient. In Fig. S8, we show how its
properties — depth, $\tilde{D}_{\text{max},0}$, width, $\tilde{W}_{0}$, non-
dimensional fluid discharge, $\tilde{Q}_{w,0}$, and aspect ratio,
$\tilde{W}_{0}/\tilde{D}_{\text{max},0}$ — depend on $\mu_{t}$ (blue lines in
Fig. S8).

We can understand how these properties generally depend on $\mu_{t}$ by
comparing them to the inert, shallow-water river whose shape can be found
analytically (black dashed lines in Fig. S8). After neglecting the momentum
diffusion, Eq. S42 becomes

$\sqrt{\tilde{D}_{\text{sw},0}^{2}+\tilde{D}_{\text{sw},0}^{\prime
2}}-\mu_{t}=0\text{ , }$ (S43)

where $\tilde{D}_{\text{sw},0}$ is the depth of an inert shallow-water river.
This equation has a simple solution (5)

$\tilde{D}_{\text{sw},0}=\mu_{t}\cos\tilde{y}\text{ . }$ (S44)

From here, we find that the shallow-water inert river has
$\tilde{D}_{\text{max},\text{sw},0}=\mu_{t}$, $\tilde{W}_{\text{sw},0}=\pi$,
$\tilde{Q}_{w,\text{sw},0}=4\mu_{t}^{3}/9$, and an aspect ratio of
$\tilde{W}_{\text{sw},0}/\tilde{D}_{\text{max},\text{sw},0}=\pi/\mu_{t}$.
These values represent bounds for our model — Fig. S8 shows that the inert
river in our model has a depth greater than $\mu_{t}$, a width greater than
$\pi$, a non-dimensional fluid discharge greater than $4\mu_{t}^{3}/9$, and an
aspect ratio greater than $\pi/\mu_{t}$. For small $\mu_{t}$, our model
approaches the shallow-water inert river, while, as we increase $\mu_{t}$, the
non-dimensional inert river in our model becomes less and less like the
shallow-water one.

###  S.4.4 Limiting river depth, $\tilde{D}_{\mathrm{max},c}$

In the main text, we noted that the value of the non-dimensional depth of an
infinite, flat river, $\tilde{D}_{\text{max},c}$, depends on parameters
$\mu_{t}$ and $\tilde{\lambda}$ (section Inert, active, and limiting river of
the main text), and we explained that for vanishing $\tilde{\lambda}$,
$\tilde{D}_{\text{max},c}$ equals the inert river depth,
$\tilde{D}_{\text{max},0}(\mu_{t})$ (section The Parker regime of the main
text). In Fig. S9, we show how
$\tilde{D}_{\text{max},c}(\tilde{\lambda},\mu_{t})$ depends on
$\tilde{\lambda}$ with $\mu_{t}$ fixed. There, we numerically show that, for
small $\tilde{\lambda}$, the depth $\tilde{D}_{\text{max},c}$ behaves as

$\tilde{D}_{\text{max},c}(\tilde{\lambda},\mu_{t})\approx\tilde{D}_{\text{max},0}(\mu_{t})+\tilde{\lambda}\quad\text{
for }\quad\tilde{\lambda}\ll 1\text{ . }$ (S45)

For experimental parameters ($\mu_{t}=0.9$ and $\tilde{\lambda}=0.02$), we
find that the relative error between $\tilde{D}_{\text{max},c}$ and Eq. S45 is
about $0.06\%$. Although we do not understand exactly why this relationship
holds, we found it to be true for all the values of $\mu_{t}$ we tested.

Figure S6:  (a) A phase portrait of Eq. 17 of the main text for $\mu_{t}=0.9$,
$\tilde{\lambda}=0.1$, and $\tilde{\xi}=1.37$. The horizontal axis is the non-
dimensional depth, $\tilde{D}$, and the vertical axis is its
$\tilde{y}$-derivative, $\tilde{D}^{\prime}$. Black dashed curves mark the
boundaries of the region of $\tilde{D}$-$\tilde{D}^{\prime}$ space on which
Eq. 17 is well-defined. Black trajectories represent solutions to Eq. 17 for
different initial conditions. The blue line represents the river solution
shown in panel (b), while the blue dot and blue squares represent the river
center and banks. The river solution is unique and acts as a separatrix
between two regions of qualitatively different solutions in the phase space
(streamwise streaks and underwater channels). Brown circle and brown star
represent the two fixed points at depths $\tilde{D}_{1}$ and $\tilde{D}_{2}$.
Brown lines are solutions ending in the second fixed point and separate
regions of the phase space with qualitatively different solutions. A constant
depth solution corresponding to the first fixed point, $\tilde{D}_{1}$, is
shown in panel (c). Different colored shadings are regions of the phase space
with qualitatively different solutions. An example of a solution from each of
these regions is marked with a colored dotted line and shown in panels
(d)-(g). Figure S7: (a)-(e) Phase portraits of Eq. 17 of the main text for
$\mu_{t}=0.9$, $\tilde{\lambda}=0.1$, and varying $\tilde{\xi}$. Notation is
the same as in Fig. S6. (a) Inert river ($\tilde{\xi}\rightarrow\infty$).
There exists only one fixed point at $\tilde{D}_{1}=\mu_{t}$. (b) Active river
($\tilde{\xi}=1.37>\tilde{\xi}_{c}$). River solution passes between the two
fixed points. (c) Infinite, limiting river
($\tilde{\xi}=1.3237\approx\tilde{\xi}_{c}$). The river solution reaches a
maximum depth at the second fixed point,
$\tilde{D}_{\mathrm{max},c}=\tilde{D}_{2}$. (d) No river solution, but fixed
points still exist ($\tilde{\xi}=1.3<\tilde{\xi}_{c}$). The blue line shows
the solutions starting at $\tilde{D}=0$. These solutions do not represent a
river since they do not join at the center. (e) No river or fixed point
solutions ($\tilde{\xi}=1.2<\xi_{\text{bif}}\approx 1.23$). The two fixed
points merge and disappear. (f) The bifurcation diagram for $\mu_{t}=0.9$ and
$\tilde{\lambda}=0.1$. The brown lines represent the depths of the fixed
points as a function of $\tilde{\xi}$ — lower branch corresponds to the first
fixed point, $\tilde{D}_{1}$ (brown circle in panels a-e), while the upper
branch corresponds to the second fixed point, $\tilde{D}_{2}$ (brown star in
panels a-e). When $\tilde{\xi}\rightarrow\infty$, the first fixed point
approaches the friction coefficient ($\tilde{D}_{1}\rightarrow\mu_{t}$), while
the second tends to infinity ($\tilde{D}_{2}\rightarrow\infty$). The two fixed
points meet at
$(\tilde{\xi},\tilde{D})=(\tilde{\xi}_{\text{bif}},\tilde{D}_{\text{bif}})$
(brown dot; $\tilde{\xi}_{\text{bif}}\approx 1.23$,
$\tilde{D}_{\text{bif}}=1$). The blue line corresponds to the river depth,
$\tilde{D}_{\text{max}}$ (blue dot in panels a-e). When
$\tilde{\xi}\rightarrow\infty$, the river approaches the inert river depth,
$\tilde{D}_{\text{max},0}$ (blue triangle). The river meets the second fixed
point at $(\tilde{\xi},\tilde{D})=(\tilde{\xi}_{c},\tilde{D}_{\text{max},c})$
(blue star; $\tilde{\xi}_{c}\approx 1.3237$, $\tilde{D}_{\text{max},c}\approx
1.21$). At this point the river is flat and infinitely wide.  Figure S8:
Properties of non-dimensional inert rivers as a function of the friction
coefficient, $\mu_{t}$. The blue lines correspond to our model (Eq. S42),
while the black dashed lines correspond to the shallow-water inert river (Eq.
S44). (a) Difference between non-dimensional depth and the friction
coefficient, $\tilde{D}_{\text{max},0}-\mu_{t}$. Black dashed line corresponds
to $\tilde{D}_{\text{max},0}-\mu_{t}=0$. (b) Non-dimensional width,
$\tilde{W}_{0}$. Black dashed line corresponds to $\tilde{W}_{0}=\pi$. (c)
Non-dimensional fluid discharge,
$\tilde{Q}_{w,0}\equiv\int\tilde{D}^{3}/3\text{d}\tilde{y}$. Black dashed line
corresponds to $\tilde{Q}_{w,0}=4\mu_{t}^{3}/9$. (d) Aspect ratio,
$\tilde{W}_{0}/\tilde{D}_{\text{max},0}$. Black dashed line corresponds to
$\tilde{W}_{0}/\tilde{D}_{\text{max},0}=\pi/\mu_{t}$.  Figure S9:  The
dependence of limiting river depth, $\tilde{D}_{\text{max},c}$, on
$\tilde{\lambda}$. The vertical axis shows the difference between the limiting
and inert river depth, $\tilde{D}_{\text{max},c}-\tilde{D}_{\text{max},0}$,
while the horizontal axis is the value of $\tilde{\lambda}$ ($\mu_{t}=0.9$ is
fixed). The blue line results from numerical solution of Eq. 17 of the main
text. The black dashed line represents the function
$\tilde{D}_{\text{max},c}-\tilde{D}_{\text{max},0}=\tilde{\lambda}$. The red
shading represents the estimate and the variability of $\tilde{\lambda}$
observed in the experiments.

##  S5 Parker regime

The Parker regime corresponds to a river with no sediment diffusion that
splits its channel into inert banks and a flat, active bottom (section The
Parker regime of the main text). Such a river has a rectangular sediment flux
profile of width $W_{T}^{(P)}$ and height $q_{s}^{(P)}$.

In Fig. S10, we show that the Parker regime is, in fact, the limit of our
model when $\lambda\rightarrow 0$. The sediment flux profile of Eq. 13 of the
main text approaches the rectangular shape of the Parker regime as $\lambda$
decreases. This is to be expected — when the sediment diffusion length scale,
$\lambda$, vanishes, the sediment flux, $q_{s}=q_{\mu}\exp[(D-\xi)/\lambda]$,
also vanishes for $D<\xi$, and becomes infinite for $D>\xi$. Therefore, to
have a meaningful solution when $\lambda\rightarrow 0$, the river bottom must
be flat with a depth, $D_{\text{max}}=\xi$. Physically, vanishing sediment
diffusion means that gravity pulls each moving grain of sediment towards the
river bottom from which it cannot escape by random collisions with the bed.
When $\lambda$ is finite, the region over which inert banks transition to the
flat bottom always has a finite size of the order $\lambda$. When the water
and sediment discharges are large, both the bank width, $W_{0}$, and the
transport width, $W_{T}$, are large too, so we can neglect the transition
region: this is the gist of the Parker regime.

The depth, sediment flux, transport width, and total width of a river in the
Parker regime are (as we explained in section The Parker regime of the main
text):

$D_{\text{max}}^{(P)}=\frac{L_{s}}{S^{(P)}}\tilde{D}_{\text{max},0}\quad\text{
, }\quad q_{s}^{(P)}=q_{\mu}(\tilde{D}_{\text{max},0}-\mu_{t})\quad\text{ ,
}\quad W_{T}^{(P)}=\frac{Q_{s}}{q_{s}^{(P)}}\quad\text{ , }\quad
W^{(P)}=\frac{L_{s}}{S^{(P)}}\tilde{W}_{0}+\frac{Q_{s}}{q_{s}^{(P)}}\text{ .
}$ (S46)

To get the above quantities, we need to find the downstream slope, $S^{(P)}$.
This follows from the fluid discharge constraint. According to Eq. 15 of the
main text, the fluid discharge is

$Q_{w}=\frac{gS^{(P)}}{\nu}\int_{-W/2}^{W/2}\frac{D^{3}}{3}\text{d}y=Q_{w,0}+\frac{gS^{(P)}}{\nu}\frac{D_{\text{max}}^{(P)3}}{3}W_{T}^{(P)}\text{
, }$ (S47)

where $Q_{w,0}\equiv\int_{-W_{0}/2}^{W_{0}/2}\frac{D_{0}^{3}}{3}\text{d}y$ is
the contribution from the inert banks. A dimensional fluid discharge, $Q_{w}$,
can be related to the non-dimensional one,
$\tilde{Q}_{w}=\int\frac{\tilde{D}^{3}}{3}\text{d}\tilde{y}$, through

$Q_{w}=\frac{gL_{s}^{4}}{\nu S^{3}}\tilde{Q}_{w}\text{ , }$ (S48)

so that Eq. S47 becomes

$Q_{w}=\frac{gL_{s}^{4}}{\nu
S^{(P)3}}\tilde{Q}_{w,0}+\frac{gS^{(P)}}{\nu}\frac{D_{\text{max}}^{(P)3}}{3}W_{T}^{(P)}\text{
, }$ (S49)

where
$\tilde{Q}_{w,0}\equiv\int_{-\tilde{W}_{0}/2}^{\tilde{W}_{0}/2}\frac{\tilde{D}_{0}^{3}}{3}\text{d}\tilde{y}$
depends only on $\mu_{t}$. Combining Eq. S49 with Eqs. S46, we retrieve Eq. 24
of the main text

$Q_{w}=\frac{gL_{s}^{4}}{\nu
S^{(P)3}}\left(\tilde{Q}_{w,0}+\frac{Q_{s}S^{(P)}\tilde{D}_{\text{max},0}^{3}}{3q_{\mu}L_{s}(\tilde{D}_{\text{max},0}-\mu_{t})}\right)\text{
. }$ (S50)

Depending on the physical parameters that enter it, Eq. S50 can have multiple
solutions for $S^{(P)}$, but there is always only one positive, real solution.

Figure S10: Sediment flux profile of a river in our model with $Q_{s}=30\text{
grains s}^{-1}$ and $Q_{w}=1\text{ l min}^{-1}$ and the varying sediment
diffusion length scale, $\lambda$. The blue lines are the numerical sediment
flux profiles for different values of $\lambda$, while the black dashed line
is the Parker regime (Eqs. S46) in the limit $\lambda\rightarrow 0$.

##  S6 Estimating sediment discharge based on river geometry

Our model predicts a link between the river geometry and its sediment load.
Here, we will use our experimental dataset to show how to estimate the
sediment discharge of a river from the shape of its cross-section.

Equation 22 of the main text relates the sediment discharge to river depth and
width, assuming we know the universal quantities, such as the limiting
sediment flux and the inert river aspect ratio ($q_{s}^{(P)}$ and
$\tilde{W}_{0}/\tilde{D}_{\mathrm{max},0}$ in Eq. 22). Although these
universal quantities follow from the theory, various approximations we made
may make them inaccurate. Moreover, for a natural river, for which the theory
is not yet available, these quantities are unknown a priori. On an ensemble of
rivers with different water and sediment discharges, $q_{s}^{(P)}$ is the
maximum allowed sediment flux while $\tilde{W}_{0}/\tilde{D}_{\mathrm{max},0}$
is the minimum allowed aspect ratio. Therefore, instead of using our
theoretical predictions, here we estimate $q_{s}^{(P)}$ and
$\tilde{W}_{0}/\tilde{D}_{\mathrm{max},0}$ as the maximum sediment flux and
minimum aspect ratio from our experimental dataset. Thus, we estimate
$q_{s}^{(P)}|_{\mathrm{exp.}}\approx 27.4\text{ grains cm}^{-1}\text{s}^{-1}$
and $\tilde{W}_{0}/\tilde{D}_{\mathrm{max},0}|_{\mathrm{exp.}}\approx 4.3$.
Using these values, we show in Fig. S11 that the sediment discharge estimated
using Eq. 22 of the main text falls within the uncertainty range of the
measurements.

Figure S11: Comparing the measured sediment discharge (x-axis) to that
estimated using Eq. 22 of the main text (y-axis). Black dashed line is a one-
to-one line.

##  S7 Effects of momentum diffusion

In section The Parker regime of the main text, we claimed that momentum
diffusion is essential for the existence of a steady-state channel with a non-
vanishing sediment discharge, as first suggested by Parker (8). In this
section, we will first illustrate this statement by comparing our model to the
shallow-water river without momentum diffusion. Then, we will investigate the
effect of momentum diffusion in detail by considering the forces acting on the
river bed.

###  S.7.1 Comparison with shallow-water rivers

The shallow-water model, which neglects the cross-stream diffusion of momentum
entirely, is defined by a first order ordinary differential equation:

$\sqrt{\tilde{D}_{\text{sw}}^{2}+\tilde{D}_{\text{sw}}^{\prime
2}}-\mu_{t}=e^{(\tilde{D}_{\text{sw}}-\tilde{\xi})/\tilde{\lambda}}\text{ . }$
(S51)

When $\tilde{\xi}\rightarrow\infty$, the right-hand side vanishes and the
model reduces to Eq. S43 for an inert shallow-water river, which has been
successfully used to predict the shape of experimental inert rivers (5). We
compare the inert shallow-water river with the experiment and our model in
Figs. S12a, c, and e. There, we can see that, while there exists some
difference between the models and the experiment, they are largely comparable
to each other.

However, when sediment transport is included ($\tilde{\xi}<\infty$), this
model is insufficient to reproduce the experiments. We compare an active
experimental river with our model and the shallow-water model in Figs. S12b,
d, and f. We can see that, while our model and the experiment seem comparable
to each other, the shallow-water model is nothing alike. This result may seem
somewhat paradoxical — although we expect that the shallow-water approximation
should work better for a wider river, we find that the wider the river, the
worse the shallow-water approximation is. We resolve this apparent paradox in
the next section.

###  S.7.2 Momentum diffusion in the Parker regime

In this section, we look closely at the components of the force ratio, $\mu$,
to understand the role of momentum diffusion. We explain that momentum
diffusion controls the sediment transport through a matching condition at the
interface of the river banks and the flat bottom — a condition that does not
depend on the width of the river.

Figure S13a shows a river in the Parker regime (in non-dimensional
coordinates, for convenience) — it is sharply split into curved, inert banks,
and a flat, active bottom. Figure S13b shows how the various components of the
force ratio, $\mu$, depend on the position within this river. These components
are

$\displaystyle\mu_{\text{sw}}\equiv\tilde{D}\text{ , }$ (S52)
$\displaystyle\mu_{\text{md}}\equiv\frac{1}{3}(\tilde{D}^{3})^{\prime\prime}\text{
, }$ (S53) $\displaystyle\mu_{\text{g}}\equiv\tilde{D}^{\prime}\text{ , }$
(S54)
$\displaystyle\mu^{2}=(\mu_{\text{sw}}+\mu_{\text{md}})^{2}+\mu_{\text{g}}^{2}\text{
, }$ (S55)

where $\mu_{\text{sw}}$ is the shallow-water contribution, $\mu_{\text{md}}$
is the momentum diffusion contribution, and $\mu_{\text{g}}$ is the gravity
contribution to the total force ratio. The banks are kept at threshold by the
combined action of the shallow-water stress, momentum diffusion, and gravity,
while the bottom only feels the shallow-water stress. At the point where banks
and flat bottom connect the banks have a depth
$\tilde{D}=\tilde{D}_{\text{max}}$, a slope $\tilde{D}^{\prime}=0$, and a
second derivative $\tilde{D}^{\prime\prime}=-\tilde{\kappa}<0$, where
$\tilde{\kappa}$ is the curvature of the banks at their deepest point. The
components of the force ratio are, correspondingly,
$\mu_{\text{sw}}=\tilde{D}_{\text{max}}$, $\mu_{\text{g}}=0$,
$\mu_{\text{md}}=-\tilde{D}_{\text{max}}^{2}\tilde{\kappa}<0$. The threshold
condition for the banks is, thus, $\mu_{\text{sw}}+\mu_{\text{md}}=\mu_{t}$,
which means that the depth must be greater than the friction coefficient,
$\tilde{D}_{\text{max}}=\mu_{\text{sw}}=\mu_{t}-\mu_{\text{md}}>\mu_{t}$. The
bottom feels only the shallow-water stress, and, therefore, must be above the
threshold,
$\mu|_{\text{bottom}}=\mu_{\text{sw}}=\tilde{D}_{\text{max}}>\mu_{t}$.

Therefore, the absence of sediment transport on the banks translates into
excess stress on the flat bottom. The sediment flux on the bottom,
$q_{s}/q_{\mu}=\mu|_{\text{bottom}}-\mu_{t}=-\mu_{\text{md}}$, corresponds to
the jump in the momentum diffusion contribution due to the sudden
disappearance of curvature. Thus, the sediment flux is driven by the boundary
condition at the junction of the curved banks with the flat bottom.

Written in terms of depth and its derivatives, the sediment flux is
$q_{s}/q_{\mu}=\tilde{D}_{\text{max}}^{2}\tilde{\kappa}$. The curvature,
$\tilde{\kappa}\equiv-\tilde{D}^{\prime\prime}$, scales roughly as
$\tilde{\kappa}\sim\tilde{D}_{\text{max},0}/\tilde{W}_{0}^{2}$, so the
sediment flux scales as
$q_{s}/q_{\mu}\sim\tilde{D}_{\text{max},0}^{3}/\tilde{W}_{0}^{2}$. Taking the
depth, $\tilde{D}_{\text{max},0}$, to be of order one, the sediment flux is
inversely proportional to the square of the aspect ratio of an inert river,
$W_{0}^{2}/D_{\text{max},0}^{2}$. Parker (8) likewise found that the distance
to threshold is inversely proportional to the aspect ratio squared — a
signature of the cross-stream diffusion of momentum.

Without cross-stream momentum diffusion, there can be no sediment transport in
the limit $\lambda\rightarrow 0$ in our model. When $\lambda$ is small but
non-vanishing, the maximal sediment flux is
$q_{s,c}/q_{\mu}\approx\tilde{D}_{\text{max},0}(\mu_{t})-\mu_{t}+\tilde{\lambda}$,
where we used Eq. S45 for the depth of the infinite, limiting river,
$\tilde{D}_{\text{max},c}$. The first term in this equation,
$\tilde{D}_{\text{max},0}(\mu_{t})-\mu_{t}$, is the contribution to sediment
transport from the momentum diffusion, while the second term,
$\tilde{\lambda}$, is the contribution from the sediment diffusion. For
experimental parameters ($\mu_{t}=0.9$ and $\tilde{\lambda}=0.02$), we thus
find that momentum diffusion is responsible for about 90% of the sediment
transport.

Figure S12: River profiles in the experiment (panels (a) and (b)), our model
(panels (c) and (d)), and shallow-water approximation (panels (e) and (f)).
Profiles in the left column (panels (a), (c), and (e)) are inert rivers, while
profiles in the right column (panels (b), (d), and (f)) are active rivers with
$Q_{s}=24\text{ grains s}^{-1}$. Red lines in panels (b), (d), and (f)
represent the sediment flux profile where the zero, $q_{s}=0$, is vertically
shifted in the plots to coincide with the river bottom. While inert rivers in
both our model and the shallow-water model are comparable with experiments,
momentum diffusion is necessary to capture the shape of active rivers.  Figure
S13: Explanation of the mechanism of sediment transport by momentum diffusion.
(a) Half a river in our model for $\lambda\rightarrow 0$ and $\mu_{t}=0.9$.
The river sharply splits into curved, inert banks and a flat, active bottom.
(b) Components of the force acting on the bed as a function of the position
along the bed. The blue line represents the total force ratio,
$\mu=\sqrt{(\mu_{\text{sw}}+\mu_{\text{md}})^{2}+\mu_{\text{g}}^{2}}$. The
dashed brown line represents the shallow-water contribution to the stress,
$\mu_{\text{sw}}=\tilde{D}$. The solid brown line represents the momentum
diffusion contribution to the stress,
$\mu_{\text{md}}=(\tilde{D}^{3})^{\prime\prime}/3$. The dotted brown line
represents the gravity contribution, $\mu_{\text{g}}=\tilde{D}^{\prime}$. On
the banks, all three contributions are non-vanishing and keep the banks at the
threshold, $\mu=\mu_{t}$. On the bottom, only the shallow-water contribution
exists and the force ratio is above threshold, $\mu>\mu_{t}$. The shallow-
water and gravity contributions transition continuously from the banks to the
bed, but the momentum diffusion contribution experiences a jump which
corresponds to a jump in the force ratio, $\mu$. This jump drives sediment
transport.

## References

  * Abramian et al. (2020) A Abramian, O Devauchelle, and E Lajeunesse.  Laboratory rivers adjust their shape to sediment transport.  _Physical Review E_ , 102(5):053101, 2020.
  * Charru et al. (2004) François Charru, Hélene Mouilleron, and Olivier Eiff.  Erosion and deposition of particles on a bed sheared by a viscous flow.  _Journal of Fluid Mechanics_ , 519:55–80, 2004.
  * Abramian (2018) Anaïs Abramian.  _Self-organization of sediment transport in alluvial rivers_.  PhD thesis, 2018.
  * Devauchelle et al. (2021) Olivier Devauchelle, Eric Lajeunesse, and Popović Predrag.  Viscous transfer of momentum across a shallow laminar flow.  _submitted to Journal of Fluid Mechanics_ , 2021.
  * Seizilles et al. (2013) G Seizilles, O Devauchelle, E Lajeunesse, and F Métivier.  Width of laminar laboratory rivers.  _Physical Review E_ , 87(5):052204, 2013.
  * Abramian et al. (2019) Anaïs Abramian, Olivier Devauchelle, and Eric Lajeunesse.  Streamwise streaks induced by bedload diffusion.  2019\.
  * Strogatz (2018) Steven H Strogatz.  _Nonlinear dynamics and chaos with student solutions manual: With applications to physics, biology, chemistry, and engineering_.  CRC press, 2018.
  * Parker (1978) Gary Parker.  Self-formed straight rivers with equilibrium banks and mobile bed. part 2. the gravel river.  _Journal of Fluid mechanics_ , 89(1):127–146, 1978.


11institutetext: University of Missouri, USA
11email: {xwf32, yl235<EMAIL_ADDRESS>

# Fast Point Voxel Convolution Neural Network with Selective Feature Fusion
for Point Cloud Semantic Segmentation

Xu Wang Equal contribution Yuyan Li* Ye Duan

###### Abstract

We present a novel lightweight convolutional neural network for point cloud
analysis. In contrast to many current CNNs which increase receptive field by
downsampling point cloud, our method directly operates on the entire point
sets without sampling and achieves good performances efficiently. Our network
consists of point voxel convolution (PVC) layer as building block. Each layer
has two parallel branches, namely the voxel branch and the point branch. For
the voxel branch specifically, we aggregate local features on non-empty voxel
centers to reduce geometric information loss caused by voxelization, then
apply volumetric convolutions to enhance local neighborhood geometry encoding.
For the point branch, we use Multi-Layer Perceptron (MLP) to extract fine-
detailed point-wise features. Outputs from these two branches are adaptively
fused via a feature selection module. Moreover, we supervise the output from
every PVC layer to learn different levels of semantic information. The final
prediction is made by averaging all intermediate predictions. We demonstrate
empirically that our method is able to achieve comparable results while being
fast and memory efficient. We evaluate our method on popular point cloud
datasets for object classification and semantic segmentation tasks.

###### Keywords:

Point Cloud Semantic Segmentation Deep Learning

##  1 Introduction

Deep learning in 3D point cloud analysis has received increasing attention
with the rising trend of Virtual Reality and 3D scene understanding
applications, etc. Existing approaches have made great progresses in tasks
such as point cloud classification [28] and point cloud semantic segmentation
[3, 1]. One fundamental issue to be tackled with in point cloud analysis is
the representation of unstructured point clouds. Some early methods discretize
point clouds into regular volumetric grids which can be directly fed into
standard 3D CNNs. However, two main problems coupled with this volumetric
representation are information loss and huge memory consumption. A high
resolution voxel grid leads to expensive computation cost, while a low
resolution inevitably suffers from information loss during voxelization
procedure.

Figure 1: Illustration of our proposed network. For an input 3D data, we pass
it through a sequence of point-voxel convolution layers(PVC). $G$ denotes grid
resolution. Outputs from each PVC layer are concatenated together to form a
global feature. MLP is used for feature dimension reduction. This global
feature is concatenated with output from each PVC layer and passed through a
channel attention module, which re-weights the features of all channels and
increases feature disciminability. The final prediction is the average of all
auxiliary predictions.

To address the problems mentioned above, another big stream is to directly
consume sparse point clouds. The pioneer work is PointNet proposed by Qi et
al. [17]. PointNet is able to process unordered point cloud inputs with
permutation invariance using a sequence of multi-layer perceptron(MLP). The
subsequent work PointNet++[18] achieves better performance by proposing a
hierarchical network that encodes local neighborhood information. Based on
PointNet++[18], a great number of networks[32, 31, 25] with more advanced
local feature aggregation techniques are introduced. Apart from MLP-based
methods, some recent works propose kernel-based approaches to mimic standard
convolution[24, 26, 30, 2]. In general, point-based approaches suffer from
point sampling scalability, neighbor point searching efficiency, point density
inconsistency issues. Most recently, Liu et al. [13] attempt to design a
point-voxel CNN that represents 3D data as points to reduce memory footprint,
and leverages voxel-based convolution to capture neighborhood features. This
network is able to achieve reasonable performance with low memory usage and
fast training/inference speed.

In this paper, we propose a novel CNN architecture that is well-balanced
between efficiency and accuracy. Inspired by point-voxel CNN[13], we construct
our network using point-voxel layer that takes advantages of both sparse point
representation and volumetric convolution. Our point-voxel layer consists of
two parallel branches, a voxel-based branch which aggregates local neighboring
features, and a point-based branch which maintains fine-grained point-wise
features. During discretization in voxel branch, we aggregate neighboring
features on non-empty voxel centers and use standard 3D convolutions to
enhance local feature encoding. Voxel features are propagated back to point
domain through devoxelization. Outputs from point and voxel branches are fused
self-adaptively via a feature selection module(FSM), which learns channel-wise
attention for both branches.

Most of the existing studies rely heavily on point sampling strategy to avoid
expensive computation cost as network goes deeper. However, point sampling
cannot always retain the fine-detailed features for every point. Details of
points are discarded as a trade-off for larger receptive field and processing
speed. In our network, we only use a small number of point-voxel layers
(default is 4 layers) that are carefully designed with effective feature
encoding modules to facilitate processing efficiency. Supervision is applied
on outputs from all layers to enforce semantic information learning. Though no
point sampling is conducted, our network is able to remain lightweight, and
effective to process large-scale point clouds. A visualization of our proposed
network is shown in Figure 1. We evaluate the performance and efficiency of
our proposed model for object classification, object part segmentation and
indoor scene semantic segmentation tasks (see Section 4).

##  2 Related Work

###  2.1 Volumentric Representation

Some early deep learning approaches transformed point clouds into 3D voxel
structure and convolve it with standard 3D kernels. VoxNet [15] and subsequent
works [28, 14, 23] discretized point cloud into a 3D binary occupancy grid.
The occupancy grid is fed to a CNN for object proposal and classification.
These voxel-based methods suffered from high memory consumption due to the
waste of computation on empty spaces. OctNet[19, 22] proposed adaptive
representation using octree structure to reduce memory consumption. Recent
researches [8, 4] introduced approaches to process high dimensional data and
apply sparse convolution only on non-empty voxels. In general, volumetric
methods preserve neighborhood information of point clouds, enable regular 3D
CNN applications, but suffer from significant discretization artifacts.

###  2.2 Point-based Representation

Point-wise models such as PointNet [17] and PointNet++ [18] directly operates
on point clouds. The former used MLPs to extract point-wise features and
permutation- invariant max pooling operation to obtain a global feature. The
latter built a hierarchical architecture that incorporates point downsampling
and local structure aggregation strategies. Inspired by PointNet [17] and
PointNet++ [18], many recent works propose advanced local feature learning
modules. For example, PointWeb [32] built a dense fully connected web to
explore local context, and used an Adaptive Feature Adjustment module for
feature refinement. GACNet [25] proposed to selectively learn distinctive
features by dynamically assigning attention weights to neighbouring points
based on spatial positions and feature differences. ShellNet [31] built a
model with several layers of ShellConv, and solved point ambiguity by
constructing concentric shells and applying 1D convolution on ordered shells.
Derived from point-based methods, some recent works define explicit kernels
for point convolution. KCNet [20] developed a kernel correlation layer to
compute affinities between each point’s K nearest neighbors and a predefined
set of kernel points. Local features are acquired by graph pooling layers.
SpiderCNN [30] designed a family of Taylor polynomial kernels to aggregate
neighbor features. PointCNN [12] introduced $\chi$-transformation to exploit
the canonical order of points. PCNN [26] built a network using parametric
continuous convolutional layers. SPH3D [11] used spherical harmonic kernels
during convolution on quantized space to identify distinctive geometric
features. KPConv [24] defined rigid and deformable kernel points for local
geometry encoding based on the Euclidean space relations between kernel point
and neighborhood supporting points.

###  2.3 Efficiency of Current Models

When processing large-scale point clouds, efficiency is one of the fundamental
measurements to evaluate models. Most of the point-based methods utilized
point sampling to improve efficiency. However, it is non-trivial to choose an
effective point sampling method. For example, Farthest Point Sampling(FPS)
which is widely adopted in [18, 32, 31], has $O(NlogN)$ computation
complexity, meaning it does not have good scalability. Random point sampling
as used in RandLA-Net[9], has $O(1)$ time complexity, but random point
sampling cannot be invariant to point densities and key information might be
discarded. Other approaches manage to incorporate hybrid representations to
avoid the redundant computing and storing of more useful spatial information.
A recent work Grid-GCN [29] proposed a novel method which facilitates grid
space structuring and provides more complete coverage of the point cloud. This
method is able to handle massive points with fast speed and good scalability.
Point-Voxel CNN [13] is most related to our method. This work combines fine-
grained point features with coarse-grained voxel features with speedup and low
memory consumption. Compared to their work, our method builds PVC layers with
multi-resolution voxels, and incorporates more accurate local feature
aggregation which reduce information loss artifacts.

##  3 Method

We build a deep architecture with a sequence of point-voxel convolution(PVC)
layers. In this section, we introduce the details of our PVC layer, including
voxelization, local feature aggregation, devoxelization, selective feature
fusion, and deep supervision.

###  3.1 Point Voxel Convolution

####  3.1.1 Voxelization and Local Aggregation

The purpose of voxel branch is to encode contextual information through
volumetric convolution. As aforementioned, information loss from the process
of discretization is inevitable. Introducing large voxel grids reduces the
loss, but burdens the network with huge computation overhead. In our design,
we opt to use low-resolution volumetric grid, and mitigate information loss by
effective local feature aggregation.

As the scale of point clouds varies, we first normalize the input data into a
bounding box. For voxelization, we quantize point cloud by calculating voxel
coordinates $(u,v,w)\in\mathbb{N}$ from point
coordinates$(x,y,z)\in\mathbb{R}$.

$u=floor((x-x_{min})/g_{x}),\;v=floor((y-y_{min})/g_{y}),\;w=floor((z-z_{min})/g_{z})$
(1)

where $g_{x},g_{y},g_{z}$ is the grid length of $x,y,z$ axis respectively:

$g_{x}=(x_{max}-x_{min})/G,\;g_{y}=(y_{max}-y_{min})/G,\;g_{z}=(z_{max}-z_{min})/G$
(2)

where $G\in\mathbb{N}$ is the grid resolution of this PVC layer.

Given the voxelized point cloud, we calculate the center location $P_{c}$ of
every voxel. To accelerate processing speed, we locate non-empty voxels and
aggregates local features only on these voxels. For a non-empty voxel
$(u,v,w)$, we use the cell center as query position and gather K neighbors
through K-nearest-neighbor(KNN). We adopt the dilated point convolution
strategy as proposed in [6], in which $K\cdot n$ nearest neighbors and every
$n-th$ neighbor is selected. The feature $f$ for voxel $(u,v,w)$ is the
weighted summation of all K neighbors.

$f_{u,v,w}=\sum_{k}^{K}l_{k}*f_{k}$ (3)

where $l_{k}\in\mathbb{R}^{1\times C}$ is the weight,
$f_{k}\in\mathbb{R}^{1\times C}$ is the feature of k-th neighbor. Inspired by
[25, 26], we use self-attention mechanism to learn weight of different
neighboring points:

$l_{k}=\sigma(\mathcal{G}(\Delta_{P_{k}},f_{k}))$ (4)

where $\Delta_{P_{k}}$ is the normalized neighbor point coordinates,
$\Delta_{P_{k}}=P_{k}-P_{c}$. $\mathcal{G}(*)$ takes the concatenation of
$\Delta_{P_{k}}$ and $f_{k}$ as input and models the attention weight.
$\sigma$ represents ReLU activation function. Through neighbor feature
aggregation, we collect useful information and store it in the voxels. Next,
we apply 3D convolutions on the voxel grid, as an enhancement of local
neighborhood learning.

####  3.1.2 Devoxelization

To allow feature fusion from two branches, we propagate voxel features back to
point domain based on their voxel coordinates. Taking efficiency into account,
we assign voxel feature $f_{(}u,v,w)$ to all points that fall into this voxel.
We observe that in our experiment since point-wise features are carried all
along by point branch, fusing voxel feature with point feature would be
effective to discriminate individual point.

####  3.1.3 Point and Voxel Feature Fusion

As illustrated in Figure 2, for point branch, we use MLPs to extract point-
wise features. While voxel features encode local neighborhood, point branch is
able to carry fine-detailed per-point features. Next, we incorporate a feature
selection module to correlate features. First, we use element-wise summation
to fuse features from point and voxel branches:

$f^{\prime}=f_{p}+f_{v}$ (5)

where $f_{p}\in\mathbb{R}^{N\times C}$ and $f_{v}\in\mathbb{R}^{N\times C}$
are features from point and voxel branch respectively. Then a global average
pooling is applied to squeeze $N$ point to one compact point feature. Fully
connected layer with non-linearity is used to provide guidance for feature
selection:

$S=\sigma(\mathcal{F}_{gp}(f^{\prime})\cdot W_{fc})$ (6)

where $\sigma$ is the ReLU activation function, $\mathcal{F}_{gp}$ is the
global average pooling, $W_{fc}\in\mathbb{R}^{C\times d}(d=C/4)$ is the
learnable weight for fully connected layer. Two separate fully connected
layers are applied to increase channel dimensions for $S$ and produce soft
attention vector $S_{p}\in\mathbb{R}^{1\times C}$ and
$S_{v}\in\mathbb{R}^{1\times C}$.

$S_{p}=S\cdot W_{1},\;\;S_{v}=S\cdot W_{2}$ (7)

where $W_{1}\in\mathbb{R}^{d\times C}$ and $W_{2}\in\mathbb{R}^{d\times C}$
are learnable weights. We adopt the softmax mechanism on $S_{p}$ and $S_{v}$
to adaptively select features.

$S_{p,c}=\frac{e^{S_{p},c}}{e^{S_{p},c}+e^{S_{v},c}},\;S_{v,c}=\frac{e^{S_{v},c}}{e^{S_{p},c}+e^{S_{v},c}}$
(8)

where $S_{p,c}$ and $S_{v,c}$ are soft attention vector for point and voxel
feature at $c^{th}$ channel. The fused feature at $c^{th}$ channel can be
calculated as follows:

$F_{fused,c}=S_{p,c}\odot F_{p,c}+S_{v,c}\odot F_{v,c}$ (9)

Therefore, FSM adjusts channel-wise weight for different branches, and outputs
the fused feature adaptively.

Figure 2: Illustration of our PVC layer (top) and Feature Selection Module
(bottom).

###  3.2 Deep Supervision

As shown in Figure 1, we build our network with several PVC layers
sequentially. Grid resolution decreases while output feature channel increases
from shallow to deep. Different PVC layers extract different levels of
semantic information. We add supervision on each PVC layer output to enforce
different levels of semantic feature learning. Similar strategy is adopted by
[21] for multi-scale medical image segmentation. In detail, we concatenate
outputs from all PVC layers and use MLPs to produce a compact feature. Serving
as a global guidance, this feature map is concatenated with the output of each
PVC layer, then pass through a channel attention module to enhance feature
representation of specific semantics. This attention module, inspired by [7],
aggregates weighted features of all the channels into the original features,
and models discriminability between channels. The output from channel
attention module produce an auxiliary loss. We add up all the losses and
average the prediction probabilities for the final prediction. We show that
incorporating channel attention module in our network boosts performances and
not necessarily slow down inference speed(see Section 4.4).

Figure 3: Channel attention module

##  4 Experiment

In this section, we evaluate the performance of our method in point cloud
learning tasks including object classification, object part segmentation, and
indoor scene segmentation. For parameter settings, our network has a total of
four PVC layers. From the first to the last layer, grid size is set to
$G_{0}=32,G_{1}=16,G_{2}=8,G_{3}=4$ respectively. Number of neighbors for KNN
search is 32. Dilation step is $n_{0}=1,n_{1}=2,n_{2}=4,n_{3}=8$. Our method
is implemented with PyTorch [16] and run on a Nvidia TitanXP GPU. Batch size
for training is set to be 16. We use the Adam optimizer [10] with default
settings. The learning rate is initialized as 0.001 and decays by a rate of
0.5 every 20 epochs. Object classification converges around 100 epochs, part
segmentation converges at 80 epochs, and indoor scene segmentation converges
at 120 epochs. The full version of our network has $C=64$ feature channels for
the first PVC layer. Feature channel for l-th layer is $C\times 2^{l-1}$.

###  4.1 Object Classification

####  4.1.1 Datasets and Implementation Details

We evaluate our network on ModelNet40 [28] for 3D object classification task.
ModelNet40 contains 12311 meshed CAD models from 40 categories. The dataset
has 9843 objects for training and 2468 objects for testing. We prepare the
dataset following PointNet [17] conventions. We random sample $N$ points and
only use normal as input feature. As illustrated in Figure 1, the output from
all PVC layers are concatenated to form a global feature. Max pooling followed
by a fully connected layer are then used to produce a classification score.

####  4.1.2 Evaluation

We train three versions of classification network with variance on number of
points and number of feature channels. Shown in Table 1, our method has a good
balance between performance and efficiency.

Table 1: Results of ModelNet40 [28] classification. $1\times C$ represents our full size model. Method | Input Data | OA | Latency(ms)
---|---|---|---
PointNet [17] | $8\times 1024$ | 89.2 | $\mathbf{15}$
Ours($0.5\times C$) | $8\times 1024$ | $\mathbf{91.7}$ | 20
PointNet++ [18] | $8\times 1024$ | 91.9 | 27
DGCNN [27] | $8\times 1024$ | 91.9 | 27
Ours($0.75\times C$) | $8\times 1024$ | $\mathbf{92.3}$ | $\mathbf{26}$
Grid-GCN(full) [29] | $16\times 1024$ | $\mathbf{93.1}$ | 42
Ours($1\times C$) | $8\times 2048$ | 92.5 | $\mathbf{35}$

###  4.2 Shape Segmentation

####  4.2.1 Data and Implementation Details

We conduct experiment on ShapeNetPart[3] for shape segmentation. ShapeNetPart
is a collection of 16681 point clouds (14006 for training, 2874 for testing)
from 16 categories, each annotated with 2-6 labels. Input features are normals
only, while point coordinates $(x,y,z)$ are incorporated in the network for
voxelization and local aggregation. We random sample 2048 points for training
and use the original points for testing.

####  4.2.2 Evaluation

A comparison of our method and previous approaches is listed in Table 2. We
report our result with mean instance IoU. We train three versions of our
method, a compact network with $0.5\times C$ feature channels, a medium size
network with $0.75\times C$ feature channels, and a full size network. Our
compact network is able to achieve the same results as PointNet++[18], with
$2\times$ speedup, and 0.7G less memory consumption. Comparing with DGCNN our
compact method is $2\times$ faster and only half of its memory usage.
Comparing with PV CNN[13], our method achieves comparable results.

Table 2: Results of object part segmentation on ShapeNetPart [3]. Our method achieves comparable performance with fast inference speed, and low GPU consumption. Method | Input Data | InstanceIoU | Latency (ms) | GPU usage (G)
---|---|---|---|---
PointNet [17] | $8\times 2048$ | 83.7 | $\mathbf{22}$ | 1.5
3D-Unet [5] |  Volume($8\times 96^{3}$) | 84.6 | 682 | 8.8
PointNet++ [18] | $8\times 2048$ | 85.1 | 78 | 2.0
DGCNN [27] | $8\times 2048$ | 85.1 | 88 | 2.4
PV CNN($0.5\times C$) [13] | $8\times 2048$ | $\mathbf{85.5}$ | 22 | $\mathbf{1.0}$
Ours($0.5\times C$) | $8\times 2048$ | $\mathbf{85.5}$ | 32 | 1.3
PointCNN [12] | $8\times 2048$ | 86.1 | 136 | 2.5
PV CNN($1\times C$) [13] | $8\times 2048$ | 86.2 | $\mathbf{51}$ | $\mathbf{1.6}$
Ours($1\times C$) | $8\times 2048$ | $\mathbf{86.3}$ | 68 | 2.3

###  4.3 Indoor Scene Segmentation

####  4.3.1 Data and Implementation Details

We conduct experiments on S3DIS [1] for large-scale indoor scene segmentation.
S3DIS [1] is a challenging dataset which consists of point clouds collected
from six areas. Following the convention [17, 12], we leave out area 5 for
testing purpose. For data preparation, we split rooms into $2m\times 2m$
blocks, with $0.5m$ padding along each side $(x,y)$. These context points do
not involve in neither loss computation nor prediction during testing. We use
color as input feature, point coordinates are incorporated in the network for
voxelization and local aggregation. At training time, we random sample 8192
points from block data, and use original points at testing time.

To demonstrate the great potential of our proposed network, we also design
experiments which replace regular 3D convolution layers with sparse 3D
convolutions. We adopt Minkowski [4] sparse convolution in our experiment.
Sparse convolution enables our network to process high-dimensional data with
further speedup and reduce computation load on GPU.

####  4.3.2 Evaluation

A list of comparison of our method and previous approaches is shown in Table
3. We also train three versions for indoor scene segmentation. Compared with
PV-CNN++ [13], our method is faster while able to achieve $2.9\%$ higher mIoU
score.

Table 3: Results of indoor scene segmentation on S3DIS[1], evaluated on Area 5. We report the result using mean Intersection-Over-Union(mIoU) metric. Compared with previous methods, our method is able to achieve top-ranking results while being lightweight and fast. Method | Input Data | mIoU | Latency (ms) | GPU usage (G)
---|---|---|---|---
PointNet [17] | $8\times 4096$ | 43.0 | $\mathbf{21}$ | 1.0
PointNet++ [18] | $8\times 4096$ | 52.3 | - | -
3D-Unet [5] |  Volume($8\times 96^{3}$) | 55.0 | 575 | 6.8
DGCNN [27] | $8\times 4096$ | 48.0 | 178 | 2.4
PointCNN [12] | $16\times 2048$ | 57.3 | 282 | 4.6
PV CNN++($0.5\times C$) [13] | $4\times 8192$ | 57.6 | 41 | $\mathbf{0.7}$
PV CNN++($1\times C$) [13] | $4\times 8192$ | 59.0 | 70 | 0.8
Grid-GCN(full) [29] | $8\times 4096$ | 57.8 | 26 | -
Ours($0.5\times C$) | $4\times 8192$ | 60.2 | 34 | 1.4
Ours($0.75\times C$) | $4\times 8192$ | 60.8 | 51 | 1.8
Ours($1\times C$) | $4\times 8192$ | $\mathbf{61.7}$ | 71 | 2.0
Ours($1\times C$, sparse) | $4\times 8192$ | 61.4 | 42 | 0.9

###  4.4 Ablation Study

To show the effectiveness of our proposed method, we gradually add a component
while keeping the rest unchanged. To see the gain of each component, we train
a baseline network consists of four PVC layers. The baseline does not use
local feature aggregation, instead an averaged feature of all points fall into
the same voxel is taken. Fusion method is summation only. And no channel
attention module(CAM) is used for prediction. Final prediction is directly
produced from global feature without deep supervision. Experiments are
conducted on ShapeNetPart[3]. The baseline model is the compact version
($0.5\times C$). From Table 4, we can see that each component is able to boost
baseline method without necessarily increase latency too much.

Table 4: Ablation studies on ShapeNetPart[3]. | mIoU | gain | Latency(ms)
---|---|---|---
Baseline | 84.6 | - | 26
w/ Local Aggregation | 84.8 | +0.2 | 26
w/ FSM | 85.2 | +0.6 | 28
w/ CAM | 84.9 | +0.3 | 30
Full(Local Aggregation + FSM + CAM) | 85.5 | +0.9 | 32

##  5 Conclusion

In this work, we propose a novel approach for fast and effective 3D point
cloud learning. We designed a lightweight network that can incorporate both
fine-grained point features and multi-scale local neighborhood information. We
introduce feature selection module and deep supervision into our network for
performance improvement. Experimental results on several point cloud datsets
demonstrate that our method achieves comparable results while being fast and
memory efficient.

## References

  * [1] Armeni, I., Sax, S., Zamir, A.R., Savarese, S.: Joint 2d-3d-semantic data for indoor scene understanding. arXiv preprint arXiv:1702.01105 (2017)
  * [2] Boulch, A.: Convpoint: Continuous convolutions for point cloud processing. Computers & Graphics (2020)
  * [3] Chang, A.X., Funkhouser, T., Guibas, L., Hanrahan, P., Huang, Q., Li, Z., Savarese, S., Savva, M., Song, S., Su, H., et al.: Shapenet: An information-rich 3d model repository. arXiv preprint arXiv:1512.03012 (2015)
  * [4] Choy, C., Gwak, J., Savarese, S.: 4d spatio-temporal convnets: Minkowski convolutional neural networks. In: Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition. pp. 3075–3084 (2019)
  * [5] Çiçek, Ö., Abdulkadir, A., Lienkamp, S.S., Brox, T., Ronneberger, O.: 3d u-net: learning dense volumetric segmentation from sparse annotation. In: International conference on medical image computing and computer-assisted intervention. pp. 424–432. Springer (2016)
  * [6] Engelmann, F., Kontogianni, T., Leibe, B.: Dilated point convolutions: On the receptive field of point convolutions. arXiv preprint arXiv:1907.12046 (2019)
  * [7] Fu, J., Liu, J., Tian, H., Li, Y., Bao, Y., Fang, Z., Lu, H.: Dual attention network for scene segmentation. In: Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition. pp. 3146–3154 (2019)
  * [8] Graham, B., Engelcke, M., van der Maaten, L.: 3d semantic segmentation with submanifold sparse convolutional networks. In: Proceedings of the IEEE conference on computer vision and pattern recognition. pp. 9224–9232 (2018)
  * [9] Hu, Q., Yang, B., Xie, L., Rosa, S., Guo, Y., Wang, Z., Trigoni, N., Markham, A.: Randla-net: Efficient semantic segmentation of large-scale point clouds. arXiv preprint arXiv:1911.11236 (2019)
  * [10] Kingma, D.P., Ba, J.: Adam: A method for stochastic optimization. arXiv preprint arXiv:1412.6980 (2014)
  * [11] Lei, H., Akhtar, N., Mian, A.: Spherical kernel for efficient graph convolution on 3d point clouds. arXiv preprint arXiv:1909.09287 (2019)
  * [12] Li, Y., Bu, R., Sun, M., Wu, W., Di, X., Chen, B.: Pointcnn: Convolution on x-transformed points. In: Advances in neural information processing systems. pp. 820–830 (2018)
  * [13] Liu, Z., Tang, H., Lin, Y., Han, S.: Point-voxel cnn for efficient 3d deep learning. In: Advances in Neural Information Processing Systems. pp. 963–973 (2019)
  * [14] Maturana, D., Scherer, S.: 3d convolutional neural networks for landing zone detection from lidar. In: 2015 IEEE international conference on robotics and automation (ICRA). pp. 3471–3478. IEEE (2015)
  * [15] Maturana, D., Scherer, S.: Voxnet: A 3d convolutional neural network for real-time object recognition. In: 2015 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS). pp. 922–928. IEEE (2015)
  * [16] Paszke, A., Gross, S., Chintala, S., Chanan, G., Yang, E., DeVito, Z., Lin, Z., Desmaison, A., Antiga, L., Lerer, A.: Automatic differentiation in pytorch (2017)
  * [17] Qi, C.R., Su, H., Mo, K., Guibas, L.J.: Pointnet: Deep learning on point sets for 3d classification and segmentation. In: Proceedings of the IEEE conference on computer vision and pattern recognition. pp. 652–660 (2017)
  * [18] Qi, C.R., Yi, L., Su, H., Guibas, L.J.: Pointnet++: Deep hierarchical feature learning on point sets in a metric space. In: Advances in neural information processing systems. pp. 5099–5108 (2017)
  * [19] Riegler, G., Osman Ulusoy, A., Geiger, A.: Octnet: Learning deep 3d representations at high resolutions. In: Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition. pp. 3577–3586 (2017)
  * [20] Shen, Y., Feng, C., Yang, Y., Tian, D.: Mining point cloud local structures by kernel correlation and graph pooling. In: Proceedings of the IEEE conference on computer vision and pattern recognition. pp. 4548–4557 (2018)
  * [21] Sinha, A., Dolz, J.: Multi-scale guided attention for medical image segmentation. arXiv preprint arXiv:1906.02849 (2019)
  * [22] Tatarchenko, M., Dosovitskiy, A., Brox, T.: Octree generating networks: Efficient convolutional architectures for high-resolution 3d outputs. In: Proceedings of the IEEE International Conference on Computer Vision. pp. 2088–2096 (2017)
  * [23] Tchapmi, L., Choy, C., Armeni, I., Gwak, J., Savarese, S.: Segcloud: Semantic segmentation of 3d point clouds. In: 2017 international conference on 3D vision (3DV). pp. 537–547. IEEE (2017)
  * [24] Thomas, H., Qi, C.R., Deschaud, J.E., Marcotegui, B., Goulette, F., Guibas, L.J.: Kpconv: Flexible and deformable convolution for point clouds. In: Proceedings of the IEEE International Conference on Computer Vision. pp. 6411–6420 (2019)
  * [25] Wang, L., Huang, Y., Hou, Y., Zhang, S., Shan, J.: Graph attention convolution for point cloud semantic segmentation. In: Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition. pp. 10296–10305 (2019)
  * [26] Wang, S., Suo, S., Ma, W.C., Pokrovsky, A., Urtasun, R.: Deep parametric continuous convolutional neural networks. In: Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition. pp. 2589–2597 (2018)
  * [27] Wang, Y., Sun, Y., Liu, Z., Sarma, S.E., Bronstein, M.M., Solomon, J.M.: Dynamic graph cnn for learning on point clouds. ACM Transactions on Graphics (TOG) 38(5), 1–12 (2019)
  * [28] Wu, Z., Song, S., Khosla, A., Yu, F., Zhang, L., Tang, X., Xiao, J.: 3d shapenets: A deep representation for volumetric shapes. In: Proceedings of the IEEE conference on computer vision and pattern recognition. pp. 1912–1920 (2015)
  * [29] Xu, Q.: Grid-gcn for fast and scalable point cloud learning. arXiv preprint arXiv:1912.02984 (2019)
  * [30] Xu, Y., Fan, T., Xu, M., Zeng, L., Qiao, Y.: Spidercnn: Deep learning on point sets with parameterized convolutional filters. In: Proceedings of the European Conference on Computer Vision (ECCV). pp. 87–102 (2018)
  * [31] Zhang, Z., Hua, B.S., Yeung, S.K.: Shellnet: Efficient point cloud convolutional neural networks using concentric shells statistics. In: Proceedings of the IEEE International Conference on Computer Vision. pp. 1607–1616 (2019)
  * [32] Zhao, H., Jiang, L., Fu, C.W., Jia, J.: Pointweb: Enhancing local neighborhood features for point cloud processing. In: Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition. pp. 5565–5573 (2019)