3D–3D Computations on Submerged-Driftwood Motions in Water Flows with Large Wood Density around Driftwood Capture Facility

: The accumulation of driftwood during heavy rainfall may block river channels and damage structures. It is necessary to mitigate such effects by periodically capturing and removing driftwood from rivers. In this study, the behavior of driftwood in open-channel ﬂows with a relatively large wood density was modeled numerically. The water ﬂow and driftwood motion were solved three-dimensionally, with an Euler-type ﬂow model coupled with a Lagrange-type driftwood motion model. A piece of driftwood was modeled as a set of connected spherical elements in a straight line for easy analysis using a discrete element method. Wood with speciﬁc gravity exceeding 1 will travel along a position near the riverbed and will be affected by bed friction. In addition, friction forces for sliding and rolling motions are considerably different. Therefore, in the numerical model, a bed friction term was introduced between the bed and driftwood considering the anisotropy of the friction force. The variation in the drag force of water ﬂow on driftwood was also considered depending on the angle between the driftwood trunkwise direction and ﬂow direction. The model was applied under the same conditions as those used in a laboratory experiment on driftwood behavior around an inlet-type driftwood capture facility. The computational results showed that the proposed model could qualitatively reproduce the driftwood behavior around the capture facility. The secondary ﬂow patterns at the approaching reach and the capture ratio were found to be strongly affected by the turbulence model and the Manning roughness coefﬁcient.


