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

    一種考慮動態(tài)磁滯效應(yīng)的高效穩(wěn)定時域有限元計(jì)算方法

    2023-11-11 06:12:42井立兵
    電工技術(shù)學(xué)報 2023年21期
    關(guān)鍵詞:磁滯回線磁阻磁通

    魏 鵬 陳 龍 賁 彤 井立兵 張 獻(xiàn)

    (1.三峽大學(xué)電氣與新能源學(xué)院 宜昌 443002 2.湖北省微電網(wǎng)工程技術(shù)研究中心(三峽大學(xué)) 宜昌 443002 3.省部共建電工裝備可靠性與智能化國家重點(diǎn)實(shí)驗(yàn)室(河北工業(yè)大學(xué)) 天津 300130)

    0 引言

    時域有限元方法廣泛應(yīng)用于電工裝備的運(yùn)行狀態(tài)預(yù)測與損耗評估之中[1],準(zhǔn)確預(yù)測電機(jī)、變壓器等電工裝備的時域暫態(tài)損耗分布特性不僅對降低設(shè)備能耗、改善局部過熱具有重要意義,同時考慮損耗對暫態(tài)磁場求解的反饋?zhàn)饔茫瑢M(jìn)一步分析電工裝備的瞬時功率平衡、暫態(tài)勵磁特征等電磁特性至關(guān)重要[2-4]。

    目前,對于電工裝備鐵心損耗的有限元計(jì)算方法可大致分為兩類:一是基于場計(jì)算結(jié)果的后處理方法,該類方法利用Steinmetz 損耗模型或Bertotti統(tǒng)計(jì)學(xué)損耗分離理論對僅單值非線性磁場的計(jì)算結(jié)果進(jìn)行后處理,得到每一個單元的損耗,最終進(jìn)行積分得到電工裝備整體的損耗[5-7]。這一類方法主要用于對電機(jī)或變壓器的損耗進(jìn)行頻域分析,并通過傅里葉變換來分析含有諧波激勵時的鐵心損耗。Ansys 公司D.Lin 等基于等效橢圓環(huán)思想,將上述方法進(jìn)行了時域瞬時功率損耗計(jì)算擴(kuò)展[8]。該方法的優(yōu)勢在于在損耗計(jì)算過程中無需迭代,在磁路簡單的結(jié)構(gòu)中具有較好的精度,但由于沒有考慮鐵心磁滯、渦流損耗、剩余損耗等效動態(tài)磁場對鐵心中磁場分布的反饋?zhàn)饔茫陔姍C(jī)等復(fù)雜磁路電工裝備的損耗計(jì)算中將帶來較大誤差。另一類鐵心損耗計(jì)算方法則是在時域有限元分析中,考慮動態(tài)磁滯效應(yīng)對每一瞬時的功率損耗進(jìn)行精確求解[9-10]。這一類方法可最大限度地模擬鐵心的真實(shí)磁化過程,在處理含有磁通密度波形畸變工況下的損耗計(jì)算問題時具有較高精度。同時,在考慮鐵心動態(tài)磁滯效應(yīng)后,鐵心中動態(tài)磁滯回線波形的精確模擬可在時域中對勵磁電流的波形進(jìn)行實(shí)時反饋,進(jìn)而可對電機(jī)或變壓器中暫態(tài)運(yùn)行特性進(jìn)行有效預(yù)測。因此,電工裝備損耗的精確預(yù)測與運(yùn)行狀態(tài)的有效評估都離不開在時域有限元分析時對動態(tài)磁滯特性的影響的考慮。

    然而,目前在有限元分析中常用的磁滯模型,如J-A 模型[11-13]、Preisach 模型[14-17]等均是與磁化速率無關(guān)的磁滯模型,并不能直接考慮電機(jī)或變壓器疊片宏觀渦流損耗與剩余損耗等動態(tài)因素對鐵心磁場分布的影響。為此,需要進(jìn)一步將二者的等效動態(tài)磁場考慮到有限元的計(jì)算方程之中。傳統(tǒng)方法是將動態(tài)效果等效為電導(dǎo)率乘以磁矢位A對時間的一階導(dǎo)數(shù),但該種等效僅適合均勻?qū)嵭膮^(qū)域渦流路徑的求解,并不適合疊片鐵心結(jié)構(gòu)[18-19]。對于鐵心疊片結(jié)構(gòu),進(jìn)行實(shí)際三維建模將會帶來無法承受的計(jì)算成本。為此,E.Dlala 等提出了一種耦合一維疊片有限元方法的二維時域有限元計(jì)算方法[20-21]。在疊片維度,忽略邊緣效應(yīng),認(rèn)為磁通密度僅是疊片厚度方向上的函數(shù),進(jìn)而僅需要求解一維有限元方程,即可得到疊片的渦流損耗等效磁場。然而,該方法實(shí)質(zhì)上是在二維有限元計(jì)算中嵌套了一維有限元程序,仍具有較高的求解難度。為了對上述求解方法進(jìn)行簡化,目前常采用Bertotti 模型的變形形式,即在耦合靜態(tài)磁滯模型的過程中,通過耦合等效微分磁阻率來考慮渦流損耗磁場與剩余損耗磁場的效果,該方法避免了一維有限元的求解,在一定程度上大大簡化了計(jì)算,但在每次迭代過程中仍需更新磁阻率和剛度矩陣,導(dǎo)致計(jì)算效率低下[22-23]。因此,在進(jìn)一步求解含有等效磁阻率的方程時,需要進(jìn)一步耦合非線性磁場計(jì)算中的固定點(diǎn)迭代技術(shù)。在僅考慮靜態(tài)磁滯模型的有限元計(jì)算中,固定點(diǎn)迭代方法已經(jīng)取得了較好的收斂效果與計(jì)算速度[24-26],并已在目前的商業(yè)軟件中取得了較好的計(jì)算效果。然而,在引入渦流損耗磁場與剩余損耗磁場后,收斂將難以得到保證,特別是在計(jì)算至磁通密度零點(diǎn)附近時間步時算法極易發(fā)散,分析其原因主要是渦流損耗、剩余損耗等效動態(tài)磁場的介入改變了固定點(diǎn)法中磁通密度與磁場強(qiáng)度的收斂過程,使得傳統(tǒng)微分磁阻率的定義無法準(zhǔn)確描述收斂路徑的斜率,計(jì)算得到微分磁阻率過低導(dǎo)致考慮動態(tài)磁滯效應(yīng)的有限元方程求解不收斂。

    為解決在考慮動態(tài)磁滯特性的時域有限元計(jì)算難以收斂的問題,本文提出一種基于等效磁阻率的固定點(diǎn)迭代求解策略,進(jìn)而提出一種耦合動態(tài)磁滯特性的時步有限元分析方法。首先,為了獲得電工鋼片的動態(tài)、靜態(tài)磁滯特性,基于愛潑斯坦方圈法搭建了一維磁特性測量系統(tǒng);其次,考慮到逆Preisach模型可以考慮磁化歷史具有較高的計(jì)算精度,同時在耦合有限元計(jì)算時易于數(shù)值實(shí)現(xiàn)的特點(diǎn),本文構(gòu)建了基于參數(shù)化解析一階回轉(zhuǎn)曲線的逆 Preisach模型,并進(jìn)行了參數(shù)辨識;再次,分析了渦流損耗和剩余損耗對算法收斂特性的影響,提出基于等效磁阻率的固定點(diǎn)迭代策略,并在有限元迭代過程中加入松弛因子保證穩(wěn)定性;最后,以愛潑斯坦方圈為例,進(jìn)行了考慮動態(tài)磁滯效應(yīng)的時域有限元分析,求解了瞬時功率損耗特性,并進(jìn)行了實(shí)驗(yàn)對比驗(yàn)證。結(jié)果驗(yàn)證了本文所提方法的準(zhǔn)確性、高效性與穩(wěn)定性。

    1 準(zhǔn)靜態(tài)與動態(tài)磁滯特性測量

    為建立磁滯模型及驗(yàn)證有限元計(jì)算準(zhǔn)確性,基于 IEC 60404—2 標(biāo)準(zhǔn)建立了一維磁性能測量平臺,對牌號為B30P105 的取向硅鋼片進(jìn)行準(zhǔn)靜態(tài)與動態(tài)磁滯特性測量。電工鋼片一維磁特性測量平臺如圖1 所示,該平臺包括寬頻功率放大器、電壓前置放大器、愛潑斯坦方圈、多功能數(shù)據(jù)采集卡和 LabVIEW 控制系統(tǒng)。本實(shí)驗(yàn)采用的硅鋼片規(guī)格及測試條件見表1。分別設(shè)定正弦交流電壓激勵的頻率為5 Hz 與50 Hz,測量得到準(zhǔn)靜態(tài)磁滯回線和50 Hz 動態(tài)磁滯回線,結(jié)果如圖2 和圖3 所示。

    表1 硅鋼片規(guī)格及測試條件Tab.1 Specifications and test conditions of silicon steel sheets

    圖1 電工鋼片一維磁特性測量平臺Fig.1 Measuring platform for 1-D magnetic propertiesof electrical steel sheets

    圖2 5 Hz 準(zhǔn)靜態(tài)磁滯回線測量結(jié)果Fig.2 Measured quasi-static hysteresis loop (5 Hz)

    圖3 50 Hz 動態(tài)磁滯回線測量結(jié)果Fig.3 Measured dynamic hysteresis loop (50 Hz)

    2 靜態(tài)逆Preisach 磁滯模型的構(gòu)建

    本文采用離散形式的靜態(tài)逆Preisach 磁滯模型描述硅鋼片的磁滯特性,該模型采用一種一階回轉(zhuǎn)曲線(First Order Reversal Curves, FORCs)的解析計(jì)算方法,并基于少量準(zhǔn)靜態(tài)磁滯回線測量數(shù)據(jù)實(shí)現(xiàn)解析式的參數(shù)辨識,從而構(gòu)造出辨識逆Everett 函數(shù)所需的一階回轉(zhuǎn)曲線,并進(jìn)一步建立靜態(tài)逆Preisach 磁滯模型。解析一階回轉(zhuǎn)曲線構(gòu)造示意圖如圖4 所示:一條起始于準(zhǔn)靜態(tài)極限磁滯回線下降支的一階回轉(zhuǎn)曲線R-P-N-T,該曲線起點(diǎn)為回轉(zhuǎn)點(diǎn)R,極限磁滯回線上升支Ha和下降支Hd交匯于頂點(diǎn)T。該一階回轉(zhuǎn)曲線構(gòu)造時需首先確定極限磁滯回線在回轉(zhuǎn)點(diǎn)R 處寬度ΔHrev和R 點(diǎn)與T 點(diǎn)的垂直高度ΔBrev為

    圖4 解析一階回轉(zhuǎn)曲線構(gòu)造示意圖Fig.4 Construction ofanalytical first order reversal curve

    圖4 中一階回轉(zhuǎn)曲線上點(diǎn)P 所在磁通密度為BP,磁場強(qiáng)度HP的計(jì)算式為

    式中,ΔHout為極限磁滯回線在高度為BP處的寬度;x表征不同的P 點(diǎn)在一階回轉(zhuǎn)曲線上的相對位置,定義為

    a、b、c為待擬合的模型參數(shù)。當(dāng)滿足a>0,0<b<1,c>0 等條件時,可使得一階回轉(zhuǎn)曲線斜率始終為正且不會超出最大極限磁滯回線。同時,以上模型參數(shù)與回轉(zhuǎn)點(diǎn)R 在下降支上的位置有關(guān),此處定義位置系數(shù)β為

    式中,ΔBout為最大環(huán)磁滯回線的高度。

    參數(shù)a、b、c可表示為β的多項(xiàng)式,即

    式中,y為多項(xiàng)式系數(shù)向量。

    以上多項(xiàng)式系數(shù)通過準(zhǔn)靜態(tài)磁滯回線測量數(shù)據(jù)進(jìn)行辨識,由于對稱性,辨識過程僅需要磁滯回線的上升支。逆Everett 函數(shù)值可利用上升支一階回轉(zhuǎn)曲線Hforc通過式(7)計(jì)算。

    磁通密度幅值為Bm的同心磁滯回線上升支可通過逆Everett 函數(shù)進(jìn)行插值運(yùn)算得到。

    辨識程序中一共使用9 條磁滯回線,將每條磁滯回線上升支根據(jù)磁通密度均勻等分為50 個點(diǎn)?;谑剑?)、式(8)可得到由一階回轉(zhuǎn)曲線計(jì)算磁滯回線的數(shù)值關(guān)系,進(jìn)而建立適應(yīng)度函數(shù)F為

    式中,Hmeas、Hcal分別為磁滯回線磁場強(qiáng)度的測量值和計(jì)算值;Bki為辨識程序中第k條磁滯回線上第i個點(diǎn)的磁通密度值。采用遺傳算法基于適應(yīng)度函數(shù)進(jìn)行參數(shù)尋優(yōu),結(jié)果見表2。獲得全部參數(shù)后,采用該解析式模擬所需一階回轉(zhuǎn)曲線數(shù)據(jù),進(jìn)一步可計(jì)算逆Everett 函數(shù)E(bu,bv),通過解析一階回轉(zhuǎn)曲線辨識的逆Everett 函數(shù)如圖5 所示。

    表2 B30P105 型電工鋼片多項(xiàng)式參數(shù)辨識結(jié)果Tab.2 The identified polynomial parameters of B30P105 elctrical steel sheet

    圖5 通過解析一階回轉(zhuǎn)曲線辨識的逆Everett 函數(shù)Fig.5 Inverted Everett functionidentified by analytical first order reversal curves

    通過上述逆Everett 函數(shù)辨識結(jié)果,可計(jì)算出一個周期內(nèi)不同磁通密度下電工鋼片的靜態(tài)磁滯回線。B30P105 型取向硅鋼的準(zhǔn)靜態(tài)磁滯回線實(shí)驗(yàn)結(jié)果和模型計(jì)算結(jié)果對比如圖6 所示,可以看出,本文所構(gòu)建的逆Preisach 磁滯模型具有一定的全局意義,可準(zhǔn)確模擬不同磁通密度幅值下的準(zhǔn)靜態(tài)磁滯回線,為進(jìn)一步考慮動態(tài)磁滯特性的有限元計(jì)算提供了可靠的模型基礎(chǔ)。

    圖6 準(zhǔn)靜態(tài)磁滯回線測量值與計(jì)算值對比Fig.6 Comparison between measured and calculated quasi-static hysteresis loops

    3 考慮動態(tài)磁滯特性的時域有限元方法

    3.1 動態(tài)磁場的表征方法

    Bertotti 基于磁疇理論和統(tǒng)計(jì)學(xué)原理提出損耗分離理論,將鐵磁材料單位體積的總損耗Ptot分為磁滯損耗Phys、渦流損耗Peddy和異常損耗Pexc三項(xiàng)。即

    式中,σ為材料的電導(dǎo)率;d為硅鋼片厚度;S為橫截面積;G為無量綱耦合常數(shù),G=0.135 6;V0為與材料內(nèi)部磁性單元有關(guān)的統(tǒng)計(jì)參數(shù),同磁通密度幅值Bm相關(guān),可通過不同頻率下的磁滯回線實(shí)驗(yàn)數(shù)據(jù)擬合得到。

    磁場損耗W與B、H之間的關(guān)系為

    將式(11)代入式(10)可得動態(tài)磁場表達(dá)式為

    式中,Htot、Hhys、Heddy、Hexc分別為動態(tài)磁場總磁場強(qiáng)度、靜態(tài)磁滯磁場強(qiáng)度、渦流磁場強(qiáng)度、異常磁場強(qiáng)度;δB為符號函數(shù),δB=sign(dB/dt)=±1。然而,由于較強(qiáng)的非線性特征,上述磁場難以直接在有限元分析中進(jìn)行迭代耦合計(jì)算,需進(jìn)一步研究其有限元迭代形式。

    3.2 動態(tài)磁滯特性在有限元中的直接耦合形式

    基于后向差分法,總磁場強(qiáng)度的三個分量的時域離散形式為

    式中,t為時間步;Δt為步長間隔;vh為靜態(tài)磁滯微分磁阻率,定義為

    相關(guān)麥克斯韋磁場方程為

    結(jié)合方程式(13)~式(19)可得基于矢量磁位A的控制方程為

    式中,ve為等效磁阻率,表達(dá)式為

    3.3 基于等效磁阻率的固定點(diǎn)迭代策略

    固定點(diǎn)法通過將非線性磁滯關(guān)系分為線性和非線性兩部分,通過在迭代過程內(nèi)不斷修改非線性部分以達(dá)到收斂,在解決磁滯非線性問題中具有較好的收斂性。為此,本文考慮將固定點(diǎn)迭代技術(shù)引入考慮動態(tài)磁滯特性的有限元迭代過程之中,在考慮靜態(tài)磁滯特性的基礎(chǔ)上,引入考慮動態(tài)磁場分量影響的等效磁阻率,采用式(21)中等效磁阻率代替?zhèn)鹘y(tǒng)固定點(diǎn)法的微分磁阻率vFP。應(yīng)用固定點(diǎn)技術(shù),總磁場強(qiáng)度Htot和磁通密度B的非線性關(guān)系可表示為

    式中,R為類磁化余量。

    結(jié)合麥克斯韋方程,可得控制方程為

    在二維求解域Ω中,將方程式(23)展開為離散形式為

    式中,N為離散單元的基函數(shù)。將方程式(24)改寫為矩陣形式,即

    式中,S為剛度矩陣;A為磁矢位向量;Y為電流區(qū)域系數(shù)向量;I為繞組電流向量;G為余量R的旋度矩陣。該方法中,渦流損耗等效磁場和剩余損耗等效磁場表達(dá)在類磁化余量R的計(jì)算中。通過這種方式,不需要大幅改動有限元方程,使得算法流程更簡潔,同時,將動態(tài)特性考慮在余量中,便于在以后的工作中對渦流和剩余損耗計(jì)算模型的改進(jìn)與修正。而且,在每一時刻的迭代過程,僅需計(jì)算一次等效磁阻率ve和剛度矩陣S,顯著提升了計(jì)算效率。

    為保證固定點(diǎn)算法收斂,磁阻率v應(yīng)滿足式(26)的必要條件[25,27]。

    在傳統(tǒng)固定點(diǎn)法中,所有網(wǎng)格單元采用統(tǒng)一的磁阻率,即

    式中,vmin、vmax分別為磁滯回線上斜率最小值與最大值。該方法稱為GCM(global-coefficient method)。顯然,采用式(27)計(jì)算的磁阻率總能滿足算法收斂的必要條件。然而,GCM 采用統(tǒng)一的磁阻率導(dǎo)致計(jì)算時間過長,需采用更具效率的微分磁阻率進(jìn)行迭代,微分磁阻率vFP定義為

    式中,dHtot、dB分別為磁場強(qiáng)度和磁通密度的差分。

    對于標(biāo)量磁滯模型,二維電磁場量間的關(guān)系可通過取模值降階為一維問題。如圖7 所示,藍(lán)色實(shí)線表示靜態(tài)磁滯回線,藍(lán)色虛線表示動態(tài)磁滯曲線,曲線L 為固定點(diǎn)迭代過程中動態(tài)磁滯磁場Htot與磁通密度B的軌跡。由于磁場分量Heddy、Hexc的計(jì)算是基于磁通密度B的后向差分對時間的偏導(dǎo)數(shù),因此,該曲線L 起始于t時刻靜態(tài)磁滯回線上的點(diǎn)P1,止于t+1 時刻動態(tài)磁滯回線上點(diǎn)P3。顯然,由于渦流損耗和剩余損耗的影響,曲線L 的斜率明顯不同于動態(tài)磁滯回線的斜率,在磁通密度的零點(diǎn)附近,兩者相差很大。傳統(tǒng)的固定點(diǎn)迭代法采用式(28)來計(jì)算微分磁阻率vFP,由于只考慮到動態(tài)磁滯磁場的軌跡,其計(jì)算結(jié)果偏小,往往不滿足收斂的必要條件,在迭代計(jì)算中算法極易發(fā)散。相反,在本文所提式(21)中等效磁阻率ve的計(jì)算考慮了渦流損耗和剩余損耗的等效磁場,反映了鐵心動態(tài)磁滯損耗對迭代過程的影響大小,可準(zhǔn)確估計(jì)動態(tài)B-Htot軌跡L 的斜率,從而保證迭代過程穩(wěn)定收斂。在有限元程序中ve根據(jù)前兩個時間點(diǎn)的結(jié)果進(jìn)行估算,即

    圖7 動態(tài)磁滯收斂路徑示意圖Fig.7 Schematic diagramof dynamic hysteresis convergence path

    3.4 考慮動態(tài)磁滯特性的場-路耦合有限元分析

    鑒于電機(jī)、變壓器等大多數(shù)電氣設(shè)備工作于50 Hz正弦交流電壓激勵下,為反映動態(tài)磁場和電壓源的相互作用,需在時域有限元分析中進(jìn)行“場-路”耦合。由于愛潑斯坦方圈實(shí)質(zhì)可看作一臺等比的空載變壓器,可以滿足對變壓器鐵心損耗分析的需要。為驗(yàn)證本文所提方法的有效性,本文以愛潑斯坦方圈為例構(gòu)建了二維有限元計(jì)算模型。

    在模型的計(jì)算過程中,以50 Hz 正弦交流電壓作為激勵,從電路角度分析,不考慮勵磁線圈漏磁和匝間電容,則端口電路方程為

    式中,U為激勵電壓;Uin為繞組上的感應(yīng)電動勢;r為線路電阻;I為勵磁電流。對于一階三角形單元,Uin的計(jì)算式為

    式中,Φ為繞組磁通;Nc為線圈匝數(shù);Sc為繞組截面積;Ωc為電流區(qū)域;se為單元面積;lz為二維求解系統(tǒng)縱向深度;Ae為單元節(jié)點(diǎn)磁矢位。將式(31)代入式(30)得到

    式中,U為激勵電壓向量;r為支路電阻矩陣。采用Crank-Nicolson 差分格式,將式(25)、式(32)聯(lián)立求解。

    通過求解方程式(33)得到磁矢勢向量,然后根據(jù)B=?×A計(jì)算出每個單元的磁通密度大小。如果計(jì)算結(jié)果未達(dá)到收斂條件eB=∣Bi+1-Bi∣<Bε(相對殘差閾值一般取εB= 1× 1 0-5),則通過磁滯模型計(jì)算靜態(tài)磁滯磁場Hhys、渦流磁場Heddy、剩余損耗磁場Hexc,進(jìn)而計(jì)算總磁場強(qiáng)度Htot,并通過松弛因子λ對Htot進(jìn)行修正。

    式中,i為迭代次數(shù)。若磁通密度B相鄰迭代誤差滿足收斂條件,則判斷是否達(dá)到最大時間步,若t=tmax,則停止迭代,輸出全部計(jì)算結(jié)果,否則跳轉(zhuǎn)至下一時間步。整個有限元迭代過程如圖8 所示。

    圖8 考慮動態(tài)磁滯特性的有限元分析流程Fig.8 Finite element iterative calculation process considering dynamic hysteresis characteristics

    需要說明的是:不同于靜態(tài)磁滯磁場,動態(tài)磁場強(qiáng)度和磁通密度間的非線性關(guān)系不是一開始就確定的,由于渦流磁場和剩余損耗磁場不是直接通過有限元磁場方程表達(dá),而是在固定點(diǎn)迭代過程加入,這意味著在算法迭代過程中,非線性磁場關(guān)系不斷改變,由動態(tài)磁場分量帶來的數(shù)值擾動使得磁矢位A和磁通密度B計(jì)算結(jié)果變化劇烈,造成收斂困難,尤其在磁通密度的零點(diǎn)附近,dB/dHtot較大,磁場強(qiáng)度的計(jì)算誤差引起的數(shù)值振蕩更大。因此松弛因子λ的引入是十分必要的。

    4 結(jié)果與討論

    4.1 愛潑斯坦方圈動態(tài)磁滯回線有限元模擬

    為說明本文所提方法的有效性,在計(jì)算模型中,選取了愛潑斯坦方圈臂上的三角形單元作為磁滯回線與損耗的觀測點(diǎn),其二維網(wǎng)格剖分情況如圖9 所示,三角網(wǎng)格M 為進(jìn)行實(shí)驗(yàn)數(shù)據(jù)與模型計(jì)算數(shù)據(jù)的對比的參考單元。在50 Hz 正弦條件下,圖10a、圖11a 分別給出了激勵電壓幅值為Umax=9 V 和Umax=13 V 時M 單元處磁滯回線測量值與計(jì)算值的對比??梢钥闯?,再考慮動態(tài)磁滯特性后,本文所提方法具有較高的磁滯回線模擬精度。進(jìn)一步地,圖10b、圖11b 給出了相應(yīng)電壓激勵條件下愛潑斯坦方圈勵磁電流波形的模擬結(jié)果。由于考慮了動態(tài)磁場與磁滯效應(yīng)對電流波形的反饋?zhàn)饔?,本文所提方法可以較好地模擬勵磁電流時域波形。

    圖9 愛潑斯坦方圈有限元剖分Fig.9 Finite element mesh of Epstein frame

    圖10 Umax=9 V 電壓激勵下測量值與計(jì)算值對比Fig.10 Comparison between measured and calculated values under voltage excitation Umax=9 V

    圖11 Umax=13 V 電壓激勵下測量值與計(jì)算值對比Fig.11 Comparison between measured and calculated values under voltage excitation Umax=13 V

    為考察松弛因子λ對程序收斂性的影響,控制電壓激勵分別為Umax=9 V,調(diào)節(jié)松弛因子大小,記錄程序是否收斂,以及程序收斂情況下一周期內(nèi)的平均迭代次數(shù)和程序總計(jì)算時間。結(jié)果見表3,當(dāng)λ=0.8時,有限元迭代過程不收斂,在確保收斂的情況下,采用較大的松弛因子可減少每時步的平均迭代次數(shù)和總計(jì)算時間,提升計(jì)算效率,最短總計(jì)算時間為228 s。通常,為同時滿足穩(wěn)定性和計(jì)算速度的要求,設(shè)置λ=0.5~0.7 即可。

    表3 9 V 電壓激勵時不同松弛因子對迭代次數(shù)和總迭代時間影響Tab.3 Effects of different relaxations on the iterations and computation time with voltage excitation amplitude of 9 V

    4.2 磁通密度分布與瞬時損耗特性

    為進(jìn)一步檢驗(yàn)算法的有效性,本文對愛潑斯坦方圈的磁通密度分布特性與瞬時功率損耗特性進(jìn)行了分析與計(jì)算。在t=0.005 s 時的磁通密度分布計(jì)算結(jié)果如圖12 所示??梢钥闯?,主磁路上磁通密度分布較為均勻,側(cè)壁上磁通密度幅值計(jì)算結(jié)果與實(shí)驗(yàn)測量值基本一致;而在疊片的內(nèi)角處,磁通密度較大,在外角部分磁通密度較小。

    圖12 磁通密度分布(Umax=13 V、t=0.005 s)Fig.12 Distribution of magnetic flux density (Umax=13 V、t=0.005 s)

    根據(jù)4.1 節(jié)中磁通密度B和動態(tài)磁場Htot計(jì)算結(jié)果,本文進(jìn)一步計(jì)算了Umax=9 V 和Umax=13 V 兩種激勵大小在參考單元M 處一周期內(nèi)的瞬時損耗功率tp并與測量值進(jìn)行對比,結(jié)果如圖13 所示。在t=0.00 2 s 和t=0.012 s 附近,測量與計(jì)算的瞬時損耗達(dá)到極大值;在t=0.005 s 和t=0.015 s 附近,瞬時損耗達(dá)到極小值。損耗負(fù)值來源于可逆磁化的貢獻(xiàn)??梢钥闯觯疚乃岱椒ㄔ谀M瞬時功率方面具有較高的精度,瞬時損耗的最大誤差不超過6%。

    圖13 動態(tài)磁滯瞬時損耗測量值與計(jì)算值對比Fig.13 Comparison between measured and calculated instantaneous loss when considering dynamic hysteresis

    完成一個周期有限元計(jì)算后,積分得到每一個單元磁滯回線的面積并求和即可得到鐵心區(qū)域整體鐵心總損耗PFe為

    式中,f為激勵的頻率;ne為鐵心離散單元總數(shù)。圖14 所示為不同磁通密度幅值時的鐵心總損耗計(jì)算值和測量值大小??梢钥闯觯ㄟ^提出的有限元算法計(jì)算得到的一個周期內(nèi)的鐵心損耗與實(shí)驗(yàn)測量結(jié)果基本一致。

    圖14 鐵心損耗測量值與計(jì)算值對比Fig.14 Comparison between measured and calculated power losses

    5 結(jié)論

    為解決時域有限元分析中考慮鐵心材料的動態(tài)磁滯特性出現(xiàn)難以收斂的問題,本文提出了一種基于等效磁阻率的高效穩(wěn)定改進(jìn)固定點(diǎn)迭代策略,并編寫了相應(yīng)的有限元程序,實(shí)現(xiàn)了硅鋼片的動態(tài)磁滯特性進(jìn)行模擬計(jì)算,并通過實(shí)驗(yàn)測量結(jié)果對比,驗(yàn)證了該有限元算法的有效性與準(zhǔn)確性,具體貢獻(xiàn)如下:

    1)采用解析方法生成一階回轉(zhuǎn)曲線,辨識逆Everett 函數(shù),建立逆Preisach 磁滯模型,可準(zhǔn)確模擬電工鋼片靜態(tài)磁滯特性。

    2)采用固定點(diǎn)迭代技術(shù)耦合動態(tài)磁滯特性,通過分析渦流損耗和剩余損耗等效磁場對有限元迭代過程的影響,提出基于等效磁阻率的固定點(diǎn)策略,并在有限元迭代過程引入松弛因子,通過實(shí)驗(yàn)數(shù)據(jù)與計(jì)算結(jié)果的對比分析,驗(yàn)證了本文提出的方法可實(shí)現(xiàn)低頻動態(tài)磁滯磁場的準(zhǔn)確計(jì)算,并具有良好的計(jì)算速度和穩(wěn)定性。

    3)通過場-路耦合有限元算法與動態(tài)磁滯模型結(jié)合,可充分考慮鐵心總損耗對電路和磁場的反饋效應(yīng),實(shí)現(xiàn)瞬時損耗的準(zhǔn)確計(jì)算,最大計(jì)算誤差不超過6%,可滿足工程計(jì)算要求。

    猜你喜歡
    磁滯回線磁阻磁通
    軸向磁通電勵磁雙凸極電機(jī)及容錯運(yùn)行控制策略
    基于MATLAB處理大學(xué)物理實(shí)驗(yàn)數(shù)據(jù)探究
    永磁磁阻電動機(jī)的研究
    磁場強(qiáng)度波形畸變對交流磁滯回線形狀的影響
    基于LabVIEW的微型磁通門磁強(qiáng)計(jì)測試系統(tǒng)搭建
    基于磁通門原理的零磁通交直流電流傳感器
    高頻脈沖激勵下磁滯回線動態(tài)測量裝置的設(shè)計(jì)及分析
    巨磁阻電渦流傳感器設(shè)計(jì)
    基于FPGA的數(shù)字磁通計(jì)設(shè)計(jì)
    電測與儀表(2015年3期)2015-04-09 11:37:52
    四相開關(guān)磁阻電機(jī)的四電平DITC調(diào)速系統(tǒng)
    看免费成人av毛片| 日韩电影二区| 久久精品亚洲av国产电影网| 国产成人精品婷婷| av在线播放精品| 97人妻天天添夜夜摸| 亚洲欧美色中文字幕在线| 久久久久精品人妻al黑| 中文欧美无线码| 国产精品欧美亚洲77777| 丰满少妇做爰视频| 如何舔出高潮| 90打野战视频偷拍视频| 少妇人妻久久综合中文| 国产av国产精品国产| 黄色视频在线播放观看不卡| 考比视频在线观看| 免费人妻精品一区二区三区视频| 国产一区二区 视频在线| 丝袜脚勾引网站| 亚洲欧美色中文字幕在线| 国产av码专区亚洲av| 久久精品国产自在天天线| 人人澡人人妻人| 国产乱来视频区| 精品人妻在线不人妻| 免费观看av网站的网址| av又黄又爽大尺度在线免费看| 一级毛片我不卡| 亚洲图色成人| 久久久国产精品麻豆| 欧美bdsm另类| 丰满少妇做爰视频| 国产乱人偷精品视频| 美女高潮到喷水免费观看| 日本午夜av视频| 国产男女内射视频| 秋霞在线观看毛片| 精品人妻熟女毛片av久久网站| 最新的欧美精品一区二区| 日韩中文字幕视频在线看片| 777久久人妻少妇嫩草av网站| 久久精品国产a三级三级三级| 精品一品国产午夜福利视频| 超色免费av| 国产亚洲最大av| 九色亚洲精品在线播放| 热99久久久久精品小说推荐| 一级,二级,三级黄色视频| av网站在线播放免费| 国产精品熟女久久久久浪| 欧美精品高潮呻吟av久久| 欧美日韩一区二区视频在线观看视频在线| 自拍欧美九色日韩亚洲蝌蚪91| 午夜老司机福利剧场| 精品亚洲乱码少妇综合久久| 你懂的网址亚洲精品在线观看| 色吧在线观看| 久久久久精品性色| 日韩在线高清观看一区二区三区| 成人黄色视频免费在线看| 黄色怎么调成土黄色| 亚洲成人一二三区av| 国产精品二区激情视频| 日韩伦理黄色片| 卡戴珊不雅视频在线播放| 中文字幕制服av| 大香蕉久久网| 久久久精品国产亚洲av高清涩受| 色婷婷久久久亚洲欧美| 成年女人在线观看亚洲视频| 免费人妻精品一区二区三区视频| 久久青草综合色| 18禁动态无遮挡网站| 亚洲精品美女久久av网站| 成人黄色视频免费在线看| 国精品久久久久久国模美| 久久人人97超碰香蕉20202| 国产免费现黄频在线看| av网站在线播放免费| 777久久人妻少妇嫩草av网站| 午夜福利一区二区在线看| a级毛片在线看网站| √禁漫天堂资源中文www| 日韩 亚洲 欧美在线| 80岁老熟妇乱子伦牲交| 欧美日韩精品成人综合77777| 久久97久久精品| 日韩制服丝袜自拍偷拍| 男女免费视频国产| 伦理电影免费视频| 中文字幕亚洲精品专区| 熟女电影av网| 丰满乱子伦码专区| 丰满饥渴人妻一区二区三| 人妻 亚洲 视频| 男女边摸边吃奶| av一本久久久久| 午夜激情av网站| a级片在线免费高清观看视频| 18在线观看网站| 精品第一国产精品| 九九爱精品视频在线观看| 国产白丝娇喘喷水9色精品| 日韩中字成人| 99久国产av精品国产电影| 久久精品熟女亚洲av麻豆精品| 人成视频在线观看免费观看| 搡老乐熟女国产| 啦啦啦啦在线视频资源| 国产免费现黄频在线看| 有码 亚洲区| 又黄又粗又硬又大视频| av视频免费观看在线观看| 久久99精品国语久久久| 亚洲欧美一区二区三区久久| av视频免费观看在线观看| 两性夫妻黄色片| 最近中文字幕高清免费大全6| www.自偷自拍.com| 美女大奶头黄色视频| av天堂久久9| 晚上一个人看的免费电影| 午夜福利在线免费观看网站| 国产精品欧美亚洲77777| 久久精品久久久久久久性| 国产欧美亚洲国产| av片东京热男人的天堂| 麻豆乱淫一区二区| 不卡视频在线观看欧美| 桃花免费在线播放| 免费黄色在线免费观看| 欧美黄色片欧美黄色片| 人妻 亚洲 视频| 两性夫妻黄色片| 午夜日韩欧美国产| 青春草视频在线免费观看| 一级,二级,三级黄色视频| 观看美女的网站| 妹子高潮喷水视频| 国产精品.久久久| 女人精品久久久久毛片| 青春草国产在线视频| 九草在线视频观看| 色哟哟·www| 国产乱来视频区| www.av在线官网国产| 免费播放大片免费观看视频在线观看| 久久久久精品性色| 国产探花极品一区二区| 日韩一卡2卡3卡4卡2021年| 最近手机中文字幕大全| 亚洲男人天堂网一区| 国产国语露脸激情在线看| 午夜福利影视在线免费观看| 亚洲精华国产精华液的使用体验| 国产精品国产三级专区第一集| 日本爱情动作片www.在线观看| 国产深夜福利视频在线观看| 午夜日本视频在线| av电影中文网址| 一区福利在线观看| 日本vs欧美在线观看视频| 黄色视频在线播放观看不卡| 久久久久久久久久久免费av| av天堂久久9| 伊人久久大香线蕉亚洲五| 免费观看av网站的网址| av在线播放精品| 午夜影院在线不卡| 日本av免费视频播放| 成人18禁高潮啪啪吃奶动态图| 9191精品国产免费久久| 国产精品久久久久久av不卡| av女优亚洲男人天堂| av线在线观看网站| 久久久久久免费高清国产稀缺| 国产成人aa在线观看| 国产精品久久久久久久久免| 欧美日韩亚洲高清精品| 啦啦啦啦在线视频资源| 国产亚洲一区二区精品| 午夜免费观看性视频| 寂寞人妻少妇视频99o| 欧美日韩一级在线毛片| 中国国产av一级| 最近的中文字幕免费完整| 日韩欧美一区视频在线观看| 久久精品国产亚洲av高清一级| 日日撸夜夜添| 国产在视频线精品| 精品亚洲成国产av| 叶爱在线成人免费视频播放| 欧美日韩精品成人综合77777| 中文字幕av电影在线播放| 久久久久久久大尺度免费视频| 十八禁网站网址无遮挡| 两性夫妻黄色片| av国产精品久久久久影院| 久久人人爽人人片av| 欧美精品高潮呻吟av久久| 亚洲第一区二区三区不卡| 日韩中文字幕视频在线看片| 国产精品嫩草影院av在线观看| 国产日韩欧美在线精品| 久久人妻熟女aⅴ| 亚洲精品久久久久久婷婷小说| 久久久久视频综合| 一区二区三区乱码不卡18| 91成人精品电影| 国产精品香港三级国产av潘金莲 | 最近中文字幕2019免费版| 国产成人91sexporn| 精品久久久精品久久久| 一个人免费看片子| 男女免费视频国产| 欧美中文综合在线视频| a级毛片黄视频| 美女中出高潮动态图| 久久免费观看电影| av天堂久久9| 成人二区视频| 美女高潮到喷水免费观看| 日韩不卡一区二区三区视频在线| 99九九在线精品视频| 国产亚洲av片在线观看秒播厂| 熟女电影av网| 日韩三级伦理在线观看| 老司机影院成人| www.av在线官网国产| 日本黄色日本黄色录像| 国产国语露脸激情在线看| 亚洲成人手机| 天天躁狠狠躁夜夜躁狠狠躁| 午夜福利在线观看免费完整高清在| 少妇的丰满在线观看| 老司机影院成人| 国产av一区二区精品久久| 亚洲美女黄色视频免费看| 黑人欧美特级aaaaaa片| 日韩一区二区三区影片| 久久毛片免费看一区二区三区| 亚洲av国产av综合av卡| 国产免费视频播放在线视频| 亚洲av在线观看美女高潮| 亚洲在久久综合| 国产国语露脸激情在线看| 欧美人与性动交α欧美软件| 久久国产精品男人的天堂亚洲| 国产成人精品无人区| 多毛熟女@视频| 人妻少妇偷人精品九色| 日韩一区二区视频免费看| 日本免费在线观看一区| 亚洲av欧美aⅴ国产| 有码 亚洲区| 亚洲国产色片| 考比视频在线观看| 男女国产视频网站| 热re99久久国产66热| 久久亚洲国产成人精品v| 免费人妻精品一区二区三区视频| 亚洲成人手机| 精品一区二区三卡| 婷婷色麻豆天堂久久| av在线播放精品| 中文字幕制服av| 你懂的网址亚洲精品在线观看| 爱豆传媒免费全集在线观看| 中文字幕精品免费在线观看视频| 日本欧美视频一区| 国产精品无大码| 一本大道久久a久久精品| 1024香蕉在线观看| 日本wwww免费看| 日韩 亚洲 欧美在线| 性高湖久久久久久久久免费观看| 国产亚洲午夜精品一区二区久久| 建设人人有责人人尽责人人享有的| 久久久久国产网址| 精品久久久久久电影网| 黄色视频在线播放观看不卡| 国产一区有黄有色的免费视频| 亚洲激情五月婷婷啪啪| 国产精品久久久久久av不卡| 在线观看一区二区三区激情| 日韩一本色道免费dvd| 国产国语露脸激情在线看| 99久久精品国产国产毛片| 精品一品国产午夜福利视频| 69精品国产乱码久久久| 咕卡用的链子| 久久毛片免费看一区二区三区| 午夜av观看不卡| 精品少妇久久久久久888优播| 韩国精品一区二区三区| 高清黄色对白视频在线免费看| 国产一区二区 视频在线| 欧美日韩亚洲高清精品| 中文字幕人妻丝袜制服| av一本久久久久| 国产日韩欧美视频二区| 午夜免费男女啪啪视频观看| 国产精品一国产av| 精品卡一卡二卡四卡免费| 日本91视频免费播放| 天堂中文最新版在线下载| 中文字幕最新亚洲高清| 三上悠亚av全集在线观看| 美女高潮到喷水免费观看| 国产97色在线日韩免费| av女优亚洲男人天堂| 精品一区二区三卡| 一级,二级,三级黄色视频| 欧美国产精品一级二级三级| 麻豆精品久久久久久蜜桃| 久久久精品免费免费高清| 欧美bdsm另类| 欧美日韩视频高清一区二区三区二| 涩涩av久久男人的天堂| 亚洲av日韩在线播放| 成人午夜精彩视频在线观看| 18禁动态无遮挡网站| 香蕉丝袜av| 有码 亚洲区| 女的被弄到高潮叫床怎么办| 老司机影院毛片| 黑人猛操日本美女一级片| 如何舔出高潮| 亚洲欧美一区二区三区黑人 | 黄色视频在线播放观看不卡| 老汉色av国产亚洲站长工具| 色吧在线观看| 久久久久久免费高清国产稀缺| 亚洲成人av在线免费| 波多野结衣一区麻豆| a级毛片黄视频| 天堂中文最新版在线下载| 日韩欧美精品免费久久| 黑丝袜美女国产一区| 久久97久久精品| 1024香蕉在线观看| 亚洲国产精品一区三区| 久久久欧美国产精品| 国产精品欧美亚洲77777| 精品视频人人做人人爽| 2018国产大陆天天弄谢| 亚洲国产看品久久| 国产精品免费大片| 精品人妻在线不人妻| 国产精品免费大片| 精品国产一区二区三区四区第35| 自拍欧美九色日韩亚洲蝌蚪91| 欧美av亚洲av综合av国产av | 色播在线永久视频| 亚洲成人av在线免费| 午夜日韩欧美国产| 国产av精品麻豆| 久久精品国产综合久久久| 成人亚洲精品一区在线观看| 亚洲成人手机| 精品国产一区二区三区久久久樱花| 王馨瑶露胸无遮挡在线观看| 色94色欧美一区二区| 亚洲色图综合在线观看| 成年人免费黄色播放视频| 久久99热这里只频精品6学生| 中文字幕人妻丝袜一区二区 | 色婷婷久久久亚洲欧美| 亚洲av欧美aⅴ国产| 色婷婷久久久亚洲欧美| 国产精品 国内视频| 国产精品99久久99久久久不卡 | 大码成人一级视频| 午夜福利视频在线观看免费| av免费在线看不卡| 亚洲第一av免费看| 国产精品亚洲av一区麻豆 | 一级毛片黄色毛片免费观看视频| 婷婷色综合www| 免费女性裸体啪啪无遮挡网站| 亚洲综合精品二区| 亚洲欧美精品自产自拍| 水蜜桃什么品种好| 只有这里有精品99| 黄色 视频免费看| av免费观看日本| 欧美亚洲日本最大视频资源| 日韩中字成人| 女人高潮潮喷娇喘18禁视频| 成年女人毛片免费观看观看9 | 亚洲一区中文字幕在线| 国产福利在线免费观看视频| 久久久久国产精品人妻一区二区| 亚洲激情五月婷婷啪啪| 日韩一区二区三区影片| 美女主播在线视频| 久久精品国产鲁丝片午夜精品| 日韩制服骚丝袜av| 巨乳人妻的诱惑在线观看| 啦啦啦在线观看免费高清www| 青春草视频在线免费观看| 色哟哟·www| 热99久久久久精品小说推荐| 校园人妻丝袜中文字幕| 亚洲成色77777| 亚洲婷婷狠狠爱综合网| 国产欧美亚洲国产| 亚洲欧美成人综合另类久久久| 免费观看在线日韩| 午夜影院在线不卡| 七月丁香在线播放| 成人黄色视频免费在线看| av免费在线看不卡| 午夜久久久在线观看| 日韩精品有码人妻一区| 国产精品秋霞免费鲁丝片| 交换朋友夫妻互换小说| 亚洲国产最新在线播放| 国产麻豆69| 欧美bdsm另类| 一边亲一边摸免费视频| 精品国产乱码久久久久久小说| 精品久久久久久电影网| 视频在线观看一区二区三区| 亚洲一码二码三码区别大吗| 久久女婷五月综合色啪小说| 午夜福利,免费看| 日本欧美视频一区| 天天躁狠狠躁夜夜躁狠狠躁| 狠狠精品人妻久久久久久综合| 国产精品久久久久久久久免| 欧美精品高潮呻吟av久久| 天天躁狠狠躁夜夜躁狠狠躁| 国产在线一区二区三区精| 国产成人欧美| 欧美 日韩 精品 国产| 国产淫语在线视频| 最近中文字幕2019免费版| 国产精品一区二区在线观看99| 麻豆乱淫一区二区| 成人影院久久| 老鸭窝网址在线观看| 视频在线观看一区二区三区| 在线观看免费日韩欧美大片| 99久久综合免费| 美女中出高潮动态图| 亚洲一码二码三码区别大吗| 亚洲精品,欧美精品| 亚洲国产欧美在线一区| 日韩制服骚丝袜av| 在线看a的网站| 热99久久久久精品小说推荐| 日日啪夜夜爽| 国产av一区二区精品久久| 日本免费在线观看一区| 中文字幕最新亚洲高清| 久久精品久久久久久久性| 岛国毛片在线播放| 国产av码专区亚洲av| 精品人妻偷拍中文字幕| 久久免费观看电影| 国产 一区精品| 日韩大片免费观看网站| 日本猛色少妇xxxxx猛交久久| 久久久久久久久久久免费av| 99精国产麻豆久久婷婷| 中文字幕人妻熟女乱码| 大话2 男鬼变身卡| 精品人妻在线不人妻| kizo精华| 99re6热这里在线精品视频| 国产精品无大码| 精品久久久精品久久久| 少妇熟女欧美另类| 伊人久久大香线蕉亚洲五| 蜜桃国产av成人99| 女的被弄到高潮叫床怎么办| 日韩电影二区| 美女高潮到喷水免费观看| 免费女性裸体啪啪无遮挡网站| 男的添女的下面高潮视频| 只有这里有精品99| 国产欧美日韩综合在线一区二区| 又粗又硬又长又爽又黄的视频| 免费看av在线观看网站| 2022亚洲国产成人精品| 国产 精品1| 日本-黄色视频高清免费观看| 久久综合国产亚洲精品| 午夜影院在线不卡| 亚洲美女黄色视频免费看| 久久精品夜色国产| 精品国产一区二区三区久久久樱花| 亚洲精品国产色婷婷电影| 中文字幕亚洲精品专区| 最近最新中文字幕大全免费视频 | 成人亚洲精品一区在线观看| 香蕉国产在线看| a级毛片黄视频| av在线老鸭窝| 十八禁网站网址无遮挡| 精品午夜福利在线看| 国产日韩欧美亚洲二区| 中文欧美无线码| 亚洲综合精品二区| 亚洲精品美女久久av网站| 日韩中文字幕视频在线看片| 久久精品久久久久久久性| 亚洲欧美一区二区三区黑人 | 国产成人欧美| 大话2 男鬼变身卡| av国产久精品久网站免费入址| 免费看不卡的av| 免费不卡的大黄色大毛片视频在线观看| 国产黄色视频一区二区在线观看| 国产欧美亚洲国产| 亚洲内射少妇av| 久久久久久久国产电影| 欧美日韩一区二区视频在线观看视频在线| 在线观看三级黄色| 国产爽快片一区二区三区| 麻豆av在线久日| 国精品久久久久久国模美| 久久国产精品男人的天堂亚洲| 国产成人91sexporn| 97精品久久久久久久久久精品| 久久精品国产亚洲av涩爱| 成人国语在线视频| 精品国产超薄肉色丝袜足j| 国产成人精品久久久久久| 色视频在线一区二区三区| 亚洲久久久国产精品| 日韩中文字幕视频在线看片| 亚洲av综合色区一区| 亚洲一码二码三码区别大吗| 中文字幕av电影在线播放| 午夜福利一区二区在线看| 日日撸夜夜添| 飞空精品影院首页| 婷婷色麻豆天堂久久| 9191精品国产免费久久| 亚洲第一av免费看| 亚洲美女黄色视频免费看| 天堂俺去俺来也www色官网| 国产精品国产三级专区第一集| 超碰成人久久| 看十八女毛片水多多多| 久久久久国产一级毛片高清牌| 欧美国产精品一级二级三级| 国产精品麻豆人妻色哟哟久久| 少妇被粗大猛烈的视频| 纯流量卡能插随身wifi吗| 极品少妇高潮喷水抽搐| 欧美精品一区二区免费开放| www日本在线高清视频| 美女脱内裤让男人舔精品视频| 国产1区2区3区精品| 午夜免费观看性视频| 久久久久国产一级毛片高清牌| 三上悠亚av全集在线观看| 国产成人精品久久二区二区91 | 精品国产乱码久久久久久男人| 王馨瑶露胸无遮挡在线观看| 色网站视频免费| 亚洲国产欧美网| 久久国产精品大桥未久av| av国产精品久久久久影院| 18禁国产床啪视频网站| 国产亚洲最大av| 少妇被粗大的猛进出69影院| 亚洲精品国产一区二区精华液| 久久久精品免费免费高清| 熟女少妇亚洲综合色aaa.| 成年女人在线观看亚洲视频| 欧美精品国产亚洲| 男人添女人高潮全过程视频| 精品久久蜜臀av无| 女人高潮潮喷娇喘18禁视频| 亚洲av.av天堂| 国产乱人偷精品视频| 欧美bdsm另类| 欧美 日韩 精品 国产| 这个男人来自地球电影免费观看 | 日日爽夜夜爽网站| 国产一区有黄有色的免费视频| 哪个播放器可以免费观看大片| 性色av一级| 2022亚洲国产成人精品| 一二三四在线观看免费中文在| 一级,二级,三级黄色视频| 中文字幕精品免费在线观看视频| 国产女主播在线喷水免费视频网站| 五月天丁香电影| 99国产综合亚洲精品| 免费在线观看视频国产中文字幕亚洲 | 国产成人a∨麻豆精品| 久久久久久久久久人人人人人人| 王馨瑶露胸无遮挡在线观看| 亚洲av电影在线进入| 99久久综合免费| 色94色欧美一区二区| 最近中文字幕高清免费大全6| 亚洲在久久综合| 最近中文字幕高清免费大全6| 色94色欧美一区二区| 99久久综合免费| 在线天堂最新版资源| 热99国产精品久久久久久7| 欧美精品人与动牲交sv欧美| 免费黄色在线免费观看| 只有这里有精品99| 如何舔出高潮| 亚洲欧美成人综合另类久久久| 91成人精品电影|