999精品在线视频,手机成人午夜在线视频,久久不卡国产精品无码,中日无码在线观看,成人av手机在线观看,日韩精品亚洲一区中文字幕,亚洲av无码人妻,四虎国产在线观看 ?

Nonlinear Problems via a Convergence Accelerated Decomposition Method of Adomian

2021-04-27 10:28:46MustafaTurkyilmazoglu

Mustafa Turkyilmazoglu

1Department of Mathematics,Hacettepe University,Ankara,06532,Turkey

2Department of Medical Research,China Medical University Hospital,China Medical University,Taichung,40447,Taiwan

ABSTRACT The present paper is devoted to the convergence control and accelerating the traditional Decomposition Method of Adomian (ADM).By means of perturbing the initial or early terms of the Adomian iterates by adding a parameterized term,containing an embedded parameter,new modified ADM is constructed.The optimal value of this parameter is later determined via squared residual minimizing the error.The failure of the classical ADM is also prevented by a suitable value of the embedded parameter,particularly beneficial for the Duan–Rach modification of the ADM incorporating all the boundaries into the formulation.With the presented squared residual error analysis,there is noneed to check out the results against the numerical ones,as usually has to be done in the traditional ADM studies to convince the readers that the results are indeed converged to the realistic solutions.Physical examples selected from the recent application of ADM demonstrate the validity,accuracy and power of the presented novel approach in this paper.Hence,the highly nonlinear equations arising from engineering applications can be safely treated by the outlined method for which the classical ADM may fail or be slow to converge.

KEYWORDS Nonlinear equations;Adomian decomposition method;modification;convergence acceleration

1 Introduction

Researchers prefer an easily accessible and user friendly method requiring less computational labor while accurately approximating highly nonlinear equations resulting from mathematical modeling of real-life phenomena.The Adomian decomposition method (ADM) is one such popular technique capable of dealing with the prevailing nonlinearities by means of Adomian polynomials [1,2].A modification of the classical ADM is proposed within the current study based on the recent publications [3,4] successfully generating fast convergent ADM series solutions with as small Adomian polynomials as possible in the solution series.

A quick literature survey exhibits that ADM has been applied to many nonlinear equations [5].To classify some of the recent bibliography,algebraic equations were contained within the references [6,7].The ordinary differential equations were dealt within the citations [8–14].The articles [15–18] covered the efforts to partial differential equations.Mathematical analysis of the convergence of ADM to certain nonlinear equations was fulfilled in the publications [19–21].It can be successfully used to gain correct physical parameters domain [22].A traffic model was also very recently treated in [23] via the Adomian method.The publications by [24,25] present investigation of some nonlinear problems via different numerical approaches.

It is now well-known that an inadequate arrangement of the classical ADM series may lead to non-convergent solutions or solutions with a poor convergence rate.To avoid these shortcomings,a parameter is generally inserted at the leading term of the Adomian series and later it is subtracted at the first order term not to break down the equation structure.This procedure was pursued by the recent publications [11,17,18].However,how a proper value of the inserted parameter will be determined was not mentioned in these references.Instead,a randomly chosen value was assigned to it.A variety of modifications were also offered in the articles [26–31].A successful formulation of the ADM was made in the recent work of [3] which was named as the optimal ADM.Further applications of the homotopy analytic approximate method may be found in the literature [32–35].

The motivation of the current work is,benefiting from the idea in [11],to devise a method that greatly improves the mathematical property of classical ADM.Within this aim,a reorganization of the ADM series is proposed by altering the early terms so that they incorporate extra controllable terms.The reason of such a treatment is to get a rapidly converging ADM solutions with the least Adomian polynomials.In place of randomly selecting,an optimum value of the introduced parameter is later determined through error on the grounds of total residual.With this value at our disposal,there is no doubt that the ADM method is convergent to the true solution in a most rapid way,not demanding a verification of the ADM solutions by numerical ones.The failure of the classical ADM in the usual form or in the Duan–Rach formulation is also prevented by a suitable value of the embedded parameter.The present approach can also extend the region of convergence of the traditional method.Examples of physical value are provided to justify and validate the given procedure.