Introduction
A common characteristic of recent disasters caused by heavy rainfall is the outflow of a large amount of driftwood and the consequent destruction of river structures such as bridges and houses or the blockage of river channels due to the accumulation of driftwood. To mitigate the effects of these disasters, measures to control the generation of driftwood and to capture and remove driftwood that has flowed out into the river channel can be considered. However, the former methods are not realistic because driftwood can be generated from a wide variety of locations. The latter methods are more realistic and economic if driftwood capture facilities are designed appropriately and installed at adequate locations in the river channels. In addition, the captured driftwood can be collected easily and used as a resource.
Several types of facilities can be employed to capture driftwood in river channels. One of the typically used facilities involves the creation of a bypass channel or cove on the side bank of the river to guide the driftwood and capture it. This type of capture facility has the advantage that it hardly obstructs the main stream of the river channel [1]. However, because the driftwood storage area is inevitably set up at a location that is divergent from the main stream area, it is difficult to efficiently guide the driftwood to the storage area. Okamoto et al. [2] performed a laboratory experiment on this type of facility and found that it is effective to set up a dead water area on the side bank of a straight channel and install a groyne on the opposite bank. Consequently, the flow diverted by the groyne can carry driftwood into the dead water area. Okamoto et al. also stated that, even if the lower part of the groyne is made water permeable, the effect of guiding the driftwood is retained. Kato et al. [3] examined the driftwood capture project planned in the Omoto River, Iwate Prefecture, Japan, by performing detailed laboratory experiments. This capture project aims to guide driftwood to a side capture area using inertial force, which can be realized by installing the facility on the outer bank side of the mountainous curved part of the river. Through a relatively large-scale hydraulic experiment (1/50 scale), Kato et al. investigated the shape that would be the most effective for the facility and found that the capture ratio could reach 65% if the facility were to be designed appropriately.
As a method for examining driftwood motion, numerical analysis is an efficient approach along with laboratory experiments. Various numerical analysis models have been developed to reproduce scenarios such as the driftwood behavior in river channels and driftwood deposition in capture facilities. In particular, two-dimensional (2D) hydrodynamic flow models based on Eulerian and Lagrangian modeling were proposed to model 2D driftwood transport. Ruiz-Villanueva et al. [4] developed the Iber-Wood model, which considers the balance of forces between the water flow and large wood. This model was verified through laboratory experiments on driftwood by considering steady flow [4] and unsteady flow (Ruiz-Villanueva et al. [5]). Kang et al. [6] performed numerical experiments to evaluate the model reproducibility of the driftwood motions using the two-dimensionally-connected spheres model, which simulates the bed deformation by driftwood. Moreover, Kang et al. [7] conducted a flume experiment and numerical simulation to study driftwood collision on the moveable bed. Persi et al. [8] also performed a flume experiment and calibrated ORSA2D_WT, which calculates the multiple local forces from the four segments of a cylindrical wooden body to simulate the wood motion. Such studies have mainly analyzed the driftwood based on a straight flume channel and two-dimensional wood motion, and have actively advanced the numerical modules for driftwood interactions, such as the movement of the riverbed with driftwood deposition, collisions between driftwoods, and water flow patterns created by the driftwood.
Shimizu and Osada [9] modeled the river flow using 2D (two-dimensional) depthintegrated equations and the driftwood using connected spherical elements. They assumed that driftwood moves two-dimensionally near the water surface and used a discrete element method (DEM) to consider the collisions between the wood pieces. Their model could simulate the capture phenomenon of driftwood by bridge piers. Hatta et al. [10] employed a similar model and showed that driftwood behavior in meandering channels can be appropriately reproduced. Because this type of model expresses both the flow and driftwood two-dimensionally, we hereafter refer to them as 2D-2D models (2D river flow and 2D driftwood model). On the other hand, Osada and Shimizu [11] built a sophisticated model that considered the elastic deformation of driftwood while extending the driftwood behavior to three dimensions under a 2D flow field (i.e., a 2D-3D model: 2D river flow and 3D driftwood model). They also used the spherical element type model for driftwood. Because the three-dimensionality of the flow is predominant at river bends, Kimura and Kitazono [12] constructed a 3D-2D model (3D river flow and 2D driftwood (spherical element model)) to study driftwood behavior at a river bend. Kato et al. [13] developed a novel model in which the flow was modeled using an advanced depthintegrated 2D model that included the effect of the secondary current of the first kind by assuming the vertical velocity profile with the Engelund model [14] and the driftwood motion was modeled using a 2D plane model. We can call such a model a q3D-2D model (quasi-3D flow and 2D driftwood model). It is noteworthy that this model uses a columnar shape directly instead of connected spheres for modeling the driftwood; therefore, the model has a smaller computational load and can be easily applied to real-scale phenomena. This model has been applied in the hydraulic experiments pertaining to the planned capture facility in the Omoto River, Iwate Prefecture, Japan, and has generally provided reasonable results. However, this model assumes a well-developed secondary flow field that occurs in an infinitely long uniform bend. Therefore, the secondary flow could possibly be overestimated in a finite length of a curved channel, which is usually the case in actual rivers. In a continuous periodic meandering channel, a secondary flow with opposite circulation to that predicted by the Engelund model may occur because of the lag between the development of the secondary flow and the curvature of a streamline. To consider the lag, it is necessary to use more advanced 2D flow models (e.g., [15,16]) or a fully 3D flow model.
Based on the existing studies on driftwood mentioned above, the models of driftwood are classified into two types-modeling of the driftwood as a columnar shape and representing a piece of driftwood as connected spheres. The models of the former type have the advantage of computational efficiency. On the other hand, the models of the latter type have the advantage that an object of arbitrary shape can be flexibly represented by altering the arrangement of the spheres. Furthermore, the collision between driftwood pieces can be easily calculated by the well-known DEM for spheres.
The objective of the present study was to numerically reproduce the driftwood capture experiment performed by Kato et al. [3] for the planned capture facility in the Omoto River. The present study used a 3D-3D model that expressed both water flow and driftwood motion three-dimensionally. The computational model was based on the model proposed by Kimura [17], which was developed to calculate driftwood stacking around bridge piers. In the experiment by Kato et al., driftwood with a specific gravity larger than 1.0 was considered. Because such driftwood primarily moves along the riverbed, it is necessary to model the friction between the riverbed and driftwood. Kang and Kimura [18] and Kang et al. [19] proposed a detailed computational model for the motion of driftwood in shallow flows; in this model, the sliding, rolling, stopping, and relocation of driftwood on a riverbed were modeled adequately. In the present study, a similar model was incorporated into the 3D-3D model for expressing the contact of driftwood with the riverbed. It should be noted that the frictional force differs greatly between the sliding movement in the driftwood axis direction and the rolling movement in the transverse direction if the shape of the driftwood is cylindrical. Therefore, an anisotropic friction model is employed. In addition, Kang et al. [19] modeled the change in drag due to variations in the collision angle between the driftwood in the trunkwise direction and that in the flow direction, which was ignored in the previous study [18]. A similar model for the drag was also incorporated into the present 3D-3D model. The computational model was applied under the same conditions as those used in the experiment by Kato et al., and the reproducibility of the model was verified.

