Typical hydrodynamic models for aquaculture nets_ A comparative study under pure current conditions

The net is regarded as the most critical component in marine aquaculture facilities as it is the only barrier which protects the environment from fish escapes. Accurate predictions of the net cage deformation and drag force on the nets are needed, both for ensuring fish welfare and for dimensioning of the mooring system. Thus, an appropriate hydrodynamic model is essential. In practice, two types of hydrodynamic force models, i.e., the Morison type and the Screen type, are commonly used to calculate the hydrodynamic forces on nets. Application of the models depends on the underlying structural model and the availability of data. A systematic review of hydrodynamic models is therefore undertaken to compare the models and various parameterisations, in aid of model selection during the design. In this study, eleven commonly used hydrodynamic models, i.e., five Morison models and six Screen models, are reviewed comprehensively, and implemented into a general finite element (FE) solver for dynamic simulations. Sensitivity studies on different current velocities, inflow angles and solidities of the nets are carried out. Moreover, different wake effects are also considered in numerical simulations. The numerical results from different models are compared against existing experimental data under pure current conditions. Suggestions for selection of suitable hydrodynamic models are provided, based on the model comparison.


Introduction
According to the Food and Agriculture Organization of the United Nations, aquaculture has been the world's fastest-growing food production method in the past 40 years (FAO, 2018). The development of high-value seafood such as Atlantic Salmon (Salmo salar) and Rainbow Trout (Oncorhynchus mykiss) has led to significant investments in the aquaculture industry. According to the Norwegian Seafood Research Fund, the seafood industry has invested more than NOK 115 billion in Norway since 2000 (Blomgre et al., 2019). The investment in the aquaculture industry upgraded conventional farming facilities and generated novel aquaculture structures such as Ocean Farm 1 and Havfarm. These innovative structures are designed to operate in the open sea and aiming to minimise the environmental impacts of aquaculture.
Moving aquaculture to the open sea can benefit the fish welfare and the ecosystem through better water exchange and dispersal of waste over a larger area (Cardia and Lovatelli, 2015). The more exposed setting implies larger waves and currents, which can increase the environmental loads on aquaculture structures. As the environmental loads on nets account for more than 85 % of total loads on a conventional fish cage (Cheng, 2017), accurate predictions on the hydrodynamic responses of nets are essential in the structural design.
The hydrodynamic forces on aquaculture nets depend on both the experienced current velocity and the hydrodynamic characteristics of the net. The hydrodynamic characteristics depend on the material, geometrical parameters (including mesh shape, mesh size, twine diameter), and production methods, i.e., twine weaving method (twisted or braided) and net weaving method (knotless or knotted). As shown in Fig. 1, the four nets which are commonly used in marine aquaculture have different hydrodynamic characteristics due to their different materials and production methods. Different twine materials and twine weaving methods make the surface roughness different. Higher surface roughness will generate larger turbulence regions; and thus, higher drag force (Balash et al., 2009). According to experimental data from Tsukrov et al. (2011), copper nets (smooth) exhibit significantly lower drag resistance in steady currents than nylon nets (rough) of the similar solidity. The experimental data from Lader et al. (2014) indicates that the drag force on the knotted net is up to 10 % higher than that of the knotless net, given the same environmental condition. Thus, in order to predict the hydrodynamic forces precisely, all the related parameters should be considered in calculations.
However, considering all the related parameters complicates the force prediction and make it impractical in numerical simulations. To make it feasible, one has to focus on key parameters while ignoring the secondary parameters. Through a large number of experiments (Tsukrov et al., 2011;Tang et al., 2018;Lader et al., 2014;Balash et al., 2009) researchers found that hydrodynamic characteristics are mainly dependent on two dimensionless variables, Reynolds number (Re) and solidity (Sn). The Reynolds number is defined as the following: where U is the undisturbed fluid velocity, is the kinematic viscosity of the fluid, d w is the twine diameter. For a typical aquaculture net, the Reynolds number is in the range of 100-10000. In some research (Kristiansen and Faltinsen, 2012;Balash et al., 2009), the Reynolds number is defined with local velocity. Models that define the Reynolds number locally lead to increased Re, although the modelled flow field does not change. In addition, the effect of the speed-up velocity on the hydrodynamic loads may need further investigations (Moe-Føre et al., 2015). Thus, the Reynolds number in the present study is calculated based on Eq. (1).
Solidity is the other key parameter for hydrodynamic characteristics. For a net panel, the drag force is mainly dependent on the value of solidity without obvious effects of twine diameter and mesh size, which both define solidity itself (Klebert et al., 2013). By definition, solidity is the ratio between the projected net area A p (area of the dark lines in Fig. 2(b)) and the total area of the net panel A t (area of the dashed box in Fig. 2(b)). For an ideal knotless square net, a mathematical expression for Sn can be formulated as: where d w is the twine diameter, L is the half mesh size. Since knotless nets are not always "mathematically perfect" and the solidity can change when the nets are submerged, some use a simplified expression for the solidity (Berstad et al., 2013): The deviation between Eqs. (2) and (3) is less than 10 % for a typical net in aquaculture where the solidity is in the range of 0.1-0.4. These two equations are obtained using the assumption of a perfect net geometry. Besides, digital image processing (DIP) techniques can also be used to estimate the solidity of an aquaculture net. The net solidity is  (Tang et al., 2018), (b) Knotted nylon net with rhombic mesh and Single English knot (Tang et al., 2018), (c) Welded silicon-bronze net (Tsukrov et al., 2011) and (d) Knotless nylon net with square mesh. H. Cheng, et al. Aquacultural Engineering 90 (2020) 102070 processed as the ratio between pixels in different colours in the DIP method (Yu, 2017). In general, the estimated solidities based on the DIP evaluations are less than 4% discrepancy compared with Eq. (2) for typical aquaculture nets (Tsukrov et al., 2011). Thus, the solidity of the net panel in the present study is calculated based on Eq. (2).
To acquire acceptable force predictions on aquaculture nets, researchers have conducted considerable amounts of both experiments and theoretical analyses. Based on the experimental results, hydrodynamic models are proposed to calculate the forces on net panels or net twines for numerical simulations. Some of them are implemented into in-house codes or commercial software (Table 1) for the fish cage simulation. In general, the hydrodynamic models can be classified as two types, the Morison model and the Screen model. Tsukrov et al. (2000) applied the Morison model extensively from fishing gear modelling to fish cage simulations and developed software, Aqua-FE. Lee et al. (2005) proposed a Screen model and implemented it into software, MPSL, for numerical simulations of both fishing gear and fish cage. In SIMA, the hydrodynamic forces on the net panels are calculated based on the Screen model proposed by Løland (1991). From the published articles (references in Table 1), the simulation results based on these software and codes agreed well with the experimental data when the velocities were smaller than 0.5 m/s. However, the agreement with the experimental data is weaker when the velocities were higher than 0.5 m/s and/or the solidities of the net panels were higher than 0.35.
Besides the hydrodynamic models for force prediction, the wake effect behind upstream net panels is also essential in the dynamic analysis of fish cages. According to the previous research (Faltinsen and Shen, 2018), the anchor force increases up to 22 % if the net-to-net wake effect was neglected in the numerical analysis.
The purpose of the present study is to summarise and compare the result of the commonly used hydrodynamic models under pure current conditions. Reviews of eleven hydrodynamic models are given in Section 2. Then, the present program for fish cage simulation, which includes the eleven hydrodynamic models, is introduced in Section 3. In Section 4, extensive validation cases are performed against experimental results. Series of experiments conducted by Lee et al. (2005); Tang et al. (2018);Tsukrov et al. (2011), and Moe-Føre et al. (2016) are reproduced by the program to validate the code and investigate the performances of different hydrodynamic models. Finally, the results of this study are discussed with concluding remarks.

