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

    Interaction tonal noise generated by contra-rotating open rotors

    2023-05-19 03:39:02WngjinSHUCongongCHENLinDUXingHEXiofengSUN
    CHINESE JOURNAL OF AERONAUTICS 2023年4期

    Wngjin SHU, Congong CHEN, Lin DU, Xing HE, Xiofeng SUN,b

    aSchool of Energy and Power Engineering, Beihang University, Beijing 100191, China

    bResearch Institute of Aero-Engine, Beihang University, Beijing 100191, China

    cAECC Hunan Aviation Powerplant Research Institute, Zhuzhou 412002, China

    KEYWORDSContra-Rotating Open Rotors (CRORs);Interaction tonal noise;Nonlinear Harmonic(NLH);Frequency domain;Acoustic analogy

    AbstractFast and accurate prediction of sound radiation of Contra-Rotating Open Rotors(CRORs) is an essential element of design methods of low-noise open rotor propulsion systems.In the present work, a previous frequency-domain model is extended to predict CROR noise.It builds explicitly the relationship between harmonic loadings and corresponding tonal noise, by which the influential parameters to noise generation can be clearly understood.The real distributions of steady and unsteady blade loadings are calculated by the Nonlinear Harmonic (NLH)method.In the present hybrid approach,both the CFD and acoustic modules are solved in the frequency domain.To assess the accuracy of the developed method, the loading noise of a CROR is calculated and compared against results by using the time-domain FW-H module of NUMECA.The predicted sound directivities by the two methods are in good agreements.The present acoustic model in the frequency domain is proven to be accurate and have high efficiency in far-field noise prediction and data processing.Furthermore, the characteristics of the CROR interaction tonal noise are analyzed and discussed.

    1.Introduction

    The configuration of Contra-Rotating Open Rotor (CROR)propulsion systems shows a preponderance of reductions in fuel burn and low pollutant emission relative to turbo-fan engines1,2.However, the control of noise generated by a CROR is a tough problem, which restricts the development and application of CROR airplanes.

    It is well known that the major noise component of rotating blades is the loading noise if the flight Mach number is small3.The loading noise contains steady and unsteady loading noises.The steady loading noise is at the Blade Passing Frequency(BPF)and its harmonics of a Single-Rotation Propeller(SRP)4, while the frequency of the unsteady loading noise is a linear (sum and difference) combination of BPFs of the front and rear rotors.However,a difference tone has poor radiation efficiency whose mode phase speed is subsonic across the whole span of a blade, and such a mode is a cut-off mode in the case of a ducted rotor5.The noise of a CROR shows significantly higher levels and richer frequency components than those of an SRP due to rotor-rotor interactions, especially in front of the plane of rotation5.Meanwhile,the circumferential mode of an SRP is multiple of the blade number increasing with the number of harmonics; however, low circumferential mode noise can be generated by a CROR due to rotor-rotor interaction.This may greatly influence the radiation efficiency and directivity of CROR noise.

    In the past decades,many experimental measurements have been set up to study the characteristics of the CROR noise.Experiments revealed that the wake introduced spikes into the pressure time history that produced higher noise levels in higher harmonics6, and the unsteady loading noise generated by rotor-rotor interactions made great contributions to the noise pattern.The noise spectrum was found to be richer in higher frequency harmonics of BPF than those of an SRP.The noise levels fore and aft of the plane of rotation were significantly higher for a CROR than those of an SRP due to aerodynamic interaction between the two rotors7.In addition,comparisons between SRP and CROR harmonic levels showed that the largest differences in noise occurred in the second and third harmonics and toward the axis of rotation8.In summary,the noise of a CROR becomes complicated due to the aerodynamic interaction between rotors, and the influential parameters to the Sound Pressure Level (SPL) and directivity pattern are of importance in the design of a CROR9.To this end, it is necessary to develop a reliable and effective acoustic model for CROR noise prediction as well as acoustic sources.

    Nowadays, the most popular methods to predict CROR noise are hybrid approaches10–14, in which the unsteady flow and blade loading are computed by CFD and the far-field sound is obtained by solving the Ffowcs-Williams Hawkings equation15.One challenge is to simulate accurately and effectively the unsteady loading due to rotor-rotor interaction by CFD.To this end,the URANS(Unsteady Reynolds Averaged Navier-Stokes) equation was used by Colin et al.10,11and Delattre and Falissard14.Meanwhile,NUMECA International developed the FW-H solver coupled with the CFD module of the Nonlinear Harmonic(NLH)method which has been extensively tested against various analytical test cases16,17.The instantaneous values of static pressure can be reconstructed by cumulating the harmonics and the mean flow given by the NLH simulation(Fourier reconstruction in time).The application of NLH results in a significant reduction in the CPU time requirements and resource compared to full unsteady simulations.The NLH method and the coupled FW-H solver were applied successfully to predict CROR noise by Ferrante12and Deconinck13,18et al.

    Hanson19,20derived the well-known formulas of tonal noise prediction of an SRP based on the acoustic analogy in the frequency domain, in which a computation of the retarded time was avoided.These formulas showed explicitly the dependence of sound pressure on design parameters and loadings.It is helpful for the understanding of the mechanism and design of low noise.Further, Hanson developed a theoretical framework for analysis and prediction of noise of CRORs.It includes analytical models for the aerodynamic unsteady interaction between rotors and for the noise from the resulting unsteady blade loadings21.However, the acoustic sources of Hanson’s method were given in terms of the lift and drag harmonic coefficients at many strips, which were not discussed in detail.This was limited by the CFD capability in that period.Recently, a series of acoustic theories dedicated to the wake interaction noise of CRORs has been presented22–24, which are based on the helicoidal surface theory of Hanson.An asymptotic analytical model for wake interaction tones was established in the frequency domain, and the effect of blade sweep was investigated22.In the solution of Kingan and Parry,the unsteady aerodynamic response was obtained from a Wiener-Hopf analysis, by which the potential interaction was neglected and the practical interaction between the wake and the rear rotor blade was difficult to predict.

    In this paper, a frequency-domain acoustic model25is extended for CROR noise prediction, coupled with the NLH method26for aerodynamic loading computation.Compared the hybrid approaches mentioned above, the present method can simplify data processing and reduce the required computational time and resource, because both the CFD and acoustic module are solved in the frequency domain.With the aerodynamic loading in the frequency domain of NLH methods, the derived acoustic formulas build the relationship between the harmonic loading and the corresponding unsteady loading noise.It is efficient for understanding of the mechanism of noise generation and useful for CROR noise evaluation in design.

    The layout of this paper is as follows.In Section 2, the acoustic model for predicting the loading noise generated by a CROR in the frequency domain is introduced.To examine the accuracy of the present model, a real CROR geometry is introduced as an example in Section 3.Then verifications and validations of the acoustic model are demonstrated in detail in Section 4.The influential parameters to the SPL and directivity pattern are analyzed and discussed in Section 5.Primary conclusions from this investigation are given in the last section.

    2.Acoustic model

    This paper extends the integral formula by Yang and Wang25to a CROR model,which is based on Goldstein’s integral version of the FW-H equation27.The details are as follows, and the loading noise part can be expressed as

    where U is the forward flight velocity, x1is the axial direction of coordinate system.

    Eq.(4) can be written as follows in the cylindrical coordinate system:

    where A is the stretched parameter of Green function Gω,with the Dirac delta function formula as

    Eq.(5) can be transformed into an inhomogeneous convected wave equation as

    The Green function of an inhomogeneous convected wave equation can be expressed as

    in the cylindrical coordinate system, and Eq.(1) can be expressed as28

    The partial derivative of the Green function of Eq.(10) in the respective direction can be expressed as

    In this paper, the acoustic model of the loading noise generated by the rear blade is developed as an example.The deduction of the loading noise generated by the front blade has the same process.

    For the rear blade,it is assumed that the unsteady loadings which are distributed over the blade surface are caused by the wake and potential interaction of the front blade, and the frequencies of the unsteady loadings are s1B1(Ω1+Ω2);s1=1;2;???.The unsteady loadings on the rear blade can be expressed as27

    Fig.1 Configuration of CROR.

    and the summation of blade loadings on all B2blades can be expressed as

    In the cylindrical coordinate system, there are forces in three directions.The acoustic pressure radiated by the axial force is taken as an example with Eq.(11) and Eq.(12) as follows:

    Letting θrot=θy-Ω2τ, the limits of integration in the surface integral over the rotor blades become independent of τ.Then, we can interchange the order of the integration, so Eq.(17) becomes

    does the interaction tonal noise caused by harmonic loadings exist.Eq.(18) is the same as the circumferential mode theory of an axial compressor defined by Tyler and Sofrin29.

    From the above deduction, the final integral formula expressing the acoustic pressure caused by the axial force component fy1can be obtained.With the same deduction procedure, the acoustic pressures radiated by the radial force component fryand the circumferential force component fθyare also obtained as follows:

    where BPF1;BPF2are the blade passing frequency of front blade and rear blade, respectively.

    In summary, the integrals of the loading noise can be written as

    The loading noise contains steady loading noise and unsteady loading noise.For the rear blade,when s1=0,noise at these frequencies is steady loading noise, and Fs1α(α=y1;ry;θy)in Eq.(19) are steady forces.When s1≠0,noise at these frequencies is unsteady loading noise, and Fs1α(α=y1;ry;θy)in Eq.(19) are harmonic loading (including real and imaginary parts).It is important that the frequencies of unsteady loadings exerting on the rear blade are s1B1(Ω1+Ω2);s1=1;2;???, while the frequencies of unsteady loading noise radiated by the rear blade are s1B1Ω1+s2B2Ω2;s1=1;2;???;s2=0;1;2;???.Naturally, the unsteady loading at a single frequency can radiate more discrete frequencies of unsteady loading noise due to the rotation of the circumferential mode m=s2B2-s1B1.For example,the 1st harmonic loading of the rear blade at the frequency of B1(Ω1+Ω2)can radiate unsteady loading noise at frequencies of B1Ω1+s2B2Ω2;s2=0;1;2;???, with Eq.(20).

    3.Aerodynamic module

    The acoustic model of the loading noise generated by a CROR has been developed in the frequency domain.Acoustic sources are required for the acoustic model to predict the noise of a CROR.In this section, an example of CRORs is introduced,and the calculation of the distributed blade loading is briefly introduced.The aerodynamic module is based on the NLH method developed by He and Ning26,and aerodynamic results are the input of the aeroacoustic module.

    The NLH method can solve a mean and finite(pre-selected)number of blade passing frequency harmonic components of the time-dependent solution of a CROR, ignoring all the higher unsteady components.The NLH method has been implemented in the Navier-Stokes solver of the FINE/Turbo which is applied in this paper and has been successfully validated and tested against basic and turbomachinery test cases30–32.This approach has also been validated in CFD simulation for aeroacoustic input of noise prediction of a CROR by Envia33.To predict the sound of a CROR,the accuracy and reliability of the CFD method to capture the aerodynamic loading are essential.It is well known that the mechanism of blade unsteady loading generation is the same for a CROR and an axial compressor.To validate the accuracy of the aerodynamic method, a reliable measurement of a rotor/stator axial compressor in annular duct is introduced for verification34.Blade parameters at design condition are shown in Table 1.

    The NLH method in FINE/Turbo is applied with the first three harmonics in the simulation of rotor/stator flow.With the time reconstruction in Fourier series,the instantaneous static pressure is obtained.The unsteady force on the mid span of the stator is extracted for each time step.Results of the NLH method are compared with those of the time marching Immersed Boundary Method(IBM)by Chen et al.35and measurements by Hsu et al.34in Fig.2.Vbis the rotor blade wheel velocity;ΔF is the difference between instantaneous static pressure and average pressure.q is the density of perfect air; T′is the period of unsteady flow.

    The black dot line represents the results of the NLH method in this paper.There is a good agreement between the results of the NLH method and experiment measurements indicated by the red point solid line.It is verified that theNLH method has the capability to capture unsteady forces on blade surface.

    Table 1 Compressor and blade parameters at design condition.

    Fig.2 Verification for NLH method.

    3.1.Test case

    A generic 6 × 6 CROR model is used in this paper, and the configuration is shown in Fig.3.

    The diameters are 3.95 m and the numbers of blades are 6,for both the front and rear blades, as shown in Table 2.The front and rear rotors operate at a 100 r/min speed difference.The BPFs of the front and rear blades are BPF1=90 Hz and BPF2= 100 Hz, respectively.The Reynolds number is about 9×105as the reference length is the chord length of the rotor.

    Fig.3 CROR geometric layout.

    Table 2 CROR characteristics.

    Fig.4 Sketch of computational domain.

    Table 3 Computational grids case and spanwise points.

    3.2.Computational domain and mesh

    The CFD mesh generation tool used in this paper is the AutoGridv5TMsoftware developed by NUMECA International, which has been introduced in detail and applied to calculate the aerodynamic input for time-domain formulation for propeller noise prediction36,37.The computational domain of CFD is shown in Fig.4, where Rtipis the tip radius of blade.The far field domain in the axial direction is about 10Rtip, and the far field domain in the radial direction is about 3Rtip.

    The non-reflecting boundary conditions based on a characteristic analysis of the linearized Euler equations38are applied at the rotor-rotor interface.

    Three mesh resolutions have been performed for mesh convergence assessment as shown in Table 3.

    The harmonic amplitudes over the blade surface are obtained by different mesh resolutions.Fig.5 shows the results extracted along the spanwise direction at the leading edge of two blade rows,and the aerodynamic results in this paper present a good mesh convergence behavior.

    There is a fictitious shroud which separates the main channel and the far-field domain.The far-field domain is extended 3 times of the radius of the blade tip.The grid points on the surface of the blade are 67 in the spanwise direction and 55 in the chordwise direction, for each side of both the front and rear blades.The topology of the mesh around the blade in Fig.6 and the values of y+above blade surfaces are almost not exceeding 10 on the first layer of cells.

    3.3.CFD settings

    Fig.5 Harmonic pressure amplitudes along spanwise direction of blade surface at leading edge.

    Fig.6 An O4H topology around rear blade and its leading edge at 50 % span.

    Table 4 Boundary conditions.

    The CFD solver used in this paper is the FINE/Turbo developed by NUMECA International.The flow model is from the NLH method, and the turbulence is modeled by the eddy-viscosity-one-equation Spalart-Allmaras model.In this paper,the first six harmonics of the blade passing frequency of the NLH method are selected.The boundary conditions of the CFD module are shown in Table 4.The blades and the hub are considered as smooth and adiabatic at the solid-wall boundaries.

    With the NLH method and Fourier reconstruction in time,one can obtain the harmonic amplitude and the instantaneous values of the flow field.Fig.7 shows the instantaneous solutions of the entropy field at 50%and 90%spans,respectively.Rotor-rotor interactions can be clearly observed at both radial locations.

    A mean axial thrust of 3.215 × 104N is produced by the CROR, while a mean power of 5.519 × 103N?m is required to drive them.Aerodynamic results used for the acoustic model are the time-average static pressure and the first six harmonic loadings.Higher-harmonic results and details of aerodynamic analysis are not discussed here for brevity.

    4.Verifications of acoustic model

    The predicted CROR sound radiation by the present frequency-domain acoustic model will be compared with the results of the time-domain model, namely the famous FW-H equation and Farassat formulation.In the rotating axial direction, 101 observers are fixed around the CROR, all of which are located 20 m (about 10 times of the blade tip radius) far away from the origin coordinates, with an interval of 1.8°,and the central observer placed off the axis above the front rotor is called Observer 1 as shown in Fig.8.

    Fig.7 Entropy field at 50% and 90% spans.

    Fig.8 Observer locations for aeroacoustic analysis.

    The periodicity of the acoustic signal is determined by the common divisor frequency: 10 Hz (the corresponding time interval is 0.1 s18), for which the number of time steps is 512.Loading noises are firstly predicted by the application of the NLH method and the FW-H solver of NUMECA as Refs.36,37.The FW-H solver can decompose the acoustic signal into thickness noise and loading noise of the front and rear blades, respectively.For Observer 1, the time signal of the loading noise by the front blade is predicted by the FW-H solver as shown in Fig.9.The corresponding spectrum is also obtained by the FFT (Fast Fourier Transform)process.

    Then the extended acoustic model in this paper, namely frequency-domain method, is applied to predict the loading noise at each interaction tonal frequency at all observer locations.The loading noises radiated by the front and rear blades are predicted respectively.Acoustic results of the two methods are compared at all tonal frequencies and 101 observer locations as below.

    Fig.9 Time signal of loading noise of front blade at Observer 1.

    Take the sound radiation of the rear blade as an example.Black square symbols are the loading noise predicted by the frequency-domain method while red circle symbols are the loading noise predicted by the time-domain method as shown in Fig.10, and we define f1=BPF1, f2=BPF2for brevity.Amplitudes below 40 dB are neglected.The comparison shows that there is a good agreement between the frequency- and time-domain methods at every discrete frequency at Observer 1.In fact, the agreement between the two methods for several dominant tones (i.e.,f1+f2, f1+2f2) is excellent, being less than 0.5 dB at Observer 1.

    Then the directivities of the loading noise at different tonal frequencies are predicted to verify the accuracy of the present method as shown in Fig.11.There are good agreements between the frequency-and time-domain methods at all observer locations.

    The present hybrid method can simplify data processing and reduce the required computational time and resource,because both the CFD and acoustic modules are solved in the frequency domain.A reconstruction of the instantaneous values of static pressure over blades is not necessary.One can directly input unsteady loadings, which are obtained in the frequency domain by using the NLH method, into Eq.(19) to calculate the noise.In this paper, the same processor(AMD Ryzen 9 3950X 16-core processor) is chosen to test the time cost by the two methods.With the same CFD process and aerodynamic results, the time cost of the aeroacoustic module is evaluated.Taking the results above as an example,the time-domain method can be used to obtain the time signal and spectrum of sound pressure for 101 observers for the first six harmonics.The frequencies in the linear combination of BPF1and BPF2are extracted as results.Contrastively the frequency-domain method can directly solve the SPL of corresponding frequencies based on the acoustic model developed in this paper.The CFD module costs about 30 CPU hours, and Table 5 shows the time costs of directivities results of the rear blade by both methods.

    5.Results and discussions

    5.1.Amplitudes of acoustic sources

    It should be noted that the rotor-rotor interaction of the front blade is dominated by the potential flow, while the harmonic loading noise of the rear blade is generated by both potential and wake interactions.The concentrated harmonic loading amplitudes by integration over the front and rear blade surfaces are expressed as

    Fig.10 Prediction of loading noise by rear blade at Observer 1 by frequency- and time-domain methods.

    Fig.11 Interaction tonal noise directivity of front blade.

    Table 5 CPU time cost.

    Fig.12 Concentrated harmonic loading amplitudes.

    The noise spectrum of the rear blade radiated by steady and unsteady loadings is expressed in different symbols as shown in Fig.13.As indicated by Eq.(19), the frequency-domain method builds the relationship between the loadings and corresponding tonal noise, which requires higher cost when using the time-domain method.The steady loading noise (black square symbols) at a frequency of s2BPF2decays rapidly with the harmonic s2increasing, and the trend is similar to the results of SRP4.The overall SPL of the harmonic loading noise shows a similar trend corresponding to the amplitude of the harmonic loading in Fig.12.

    Fig.13 Spectrum of loading noise by rear blade at Observer 1.

    It should be emphasized that the SPL of the harmonic loading noise shows various trends even with the same harmonic loading as acoustic source, as shown with different symbols in Fig.13.The circumferential mode m plays a key role in determining the amplitudes and phases of the interaction tones.Variation of the circumferential mode m leads to great changes of the harmonic amplitude of the unsteady loading noise even radiated by the same harmonic loading.For example, the amplitude of the unsteady loading noise radiated by the 1st harmonic loading (red point symbols) reaches a peak at BPF1+2BPF2,m=6,and then falls rapidly as m increases.The trend of higher frequencies radiated by a higher harmonic loading becomes more complicated.

    The acoustic source over the blade surface contains both the amplitude and phase of the harmonic loading.The 1st harmonic loading amplitude and phase of the suction side on both the front and rear blades are shown in Fig.14.The amplitudes on both the front and rear blades are in the same range; however, the phases of the suction side on the front blade are well distributed, while those of the rear blade are not well distributed, which leads to constructive interference of the front blade and destructive interference of the rear blade.Therefore,the concentrated harmonic loading amplitude of the front blade is evidently higher than that of the rear blade as shown in Fig.12.

    The 1st harmonic loading of the front blade is caused by the potential flow effect of the rear blade.The amplitude is concentrated on the leading edge and top region of the pressure side, while the corresponding phase in the same region is well distributed as labeled in Fig.14(a).However,the 1st harmonic loading of the rear blade is mostly caused by the incidence wake of the front blade.The amplitude of the suction side on the rear blade is concentrated on the leading edge and blade tip, and in contrast, the corresponding phase is not well distributed as labeled in Fig.14(b).In Fig.14(c), the amplitudes and corresponding phases of the dash line on both the front and rear blades are presented in detail.Variations of the amplitudes and phases of the front blade are slight; however, the amplitudes of the rear blade are in the range from 100 Pa to 300 Pa and the corresponding phases vary from 0.75 rad to 3.5 rad.As a result,they are in phase of the harmonic loading on the front blade and out of phase on the rear blade.Even the amplitude of element of the rear blade is higher than that of the front blade in Fig.14, the concentrated harmonic loading amplitude, defined by Eq.(22), of the front blade is evidently higher than that of the rear blade as shown in Fig.12.

    For the unsteady loading noises generated by the front and rear blades with the same harmonic number and circumferential mode,the sound radiation is mostly dominated by acoustic source.The sound directivities at 190 Hz (s1=1;s2=1) of both the front and rear blades are shown in Fig.15.The SPL of the front blade is higher than that of the rear blade within the range of 0°-55° and 125°-180° in polar angle.The sound radiation of the front blade is also significant.

    Fig.15 Sound directivities at 190 Hz (BPF1+BPF2)of both frontandrearblades.

    From Figs.12-15, it is emphasized that the SPL of the CROR is not only simply dominated by the strength of the harmonic loadings as indicated by Eq.(19), but also significantly influenced by the circumferential mode and loading distribution.

    5.2.Distribution of harmonic loadings

    To investigate the effect of the distribution of harmonic loadings, the blade surface is divided into strips in the spanwise direction in this paper.Distributed loadings in each strip are concentrated on a point in the center of element, and each point forms the line source model.For the rear blade, when s1= 1 and s2= 1, an unsteady loading noise is generated by the 1st harmonic pressure, and the frequency is f1= 190 Hz.The wavelength λ1is longer than the spanwise length, and the distribution of the 1st harmonic pressure amplitude is shown in Fig.16.

    Results of the line source model and the blade source are compared in Fig.17.There is a good agreement between the line and blade sources.Actually, the wavelength is much longer than the chord length, and the blade source is acoustically compact.The line source model is reasonable in this case.

    Nevertheless,when s1=4 and s2=4,an unsteady loading noise is generated by the 4th harmonic pressure, and the frequency is f2= 760 Hz.The distributions of the 4th harmonic pressure amplitude and wavelength λ2are shown in Fig.18.

    Fig.14 1st harmonic pressure amplitude and corresponding phase on blade surface.

    Fig.16 1st harmonic amplitudes of rear blade and line source model.

    Fig.17 Comparison of noise directivities at 190 Hz.

    The directivity of the unsteady loading noise at 760 Hz is shown in Fig.19.There is an evident difference between the line and blade sources at most observer locations.In this case,the wavelength is not long enough compared to the chord length, so the blade is not acoustically compact in the chordwise direction.The line source model is not reasonable.

    Fig.18 4th harmonic amplitudes of rear blade and line source model.

    Fig.19 Comparison of noise directivities at 760 Hz(4BPF1+4BPF2).

    More interaction tonal noises of the rear blade at three observer locations are shown in Fig.20.The amplitudes below 50 dB are neglected for brevity.A comparison between the blade and line source models shows a big difference, even larger than 5 dB at several interaction tonal frequencies.The drawback of the line source model reveals that it must consider the distribution of harmonic loadings for interaction tonal noise at all observer locations.A blade cannot be treated simply by strips assumption.

    5.3.Sound interference between harmonic loadings in three directions

    The sound pressure can be easily decomposed into contributions of radial, circumferential, and axial directions by using the developed model in Eq.(19).The loading noise at a frequency can be separately obtained for the three components of forces in three directions without additional effort and cost, which is plotted in Fig.21 for the front blade at 100 Hz.It is found that the SPL radiated by the radial loading is the smallest while the circumferential force-related noise dominates the SPL at this frequency.It should be noted that the total SPL is lower than the noises generated by the circumferential loading in the range of 0°-95° and the axial loading in the range of 0°-80°.This is because the noise signal radiated by the loading from each direction has a different phase and amplitude, which may lead to destructive interference in the superposition of sound.Naturally, this mechanism can be applied to reduce the SPL in the design of a CROR.

    The decomposition of the unsteady loading noise of the front blade at 380 Hz is shown in Fig.22.At this frequency,s1=1,s2=1, and the circumferential mode m=s2B2-s1B1=0 in Eq.(18), which leads to an acoustic pressure radiated by the circumferential force pθy=0 in Eq.(19).This means that the circumferential mode propagating in free space is zero, and observations at 0° and 180° polar angles have strong amplitudes.This phenomenon was observed in Block et al.’s experiment9.

    The SPL is dominated by the axial loading when m=0,which shows a typical sound directivity of a dipole in the axial direction.The noise radiated by the radial loading can be negligible except in the vicinity of the planes of rotation.This is partly due to the sharp fall of noise generated by the axial loading in this range.The directivity at the circumferential mode m=0 reveals that there is a clear difference on the SPL between the axis direction and plane of rotation.It is important to choose suitable observer locations for an estimation of the Overall Sound Pressure Level (OASPL), especially in experimental measurements.

    Fig.20 Comparison between models at different observer locations.

    Fig.21 Directivity of unsteady loading noise and contributions of three directional loadings at 100 Hz of front blade.

    Fig.22 Decomposition of sound pressure at three directions at 380 Hz of front blade.

    5.4.Interference of two blade rows

    The interaction tonal noise by the CROR is the sum of the harmonic sound pressures of the front and rear blades at the same frequency.However, the interaction tonal noise at frequency s1B1Ω1+s2B2Ω2of the front and rear blades is radiated from different harmonic loadings.The noise at frequency s1B1Ω1+s2B2Ω2of the front blade is radiated by the s2harmonic loading, while that of the rear blade is radiated by the s1harmonic loading.With an increase of the harmonic number,the amplitude of the harmonic loading decreases as shown in Fig.12.At most tonal frequencies, the total SPL is dominated by the rear blade if s1

    For a tone with the same harmonic number, for example,s1=s2=1, the amplitudes of the harmonic loadings exerted on the front and rear blades are in the same order as shown in Fig.12.This results in a close SPL in the sound directivities of the front and rear blades as shown in Fig.23(a).There exists a constructive interference at a polar angle of 45°-135°.Meanwhile, the total SPL is smaller than that of the front blade at polar angles of 0°-45° and 150°-180°, for which the design of the CROR is beneficial for interference of two blade rows.

    Fig.23 Interferences at several tonal noises.

    Results in Fig.23 show that the interaction tonal noise of the CROR might be interfered constructively or destructively at different polar angles, which is essentially determined by Eq.(19).This is much more complicated than noise of an SRP.The integral of harmonic loading amplitudes, as shown in Fig.12,is helpful in estimating the contribution of the front and rear blades to a tone noise with a combination of s1and s2,although the sound signal is also influenced by the interference among strips of blades.

    6.Conclusions

    In this paper, a hybrid CFD-acoustic method is extended for prediction of loading noise generated by a CROR in the frequency domain.The real distributions of steady and unsteady loadings are calculated by the NLH method.In the present hybrid approach, both the CFD and acoustic modules are solved in the frequency domain.Data processing is simplified,and the required computational time and resource are reduced.The formulas build an explicit relationship between harmonic loading and corresponding tonal noise in the frequency domain of a CROR.With the same acoustic sources by the NLH method, the acoustic model shows good agreements on the sound directivities of all tonal loading noises compared with the FW-H solver in the time domain.

    Then, the characteristics of CROR interaction tonal noise are discussed.The amplitudes of acoustic sources are analyzed.The integrals of harmonic loadings show decreasing trends with an increase of the harmonic number, but the trends of the corresponding noise spectrum are more complex.The potential flow effect on the front blade proves to be important.Results show that the simplification of acoustic sources with strips is unreasonable for high frequencies of interaction tonal noise.The blades of a CROR cannot be treated as acoustically compact.The contribution of harmonic loadings in each direction to the sound directivity is analyzed in detail with different combinations of s1and s2.The characteristics of sound radiation,especially when the circumferential mode m=0,are analyzed by the integrals of the acoustic model.The noise generated by the radial force cannot be completely ignored for m=0.The complexity of the CROR noise is not only in spectrum but also in sound directivity.Compared with an SRP, more attentions should be paid to the sound radiation along the axial direction for a CROR.Finally, the contribution and interference of the front and rear blades are analyzed.The plots of sound directivities show that, at most tonal frequencies, the total SPL is dominated by the rear blade if s1

    The preponderance of this approach is that the model explains the relationship between the harmonic loading and the corresponding tonal noise directly, and data processing is simple.This approach shows a great efficiency and accuracy,and it can be applied to prediction of loading noise with an accurate input.For further investigation, experimental results are required for further verification of this approach.

    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.

    Acknowledgements

    This study was co-supported by the National Natural Science Foundation of China (Nos.52022009 and 51790514), the National Science and Technology Major Project, China (No.2017-II-003-0015), and the Key Laboratory Foundation,China (No.2021-JCJQ-LB-062-0102).

    av国产精品久久久久影院| 精品人妻在线不人妻| 亚洲av二区三区四区| 国产午夜精品久久久久久一区二区三区| 日韩av免费高清视频| 久久97久久精品| 99热这里只有是精品在线观看| 亚洲第一av免费看| 亚洲第一区二区三区不卡| 成人国产麻豆网| 夫妻午夜视频| 一边摸一边做爽爽视频免费| 国产不卡av网站在线观看| 人成视频在线观看免费观看| 爱豆传媒免费全集在线观看| a级毛片黄视频| 国产精品麻豆人妻色哟哟久久| 免费av不卡在线播放| av福利片在线| 国产精品欧美亚洲77777| 国产精品嫩草影院av在线观看| 国产日韩欧美亚洲二区| 久久人妻熟女aⅴ| av免费在线看不卡| 美女cb高潮喷水在线观看| 18禁观看日本| 国产成人精品在线电影| av免费在线看不卡| 精品一区二区三区视频在线| 免费人成在线观看视频色| 免费观看在线日韩| 国产成人免费观看mmmm| 国产精品嫩草影院av在线观看| 国产一区有黄有色的免费视频| 国产精品一区二区在线观看99| 观看av在线不卡| 日韩成人av中文字幕在线观看| 啦啦啦啦在线视频资源| 秋霞伦理黄片| 一区二区三区精品91| 精品一区在线观看国产| 热re99久久精品国产66热6| 日韩在线高清观看一区二区三区| 成人黄色视频免费在线看| 久久人人爽av亚洲精品天堂| 亚洲av不卡在线观看| 欧美激情国产日韩精品一区| 精品国产一区二区久久| 99久久精品国产国产毛片| av在线app专区| 免费高清在线观看日韩| 精品久久蜜臀av无| 亚洲人与动物交配视频| 欧美+日韩+精品| 黄片无遮挡物在线观看| 国产高清有码在线观看视频| 亚洲国产精品国产精品| 久久女婷五月综合色啪小说| 精品久久蜜臀av无| 午夜激情福利司机影院| 国产成人freesex在线| 91aial.com中文字幕在线观看| 99九九线精品视频在线观看视频| 人妻系列 视频| 亚洲精华国产精华液的使用体验| av免费观看日本| 91在线精品国自产拍蜜月| 国产一级毛片在线| 色网站视频免费| 国产成人免费观看mmmm| 亚洲欧美一区二区三区国产| 日韩,欧美,国产一区二区三区| 免费观看a级毛片全部| 啦啦啦中文免费视频观看日本| 国产一级毛片在线| 精品国产一区二区三区久久久樱花| 久久久久久伊人网av| 国语对白做爰xxxⅹ性视频网站| 国产日韩欧美亚洲二区| 日韩一本色道免费dvd| 亚洲欧洲精品一区二区精品久久久 | 啦啦啦中文免费视频观看日本| 中文字幕亚洲精品专区| 国产免费现黄频在线看| 亚洲精品乱久久久久久| 九色亚洲精品在线播放| 人人妻人人爽人人添夜夜欢视频| 亚洲色图 男人天堂 中文字幕 | 我要看黄色一级片免费的| 国产乱人偷精品视频| 国产av码专区亚洲av| 亚洲精品乱久久久久久| 亚洲精品第二区| a级片在线免费高清观看视频| 日韩av在线免费看完整版不卡| 中文精品一卡2卡3卡4更新| 亚洲国产精品一区三区| 精品一区二区三区视频在线| 新久久久久国产一级毛片| 九九在线视频观看精品| 久久久精品94久久精品| 晚上一个人看的免费电影| 少妇高潮的动态图| 国产成人精品婷婷| 婷婷成人精品国产| 人人妻人人澡人人爽人人夜夜| 亚洲国产精品专区欧美| 国内精品宾馆在线| 色5月婷婷丁香| av免费观看日本| 美女内射精品一级片tv| 午夜免费观看性视频| 亚洲精品乱码久久久久久按摩| 成人漫画全彩无遮挡| 国产一级毛片在线| 日本黄大片高清| 男女边吃奶边做爰视频| 婷婷成人精品国产| 制服丝袜香蕉在线| 黄片播放在线免费| 乱人伦中国视频| 午夜福利在线观看免费完整高清在| 国产亚洲欧美精品永久| 久久99热6这里只有精品| 亚洲国产日韩一区二区| 大码成人一级视频| 国产69精品久久久久777片| 黄色配什么色好看| 免费不卡的大黄色大毛片视频在线观看| 中文乱码字字幕精品一区二区三区| 免费看不卡的av| 69精品国产乱码久久久| 亚洲四区av| 人妻 亚洲 视频| 日本与韩国留学比较| 亚洲色图综合在线观看| 久久精品熟女亚洲av麻豆精品| 午夜日本视频在线| 成人无遮挡网站| 国产免费又黄又爽又色| 国产 一区精品| 丝袜美足系列| 另类亚洲欧美激情| 最新中文字幕久久久久| 狂野欧美激情性bbbbbb| 国产免费现黄频在线看| av一本久久久久| 亚洲精品中文字幕在线视频| 天堂中文最新版在线下载| 久久久久人妻精品一区果冻| 久久精品国产a三级三级三级| 国产成人精品久久久久久| 免费av中文字幕在线| 999精品在线视频| 免费av不卡在线播放| 又粗又硬又长又爽又黄的视频| 国产日韩欧美视频二区| 少妇高潮的动态图| 亚洲综合色网址| 国产日韩欧美视频二区| 欧美xxⅹ黑人| 亚洲欧美日韩卡通动漫| 亚洲精品久久成人aⅴ小说 | 建设人人有责人人尽责人人享有的| 亚洲在久久综合| 精品久久久久久电影网| 国产精品99久久99久久久不卡 | 我要看黄色一级片免费的| av视频免费观看在线观看| 一区二区三区精品91| 亚洲婷婷狠狠爱综合网| 欧美激情 高清一区二区三区| 免费播放大片免费观看视频在线观看| 熟女人妻精品中文字幕| 免费高清在线观看日韩| 最新的欧美精品一区二区| 久久99热这里只频精品6学生| 99久久人妻综合| 亚洲av欧美aⅴ国产| 久久精品久久精品一区二区三区| 久久国产精品男人的天堂亚洲 | 久久99热这里只频精品6学生| 丁香六月天网| 一区二区三区精品91| 丰满少妇做爰视频| 久久亚洲国产成人精品v| 日本wwww免费看| 曰老女人黄片| 国产一区二区三区综合在线观看 | 一本一本综合久久| 伦精品一区二区三区| 少妇人妻久久综合中文| 亚洲av综合色区一区| 精品国产一区二区三区久久久樱花| 久久久精品94久久精品| 国产精品麻豆人妻色哟哟久久| 蜜臀久久99精品久久宅男| 亚洲av片天天在线观看| 高清视频免费观看一区二区| 国产成人免费观看mmmm| 老熟妇仑乱视频hdxx| cao死你这个sao货| 一级毛片精品| 一区二区日韩欧美中文字幕| 99久久精品国产亚洲精品| 少妇 在线观看| 成人免费观看视频高清| 黑人欧美特级aaaaaa片| av在线播放免费不卡| 午夜精品久久久久久毛片777| 欧美性长视频在线观看| 午夜福利免费观看在线| 国产精品一区二区在线不卡| 免费少妇av软件| 男人操女人黄网站| 一级a爱视频在线免费观看| 五月开心婷婷网| 色尼玛亚洲综合影院| 日韩制服丝袜自拍偷拍| 欧美成狂野欧美在线观看| 又黄又粗又硬又大视频| 精品少妇黑人巨大在线播放| 亚洲成国产人片在线观看| 国产亚洲一区二区精品| 18禁美女被吸乳视频| 精品一区二区三区av网在线观看 | 国产成人一区二区三区免费视频网站| 精品少妇一区二区三区视频日本电影| 黄频高清免费视频| 中文字幕av电影在线播放| 精品一品国产午夜福利视频| 日韩一卡2卡3卡4卡2021年| 午夜两性在线视频| 国产有黄有色有爽视频| 国产精品一区二区免费欧美| 欧美 日韩 精品 国产| 色94色欧美一区二区| 女人精品久久久久毛片| 一本大道久久a久久精品| 亚洲国产看品久久| 精品福利永久在线观看| 国产激情久久老熟女| 午夜激情av网站| 欧美变态另类bdsm刘玥| 中国美女看黄片| 中文字幕最新亚洲高清| av免费在线观看网站| 蜜桃在线观看..| 青草久久国产| 人人妻人人添人人爽欧美一区卜| 亚洲精品中文字幕在线视频| 亚洲精品粉嫩美女一区| av福利片在线| 黑人操中国人逼视频| 久久久国产欧美日韩av| 亚洲一码二码三码区别大吗| 成人18禁在线播放| 黄网站色视频无遮挡免费观看| 美女高潮喷水抽搐中文字幕| 国产av又大| 一区二区三区国产精品乱码| av线在线观看网站| 国产精品成人在线| 老司机亚洲免费影院| 午夜福利免费观看在线| 免费一级毛片在线播放高清视频 | 中文字幕av电影在线播放| 他把我摸到了高潮在线观看 | 蜜桃在线观看..| 中文字幕色久视频| 国产深夜福利视频在线观看| 欧美人与性动交α欧美精品济南到| 九色亚洲精品在线播放| 人成视频在线观看免费观看| 女性被躁到高潮视频| 久久精品人人爽人人爽视色| 久久国产亚洲av麻豆专区| a级片在线免费高清观看视频| 成人三级做爰电影| 精品亚洲成国产av| 亚洲成a人片在线一区二区| 亚洲国产欧美日韩在线播放| 制服诱惑二区| 777久久人妻少妇嫩草av网站| 香蕉丝袜av| 久久国产亚洲av麻豆专区| 日韩欧美三级三区| 黄色丝袜av网址大全| 丁香欧美五月| 日韩中文字幕视频在线看片| 999久久久国产精品视频| 国产成人影院久久av| 免费日韩欧美在线观看| 成人免费观看视频高清| 91精品国产国语对白视频| 国产精品免费一区二区三区在线 | 中文字幕av电影在线播放| 亚洲成人国产一区在线观看| 黄色丝袜av网址大全| 国产xxxxx性猛交| 色婷婷av一区二区三区视频| 日韩欧美国产一区二区入口| 视频区图区小说| 国产一区二区三区视频了| 一级片免费观看大全| 夜夜骑夜夜射夜夜干| 美女高潮到喷水免费观看| 黄片播放在线免费| 国产精品二区激情视频| 亚洲欧美日韩高清在线视频 | 午夜日韩欧美国产| 考比视频在线观看| 精品国产乱码久久久久久小说| 国产精品久久久久久精品电影小说| 欧美日韩黄片免| 久久久久久人人人人人| 午夜福利欧美成人| 成人手机av| 亚洲国产欧美网| 男女无遮挡免费网站观看| 99热国产这里只有精品6| 在线 av 中文字幕| 女性生殖器流出的白浆| 免费久久久久久久精品成人欧美视频| 一级片'在线观看视频| 亚洲专区国产一区二区| 久久精品熟女亚洲av麻豆精品| 国产av又大| 国产精品99久久99久久久不卡| 亚洲欧美一区二区三区黑人| 操出白浆在线播放| 午夜日韩欧美国产| 一区福利在线观看| 别揉我奶头~嗯~啊~动态视频| 午夜91福利影院| 久久精品国产亚洲av香蕉五月 | 中国美女看黄片| 美女福利国产在线| tube8黄色片| 十八禁人妻一区二区| 窝窝影院91人妻| 搡老岳熟女国产| 如日韩欧美国产精品一区二区三区| 欧美乱码精品一区二区三区| 国产不卡av网站在线观看| 精品国产一区二区久久| 9热在线视频观看99| 岛国毛片在线播放| 自拍欧美九色日韩亚洲蝌蚪91| 中文亚洲av片在线观看爽 | 亚洲精华国产精华精| 99国产精品99久久久久| 久久亚洲精品不卡| 一区二区av电影网| 久热爱精品视频在线9| 午夜福利在线免费观看网站| 欧美一级毛片孕妇| svipshipincom国产片| videos熟女内射| 在线观看免费午夜福利视频| 免费日韩欧美在线观看| 国产精品av久久久久免费| 成人影院久久| 亚洲国产毛片av蜜桃av| 国产色视频综合| 亚洲精品久久午夜乱码| 欧美在线黄色| 欧美日韩一级在线毛片| 黄色视频在线播放观看不卡| 99精品欧美一区二区三区四区| 99久久精品国产亚洲精品| 在线观看66精品国产| 亚洲人成77777在线视频| 最新美女视频免费是黄的| 亚洲国产成人一精品久久久| 纵有疾风起免费观看全集完整版| 亚洲精品久久午夜乱码| 一区福利在线观看| 无人区码免费观看不卡 | 免费在线观看完整版高清| 亚洲国产欧美在线一区| 欧美 日韩 精品 国产| 久久中文字幕人妻熟女| 国产人伦9x9x在线观看| 男女午夜视频在线观看| 嫁个100分男人电影在线观看| 日韩免费高清中文字幕av| www.999成人在线观看| 亚洲精品国产一区二区精华液| 中文字幕制服av| a级毛片在线看网站| 色综合婷婷激情| 两个人免费观看高清视频| 日本撒尿小便嘘嘘汇集6| 在线观看www视频免费| 久久99一区二区三区| 咕卡用的链子| 日韩熟女老妇一区二区性免费视频| 男女高潮啪啪啪动态图| 精品亚洲成国产av| 亚洲av欧美aⅴ国产| 最近最新免费中文字幕在线| 窝窝影院91人妻| 国产精品国产av在线观看| 亚洲国产av新网站| 日韩中文字幕欧美一区二区| 国产有黄有色有爽视频| 色视频在线一区二区三区| 一区福利在线观看| 性少妇av在线| 999精品在线视频| av网站免费在线观看视频| 热re99久久国产66热| 欧美老熟妇乱子伦牲交| 成人三级做爰电影| www.自偷自拍.com| 精品国产乱码久久久久久男人| 91老司机精品| 精品人妻在线不人妻| 色婷婷av一区二区三区视频| 热re99久久精品国产66热6| 99国产精品一区二区蜜桃av | 黄色a级毛片大全视频| 国产真人三级小视频在线观看| 中文字幕人妻熟女乱码| 亚洲五月色婷婷综合| 男女下面插进去视频免费观看| 久久精品亚洲精品国产色婷小说| 亚洲国产精品一区二区三区在线| 黑人巨大精品欧美一区二区mp4| 日本av免费视频播放| 亚洲 国产 在线| 欧美黄色淫秽网站| 欧美变态另类bdsm刘玥| 99精品在免费线老司机午夜| 国产av精品麻豆| 欧美激情久久久久久爽电影 | 国产在线一区二区三区精| 桃花免费在线播放| 久久亚洲精品不卡| 精品少妇久久久久久888优播| 亚洲中文av在线| 天天躁日日躁夜夜躁夜夜| 亚洲av第一区精品v没综合| 国产精品1区2区在线观看. | 一区二区三区精品91| 亚洲中文字幕日韩| 久久久精品国产亚洲av高清涩受| 天天影视国产精品| 久久人妻av系列| 国产精品.久久久| 中文字幕精品免费在线观看视频| 欧美精品人与动牲交sv欧美| 欧美亚洲日本最大视频资源| 女人精品久久久久毛片| 考比视频在线观看| 一级毛片精品| 国产免费福利视频在线观看| 女同久久另类99精品国产91| 国产精品美女特级片免费视频播放器 | 在线观看免费日韩欧美大片| 18禁观看日本| 人妻久久中文字幕网| 欧美激情 高清一区二区三区| 久久 成人 亚洲| 一夜夜www| 一边摸一边做爽爽视频免费| 午夜老司机福利片| av欧美777| 老司机午夜福利在线观看视频 | 天天添夜夜摸| 欧美日韩国产mv在线观看视频| 国产福利在线免费观看视频| 一边摸一边做爽爽视频免费| 久久久精品国产亚洲av高清涩受| 日韩精品免费视频一区二区三区| 日本av手机在线免费观看| 天天添夜夜摸| 少妇 在线观看| 最黄视频免费看| 欧美av亚洲av综合av国产av| 久久人妻av系列| 国产精品一区二区在线观看99| 亚洲午夜精品一区,二区,三区| 国产色视频综合| 亚洲少妇的诱惑av| 国产高清videossex| av有码第一页| 手机成人av网站| 制服诱惑二区| 国产一区有黄有色的免费视频| av天堂久久9| 久久精品人人爽人人爽视色| 成年版毛片免费区| 欧美人与性动交α欧美精品济南到| 色在线成人网| 精品第一国产精品| 久久国产精品大桥未久av| 九色亚洲精品在线播放| bbb黄色大片| 男女无遮挡免费网站观看| 亚洲欧美日韩高清在线视频 | 国产日韩欧美视频二区| 亚洲伊人久久精品综合| h视频一区二区三区| 黄色a级毛片大全视频| 国内毛片毛片毛片毛片毛片| 亚洲中文字幕日韩| 最新在线观看一区二区三区| 亚洲第一青青草原| 欧美乱码精品一区二区三区| 亚洲av国产av综合av卡| 亚洲伊人色综图| 久久久精品94久久精品| 国产成人精品在线电影| 欧美精品亚洲一区二区| 99精品久久久久人妻精品| 日韩视频在线欧美| 亚洲色图综合在线观看| 欧美在线一区亚洲| 狠狠精品人妻久久久久久综合| 大码成人一级视频| 久久久久精品人妻al黑| 在线 av 中文字幕| 欧美日韩福利视频一区二区| 99re在线观看精品视频| 亚洲精品av麻豆狂野| 纵有疾风起免费观看全集完整版| 欧美日韩亚洲高清精品| 国产男女内射视频| 欧美日韩中文字幕国产精品一区二区三区 | 亚洲第一av免费看| 看免费av毛片| 精品福利永久在线观看| 久久人人97超碰香蕉20202| 国产男女内射视频| 成人免费观看视频高清| 久热爱精品视频在线9| 久久这里只有精品19| 国产欧美日韩一区二区三| 王馨瑶露胸无遮挡在线观看| a在线观看视频网站| 精品人妻1区二区| 国产精品一区二区在线不卡| 国产精品自产拍在线观看55亚洲 | a级片在线免费高清观看视频| 91九色精品人成在线观看| 黄片大片在线免费观看| 国产成人系列免费观看| 怎么达到女性高潮| 久久久久久久精品吃奶| 久久久久久人人人人人| 国产精品 欧美亚洲| 色94色欧美一区二区| 麻豆成人av在线观看| e午夜精品久久久久久久| 日韩欧美三级三区| 久久精品亚洲熟妇少妇任你| 欧美日韩福利视频一区二区| 精品人妻在线不人妻| 亚洲av成人一区二区三| 欧美亚洲 丝袜 人妻 在线| 新久久久久国产一级毛片| 久久人人爽av亚洲精品天堂| 在线看a的网站| 18禁观看日本| 欧美精品av麻豆av| 俄罗斯特黄特色一大片| 亚洲精品在线美女| 色综合欧美亚洲国产小说| 99在线人妻在线中文字幕 | 久久久欧美国产精品| 在线观看舔阴道视频| 中文欧美无线码| 99re在线观看精品视频| 热re99久久国产66热| 亚洲人成电影观看| 欧美日韩亚洲综合一区二区三区_| 色老头精品视频在线观看| 十八禁人妻一区二区| 一级a爱视频在线免费观看| 一区二区日韩欧美中文字幕| 国产精品欧美亚洲77777| 一本一本久久a久久精品综合妖精| 国产免费视频播放在线视频| 90打野战视频偷拍视频| 老司机在亚洲福利影院| 国产熟女午夜一区二区三区| 亚洲成人免费av在线播放| 国产精品一区二区在线不卡| 人人妻人人澡人人爽人人夜夜| 五月天丁香电影| 日韩大码丰满熟妇| 每晚都被弄得嗷嗷叫到高潮| 欧美黄色淫秽网站| 国产精品久久久人人做人人爽| 天堂8中文在线网| 俄罗斯特黄特色一大片| 在线观看免费视频日本深夜| 精品国产一区二区久久| 国产不卡一卡二| xxxhd国产人妻xxx| 黄频高清免费视频| 天堂俺去俺来也www色官网| 操美女的视频在线观看| 久久亚洲精品不卡| 免费一级毛片在线播放高清视频 | 日本一区二区免费在线视频| 国产又色又爽无遮挡免费看| 久久ye,这里只有精品| 国产一区二区三区在线臀色熟女 | 十分钟在线观看高清视频www| 蜜桃在线观看..| 热re99久久国产66热| 国产福利在线免费观看视频|