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

    Modeling of microstructure evolution of AZ80 magnesium alloy duringhot working process using a unified internal state variable method

    2022-07-14 08:56:16ZexingSuChaoyangSunMingjiaWangLingyunQianXintongLi
    Journal of Magnesium and Alloys 2022年1期

    Zexing Su, Chaoyang Sun,d,?, Mingjia Wang, Lingyun Qian,?, Xintong Li

    b Beijing Key Laboratory of Lightweight Metal Forming, Beijing 100083, China

    c College of Nuclear Equipment and Nuclear Engineering, Yantai University, Yantai 264005, China

    d Shunde Graduate School of University of Science and Technology Beijing, Foshan, 528000, China

    Abstract In this paper, a unified internal state variable (ISV) model for predicting microstructure evolution during hot working process of AZ80 magnesium alloy was developed.A novel aspect of the proposed model is that the interactive effects of material hardening, recovery and dynamic recrystallization(DRX)on the characteristic deformation behavior were considered by incorporating the evolution laws of viscoplastic flow, dislocation activities, DRX nucleation and boundary migration in a coupled manner.The model parameters were calibrated based on the experimental data analysis and genetic algorithm (GA) based objective optimization.The predicted flow stress, DRX fraction and average grain size match well with experimental results.The proposed model was embedded in the finite element(FE)software DEFORM-3D via user defined subroutine to simulate the hot compression and equal channel angular extrusion (ECAE) processes.The heterogeneous microstructure distributions at different deformation zones and the dislocation density evolution with competitive deformation mechanisms were captured.This study can provide a theoretical solution for the hot working problems of magnesium alloy.

    Keywords: AZ80 magnesium alloy; Internal state variable model; Microstructure evolution; Dynamic recrystallization; Hot working process; Finite element simulation.

    1.Introduction

    Magnesium alloys have exhibited considerable application prospects in the aerospace and automobile lightweight fields owing to the high specific strength, excellent damping characteristics and favorable electromagnetic shielding capability[1,2].However, due to the few slip systems and poor formability, extensive applications of these alloys are restricted to die-cast components [3,4].Bulk formation methods including extrusion and forging are carried out at elevated temperature to activate more slip systems and improve workability [51].During hot working process, the dynamic recrystallization (DRX) mechanism is beneficial to material softening and grain refinement [5][52].Unfortunately, the microstructure evolution is sensitive to deformation parameters [6][53],which increases the control difficulties in mechanical properties and forming quality.It is critical to clarify the microstructure evolution affected by processing parameters, which is the premise for the workability optimization and the production of high-quality magnesium alloy parts.

    At present, plentiful constitutive modeling research has been conducted to investigate the material deformation characteristics.For the understanding of flow behavior,the empirical approaches including Johnson-Cook and Arrhenius models show simplicity in the determination of material constants,which have been widely used [7-9].However, the nonlin-ear relationship among material parameters brings difficult to the prediction of hot deformation characteristic beyond the experimental conditions.In contrast, the advanced statistical method such as neural network model shows considerable expertise in self-learning, which can accurately predict the flow stresses within wide deformation conditions [10,11].It is a pity that the prediction capability of neural network model largely depends on the characteristic variables and high quality data.It should be noted that the empirical and advanced statistical methods mainly select strain rate, temperature and plastic strain as macroscopic state variables, which offers less physical insight in metallurgical phenomena.To understand the microstructure evolution of alloys with low stacking fault energy (SFE) including magnesium alloy [12], titanium alloy [13]and nickel-based alloy [14], the Avrami equation is used and contains the DRX fraction as a microstructural variable [15].However, the separating application of macroscopic and microcosmic models ignores the relationship between microstructure evolution and mechanical response.Therefore, a physically-based constitutive model for hot working process,which can relate the deformation behavior to microstructure characteristic, is urgently needed.

    Compared with phenomenological modeling approaches,physically-based model involves internal state variables(ISVs) including dislocation density, recrystallized fraction and average grain size (AGS), which can describe the hot deformation characteristic coupled with different microstructure mechanisms.The evolution of physical variables at different scales has a clear correlation with time and strain, which can easily achieve the interactivity and integration with finite element (FE) platform [16].Lin et al.studied the grain evolution with recrystallization effect during the rolling processes of C-Mn steel through a dislocation density based ISV model [17].Li et al.set up a physically-based model for titanium alloy, which involves the physical variables including dislocation density, damage and globularization [18].Wang et al.established the constitutive equations for Inconel 740 superalloy, which includes the grain evolution, critical strain and volume fraction of DRX behavior during hot deformation process [19].Mohamed et al.developed a viscoplastic constitutive framework coupling dislocation density and damage to describe the deformation and failure characteristics in hot stamping of aluminum alloy [20].

    At present, most constitutive models for magnesium alloys use classical viscoplastic equations based on work hardening (WH) and dynamic softening mechanisms, which lacks an in-depth insight of microstructure evolution.Although the flow behavior coupling microstructure mechanisms has been widely studied via the unified ISV approach among different materials such as titanium, aluminum and nickel-based alloys.An accurate ISV model for hot working of magnesium alloys is rarely reported, and is urgently needed for the material processing applications.Xu et al.put forward a unified ISV model of Mg-RE alloy by incorporating dislocation density, recrystallized grain size, un-recrystallized grain size and DRX fraction [21], which provides theory guidance for the hot working process of magnesium alloy.However, the characteristic flow behavior of magnesium alloy was not revealed, especially the dynamic softening with the DRX sensitivity at different deformation conditions.In our previous research, a physically-based model involving dislocation density and DRX fraction in hot working of AZ80 magnesium alloy was developed [22].It is a pity that the mechanical response with the interactions of different deformation mechanisms and grain size evolution has not been involved.And the modeling applicability for an actual processing problem combined with experimental validation is not considered.

    Fig.1.A typical flow stress curve coupling various microstructural mechanisms in hot compression process of AZ80 magnesium alloy.

    In the current study,firstly a unified ISV model for predicting microstructure evolution during hot working process of AZ80 magnesium alloy was established.Especially, the flow stress characteristic with microstructural effect was reasonably described by incorporating the evolution laws of viscoplastic flow, dislocation activities, DRX nucleation and boundary migration.The experimental flow stress, DRX fraction and AGS matched well with calculated results.Secondly, the FE applicability for hot compression process based on the secondary development technology was analyzed.The heterogeneous distribution of ISVs in the deformed specimen and the dislocation density evolution with competitive deformation mechanisms were well captured.Finally, the ISV modeling approach was used to simulate equal channel angular extrusion (ECAE) process.The microstructure evolution under complex thermomechanical boundary condition at characteristic zones was revealed.This study can provide a solution for the hot working problems of magnesium alloy.

    2.The unified ISV model

    2.1.Modeling viscoplastic flow coupling microstructural mechanisms

    Fig.1 shows a typical flow stress curve coupling various microstructural mechanisms in hot compression process of AZ80 magnesium alloy.The overall stressσis composedof the viscoplastic flowσv, isotropic hardeningHand initial yieldk[23].The general flow characteristic shows an initial WH to the peak and a subsequent softening, which is dominated by the occurrence of DRX [24].Generally, DRX process needs to achieve a critical strainεc[25].Dynamic recovery (DRV) appears in the whole deformation stage, and static recovery (SRV) occurs at elevated temperature in the heating and insulation stages, both of which play a competitive role to DRX in the consumption process of stored energy[26].At the grain level, the grain growth in atomic diffusion process and the grain refinement in DRX process result in the evolution of AGS, which can also affect the material behavior.Therefore,the microstructural effects including SRV,WH,DRV, DRX and grain size should be considered in modeling flow behavior.

    The stress-strain relation is given by Hooke’s law [27]:

    whereEdenotes the Young’s modulus, and its values at different temperatures can be calculated by JmatPro software;εTandεpdenote the total strain and plastic strain, respectively.

    To describe the characteristic metal flow within wide deformation conditions and relate it to microstructure evolution,a hyperbolic-sine type viscoplastic flow equation is introduced[21]:

    The evolution rate of isotropic hardening stress in Eq.(2)is defined as [28]:

    whereBdenotes the material hardening constant, which is temperature-dependent.denotes a normalized concept for dislocation density [29]:

    whereρ0andρrepresent the initial and deformed dislocation densities, respectively.With the development of hot working process from the initial state to saturated state,the normalized dislocation density varies from 0 to 1.

    The evolution rate of normalized dislocation density is expressed as [30]:

    whereA4,A3,δ1,δ2,δ3andδ4are material constants.The first term represents the storage and DRV of dislocations,and(d/d0)δ1describes the grain size sensitivity on the dislocation multiplication.The second term denotes the SRV of dislocations, andCdenotes the SRV coefficient, which is temperature-dependent.The third term represents the DRX,andSdenotes the DRX fraction relevant to the evolution of new grains, which is clarified in the next section.

    2.2.Modeling DRX nucleation and grain boundary migration

    During the hot deformation of low SFE metals, the formation and migration of high angle grain boundaries are driven by the stored energy.Then the deformed microstructure is replaced by the recrystallized grains and results in a decrease of dislocation density.The DRX nucleation rate is sensitive to temperature and strain rate, which can be designated as[31,32]:

    whereN0andbare material constants;Qndenotes the recrystallization activation energy.

    Based on the DRX kinetics equation in the literature [33],by coupling the nucleation rate and critical strain of DRX process with normalized dislocation density, the evolution of recrystallized fraction is designated as:

    whereNis a parameter related to temperature;λ1is a material constant.εcdenotes the DRX critical strain relevant to deformation conditions, which can be determined by Zener-Hollomon (Z) parameter [34]:

    xin Eq.(7) represents the incubation period of DRX initiating, which is given by [35]:

    whereA0is a material parameter related to temperature.

    The grain boundary migration in hot working process contains the grain growth with thermal effect and the grain refinement with DRX, which affects the variation of AGS [36].The grain growth associated with grain boundary migration is expressed as [37]:

    whereσsurfandMdenote the energy density and mobility of grain boundaries, respectively.TheMσsurfterm is affected by temperature, thus Eq.(10) can be designated as:

    whereG1is a material constant relevant to temperature;φ1is a material constant.

    The grain refining rate coupling DRX fraction is expressed as [38]:

    whereλ2andφ2are material constants.G2is a material parameter related to temperature.

    By taking the relative grain size ˉd=d/d0into account and combining Eq.(11) with Eq.(12), the AGS evolves as:

    The two terms represent the grain growth and grain refinement, respectively.

    3.Model calibration and validation

    3.1.Calibration of model parameters

    The developed constitutive equations can be formulated as:

    B, k, C, N, A0,G1andG2are the temperature-dependent variables with thermal activation effect, which can be formulated as Arrhenius-type relations [19,26]:

    whereB0,k0,C0,N0,A00,G10andG20are material constants;QB, Qk, QC, QN, QA0,QG1andQG2are the corresponding activation energies.

    For the calibration of model parameters, firstly the experimental verifiability and physical meaning need to be considered in the corresponding equations.Especially, the critical strain for DRX initiation in Eq.(7) was obtained according to the criterion of irreversible thermodynamics proposed by Poliak and Jonas [39,40].The inflection points in the lnθ-εcurves (θdenotes the strain hardening rate,θ= dσ/dε) indicate the DRX initiation.On the basis of the literature [21],the function betweenεcandZparameter in Eq.(8) was determined by the linear regression analysis.The initial values of AGSd0, initial yield stressk0and Young’s modulusE0were obtained from the experimental data.

    Fig.2.Cubic polynomial fitting between SRV coefficient (C0) and temperature (T).

    For the material parameters that are difficult to determine by the experimental data, the initial values were selected based on the reference papers [33,41], and then the lower and upper limits of 25%?200% were defined for the initial values.Since many ISVs in the equations(such as normalized dislocation density, AGS and DRX fraction) are expressed in differential forms.An optimization solution based on genetic algorithm (GA) [42]was adopted to solve the nonlinear relationship among different material constants.The MATLAB GA toolbox was applied to the minimization of residuals between calculated and experimental flow stress, DRX fraction and AGS.The detailed optimization process has been reported elsewhere [38,43].It should be noted from Fig.2 that the stress softening effect at low temperature is much more significant compared with that at high temperature.The stress softening differences at various temperatures can be attributed to the SRV effect in the heating and insulation process, which further affects the microstructure evolution as well as flow characteristics during the hot working process.Taking this physical meaning into account, the SRV coefficientC0at a fixed temperature was firstly determined.Then a cubic polynomial fitting betweenC0andTwas employed,as shown in Fig.2.The empirical equation was given by:C0= ?157.818 + 1.569T?0.0053T2+ 6.08 × 10?6T3,and the determination coefficient(R-square)was calculated as 0.98841, which indicates a good fitting effectiveness.The calibrated model parameters with the physical significance and feasibility are listed in Table 1.

    Table 1The calibrated model parameters with the physical significance and feasibility.

    3.2.Validation of flow stress

    Fig.3(a)-(c) show the experimental (symbols) and calculated (lines) true stresses during hot compression process at 0.1 s?1, 0.01 s?1and 0.001 s?1, respectively.The flow characteristics including WH and dynamic softening can be well reflected, and the calculated results match well with experimental data.

    To quantitatively evaluate the modeling accuracy, the root mean square error (RMSE), average absolute relative error(AARE) and correlation coefficient (R) are introduced as:

    Fig.3.Comparison between experimental (symbols) and calculated (lines) true stresses during hot compression process: (a) 0.1 s?1, (b) 0.01 s?1, (c) 0.001 s?1 and (d) correlation analysis.

    whereNdenotes number of sampling point;EiandCirepresent experimental and calculated true stress, respectively; ˉEand ˉCdenote the averages ofEiandCi,respectively.Fig.3(d)shows the correlation analysis between experimental and calculated true stress.Most of stress data locate near the best fitting line,while a few stress data deviate slightly from the best fitting line at high stress level.The calculatedRMSE, AAREandRare 2.38 MPa, 3.23% and 0.994, respectively, indicating a low deviation and high correlation between experimental and calculated true stresses.Therefore, the characteristic flow behavior with the WH and dynamic softening effects can be well presented by the unified constitutive model.

    3.3.Validation of microstructure evolution

    To verify the microstructure evolution under different deformation conditions, the microstructure characterization was conducted at the center zone of specimen by electron backscatter diffraction (EBSD) and optical microscopy(OM) tests.The as-extruded billet was annealed before the hot compression process.The initial microstructure with near-equiaxial static recrystallized grains can be observed in Fig.4(a).The AGS was determined to be 20.8 μm.The deformed microstructures at 300 °C and 0.001 s?1with true strains of 0.1, 0.3 and 0.916 are shown in Fig.4(b)-(d), respectively.At lower true strains, the necklace structure composed of DRX grains appears at initial grain boundaries, as shown in Fig.4(b) and (c).As true strain increases to 0.916,DRX fraction increases and AGS decreases, resulting in a refined microstructure with full DRX grains, as shown in Fig.4(d).Increasing deformation degree can accumulate deformation stored energy and provide driving force for grain boundary migration [44].Fig.4(e) shows the comparisonsbetween experimental and calculated results of DRX fraction and AGS.As deformation degree reaches DRX critical strainεc, the AGS curve significantly decreases with the increasing DRX fraction curve.The experimental DRX fraction and AGS match well with the corresponding calculated curves.

    Fig.4.The microstructure evolution at different deformation conditions: initial state (a) and deformed state (300 °C and 0.001 s?1) with true strains of 0.1(b), 0.3 (c) and 0.916 (d); (e) Comparisons between experimental (symbols) and calculated (lines) results of DRX fraction and AGS.

    To verify the AGS under different strain rates and temperatures, the deformed specimens with true strain of 0.916 were characterized.Fig.5(a) shows an inhomogeneous microstructure composed of small and large grains at 300 °C and 0.01 s?1.As strain rate decreases to 0.001 s?1, a refined and homogeneous microstructure containing complete DRX grains can be observed in Fig.5(b).Low strain rate provides sufficient time for DRX initiating, and low temperature can inhibit the grain growth.As temperature increases to 400 °C,grains grow significantly with enhanced mobility, as shown in Fig.5(c).In Fig.5(d), a necklace structure composed of numerous DRX grains appears at grain boundaries at 400 °C and 0.1 s?1, which is due to the incomplete recrystallization kinetics at high strain rate.Fig.5(e) shows the comparison between experimental and calculated AGS.With the help of Image Pro-Plus software 6.0, three representative locations were selected, and the linear intercept method was used to measure the AGS of each location.Then the standard deviation range and the average value of the experimental AGS were obtained.The calculated AGS locates in the range of experimental results with relatively low standard deviation.The proposed constitutive model can capture the microstructureevolutions including DRX fraction and AGS under various deformation conditions.

    Fig.5.The microstructure evolution with a true strain of 0.916: (a) 300 °C, 0.01 s?1, (b) 300 °C, 0.001 s?1, (c) 400 °C, 0.001 s?1, (d) 400 °C, 0.1 s?1 and(e) Comparison between experimental (black symbol) and calculated (red symbol) AGS.

    4.Results and discussion

    4.1.FE simulation of hot compression process

    Based on the secondary development technique,the constitutive equations were embedded in FE software DEFORM-3D via the user defined subroutine to simulate the hot compression experiment.The cylindrical specimen adopted the dimension ofΦ8 × 12 mm, and the tetrahedron element mesh was used.The shear friction mode was applied between anvil and specimen, and a friction coefficient of 0.3 was adopted.The deformation conditions were set as 300 °C and 0.001 s?1,with a true strain of 0.916 (corresponding to a reduction of 7.2 mm).

    Fig.6(a) and (b) show the experimental and simulated compressed specimens, respectively.The obvious bulging phenomenon can be observed due to the effect of friction.The experimental and simulated section sizes at different zones were compared, as shown in Table 2.The experimental section sizes at the end and center zones match well with the simulated results.The maximum relative error was calculated to be 1.30%.The established FE model can present the hot compression process of AZ80 magnesium alloy.

    Table 2Experimental and simulated section sizes at different locations.

    Table 3The decomposition of dislocation density evolution.

    The comparison of the experimental and simulated loads during hot compression process at 300 °C and 0.001 s?1is shown in Fig.7.The experimental load shows an initial increasing stage and remains stable at medium strain level,followed by a continuous increasing stage with the increase of deformation degree.The experimental load accords with the simulated result, while a few experimental points deviate slightly at high strain level, which can be attributed to the increasing friction effect at the specimen-anvil interface.

    Fig.6.Experimental and simulated results of compressed specimen.

    Fig.7.Experimental and simulated loads during hot compression process at 300 °C and 0.001 s?1.

    4.2.Heterogeneous ISV distribution

    Fig.8 shows the simulated distributions of various ISVs in compressed specimen.In Fig.8(a), (b), effective strain locally accumulates at center zone P1, and then decreases from bulging zone P2 to end zone P3.With the increase of stroke to 7.2 mm, effective strains at different deformation zones reach 1.4,0.7 and 0.3,respectively.The inhomogeneous distribution of effective strain is caused by the different stress states with friction at specimen-anvil interface [45].The distributions of DRX fraction and AGS are shown in Fig.8(c),(d)and(e),(f),respectively, which accords with the distribution of effective strain.The zone with larger deformation exhibits higher DRX fraction and lower AGS.Maximum DRX fraction of 0.98 and minimum AGS of 3.9 μm appear at center zone P1, corre-sponding to a maximum effective strain.High plastic strain provides much stored energy and promotes the activation of DRX process, which results in refined microstructure with high DRX fraction [46].The evolving forms of microstructural ISVs simulated by FE model at P1 and P2 accord with the results obtained by MATLAB in Fig.4.An interesting phenomenon is that the state variable evolution at end zone P3 differs from the other two zones.As the stroke increases to~5 mm, the state variables at P3 accelerate to rise.This can be attributed to the increasing friction effect on deformation and microstructure evolution at high strain level.

    To verify the FE model during hot compression process,the EBSD characterization was supplemented at the bulging and end zones, as shown in Fig.9.A wide observation area was adopted to obtain the statistical accuracy of heterogeneous microstructure with relatively low DRX degree.Combined with the microstructure at the center zone in Fig.4, the statistical AGS and DRX fraction at different zones were plotted in Fig.8.In the zones of P1, P2 and P3, the AGSs were 3.9, 4.9 and 6.3 μm, respectively; the corresponding DRX fractions were 0.98, 0.88 and 0.68, respectively.The heterogeneous microstructure distributions at different deformation zones are in good agreement with the predicted results.

    4.3.Dislocation activities with microstructural mechanisms

    The normalized dislocation density provides a physical insight for microstructure evolution in hot working.In this section, the dislocation density evolution is discussed based on the interactions among various deformation mechanisms.The evolution rate of normalized dislocation density in Eq.(5) is divided into three parts, including the WH as well as DRV,SRV and DRX.The decomposition of dislocation density evolution is shown in Table 3.

    Based on the FE simulation of hot compression, the dislocation density evolution at 350 °C and 0.1 s?1are obtained,as shown in Fig.10.The contributed proportions of different microstructure mechanisms are calculated in percentage form.The normalized dislocation density ˉρrises rapidly to a peak and decreases gradually,which is in consistence with the variation of true stress in Fig.3.During the initial deformation stage(ε<εc),dislocations proliferate rapidly,and DRV occurs to annihilate the entangled dislocation network, which results in an initial decreasing proportion of WH+DRV term.The increasing normalized dislocation density ˉρa(bǔ)t this stage provides the deformation storage energy for the subsequent DRX process.As the deformation reaches the DRX critical strainεc, the DRX initiates and then evolves sufficiently, which leads to the annihilation of abundant dislocations.Thus the proportion of WH+DRV term further decreases significantly,and the proportion of DRX term increases rapidly.When the evolution rate of DRX term exceeds WH+DRV term, the normalized dislocation density continuously decreases.By involving the effects of different deformation mechanisms including WH+DRV and DRX, the proportion of SRV term shows a similar evolution characteristic to the normalized dislocation density.

    Fig.10.The simulated dislocation density evolution with contributed proportions of different microstructure mechanisms including WH+DRV, SRV and DRX at 350 °C and 0.1 s?1.

    Fig.8.The simulated distributions of effective strain (a) (b), DRX fraction (c) (d) and AGS (e) (f) in compressed specimen at 300 °C and 0.001 s?1.

    Fig.9.Microstructures in bulging zone P2 (a) and end zone P3 (b) at 300 °C and 0.001 s?1.

    Fig.11 shows the simulated dislocation density evolution with contributed proportions of different microstructure mechanisms at different deformation conditions.As temperature decreases and strain rate increases, normalized dislocation density and the contributed proportion of WH+DRV both increase, as shown in Fig.11(a) and (b), respectively.High strain rate and low temperature contribute to increasing critical shear stress and decreasing fraction of mobile dislocations,which is beneficial to the multiplication and entanglement of dislocations.In Fig.11(c), the contributed proportion of SRV increases with increasing temperatures and decreasing strain rates.High temperature and low strain rate increase SFE and provide a stronger thermal fluctuation for activation, which decreases the energy barrier of dislocation annihilation [47].With deformation storage energy consumed by recovery at high temperature and low strain rate, the contributed proportion of DRX decreases, as shown in Fig.11(d).At one giventemperature (such as 350 °C), the decrease of strain rate results in decreasing DRX critical strain and contributes to the activation of DRX process, which leads to the intersections among DRX curves with various strain rates.The competition between recovery and DRX process with effects of various deformation conditions can also be intuitively explained via true stress curves in Fig.3.The small difference between steady stress and peak stress exists at high temperature and low strain rate (such as 400 °C and 0.001 s?1), which is similar to a recovery type curve.

    Fig.11.The simulated dislocation density evolution (a) with the contributed proportions of WH+DRV (b), SRV (c) and DRX (d) at various deformation conditions.

    Fig.12.Experimental setup (a) and FE model (b) of ECAE process.

    5.Application in ECAE process

    5.1.Experimental and numerical scheme

    During the complicated hot working processes of magnesium alloy such as forging and extrusion, the material behavior and microstructure evolution are sensitive to the nonuniform deformation, stress states and thermal conditions.Thus the reliability of the proposed constitutive model need to be validated via the complex hot working problems.In this section, the FE modeling approach was used to simulate the ECAE process.Meanwhile, the extrusion experiment was carried out for validation.

    Fig.12(a) shows the experimental setup of ECAE process.The ECAE moulds including the top and bottom dies were assembled and then fixed on a 600 ?N hydraulic press.The channel in the bottom die bents through an angular of 120°The channel diameter is 15 mm, and the dimension of workpiece isΦ15 × 40 mm.The inlet and outlet diameters of bottom die were determined to be 20 mm to facilitate the entrance and exit of the workpiece.Before the ECAE process,the die-sets were heated by the heating coil and held for 2 h,and the insulation was used for the reduction of heat dissipation.The workpiece was heated by the heating furnace and held for 30 min.The mold cavity and workpiece were lubricated using molybdenum disulfide (MoS2).The workpiece was placed in the cavity and held for ~10 min to maintain stable temperature distribution.The extrusion parameters were set to be 380 °C and 3 mm/s.The stroke of top die was set to be 40 mm, and the ECAE process lasted for ~40 s.After ECAE process, the workpiece was water-quenched.

    Fig.13.Experimental and simulated loads during ECAE process.

    Fig.12(b) shows the FE model of ECAE process.The dies and workpiece were set as rigid bodies and plastic body,respectively.The die material adopted H-13 steel.The boundary conditions of heat transfer (including convection and radiation) among the workpiece, dies and ambient atmosphere were considered.A heat transfer coefficient of 11 N/s/mm/°C was set between workpiece and dies, while the heat transfer coefficient to atmosphere was 0.02 N/s/mm/ °C.The workpiece was meshed via 20,000 tetrahedron elements.A shear friction coefficient of 0.2 was applied between dies and workpiece.

    5.2.Deformation and microstructure evolution analysis

    The comparison of the experimental and simulated loads during ECAE process is shown in Fig.13.The experimental load increases to the peak, followed by a relatively stable stage and decreases ultimately.The variation trend of sim-ulated load matches well with the experimental result.The experimental load is slightly higher than the predicted curve,which is due to the increasing deformation resistance of alloy with heat dissipation effect.The established FE model based on the unified ISV approach can reflect the plastic deformation feature of AZ80 magnesium alloy during ECAE process.

    Fig.14.Simulated distributions of normalized dislocation density (a), DRX fraction (b) and AGS (c).

    The simulated distributions of microstructural ISVs including normalized dislocation density, DRX fraction and AGS are shown in Fig.14(a)-(c), respectively.The shear band at the channel angular can be observed by the distribution of normalized dislocation density in Fig.14(a).As workpiece entering the shear deformation zone, the proliferation and entanglement of dislocations result in high dislocation density, which provides storage energy for the subsequent DRX [48].When the workpiece flows out of the shear deformation zone, the entangled dislocation network annihilatesby the DRX effect, and the dislocation density significantly decreases.Fig.14(b) and (c) show the similar distribution characteristics of DRX fraction and AGS, respectively.With the shear deformation effect at the channel angular, the DRX fraction increases and AGS decreases.At the outlet of extrusion channel, the distributions of DRX fraction and AGS are different to the distributions at the stable deformation zone,which can attributed to the uneven metal flow with the variation of channel diameter.

    Fig.16.Microstructures in different deformation zones at A1 (a), A2 (b) and A3 (c).

    Fig.15 shows the evolutions of normalized dislocation density, DRX fraction and AGS with the variation of deformation pattern.The distance at horizontal coordinate is defined as the center axis started at the upper surface of channel in Fig.14(c).With the development of ECAE process, the normalized dislocation density reaches a peak and then decreases gradually, which indicates the hot deformation mechanisms including WH and DRX softening.The AGS significantly decreases with the increase of DRX fraction.

    To validate the reliability of FE model in predicting microstructure evolution, three representative locations in the center of deformed workpiece in Fig.14(c) were selected.The microstructures in different deformation zones of A1, A2 and A3 are shown in Fig.16.In Fig.16(a), the workpiece undergoes initial upsetting deformation with large extrusion force before entering the shearing zone, which results in the equiaxed grains with a few DRX grains nucleating at grain boundaries.In Fig.16(b), with the effect of shear deformation at the channel angular, the recrystallization process is developed sufficiently, and the AGS decreases rapidly.As the shearing deformation process finished, the proportion of DRX grains increases, and the microstructure is further refined, as shown in Fig.16(c).The statistical AGS and DRX fraction in different deformation zones were plotted in Fig.15.The experimental microstructure data matches well with the predicted curves.The ECAE process contributes to the refinement of microstructure and the improvement of comprehensive mechanical properties, which has been proved by many literatures [49,50].It can be concluded that the ISV modeling approach can describe the microstructure evolution coupling deformation mechanisms during the ECAE process of AZ80 magnesium alloy.

    Fig.15.The evolutions of normalized dislocation density, DRX fraction and AGS during ECAE process.

    6.Conclusions

    In this study, a unified ISV model for predicting microstructure evolution during hot working process of AZ80 magnesium alloy was established.The macroscopic and microscopic kinetics including viscoplastic flow, dislocation activities, DRX nucleation and boundary migration were introduced.The microstructure evolution features during hot compression and ECAE processes were investigated.The following conclusions can be summarized:

    (1) The characteristic deformation behavior including work hardening and dynamic softening was dominated by the interactive effects among dislocation multiplication, recovery, recrystallization and grain size.The recrystallized microstructure evolution was modeled by incorporating the evolution laws of critical strain, nucleation rate and boundary migration rate.The model parameters were calibrated based on experimental data analysis and GA-based objective optimization.Predicted flow stress, DRX fraction and AGS match well with experimental results.

    (2) The set of constitutive equations was integrated into FE model via the user defined subroutine to simulate hot compression process.The friction at specimenanvil interface leads to the heterogeneous microstructure distribution.Maximum DRX fraction and minimum AGS locate at the center zone of deformed specimen.The experimental microstructure distributions at center,bulging and end zones agree well with the simulated results.The effect of deformation conditions on dislocation density evolution was explained by competitive processes between different microstructural mechanisms.

    (3) Using the calibrated FE modeling approach, the heterogeneous microstructure features during ECAE process were captured.With the shear deformation effect, high dislocation density accumulates at the channel angular.DRX fraction increases and AGS decreases with the development of ECAE process.The unified ISV model can provide scientific understanding and reasonable prediction for the hot working process of magnesium alloy.

    Declaration of Competing Interest

    None.

    Acknowledgments

    Authors acknowledge the funding supported by National Natural Science Foundation of China(No.52175285),Beijing Municipal Natural Science Foundation (No.3182025), National Defense Science and Technology Rapid support Project(No.61409230113), Scientific and Technological Innovation Foundation of Shunde Graduate School, USTB and Fundamental Research Funds for the Central Universities(No.FRFBD-20-08A, FRF-TP-20-009A2).

    叶爱在线成人免费视频播放| 热re99久久国产66热| 国产精品一区二区三区四区久久 | 精品久久久久久,| 两人在一起打扑克的视频| 18禁观看日本| 女生性感内裤真人,穿戴方法视频| 在线av久久热| 老汉色∧v一级毛片| 国产成人精品久久二区二区91| 免费无遮挡裸体视频| av超薄肉色丝袜交足视频| 夜夜躁狠狠躁天天躁| 麻豆成人av在线观看| 中国美女看黄片| 久热爱精品视频在线9| 亚洲成国产人片在线观看| 国产xxxxx性猛交| 91老司机精品| 久久香蕉国产精品| 好看av亚洲va欧美ⅴa在| 一边摸一边抽搐一进一小说| 婷婷精品国产亚洲av在线| 九色国产91popny在线| 俄罗斯特黄特色一大片| 好看av亚洲va欧美ⅴa在| 日韩视频一区二区在线观看| 国产成人欧美| 高清在线国产一区| 久久性视频一级片| 老司机靠b影院| 中文字幕最新亚洲高清| 人成视频在线观看免费观看| 欧美在线黄色| av电影中文网址| 777久久人妻少妇嫩草av网站| 久久久久久国产a免费观看| 9热在线视频观看99| 国产97色在线日韩免费| 午夜福利免费观看在线| 精品人妻1区二区| 咕卡用的链子| 桃色一区二区三区在线观看| 极品教师在线免费播放| 99re在线观看精品视频| 亚洲国产高清在线一区二区三 | 国产精品免费一区二区三区在线| 一边摸一边做爽爽视频免费| 一区二区三区国产精品乱码| 操美女的视频在线观看| 亚洲中文av在线| 国产又色又爽无遮挡免费看| 亚洲国产精品999在线| 一区二区三区激情视频| 免费搜索国产男女视频| 国产成人系列免费观看| 欧美精品啪啪一区二区三区| 色av中文字幕| 99国产综合亚洲精品| 97人妻精品一区二区三区麻豆 | 久久精品国产清高在天天线| 一区二区三区激情视频| 中文字幕av电影在线播放| 又大又爽又粗| 黑人巨大精品欧美一区二区蜜桃| 国产人伦9x9x在线观看| 一级a爱视频在线免费观看| 久久伊人香网站| 亚洲黑人精品在线| 久久中文看片网| 国产精品综合久久久久久久免费 | 亚洲五月天丁香| 一区二区日韩欧美中文字幕| 黄网站色视频无遮挡免费观看| 多毛熟女@视频| 午夜亚洲福利在线播放| 男人舔女人下体高潮全视频| 男女下面进入的视频免费午夜 | 亚洲成av人片免费观看| 99国产精品一区二区蜜桃av| 黄频高清免费视频| 视频区欧美日本亚洲| 十八禁网站免费在线| 国产精品日韩av在线免费观看 | 人妻久久中文字幕网| 国产97色在线日韩免费| 精品久久久久久成人av| av片东京热男人的天堂| 国产欧美日韩一区二区精品| 国产欧美日韩一区二区精品| 久久久久久亚洲精品国产蜜桃av| 国产成人欧美在线观看| 午夜福利视频1000在线观看 | 亚洲av片天天在线观看| 国产亚洲av高清不卡| 大型黄色视频在线免费观看| 色综合欧美亚洲国产小说| 日韩av在线大香蕉| 12—13女人毛片做爰片一| 熟女少妇亚洲综合色aaa.| 老鸭窝网址在线观看| 97人妻精品一区二区三区麻豆 | 国产伦人伦偷精品视频| 欧美日韩黄片免| 91大片在线观看| 我的亚洲天堂| 一级作爱视频免费观看| 午夜福利视频1000在线观看 | 91大片在线观看| 69精品国产乱码久久久| 两性夫妻黄色片| 久久精品成人免费网站| 亚洲国产精品成人综合色| 国产三级黄色录像| 国产精品二区激情视频| 色综合亚洲欧美另类图片| 99国产精品免费福利视频| 欧美老熟妇乱子伦牲交| 精品国产美女av久久久久小说| 国产精品久久久久久精品电影 | 巨乳人妻的诱惑在线观看| av天堂在线播放| 欧美日韩黄片免| 久久人人精品亚洲av| 成人三级黄色视频| av欧美777| 久久亚洲精品不卡| 在线视频色国产色| 日韩精品青青久久久久久| 男女午夜视频在线观看| 亚洲成人久久性| 每晚都被弄得嗷嗷叫到高潮| 麻豆成人av在线观看| 国产免费男女视频| 九色亚洲精品在线播放| 午夜福利影视在线免费观看| 一本综合久久免费| 色综合婷婷激情| 丁香六月欧美| netflix在线观看网站| 欧美乱码精品一区二区三区| 国产成+人综合+亚洲专区| 亚洲午夜精品一区,二区,三区| 成人永久免费在线观看视频| 午夜日韩欧美国产| 欧美国产精品va在线观看不卡| 亚洲成人免费电影在线观看| 大型黄色视频在线免费观看| 午夜福利18| 日本五十路高清| 国产男靠女视频免费网站| 国产男靠女视频免费网站| 最好的美女福利视频网| 涩涩av久久男人的天堂| 亚洲视频免费观看视频| 国产成年人精品一区二区| 成在线人永久免费视频| 国产精品久久视频播放| 亚洲一卡2卡3卡4卡5卡精品中文| tocl精华| 三级毛片av免费| 人妻久久中文字幕网| 制服人妻中文乱码| 中文字幕久久专区| 国产成人av激情在线播放| 欧美乱色亚洲激情| 国产视频一区二区在线看| 亚洲精品中文字幕在线视频| 国产精品影院久久| av视频免费观看在线观看| 99re在线观看精品视频| 国产蜜桃级精品一区二区三区| 宅男免费午夜| 男女之事视频高清在线观看| 高清黄色对白视频在线免费看| 黄片小视频在线播放| 日本精品一区二区三区蜜桃| 琪琪午夜伦伦电影理论片6080| 成人免费观看视频高清| 亚洲性夜色夜夜综合| 欧美激情极品国产一区二区三区| 久久中文看片网| 国产午夜精品久久久久久| 免费高清在线观看日韩| 亚洲人成77777在线视频| 性色av乱码一区二区三区2| av在线天堂中文字幕| 黄色视频,在线免费观看| www.www免费av| 成年女人毛片免费观看观看9| 夜夜爽天天搞| av超薄肉色丝袜交足视频| 又大又爽又粗| 视频在线观看一区二区三区| 久久精品国产清高在天天线| 国产精品香港三级国产av潘金莲| 国产不卡一卡二| 中文字幕人妻熟女乱码| 一区二区日韩欧美中文字幕| 国产蜜桃级精品一区二区三区| 日韩欧美免费精品| 男女床上黄色一级片免费看| 亚洲精品中文字幕一二三四区| 久久天躁狠狠躁夜夜2o2o| 69精品国产乱码久久久| 色av中文字幕| 精品国内亚洲2022精品成人| 中文亚洲av片在线观看爽| 国产精品久久久av美女十八| 亚洲国产高清在线一区二区三 | 亚洲天堂国产精品一区在线| 中文字幕精品免费在线观看视频| 少妇粗大呻吟视频| 亚洲片人在线观看| 国产精品爽爽va在线观看网站 | 99久久久亚洲精品蜜臀av| 久热这里只有精品99| www.999成人在线观看| 日日干狠狠操夜夜爽| 欧美成狂野欧美在线观看| 性欧美人与动物交配| 国产精品 欧美亚洲| 热re99久久国产66热| 久久香蕉精品热| 伊人久久大香线蕉亚洲五| 无遮挡黄片免费观看| 在线免费观看的www视频| 免费女性裸体啪啪无遮挡网站| 久久香蕉激情| 人妻久久中文字幕网| 精品高清国产在线一区| 亚洲欧美激情在线| 国产亚洲精品第一综合不卡| 19禁男女啪啪无遮挡网站| 视频在线观看一区二区三区| 亚洲欧美激情综合另类| 日日夜夜操网爽| 亚洲中文字幕日韩| 成人三级黄色视频| 黑丝袜美女国产一区| 在线观看免费日韩欧美大片| 欧美成人免费av一区二区三区| 国产成人精品久久二区二区免费| 亚洲中文字幕一区二区三区有码在线看 | 久久久国产成人精品二区| 久久伊人香网站| 欧美性长视频在线观看| 婷婷精品国产亚洲av在线| 欧美激情 高清一区二区三区| 中文字幕另类日韩欧美亚洲嫩草| 成人手机av| 婷婷丁香在线五月| 精品国产国语对白av| 亚洲专区中文字幕在线| 岛国视频午夜一区免费看| cao死你这个sao货| 久久久国产成人精品二区| 国产亚洲精品综合一区在线观看 | 久久久久久久久久久久大奶| 很黄的视频免费| 久久国产亚洲av麻豆专区| 国产视频一区二区在线看| 精品国产一区二区三区四区第35| 日韩视频一区二区在线观看| 精品熟女少妇八av免费久了| xxx96com| 不卡一级毛片| 国产欧美日韩一区二区三| 免费搜索国产男女视频| 国产单亲对白刺激| 欧美日本视频| 神马国产精品三级电影在线观看 | 自拍欧美九色日韩亚洲蝌蚪91| 亚洲av成人一区二区三| 此物有八面人人有两片| 亚洲av日韩精品久久久久久密| 波多野结衣一区麻豆| 国产精品一区二区在线不卡| 亚洲午夜理论影院| 麻豆久久精品国产亚洲av| 黄色女人牲交| 久久亚洲真实| 美女扒开内裤让男人捅视频| 国产成+人综合+亚洲专区| 人人妻人人澡欧美一区二区 | 久久人妻av系列| 精品不卡国产一区二区三区| 少妇被粗大的猛进出69影院| 中国美女看黄片| 在线播放国产精品三级| 国产精品乱码一区二三区的特点 | 丝袜人妻中文字幕| 国产高清videossex| 免费无遮挡裸体视频| 亚洲国产看品久久| 午夜福利在线观看吧| 在线观看免费午夜福利视频| 男人的好看免费观看在线视频 | 亚洲人成网站在线播放欧美日韩| 久久亚洲精品不卡| 国产高清有码在线观看视频 | 天天添夜夜摸| 欧美成狂野欧美在线观看| 免费在线观看日本一区| 精品国产超薄肉色丝袜足j| 岛国视频午夜一区免费看| 亚洲精品久久国产高清桃花| 久久国产精品人妻蜜桃| 99久久久亚洲精品蜜臀av| 免费高清在线观看日韩| 首页视频小说图片口味搜索| 亚洲专区国产一区二区| 国产一区二区三区综合在线观看| 热99re8久久精品国产| 又黄又粗又硬又大视频| 久热这里只有精品99| 中亚洲国语对白在线视频| 成熟少妇高潮喷水视频| 欧美精品亚洲一区二区| 国产91精品成人一区二区三区| 一区二区三区精品91| 91大片在线观看| 999久久久精品免费观看国产| 午夜福利,免费看| av视频在线观看入口| 久久久久国产一级毛片高清牌| 麻豆一二三区av精品| 波多野结衣高清无吗| 亚洲熟妇中文字幕五十中出| 午夜福利视频1000在线观看 | 国产不卡一卡二| 久久精品成人免费网站| 正在播放国产对白刺激| 日韩中文字幕欧美一区二区| 午夜免费观看网址| 桃红色精品国产亚洲av| 亚洲国产精品久久男人天堂| 国产精品98久久久久久宅男小说| 一a级毛片在线观看| 麻豆av在线久日| 亚洲精品国产一区二区精华液| 两性午夜刺激爽爽歪歪视频在线观看 | 在线观看一区二区三区| 久久久久久大精品| 99精品久久久久人妻精品| 午夜久久久久精精品| 999久久久精品免费观看国产| 黄片小视频在线播放| a级毛片在线看网站| 黄色丝袜av网址大全| 身体一侧抽搐| 久久久久久大精品| 亚洲成国产人片在线观看| 欧美激情久久久久久爽电影 | 可以免费在线观看a视频的电影网站| 国产伦一二天堂av在线观看| 亚洲少妇的诱惑av| 欧美日韩乱码在线| 国产高清videossex| 亚洲一区二区三区色噜噜| 久久精品人人爽人人爽视色| 69精品国产乱码久久久| 夜夜躁狠狠躁天天躁| 亚洲av成人不卡在线观看播放网| 亚洲精品国产精品久久久不卡| 免费不卡黄色视频| 国产熟女xx| 丁香欧美五月| 色精品久久人妻99蜜桃| 成人18禁在线播放| 国产午夜福利久久久久久| 欧美黄色片欧美黄色片| 可以在线观看的亚洲视频| 午夜日韩欧美国产| 免费少妇av软件| 夜夜爽天天搞| 精品久久久久久久久久免费视频| 国产精品综合久久久久久久免费 | 热99re8久久精品国产| 久久精品国产综合久久久| 国产av在哪里看| 男女下面插进去视频免费观看| 欧美久久黑人一区二区| 黄色 视频免费看| 国产一区二区三区在线臀色熟女| 一卡2卡三卡四卡精品乱码亚洲| 午夜激情av网站| 日韩欧美三级三区| 欧美人与性动交α欧美精品济南到| 国产亚洲精品一区二区www| 欧美黑人精品巨大| 久久午夜亚洲精品久久| 久久 成人 亚洲| 亚洲国产精品sss在线观看| 一级片免费观看大全| 欧美久久黑人一区二区| 欧美亚洲日本最大视频资源| or卡值多少钱| 大码成人一级视频| 激情视频va一区二区三区| 90打野战视频偷拍视频| 老司机午夜福利在线观看视频| 日韩中文字幕欧美一区二区| 久久久水蜜桃国产精品网| 国产成人欧美在线观看| 欧美日韩精品网址| 91麻豆精品激情在线观看国产| 无人区码免费观看不卡| 色老头精品视频在线观看| 欧美黄色片欧美黄色片| 免费在线观看黄色视频的| 久久久久国产精品人妻aⅴ院| 人人妻人人爽人人添夜夜欢视频| 性少妇av在线| 成年版毛片免费区| 欧美黑人欧美精品刺激| e午夜精品久久久久久久| 成人国产综合亚洲| 欧美黑人欧美精品刺激| 国产97色在线日韩免费| 亚洲欧美激情在线| 国产亚洲欧美在线一区二区| www.精华液| 久久人人97超碰香蕉20202| 日韩av在线大香蕉| 91精品国产国语对白视频| 91成人精品电影| 色综合站精品国产| 人妻久久中文字幕网| 亚洲专区中文字幕在线| 成在线人永久免费视频| 99久久综合精品五月天人人| 每晚都被弄得嗷嗷叫到高潮| 国产成人精品在线电影| 欧洲精品卡2卡3卡4卡5卡区| 啦啦啦观看免费观看视频高清 | 夜夜爽天天搞| 久久青草综合色| 亚洲欧美日韩无卡精品| 午夜精品久久久久久毛片777| 搞女人的毛片| 免费在线观看完整版高清| 亚洲一区高清亚洲精品| 亚洲成人免费电影在线观看| 大陆偷拍与自拍| av天堂久久9| 丁香欧美五月| tocl精华| 国产精品99久久99久久久不卡| 视频区欧美日本亚洲| 国内久久婷婷六月综合欲色啪| 欧美日本亚洲视频在线播放| 91国产中文字幕| 一边摸一边抽搐一进一出视频| 精品久久久久久久毛片微露脸| 级片在线观看| 精品人妻1区二区| 91国产中文字幕| 搞女人的毛片| 啦啦啦 在线观看视频| 法律面前人人平等表现在哪些方面| 99国产极品粉嫩在线观看| 亚洲人成电影观看| 757午夜福利合集在线观看| 欧美黑人欧美精品刺激| 午夜福利成人在线免费观看| 亚洲成人久久性| 一级,二级,三级黄色视频| 熟妇人妻久久中文字幕3abv| 国产精品九九99| 校园春色视频在线观看| 欧美+亚洲+日韩+国产| 午夜福利视频1000在线观看 | 久久久久久大精品| 久久午夜综合久久蜜桃| 成人国产一区最新在线观看| 午夜福利视频1000在线观看 | 涩涩av久久男人的天堂| 好男人在线观看高清免费视频 | 欧美日韩亚洲国产一区二区在线观看| 国产高清视频在线播放一区| 这个男人来自地球电影免费观看| 香蕉久久夜色| 老汉色av国产亚洲站长工具| 欧洲精品卡2卡3卡4卡5卡区| 午夜免费成人在线视频| 99香蕉大伊视频| 国产精品免费视频内射| 亚洲精品中文字幕一二三四区| 看黄色毛片网站| 黄片小视频在线播放| 国内久久婷婷六月综合欲色啪| 亚洲国产毛片av蜜桃av| 国产成人av教育| 午夜两性在线视频| 一级,二级,三级黄色视频| 无限看片的www在线观看| 19禁男女啪啪无遮挡网站| 久久午夜综合久久蜜桃| 久久国产乱子伦精品免费另类| 国产av又大| 亚洲专区国产一区二区| 久久人妻熟女aⅴ| 亚洲av成人一区二区三| 亚洲熟妇熟女久久| 午夜成年电影在线免费观看| 18禁观看日本| 久久人妻熟女aⅴ| 国产成人欧美在线观看| 香蕉久久夜色| 久久久久久免费高清国产稀缺| 999精品在线视频| 日韩欧美一区视频在线观看| 99精品欧美一区二区三区四区| 日韩欧美免费精品| 日韩欧美国产在线观看| 波多野结衣av一区二区av| 热99re8久久精品国产| 亚洲欧美激情在线| 亚洲欧美一区二区三区黑人| 亚洲精品av麻豆狂野| 十八禁网站免费在线| svipshipincom国产片| 天天躁狠狠躁夜夜躁狠狠躁| 可以在线观看的亚洲视频| 欧美亚洲日本最大视频资源| 中国美女看黄片| 手机成人av网站| 日韩av在线大香蕉| 欧美乱色亚洲激情| 国产精品爽爽va在线观看网站 | 高清毛片免费观看视频网站| 亚洲欧美激情在线| 亚洲国产欧美一区二区综合| 超碰成人久久| 免费女性裸体啪啪无遮挡网站| 看黄色毛片网站| 美女高潮到喷水免费观看| 看免费av毛片| 亚洲,欧美精品.| 高清在线国产一区| 国产主播在线观看一区二区| 可以免费在线观看a视频的电影网站| 韩国精品一区二区三区| 亚洲人成伊人成综合网2020| 亚洲 国产 在线| 精品国产一区二区久久| 一区二区三区高清视频在线| 波多野结衣一区麻豆| 亚洲精品在线美女| 久久午夜亚洲精品久久| 激情视频va一区二区三区| 18美女黄网站色大片免费观看| 亚洲色图av天堂| 人成视频在线观看免费观看| 欧美国产精品va在线观看不卡| 亚洲一码二码三码区别大吗| 91九色精品人成在线观看| 女性生殖器流出的白浆| 一区二区三区国产精品乱码| 国产精品日韩av在线免费观看 | 男女床上黄色一级片免费看| 日韩欧美国产在线观看| 国产精品免费一区二区三区在线| 国产精品 欧美亚洲| 美女免费视频网站| 99久久综合精品五月天人人| 在线视频色国产色| 夜夜夜夜夜久久久久| 久久欧美精品欧美久久欧美| 人人妻人人澡人人看| 12—13女人毛片做爰片一| 99精品久久久久人妻精品| 中文亚洲av片在线观看爽| www.www免费av| 亚洲精品在线美女| 给我免费播放毛片高清在线观看| 亚洲国产欧美日韩在线播放| 久久久久九九精品影院| 夜夜看夜夜爽夜夜摸| 午夜福利视频1000在线观看 | a级毛片在线看网站| 亚洲欧美激情在线| 亚洲五月婷婷丁香| 国产99白浆流出| 亚洲人成电影免费在线| 天天躁夜夜躁狠狠躁躁| 日韩有码中文字幕| 日韩欧美国产一区二区入口| av片东京热男人的天堂| 亚洲精品一区av在线观看| 亚洲欧美激情综合另类| 人人妻人人澡欧美一区二区 | 午夜福利视频1000在线观看 | 精品国产一区二区久久| 无遮挡黄片免费观看| 婷婷丁香在线五月| 18禁国产床啪视频网站| 亚洲欧洲精品一区二区精品久久久| 18禁国产床啪视频网站| 一级毛片女人18水好多| АⅤ资源中文在线天堂| 精品国内亚洲2022精品成人| 18禁国产床啪视频网站| 精品一区二区三区av网在线观看| 中文字幕人妻丝袜一区二区| 亚洲色图综合在线观看| 黄色视频,在线免费观看| 日韩国内少妇激情av| 亚洲av熟女| 亚洲熟女毛片儿| 欧美色视频一区免费| 日韩欧美国产一区二区入口| 国产成人欧美在线观看| av在线播放免费不卡|