2 Traditional Decomposition Method of Adomian

The usual steps of traditional ADM can be inferred from the aforementioned citations.The methodology in brief is such that under an invertible linear operatorLand a forcing functionf,it is desired to approximate the functionuhaving the nonlinearityN(u)and satisfying the general nonlinear equation

with the initial and/or boundary restrictions

Having inverted (1) under the restrictions (2) generally leads to

wheregis due to the conditions in (2).Ifuis a single scalar parameter like (1) representing an algebraic equation,then there is no suchgin (3),whereas,in the case of a variableu,L?1denotes an integral operator giving rise togin (3).Then (3) is a mixed Volterra-Fredholm type equation so-called as the Duan–Rach formulation in the recent literature,see for instance [13,14].

The subsequent series decompositions ofuandN(u)

in whichAn’s are the classical Adomian polynomials,are later substituted into (3).The solutionuof (1) is finally generated from the recurrence relation

As a result,by means of the relations from (5),an approximate series solution of orderMis obtained as

which serves for practical purposes.

In general,the procedure in (5) yields convergent ADM series solutions,see for instance [19–21].If not,to achieve convergent solutions or for computational conveniences some modifications in the termsuiin (5) are implemented as in the articles [11,12],without a proper mathematical evidence and support.

3 A Modified Decomposition Method of Adomian

To overcome the divergence of classical ADM or to speed up the convergence rate of the ADM series,the leading order termu0in (5) (which is in compliance with the previous implementations,in for instance [11]) or some of the early terms,callue,0 ≤e≤Me,in (5) will be modified so that the modification will account for the change in the ADM by incorporating parameterized terms(h)withhan embedded unknown parameter [3].To exemplify,consider the modified version of (5) in the manner

The following conditions for parameterized terms in the new algorithm (7) should be added

so that it can be reduced to the traditional ADM forh=0.

It is remarked that there is no a unique way of selecting the(h)terms,but they must be as simple as possible in line with the simplicity inspiration of the ADM method itself.For instance,constants depending onhor simple power functions (or simple integrable functions)in combination withhcan be used,refer to the illustrations belove for more details.We should remark that in the absence of tilde terms (orh) (7) dully conforms with (5).Having inserted such tilde terms in the modified ADM (7),the question now arises,how to determine the correct or proper value ofh? So far,in ADM applications of this kind,only random or trial values are preassigned,see for instance [11] amongst others.On the other hand,we may either observe the least change in the plot of some physical quantities,the method is so-called as the constanth-level curves [36] giving rise to rough estimates forh,or we outline the following rigorous algorithm to determine the best or optimalh.

Algorithm.Consider the squared residual error corresponding to (1) defined by

where eitherL2(Ω)orL(Ω)norms are employed.Moreover,uin (8) is owing to the modified ADM from (7).The optimum value ofhis such that it minimizes the squared residual error (8)at the approximation levelM.

As a consequence,the above Algorithm will generate the best value ofhwhich will ensure the convergence of ADM series solution (7) in a fastest rate of convergence.The minimization task of(8) may be fulfilled by means of contemporary softwares,such as MAPPLE or MATHEMATICA.

4 Applications

Potential applications of the introduced ADM in (7) are given here.To control the error,we use the norm

with the exactueand ADM solutionu.

4.1 An Algebraic Equation

As stated by Adomian [1] the classical ADM method (5) fails to result in a convergent solution of

for the solutionu=?0.73205080757.On the other hand,when the new ADM is built via

Fig.1 displaysh-level curves at selected truncation ordersM.The intervalh∈(?0.35,?0.25)is observed to yield convergent ADM solutions,which excludes the traditional ADM withh=0 [1].