Flow Model
In this study, the river flow was calculated using a 3D flow model. The model is similar to the 3D solver NaysCUBE [20] available on the river analysis platform iRIC [21]. This model solves the Reynolds-averaged Navier-Stokes equations and continuity equation on a moving-boundary adaptive-mesh system. Hydraulic variables were defined on a fully staggered grid using contravariant components for vectors and tensors. The TVD MUSCL method with spatial third-order accuracy was used for the discretization of the advection terms; the quadratic Adams-Bashforth method was used for time integration; and the nonlinear k-ε model by Kimura et al. [22,23] was employed for turbulence closure.
To check the effect of the turbulence model, we also performed computations with the linear standard k-ε model. The model was used to simulate the secondary current of the first kind [24] and the secondary current of the second kind [25], and showed excellent reproducibility in each case. Please refer to the literature on the iRIC website for details on basic equations, discretization methods, boundary conditions, and so on [21].
For the interaction between driftwood and flow, Kimura [17] adopted a one-way model that ignored the effect of driftwood on the flow for the sake of simplicity. In the present study, we used a two-way model that considered the effect of driftwood on the flow as drag. The drag force on the flow in a computational cell was evaluated by the following formula: where F dr i is the drag force that the driftwood exerts on the flow in a certain calculation cell, ρ is the density of water, C D is the drag coefficient, √ g is the Jacobian of the cell, A k is the flow direction projection area existing in the cell, u is the flow velocity vector, u p is the driftwood motion velocity vector, U i is the flow velocity vector component in the i direction, U p i is the driftwood motion velocity component in the i direction, i is the generalized curvilinear coordinate system (ξ, η, ζ), and N cell is the number of spherical elements that constitute the driftwood existing in the calculation cell.
The drag coefficient C D was considered as a function of the particle Reynolds number Re d as follows [26]: In Equation (2), the particle Reynolds number was defined as where d is the diameter of a spherical element and ν is the kinematic viscosity coefficient of water.

Driftwood Motion Model
•

Outline of Driftwood Motion Model
The driftwood motion model used in this study expressed a single piece of driftwood using connecting spheres. Figure 1 shows the basic concept of the model. To calculate the motion of the model, first, each sphere was advected without considering the connections in one time step. Then, the spheres were rearranged again in a straight line. By repeating these processes, the motion of the driftwood was calculated continuously. This method is similar to the models proposed by Shimizu and Osada [9], Hatta et al. [10], and Kimura and Kitazono [12], which applied the calculation method originally proposed by Koshizuka [27] based on the moving particle semi-implicit method to 2D driftwood calculations. First, the movement of each spherical element in one time step (∆t) is computed independently (using Equation (4)), and the center of gravity is obtained by averaging the values of the coordinates at the sphere centers. Next, the average linear velocity and average angular velocity around the center of gravity of the driftwood piece are calculated. Then, the spheres constituting the piece of driftwood are rearranged again in a straight line. The driftwood behavior is tracked by repeating this process step by step. In this method, the motion of the driftwood can be approximated as a rigid body movement when the calculation time interval is set sufficiently small. We extended this method to a 3D space. The details of the calculation method can be found in the paper by Kimura [17].