Hydrodynamic models
According to the summary in Table 1, there are two types of hydrodynamic models for aquaculture nets: Morison model and Screen model. In Morison models, the forces on the net panel are treated as the sum of forces on individual twines. In Screen models, the forces are calculated based on considering the net as a panel. These models were initially proposed to calculate the hydrodynamic forces on fishing nets, especially trawl nets. With the slowdown of the fishing industry and booming of aquaculture, researchers have extended these models for aquaculture nets. In general, the hydrodynamic forces on net panels in an oscillatory flow can be written in Eq. (4).
where is the fluid density, ∇ is the submerged volume of the net structure, A is the reference area (the projected area of a net twine or a net panel), u is the fluid velocity vector, v is the structure velocity vector. The equation contains two empirical hydrodynamic coefficients-the added mass coefficient, C a and the drag coefficient, C d . The coefficients are determined from experimental data and dependent on the Keulegan-Carpenter number, Reynolds number, twine's surface roughness and the solidity of the net (Kristiansen and Faltinsen, 2012). The first three terms in the right side of Eq. (4) are Froud-Krylov force, diffraction force and radiation force, respectively. Combining these three terms will get the inertial force, F I in Eq. (5). The last term in Eq.
(4) is also called viscous force. Thus, the total hydrodynamic forces on the net are the sum of the inertial force and the viscous force. The inertial force can be rewritten to the following format: where C m is the inertia coefficient (C m = C a +1). For aquaculture nets, a typical value of C m = 2 is used for the inertia coefficient (Lader et al., 2007). The inertial force is experienced either during the transition to steady-state or when the net is experiencing wave loads. According to the experimental data and the dimensional analysis by Zhao et al. (2008), the inertial force on commonly used nets is in the order of O [d w /H] times the viscous force, where the twine diameter d w is far smaller than the wave height H. Since the twine diameter d w is far smaller than the wave height H, the inertial force can be negligible compared to the viscous force, e.g. the inertial force is approximately 0.1 % of the viscous force when d w = 1 mm and H = 1 m. Thus, it is reasonable to ignore the inertial force when one calculates the hydrodynamic forces on an aquaculture net. Due to this reason, the inertial force on aquaculture nets is not discussed in this study.