Through the residual minimization at the approximation levelM=8,h=?0.2679492 is obtained as the optimum.The history and why this value is the best for the convergence control,as compared to the failure of classical ADM can be visualized in Tab.1.

Figure 1:Convergence control parameter h regarding (11)

Table 1:Convergence of modified ADM (12) for Eq.(11) with different h.Paranthesis is for the absolute error (10)

The convergent solution of (11) with the new modified ADM (12) at the approximation levelM=8 is found to be

Steve went home from school, thoughtful, that afternoon. Walking into the house, he took one look around. Both parents were passed out, in various stages of undress, and the stench was overpowering! He, quickly, gathered up his camping gear, a jar of peanut butter, a loaf of bread, a bottle of water, and this time...his schoolbooks. Grim faced and determined7, he headed for the woods.

for which the optimumhis tabulated in Tab.1.

4.2 Equation Involving Integral

Consider the equation given in [30]

The modified ADM method (7) for the current integral problem is adopted as

We findh∈(?1500,1500)for the convergence interval in Fig.2,refer also to [3].A comparison and convergence accelerating feature of the present ADM can be visualized from the Tab.2.It is noticed from Tab.2 that even the 4th-order modified ADM (16) is able to produce the solution accurate of order 10?10,which is adequate in practical purposes.

Figure 2:Convergence control parameter h regarding (15)

Table 2:Errors in (15) regarding (14).Parenthesis denotes the optimum values of h

4.3 A Fin with Porosity Feature

As taken from [14],a porous fin can be modelled via

In (17),temperature along the fin isu,andsandβare physical parameters.

The modified ADM algorithm (7) here is

At the selected valuess=5 andβ=1,the history of convergence is tabulated in Tab.3.It is seen how useful the modified ADM over the classical one by more than doubling the convergence rate at the same number of iterations.The CPU times evaluated by MATHEMATICA through the residual are shown in Tab.3.The list proves the improved efficiency of the current method over the classical ADM.

Table 3:Residual errors of classical ADM and modified ADM (18) for Eq.(17) with s=5 and β=1.Parenthesis denote the optimum values of h

Table 4:Numerical and modified ADM results regarding (17) [3]

Instead of the modification of ADM in (18),we may use the Duan–Rach formulation involving no unknown parameters within it except the embedded parameterh

Choosings=β=1,Fig.3 showsh?level curves,indicating that the prescription ofhmust be made within the range(0.39,0.49)to ensure the convergence of the modified ADM 19.The exact values ofu(0)andu′(1)are respectively,0.522738093570 and 1.13937891581.

Figure 3:Convergence control parameter h regarding (17)

To demonstrate the power of the modified ADM (19),Tab.5 shows the squared residual error(9) from both the novel and classical ADM.It is unfortunate to observe that the Duan–Rach formulation (19) withh=0 fails to converge,however,the optimum embedded parameterhinsures that the modified ADM is convergent for the present physical problem,even if the convergence is not as fast as the modified formulation in (18).

Table 5:Values offor (17) from the classical and modified ADM methods

Table 5:Values offor (17) from the classical and modified ADM methods

M √Res(0) √Res(h) h 4 0.91×101 9.379592752×10?3 0.5717083134 8 5.90×102 8.485113566×10?5 0.4466738418 12 1.90×105 1.712138085×10?5 0.4205569783 16 9.50×107 8.888688579×10?6 0.4108495377 20 5.68×1010 3.872408211×10?6 0.4061899784

4.4 Gelfand Equation

The Gelfand equation [5] involves exponential nonlinearity [8]

withu′(0)=?0.4636325917 [3].

In line with the publications [5,8] whenh=0,the modified ADM is

u0(η)=0,

whereAn(h,η)are the polynomials of Adomian foreu(η),encompassing all boundary conditions in (20) within the Duan–Rach approach,and we have

