• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    A study of hybrid prediction method for ship parametric rolling*

    2016-10-18 01:45:29YaohuaZHOU周耀華NingMA馬寧JiangLU魯江MinGU顧民
    關(guān)鍵詞:馬寧耀華

    Yao-hua ZHOU (周耀華), Ning MA (馬寧), Jiang LU (魯江), Min GU (顧民)

    1. School of Naval Architect, Ocean and Civil Engineering, Shanghai Jiao Tong University, Shanghai 200240,China

    2. Shanghai Rules and Research Institute, China Classification Society, Shanghai 200135, China,E-mail:yhzhou@ccs.org.cn

    3. State Key Laboratory of Ocean Engineering, Shanghai Jiao Tong University, Shanghai 200240, China

    4. China Ship Scientific Research Center, Wuxi 214082, China

    ?

    A study of hybrid prediction method for ship parametric rolling*

    Yao-hua ZHOU (周耀華)1,2, Ning MA (馬寧)1,3, Jiang LU (魯江)4, Min GU (顧民)4

    1. School of Naval Architect, Ocean and Civil Engineering, Shanghai Jiao Tong University, Shanghai 200240,China

    2. Shanghai Rules and Research Institute, China Classification Society, Shanghai 200135, China,E-mail:yhzhou@ccs.org.cn

    3. State Key Laboratory of Ocean Engineering, Shanghai Jiao Tong University, Shanghai 200240, China

    4. China Ship Scientific Research Center, Wuxi 214082, China

    The parametric rolling (PR) in the head or following waves has been considered as one of the main stability failure modes in the development of the 2nd generation Intact Stability criterion by the International Maritime Organization (IMO). According to previous studies, the estimation methods of the roll damping affect the prediction of the PR significantly, and most of them are based on experiment data or Ikeda's empirical formula. The accuracy of the estimation method for the roll damping could be a key aspect for the validity of its prediction for the full scale ship. In this research, a hybrid prediction method is developed for the numerical prediction of the parametric rolling when experiment data are not available for the roll damping. Comparison study is also carried out between the hybrid method and a nonlinear dynamics method, where the roll damping is estimated by the simplified Ikeda's method and the direct CFD prediction method in a direct non-linear simulation based on the 3-D CFD approach in the model scale. It is shown that the results of the hybrid method are in satisfactory agreements with the model experiment results, and the method can be used for analysis especially at the early design stage where experiment data are often not available.

    hybrid prediction, parametric rolling (PR), roll damping, CFD

    Introduction

    Since the parametric rolling is a remarkably nonlinear dangerous phenomenon, a new intact stability criterion on the parametric rolling is developed at the IMO and will be introduced to the engineering application in a near future. This phenomenon is not caused by the external excitation, but by the time-varying changes in the parameters[1]. The prediction of the parametric rolling (PR) was widely studied, especially with respect to the influence of the estimation method of the roll damping[2]or the improvement of the design by increasing the roll damping[1,3], but the roll damping is still usually obtained by a semi-empirical method or the conventional model test. Due to the lack of experiment data in the initial design stage, the prediction of the PR could be unreliable, because the roll damping is not known a priori. Therefore, the estimation method of the roll damping should be further studied.

    The estimation methods based on the 3-D CFD approach have made a great progress in recent years. Yang et al.[4]carried out a simulation of the forced roll of a ship moving with or without the forward speed. Wilson et al.[5]carried out a simulation of the roll decay for a model with/without bilge keels of a combatant. Furthermore, Sadat-Hosseini et al.[6]even successfully carried out a direct simulation of the PR for a surface combatant based on the 3-D CFD approach. However in view of the computational efficiency as well as the applicability for the ship design, it is preferable to directly obtain the roll damping by the CFDapproach for the numerical prediction model of the PR.

    For this purpose, in this study, a hybrid method is developed based on the 3-D CFD approach and the 3-D potential method. By estimating the roll damping with the CFD approach and solving for the 3 DOF motions, the parametric rolling is simulated and validated for a containership. Good agreement is achieved between the numerical simulation results and the experiment data. Comparison with other two methods,including the direct CFD prediction method and the nonlinear dynamics method is made to see the applicability of this method.

    1. Theoretical formulation and numerical method

    The hybrid method is developed based on a nonlinear 3-D CFD approach and the 3-D potential method. The 3-D CFD approach is used for the calculation of the roll damping moment and the 3-D potential method is used for the calculation of the radiation and diffraction forces. The method follows a process shown in Fig.1.

    Fig.1 Conceptual scheme of hybrid method for PR assessment

    1.1 Theoretical and numerical models

    The hybrid method adopts a 3 DOF weakly nonlinear model (roll, heave and pitch) for the simulation of the ship motion. Such kind of models with consideration of the time delay effect and the nonlinearity of Froude-Krylov forces were successfully applied for the simulation of the parametric rolling[7-9]. Motion equations[10]are as follows:

    where FISis the composition force of the Froude-Krylov force and the restore force, which are calculated based on the 3-D pressure integration method for the instantaneous wetted hull, the diffraction forces FDare predicted by the 3-D frequency domain potential method, the radiation forces are calculated based on the impulse response theory in the motion equation with consideration of the memory effect as

    The roll damping coefficients with/without the forward speed are calculated based on the motion or moment data of the numerical simulations for the free decay or the forced roll of the scaled model. The simulation is carried out by the 3-D CFD solver ISISCFD[11,12]. This flow solver uses the incompressible unsteady Reynolds-averaged Navier Stokes (RANS)equations, based on the finite volume method to build the spatial discretization of the transport equations. The face-based method is generalized to 3-D unstructured meshes, where the non-overlapping control volumes are bounded by an arbitrary number of constitutive faces. The flow solver deals with the multi-phase flows and the moving grids.

    The k-ω(SST-Menter) turbulence model is used by the RANS solver, and the adaptive grid refinement method is adopted for capturing the free surface. The 6 DOF motion equations are directly based on the Newton's second law. The forces due to the interaction of the ship and the fluids (water and air) are calculated by integrating the 3-D pressures of the fluids on the hull surface of the ship and updating every time step, then the motions are solved. For the Newton's law, the RANS equations solver is used and the flow/motion coupling is updated each non-linear iteration until a converged solution at each time step is obtained. Then the pressure of the water/air will be calculated for each cell and integrated on the hull of ship.

    1.2 Estimation method of roll damping

    The free decay and forced roll experiments are usually carried out in the calm water to estimate the roll damping of the ship. In this study, numerical simulations for these two experiments are carried out by the CFD solver for the 3 DOF motions of the 3-D geometry at the model scale. And the roll damping coefficients are obtained by the conventional experimental data processing method.

    For simulation of the free decay, the validation for a ship at a forward speed was made in the Gothenburg 2010 Workshop for a bare hull with bilge keels[13]. Zhou et al.[14]validated the simulations for a ship at zero speed with experiments of four test cases at the model scale, including mono-hull ships with/ without complex appendages and a trimaran without appendages. In this study, the simulations for the free decay and the calculations of the roll damping coefficients are performed based on the method developed by Zhou et al.[14].

    For the simulation of the forced roll, the forward speed is accelerated from zero to the target velocity. Then the ship is forced to roll following the sinus law,and the roll amplitude slowly reaches the target amplitude to enter the stable rolling stage. Finally, the roll damping coefficient could be calculated by analyzing the moment of roll. The simulation of the forced roll is based on the 3 DOF motions, including roll, pitch and heave. Pitch and heave follows the free trim principle. This means that trim and drafts are not restrained and the ship is to “balance itself” in the water.

    The forced roll motion in the stable stage can be written as

    According to the Newton second law, the motion equation and the moment of roll are as follows:

    where I44,?I44are the roll moment of inertia and the added roll moment of inertia,A,C are the roll damping coefficients,C44is the restoring coefficient,M( t)=M0sin(ωt+ε1)is the exciting moment.

    Then the moment of roll is calculated by the CFD solver and can be written as

    The roll damping coefficients could be calculated by fitting the following formula:

    Based on the principle of conservation of energy,the equivalent linear damping could be calculated as

    2. Validations and discussions

    2.1 Ship model

    The well-known container ship C11[15]is adopted for the numerical simulation to validate the hybrid method. Table 1 shows the principal dimensions of C11. Figure 2 shows the 3-D numerical scaled model of C11. The scale factor is 0.01.

    Table1 Principal dimensions of C11 (full scale)

    Fig.2 The ship model of C11 containership

    Fig.3 Meshes of typical sections of C11

    2.2 Parameter settings and mesh generation

    (1) Mesh generation

    The simulations of the free decay and the forced roll with/without a forward speed are based on the same mesh. Figure 3 shows the geometries and the meshes of typical sections.

    The CFD computation is performed with the samescale factor as in the model test. The computational domain is -4<x<3,-4<y<4,-1.5<z<0.5for the free decay and the forced roll cases, in dimensionless coordinates which are non-dimensionalized by usingLpp. With the scale of the domain, the reflection of the wave-making caused by the forward speed of the ship at domain boundaries can be prevented, and the computational resource can be saved for the simulation of air with little effect on the roll damping. The ship axis is aligned with the x-axis with the bow at x =0.52and the stern at x=-0.48. The free surface lies atz=0. Enough grid points are generated near the free surface to capture the radiation waves and the6extreme motion. All predictions are made on 1.44×10 grids. The size of the grid near the free-surface is taken asLpp/500in the vertical direction,Lpp/50alongx axis, and 4Lpp/10alongyaxis. If the initial roll amplitude or the stable roll amplitude is very large,the sliding grid is essential for simulating the large roll motion. In this simulation, the sliding grid is not used for the simulation of the roll motion.

    (2) Parameter settings

    The numerical simulation is based on solving the 3 DOF motions (roll, pitch and heave) for the ship in calm water. The forward speed accelerates from zero to the target value following the 1/2sin law, and then keeps as a constant value. The heading is maintained,so the influences of sway and yaw are small and could be neglected.

    Fig.4 Errors of natural roll period of C11

    Fig.5 Forced roll simulation of C11

    For a blind simulation of the parametric rolling,the initial heel angle or the forced roll amplitude is difficult to determine for the estimation of the roll damping. Therefore, these two values are taken as 20ofor both CFD simulations. This is consistent with the Level 2 criteria of the parametric rolling developed by the IMO[16].

    The time step is a key parameter for the accuracy of the numerical simulation of the free decay[14]. In order to achieve an acceptable accuracy and avoid the divergence of the motion solving, 350 time steps per roll period (natural roll period and forced roll period)is adopted. The maximum number of the non-linear iterations is taken as 8, and the convergence criteria of the solver are taken as 2 orders of magnitude. For about 98.37% of cells, the courant number is from 0 to 1.1132.

    (3) Estimation of roll moment of inertia

    The predictions of the parametric rolling are blind in the sense that the actual roll moment of inertia(Ixx)of the scaled experiment model is not known a priori. This value depends not only on the geometry, but also on the mass distribution inside the geometry. If the exact value of this parameter is unknown, the validation with the measurement data will be problematical, since the predictions of the roll damping and the parametric rolling depend strongly on the accuracy of Ixx[9,14]. Fortunately, the natural roll period is measured during the model test. Since the natural roll period depends on Ixx, it can be used to determine the exact value ofxxI iteratively using the procedure described below.

    First, the initial value of Ixxis calculated with Eq.(16)[17]. Then the value of Ixxcould be adjusted based on the error of the natural roll period by simulating the free decay in calm water at zero speed with the CFD approach. Once the iteration ends with a satisfactory accuracy, the appropriate value ofIxxcould be taken for the further prediction of the roll damping and the parametric rolling.

    whereT is the natural period of the roll,Ixx-iis the initial value of the roll moment of inertia,?Ixxis the added roll moment of inertia,M is the mass of the ship,GM is the metacentric height,gis the gravity acceleration. The value of?Ixxis estimated by the 3-D frequency domain potential method at zero speed.

    Figure 4 shows the errors(ε)of the natural roll periods of C11 in the iterative process. The accuracy of Ixxcould reach an acceptable level after 5 iterations.

    2.3 Estimation of roll damping

    The roll damping coefficients are estimated for Fr=0, 0.05 and 0.10 by simulations of the free decay and the forced roll. Fig.5 and Fig.6 show the simulations of the forced roll and the free decay for C11.φ,θis the amplitude of roll and pitch,Mxis the roll moment.

    Fig.6 Time history of free decay simulation of C11

    Unlike the traditional simulation in 2-D cases,the heave and the pitch are unrestrained in the forced roll motions, thus the simulation results could be in better agreement. The time-domain response of the forced roll (Fig.5(b)) shows that the period of the roll is approximately twice the period of the pitch. This is very similar to the response period relation of the parametric rolling, which is the so-called 2:1 response. This phenomenon indicates that the roll damping characteristics simulated by the forced roll may be closer to the actual parametric rolling.

    Fig.7 Roll amplitudes predicted for PR

    2.4 Validation and discussions of parametric rolling results

    Figure 7 shows the predictions of the parametric rolling for C11 by the hybrid method, and Fig.8 shows the time history. These cases include 3 Froude numbers and the range of wave-steepness is from 0.0025 to 0.0400. Both of two options of the hybrid method (free decay and forced roll) are compared with the experiment data.

    According to the prediction results of C11 (as shown in Fig.7), with the hybrid method (the forced roll), the parametric roll could be successfully predicted for different speeds and wave-steepness. For cases when those amplitudes are around and less than 27o,a satisfactory accuracy is achieved.

    Fig.8 Time history of PR (Fr=0.05, wave-steepness: 0.02)

    Small refinements are made for meshes of each simulation of the roll damping with different forward speeds. Thus, the difference of the mesh quality may slightly affect the accuracies in cases at Fr =0.05.

    Fig.9 Comparison of effects of initial heel angles on equivalent linear damping coefficients

    (1) Influence of damping coefficient with different initial roll angles

    The initial heel angle or the forced roll amplitude plays an important role in obtaining the roll damping characteristics such as the equivalent damping coefficients[1]. Figure 9 shows the comparison of the effect of the initial heel angle on the damping coefficients estimated by the simulation of the free decay. In this study, the initial heel angle and the forced roll amplitude are both taken as 20ofor the estimation of the roll damping coefficients. Thus, the agreement is not good for cases with large or small amplitudes as shown in Fig.9. Therefore, it still needs further study how to determine the initial heel angle or the forced roll amplitude for the blind simulation of the PR by the hybrid method. These values could be taken as20otemporarily as is consistent with the 2nd generation criteria of IMO[16].

    (2) Influence of roll period

    On the other hand, the roll periods of the free decay and the forced roll also affect the estimated roll damping coefficients and the results of the hybrid method. As shown in Fig.7, if the roll damping coefficients estimated by simulating the forced roll, the results of the hybrid method could be in better agreement with the experiment data. This is because the roll period of the forced roll is taken as the roll period of the PR, but the period of the free decay is the natural roll period. So the roll damping characteristics simulated by the forced roll is more close to the actual parametric rolling. Figure 10 shows the influence of the roll period. In this comparison, when simulating the free decay to estimate the roll damping coefficients, the natural roll period of the “free decay-M ” is taken as the roll period of the PR by adjusting the value of Ixx. The results indicate that the roll period could have a significant effect on the prediction of the PR by the hybrid method.

    Fig.10 Comparison of effect of roll period on prediction of PR(F r=0)

    However, according to the basic assumption of the 2nd generation criteria of the IMO[16], this case is unacceptable to this failure model only if the roll amplitude of the PR is larger than 25o. And the vulnerability of a loading condition is judged by an index weighted according to these critical wave-steepness(for different Froude numbers and wave headings) and the wave scatter diagram of North Atlantic. Therefore,from the view point of engineering application, it is more important to ensure the prediction accuracy of the critical wave-steepness that the roll amplitude of the PR is 25o. According to the present result of C11,it is appropriate to adopt the hybrid method (forced roll) for this purpose, and this method could bring a great advantage in the initial design stage especially in the lack of experiment data for a new design. Of course, further studies of the scale effect are needed to predict the roll damping of the full scale ship properly.

    3. Comparison of hybrid, nonlinear dynamics, and direct CFD prediction methods

    In this study, in order to investigate the applicability of the proposed hybrid method, two other methods including a non-linear dynamics method and the direct CFD prediction method are utilized to compare with the hybrid method.

    3.1 Nonlinear dynamics method

    The non-linear dynamics method is adopted as the 2nd check of Level 2 criteria by the IMO[16]for vulnerability judgment of the full scale ships. Due to the fact that the Level 2 criteria are designed for the “early warning” of vulnerability in the design stage, this method is supposed to be simple and conservative predictions are to be obtained. The roll amplitude of the PR is obtained from the analytical solutions of the 1 DOF non-linear dynamics method. So it is very easy to use for engineering applications.

    The 1 DOF motion equation is written as:

    where GZwis the restoring variation due to waves,GMmeanis the mean of the metacentric height variation,GMampis the amplitude of the metacentric height variation,ωφis the natural roll circular frequency,ωeis the encounter frequency; the l2k+1coefficients shall be determined with a least square fit to theGZ curve in calm water,αand γare the non-dimensional damping coefficients that could be estimated by the simplified Ikeda's method[16,18], in the absence of the free decay data. The stable solution could be assumed asis the roll period of the PR,φis the roll angle,GMis the metacentric height in calm water andk!=1· 2…k.

    Applying the averaging method, the roll amplitudeAand the phaseεcould be calculated separately by analytical solutions[16,19], as follows

    After obtaining the complex solutions in the polynomial form of Eq.(18), the stability of the solutions need to be examined using the Jacobian with respect to the averaged system for the eigen analysis. Finally the effective amplitude of the PR could be obtained.

    3.2 Direct CFD prediction method

    In order to improve the forecasting precision of the PR for the optimal design, the best way, in theory,is to carry out a good simulation for the encountered wave surface accounting for the action of the ship, the highly nonlinear restoring forces and the hydrodynamic forces, and the large roll-heave-pitch resonance. Different from the traditional potential methods, the fully nonlinear CFD approach could be a good choice for the art of state method for this purpose. Thus, it is also necessary to carry out a comparison study for the hybrid method with the “state of art” methods such as the 3-D direct CFD approach.

    Fig.11 Computation domain for prediction of PR (F r=0)

    Fig.12 Simulations of PR (Fr=0, wave-steepness: 0.03)

    In the direct CFD prediction method, the same RANS solver is used for the estimation of the roll damping. Most of the parameter settings and the mesh generation also follow the same principles. The computational domain is -4<x<2,-1.5<z<0.5, in the dimensionless coordinates non-dimensionalized by usingLpp. All predictions are made on 2.77×106grids. Since the simulations are carried out for the large roll response of the PR, a cylindrical computational domain is created with the sliding grid for simulating the near field of the ship (see Fig.11(a)). We have -0.6Lpp<x<0.6Lpp, and the radius is 1.2 of the depth of the ship. And the grid in the cylindrical domain is refined from-35oto 35oand 145oto 215o(see Fig.11(b)).

    The 3 DOF motions (roll, pitch and heave) are the free motions for the simulation of the PR. The forward speed is forced to accelerate from zero to the target value following the 1/2sin law, and then keeps as a constant value sailing in head waves. Sway and yaw are assumed to be neglectable.

    Fig.13 Simulations of PR (Fr=0.05, wave-steepness: 0.02)

    Figure 12 shows the time histories ofthe motionresponses and the interactions of the fluid field and the ship at Fr =0.zis the heave amplitude. Figure 13 shows the motion responses atFr =0.05.

    Fig.14 Comparisons of roll amplitudes for PR

    3.3 Comparison of hybrid, nonlinear dynamics and direct CFD prediction methods

    Figure 14 shows the comparisons of the roll amplitudes predicted by different methods. These methods include the 1 DOF nonlinear dynamics method(Method 1-A, 1-B in Table 2), the 3 DOF potential method with semi-empirical or experimental roll damping (Method 2-A, 2-B in Table 2), the hybrid method(Method 2-C, 2-D in Table 2) and the direct CFD prediction method (Method 3 in Table 2). The roll damping coefficients utilized in the nonlinear dynamics method and the 3 DOF potential method are estimated by the free decay data from the experiment and the simplified Ikeda's method separately. The differences of the results obtained by these methods are shown in Table 2. There are four options for Method 2 that adopt the same 3 DOF potential method for simulating the motion response, but four different methods for estimating the roll damping.

    As shown in the results, the 1 DOF nonlinear dynamics method is more conservative than other methods. This is because the nonlinear dynamics method is a highly simplified 1 DOF method and the motion response is obtained by the analytical method. On the contrary, the 3 DOF method and the hybrid method could take account of the nonlinear effect of the Froude-Krylov forces and the effect of heave and pitch in the time domain. Furthermore, the direct CFD prediction method is a fully nonlinear method developed based on the viscous flow theory. This means that in the direct CFD prediction, the viscous effect and the complex nonlinear action between the incident wave and 3 DOF motion of ship are considered.

    On the other hand, the predictions of the PR based on the 1 DOF method and the 3 DOF method for the comparison are both very sensitive to the estimation methods of the roll damping coefficients. Thus the calculation method of the roll damping is a key aspect for the prediction of the PR, and this result agrees with the previous studies.

    The simplified Ikeda's method is very easy to use,but could also cause uncertainties for the damping prediction as compared with the experiment. These uncertainties make the predicted amplitude of the PR deviate from the results predicted based on the free decay data of the experiment, especially for the nonlinear dynamics method (as shown for Method 1-A and Method 1-B).

    In Fig.14, the results of 4 options of Method 2 show the influence of the four estimating methods of the roll damping. The roll amplitudes obtained by Method 2-C (hybrid method) is very close to that obtained by Method 2-B, It indicates that the roll damping estimated by the free decay based on the CFD simulation enjoys a good accuracy for the prediction of the PR, and can be a good option to replace the free decay tests which are currently being carried out.

    On the whole, Method 2-D (hybrid method) and Method 3 (direct CFD prediction method) enjoy the best accuracy for the prediction of the PR. Both of them are much more complicated than the nonlinear dynamics method. They aim to solve for the 3 DOF motion responses, taking account the effect of the roll period on the roll damping. They are considered to be more appropriate as the numerical models for the direct stability assessment of the PR as suggested by the IMO[16].

    Table 2 Methods for predicting PR

    The hybrid method is only based on a weakly nonlinear model, and for the direct CFD prediction method, it is necessary to pay a special attention to the mesh generation and it takes more time to simulate the motion response in regular waves, which is even more complicated than simulating the free decay or the forced roll. It is also necessary to develop a standard procedure for the CFD approach to calculate the roll damping or the direct prediction of the PR. However,this means that sufficient experimental data are needed for validation. For the former, the free decay and forced roll tests have been conducted for many years, no matter whether these ships are vulnerable to the PR. Therefore, we will suffer a more serious data shortage problem, if we adopt the direct CFD prediction method. Therefore, the hybrid method is a better option for the direct stability assessment of the PR.

    4. Conclusions

    A hybrid method is developed and validated for predicting the parametric rolling in this study. Comparison study is also carried out with two methods, and their advantages are shown. The following conclusions are drawn:

    With the hybrid method, the roll damping efficiently can be estimated efficiently and the PR be predicted with a good accuracy. The roll damping estimated by simulating the forced roll is more flexible in simulating the damping characteristics of the PR, thus with the hybrid method, better agreement could be achieved. In addition, with the hybrid method, similar accuracy could be achieved but with far less computational time than the direct CFD prediction method.

    Therefore, the hybrid method could be a good option for the direct stability assessment of the PR in the near future, and could be utilized for the optimization of the ship design, especially in the initial design stage.

    References

    [1] FOSSEN T., NIJMEIJER I. H. Parametric rsonance in dynamical systems[M]. New York, USA: Springer, 2012.

    [2] HASHIMOTO H., UMEDA N. A study on quantitative prediction of parametric roll in regular waves[C]. Proceedings of the 11th International Ship Stability Workshop. Wageningen, The Netherlands, 2010.

    [3] UMEDA N., HASHIMOTO H. and MINEGAKI S. et al. An investigation of different methods for the prevention of parametric rolling[J]. Journal of Marine Science and Technology, 2008, 13 (1): 16-23.

    [4] YANG Chun-lei, ZHU Ren-chuan and MIAO Guo-ping et al. A roll damping prediction method of three-dimensional ships based on CFD computation[J]. Journal of Shanghai Jiaotong University, 2012. 46(8): 1190-1202(in Chinese).

    [5] WILSON R. V., CARRICA P. M. and STERN F. Unsteady RANS method for ship motions with application to roll for a surface combatant[J]. Computers and Fluids,2006, 35(5): 501-524.

    [6] SADAT-HOSSEINI H., STERN F. and OLIVIERI A. et al. Head-wave parametric rolling of a surface combatant[J]. Ocean Engineering, 2010, 37(10): 859-878.

    [7] TURAN O., AYAZ Z. Parametric rolling behaviour of azimuthing propulsion-driven ships[J]. Ocean Engineering, 2008, 35(13): 1339-1356.

    [8] CHANG Yong-quan, FAN Ju and ZHU Ren-chuan et al. Analysis of ship parametric rolling in head sea[J]. Chinese Journal of Hydrodynamics, 2008, 23(2): 204-211(in Chinese).

    [9] ZHOU Yao-hua, MA Ning and GU Jie-chong. Effect of calculation methods of roll moment of inertia on sensitivity prediction of parametric rolling for container ships[J]. Shipbuilding of China, 2013, 54(3): 11-20(in Chinese).

    [10] ZHOU Yao-hua. The prediction of roll damping and nonlinear motion of trimaran[D]. Master Thesis. Harbin, China: Harbin Engineering University, 2010(in Chinese).

    [11] Deng G. B., QUEUTEY P. and VISONNEAU M. RANS prediction of the KVLCC2 tanker in head waves[C]. Proceedings of the 9th International Conference on Hydrodynamics. Shanghai, China, 2010.

    [12] DENG G. B., QUEUTEY P. and VISONNEAU M. Seakeeping prediction for a container ship with RANS computation[C]. Proceeding of the 9th National congress on Hydrodynamics and 22nd National conference on Hydrodynamics. Chengdu, China, 2009.

    [13] LARSSON L., STERN F. and VISONNEAU M. CFD in Ship hydrodynamics-Results of the Gothenburg 2010 Workshop[C]. The Gothenburg 2010 Workshop. Gothenburg, Sweden, 2010.

    [14] ZHOU Yao-hua, MA Ning and SHI Xun et al. Direct calculation method of roll damping based on three dimensional CFD approach[J]. Journal of Hydrodynamics,2015, 27(2): 176-186.

    [15] LU J., UMEDA N. and MA K. Predicting parametric rolling in irregular head seas with added resistance taken into account[J]. Journal of Marine Science and Technology, 2011, 16(2): 462-471.

    [16] IMO. Information collected by the correspondence group on intact stability regarding the second generation intact stability criteria development[C]. Sub-Committee on Ship Design and Construction (SDC) 1st session. London,UK, 2014.

    [17] LI Ji-de. Ship seakeeping performance[M]. Harbin,China: Harbin Engineering University Press, 2007(in Chinese).

    [18] KAWAHARA Y., MAEKAWA K. and IKEDA Y. A simple prediction formula of roll damping of conventional cargo ships on the basis of Ikeda's method and its limitation[C]. Proceedings of the 10th International Conference on Stability of Ships and Ocean Vehicles. St Petersburg, Russia, 2009, 387-398.

    [19] MAKI A., UMEDA N. and SHIOTANI S. et al. Parametric rolling prediction in irregular seas using combination of deterministic ship dynamics and probabilistic wave theory[J]. Journal of Marine Science and Technology,2011, 16(3): 294-310.

    10.1016/S1001-6058(16)60666-2

    December 9, 2014, Revised August 17, 2015)

    * Project supported by the Ministry of Industry and Information Technology of China (Grant No. [2012]533), the Chinese Government Key Research Knowledge-based Ship Design Hyper-Integrated Platform (Grant No. 201335).

    Biography: Yao-hua ZHOU (1984-), Male, Ph. D., Engineer

    Ning MA,

    E-mail: ningma@sjtu.edu.cn

    2016,28(4):617-628

    猜你喜歡
    馬寧耀華
    滿江紅·七月半到韶山
    最高樓·梅花夢
    沁園春·端陽節(jié)祭屈原
    奮楫揚帆 十年輝煌
    星環(huán)上的雜貨鋪(下)
    星環(huán)上的雜貨鋪(中)
    星環(huán)上的雜貨鋪(上)
    常耀華
    流行色(2019年7期)2019-09-27 09:32:58
    我出高考數(shù)學(xué)題(十三)
    翻番的1季度與瘋狂的3月
    成年女人看的毛片在线观看| 草草在线视频免费看| 国产黄片美女视频| 观看免费一级毛片| 禁无遮挡网站| 97超视频在线观看视频| 激情在线观看视频在线高清| 美女扒开内裤让男人捅视频| 亚洲五月婷婷丁香| 亚洲精品粉嫩美女一区| 国产精品影院久久| 亚洲午夜理论影院| 久久这里只有精品19| 国内毛片毛片毛片毛片毛片| 男女下面进入的视频免费午夜| 亚洲人成网站在线播放欧美日韩| 精品一区二区三区av网在线观看| 亚洲色图 男人天堂 中文字幕| 久久久水蜜桃国产精品网| 久久欧美精品欧美久久欧美| 国产成+人综合+亚洲专区| 成人三级做爰电影| 欧美xxxx黑人xx丫x性爽| 午夜福利欧美成人| 九九在线视频观看精品| 精品久久久久久久人妻蜜臀av| 成熟少妇高潮喷水视频| 毛片女人毛片| 亚洲av日韩精品久久久久久密| 男女午夜视频在线观看| 国产成人精品久久二区二区免费| 老司机午夜福利在线观看视频| 久久精品91无色码中文字幕| 琪琪午夜伦伦电影理论片6080| 日韩欧美精品v在线| 欧美日韩瑟瑟在线播放| 国产aⅴ精品一区二区三区波| 非洲黑人性xxxx精品又粗又长| 啪啪无遮挡十八禁网站| 久久久久久久久久黄片| 免费看a级黄色片| 亚洲欧美日韩东京热| 91老司机精品| 精品国内亚洲2022精品成人| 看免费av毛片| 9191精品国产免费久久| 成人18禁在线播放| 午夜激情福利司机影院| 日本a在线网址| 国产精品一及| 午夜福利高清视频| 久久精品国产综合久久久| 老汉色∧v一级毛片| 精品国产超薄肉色丝袜足j| 国产极品精品免费视频能看的| 神马国产精品三级电影在线观看| 一个人看的www免费观看视频| 精品久久久久久久久久免费视频| 操出白浆在线播放| 亚洲av熟女| 亚洲国产欧美人成| 男女床上黄色一级片免费看| 欧美又色又爽又黄视频| 久久国产精品影院| 国产亚洲精品综合一区在线观看| 国产精品久久电影中文字幕| 给我免费播放毛片高清在线观看| 亚洲九九香蕉| 波多野结衣高清作品| 亚洲国产精品sss在线观看| 亚洲中文字幕一区二区三区有码在线看 | 黄片大片在线免费观看| 我要搜黄色片| 美女黄网站色视频| 97碰自拍视频| 女同久久另类99精品国产91| 国产精品99久久久久久久久| 夜夜躁狠狠躁天天躁| 午夜福利成人在线免费观看| 久久性视频一级片| 久久精品综合一区二区三区| 国产精品久久久久久精品电影| 精华霜和精华液先用哪个| 亚洲专区中文字幕在线| 欧美日韩国产亚洲二区| 国产精品一区二区免费欧美| 亚洲av电影不卡..在线观看| 久久这里只有精品19| 男插女下体视频免费在线播放| 亚洲九九香蕉| 88av欧美| 色噜噜av男人的天堂激情| 神马国产精品三级电影在线观看| 成年女人看的毛片在线观看| 又黄又爽又免费观看的视频| 国产探花在线观看一区二区| 午夜成年电影在线免费观看| 最新在线观看一区二区三区| 真人一进一出gif抽搐免费| 男女视频在线观看网站免费| 亚洲精品国产精品久久久不卡| 九九久久精品国产亚洲av麻豆 | 在线永久观看黄色视频| 男女午夜视频在线观看| 999精品在线视频| 亚洲无线观看免费| 成人欧美大片| 久久久久久久久免费视频了| svipshipincom国产片| 亚洲精品456在线播放app | 制服丝袜大香蕉在线| 视频区欧美日本亚洲| 色综合亚洲欧美另类图片| 久久久国产欧美日韩av| 又粗又爽又猛毛片免费看| 国产亚洲精品久久久久久毛片| 超碰成人久久| 一级毛片高清免费大全| 国产精品一区二区免费欧美| svipshipincom国产片| 欧美日韩亚洲国产一区二区在线观看| 女人被狂操c到高潮| 一本一本综合久久| 亚洲av熟女| 90打野战视频偷拍视频| АⅤ资源中文在线天堂| 88av欧美| 国产精品久久久久久精品电影| 久久久国产欧美日韩av| www.自偷自拍.com| 亚洲成人久久性| 国产精品综合久久久久久久免费| av天堂中文字幕网| 熟妇人妻久久中文字幕3abv| 一区二区三区国产精品乱码| 日韩大尺度精品在线看网址| 成人亚洲精品av一区二区| 美女高潮的动态| а√天堂www在线а√下载| 老汉色av国产亚洲站长工具| 一卡2卡三卡四卡精品乱码亚洲| 男人和女人高潮做爰伦理| 香蕉国产在线看| bbb黄色大片| 久久国产精品人妻蜜桃| 淫妇啪啪啪对白视频| 亚洲欧美日韩高清在线视频| 真实男女啪啪啪动态图| 日韩欧美在线二视频| 免费看日本二区| 国产精品亚洲av一区麻豆| 欧美丝袜亚洲另类 | 极品教师在线免费播放| 一夜夜www| 亚洲欧美日韩卡通动漫| 两个人视频免费观看高清| 精品久久久久久久毛片微露脸| 日韩三级视频一区二区三区| 国产一区二区激情短视频| 琪琪午夜伦伦电影理论片6080| 成年女人永久免费观看视频| 黄片小视频在线播放| 中文字幕人妻丝袜一区二区| 国产精品一及| 天堂动漫精品| 亚洲精品粉嫩美女一区| 露出奶头的视频| 19禁男女啪啪无遮挡网站| 在线免费观看的www视频| 国产蜜桃级精品一区二区三区| 他把我摸到了高潮在线观看| 1024香蕉在线观看| 国产精品一及| 丰满的人妻完整版| 国产又色又爽无遮挡免费看| 国产精品香港三级国产av潘金莲| av欧美777| 国产精品久久久久久精品电影| 午夜a级毛片| 午夜福利欧美成人| 国产高清激情床上av| 黄色女人牲交| 他把我摸到了高潮在线观看| 精品国产美女av久久久久小说| 久久精品国产清高在天天线| 精品乱码久久久久久99久播| 亚洲国产欧洲综合997久久,| 久久久久久国产a免费观看| 亚洲狠狠婷婷综合久久图片| 国产一区二区三区视频了| 亚洲精品456在线播放app | 国产蜜桃级精品一区二区三区| 18禁国产床啪视频网站| 男女床上黄色一级片免费看| 久久香蕉国产精品| 香蕉丝袜av| 日日干狠狠操夜夜爽| 日韩av在线大香蕉| 麻豆av在线久日| tocl精华| 久久久久性生活片| 国产97色在线日韩免费| 亚洲人成电影免费在线| 精品无人区乱码1区二区| 首页视频小说图片口味搜索| 久久精品国产99精品国产亚洲性色| 欧美丝袜亚洲另类 | 久久久久久久精品吃奶| 18禁裸乳无遮挡免费网站照片| 一级作爱视频免费观看| 97超级碰碰碰精品色视频在线观看| 国产成人欧美在线观看| 校园春色视频在线观看| 久久久久精品国产欧美久久久| 99热精品在线国产| 亚洲成av人片在线播放无| 国产成+人综合+亚洲专区| 国产乱人伦免费视频| 国产一级毛片七仙女欲春2| 午夜视频精品福利| 嫩草影视91久久| 美女午夜性视频免费| 最近视频中文字幕2019在线8| 两性夫妻黄色片| 99久久综合精品五月天人人| 在线播放国产精品三级| 成人欧美大片| 色综合欧美亚洲国产小说| 免费看日本二区| 国产av麻豆久久久久久久| 国产精品电影一区二区三区| 又黄又粗又硬又大视频| 一二三四在线观看免费中文在| 手机成人av网站| 一区二区三区高清视频在线| 日韩成人在线观看一区二区三区| 国产精品香港三级国产av潘金莲| 在线观看一区二区三区| 变态另类丝袜制服| 国产成人精品久久二区二区91| 嫁个100分男人电影在线观看| 亚洲成av人片免费观看| 两人在一起打扑克的视频| 久久伊人香网站| 精品国产乱子伦一区二区三区| 久久久久免费精品人妻一区二区| 又紧又爽又黄一区二区| 嫩草影视91久久| 色噜噜av男人的天堂激情| 欧美午夜高清在线| 国产精品精品国产色婷婷| 亚洲精华国产精华精| 精品一区二区三区视频在线 | 久久久久九九精品影院| 午夜精品久久久久久毛片777| 亚洲第一欧美日韩一区二区三区| 久久久久久久久免费视频了| 日本熟妇午夜| 999久久久国产精品视频| 99re在线观看精品视频| 免费高清视频大片| 我要搜黄色片| 午夜免费激情av| 欧美3d第一页| 精品国产乱子伦一区二区三区| 国产伦一二天堂av在线观看| 欧美不卡视频在线免费观看| 91字幕亚洲| 日韩欧美一区二区三区在线观看| 老汉色av国产亚洲站长工具| 国产成人精品久久二区二区91| 给我免费播放毛片高清在线观看| 亚洲欧美激情综合另类| 国产精品一区二区免费欧美| 欧美不卡视频在线免费观看| 国产精品一及| 麻豆av在线久日| 精品午夜福利视频在线观看一区| 99久久精品一区二区三区| 变态另类成人亚洲欧美熟女| 久久人人精品亚洲av| 九色成人免费人妻av| 久久天躁狠狠躁夜夜2o2o| 桃色一区二区三区在线观看| 老司机在亚洲福利影院| 岛国视频午夜一区免费看| 中亚洲国语对白在线视频| 久99久视频精品免费| 法律面前人人平等表现在哪些方面| 99久久精品一区二区三区| 99视频精品全部免费 在线 | 久久久久久久精品吃奶| 国产成人aa在线观看| 国产人伦9x9x在线观看| 亚洲一区二区三区不卡视频| 久久国产精品影院| 久久人人精品亚洲av| 九色成人免费人妻av| 欧美在线一区亚洲| 午夜激情福利司机影院| 精品久久久久久久久久免费视频| 曰老女人黄片| 波多野结衣高清无吗| 99热只有精品国产| 嫩草影院入口| 成人三级做爰电影| 99久久精品一区二区三区| 色老头精品视频在线观看| 亚洲精品在线观看二区| 国产欧美日韩精品亚洲av| 亚洲aⅴ乱码一区二区在线播放| 狂野欧美白嫩少妇大欣赏| 久久热在线av| 手机成人av网站| 成人三级黄色视频| 嫩草影院精品99| 欧美成狂野欧美在线观看| 午夜福利成人在线免费观看| 久久精品国产综合久久久| 国产一区二区三区视频了| 国产男靠女视频免费网站| 日韩三级视频一区二区三区| 免费看十八禁软件| 国产亚洲欧美98| 99国产精品99久久久久| 国内精品美女久久久久久| 国内揄拍国产精品人妻在线| 欧美成人性av电影在线观看| 两人在一起打扑克的视频| 看片在线看免费视频| 小说图片视频综合网站| 日韩 欧美 亚洲 中文字幕| 国产一区二区三区在线臀色熟女| 亚洲精华国产精华精| 午夜福利成人在线免费观看| 在线国产一区二区在线| 精品一区二区三区视频在线观看免费| 九九久久精品国产亚洲av麻豆 | 欧美激情在线99| 午夜激情福利司机影院| 淫妇啪啪啪对白视频| 又紧又爽又黄一区二区| bbb黄色大片| xxx96com| 日韩国内少妇激情av| 亚洲精品在线观看二区| 亚洲精品粉嫩美女一区| 亚洲成av人片免费观看| 九色国产91popny在线| 国产黄片美女视频| 午夜福利欧美成人| 亚洲欧美一区二区三区黑人| 精品无人区乱码1区二区| 99久久精品一区二区三区| 色哟哟哟哟哟哟| 日本 av在线| 日韩av在线大香蕉| 成年免费大片在线观看| 国产成人精品久久二区二区免费| 欧美又色又爽又黄视频| 久久香蕉国产精品| 美女cb高潮喷水在线观看 | 一区二区三区高清视频在线| 亚洲精品色激情综合| 午夜精品在线福利| 九色国产91popny在线| 亚洲欧美精品综合一区二区三区| 国产精品av视频在线免费观看| 日韩欧美免费精品| 久久精品国产99精品国产亚洲性色| 国产亚洲精品一区二区www| 国产高清视频在线播放一区| АⅤ资源中文在线天堂| 1024手机看黄色片| 夜夜躁狠狠躁天天躁| 三级国产精品欧美在线观看 | 亚洲aⅴ乱码一区二区在线播放| 日本黄色片子视频| 少妇熟女aⅴ在线视频| 欧美黑人欧美精品刺激| 青草久久国产| 亚洲在线观看片| 我的老师免费观看完整版| 国产精品九九99| 国产美女午夜福利| 中文字幕av在线有码专区| 亚洲av成人不卡在线观看播放网| 1024香蕉在线观看| 美女cb高潮喷水在线观看 | 亚洲欧美日韩无卡精品| 亚洲av第一区精品v没综合| 性欧美人与动物交配| 757午夜福利合集在线观看| 一a级毛片在线观看| 精品久久久久久久末码| 午夜久久久久精精品| 亚洲精品在线美女| 亚洲av成人av| 狂野欧美白嫩少妇大欣赏| netflix在线观看网站| 日韩欧美精品v在线| 欧美激情在线99| 国产精品 欧美亚洲| 九九在线视频观看精品| 美女高潮喷水抽搐中文字幕| 麻豆国产av国片精品| 免费一级毛片在线播放高清视频| 久久亚洲精品不卡| 黄色片一级片一级黄色片| 中文字幕av在线有码专区| 日本与韩国留学比较| 亚洲七黄色美女视频| 欧洲精品卡2卡3卡4卡5卡区| 在线观看66精品国产| 久久久久久久久久黄片| 久久久国产成人免费| 色综合婷婷激情| 天天躁狠狠躁夜夜躁狠狠躁| 亚洲aⅴ乱码一区二区在线播放| 亚洲国产精品久久男人天堂| 国产亚洲精品久久久久久毛片| 五月玫瑰六月丁香| 亚洲成人久久爱视频| 欧美乱码精品一区二区三区| www日本黄色视频网| 久久精品91无色码中文字幕| 男女那种视频在线观看| 国产精品综合久久久久久久免费| 精品久久久久久久毛片微露脸| 久久精品国产99精品国产亚洲性色| 精品国产乱子伦一区二区三区| 三级国产精品欧美在线观看 | 国产又黄又爽又无遮挡在线| 国产一区二区三区视频了| 久久久国产精品麻豆| www国产在线视频色| 久久国产精品影院| 波多野结衣高清作品| 亚洲人成伊人成综合网2020| 成人高潮视频无遮挡免费网站| 国产欧美日韩一区二区三| 精品一区二区三区视频在线 | 日日摸夜夜添夜夜添小说| 亚洲人成网站高清观看| 亚洲美女黄片视频| 亚洲成人久久爱视频| 欧美乱码精品一区二区三区| 国产又色又爽无遮挡免费看| 男插女下体视频免费在线播放| 成人欧美大片| 精品久久久久久,| 免费看十八禁软件| 国产av麻豆久久久久久久| 欧美av亚洲av综合av国产av| 国产精品久久视频播放| 可以在线观看毛片的网站| 亚洲欧美日韩无卡精品| 亚洲精品久久国产高清桃花| 免费一级毛片在线播放高清视频| 精品久久久久久,| 18禁裸乳无遮挡免费网站照片| 国产av麻豆久久久久久久| 搡老岳熟女国产| 青草久久国产| 91九色精品人成在线观看| 国产伦精品一区二区三区视频9 | 国产私拍福利视频在线观看| 亚洲人与动物交配视频| 色老头精品视频在线观看| 亚洲精品美女久久久久99蜜臀| 日本黄色片子视频| 中文亚洲av片在线观看爽| 精品久久久久久久毛片微露脸| 国产成人啪精品午夜网站| 日韩欧美 国产精品| 狠狠狠狠99中文字幕| 午夜精品一区二区三区免费看| 久久人人精品亚洲av| 久久久久九九精品影院| 欧美乱码精品一区二区三区| 色综合亚洲欧美另类图片| 久久这里只有精品19| av天堂中文字幕网| 亚洲无线在线观看| 成人18禁在线播放| 亚洲av成人不卡在线观看播放网| 男女视频在线观看网站免费| 免费一级毛片在线播放高清视频| 一a级毛片在线观看| 精品熟女少妇八av免费久了| 又黄又爽又免费观看的视频| 九色成人免费人妻av| 欧美日韩亚洲国产一区二区在线观看| 1000部很黄的大片| 亚洲aⅴ乱码一区二区在线播放| 成人一区二区视频在线观看| 在线观看一区二区三区| 熟女人妻精品中文字幕| 国产乱人伦免费视频| 国产精品爽爽va在线观看网站| www.熟女人妻精品国产| 99国产极品粉嫩在线观看| 俺也久久电影网| 特级一级黄色大片| 天堂√8在线中文| 美女扒开内裤让男人捅视频| 久久精品人妻少妇| 在线观看午夜福利视频| 香蕉久久夜色| 一二三四社区在线视频社区8| 色av中文字幕| 又粗又爽又猛毛片免费看| 99riav亚洲国产免费| 国产免费男女视频| 床上黄色一级片| 丰满人妻熟妇乱又伦精品不卡| 最好的美女福利视频网| 成人av在线播放网站| av福利片在线观看| 在线观看美女被高潮喷水网站 | 亚洲成人中文字幕在线播放| 18禁黄网站禁片午夜丰满| 动漫黄色视频在线观看| 亚洲国产欧美一区二区综合| 亚洲自偷自拍图片 自拍| 露出奶头的视频| 成人欧美大片| 亚洲18禁久久av| 国产人伦9x9x在线观看| 国产不卡一卡二| 又爽又黄无遮挡网站| 午夜免费成人在线视频| 可以在线观看的亚洲视频| 亚洲一区二区三区色噜噜| 99riav亚洲国产免费| 国产精品1区2区在线观看.| 亚洲av熟女| 日韩国内少妇激情av| 听说在线观看完整版免费高清| 国产1区2区3区精品| 国产一区二区在线av高清观看| 啦啦啦观看免费观看视频高清| 1024香蕉在线观看| 欧美国产日韩亚洲一区| 久久精品夜夜夜夜夜久久蜜豆| 老司机午夜福利在线观看视频| 最好的美女福利视频网| x7x7x7水蜜桃| 97碰自拍视频| 99久久国产精品久久久| 999久久久国产精品视频| av女优亚洲男人天堂 | 一本久久中文字幕| 小说图片视频综合网站| 久久天躁狠狠躁夜夜2o2o| 网址你懂的国产日韩在线| 一本久久中文字幕| 91字幕亚洲| 成人鲁丝片一二三区免费| 日韩欧美精品v在线| 99久国产av精品| 97人妻精品一区二区三区麻豆| 在线看三级毛片| 国产精品久久久久久亚洲av鲁大| 别揉我奶头~嗯~啊~动态视频| 999久久久精品免费观看国产| 日本一二三区视频观看| 亚洲成人免费电影在线观看| 一本久久中文字幕| 夜夜看夜夜爽夜夜摸| 两性夫妻黄色片| 身体一侧抽搐| 国产真实乱freesex| 国内少妇人妻偷人精品xxx网站 | 老汉色∧v一级毛片| 精品一区二区三区视频在线 | 日韩免费av在线播放| 麻豆成人午夜福利视频| 亚洲一区二区三区色噜噜| 成年人黄色毛片网站| 精品国产超薄肉色丝袜足j| 成年人黄色毛片网站| 色哟哟哟哟哟哟| 国产淫片久久久久久久久 | 国产精品自产拍在线观看55亚洲| 一个人看视频在线观看www免费 | 久99久视频精品免费| 色在线成人网| 天堂影院成人在线观看| 国产黄a三级三级三级人| 特级一级黄色大片| 日本 欧美在线| 欧美三级亚洲精品| 女警被强在线播放| 村上凉子中文字幕在线| 成年女人毛片免费观看观看9| 亚洲国产看品久久| 久久伊人香网站| 国产精品免费一区二区三区在线| 色尼玛亚洲综合影院| 亚洲av美国av| 亚洲中文日韩欧美视频| 成在线人永久免费视频| 一本综合久久免费| 精品免费久久久久久久清纯| 天天一区二区日本电影三级| 90打野战视频偷拍视频| av黄色大香蕉| 久久亚洲精品不卡| 制服丝袜大香蕉在线| 亚洲av成人一区二区三| 国产精品av久久久久免费| 无限看片的www在线观看| 黄频高清免费视频| 色老头精品视频在线观看|