Morison model
In the Morison model, the forces on nets are calculated based on individual twines and knots. The twines are taken as cylindrical elements, and knots are taken as spheres. In practice, the viscous force (the last term in Eq. (4)) is usually decomposed into the two components: normal drag force (F n ) and tangential drag force (F t ): where L is the twine length, d w is the twine diameter, is the fluid  (Shainee et al., 2013;Decew, 2011;DeCew et al., 2010;Tsukrov et al., 2003) SIMA Screen type Truss (Li and Ong, 2017;Faltinsen and Shen, 2018;Li et al., 2018) FhSim Morison type / Screen type Triangles / Mass-spring (Reite et al., 2014;Endresen et al., 2013) AquaSim Morison type Truss (Berstad and Aarsnes, 2018;Berstad and Heimstad, 2017;Reichert, 1994 H. Cheng, et al. Aquacultural Engineering 90 (2020) 102070 density. u r n and u r t are the normal and tangential velocity of fluid relative to the twine. C n and C t are the normal and tangential drag coefficients. A 2D illustration of the force directions is given in Fig. 3.
C n and C t are crucial for force predictions as they determine how much force is generated in numerical simulations. In most of the cases, the normal and tangential force coefficients are functions of Reynolds number. Table 2 summaries the two coefficients for the twines when 100 < Re < 10 000 based on literature. From M1 to M5, the expressions of normal drag coefficient increase in complexity. M1 and M2 treat the normal drag coefficient as a constant value independent of Reynolds numbers. M3 and M4 include the variable Re to reflect the different fluid flow regimes. M5 adds another variable, Sn, to include the effect of net solidity. The details of these models and their applicable regions are given in Appendix C.
According to Table 2, C t is much smaller as compared to C n . Therefore, the simulations results are still acceptable, even C t is ignored in M2 and M5 (Cifuentes and Kim, 2017a;Wan et al., 2002). The Morison model was originally used in the fisheries research to calculate the forces and deformations of the fishing gears, especially the trawl net. The solidity has a small effect on the drag coefficients due to the large ratio of mesh size to twine diameter for typical trawl nets. Thus, the effect of solidity on the drag coefficients has not been included in Morrison type models. M5 (Cifuentes and Kim, 2017a) is the first Morison model to considers the solidity. However, one should note the strict application area of M5, since the negative quadratic term in C n can result in unrealistic values for large values of the solidity or/and the Reynolds number. Fig. 4 shows the normal drag coefficients of twines in the five models together with the normal drag coefficient of a smooth cylinder. It indicates that when 100 < Re < 10 000, all the C n of twines in the different hydrodynamic models are similar with values between 1.1 and 1.3. These values are close to C n of a smooth cylinder in the subcritical Reynolds number region.
The advantage of Morison models lies in their format. Since the formulation of Morison models is coincident to the line-type elements in structural models, application of a Morison model is directly compatible with the structural model. Thus, it is easy to implement Morison models into FE solvers to calculate the hydrodynamic forces. As shown in Table 1, Morison model is dominated among the software and codes. However, there are some drawbacks in the Morison models, i.e., (1) The velocity decomposing follows the independence principle, while this principle is only partially successful in correlating measured force data (Zdravkovich, 2003); (2) The Morison models can overestimate the drag force when the angle of attack (α) is small as it is not able to capture the interaction between the twines (Kristiansen and Faltinsen, 2012); (3) It is impractical to build a numerical mesh for a fish cage twine by twine because it requires millions of truss or spring elements to represent the net. A large number of elements will slow down the simulation significantly and make the problem numerically stiff and impossible to solve in a reasonable amount of time; (4) Although one can reduce the number of elements with assumptions about the deformation of the geometry, the assumptions can be incorrect under some circumstances and challenging to be implemented in existing structural analysis software. In order to mitigate the defects of Morrison models, Screen models are formulated. A detailed comparison between Morison and Screen models is shown in Section 4.2.

Screen model
In Screen models, the hydrodynamic forces are calculated based on a panel section of the net. The twines and knots in the net panel are considered as an integrated structure. In practice, the viscous force on the net panel is decomposed into components, either relative to the panel or relative to the flow. In some screen models (Fridman, 1973), the viscous force is decomposed into normal drag force (F N ) and tangential drag force (F T ), which are related to the orientation of the net panel. The expressions of these two components (Eqs. (8) and (9)) have a similar form with Morison models (Eqs. (6) and (7)), except that the reference area is changed from the projected area of a net twine d L w to the total area of a net panel A t .
where u r n and u r t are the normal and tangential components of the fluid velocity relative to the net panel. C N and C T are the normal and tangential drag coefficients of the net panel, which are dependent on the Reynolds number and the solidity. This Screen model formulation decomposes the ambient velocity into tangential and normal velocities similar to the Morison model.
Other Screen models decompose the viscous force on the net plane into drag and lift forces (F D and F L in Eqs. (10) and (11)) relative to the direction of the ambient stream velocity. Fig. 3. A 2D illustration of the hydrodynamic forces on twines. F n and F t are the normal and tangential drag forces on twines, respectively. The angle of attack α is the angle between the current direction and the axis of net twine. H. Cheng, et al. Aquacultural Engineering 90 (2020) 102070 where A t is the area of the net panel element, u r is the fluid velocity relative to the net panel, i D and i L are unit force vectors to indicate the directions of drag and lift forces. Procedures to calculate A t and unit force vectors (i D and i L ) are presented in Appendix A. C D and C L are the drag and lift force coefficients, respectively. The coefficients are determined from experimental data and depend on the Reynolds number, solidity and inflow angle (Fig. 5). The relationships of F N , F T , F D , F L and the inflow angle (θ) are shown in Fig. 5. The relationships of C D , C L , C N , C T are given in Eqs. (12) and (13).
Six Screen models which are the most commonly used in fish cage simulations are compared in this section and C N and C T are converted to C D and C L using Eqs. (12) and (13) when required. Figs. 6 and 7 present the drag and lift force coefficients of the six Screen models within the applicable solidity range of corresponding models. S6 (Balash et al., 2009) is not shown in Fig. 7 since no expression for the lift coefficient is given. The Reynolds number is assumed as 1000 for S3, S4 and S6 models in Figs. 6 and 7. For S3, the harmonic terms (a 3 and b 4 ) should increase with the increasing solidity, but no quantitative relationship is given in Kristiansen and Faltinsen (2012). Thus, the harmonic terms are set according to the experimental data reproduced by Shimizu et al. (2018), whereas Sn = 0.29 and a 3 = 0.15. Since S5 applied in MPSL (Table 1) did not disclose its formulation in the published article (Lee et al., 2005), the coefficients are assumed linearly proportional to Sn and independent of Re.
As shown in Fig. 6, the values of C D decay as expected with increasing inflow angle. The first two models (S1 and S2) have a similar shape and follow the dashed line, which indicates their C D decay following the cosine function. The values of C D in S3, S4 and S6 decay faster than the cosine function with the increasing inflow angle. According to the expressions of S3 and S6, the drag coefficient is zero when θ = 90˚. It means that the drag force is zero when the flow is parallel to a net panel, which is incorrect as there must be a drag force, although very small. If one used these two models to design a squared fish cage, the drag force can be underestimated when half of the nets are parallel to the ambient flow. Among the six Screen model, only the value of C D in S5 decays slower than the cosine function with the increasing inflow angle. A detailed discussion on the drag force on a net panel with experimental data is given in Section 4.2.2.
According to Fig. 7, the values of C L first increase, and then decrease with the increasing inflow angle, similar to an airfoil. The value of C L is zero when the flow is parallel or perpendicular to the net panel. The curves of S1 and S2 are similar to the sine function of 2 . While for S3, S4 and S5, the crests of curves are located between 30°and 45°. Compared to the values of C D , the values of C L is relatively small when θ < 30˚. It means the drag force is the dominated force when the inflow angle is small. However, the ratio between lift and drag forces is changing with different inflow angles. The appropriate hydrodynamic model should represent the observed ratio of lift/drag in experiments.
Solidity is an essential parameter in Screen models. Through a large number of experiments (Klebert et al., 2013;Zhou et al., 2015;Tang et al., 2018), researchers found that the hydrodynamic coefficients are highly dependent on the solidity of the net panel. Thus, all the Screen models in this paper take the solidity, Sn, into their expressions (the detailed expressions and its applicable region are given in Appendix C). In general, the values of C D and C L increase with increasing solidity, which indicates that nets with larger solidity have higher hydrodynamic forces when the other conditions are the same. Fig. 8 shows the drag coefficients of knotless nylon net panels when θ = 0°with different solidities from the available experimental data (Zhou et al., 2015;Tsukrov et al., 2011;Gansel et al., 2015). The regression curves in the figure are fitted using the ordinary least squares methods. The coefficients of determination (R 2 ) show that the cubic regression fits the data better than the simple linear regression. This observation complies with the expressions in S1 (Aarsnes et al., 1990) and S2 (Løland, 1991).
The flow patterns around nets should change with the Reynolds numbers, and thus influence the hydrodynamic coefficients. The   Fig. 4. Normal drag coefficients (C n ) versus Reynolds number (Re) according to different hydrodynamic models. Because C n in M5 model changes with different solidities (Sn), the filled area represents its variation range for its applicable solidities (0.172 < Sn < 0.208). The normal drag coefficient of the twine in S6 is a polynomial function according to its formulae (Kristiansen et al. 2012). The normal drag coefficients for a smooth cylinder and the typical Re region (100-10 000) for twines are also shown in the figure. H. Cheng, et al. Aquacultural Engineering 90 (2020) 102070 hydrodynamic coefficients in S1, S2 and S5 are constant with changing Re, as they do not include the Reynolds number in the expressions. While in S3, S4 and S6, the hydrodynamic coefficients are changing with different Reynolds numbers. According to Fig. 1, the effect of Re might be negligible since the drag coefficient of a net twine is almost unchanged when 100 < Re < 10 000. A detailed comparison of the drag force on net panels with experimental data under different Reynolds numbers is presented in Section 4.2.1.
Screen models are seldom used in the commercial FE solvers (see Table 1) for fish cage simulations, due to the complexity of implementation. The structural solver usually calculates the motion and deformation of aquaculture nets based on the line-type elements (bar, pipe or beam). In order to implement Screen models into the existing FE solver, other types of element (shell or plane) must be introduced to calculate the hydrodynamic forces and extra steps are required to map the hydrodynamic forces to the line-type structural elements. From a . S6 is omited due to the lack of formulas for the lift coefficient.
H. Cheng, et al. Aquacultural Engineering 90 (2020) 102070 programmer's point of view, Screen models require more algorithms than Morison models to fulfil dynamic analyses of fish cages. Thus, Screen models are not commonly used in the software and codes, referring to Table 1. Although more algorithms are required when integrating Screen models with a general FE solver, some codes provide Screen models to mitigate the defects of Morison models. In this study, the authors build a new module in a general FE solver to calculate the hydrodynamic forces on aquaculture nets using both Screen models and Morison models. The aforementioned models (five Morison models and six Screen models) are integrated into the new module to compare their performance with experimental data. The detail description of the new module is given in Section 3.

Wake effect
Wake effect is an essential and complex mechanism in analyses of permeable structures, such as the nets in fish cage and fishing gear. The wake is the region of disturbed flow (often turbulent) downstream a structure, caused by the viscosity of the fluid (Zhao et al., 2013a,b). In a fish cage, the wake effect means that the presence of upstream nets modifies the incoming flow velocity for downstream nets. In structural analyses of fish cages, the solver calculates the equilibrium between the  (1) Twine-to-twine wake effect, where a grid of i+1 cylinders (cross-section of a net panel) are exposed to an incident current velocity U. The U i (i = 0, 1 …) denotes the velocity experienced by cylinder i, which is modified due to the presence of upstream cylinders. (2) Net-to-net wake effect, where the upstream (left) net panel is exposed to an incoming current velocity U. The net-to-net wake effects from the upstream net panel result in a reduced flow (rU) at the downstream net. (3) Cage-to-cage wake effect, where the incoming flow for the downstream (right) fish cage is anisotropic and might be smaller than the incoming flow for the upstream (left) fish cage.
H. Cheng, et al. Aquacultural Engineering 90 (2020) 102070 external and internal forces of the structure and neglects the fluid mechanics. Therefore, the hydrodynamic forces on the downstream nets can be overestimated if no special precautions are taken to include the wake effect. In general, there are two ways to include the wake effect in the fish cage simulation. One can couple a structural solver with a fluid solver to include the fluid-structure interaction (Chen and Christensen, 2017). Alternatively, one can use the quasi-static assumption to "register" a wake region in the structural solver (Endresen et al., 2013). For the latter method, the wake effect for a conventional fish farm can be subdivided into three scales (Fig. 9) in the FE solver for implementation: (1) twine-to-twine wake effect; (2) net-to-net wake effect; (3) cage-to-cage wake effect.

Twine-to-twine wake effect
The twine-to-twine wake effect represents the interactions between net twines (the influence region is in the order of centimetres). In a net panel, the velocity of the downstream twine is smaller than that of the upstream twine when the inflow angle of the net plane is larger than 70°. According to (Endresen et al., 2013), when the inflow angle is 90°, the drag force on a net panel without twine-to-twine wake effect can be maximumly eight times larger than that with twine-to-twine wake effect.
The Morison model has a natural drawback on the implementation of twine-to-twine wake effect. To include this effect in the Morison model, one needs to make a function to describe the flow pattern behind a cylinder. For example, the wake shape around a 2D circular cylinder in an infinite fluid can be calculated based on Blevins formula (Eq. (14)): where U downstream is the velocity for the downstream cylinder at coordinate (x, y). In Morison models, the wake effect experienced by a single twine depends on contributions from all other twines in the model. Excessive numerical calculations and sophisticated algorithms are required to determine the spatial relationships among a large number of twines in the numerical model. On the other hand, the twineto-twine effect is naturally included in Screen models, since the hydrodynamic coefficients of net panels consider the interactions between twines implicitly. A detail demonstration on the twine-to-twine wake effect is shown in Section 4.2.2.

Net-to-net wake effect
The net-to-net wake effect is used to represent the interaction between nets inside a single fish cage (the influence region is in the order of tens of metres). Approximately half of the nets in a cylindrical fish cage will experience the net-to-net wake effect. The mooring force can be overestimated up to 22 % if the net-to-net wake effect is neglected in fish cage dynamic analyses (Faltinsen and Shen, 2018). In practice, a flow reduction factor (r) is adopted in software and codes to represent the net-to-net wake effect. Eq. (15) is a typical expression for the net-tonet wake effect, where r is the flow reduction factor (0 < r < 1), U is the ambient current velocity. According to this equation, the downstream nets experience smaller current velocity compared to the upstream nets.
In numerical simulations, whether a net is located in the wake can be determined based on its position, the centre of the fish cage and the incoming current direction (see Fig. 10). The flow reduction factor can be set as an attribute of the downstream nets to reduce the ambient velocity numerically.
An accurate flow reduction factor is critical for predicting the hydrodynamic forces on and the deformation of the nets in the wake (the blue part in Fig. 10). Table 3 shows flow reduction factors from theoretical analyses and experimental results. According to Table 3, the flow reduction factor should be a function of Reynold number, solidity and inflow angle. The most commonly used flow reduction factor, r = 1-0.46 = C D ( 0 ) , is consistent with different inflow angles. That means all the downstream nets in the rear half of a fish cage experience the same reduced current velocity, which is unphysical and contrary to the experimental results reported by Bi et al. (2013). Fig. 11 shows the flow reduction factor with respect to the inflow angle for downstream net panels in a cylindrical fish cage. In this figure, the data from experiment by Bi et al. (2013)  can be used to represent the equivalent drag coefficient of a net panel in the wake. In this figure, the drag coefficients C ( D ) are calculated based on the S1 model for both curves. For the dashed line, r is calculated according to Endresen et al. (2013). For the solid line, r is calculated by 1-0.46 = C D ( 0 ) according to Aarsnes et al. (1990). As shown in this figure, the equivalent drag coefficients of downstream net panels based on the two flow reduction factors are similar when θ < 30°. However, with increasing inflow angle, the equivalent drag coefficient based on the constant flow reduction factor is larger than the one based on the variable flow reduction factor. It means around 2/3 of the downstream net panels in a cylindrical fish cage contribute larger hydrodynamic forces if the constant flow reduction factor is applied. A detailed study of the two forms of flow reduction factors is given in Section 4.3.

Cage-to-cage wake effect
The cage-to-cage wake effect is used to represent the interaction between cages (the influence region is in the order of a few hundred metres). In the marine aquaculture industry, fish cages are usually grouped in arrays as a fish farm. Due to the block effect of the upstream cage, the flow for the downstream cage is affected by the existence of the upstream cage. In the previous study, a theoretical expression for the velocity reduction behind a net panel (r = 1-0.46 C D , r is the Fig. 10. Illustration of the method to identify the nets which experience the netto-net wake effect caused by upstream nets in a cylindrical fish cage. The fish cage is shown from the top, and the blue part is the rear half of a cage where the nets will experience a reduced flow. The inflow angle θ of a net panel is the angle between the normal vector of the net panel and the incoming flow. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article).
H. Cheng, et al. Aquacultural Engineering 90 (2020) 102070 velocity reduction factor, C D is the drag coefficient of a net panel) was used to represent the wake effect between cages. The theoretical expression gives a uniform reduced flow throughout the entire wake. However, according to the numerical simulations reported by Bi and Xu (2018), the current velocity around a fish cage is reduced by 38.3 % at the back and increased by 14.4 % at the two sides. According to the experiment by Turner et al. (2016), the wake after a fish cage is nonuniform, and the current velocity was reduced up to 62 % behind the fish cage and increased 19 % underneath the fish cage. Compared to the theoretical expression, the nonuniform wake flow is more physical and realistic. In addition, high levels of large scale turbulence were also observed behind a fish cage (Turner et al., 2016). Therefore, the theoretical expression cannot sufficiently describe the wake behind a fish cage.
The cage-to-cage wake has not been fully implemented into any FE solver or codes now due to its complexity. The wake topology is dependent on the environmental conditions, the status of the upstream fish cage and the spatial relationships among the fish cages. Although the wake shape and velocity profiles can be pre-predicted through accurate computational fluid dynamics (CFD) simulations, complex and verified algorithms are still needed to implement such pre-predictions into a FE solver for fish farm analyses.