Fig.4 shows the predicted convergence control parameters.WithM=12,an optimum value for the embedding parameterhis found to be?0.01274 from the Algorithm in (9).We find that the residual error is 1.6835117258×10?15withu′(0)=?0.4636325917 from the present approach.On the other hand,and error of 4.813740934×10?9withu′(0)=?0.4636325899 is calculated from the traditional ADM.Our value is given by

The success of the present modified ADM (21) is thus obvious.

Figure 4:Convergence control parameter h regarding (20)

4.5 Electrostatic Cantilever Micro-Electromechanical System

The beam-type electrostatic actuators for the nonlinear cantilever micro-electro mechanical systems are modelled by the fourth-order boundary value problem from [11]

To comply with the Duan–Rach Adomian decomposition method in [11],the present modified ADM is

where the Adomian polynomialsAn(h,η)are due to the negative-power nonlinearities in (23),see[11] for more details.

For the fixed parametersK=3,αK=0.2,β=0.5 andγ=0.25,Fig.5 produced from (24)helps us guess the proper values of embedded parameterh.Small values ofhare seen to be adequate to get convergent ADM solutions.

Figure 5:Convergence control parameter h regarding (23)

Fig.6 demonstrates different approximation levelsM,and it signifies toh=0.1045421730 as the optimumhwhenM=10.With this optimum value of the embedding parameter,the squared residual error for the current problem is=3.2727623195×10?9,whereas the classical residual error withh=0 is=8.7520447833×10?5.This implies that much more Adomian series terms are required for the classical ADM to reach the accuracy of the modified ADM here.Moreover,from [11] the physical values areu′′(0)=?0.56764 andu′′′(0)=1.07585 evaluated with 12 Adomian polynomials.The present values areu′′(0)=?0.5676684138 andu′′′(0)=1.0758898885 evaluated with only 10 Adomian polynomials,correct to 9 decimal places.Hence,the advantage of the present modification is clear for the present highly nonlinear physical problem.We should remark that with the present Algorithm,the validation of the ADM results against the numerical ones is no longer a prerequisite.Also,there is no need to evaluate error remainder functions,nor to consult to Pade-approximates to increase the accuracy,both of which as implemented in [11].

Figure 6:Error regarding (23) for various M

The convergence accelerating feature of the present modified ADM (24) as compared to the classical ADM is better visualized from the Tab.6.Table also shows the comparable CPU times.

Table 6:Convergence history of modified ADM vs.classical ADM for Eq.(23)

4.6 Electrostatic Cantilever Nano-Electromechanical System

Nonlinear model for the electrostatic double cantilever nano-electromechanical system in the case of Casimir force (K=4) is given by [11]

We adopt the subsequent modified ADM,that conforms to the classical ADM (h=0) given in [15]

where the Adomian polynomialsAn(h,η)are due to the negative-power nonlinearities in (25),see[11] for more details.

For the specific parametersαK=1,β=1.5 andγ=0.5,Fig.7 displays the constanth-level curves drawn at the approximation levelM=8.Similar to the previous example,very small values ofhare seen to be adequate to get convergent ADM.

In order to evaluate the performance of modified ADM over the classical one,Tab.7 shows the unknown physical quantitiesu′′(0)=u′′(1)andu′′′(1)at several truncation ordersM.The faster convergence rate of the present ADM is apparent.

To illustrate,the analytical formula computed via the present algorithm (26) atM=4 for the value ofu′′(0)is given by

which is of almost nine degree of accuracy as seen from Tab.7.

Table 7:Values of u′′(0) and u′′′(1) for (25) from the classical and modified ADM methods

4.7 Lane–Emden Equation

We consider the Lane–Emden type boundary value problem from [9]

that models the oxygen diffusion in a spherical cell with Michaelis–Menten uptake kinetics.We take into account the subsequent constants to comply with the literature [9]

k=4,r=0.76129,s=0.03119.

The modified ADM that is offered for the present problem is then

which conforms with the classical ADM of [9] in the limith→0 and the Adomian polynomialsAn(h,η)are to account for the nonlinearity.The interest is to determine the physical value ofu(0)=a.