• Bed Friction Terms
This study targeted driftwood with a specific gravity greater than 1. Therefore, driftwood movement occurred primarily near the bottom surface and was accompanied by friction with the riverbed. In the previous model [17], driftwood with a specific gravity less than 1 was assumed; therefore, the bottom friction was not taken into consideration. Bed friction terms were newly introduced in the present model as follows.
First, each sphere constituting a piece of driftwood was moved during a calculation step ∆t according to the Lagrange equation while ignoring the constraint. The Lagrange equation considering the bottom friction term is as follows: where σ is the sphere density, ρ is the water density, C D is the drag coefficient, u is the water flow velocity vector, u p is the moving velocity vector of a sphere, g is the gravity acceleration vector, F p is the interparticle collision force vector, t is the time, C M is the additional mass coefficient (C M = 0.5), A 2 and A 3 are the 2D and 3D shape coefficients (A 2 = π/4, A 3 = π/6), λ A-sub is the ratio of the submerged part of the projected area of the sphere in the flowing water direction, λ V-sub is the ratio of the submerged part of the sphere volume, and F bed is the bottom friction force vector.
The diameter of the sphere, d, was calculated using the following equation, assuming that the densities of the sphere and driftwood are equal, and that the total mass of the spheres constituting one piece of driftwood and the mass of the driftwood are equal.
where d is the diameter of the sphere, D is the diameter of the driftwood, L is the length of the driftwood, and N is the number of spherical elements constituting a piece of driftwood. The bottom friction vector F bed is expressed as where z p is the height of the sphere center from the riverbed, F b is the bottom friction vector, and F s is the effect of gravity associated with the riverbed gradient. F b is expressed by the following equation: where µ p is the friction coefficient between the sphere and riverbed and is described in detail later. In addition, N is the normal force at the contact point between the bottom and sphere. If the bottom surface gradient is small, it is approximately expressed by the following equation: where W and B are the gravity and buoyancy forces acting on the sphere, respectively, and are expressed as As mentioned above, λ v-sub in Equation (10) is a coefficient representing the volume ratio of the submerged part of the sphere, which is usually 1 for the driftwood near the bottom, but this value becomes less than 1 when the water depth is small. Regarding the friction coefficient µ p , when the driftwood has a cylindrical shape, the friction coefficient of sliding in the driftwood axis direction and that of rolling in the transverse direction are considerably different. Thus, the friction coefficient was determined while considering the anisotropy of the frictional force based on the results of Kang and Kimura [18] and Kang et al. [19]. The profile of the friction force in such models is illustrated in Figure 2, which shows a 2D view of the presence of driftwood near the bottom. The distribution of friction coefficients in the driftwood axis direction (t direction) and transverse direction (n direction) has the maximum and minimum values, and the values in the other directions are assumed to be elliptical, as shown in the figure. Assuming that the friction coefficient of sliding in the driftwood axis (trunkwise) direction is µ t and the friction coefficient of rolling in the lateral direction is µ n , the friction coefficient µ p in the moving direction of the driftwood is expressed by the following equation: where ψ t is the angle formed by the driftwood movement direction and driftwood axis direction. In Equation (11), cos ψ t is obtained by the following equation: where θ is the angle formed by the x-axis and driftwood axis. With regard to the sliding friction µ t , we considered the difference between the static friction coefficient µ ts and dynamic friction coefficient µ tk (µ ts ≥ µ tk ). Because it is impossible to theoretically determine the values of these coefficients, trial-and-error procedures were performed to more faithfully reproduce the driftwood behavior of the experiment. Table 1 lists the values of the friction coefficients used in the present study. Driftwood movement that involves sliding and rolling contact with the bottom surface is affected by the bottom surface gradient. To reproduce this effect in the model, we considered the local x -y plane along the bottom surface. The bottom slopes in the x and y directions were approximated as follows using the slopes in the ξ and η directions in a generalized curvilinear coordinate system under the assumption of a small bottom slope: where z b is the bed elevation. The components in the x and y directions of the gravity forces due to the slope effect F s in Equation (6) were modeled as follows: Profile of bed friction using anisotropic friction model. In this study, driftwood was modeled by connecting the spherical elements. The drag coefficient C D for a sphere can be evaluated using Equations (2) and (3). However, in the case of driftwood, the drag force changes depending on the directions of the flow and driftwood axis. The variation in the drag coefficient was recently studied by Persi et al. [28], who considered the relationship between the angle of attack and the drag in the case of a columnar object. On the other hand, Kang et al. [19] assumed that the drag coefficient C D was constant, but the projected area for each spherical element changes depending on the angle between the flow direction and driftwood axis direction. The difference between these two studies is whether the variable drag force was included in the drag coefficient C D or projected area A. The latter approach was adopted in this study. Based on the study by Kang et al. [19], the variable projected area was described as follows: where A is the projected area of an isolated sphere (=π/4 d 2 ), A is the modified projected area considering the angle between the driftwood and flow, φ is the angle between the driftwood axis direction and flow direction (Figure 3), and M is the number of constituent spheres of the driftwood. The direction of the driftwood axis is defined to be positive from the sphere with m = 1 to the sphere with m = M. Figure 4 shows the relationship between the projection area and flow direction. From this figure, cosφ and |sinφ| are calculated as follows ( Figure 4): where u x , u y , u z are the relative velocity vectors of the flow with respect to the sphere and t x , t y , t z are the unit direction vectors along the driftwood axis.  Kato et al. [3] performed a laboratory experiment on a driftwood capture facility planned on a mountainous curved part of the Omoto River in Iwate Prefecture. Figure 5 (left) shows the aerial photograph of the planed site, and Figure 5 (right) shows the schematic diagram of the planned capture facility. The structure of the driftwood capture facility is planned such that a new river channel is dug on the inner bank side of the originally curved river channel, and the original river channel on the outer bank side is used as a driftwood capture area. Therefore, it is considered that the proposed shape will make it easy for the driftwood to be guided to the capture area by natural inertial force. Kato et al. [3] performed experiments by changing the position and size of the entrance of the facility in several configurations and found the shape with the highest capture performance. Figure 6 shows the shape optimized in this manner. The scale of the experimental equipment was 1/50. The river kilometer (KP: approximate distance from the river mouth along the river, unit is "km") ranges from 33.7 km to 34.3 km. The capture pond exists approximately between 33.9 km and 34.05 km.  The discharges used in the experiment were for 10-year return period flood (Q 10 = 469.7 m 3 /s) and 30-year return period flood (Q 30 = 680 m 3 /s). Two types of driftwood models with different lengths (6 and 12 m) and a common diameter (30 cm) were used. All the values were real-scale values. In each case, the specific gravity of driftwood was set to 1.1, considering the dominant species of trees (hardwood) at the site. The number of driftwood pieces was 3600 for 6 m long pieces and 1800 for 12 m long pieces.
The Manning roughness coefficient in the experiment was calculated using the velocities and depths observed in the experiment. The estimated values of the Manning roughness along the streamwise direction are shown in Figure 7. The averaged value of the Manning roughness was approximately n = 0.03.