Structural solver
A general FE solver, Code_Aster, is selected as the structural solver in this study. Code_Aster is EDF R&D's open-source FE solver for the thermo-mechanical study of structures (Electricité de France (EDF), 1989-2017). Since the code is open-source, it can be extended with additional functionality. With over 20 years' development, this software offers 400 finite element typologies for the discretisation of solids and a broad range of solvers. It enables the static, dynamic, vibrational analysis as well as modal analysis. Thus, it satisfies the requirements for structural analysis of fish cages.
The finite element discretisation of nets in a moving fluid environment results in the same differential equations as Tsukrov et al. (2003), shown as below: where x is the time-dependent vector of nodal displacements, M is the mass matrix, K is the stiffness matrix, F s is the nodal force vector due to gravity and buoyancy forces and F h is the nodal force vector for the hydrodynamic forces. The equation is highly nonlinear because the terms in the right-hand side (F s and F h ) depend on time, motion and deformation. The nonlinearity makes the classic linearized methods impossible to accurately describe the displaced configuration of the structure (Aubry, 2019). In the present FE solver, the solution technique for Eq. (16) is based on the unconditionally stable HHT-α method to integrate the equations in time. The discretised form of Eq. (16) in time is: Together with the displacements and velocities discretions (Eqs. (18) and (19)), the recurrence relation for the HHT-α method is obtained.

