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

    具有恐懼效應(yīng)及修正的Holling-Ⅱ捕食者-食餌模型動(dòng)力學(xué)分析

    2024-01-01 00:00:00劉宇鵬石垚
    關(guān)鍵詞:捕食者穩(wěn)定性

    摘要: 利用微分方程的特征值理論、 Poincare-Bendixson環(huán)域定理和Hopf分支理論分析具有恐懼效應(yīng)及修正的Holling-Ⅱ捕食者-食餌模型, 給出該模型平衡點(diǎn)的穩(wěn)定性, 并證明該模型具有穩(wěn)定的極限環(huán)以及在共存平衡點(diǎn)處會(huì)出現(xiàn)Hopf分支. 結(jié)果表明, 恐懼效應(yīng)和修正的Holling-Ⅱ函數(shù)對(duì)系統(tǒng)穩(wěn)定性有顯著影響.

    關(guān)鍵詞: 捕食者-食餌模型; 恐懼效應(yīng); 修正的Holling-Ⅱ功能反應(yīng)函數(shù); 穩(wěn)定性; Hopf分支

    中圖分類(lèi)號(hào): O175.26" 文獻(xiàn)標(biāo)志碼: A" 文章編號(hào): 1671-5489(2024)04-0800-09

    Dynamic Analysis of a Predator-Prey Model ofHolling-Ⅱ with Fear Effect and Modification

    LIU Yupeng, SHI Yao

    (School of Science, Chang’an University, Xi’an 710064, China)

    Abstract: By using" the eigenvalue theory of differential equations, Poincare-Bendixson ring theorem and Hopf bifurcation theory,

    we" analyzed the predator-prey model of Holling-Ⅱ with fear effect and modification, gave the stability of the equilibrium point of the model, and proved that the model had stable

    limit cycles and Hopf bifurcations appeared at coexistence equilibrium points. The results show that the fear effect and the modified Holling-Ⅱ function have si

    gnificant effects on the stability of the system.

    Keywords: predator-prey model; fear effect; modified Holling-Ⅱ functional response function; stability; Hopf bifurcation

    0 引 言

    目前, 對(duì)捕食者-食餌模型的相關(guān)研究及改進(jìn)備受關(guān)注[1-4]. Holling[5]在大量實(shí)驗(yàn)和分析的基礎(chǔ)上, 提出了3種不同類(lèi)型的功能反應(yīng)函數(shù): Holling-Ⅰ,Holling-Ⅱ,Holling-Ⅲ, 且這些功能反應(yīng)函數(shù)只依賴于食餌的種群密度. 文獻(xiàn)[6-10]提出了其他類(lèi)型的功能反應(yīng)函數(shù). Dalziel等[11]在研究可變搜索率的捕食者-食餌模型時(shí), 提出了修正的Holling-Ⅱ功能反應(yīng)函數(shù), 研究結(jié)果表明, 與經(jīng)典Holling-Ⅱ模型相比, 該模型不總出現(xiàn)富集悖論, 即使出現(xiàn)富集悖論, 捕食者也能通過(guò)降低搜索速度進(jìn)行調(diào)整, 從而使系統(tǒng)穩(wěn)定.

    在一些生態(tài)系統(tǒng)中, 食餌可能會(huì)對(duì)捕食者感到恐懼, 從而使捕食者的捕獵更困難. Zanette等[12]在整個(gè)繁殖季節(jié), 利用電籬笆對(duì)歌雀進(jìn)行了田間實(shí)驗(yàn), 結(jié)果表明, 歌雀在感知到捕食風(fēng)險(xiǎn)后, 其繁殖數(shù)量下降40%. 文獻(xiàn)[13-14]對(duì)其他鳥(niǎo)類(lèi)和脊椎動(dòng)物進(jìn)行了類(lèi)似實(shí)驗(yàn), 也得出了同樣的結(jié)論: 即使捕食者和食餌之間沒(méi)有直接捕殺, 但捕食者的存在會(huì)由于反捕食者行為而導(dǎo)致食餌數(shù)量減少. Wang等[15]首次提出了恐懼因子, 并將恐懼因子分別與線性功能反應(yīng)、 Holling-Ⅱ功能反應(yīng)結(jié)合, 建立了捕食者-食餌相互作用中的恐懼效應(yīng)模型, 通過(guò)數(shù)學(xué)分析, 得出無(wú)論是高水平, 還是低水平的恐懼效應(yīng), 都可以使振蕩的系統(tǒng)穩(wěn)定. 此外, Pal等[16\|17]分別研究了恐懼對(duì)帶有狩獵合作的捕食者-食餌模型和恐懼對(duì)帶有狩獵合作的Leslie-Gower模型. 文獻(xiàn)[18-19]分別將具有加法的Allee效應(yīng)、 具有乘法的Allee效應(yīng)和恐懼效應(yīng)結(jié)合, 建立了捕食模型, 并研究了其動(dòng)力學(xué)性質(zhì).

    本文提出將修正的Holling-Ⅱ功能反應(yīng)函數(shù)[11]和Wang等[15]提出的恐懼效應(yīng)因子引入捕食者-食餌模型, 建立如下具有恐懼效應(yīng)及修正的Holling-Ⅱ捕食者-食餌模型:

    dudt=ru1+kv-du-au2-bu2vbHu2+u+g=P(u,v),

    dvdt=cbu2vbHu2+u+g-mv=Q(u,v),(1)

    其中u表示食餌的數(shù)量, v表示捕食者的數(shù)量, r表示食餌的內(nèi)稟增長(zhǎng)率, d和m分別表示食餌和捕食者的死亡率, 參數(shù)a表示食餌在種群內(nèi)部直接的競(jìng)爭(zhēng)強(qiáng)度, 參數(shù)k用來(lái)刻畫(huà)食餌見(jiàn)到捕食者時(shí)的恐懼程度, 參數(shù)c刻畫(huà)捕食轉(zhuǎn)化程度. 函數(shù)bu2vbHu2+u+g表示修正的Holling-Ⅱ型功能反應(yīng)函數(shù), 其中b表示捕食者的最大搜索速度, H表示捕食者處理一個(gè)食餌所需的時(shí)間, g表示半飽和常數(shù), 對(duì)應(yīng)于搜索速率等于最大值b的一半時(shí)的食餌數(shù)量.

    1 平衡點(diǎn)的存在性和穩(wěn)定性

    定義瘙綆2+={(u,v)u≥0, v≥0}. 在初始條件u≥0, v≥0下, 系統(tǒng)(1)對(duì)應(yīng)的非線性項(xiàng)滿足局部Lipschitz條件且連續(xù)可微, 因此系統(tǒng)(1)存在局部解.

    定理1 定義Ω=(u(t),v(t))cu(t)+v(t)≤c(r-d+m)24am, 則系統(tǒng)(1)的解一致最終有界.

    證明: 令N(t)=cu(t)+v(t), 將N(t)沿系統(tǒng)(1)的軌線求導(dǎo), 得

    dN(t)dt= "cu′(t)+v′(t)=cru1+kv-cdu-cau2-mv≤ "cru-cdu-cau2+cmu-mN

    = "c(r-d+m)u-cau2-mN≤c(r-d+m)24a-mN,

    N(t)≤c(r-d+m)24am+N(0)-c(r-d+m)24ame-mt.

    則當(dāng)t→∞時(shí), 有N(t)≤c(r-d+m)24am. 從而任給系統(tǒng)(1)一個(gè)初值, 系統(tǒng)(1)的所有解最終進(jìn)入?yún)^(qū)域Ω=(u(t),v(t))cu(t)

    +v(t)≤c(r-d+m)24am, 所以Ω是系統(tǒng)(1)的正不變集, 吸引瘙綆2+中的所有正解, 即系統(tǒng)(1)的解是滿足一致有界的.

    定理2 1) 系統(tǒng)(1)一直存在一個(gè)零平衡點(diǎn)E0=(0,0);

    2) 當(dāng)rgt;d時(shí), 系統(tǒng)(1)存在一個(gè)邊界平衡點(diǎn)E1=r-da,0;

    3) 當(dāng)rgt;d且[b(r-d)(c-mH)-ma](r-d)gt;ma2g時(shí), 系統(tǒng)(1)存在共存平衡點(diǎn)E2=(u*,v*), 其中u*=m+m2+4mbg(c-mH)2b(c-mH), v*在證明中給出.

    證明: 系統(tǒng)(1)的所有平衡點(diǎn)都滿足:

    ur1+kv-d-au-buvbHu2+u+g=0,

    vcbu2bHu2+u+g-m=0.(2)

    顯然, 滅絕平衡點(diǎn)E0=(0,0)總存在. 當(dāng)rgt;d時(shí), 邊界平衡點(diǎn)E1=r-da,0. 下面考慮共存平衡點(diǎn)E2=(u*,v*)的存在性, 由式(2)的第二個(gè)方程得

    (cb-mbH)u2-mu-mg=0,

    解得u=m±m(xù)2+4mgb(c-mH)2b(c-mH).

    因?yàn)閏-mHgt;0, 所以

    u=m+m2+4mgb(c-mH)2b(c-mH)u*.(3)

    將式(3)代入式(2)中的第一式, 則v*滿足方程

    M1v*2+M2v*+M3=0,(4)

    其中

    M1=bku*gt;0,M2=k(d+au*)(bHu*2+u*+g)+bu*gt;0,M3=(d+au*-r)(bHu*2+u*+g).

    下面分兩種情形討論:

    情形1) 當(dāng)M3lt;0時(shí), 式(4)有一正根

    v*=-M2+M22-4M1M32M1;

    情形2) 當(dāng)M3≥0時(shí), 式(4)無(wú)正根. 因?yàn)镸3lt;0, 所以d+au*-rlt;0, 從而

    u*lt;r-da.(5)

    將式(3)代入式(5)得

    [b(r-d)(c-mH)-ma](r-d)gt;ma2g.

    從而結(jié)論得證.

    定理3 1) 當(dāng)r≤d時(shí), 滅絕平衡點(diǎn)E0=(0,0)全局漸近穩(wěn)定;

    2) 當(dāng)rgt;d時(shí), 滅絕平衡點(diǎn)E0=(0,0)不穩(wěn)定.

    證明: 構(gòu)造Lyapunov函數(shù)

    V(t)=cu(t)+v(t),

    其中c為正常數(shù). 顯然, V(t)是在原點(diǎn)鄰域內(nèi)的正定函數(shù). 從而V(t)沿著系統(tǒng)(1)軌線的全導(dǎo)數(shù)為

    V′(t)=cu(r-d)1+kv-cdkuv1+kv-cau2-mv.

    當(dāng)r≤d時(shí), 對(duì)任意的u≥0和v≥0, 有V′(t)≤0, 則V′(t)是半負(fù)定的. 又因?yàn)榧?/p>

    D={(u,v)V′(t)=0}={(0,0)},

    而集合D內(nèi)除(0,0)外不再包含系統(tǒng)(1)的其他軌線. 由Lyapunov-LaSalle不變集原理知, 當(dāng)r≤d時(shí), 滅絕平衡點(diǎn)E0=(0,0)全局漸近穩(wěn)定. 此外, 系統(tǒng)(1)在E0=(0,0)處的Jacobi矩陣為

    JE0=r-d00-m,

    JE0的特征值為λ1=r-d和λ2=-m. 從而當(dāng)rgt;d時(shí), λ1gt;0, 因此滅絕平衡點(diǎn)E0不穩(wěn)定.

    定理4 若rgt;d, 則:

    1) 當(dāng)clt;mH時(shí), 邊界平衡點(diǎn)E1=r-da,0是局部漸近穩(wěn)定的;

    2) 當(dāng)cgt;mH且[b(r-d)(c-mH)-ma](r-d)lt;ma2g時(shí), E1是全局漸近穩(wěn)定的;

    3) 當(dāng)cgt;mH且[b(r-d)(c-mH)-ma](r-d)gt;ma2g時(shí), E1不穩(wěn)定.

    證明: 在邊界平衡點(diǎn)E1=r-da,0處, 系統(tǒng)(1)的Jacobi矩陣為

    JE1=d-r-kr(r-d)a-b(r-d)2bH(r-d)2+a(

    r-d)+a2g0bc(r-d)2bH(r-d)2+a(r-d)+a2g-m,

    則求得JE1的特征值為

    λ1=d-rlt;0," λ2=bc(r-d)2bH(r-d)2+a(r-d)+a2g-m.

    因?yàn)棣?lt;0等價(jià)于bc(r-d)2bH(r-d)2+a(r-d)+a2g-mlt;0, 所以[b(r-d)(c-mH)-ma](r-d)lt;ma2g.

    從而當(dāng)[b(r-d)(c-mH)-ma](r-d)lt;ma2g時(shí), 邊界平衡點(diǎn)E1是局部漸近穩(wěn)定的結(jié)點(diǎn); 當(dāng)cgt;mH且[b(r-d)(c-mH)-ma](r-d)gt;ma2g時(shí), 邊界平衡點(diǎn)E1是不穩(wěn)定的鞍點(diǎn).

    由定理2可知, 系統(tǒng)(1)除平衡點(diǎn)E0和E1外沒(méi)有其他的平衡點(diǎn). 由于E0是不穩(wěn)定的平衡點(diǎn), E1是局部漸近穩(wěn)定的平衡點(diǎn), 因此系統(tǒng)(1)在瘙綆2

    +內(nèi)不存在周期解, 從而可知E1是全局漸近穩(wěn)定的平衡點(diǎn).

    定理5 若當(dāng)rgt;d且[b(r-d)(c-mH)-ma](r-d)gt;ma2g時(shí), 共存在平衡點(diǎn)E2=(u,v)存在, 則:

    1) 當(dāng)ac2bu*3-m2bHu*2v*+m2gv*gt;0時(shí), 共存在平衡點(diǎn)E2是局部漸近穩(wěn)定的;

    2) 當(dāng)ac2bu*3-m2bHu*2v*+m2gv*lt;0時(shí), 共存平衡點(diǎn)E2是不穩(wěn)定的.

    證明: 系統(tǒng)(1)在E2=(u,v)處的Jacobi矩陣為

    JE2=-au+b2Hu3v-bguv(bHu2+u+g)2-kru(1+kv)2-bu2bHu2+u+gbcu2v+2bcguv(bHu2+u+g)20,

    對(duì)應(yīng)的特征方程為

    λ2-tr(JE2)λ+det(JE2)=0,(6)

    其中

    tr(JE2)=-au+b2Hu3v-bguv(bHu2+u+g)2,

    det(JE2)=bcu2v+2bcguv(bHu2+u+g)2kru(

    1+kv)2+bu2bHu2+u+g.

    當(dāng)ac2bu*4-m2bHu*2v*+m2gv*gt;0時(shí), 可計(jì)算方程(6)對(duì)應(yīng)特征值λ1,λ2的實(shí)部都小于零. 由Routh-Hurwitz穩(wěn)定性判據(jù)

    [20], E2是局部漸近穩(wěn)定平衡點(diǎn). 當(dāng)ac2bu*4-m2bHu*2v*+m2gv*lt;0時(shí), 方程(6)的特征值λ1,λ2的實(shí)部均大于零, 則E2是不穩(wěn)定的.

    定理6 當(dāng)[b(r-d)(c-mH)-ma](r-d)gt;ma2g, k≥1且3a2bHg≥b2H2(r-d)2+abH(r-d)+a2時(shí), 共存平衡點(diǎn)E2是全局漸近穩(wěn)定的.

    證明: 設(shè)Dulac函數(shù)為B(u,v)=(1+kv)(bHu2+u+g)u-1vβ-1, 其中參數(shù)β待定, 則

    D=(P(u,v)B(u,v))u+(Q(u,v)B(u,v))v=u-1vβ-1[f1(u,β)v2+f2(u,β)v+f3(u,β)],

    其中,

    f1(u,β)=-bku,

    f2(u,β)=-dk(2bHu2+u)-ak(3bHu3+2u2+gu)-bu+cbk(β+1)u2-mk(β+1)(bHu2+u+g)=-3abkHu3-(2dbkH+2ak)u2-(d

    k+akg+b)u,f3(u,β)= "(r-d)(2bHu2+u)-a(3bHu3+2u2+gu)+cbβu2-mβ(bHu2+u+g)= "-3abHu3+[2bH(r-d)-2a]u2+(r-d-ag)u.

    易知, f1(u,β)lt;0, 當(dāng)k≥1時(shí), 有

    f2(u,β)-f3(u,β)=-r(2bHu2+u)+(2bHu2+u)(1-k)+a(3bHu3+2u2+gu)(1-k)-bu≤0,

    故f2(u,β)≤f3(u,β).

    因?yàn)楫?dāng)f3(u,β)≤0時(shí), f2(u,β)lt;0, 所以D(v)在[0,+∞)上單調(diào)遞減, 而最大值D(0)=f3(u,β). 因此要使D≤0對(duì)任意的(u,v)∈瘙綆2+成立, 只需

    f3(u,β)≤0," u∈[0,+∞),

    只需證

    -3abHu3+[2bH(r-d)-2a]u2+(r-d-ag)u≤0.

    因?yàn)閡gt;0, 所以只需證

    -3abHu2+[2bH(r-d)-2a]u+r-d-ag≤0,(7)

    又因?yàn)?a2bHg≥b2H2(r-d)2+abH(r-d)+a2, 所以式(7)得證, 從而f3(u,β)≤0得證.

    進(jìn)一步, 利用Bendixson-Dulac定理可得系統(tǒng)(1)不存在周期軌道. 因此兩個(gè)平衡點(diǎn)E0,E1是不穩(wěn)定的, 但E2為唯一的局部穩(wěn)定正平衡點(diǎn). 從而E2是全局漸近穩(wěn)定平衡點(diǎn).

    2 極限環(huán)的存在性和Hopf分支

    定理7 若ac2bu*3-m2bHu*2v*+m2gv*lt;0, 則系統(tǒng)(1)至少存在一個(gè)包含共存平衡點(diǎn)E2的穩(wěn)定極限環(huán).

    證明: 由定理5知, 當(dāng)ac2bu*3-m2bHu*2v*+m2gv*lt;0時(shí), E2是不穩(wěn)定的. 為證明極限環(huán)的存在性, 需構(gòu)造Poincare-Bendixson環(huán)域的外境線L.

    首先, 考慮直線

    L1u-u=0,

    其中u=r-da. 當(dāng)vgt;0時(shí), 有

    dL1dtL1=0=ur1+kv-d-au-buvbHu2+u+g≤-bu2vbHu2+u+glt;0,

    所以當(dāng)軌線與直線L1=0相遇時(shí), 均從直線L1=0的右方穿入左方.

    其次, 考慮直線L2cu+v-μ=0,

    其中μgt;0, 0lt;u≤u, 則

    dL2dtL2=0= "cru1+kv-cdu-cau2-mv≤cru-cdu-cau2+cmu-mμ

    = "-cau2+c(r-d+m)u-mμF(u).

    此時(shí), F(u)是關(guān)于u的一個(gè)開(kāi)口向下的二次函數(shù), 所以對(duì)于足夠大的μ, 有F(u)lt;0, 即dL2dtL2=0lt;0, 所以當(dāng)軌線與直線L2=0相遇時(shí), 均從直線L2=0的右上方穿入左下方.

    因?yàn)橹本€u=0和v=0都為系統(tǒng)(1)的軌線, 所以直線L1=0、 L2=0、 u軸和v軸圍成了Poincare-Bendixson環(huán)域的外境線L, 而E2是不穩(wěn)定的

    奇點(diǎn), 邊界上的奇點(diǎn)E0和E1都是鞍點(diǎn), 由Poincare-Bendixson環(huán)域定理[21]知, 系統(tǒng)(1)至少存在一個(gè)包含共存平衡點(diǎn)E2的穩(wěn)定極限環(huán).

    當(dāng)特征方程(6)中tr(JE2)=0時(shí), 此時(shí)方程有一對(duì)純虛特征根±iβ0, 其中β0=-det(JE2)

    . 而由tr(JE2)=0可得

    -au+b2Hu3v-bguv(bHu2+u+g)2=0.(8)

    聯(lián)立式(2),(3),(8)可解得

    k=-m3(d+au*-r)(bHu*2-g)2+abcm2u*3(bHu*2-g)abc2mu*4(d+au*)(bHu*2-g)+a2b2c3u*7k*.

    下面討論當(dāng)恐懼水平k=k*為分支參數(shù), 其余參數(shù)保持不變時(shí), 系統(tǒng)(1)在E2=(u,v)處出現(xiàn)Hopf分支的可能性.

    定理8 若k=k*, 則系統(tǒng)(1)在E2=(u,v)處出現(xiàn)Hopf分支.

    證明: 設(shè)特征方程(6)的特征根λ1,2(k)=α(k)±iβ(k), 代入方程(4)得

    [α(k)±iβ(k)]2-tr(JE2)[α(k)±iβ(k)]+det(JE2)=0,

    分離實(shí)部和虛部得

    α2(k)-β2(k)-tr(JE2)α(k)+det(JE2)=0,

    ±2α(k)tr(JE2)=0,

    解得

    α(k)=12tr(JE2),β(k)=124det(J

    E2)-tr2(JE2),

    當(dāng)k=k*時(shí), α(k*)=12tr(JE2)k=k*=0, β(k*)=det(JE2)k=k*gt;0.

    通過(guò)計(jì)算, 橫截條件為

    dtr(JE2)dk

    k=k*=b2Hu3-bgu(bHu2+u+g)2dvdkk=k*lt;0,

    此時(shí), 說(shuō)明系統(tǒng)(1)滿足Poincare-Andronow-Hopf分支定理[22], 因此系統(tǒng)(1)在E2=(u,v)處出現(xiàn)Hopf分支.

    定理9 設(shè)L=116(pxxxp2y+qxxyp2y-pxypxxpy+pxypyyqx)-116

    qxp2yqxx(qxy+pxx), 當(dāng)Llt;0時(shí), 系統(tǒng)(1)在共存平衡點(diǎn)E2=(u,v)處產(chǎn)生超臨界Hopf分支; 當(dāng)Lgt;0時(shí), 其為亞臨界Hopf分支.

    證明: 令x=u-u, y=v-v, E2=(u,v), 代入系統(tǒng)(1)得

    dxdt= "r(x+u)1+k(y+v)-d(x+u)-a(x+u)2- "b(x+u)2(y+v)

    bH(x+u)2+(x+u)+g=p(u,v),dydt= "cb(x+u)

    2(y+v)bH(x+u)2+(x+u)+g-m(y+v)=q(u,v).(9)

    在(x,y)=(0,0)處分別利用Taylor級(jí)數(shù)將p(u,v),q(u,v)展開(kāi)至3階, 則系統(tǒng)(9)轉(zhuǎn)化為

    dxdt=px(0,0)x+py(0,0)y+12pxx(0,0)x2+pxy(

    0,0)xy+12pyy(0,0)y2+" 16pxxx(0,0)x3+12

    pxxy(0,0)x2y+12pxyy(0,0)xy2+16pyyy(0,0)y3+…,

    dydt=qx(0,0)x+qy(0,0)y+12qxx(0,0)x2+qxy(0,0)xy+12

    qyy(0,0)y2+" 16qxxx(0,0)x3+12qxxy(0,0)x

    2y+12qxyy(0,0)xy2+16qyyy(0,0)y3+…,(10)

    其中,

    px(0,0)=-au+b2Hu3v-bguv(bHu2+u+g)2," py

    (0,0)=-kru(1+kv)2-bu2bHu2+u+g,

    pxx(0,0)=-a+3b2Hu2v-bgv(bHu2+u+g)2-2(2bHu+1)(b2Hu3v

    -bguv)(bHu2+u+g)2,pxx(0,0)=-a+3b2Hu2v-bgv(

    bHu2+u+g)2-2(2bHu+1)(b2Hu3v-bguv)(bHu2+u+g)2,

    pxy(0,0)=b2Hu3-bgu(bHu2+u+g)2," pyy(0,0)=2k2ru(1+kv)3,

    pxxx(0,0)=6b2Hu(bHu2+u+g)2-"""" 4b3H3u4v+

    16b3H3u3v+2b2Hu3v-8b2Hguv-2bguv-2bgv(bHu2+u+g)3

    +"""" 6(2bHu+1)2(b2Hu3v-bguv)(bHu2+u+g)4,

    pxxy(0,0)=3b2Hu2-bg(bHu2+u+g)2-2(2bHu+1)

    (b2Hu3-bgu)(bHu2+u+g)2,pxyy(0,0)=0," pyyy(0,0)=-6k3ru(1+kv)4,

    qx(0,0)=bcu2v+2bcguv(bHu2+u+g)2," qy(0,0)=0,

    qxx(0,0)=2bcuv+2bcgv(bHu2+u+g)2-2(2bHu+1)(bcu2v+2bcguv)(bHu2+u+g)3,

    qxy(0,0)=bcu2+2bcgu(bHu2+u+g)2," qyy(0,0)=0,

    qxyy(0,0)=0," qyyy(0,0)=0,qxxy(0,0)=2bcu+2bcg(bHu2+u+g)2

    -2(2bHu+1)(bcu2+2bcgu)(bHu2+u+g)2,

    qxxx(0,0)=2bcv(bHu2+u+g)2-"""" 8b2cHu3v+8b2cHu2v+12b2cu3v+16

    b2cHguv+8bcuv+8bcgv(bHu2+u+g)3+"""" 6(2bHu+1)2(bcu2v+2bcguv)(bHu2+u+g)4.

    去掉式(8)的高階4次項(xiàng), 然后將式(8)改寫(xiě)成

    =JE2X+G(X),

    其中,

    X=xy, G=G1G2=12pxxx

    2+pxyxy+12pyyy2+16pxxxx3+12pxxyx2y+16pyyyy312qxxx2+qxyxy+16qxxxx3+12qxxyx2y.

    當(dāng)k=k*時(shí), px=0, Jacobi矩陣JE2的一個(gè)特征值是純虛數(shù)iβ0, 其中β0=i-pyqx, 對(duì)應(yīng)的特征向量v=(py,i-pyqx)T, 不妨設(shè)

    Y=(Re v,Im v)=py00--pyqx, JE2=0p

    yqx0, Y-1=1py00-1-pyqx.

    令X=YW, 則W=Y-1X, 其中W=(w1,w2)T, 得

    =(Y-1JE2Y)W+Y-1

    G(YW),

    12=0--pyqx-pyqx0w1w2+G1(w1,w2)G2(w1,w2),

    其中,

    G1(w1,w2)= "1py12pxxp2yw21-pxypy-pyqxw

    1w2-12pyypyqxw22+16pxxxp3yw31- "12pxxyp2y-pyqxw21w2+16pyyypyqx-pyqxw32,

    G2(w1,w2)= "1-pyqx12qxxp2yw21-qxypy-pyqxw1w2+16pxxxp3yw31-12qxxyp2y-pyqxw21w2.

    由于Hopf分支的方向由第一Lyapunov系數(shù)的符號(hào)決定, 所以下面計(jì)算第一Lyapunov系數(shù):

    L= "1163G1w31+3G1w1w22+3G2w21w

    2+3G2w32+116-pyqx2G1w1w22G1w21+2G1w22- "2G2w1w22G2w21+2G2w22-2G1w21·2G2w2

    1+2G1w22·2G2w22,

    化簡(jiǎn)得

    L=116(pxxxp2y+qxxyp2y-pxypxxpy+pxypyyqx)-116qxqxxp

    2y(qxy+qxx).

    由Poincare-Andronow-Hopf分支定理知: 當(dāng)Llt;0時(shí), 系統(tǒng)(1)在共存平衡點(diǎn)E2=(u,v)處產(chǎn)生Hopf分支為超臨界分支; 當(dāng)Lgt;0時(shí), 該Hopf分支為亞臨界分支.

    綜上所述, 本文在均勻空間分布下, 建立了一個(gè)具有恐懼效應(yīng)及修正的Holling-Ⅱ捕食者-食餌模型, 并研究了恐懼因子k、 捕食者的最大搜索速度b和捕食者處理一個(gè)食餌

    所需的時(shí)間H對(duì)系統(tǒng)(1)動(dòng)力學(xué)行為的影響. 理論分析和計(jì)算結(jié)果表明: 1) 恐懼因子k對(duì)滅絕平衡點(diǎn)和邊界平衡點(diǎn)的穩(wěn)定性沒(méi)有影響, 但當(dāng)k發(fā)生變化時(shí), 對(duì)共存平衡點(diǎn)

    有影響; 2) 捕食者的最大搜索速度b和捕食者處理一個(gè)食餌所需的時(shí)間H對(duì)滅絕平衡點(diǎn)的穩(wěn)定性沒(méi)有影響, 但對(duì)邊界平衡點(diǎn)和共存平衡點(diǎn)的穩(wěn)定性都有影響, 并且當(dāng)H充分大

    時(shí), 系統(tǒng)(1)不存在共存平衡點(diǎn); 3) 當(dāng)滿足定理8的條件時(shí), 系統(tǒng)(1)存在一個(gè)包含共存平衡點(diǎn)的穩(wěn)定極限環(huán); 4) 以恐懼因子k=k*為分支參數(shù), 系統(tǒng)(1)在共存平衡點(diǎn)處出現(xiàn)Hopf分支.

    參考文獻(xiàn)

    [1] FU S M, ZHANG H S. Effect of Hunting Cooperation on the Dynamic Behavior for

    a Diffusive Holling Type Ⅱ Predator-Prey Model [J]. Communications in Nonlinear Science and Numerical Simulation, 2021, 99: 105807\|1\|105807\|23.

    [2] 張萌萌, 李善兵. 具有恐懼效應(yīng)和空間異質(zhì)捕食-食餌模型的穩(wěn)態(tài)解 [J]. 吉林大學(xué)學(xué)報(bào)(理學(xué)版), 2022, 60(4): 775-783. (ZHANG M M, LI S B

    . Steady State Solutions of Predator-Prey Model with Fear Effect and Spatial He

    terogeneity [J]. Journal of Jilin University (Science Edition), 2022, 60(4): 775-783.)

    [3] QI H K, MENG X Z, HAYAT T, et al. Bifurcation Dynamics of a Reaction-Diffusion P

    redator-Prey Model with Fear Effect in a Predator-Poisoned Environment [J]. Mathematical Methods in the Applied Sciences, 2022, 45(10): 6217-6254.

    [4] 高鶴, 李秀玲. 具時(shí)滯的捕食-食餌共生模型的Hopf分支 [J]. 吉林大學(xué)學(xué)報(bào)(理學(xué)版), 2023, 61(6): 1339-1350. (GAO H, LI X L. Hopf

    Bifurcation of Predator-Prey Symbiotic Model with Time Delay [J]. Journal of Jilin University (Science Edition), 2023, 61(6): 1339-1350.)

    [5] HOLLING C S. The Functional Response of Predators to

    Prey Density and Its Role in Mimicry and Population Regulation [J]. The Memoirs of the Entomological Society of Canada, 1965, 97(S45): 5-60.

    [6] ANDREWS J F. A Mathematical Model for the Continuous

    Culture of Microorganisms Utilizing Inhibitory Substrates [J]. Biotechnology and Bioengineering, 1968, 10(6): 707-723.

    [7] ARDITI R, GINZBURG L R. Coupling in Predator-Prey D

    ynamics: Ratio-Dependence [J]. Journal of Theoretical Biology, 1989, 139(3): 311-326.

    [8] BEDDINGTON J R. Mutual Interference between Parasites or Predators and Its Ef

    fect on Searching Efficiency [J]. The Journal of Animal Ecology, 1975, 44(1): 331-340.

    [9] DEANGELIS D L, GOLDSTEIN R A, O’NEILL R V. A Model for Tropic Interaction [J]. Ecology, 1975, 56(2): 881-892.

    [10] CROWLEY P H, MARTIN E K. Functional Responses and I

    nterference within and between Year Classes of a Dragonfly Population [J]. Journal of the North American Benthological Society, 1989, 8(3): 211-221.

    [11] DALZIEL B D, THOMANN E, MEDLOCK J, et al. Global An

    alysis of a Predator-Prey Model with Variable Predator Search Rate [J]. Journal of Mathematical Biology, 2020, 81(1): 159-183.

    [12] ZANETTE L Y, WHITE A F, ALLEN M C, et al. Perceived Predation Risk Reduces t

    he Number of Offspring Songbirds Produce per Year [J]. Science, 2011, 334: 1398-1401.

    [13] CREEL S, CHRISTIANSON D, LILEY S, et al. Predation Risk Affects Reproductive

    Physiology and Demography of Elk [J]. Science, 2007, 315: 960.

    [14] SHERIFF M J, KREBS C J, BOONSTRA R. The Sensitive Hare: Sublethal Effects of P

    redator Stress on Reproduction in Snowshoe Hares [J]. Journal of Animal Ecology, 2009, 78(6): 1249-1258.

    [15] WANG X Y, ZANETTE L, ZOU X F. Modelling the Fear Effect in Predator-Prey Intera

    ctions [J]. Journal of Mathematical Biology, 2016, 73(5): 1179-1204.

    [16] PAL S, PAL N, SAMANTA S, et al. Effect of Hunting C

    ooperation and Fear in a Predator-Prey Model [J]. Ecological Complexity, 2019, 39: 100770\|1\|100770\|18.

    [17] PAL S, PAL N, SAMANTA S, et al. Fear Effect in Prey

    and Hunting Cooperation among Predators in a Leslie-Gower Model [J]. Mathematical Biosciences and Engineering, 2019, 16(5): 5146-5179.

    [18] LAI L Y, ZHU Z L, CHEN F D. Stability and Bifurcatio

    n in a Predator-Prey Model with the Additive Allee Effect and the Fear Effect [J]. Mathematics, 2020, 8(8): 1280-1-1280-21.

    [19] SASMAL S K. Population Dynamics with Multiple Allee Effects Induced by Fear F

    actors—A Mathematical Study on Prey-Predator Interactions [J]. Applied Mathematical Modelling, 2018, 64: 1-14.

    [20] 王靈芝. 具有恐懼效應(yīng)的時(shí)滯捕食者-食餌模型 [J]. 吉林大學(xué)學(xué)報(bào)(理學(xué)版), 2023,

    61(3): 449-458. (WANG L Z. Delayed Predator-Prey Model with Fear Effect [J]. Journal of Jilin University (Science Edition), 2023, 61(3): 449-458.)

    [21] 馬知恩, 周義倉(cāng), 李承治. 常微分方程定性與穩(wěn)定性方法 [M]. 2版. 北京: 科學(xué)出版社, 20

    15: 187-194. (MA Z E, ZHOU Y C, LI C Z. Qualitative and Stability Methods for Ordinary Differential Equations [M]. 2nd ed. Beijing: Science Press, 2015: 187-194.)

    [22] HASSARD B D, KAZARINOFF N D, WAN Y H. Theory and Applications of Hopf Bifurca

    tion [M]. Cambridge: Cambridge University Press, 1981: 10-100.

    (責(zé)任編輯: 趙立芹)

    猜你喜歡
    捕食者穩(wěn)定性
    羊的“視界”大不同
    一類(lèi)k-Hessian方程解的存在性和漸近穩(wěn)定性
    SBR改性瀝青的穩(wěn)定性評(píng)價(jià)
    石油瀝青(2021年4期)2021-10-14 08:50:44
    交錯(cuò)擴(kuò)散對(duì)具有Ivlev型功能反應(yīng)的捕食模型共存解存在性的作用
    一類(lèi)具有時(shí)滯的捕食者-獵物-共生者系統(tǒng)的研究
    具有Allee效應(yīng)隨機(jī)追捕模型的滅絕性
    一類(lèi)隨機(jī)食餌-捕食者模型的參數(shù)估計(jì)
    非線性中立型變延遲微分方程的長(zhǎng)時(shí)間穩(wěn)定性
    半動(dòng)力系統(tǒng)中閉集的穩(wěn)定性和極限集映射的連續(xù)性
    瘋狂的捕食者
    中外文摘(2016年13期)2016-08-29 08:53:27
    日韩精品有码人妻一区| av专区在线播放| 亚洲va日本ⅴa欧美va伊人久久| 久久久久国产精品人妻aⅴ院| 一卡2卡三卡四卡精品乱码亚洲| 九九在线视频观看精品| 天堂动漫精品| 嫩草影视91久久| 精品久久久久久久人妻蜜臀av| 亚洲欧美精品综合久久99| 国产精品爽爽va在线观看网站| 国产精品一及| 亚洲一级一片aⅴ在线观看| 91午夜精品亚洲一区二区三区 | 国产精品女同一区二区软件 | 亚洲真实伦在线观看| 国产精品免费一区二区三区在线| 久久国内精品自在自线图片| 在线播放国产精品三级| 亚洲精品一卡2卡三卡4卡5卡| 97人妻精品一区二区三区麻豆| 国产精品自产拍在线观看55亚洲| 国产亚洲欧美98| 在线国产一区二区在线| 色在线成人网| 亚洲av第一区精品v没综合| 亚洲成a人片在线一区二区| 成人欧美大片| 22中文网久久字幕| 久久精品综合一区二区三区| 国产美女午夜福利| 色综合婷婷激情| 日韩一本色道免费dvd| 亚洲内射少妇av| 精品欧美国产一区二区三| 高清在线国产一区| 午夜福利在线观看免费完整高清在 | 欧美国产日韩亚洲一区| 欧美日韩综合久久久久久 | 在线看三级毛片| av专区在线播放| 国产高清视频在线播放一区| 亚洲无线在线观看| 性欧美人与动物交配| 在线播放无遮挡| 我的老师免费观看完整版| 亚洲av免费在线观看| 成人无遮挡网站| 国产av一区在线观看免费| 网址你懂的国产日韩在线| 最新中文字幕久久久久| 国产精品久久视频播放| 亚洲在线观看片| 动漫黄色视频在线观看| 精品久久国产蜜桃| 欧美丝袜亚洲另类 | 欧美性猛交黑人性爽| 欧美xxxx性猛交bbbb| 欧美最黄视频在线播放免费| 国产精品一区二区性色av| 国产男人的电影天堂91| 午夜福利视频1000在线观看| 99热这里只有是精品50| 国内精品美女久久久久久| 亚洲国产精品sss在线观看| 美女高潮喷水抽搐中文字幕| 一级毛片久久久久久久久女| 99riav亚洲国产免费| 中出人妻视频一区二区| 成人特级黄色片久久久久久久| 男女下面进入的视频免费午夜| 国产高清有码在线观看视频| av黄色大香蕉| 国产视频一区二区在线看| 精品一区二区免费观看| 欧美xxxx黑人xx丫x性爽| 国产v大片淫在线免费观看| 伊人久久精品亚洲午夜| 久久久久久久久大av| 亚洲欧美清纯卡通| 免费一级毛片在线播放高清视频| 日韩欧美 国产精品| 亚洲精品成人久久久久久| 国产蜜桃级精品一区二区三区| 久久精品国产亚洲av香蕉五月| 亚洲18禁久久av| 成年女人永久免费观看视频| 亚洲无线在线观看| 一区二区三区四区激情视频 | 成人国产综合亚洲| 亚洲精品亚洲一区二区| 三级国产精品欧美在线观看| 一区二区三区激情视频| 淫妇啪啪啪对白视频| 伊人久久精品亚洲午夜| 老熟妇仑乱视频hdxx| 国产午夜精品久久久久久一区二区三区 | 99热这里只有是精品50| 日本熟妇午夜| 亚洲国产高清在线一区二区三| 男人舔女人下体高潮全视频| 久久久久九九精品影院| 在线播放无遮挡| 亚洲精品日韩av片在线观看| 久久久久九九精品影院| 尾随美女入室| www.色视频.com| 亚洲va在线va天堂va国产| 久久久久久久亚洲中文字幕| 久久精品人妻少妇| 91久久精品国产一区二区成人| 亚洲图色成人| 国产一级毛片七仙女欲春2| 少妇的逼好多水| 狂野欧美白嫩少妇大欣赏| 午夜日韩欧美国产| 国产真实伦视频高清在线观看 | 欧美日韩综合久久久久久 | 中文字幕免费在线视频6| 亚洲va日本ⅴa欧美va伊人久久| 18禁黄网站禁片免费观看直播| 啪啪无遮挡十八禁网站| 中文字幕人妻熟人妻熟丝袜美| 国产一级毛片七仙女欲春2| 久久久久久国产a免费观看| 美女 人体艺术 gogo| 中文字幕高清在线视频| 欧美日韩精品成人综合77777| 国产精品久久久久久av不卡| 啦啦啦韩国在线观看视频| 日韩精品青青久久久久久| 超碰av人人做人人爽久久| 欧美另类亚洲清纯唯美| 国产一级毛片七仙女欲春2| 久久久久久大精品| 色尼玛亚洲综合影院| 精品一区二区三区人妻视频| 国产高清有码在线观看视频| 欧美最黄视频在线播放免费| 一a级毛片在线观看| 最近最新免费中文字幕在线| 天堂av国产一区二区熟女人妻| 亚洲 国产 在线| 国产伦人伦偷精品视频| av在线老鸭窝| 国产精品无大码| 天天躁日日操中文字幕| 午夜视频国产福利| 日韩精品有码人妻一区| bbb黄色大片| 日本熟妇午夜| 天堂动漫精品| 桃色一区二区三区在线观看| 亚洲av五月六月丁香网| 欧美日韩中文字幕国产精品一区二区三区| 欧美另类亚洲清纯唯美| 一区二区三区高清视频在线| 女人十人毛片免费观看3o分钟| 桃色一区二区三区在线观看| 亚洲美女搞黄在线观看 | 欧美成人免费av一区二区三区| 欧美+日韩+精品| 少妇被粗大猛烈的视频| 一区福利在线观看| 热99在线观看视频| 久久久久久伊人网av| 亚洲无线在线观看| 欧美一区二区亚洲| 午夜精品在线福利| 男人狂女人下面高潮的视频| 窝窝影院91人妻| 网址你懂的国产日韩在线| 97碰自拍视频| 久久久成人免费电影| 亚洲性夜色夜夜综合| 久久久久国内视频| 免费人成在线观看视频色| 国产男靠女视频免费网站| 丰满的人妻完整版| 欧美性感艳星| 日韩欧美精品免费久久| 干丝袜人妻中文字幕| 午夜精品一区二区三区免费看| 特级一级黄色大片| 日韩欧美国产一区二区入口| 制服丝袜大香蕉在线| 国产大屁股一区二区在线视频| ponron亚洲| 狠狠狠狠99中文字幕| 国产 一区 欧美 日韩| 国产免费一级a男人的天堂| 91在线观看av| 高清日韩中文字幕在线| 91精品国产九色| 亚洲中文字幕一区二区三区有码在线看| 狂野欧美激情性xxxx在线观看| 免费在线观看成人毛片| 日本黄色视频三级网站网址| 有码 亚洲区| 午夜精品在线福利| 老师上课跳d突然被开到最大视频| 国语自产精品视频在线第100页| 天堂√8在线中文| 国产成人av教育| 国产精品野战在线观看| 欧美精品啪啪一区二区三区| 精品久久久久久,| 欧美一区二区亚洲| 干丝袜人妻中文字幕| 国产精品99久久久久久久久| 国产精品福利在线免费观看| 日本爱情动作片www.在线观看 | 久久99热6这里只有精品| 99精品在免费线老司机午夜| 22中文网久久字幕| 亚洲狠狠婷婷综合久久图片| 婷婷丁香在线五月| 国产免费男女视频| 日韩精品有码人妻一区| 在线观看午夜福利视频| 黄色丝袜av网址大全| 熟女电影av网| 精品午夜福利在线看| 给我免费播放毛片高清在线观看| 久久人妻av系列| 成人毛片a级毛片在线播放| 国产成年人精品一区二区| 亚洲av熟女| 欧美中文日本在线观看视频| 中文字幕高清在线视频| 日本五十路高清| 欧美性感艳星| 男人狂女人下面高潮的视频| 欧美精品啪啪一区二区三区| 午夜福利欧美成人| 亚洲精品久久国产高清桃花| 欧美潮喷喷水| 午夜福利18| 在现免费观看毛片| 一个人看的www免费观看视频| 久久久久久国产a免费观看| 99久久成人亚洲精品观看| 男人舔奶头视频| 午夜福利在线观看吧| 日日摸夜夜添夜夜添av毛片 | 国产激情偷乱视频一区二区| 亚洲一级一片aⅴ在线观看| 免费看日本二区| 午夜福利高清视频| 毛片一级片免费看久久久久 | 无遮挡黄片免费观看| 99国产精品一区二区蜜桃av| 麻豆成人午夜福利视频| 国产精品福利在线免费观看| 一进一出好大好爽视频| 天美传媒精品一区二区| 在线播放国产精品三级| 不卡一级毛片| 亚洲狠狠婷婷综合久久图片| 国产人妻一区二区三区在| 成人欧美大片| 麻豆久久精品国产亚洲av| 床上黄色一级片| 99热6这里只有精品| 搞女人的毛片| 国产蜜桃级精品一区二区三区| 国产精品嫩草影院av在线观看 | 最新中文字幕久久久久| 在线看三级毛片| 麻豆久久精品国产亚洲av| 最近最新免费中文字幕在线| 伦理电影大哥的女人| 国产精品日韩av在线免费观看| 99热这里只有是精品在线观看| netflix在线观看网站| 91av网一区二区| 国产人妻一区二区三区在| 99视频精品全部免费 在线| 欧美性感艳星| 在线观看舔阴道视频| 国产国拍精品亚洲av在线观看| 免费看日本二区| 91午夜精品亚洲一区二区三区 | 免费看a级黄色片| 国模一区二区三区四区视频| a级一级毛片免费在线观看| 色在线成人网| 久久天躁狠狠躁夜夜2o2o| 色综合婷婷激情| 国产亚洲91精品色在线| 国产 一区 欧美 日韩| 少妇丰满av| 亚洲成a人片在线一区二区| 国产精品综合久久久久久久免费| 99riav亚洲国产免费| 色综合色国产| 男女之事视频高清在线观看| 搞女人的毛片| 日本免费a在线| 赤兔流量卡办理| 国产欧美日韩精品亚洲av| 欧美丝袜亚洲另类 | 午夜福利在线观看吧| 午夜免费男女啪啪视频观看 | 国产亚洲精品av在线| 村上凉子中文字幕在线| 久久久久国内视频| 国产精品亚洲美女久久久| 国产高清激情床上av| 中文字幕精品亚洲无线码一区| 国产精品亚洲美女久久久| 日日撸夜夜添| 一级黄片播放器| 亚洲色图av天堂| 亚洲内射少妇av| 欧美另类亚洲清纯唯美| 亚洲专区国产一区二区| 一进一出好大好爽视频| 国产人妻一区二区三区在| 亚洲va日本ⅴa欧美va伊人久久| 国产精品,欧美在线| 午夜激情福利司机影院| 伦理电影大哥的女人| 久久精品综合一区二区三区| 欧美高清性xxxxhd video| 无遮挡黄片免费观看| 久久久久性生活片| 亚洲精品在线观看二区| 桃色一区二区三区在线观看| 午夜免费成人在线视频| 日本黄色片子视频| 自拍偷自拍亚洲精品老妇| 在线免费观看不下载黄p国产 | 久久精品国产亚洲网站| av天堂在线播放| ponron亚洲| 男插女下体视频免费在线播放| 色综合色国产| 最近最新免费中文字幕在线| 女同久久另类99精品国产91| 俺也久久电影网| 亚洲色图av天堂| 欧美性猛交╳xxx乱大交人| 国产 一区 欧美 日韩| 国产三级中文精品| 欧美性猛交黑人性爽| 干丝袜人妻中文字幕| 老熟妇乱子伦视频在线观看| 成人国产麻豆网| 国产亚洲av嫩草精品影院| 亚洲精华国产精华液的使用体验 | 久久草成人影院| 小说图片视频综合网站| 深夜精品福利| 99热这里只有是精品在线观看| 国产黄片美女视频| 变态另类丝袜制服| 日本黄色片子视频| 三级国产精品欧美在线观看| 69人妻影院| 国产精品电影一区二区三区| 男人舔女人下体高潮全视频| 午夜精品久久久久久毛片777| 国产亚洲精品av在线| 日日啪夜夜撸| 亚洲五月天丁香| 最近中文字幕高清免费大全6 | 日韩欧美在线乱码| 成人国产麻豆网| 51国产日韩欧美| 亚洲专区国产一区二区| 亚洲无线在线观看| 中文字幕人妻熟人妻熟丝袜美| 精品人妻视频免费看| 成人性生交大片免费视频hd| 国产免费一级a男人的天堂| 亚洲欧美日韩高清专用| 国产老妇女一区| 俺也久久电影网| 麻豆成人av在线观看| 最近最新免费中文字幕在线| 一区二区三区激情视频| 午夜精品在线福利| 成人毛片a级毛片在线播放| 人妻少妇偷人精品九色| 亚洲精品亚洲一区二区| 啦啦啦啦在线视频资源| 午夜免费成人在线视频| 桃色一区二区三区在线观看| 午夜福利18| 极品教师在线视频| 日韩精品有码人妻一区| 在线播放国产精品三级| 91麻豆精品激情在线观看国产| 真实男女啪啪啪动态图| 亚洲自偷自拍三级| 可以在线观看毛片的网站| 国产精品一区二区三区四区久久| 国产在线男女| 我要看日韩黄色一级片| 成人高潮视频无遮挡免费网站| 欧美在线一区亚洲| 国产精品亚洲美女久久久| 久久人人爽人人爽人人片va| 国产av一区在线观看免费| 国产一区二区在线观看日韩| 99久国产av精品| 午夜影院日韩av| 国产精华一区二区三区| 亚洲午夜理论影院| 国产精品三级大全| 国产高清三级在线| 在线观看一区二区三区| 国内精品久久久久久久电影| 成人特级黄色片久久久久久久| 亚洲av第一区精品v没综合| 色综合婷婷激情| 亚洲熟妇中文字幕五十中出| 久久精品国产亚洲av天美| 午夜免费成人在线视频| 男人和女人高潮做爰伦理| 欧美不卡视频在线免费观看| 国产精品爽爽va在线观看网站| 日日夜夜操网爽| 一a级毛片在线观看| 97热精品久久久久久| 搡老熟女国产l中国老女人| 久久人人爽人人爽人人片va| 日本三级黄在线观看| 日韩中文字幕欧美一区二区| 亚洲av一区综合| netflix在线观看网站| 免费电影在线观看免费观看| 蜜桃亚洲精品一区二区三区| 18禁黄网站禁片免费观看直播| 日韩欧美在线二视频| 国产精品国产三级国产av玫瑰| 联通29元200g的流量卡| 国内毛片毛片毛片毛片毛片| 婷婷亚洲欧美| www.www免费av| 一区二区三区高清视频在线| 美女高潮的动态| 国内揄拍国产精品人妻在线| 精品一区二区三区av网在线观看| 亚洲精品国产成人久久av| 国产伦在线观看视频一区| 岛国在线免费视频观看| 久久精品91蜜桃| 午夜爱爱视频在线播放| 国产 一区精品| 九九爱精品视频在线观看| 一本久久中文字幕| 精品无人区乱码1区二区| 午夜福利高清视频| 熟女人妻精品中文字幕| 99久久久亚洲精品蜜臀av| 国产精华一区二区三区| 亚洲国产高清在线一区二区三| 亚洲va在线va天堂va国产| 夜夜看夜夜爽夜夜摸| 久久久色成人| 成人毛片a级毛片在线播放| 亚洲欧美日韩高清专用| 亚洲中文字幕一区二区三区有码在线看| 男人的好看免费观看在线视频| 一个人看的www免费观看视频| 日韩av在线大香蕉| 国产精品不卡视频一区二区| 国产成人福利小说| 免费一级毛片在线播放高清视频| 欧美日韩中文字幕国产精品一区二区三区| 在线国产一区二区在线| 欧美性感艳星| 日本熟妇午夜| 性插视频无遮挡在线免费观看| 亚洲精品一卡2卡三卡4卡5卡| 最后的刺客免费高清国语| 亚洲国产高清在线一区二区三| 日韩亚洲欧美综合| 欧美一区二区精品小视频在线| 狂野欧美白嫩少妇大欣赏| 国产三级中文精品| 午夜激情福利司机影院| 亚洲乱码一区二区免费版| 蜜桃亚洲精品一区二区三区| eeuss影院久久| 国产久久久一区二区三区| 在线免费观看的www视频| 内地一区二区视频在线| 国产不卡一卡二| 极品教师在线视频| www日本黄色视频网| 国产v大片淫在线免费观看| 美女大奶头视频| 最新在线观看一区二区三区| 日日啪夜夜撸| 亚洲欧美激情综合另类| 神马国产精品三级电影在线观看| 亚洲精品乱码久久久v下载方式| 国产亚洲欧美精品永久| 日本黄色日本黄色录像| 久久久色成人| 亚洲成人中文字幕在线播放| 亚洲人成网站高清观看| 免费大片18禁| 欧美三级亚洲精品| 欧美丝袜亚洲另类| 久久久久久伊人网av| 少妇的逼水好多| 美女福利国产在线 | 国产精品欧美亚洲77777| 久久鲁丝午夜福利片| 国产精品熟女久久久久浪| 91aial.com中文字幕在线观看| 中文字幕人妻熟人妻熟丝袜美| 国产亚洲最大av| 九草在线视频观看| 丝袜脚勾引网站| 久久精品久久久久久噜噜老黄| 国产色婷婷99| 嫩草影院新地址| 日本黄大片高清| 亚洲欧美日韩另类电影网站 | 久久99热6这里只有精品| 99热网站在线观看| 亚洲精品一二三| 国产免费一区二区三区四区乱码| 婷婷色麻豆天堂久久| 亚洲国产精品国产精品| 搡老乐熟女国产| 欧美区成人在线视频| 99热6这里只有精品| www.色视频.com| 日韩在线高清观看一区二区三区| 精品久久久精品久久久| 只有这里有精品99| 国产一区有黄有色的免费视频| 在线观看美女被高潮喷水网站| 免费观看av网站的网址| 91精品国产九色| 高清av免费在线| 精品酒店卫生间| 国产人妻一区二区三区在| 久久久久精品久久久久真实原创| 人妻 亚洲 视频| 18禁在线无遮挡免费观看视频| 久久久久久九九精品二区国产| 熟妇人妻不卡中文字幕| av一本久久久久| 美女脱内裤让男人舔精品视频| 韩国高清视频一区二区三区| 国产精品人妻久久久影院| 免费看光身美女| 在线免费十八禁| 搡女人真爽免费视频火全软件| 亚洲人成网站在线观看播放| 午夜老司机福利剧场| 街头女战士在线观看网站| 精品久久久噜噜| 日日摸夜夜添夜夜爱| 成人毛片a级毛片在线播放| 欧美成人午夜免费资源| 搡老乐熟女国产| 肉色欧美久久久久久久蜜桃| 麻豆精品久久久久久蜜桃| 日韩中字成人| 国产高清三级在线| 久久ye,这里只有精品| 亚洲性久久影院| 人人妻人人澡人人爽人人夜夜| 观看美女的网站| 伦理电影大哥的女人| 最近手机中文字幕大全| 亚洲av欧美aⅴ国产| 日本黄色片子视频| 成人国产麻豆网| 国产色婷婷99| 99九九线精品视频在线观看视频| 亚洲成色77777| 国产人妻一区二区三区在| 亚洲精品视频女| 超碰97精品在线观看| 亚洲精品,欧美精品| 三级国产精品片| 在线观看一区二区三区| 伊人久久国产一区二区| 一级片'在线观看视频| 18+在线观看网站| 91久久精品国产一区二区成人| 国产精品久久久久久久电影| 国产人妻一区二区三区在| 91aial.com中文字幕在线观看| 欧美成人午夜免费资源| 菩萨蛮人人尽说江南好唐韦庄| 成人一区二区视频在线观看| 男女国产视频网站| 黄色配什么色好看| 一区二区三区乱码不卡18| 国产精品偷伦视频观看了| 日韩电影二区| 丰满少妇做爰视频| 精华霜和精华液先用哪个| 亚洲人成网站在线观看播放| 99久久精品国产国产毛片| 久久久久精品久久久久真实原创| 亚洲精品一二三| 在线 av 中文字幕| 日本wwww免费看| 欧美老熟妇乱子伦牲交| 成年人午夜在线观看视频| 久久99热这里只有精品18| 一个人看视频在线观看www免费| 国产精品久久久久久精品电影小说 | 91在线精品国自产拍蜜月|