Computational Conditions
In the numerical analysis, calculations were performed for the case with a driftwood length of 6 m from among the cases considered by Kato et al. In the experiment, 3600 driftwood pieces were used, but if all of them were to be considered in the analysis, the computational load would become high. Therefore, in the present computations, the number of driftwood pieces was set to 200 or 800.
Regarding the number of spherical elements constituting a piece of driftwood, we considered two cases of computations with 5 and 10 spheres for each piece of driftwood. Because we confirmed that the results in the two cases were almost the same, we employed the model with five spheres considering computational efficiency.
The calculation was performed on an actual scale. The planar shape of the calculation grid is shown in Figure 8. The river kilometer (KP) used in the computation ranged from 33.9 km to 34.3 km. The number of grid cells was 112 × 40. The resolution of the planar calculation grid was set through trial calculations so as to sufficiently reflect the basic topographical shape of the river channel without excessive calculation. However, the effect of the planar grid resolution should be further investigated in the future. The slit of the capture facility was expressed by setting the mesh cells at the slit bars as impermeable obstacles. However, because the grid cell size was larger than the slit interval (2 m), it was unavoidably set as an obstacle with a coarser interval (6 m). Driftwood pieces were supplied by dropping them from the air at the position marked with a circle in Figure 8 by considering the flow distance at which the effects of turbulence caused by the drop almost disappeared at the entrance of the capture facility.
We performed eight cases of computations listed in Table 2. Regarding the vertical grid division, the cases with 10 and 20 layers were examined. We used not only the second-order nonlinear k-ε model, but also the linear standard k-ε model for comparison of turbulence closure. From Figure 7, the average Manning roughness coefficient in the test reach is deduced as n = 0.03. However, the approaching area at the upstream capture facility has a lower roughness. Therefore, we tested n = 0.025, 0.02, and 0.015 in addition to the average value of n = 0.03. The specific density of wood used in the experiment was 1.1. We also considered cases with specific gravity of 1.0 and 0.9 to determine the sensitivity of driftwood behavior to the specific gravity.