Structural element
The structural element used in the present study is the homogeneous, one-dimensional finite element denoted 'CABLE' in Code_Aster, which was initially developed to calculate the mechanical behaviour of overhead electrical lines. This two-node element is a version of the classic 'bar' element, adapted to the large displacement problem. It is suitable for representing highly flexible line-like structures (Antonutti et al., 2018). The element section is constant and maintains continuous orientation in the local frame. In the 3D space, the 'CABLE' element has six nodal degrees of freedom in the global coordinate system, which correspond to its nodal translations. Flexible nets can be achieved with less 'CABLE' elements than with bar or truss elements since the catenary shape is already taken into consideration by the element's shape function (Stengel and Mehdianpour, 2014).

Hydrodynamic elements
Two types of hydrodynamic elements were developed in the present study to represent the two types of hydrodynamic models (i.e., Morison models and Screen models). All the eleven aforementioned hydrodynamic models in Section 2 are compiled together with the hydrodynamic elements as an external module to Code_Aster. Fig. 13 shows the numerical simulation procedure with the external module which is highlighted by the red dashed box. The external module is invoked at each time step to calculate the hydrodynamic forces on the nets and maps the forces onto corresponding nodes in the structural elements. H. Cheng, et al. Aquacultural Engineering 90 (2020)

The hydrodynamic element for Morison models
The hydrodynamic element for the Morison models is a one-dimensional line. Since they are compatible with the existing finite element codes, the mesh for the hydrodynamic model is inherited directly from the structural model. In the external module, the hydrodynamic forces on net structure are first calculated by Eqs. (6) and (7), and then mapped to their corresponding nodes in the structural elements.

The hydrodynamic element for Screen models
It is difficult to integrate Screen models to a general FE solver (Bore et al., 2017), as there is no direct way to build the "hydrodynamic panel element" from the nodes and line-elements in the structural model. In this study, an automatic meshing function was developed to generate the two-dimensional virtual panel elements for Screen models. The functionalities in the automatic meshing module included: (1) create the "hydrodynamic panel element" automatically based on the position of nodes; (2) calculate the normal vector and the area of the "hydrodynamic panel element"; (3) calculate the hydrodynamic forces on the net panels based on the formulas of the different Screen models and map the forces to the nodes in corresponding structural elements.

Mesh grouping method
In a full-scale fish cage, the cage net is composed of thousands of small twines. It is impractical to build a numerical model twine by twine. A mesh grouping method is used in Section 4.3 to reduce the computational effort. In the present method, the material properties of the numerical model are assumed consistent with that of the prototype net. In order to acquire the right solutions, the M, K, F s and F h in Eq. (16) should be consistent between the physical and numerical nets. To satisfy the consistency of the aforementioned variables, three derived diameters, i.e., structural diameter (d ws ), elastic diameter (d we ) and hydrodynamic diameter (d wh ), are applied during the model building through the external module (the red dashed box in Fig. 13). The detailed derivation is illustrated in Appendix B, and only the final relationships between the three numerical diameters and the physical twine diameter (d w0 ) are presented here: ; where λ is the ratio between the half mesh size of the numerical net (Ln) and the half mesh size of the physical net (Lp).

Results and discussion
In this section, four experiments (Lee et al., 2005;Tang et al., 2018;Tsukrov et al., 2011;Moe-Føre et al., 2016) are selected to validate the newly developed program and discuss the aforementioned hydrodynamic models. The first case demonstrates the feasibility of the FE solver through a comparison with an experiment in Section 4.1. The second case compares the discrepancies among the different hydrodynamic models for net panels with respect to different current velocities, inflow angles and solidities. Based on the comparison in Section 4.2, the appropriate hydrodynamic model is selected for the fish cage simulations. In Section 4.3, a comparative study on wake effects is carried out based on fish cages under pure current conditions by using the verified solver.

Feasibility of the structural solver
As it is the first time that the open-source FE solver, Code_Aster, is used to simulate aquaculture nets, the feasibility of the solver should be assessed at the very beginning. In order to assess the FE solver itself, the external module (highlighted by the red dashed box in Fig. 13) and the mesh grouping method are not applied to the numerical model for avoiding interferences.
The feasibility is assessed through (1) examining iterative convergence and (2) examining consistency in the solution. The numerical model is set up based on the experiment by Lee et al. (2005). In the experiment, the net is 12 × 18 meshes in squared shape with twine diameter d w = 0.4 mm and half mesh size L = 100 mm. The Young's modulus of the twine is 119.37 MPa. The four corners are fixed, and three sinkers, whose masses are 0.5 kg, 1.5 kg and 0.7 kg from left to right in Fig. 14(a), are hung in the middle of the net. In the numerical simulation, the characteristics and configuration of the net are the same as those in the experiment. The net is modelled by 424 elements and 213 nodes. The three hung sinkers are represented by three vertical concentrated forces, which are 5 N, 15 N and 7 N from left to right. The density of the twine is assigned 1125 kg/m 3 by assuming the material is Nylon (Moe et al., 2010). The Newton-Raphson method is adopted to solve the FE equations iteratively.
The final shape of the net from the numerical simulation is shown in Fig. 14 and compared with the experimental results. Regarding the iterative convergence, the criterion is that the maximum force residue is less than 2e-5. The present simulation converges after 200 iterations by using 25.7 s. Regarding the consistency in the solution, the balance of forces in the vertical direction is checked. The total reaction force on the four fixed nodes in the vertical direction is 27.06 N which is equal and opposite with the to the sum of the wight of the net 0.06 N and the three concentrated forces 27 N. Through the two examinations, the present FE solver is proved feasible to simulate the flexible net.

Comparative study on the hydrodynamic models
The second study is based on net panels. Four net panels with different parameters which are wildly used in the aquaculture industry are selected from limited available experimental data (Tang et al., 2018;Tsukrov et al., 2011) to study the applicability and accuracy of the aforementioned hydrodynamic models with respect to net structures, ambient velocities and inflow angles. The mesh grouping method is not applied to the four net panels to ensure that the differences of the predicted forces only come from the hydrodynamic models themselves. The parameters of the four net panels are given in Table 4, and the structure of the four nets are shown in Fig. 1. The size of the net planes in the numerical simulation is 1 m × 1 m, and the four edges of the net planes are fixed in the simulation. Fig. 15 shows the drag forces on net planes for different current velocities with θ = 0˚using the eleven hydrodynamic models. In general, the drag forces increase with the increasing current velocity but with different increasing rates which depend on the expressions in hydrodynamic models.