We present Tab.8 to demonstrate the performance of the modified ADM (29) versus the classical ADM.The expected practical accuracy is met at lower Adomian series approximations via the modified method.

Table 8:Values of u(0)=a and for (28) from the classical and modified ADM methods

Table 8:Values of u(0)=a and for (28) from the classical and modified ADM methods

M √Res(0) √Res(h) a (h=0) a h 2 5.05920434×10?5 1.82258949×10?5 0.896870943665 0.896877048532?0.0289267765 4 2.24898539×10?7 2.18219055×10?8 0.896877026028 0.896877046709?0.0321993310 6 1.12229262×10?9 2.83931728×10?11 0.896877046624 0.896877046709?0.0335012655

4.8 The Fluid Flow of Jeffery–Hamel

The Jeffery–Hamel fluid flow problem is modelled via [22]

withα,ReandHaare the physical parameters [13].

Following the successful Duan–Rach ADM formulation of the problem (30) in [17],we propose the following modified version

whereAn(h,η)=?(4?Ha)α2u′n(η)?are the Adomian polynomials.

For the diverging channel,considering the specific parametersα=50,Re=50 andHa=1000 to be in line with [13],the constanth-level curves are depicted in Fig.8 at the approximation orderM=8.Accordingly,we expect the convergency of the modified ADM (31) in the vicinity of zero.

The performance of modified ADM (31) is next measured by computing the centerline velocityu(0.5)(numerical value is 0.764064240111) at different approximation levelsMas shown in Tab.9.It is observed that 10 digits of accuracy is quickly reached by the present ADM,whereas the classical ADM falls behind.Hence,even though it was not clearly mentioned in [13] (see Tab.1 therein),the accuracy of order 10?8as obtained via the classical ADM demands at least 15–20 Adomian polynomials,whereas only 6 Adomian polynomials are sufficient to gain the same accuracy with the present modification.

Figure 8:Convergence control parameter h regarding (30)

Table 9:Values of u(0.5) and for (30) from the classical and modified ADM methods

Table 9:Values of u(0.5) and for (30) from the classical and modified ADM methods

M √Res(0) √Res(h) u(0.5) (h=0) u(0.5) h 2 6.92050730×10?2 1.13520536×10?2 0.763699786377 0.764097626541?0.0180187189 4 1.94819917×10?3 1.87302552×10?5 0.764051346961 0.764064279681?0.0179816107 6 7.53651632×10?5 3.95453415×10?8 0.764063696627 0.764064239975?0.0178807933 8 3.34939742×10?6 9.72265613×10?11 0.764064214850 0.764064240111?0.0178880579

4.9 Squeezing Two Parallel Plates

The flow squeezed between two parallel plates are modelled by the nonlinear equations [22]

see [14] for the flow parameters.

In accordance with the Duan–Rach ADM formulation of the physical problem (32) in [14],we set the modified ADM in the form

where

are the Adomian polynomials.

To make a comparison with the classical ADM in [14],we set the parametersS=1,φ=0.02 and concentrate on Cu-Water nanofluid withρs=8933 andρf=997.1.The corresponding constanth-level curves are plotted in Fig.9.It appears that the convergence of the modified ADM(32) is guaranteed ifhis selected in the overlapping interval(?1,1).

Figure 9:Convergence control parameter h regarding (32)

The effects of iterative numberMon the skin frictionu′′(1)are next demonstrated in Tab.10.The advantage of the modified ADM (33) with optimum values of embedding parameterhis clearly observed against the classical ADM withh=0.It appears that the accuracy of classical ADM of O(10?8)as displayed in Tab.2 of [18] clearly requires more Adomian series terms (seeM=8 in Tab.10) as compared to the less terms needed in the modified ADM here.

Table 10:Values of u′′(1) and for (32) from the classical and modified ADM methods

Table 10:Values of u′′(1) and for (32) from the classical and modified ADM methods