Cross-Sectional Flow Pattern
First, we examined the structure of the secondary current of the first kind induced by an imbalance in the centrifugal force, which is considered to have a major effect on the driftwood behavior in the curved part. Figure 9 shows the flow velocity vectors in Run 1 and 2 along the cross-sections A-A' and B-B' (Figure 8) located in the upstream side of the capture facility and inside the capture area. The distribution of the vorticity in the mainstream direction is also shown as a color contour. In Figure 8, the upper and lower panels show the results of the cases with vertical grid divisions of 20 (Run 2) and 10 layers (Run 1), respectively. To make the cross-sectional flow patterns easier to understand, the scale in the vertical direction is magnified thrice when compared with that in the horizontal direction.
Regardless of the grid resolution in the vertical direction, the secondary current of the first kind was reproduced in the section A-A', and the flow from the outer bank (right bank) to the left bank was confirmed near the channel bed. In addition, in the case where the vertical grid division was 20 layers, a smaller vortex cell in the counter-clockwise direction was observed on the right bank side near the water surface. Hereafter, we call this vortex cell near the surface the "surface vortex". In the case where the division was 10 layers, no clear surface vortex was observed. In the section B-B', the secondary current was still dominant in the capture area, but the surface vortex was not observed even in the case of 20 layers.  Figure 10 shows the comparison of the cross-sectional velocity vectors and the color contour of the streamwise vorticity along the A-A' section with different Manning roughness coefficients. We can see that the secondary flow patterns are strongly affected by the magnitude of the roughness. It is clear that the vortex cell is in the clockwise direction (blue color of the streamwise vorticity denotes the secondary current of the first kind at the bend). However, the surface vortex (the vortex cell in the counter-clockwise direction near the surface, shown as the red colored area of the streamwise vorticity contour) has the opposite rotating direction. The plan view of the computational reach in Figure 8 indicates that a bend in the opposite direction exists at the upstream reach of the A-A' section (at KP34.2-KP34.3), which induces the secondary current of the first kind in the counter-clockwise direction. Because the development of the secondary current lags behind the shape of the bend, the secondary current of the first kind generated at the KP34.2-KP34.3 sections persists at the A-A' section. This is the cause of the surface vortex observed at the A-A' section.    Figure 12 shows the profile of Un/Us along the A-A' section in the results obtained with different Manning roughness coefficients. We can see clearly from this figure that the ratio of the secondary flow to the main flow becomes larger as the roughness of the bed increases. Figure 13 shows the profile of Un/Us along the A-A' section in the results obtained with different turbulence models. We can see that the result with the nonlinear model yields stronger secondary flow when compared with the result obtained with the linear model.   Figure 14 shows the simulated movement of driftwood and the flow velocity vector near the bed in the downstream area from the driftwood input point in the case of a vertical grid division of 10 layers. The time shown in the figure indicates the elapsed time after the driftwood was thrown in. The figure indicates that the pieces of driftwood travel near the channel bed and are divided into two groups: one group flows down along the main channel, and the other enters the capture area. The movement of the latter group is presumed to be due to the inertia force of driftwood resulting from the flow in the curved part. A clockwise recirculating flow was generated inside the capture area. Therefore, the driftwood that entered the capture area traveled along the left bank side with the recirculation flow, and a part of it accumulated and stopped because of the small flow velocity in and around the circulation. This simulated behavior is in qualitative agreement with the results observed in the laboratory experiment.

Driftwood Behavior
However, after t = 120 s, most of the driftwood pieces that entered the capture area passed through the slit framework on the downstream side and returned to the main stream in the computation. On the other hand, in the laboratory experiment, the ratio of driftwood passing through the slit framework was less than 10%. That is, this behavior in the laboratory experiment was considerably different from that in the computation. The reason for this difference is that the slit interval in the computation is considerably coarser than that in the experiment because the slit is expressed as an obstacle cell in the numerical grid, as previously explained. Figure 15 shows the velocity vectors near the bottom and driftwood distribution at t = 180 s. The green dashed line indicates the area where the driftwood stopped and were deposited. On the other hand, Figure 16 shows the deposition area in the laboratory experiment. Because the number of inputs was different, the deposition density was also considerably different, but the deposition location and its shape were almost the same.  In the laboratory experiment performed by Kato et al. [3], the specific gravity of driftwood was set to 1.1 considering the dominant species of trees at the site, but the density is lower for coniferous trees and dry trees immediately after runoff. To check the effect of the density of trees, we compared the results of Run 8 (specific gravity = 0.9), Run 7 (specific gravity = 1.0), and Run 1 (specific gravity = 1.1). All the other conditions (except the specific gravity) were maintained to be the same as those in the previous computation. Figure 17 shows the results at t = 60 s for three different specific gravities (0.9, 1.0, and 1.1). The flow velocity vectors show the surface flow condition at a specific gravity of 0.9 and the bottom flow condition at specific gravities of 1.0 and 1.1. When the specific gravity was 0.9, the driftwood was strongly influenced by the outward flow resulting from the secondary current (Figure 9 (lower panel)) and the driftwood headed toward the outer bank. Thus, the pieces of driftwood adhered to and accumulated near the side wall. Even when the specific gravity was 1.0, the movement toward the outer bank was very large. This behavior may be caused by the inertial force of the driftwood motion in the curved part in addition to the direction of the secondary flow near the surface.

Driftwood Capture Ratio
The driftwood capture ratios (number of captured pieces of wood/total supplied number of pieces × 100%) in each run are listed in Table 3. In the experiment, most of the driftwood pieces entering the capture area were trapped by the slit framework at the downstream of the capture area. However, in the present computations, most of the driftwood pieces passed through the slit framework because of the coarse grid resolution as explained above. Therefore, we calculated the capture ratio, which included the number of wood pieces that passed through the slit framework from the capture area. Table 3. Driftwood capture ratios in each run and in the experiment.