Drag forces under different current velocities
Solidity is an important factor for force predictions. In general, the higher solidity can induce larger drag force. The predicted drag forces for N4 (higher solidity net) using Morison models can fit the experimental data well, except for M5 when the velocity is 1 m/s. According to formulas in M5, the drag coefficient could be negative when Re Sn * 2 > 218. That means when the solidity is 0.3 and the Reynold number is higher than 2 400, the drag coefficient can be negative. Thus, one should notice this strict applicable region when using this model.
Knots can increase the hydrodynamic forces on nets. Compared to N4, N2 has smaller solidity, which means the drag force on N2 should be smaller than N4 when both net panels were experiencing the same current velocity. However, due to the additional drag force from the knots on N2, the two net panels have similar total drag force under the same environmental conditions. For the knotted net (N2), Morison models underestimate the drag forces as they ignore the effect from H. Cheng, et al. Aquacultural Engineering 90 (2020) 102070 knots. These phenomena are consistent with the findings from Lader et al. (2014), in which the drag force on the knotless net is up to 10 % less than that of the knotted net. For the knotted net (N2), the predicted forces based on Screen models fit the experimental data better than that based on Morison models. In particular, the predicted forces based on S4 and S6 are very close to the experimental data, because these two models have included the effect from knots. However, it is also observed that S4 model always over predicts the drag force on knotless net panels (N1, N3 and N4). One should notice that this model was proposed more than 40 years ago. At that time, the aquaculture was only a small industry compared to the fishing industry. The researchers used fishing nets, most likely knotted nets, to generate this hydrodynamic model. It can also be proved by the fact that the predicted drag forces based on S4 model fit the N2 best compared to the other models. Thus, S4 is suitable to predict the hydrodynamic force for the knotted nets. Different materials can make the twine surface roughness different, and the smooth surface can reduce the drag force. For the silicon-bronze net (N3), all the hydrodynamic models overestimate the drag force compared to the experimental data, especially when the velocity is higher than 0.5 m/s. When the velocity is 1 m/s, the discrepancies between the experimental data and the predicted forces are varied from 43 % to 113 %. However, for the Nylon nets (N1 and N4), the discrepancies between experimental data and the predicted forces can be as low as 0.4 %. Since all the eleven hydrodynamic models in the present study were developed based on fibred nets whose surface is rougher than that of metal nets, they are unsuitable for the smooth metal nets. Moreover, the experimental results reported by Cha et al. (2013) revealed that the coefficients of chain-link woven copper alloy nets are smaller than that of the fabric nets with similar solidity only for inflow angles under 60°. Additional research work is necessary to have a better understanding of the hydrodynamic differences between fabric nets and copper alloy nets.
Although the solidity has the same physical meaning, the formulae to estimate solidity are different in the six Screen models. Table 5 compares the solidities of N1 and N4 based on the analytical formulae in S1-S4 Screen models with that from experimental data. The analytical estimations of N1′s solidity are within 5% of the experimental data. However, for the solidity of N4, the relative difference between the estimated and experimental value can be as large as 10.1 %. That large difference in the solidity can affect the accuracy of the predicted force. Thus, the predicted drag forces on N4 (high solidity net) have large deviations than that for N1 (low solidity net) when using the Screen models. In addition, when attaching the net panels to the frame, a pre-tension is usually needed to keep the net in the desired shape in the experiments. The different pre-tensions can make the solidity of the net panel various. Thus, it would be better to use the digital image processing (DIP) technique to estimate the solidity when testing the hydrodynamic properties of nets.
The mesh orientations have negligible effects on drag force in numerical simulations when inflow angles θ = 0°. For Morison models, the drag forces on a net panel are the sum of the force on each twine. The sum of the projected area of the twines does not change with the changing orientation. For Screen models, the drag forces are calculated based on the total area of a net panel whose area is also unchangeable with the changing orientation. Thus, the predicted forces based on both types of hydrodynamic forces models are independent of orientations when the flow is perpendicular to the net panel. For instance, the two net panels in Fig. 16 should have the same drag force when experiencing the same perpendicular flow. However, the drag forces on the two net panels can be different when θ ≠ 0° (Balash et al., 2015). A numerical study indicates that when θ > 45°, the drag force on net (b) is larger than that on net (a) given the same flow velocity and direction (Bi et al., 2017).

Drag and lift forces under different inflow angles
In practice, most of the nets in a cylindrical fish cage are not perpendicular to the incoming current. Thus, it is essential to compare the hydrodynamic forces under different inflow angles. Due to limited experimental data, only N1 and N2 have the experimental results under different inflow angles. Since most of the hydrodynamic models are not applicable to the knotted net (N2), we only discuss N1 in this section. Fig. 17 shows the drag forces, lift forces, and lift to drag ratios of the net plane for different inflow angles θ. The current velocity in this figure is fixed to 0.6 m/s to be consistent with the experiments reported by Tang et al. (2018).
In general, the drag force decreases with the increasing inflow angle, but the decreasing ratios are various due to the different formulas in hydrodynamic models. As shown in Fig. 17, the calculated drag forces based on the five Morison models have similar values, and all fit the experimental data well when the inflow angle is smaller than 70˚. While θ > 70˚, all the Morison models overestimate the drag force due to the absence of the twine-to-twine wake effect. That means the drag force on 22 % nets of a cylindrical fish cage can be overestimated. The overestimated drag forces could lead to inaccuracy in the prediction of displacements and cultivation volumes. Thus, the accuracy of the simulations based on Morison models is low when the fish cage has large deformation (Moe-Føre et al., 2016). For Screen models, not all the models agree with the experimental results well. S3, S5 and S6 underestimate the drag forces when θ > 30˚; S4 overestimates the drag forces when θ < 30˚and underestimates the drag forces when θ > 30˚; S1 and S2 agree with the experimental data quite well for all inflow angles. According to the drag coefficients in Fig. 6, only S1 and S2 decrease along with the cosine function with the increasing inflow angle. The drag coefficients of the rest Screen models decrease much faster than the cosine function. That is the reason why the drag forces are underestimated by S3-S6 when θ > 30˚. It is observed that the drag forces based on S3 and S6 models are zero when θ = 90˚. That unphysical value contradicts with the experiment data by Zhou et al. (2015). Therefore, S1 and S2 are more accurate than the other Screen models in the drag force prediction for N1 net panel with different inflow angles. The relative differences between the two models and experimental results are shown in Fig. 18 and discussed later.
Regarding the lift force, it first increases and then decreases with the increasing inflow angle according to the experimental data in Fig. 17. In H. Cheng, et al. Aquacultural Engineering 90 (2020) 102070 general, curves of the lift forces are similar to the shape of the sine function. For Morison models, the predicted lift forces are similar to each other, and all of them are smaller than the experimental results when θ > 30˚. The underestimations of lift force might lead to underestimations on fish cage deformations. For the Screen models, S3, S4 and S5 overestimate the lift force when 15˚< θ < 45˚due to their large lift force coefficients; S1 and S2 slightly underestimate the lift force when 30˚< θ < 60˚; S6 has zero lift force based on its formulas. The lift-to-drag ratio is a dimensionless parameter to express the relationship between lift and drag forces. The experimental data indicate that the maximum lift-to-drag ratio is 0.5. For the Morison models, the curves are close to each other, and all of them can fit the experimental results well when θ < 30˚. All the Morison models underestimate the lift-to-drag ratio when θ > 30˚, because of the overestimated drag forces and underestimated lift forces. For Screen models, the curves are distinct among the different models. Only S1 and S2 fit the experimental results. It can be observed that when θ > 45˚, the lift-to-drag ratios based on S3 and S5 are larger than one, and are two times higher than the experimental results. That irrational relationship could lead to incorrect simulations where N1 is used as the cage net.
Based on the aforementioned discussion on Fig. 17, four models, i.e. M4, M5, S1 and S2, are chosen to calculate the relative difference between their predicted results to the experimental results. As shown in Fig. 18(a), the drag forces predicted by the four models are within 5% of the experimental results when θ < 70˚. However, the drag forces  H. Cheng, et al. Aquacultural Engineering 90 (2020) 102070 predicted by Morison models are more than twice of the experimental results when θ > 70˚, due to the lack of the twine-to-twine wake effect. These overestimations in large inflow angle were also observed by Endresen et al. (2013). According to Fig. 18(b), the lift forces predicted by Morison models are less than half of the experimental results when θ > 45˚. The underestimated lift forces together with the overestimated drag forces might lead to incorrect results when θ > 45˚. The lift forces and lift-todrag ratios predicted by Screen models are better than those predicted by Morison models, especially when θ > 45˚. In particular, the relative difference of S1 is as high as 10 % compared to the experimental results in Fig. 18(c).
According to the comparisons and discussions on the four net panels, S1 is more suitable than the others to predict the hydrodynamic forces on net panels for different velocities and inflow angles. Thus, S1 is chosen to apply to a cylindrical fish cage in the nest section to study the wake effect.