M √Res(0) √Res(h) u′′(1) (h=0) u′′(1) h 2 2.67740258×10?1 6.68835495×10?2?3.6764267174?3.6917310584?0.0991340307 4 1.61057990×10?2 7.92481375×10?4?3.6952425904?3.6965216730?0.1097576248 6 1.09939458×10?3 6.61056361×10?6?3.6964553942?3.6965506329?0.1111585030 8 8.10446170×10?5 7.71479681×10?8?3.6965431988?3.6965505168?0.1112323962 10 6.28790630×10?6 1.09803715×10?9?3.6965499276?3.6965505117?0.1112270087

The following fourth-order modified ADM series solution for the skin friction may serve good to the purpose of engineering applications if not high accuracy is required

4.10 Nonlinear Oscillator Problem

Let us consider the nonlinear oscillator Duffing problem (see [36] (Chapter 5) and [4])

which involves a cubic nonlinearity.

The improved ADM can be given via

with the Adomian termsAn(t)in (35).

The classical Adomian method withh=0 in (36) is not convergent,whereas withh=0.68981924,the residual error becomes=0.02589439,and the convergency is satisfied as revealed in Fig.10 for the domain of definition.

Figure 10:Convergence control parameter h regarding (35)

4.11 Diffusion Equation

Let us consider the nonlinear diffusion equation,see [18] and [3]

for which [18] presents an exact solution.

The form of modified ADM for the partial differential equation (37) is

where

An(h,x,t)

are the Adomian polynomials corresponding to the function?u?u3+(u2ux)x.

Fig.11 shows the constanth-level curves at the approximation levelM=10,indicating a very large range of embedding parameterh.

Figure 11:Convergence control parameter h regarding (37)

Actually at this truncation of the modified ADM series,it is obtained

Defining the squared residual error for (37) as

Tab.11 tabulates how the modified ADM has smaller residual errors.

Table 11:The residual errors for (37) from the classical and modified ADM methods

Table 11:The residual errors for (37) from the classical and modified ADM methods

M √Res(0) √Res(h) h 2 3.067520559×10?1 1.126425884×10?2 0.78155339805 4 1.905329025×10?2 2.920732957×10?4 0.84725848564 6 5.284431922×10?4 4.431860700×10?6 0.88282078470 8 8.251967252×10?6 4.360598101×10?8 0.90500266099 10 8.249541273×10?8 2.997208214×10?10 0.92013994524

4.12 Burger’s Equation

The final example is known as Burger’s equation [3]

with an exact solutionu(x,t)=.

Figure 12:Convergence control parameter h regarding (41)

The form of modified ADM for the partial differential equation (40) is

where

An(h,x,t)are the Adomian polynomials corresponding to the functionuxx?uux.

The traditional ADM withh=0 in the domain 0 ≤x≤1,0 ≤t≤1 is noticed to fail to produce a convergent series as also clear in Fig.12,becauseh=0 is not lying in the overlapping region.

With the definition

Tab.12 justifies the success of the present modified ADM over the classical divergent one,both in terms of accuracy and computational cost.

Table 12:The residual errorsfor (41) from the classical and modified ADM methods

Table 12:The residual errorsfor (41) from the classical and modified ADM methods

M √Res(0) CPU √Res(h) (hopt) CPU 5 1.75×102 0.012 4.674547981×10?2 (2.2068296567) 0.039 10 1.67×105 0.053 4.082713247×10?3 (2.3147042279) 0.078 15 1.39×108 0.101 6.755704561×10?4 (2.3159080855) 0.145 20 1.24×1011 0.242 2.741368603×10?5 (2.3820449852) 0.289

5 Concluding Remarks

The aim of the present work is to present superiority over the well-known Adomian decomposition method (ADM) often employed in the recent literature to analytically approximate solutions to highly nonlinear algebraic and differential equations of some real physical motions.Within this aim,a reformulation of the ADM is targeted to prevent first the failure and then convergence acceleration of the classical Adomian polynomials.