Run
Capture Ratio (%) The results show that the capture ratio increases as the value of Manning roughness coefficient decreases (comparison of the results in Run 1 (n = 0.03), Run 3 (n = 0.025), Run 4 (n = 0.02), and Run 5 (n = 0.015)). It is also clear that the capture ratio increases as the specific gravity of wood decreases (comparison of the results in Run 1 (specific gravity = 1.1), Run 7 (specific gravity = 1.0), and Run 8 (specific gravity = 0.9)). From the comparison of the results of Run 3 (nonlinear k-ε model) and Run 6 (linear standard k-ε model), we can see the effect of the turbulence model on the capture ratio. The computational result with the standard linear k-ε model shows a higher capture ratio than that with the linear model. In the laboratory experiment, the capture ratio was 65-87.5 % [3].

Reproducibility of Basic Characteristics of Driftwood Behavior
In terms of the movement, deposition, and capture behavior of driftwood around the inlet-type driftwood capture facility, the computational results showed that the proposed 3D-3D numerical analysis model was able to reproduce the basic characteristics of driftwood dynamics to some extent. However, the numerical analysis results had some differences when compared with the experimental results in terms of quantitative aspects such as the capture ratio.
One of the major differences between the numerical analysis results and experimental results was that the phenomenon of driftwood capturing using slits could not be reproduced by the numerical analysis. The reason for this is that, as described above, the slit framework was considered in the computation by setting the calculation grid cells as impermeable obstacles. Because the calculation grid was coarse, the slit spacing in the computation became larger than that in the laboratory experiment. There are two possible approaches for solving this problem. One is to increase the resolution in the computational grid and make the grid cells sufficiently small so that they can represent the real slit spacing. However, this approach will increase the computational load because the flow computation will have to be performed on a large number of computational cells with a smaller time increment because of the Courant-Friedrichs-Lewy (CFL) condition. The other approach is to place a fixed group of spheres at the positions of the slit rods and calculate the collision force between the spheres constituting the driftwood and those constituting the slit by a DEM. In this case, it is assumed that the spheres that constitute the slit do not move even when a collision occurs. As for the calculation grid, it is possible to use one with a resolution larger than the slit spacing. In general, the latter approach appears to be more advantageous in terms of computational efficiency. Thus, we intend to modify the model by this method in the future.
A comparison of the computational results of the eight cases listed in Table 2 showed that the simulated results are very sensitive to the Manning roughness coefficient, turbulence model, and wood density. We discuss these sensitivities below.

Effect of Manning Roughness Coefficient
We performed the computations with four different values of the Manning roughness coefficient (Run 1, 3, 4, and 5 with n = 0.03, 0.025, 0.02, and 0.03, respectively). The capture ratio becomes higher as the roughness becomes larger. The reason is obvious, namely, the change in the secondary flow pattern. Figure 10 shows that the secondary current of the first kind (clockwise vortex cell near the bottom) becomes smaller and the surface vortex (counter-clockwise cell near the surface) becomes larger as the roughness becomes smaller. In the case of submerged driftwood with specific gravity greater than 1, the flow velocity near the riverbed has a considerable influence on the behavior of the driftwood. Figure 12 demonstrates that the ratio of lateral velocity due to the secondary current, which guides the driftwood toward the inner bank side and hinders its capture, becomes smaller if the Manning roughness becomes smaller. In other words, if the roughness becomes smaller, the inertia force that pushes the driftwood toward the outer bank side becomes dominant. Thus, the capture ratio increases if the value of the Manning roughness coefficient increases (see Table 3 (Run 1, 3, 4, 5)). Figure 7 shows that the average value of the Manning roughness in the experimental channel reach is approximately 0.03. However, the computation with n = 0.03 (Run 1) considerably underestimated the capture ratio (see Table 3 (Run 1)). Figure 7 also shows that the Manning roughness at the area approaching the capture facility is smaller than the average value (= 0.03). Considering this feature, if the Manning roughness is set to be smaller than the average value (Run 3, n = 0.025), the simulated capture ratio becomes closer to the experimental result. These results imply that it is important to set the value of the Manning roughness coefficient carefully by considering its value at the approaching area.