Comparative study on the wake effect
In this section, an appropriate hydrodynamic model, S1, is used to study the wake effect with the corresponding experiment by Moe-Føre et al. (2016). Since the twine-to-twine wake effect is already included in S1 implicitly, and its effect has been discussed in Section 4.2.2, this section is focused on the net-to-net wake effect. The two forms of flow reduction factor which are discussed in Section 2.3.2, are applied to numerical models to study the wake effect. The expressions of the flow reduction factor are shown as below: Eq. (21) is the most commonly used formula in the dynamic analysis of fish cages (Løland, 1991;Aarsnes et al., 1990;Kristiansen and Faltinsen, 2012 , the nets with larger inflow angle ( ) induce smaller flow velocity in its wake region. A comparison between the two formulas with experimental results by Bi et al. (2013) and Patursson (2008) are shown in Fig. 19. Based on the experimental results, the flow reduction factor decreases with the increasing inflow angle. Compared to the commonly used formula, , shows better agreement with the experimental results.

Numerical models and loads
Two fish cages are used in this section. The two fish cages have the same diameter and height, but different solidities for the nets. The parameters for the numerical models and corresponding experimental  H. Cheng, et al. Aquacultural Engineering 90 (2020) 102070 models are listed in Table 6. The nodes in the upper circumference of the numerical model are restricted from translational motion, representing the rigid and fixed steel ring in the physical model. In the experiments by Moe-Føre et al. (2016), each sinker is a circular steel cylinder with a diameter of 4 cm, a length of 6 cm, and a submerged weight of 4.48 N, as given in Table 6. In the numerical model, the 16 sinkers are represented by 16 vertical concentrated forces corresponding to the submerged weight. Fig. 20 shows the physical and numerical fish cage models in still water. It can be observed that both the physical and numerical fish cages are stretched in Z-direction due to the weights.

Convergence studies
In order to demonstrate the reliability of the present numerical results, convergence studies on both computational mesh and time step are performed at first.
In the convergence study of computational meshes, five different sets of computational meshes shown in Table 7, are created for the fish cage with high solidity net in Table 6. Drag forces on the fish cage are estimated by using the five sets of computational meshes under flow velocity of 1 m/s. As shown in Fig. 21(a), the relative differences of drag forces are less than 3%, which demonstrates that the present mesh grouping method, as discussed in Section 3.2 and Appendix B, is workable for aquaculture nets. As shown in Table 7, with the increasing number of nodes (elements), the computer memory and computational time are increased, and the difference of the drag force compared to the finest computational mesh (Mesh 5) is reduced. In order to achieve the results within 1% difference compared to the finest computational mesh (Mesh 5) and keep the computational costs low, Mesh 3 is chosen for the subsequent simulations. By using Mesh 3, the numerical model consists of 64 elements (64 nodes) with 85.90 mm length around the circumference and 16 elements (17 nodes) with 93.75 mm length over the depth. The total numbers of elements and nodes are 2112 and 1088 in the numerical model, respectively. In addition, by assigning different diameters in the present solver according to the relationships in Eq. (21), the two physical fish cage models with different solidities in Table 6 are modelled by using the same computational mesh resolution, i.e., Mesh 3.
In the convergence study of time steps, four different time steps listed in Table 8 are applied in the simulations by using Mesh 3. Drag forces on the fish cage under different time steps are calculated under flow velocity of 1 m/s. As shown in Fig. 21(b), the drag forces first increase then decay fast with oscillations as the time increases; After 6 s, all the simulations reach equilibrium. As shown in Table 8 and Fig. 21(b), the drag forces on the fish cage calculated with the four different time steps reach the same value at the end of simulations. Increasing the time step can reduce the computational time, significantly. Since the simulations are calculated under pure current conditions without any oscillating load, the studied time steps have neglectable influences on the final results as long as the simulation is converged. Therefore, the subsequent simulations are calculated by using Mesh 3 with a time step of 0.2 s and a duration of 10 s. for the net-to-net wake. From the side view, the two models with different net-to-net wake effects have significant distinctions in the deformation, especially at the rear part. The model using f Sn ( ) 1 has larger deformation at the rear half of the fish cage. According to Figs. 12,19 and discussions in Section 2.3, the equivalent drag coefficients of the downstream nets with a constant flow reduction factor is much larger than the one with the variable flow reduction factor, especially when θ > 30˚. Therefore, the rear half of the cage experiences smaller drag forces and has less deformation when f Sn ( , ) 2 is applied. In addition, the deformations at the frontal half of the fish cage are similar in the two numerical models, because the frontal nets experience the same current velocity in both models. Fig. 24 shows the normalised height of fish cage in numerical simulations with the two net-to-net wake effects. The normalised height is calculated as the height of fish cages at a given current velocity divided by the initial height of the fish cage (1.53 m). Since the bottom nodes of the fish cage are not in a horizontal plane, the height of the fish cage is calculated as the vertical distance between the lowest node and highest node. It can be observed that the height decreases with increasing current velocity. The height of the model using = r f Sn ( , ) 2 is clearly larger than that using = r f Sn ( ) 1 , and the distinctions become significant with increasing current velocity. In particular, when the Sn = 0.347, the normalised height of the fish cage is 0.26 for the model using . The distinction in the height of fish cage can influence the design of feeding system and on-site operations related to nets, as the height of fish cage should be provided to make a precise decision. Fig. 25 shows the comparison of the drag forces using the two forms of flow reduction factor with experimental data by Moe-Føre et al. (2016). According to the experimental results: (1) the drag force on the low solidity (Sn = 0.194) fish cage is nearly proportional to the current velocity; (2) the drag force on the high solidity (Sn = 0.347) fish cage increases slower with increasing current velocity when the velocities are above 0.5 m/s than that at lower velocities.

Comparison of the drag force
The calculated drag forces using both net-to-net wake models increase with increasing current velocity and are close to the experimental results when the current velocity is less than 0.5 m/s. Compared to the experiments conducted by Moe-Føre et al. (2016), the model using = r f Sn ( ) 1 overestimates the drag force, especially when the current velocity is high, and the overestimations are more evident for the higher solidity fish cage. For the model using = r f Sn ( , ) 2 , the slope of the drag force curve decreases when the current velocity exceeds 0.5 m/s, and the predicted drag forces agree with the experimental results quite well. In particular, the maximum difference between the numerical and experimental results is only 5% when using = r f Sn ( , ) 2 . And the drag force on the fish cage can be as large as 30 % higher than the experimental results when applying = r f Sn ( ) 1 . According to the experimental photos in Fig. 23, the fish cage has large deformation, i.e., the nets have large inflow angles, when the current velocity is high. Together with the comparison in Fig. 19 which indicates that f Sn ( ) 1 highly overestimates the flow reduction factor when > 70°, the drag force on the downstream nets can be overestimated when applying = r f Sn ( ) 1 . Therefore, the total drag force on the fish cage is overestimated when using f Sn ( ) 1 . The comparison of the two net-to-net wake effects indicates that the H. Cheng, et al. Aquacultural Engineering 90 (2020) 102070 commonly used expression, = r f Sn ( ) 1 , is not sufficient to model the interaction between the fluid flow and nets (hydro-elasticity). The variable flow reduction factor, = r f Sn ( , ) 2 , is recommended for numerical simulations of the fish cage with the high solidity nets and subjected to high current velocities.

Conclusion
In the present study, five Morison models and six Screen models for calculating hydrodynamic forces on aquaculture nets are reviewed and implemented to a general FE solver, Code_Aster. The main contributions of the present study are listed as follows: 1 It is the first time that Code_Aster, the open-source FE solver developed by EDF R&D, has been used to simulate nets in marine aquaculture facilities. The successful application is fulfilled through the external module in the present work. 2 Verifications based on computational mesh and time step convergences and validations with experimental results are achieved. It is shown that by employing the newly developed external module in Code_Aster, the present numerical model can predict the response of a flexible fish cage under pure current conditions with satisfactory accuracy. In particular, the maximum difference between the numerical and experimental results is only 5% in drag force prediction. 3 The new formula proposed in the present work can fix the evident defect in the previous formula for the net-to-net wake effect. With the help of the new formula, the discrepancy between the predicted and experimental drag force for a fish cage can be reduced from 30 % to 5%.
Besides, comparative studies on the hydrodynamic models and wake effects are conducted using the developed model and workflow. The following conclusions are drawn from this study: 1 When the incoming flow is perpendicular to the net panel (θ = 0°), the drag force on knotless nylon nets can be well predicted by all the hydrodynamic models except for S4 which is originally for knotted nets. The discrepancies between experimental data and the predicted forces can be as low as 0.4 % when the current velocity is 1 m/s. 2 For the metal net with smooth surface, all the eleven hydrodynamic models overestimate the drag force. That is because all these models are initially developed for twisted or braided nets with rough surfaces. Further studies are needed to develop a new hydrodynamic model for metal nets.   H. Cheng, et al. Aquacultural Engineering 90 (2020) 102070 3 Knots can bring additional drag force on nets. Morison models underestimate the drag force on knotted nets if the effects from knots are not considered. As for Screen models, the drag forces on the knotted net can only be well predicted by S4 and S6. 4 When the incoming flow is not perpendicular to the net panel (0°< θ ≤ 90˚), drag forces predicted by the five Morison models are within 5% of the experimental results if θ < 70˚. However, the predicted drag forces are two times higher than the experimental results when θ > 70˚, due to the lack of the twine-to-twine wake effect. As for Screen models which include the twine-to-twine wake effect implicitly, the predicted forces are within 10 % of the experimental results for all inflow angles. 5 The drag force on a single fish cage is overestimated by the existing Screen models, especially when high-solidity nets experience large   H. Cheng, et al. Aquacultural Engineering 90 (2020) 102070 deformation. This is due to the inappropriate net-to-net wake effect. The consistent flow reduction factor, which is commonly used in the fish cage simulation, can overestimate the current velocities on downstream nets. Thus, the total drag force on a fish cage can be as large as 30 % higher than the experimental results.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.
In the present study, the "hydrodynamic panel element" is a triangle (e.g. P1-P2-P4 in Fig. A1). The triangular net panel elements are generated automatically by a mesh-generating function. The initial area of the triangular net panel is approximately equal to L 2 /2, where L is the half mesh size. The area of the "hydrodynamic panel element" can change during the simulations and is dependent on the position of its vertexes. In the simulation, the coordinates of every node are known. The area of the "hydrodynamic panel element" is calculated according to the following equation: The unit normal vector of the net panel is calculated by the following equation: The unit force vectors (i D and i L ) which are used to indicate the directions of drag and lift forces are calculated as follows:

Appendix B. Derivation of the mesh grouping method
This section gives an introduction on how to reduce the number of elements using the present mesh grouping method. As mentioned in Section 3.4, the principle is to keep M, K, F s and F h in Eq. (16) consistent between the physical net and numerical net. In the numerical model, the fluid density ( fluid ), the fluid velocity (u), the density of twine ( twine ) and Young's modulus of twine (E) are consistent with the physical value, and λ is the ratio between the half mesh size of the numerical net (Ln) and the half mesh size of the physical net (Lp). As shown in Fig. B1, the nets in the two blue dashed boxes should have the same mass (M), stiffness (K) and environmental loads (F s +F h ). To satisfy the consistency, three derived diameters, i.e., structural diameter (d ws ), elastic diameter (d we ) and hydrodynamic diameter (d wh ), are used in the numerical model building. Below, the relationships between the three diameters and the physical twine diameter (d w0 ) are derived.
Mass conservation (M) As shown in Fig. B1, the mass of the physical net in the blue dashed box is: Fig. 25. Measured and calculated drag force in different current velocities using the two net-to-net wake effects.
H. Cheng, et al. Aquacultural Engineering 90 (2020)   H. Cheng, et al. Aquacultural Engineering 90 (2020)  Similar to the derivation for the mass conservation, the diameter for the gravity and buoyancy forces is the same as the one for the mass. Therefore, it can use the same parameter, d ws , to calculate the gravity and buoyancy forces.
Hydrodynamic forces (F h ) For both Morison and Screen models, the hydrodynamic forces are calculated based on the following equation: In Morison models, the hydrodynamic coefficients depend on the physical twine diameters (d w0 ) and the reference area A is the projected area of twines. The projected area of twines in the physical net is: The projected area of the twines in the numerical net is: (B.14) In Screen models, the hydrodynamic coefficients depend on the solidity or the twine diameters (d w0 ) of the physical net. The reference area A in the numerical model is the net panel area, which is the same as the physical net. The solidity of the physical net is: Base on = Sn Sn p n , the derived hydrodynamic diameter satisfies the same relationship in Eq. (B.14). In summary, Eqs. (B.5),(B.8) and (B.14) have been used in the present mesh grouping method to reduce the computational effort. ; D = knot diameter 0.051<Sn<0.235 -