To accomplish the objective,the classical ADM is modified by inserting some simple parameterized terms into the early iterates involving an embedded parameter to control and pacing the convergence of the generated ADM series.In order to determine the best suitable value or the optimum value of this parameter,squared residual minimizing of the governing equation is proposed.This enables us to overcome the divergence of the classical ADM,and more importantly,there is no need to check out the results against the numerical ones,as usually has to be done in traditional ADM studies,since the optimum embedded parameter obtained is an insurance for ADM series convergence in a most rapid manner.

Physical examples selected from the recent application of ADM demonstrate the validity,accuracy and power of the present approach in terms of generating the convergent solution within the least number of iterations.In particular,the Duan-Rach modification of the ADM incorporating all the boundaries mostly used in the recent ADM applications takes great benefit of the present proposal,otherwise there is always the inevitable danger that it may lead to non physical solutions.The present approach successfully extends the convergence interval of the studied problem.In conclusion,the present formulation of ADM offers a promising tool to treat more strongly nonlinear equations/systems of real life phenomena.

Funding Statement:The author received no specific funding for this study.

Conflicts of Interest:The author declares that he has no conflicts of interest to report regarding the present study.

主站蜘蛛池模板: 亚洲欧美不卡中文字幕| 中文字幕1区2区| 久久久久久久久18禁秘| 操美女免费网站| 91九色国产在线| 久久黄色视频影| 18黑白丝水手服自慰喷水网站| 中文毛片无遮挡播放免费| 亚洲日韩第九十九页| 国产欧美日韩综合在线第一| 精品视频免费在线| 97se亚洲综合在线天天| 成人av手机在线观看| 亚洲欧洲天堂色AV| 最近最新中文字幕在线第一页| 毛片免费视频| 精品国产电影久久九九| 亚洲AV无码不卡无码| 久久亚洲高清国产| 亚洲精品va| 狠狠ⅴ日韩v欧美v天堂| 国产成人高清在线精品| 高清无码一本到东京热| 华人在线亚洲欧美精品| 亚洲品质国产精品无码| 久久精品女人天堂aaa| 中文字幕在线观| 色综合狠狠操| 亚洲日韩国产精品综合在线观看| 午夜国产大片免费观看| 亚洲天堂区| 伦伦影院精品一区| 日本午夜影院| 老司机午夜精品视频你懂的| 亚洲香蕉久久| 婷婷六月激情综合一区| 国产成人无码Av在线播放无广告| 国产簧片免费在线播放| 影音先锋丝袜制服| 黄色网址手机国内免费在线观看| 国产视频入口| 中文字幕亚洲综久久2021| 网友自拍视频精品区| 久久免费成人| 国产精品香蕉在线| 亚洲中文精品久久久久久不卡| 99伊人精品| 欧美亚洲国产精品久久蜜芽| 国产幂在线无码精品| 国产又爽又黄无遮挡免费观看| 在线免费观看a视频| 五月丁香在线视频| 国产欧美网站| 麻豆精品久久久久久久99蜜桃| 日韩第九页| 国产精品30p| 午夜国产精品视频| 亚洲中文字幕在线观看| 成人亚洲国产| 亚洲v日韩v欧美在线观看| 国产噜噜噜| 四虎永久免费网站| 99这里只有精品在线| 2021国产v亚洲v天堂无码| 亚洲永久色| 久久夜夜视频| 波多野结衣无码视频在线观看| 97色伦色在线综合视频| 毛片网站免费在线观看| 国产午夜无码片在线观看网站 | 亚洲天堂免费在线视频| 人妻中文字幕无码久久一区| 亚洲天堂视频网站| 亚洲美女一区二区三区| 成人av专区精品无码国产| 不卡无码网| 美女扒开下面流白浆在线试听| 日韩一区二区三免费高清| 日韩在线第三页| 网友自拍视频精品区| 久久久久久午夜精品| 国产美女无遮挡免费视频网站 |