Effect of Turbulence Model
We applied two types of turbulence models-a second-order nonlinear k-ε model and the linear standard k-ε model. The comparison of the results of Run 3 and Run 6 reveals the effects of the turbulence model on the flow and driftwood behavior. Figure 11 shows that the turbulence model affects the secondary flow pattern strongly at a cross section. The simulation of the surface vortex, which is induced at the previous bend in the opposite direction and persists at the next downstream bend, is larger for the linear model than for the nonlinear model. Thus, the value of Un/Us in Figure 13 is smaller in the result obtained with the linear model. Consequently, the pushing of the driftwood to the capture area located at the outer-bank side becomes dominant and the capture ratio becomes larger than in the case with the nonlinear model.
Because the secondary flow pattern was not measured in the laboratory test, we could not assess directly the accuracy of the result obtained with each turbulence model. However, the present computational results imply that the choice of the turbulence model is very important for the quantitative assessment of the driftwood capture facilities.

Effect of Vertical Grid Division
Regarding the grid division in the vertical direction, two types of vertical grid resolutions (10 layers and 20 layers) were tested. Figure 9 shows that the grid division affects the secondary flow pattern at the cross sections. The reproduction of the surface vortex was smaller in the result of Run 1 (division into 10 layers) when compared with the result of Run 2 (division into 20 layers). This difference will affect the advection of driftwood, especially if the specific gravity of the wood is less than 1. However, in the present cases with specific gravity of 1.1, the influence on the capture ratio was not large (capture ratios in Run 1 and Run 3 were 32% and 39%, respectively).

Effect of Wood Density
In the laboratory experiment by Kato et al. [3], the specific gravity of the driftwood was set to 1.1, which is the average value for hardwood at the planned site. However, dry wood or softwood (conifer) usually has a smaller value of the specific gravity. If secondary current exists at a river bend, the lateral velocity components near the bottom and near the surface are in the opposite directions (see Figure 9). Figure 17 clearly demonstrates that the behavior of driftwood at a bend is strongly influenced by its specific gravity. The capture ratios in Run 1, 7, and 8 in Table 3 are very different, and the case with smaller wood density has a larger capture ratio because of the lateral velocity component toward the inner bank near the surface. It is interesting that the secondary current of the first kind has an opposite effect on driftwood behavior regardless of whether the specific gravity of driftwood is greater than 1. If the specific gravity of driftwood is greater than 1, as in the present case, the secondary current hinders the driftwood from entering the capture area. In other words, when the specific gravity of driftwood is greater than 1, the inertial force, which pushes the driftwood toward the outer-bank side, and secondary flow have opposite effects on the driftwood around the capture facility. Therefore, the design of a driftwood capture facility should be performed carefully if the specific density of driftwood is expected to be greater than 1.

Conclusions
To simulate the behavior of driftwood motion and deposition around an inlet-type driftwood capture facility, we constructed a 3D-3D (3D flow and 3D driftwood) computa-tional model. Considering the dominant spices of trees at the planned site, we assumed that the specific gravity of the driftwood was greater than 1. Because such driftwood pieces travel near the channel bed, they are constantly affected by bed friction. Therefore, the modeling of the bed friction is important. We proposed an anisotropic bed friction model considering the different friction forces between the sliding (driftwood axis direction) and rolling (lateral direction) motions. We also modeled the variation in the drag force depending on the angle between the flow direction and driftwood axis direction by considering the change in the projection area of the flow with respect to each spherical element.
Although the numerical analysis results could generally reproduce the experimental results, the driftwood capture by the slit could not be reproduced because of insufficient grid resolution for precisely expressing the slit intervals. This problem was also pointed out by Kato et al. [13] in their computation using a 2D-2D model.
The computational results under different hydraulic conditions showed that the behavior of the driftwood and the capture ratio at the capture facility along the outerbank side of a river bend are strongly affected by the Manning roughness coefficient and turbulence model. The mechanism was explained as the variation in the secondary current at a cross section according to the roughness and turbulence model.
It was also confirmed that the driftwood behavior around the capture facility in the curved part was strongly influenced by the driftwood density because the 3D flow pattern was governed by the secondary current of the first kind and the lateral flow directions at the bottom and surface of the water were opposite because of spiral flow. When the specific gravity of driftwood is greater than 1, the inertial force of the main flow, which pushes the driftwood toward the outer-bank side, and secondary flow have opposite effects on the driftwood motion. This implies that more careful attention is required while designing a driftwood capture facility if the specific density of wood is greater than 1.