別全解析:從速度梯度張量到Matlab實(shí)現(xiàn)與Tecplot可視化)
簡(jiǎn)介這套Liutex第三代渦識(shí)別方法的Matlab-Tecplot實(shí)現(xiàn)資料面向流體力學(xué)、CFD后處理及渦識(shí)別方向的研究者和工程師幫助理解Liutex理論并快速上手代碼實(shí)現(xiàn)。資源共6個(gè)文件以3個(gè)Matlab腳本負(fù)責(zé)核心計(jì)算、矩陣處理與Liutex導(dǎo)出為主另含PDF文獻(xiàn)、PPT方法介紹及eqn公式文件壓縮包僅6.67MB小巧但知識(shí)點(diǎn)密集。內(nèi)容從渦定義講到實(shí)際計(jì)算流程讀者可對(duì)照文檔理解Liutex的張量運(yùn)算再利用現(xiàn)成腳本完成數(shù)據(jù)輸出與Tecplot可視化節(jié)省從零寫(xiě)代碼和調(diào)試的時(shí)間。附帶的PDF與PPT梳理了Liutex的數(shù)學(xué)背景和渦識(shí)別發(fā)展脈絡(luò)eqn公式文件便于對(duì)照關(guān)鍵推導(dǎo)使這套包不僅是代碼工具也是很好的理論入門(mén)材料。目前已有3582人學(xué)習(xí)適合需要系統(tǒng)掌握第三代渦識(shí)別方法、做科研復(fù)現(xiàn)或課程教學(xué)的高年級(jí)本科生、研究生與科研人員資源免費(fèi)開(kāi)放無(wú)需積分即可下載進(jìn)一步降低了學(xué)習(xí)門(mén)檻。1. Liutex 不是又一種渦判據(jù)先搞清楚它解決什么問(wèn)題做流體后處理時(shí)最頭疼的從來(lái)都不是算得慢而是算出來(lái)的東西說(shuō)不清。用渦量顯示旋渦剪切層里全是假陽(yáng)性用 Q 判據(jù)畫(huà)等值面閾值要靠肉眼調(diào)半天而且不同工況之間沒(méi)法橫向比。這些傳統(tǒng)判據(jù)都卡在同一個(gè)問(wèn)題上它們描述的“渦”不是直觀意義上的旋轉(zhuǎn)而是數(shù)學(xué)上的某個(gè)張量特性。Liutex 方法把這個(gè)問(wèn)題換了個(gè)問(wèn)法——不再問(wèn)“這里像不像有渦”而是直接算“這個(gè)點(diǎn)的局部旋轉(zhuǎn)軸指向哪個(gè)方向、旋轉(zhuǎn)強(qiáng)度是多少”輸出的是一個(gè)三維向量場(chǎng)后續(xù)無(wú)論是切片的渦核線還是三維等值面物理意義都更明確。這份資源給的是 Liutex 理論在 Matlab 端從速度場(chǎng)到 Liutex 向量的完整實(shí)現(xiàn)以及把結(jié)果導(dǎo)入 Tecplot 做可視化的配套腳本適合正在寫(xiě)論文需要一種可靠渦識(shí)別方法、或者做 CFX/開(kāi)源求解器后處理想換個(gè)判據(jù)驗(yàn)證結(jié)論的人。2. Liutex 的數(shù)學(xué)骨架旋轉(zhuǎn)軸、旋轉(zhuǎn)強(qiáng)度與特征分解2.1 為什么傳統(tǒng)判據(jù)說(shuō)不清“旋轉(zhuǎn)軸”方向先對(duì)齊一個(gè)基本概念判斷一個(gè)局部流動(dòng)是否在旋轉(zhuǎn)不能只看速度本身而是看速度的空間變化率。這個(gè)變化率由速度梯度張量描述九個(gè)分量排成一個(gè) 3×3 矩陣?V [du/dx du/dy du/dz] [dv/dx dv/dy dv/dz] [dw/dx dw/dy dw/dz]它可以分解成對(duì)稱部分應(yīng)變率張量 S和反對(duì)稱部分旋轉(zhuǎn)張量 Ω。Q 判據(jù)看的是 ‖Ω‖2 - ‖S‖2 的正負(fù)λ2 判據(jù)看的是 S2 Ω2 的特征值這兩種方法本質(zhì)上都在做同一件事拿“剪切/拉伸”和“旋轉(zhuǎn)”比誰(shuí)更占優(yōu)而不是把旋轉(zhuǎn)單獨(dú)提取出來(lái)。所以同樣的流動(dòng)閾值取 0.5 和取 1 畫(huà)出來(lái)的拓?fù)浣Y(jié)構(gòu)完全不一樣這就是玄學(xué)來(lái)源。Liutex 的關(guān)鍵變化是放棄這種比例比較直接對(duì) ?V 做特征分解找出張量中代表真實(shí)旋轉(zhuǎn)的分量。數(shù)學(xué)上有三種情況一種是三個(gè)實(shí)特征值代表該點(diǎn)沒(méi)有局部旋轉(zhuǎn)純變形另一種是一個(gè)實(shí)特征值加一對(duì)共軛復(fù)特征值這對(duì)復(fù)特征值的存在就意味著該點(diǎn)存在局部旋轉(zhuǎn)。Liutex 認(rèn)為復(fù)特征值的虛部決定了旋轉(zhuǎn)的強(qiáng)度而實(shí)特征值對(duì)應(yīng)的特征向量方向就是局部旋轉(zhuǎn)軸的方向。這個(gè)思路比在一堆張量范數(shù)里找閾值更貼近物理直覺(jué)。2.2 從特征值到旋轉(zhuǎn)強(qiáng)度的標(biāo)定關(guān)系具體算的時(shí)候會(huì)得到如下關(guān)系。設(shè) ?V 的特征值為 λr 和 λcr ± iλci其中 λr 是實(shí)數(shù)λcr 和 λci 是共軛復(fù)特征值的實(shí)部和虛部。Liutex 向量定義為R R * r其中 r 是實(shí)特征值 λr 對(duì)應(yīng)的單位特征向量也就是旋轉(zhuǎn)軸方向而旋轉(zhuǎn)強(qiáng)度 R 的大小取 2|λci|。為什么是 2 倍因?yàn)樵诩冃D(zhuǎn)流動(dòng)里速度梯度張量的復(fù)特征值虛部正好等于角速度的一半乘 2 之后Liutex 強(qiáng)度的量綱才能和渦量 ω 對(duì)齊單位都是 1/s這樣方便和后處理里的其他物理量橫向比較。這里有一個(gè)值得注意的細(xì)節(jié)當(dāng)三個(gè)特征值全是實(shí)數(shù)時(shí)表示該點(diǎn)沒(méi)有局部旋轉(zhuǎn)直接把 Liutex 設(shè)為零向量。實(shí)際計(jì)算時(shí)還會(huì)遇到數(shù)值誤差導(dǎo)致虛部為微小值的偽復(fù)特征根需要設(shè)置一個(gè)虛部閾值小于該閾值的按實(shí)特征值處理。我在實(shí)現(xiàn)時(shí)用的閾值是 1e-10但這取決于數(shù)據(jù)精度如果是單精度浮點(diǎn)場(chǎng)可能需要放大到 1e-6 級(jí)別。對(duì)比起來(lái)Q 判據(jù)需要取閾值λ2 判據(jù)也需要取閾值而 Liutex 的旋轉(zhuǎn)強(qiáng)度是逐點(diǎn)計(jì)算出來(lái)的不需要人為設(shè)閾值。這并不代表 Liutex 在可視化階段不需要選等值面但至少“旋轉(zhuǎn)強(qiáng)度”本身是客觀的不是通過(guò)一個(gè)方程左右兩邊誰(shuí)大的不等式判出來(lái)的。這也是這個(gè)資源里最值錢(qián)的部分它把計(jì)算核心從“猜閾值”變成了“解特征值”。3. Matlab 端實(shí)現(xiàn)從速度場(chǎng)到 Liutex 向量的完整流程3.1 數(shù)據(jù)準(zhǔn)備與速度梯度計(jì)算這套實(shí)現(xiàn)接收的是結(jié)構(gòu)化網(wǎng)格上的速度場(chǎng)三個(gè)速度分量分別存成三維數(shù)組 U、V、W維度順序建議是 (nx, ny, nz)對(duì)應(yīng)物理坐標(biāo) (x, y, z)。如果手里是 CFD 軟件導(dǎo)出的非結(jié)構(gòu)化網(wǎng)格需要先插值到矩形網(wǎng)格上再喂給腳本這一點(diǎn)后面避坑章還會(huì)專(zhuān)門(mén)說(shuō)。速度梯度的計(jì)算是整套流程的精度基礎(chǔ)。內(nèi)部點(diǎn)用中心差分邊界用單側(cè)差分。我這里給出一個(gè)可用的實(shí)現(xiàn)片段function [dudx, dudy, dudz] calc_gradient(U, dx, dy, dz) [nx, ny, nz] size(U); dudx zeros(nx, ny, nz); dudy zeros(nx, ny, nz); dudz zeros(nx, ny, nz); % 內(nèi)部中心差分 dudx(2:end-1, :, :) (U(3:end, :, :) - U(1:end-2, :, :)) / (2*dx); dudy(:, 2:end-1, :) (U(:, 3:end, :) - U(:, 1:end-2, :)) / (2*dy); dudz(:, :, 2:end-1) (U(:, :, 3:end) - U(:, :, 1:end-2)) / (2*dz); % x 方向邊界一階單側(cè)差分 dudx(1, :, :) (-3*U(1,:,:) 4*U(2,:,:) - U(3,:,:)) / (2*dx); dudx(end, :, :) (3*U(end,:,:) - 4*U(end-1,:,:) U(end-2,:,:)) / (2*dx); % y、z 邊界同理此處省略 enddy 和 dz 方向的邊界處理與 x 方向完全一致不再重復(fù)展開(kāi)。中心差分的截?cái)嗾`差是二階對(duì)大多數(shù)仿真數(shù)據(jù)夠用如果網(wǎng)格長(zhǎng)寬比很大比如邊界層里法向網(wǎng)格極密建議對(duì)法向梯度改用加權(quán)差分否則切向梯度會(huì)被噪聲放大。這個(gè)函數(shù)返回的是 U 分量三個(gè)方向的偏導(dǎo)。為了組裝完整的梯度張量需要分別對(duì) V 和 W 分量也調(diào)用一次然后把九個(gè)偏導(dǎo)按矩陣排列。檢查梯度的正確性有一個(gè)土辦法取一個(gè)已知旋轉(zhuǎn)速度場(chǎng)后面驗(yàn)證章會(huì)給算出來(lái)的渦量 ω ?×V 應(yīng)該和分析值吻合誤差在 1% 以內(nèi)再繼續(xù)不然后面全白做。3.2 逐點(diǎn)特征分解與 Liutex 向量組裝得到九個(gè)梯度分量之后在每個(gè)網(wǎng)格點(diǎn)上組裝梯度張量然后做特征分解。這里我給出核心循環(huán)的完整邏輯% 預(yù)分配輸出數(shù)組 LiutexX zeros(nx, ny, nz); LiutexY zeros(nx, ny, nz); LiutexZ zeros(nx, ny, nz); omega_x dWdy - dVdz; % 渦量分量用于方向矯正 omega_y dUdz - dWdx; omega_z dVdx - dUdy; for k 2:nz-1 for j 2:ny-1 for i 2:nx-1 % 組裝速度梯度張量 GradV [ dUdX(i,j,k), dUdY(i,j,k), dUdZ(i,j,k); dVdX(i,j,k), dVdY(i,j,k), dVdZ(i,j,k); dWdX(i,j,k), dWdY(i,j,k), dWdZ(i,j,k) ]; % 特征值分解 [eigVec, eigVal] eig(GradV); lambda diag(eigVal); % 找實(shí)特征值虛部絕對(duì)值小于閾值 isReal abs(imag(lambda)) 1e-10; if sum(isReal) 1 realIdx find(isReal, 1); r real(eigVec(:, realIdx)); r r / norm(r); % 歸一化旋轉(zhuǎn)軸方向 % 找到一對(duì)共軛復(fù)特征值取虛部 complexIdx find(~isReal); R_mag 2 * abs(imag(lambda(complexIdx(1)))); % 旋轉(zhuǎn)軸方向與渦量方向?qū)R保證右手定則 omega [omega_x(i,j,k), omega_y(i,j,k), omega_z(i,j,k)]; if dot(r, omega) 0 r -r; end LiutexX(i,j,k) R_mag * r(1); LiutexY(i,j,k) R_mag * r(2); LiutexZ(i,j,k) R_mag * r(3); end % 三個(gè)實(shí)特征值無(wú)旋轉(zhuǎn)保持零向量 end end end特征分解用的是 Matlab 自帶的 eig 函數(shù)對(duì)小規(guī)模網(wǎng)格速度足夠快但要注意如果網(wǎng)格是千萬(wàn)級(jí)體量逐點(diǎn)用 eig 在循環(huán)里跑會(huì)非常慢一個(gè)常見(jiàn)的優(yōu)化是改用 LAPACK 的 dgeev 接口或者把張量批量展開(kāi)后向量化速度能提一個(gè)數(shù)量級(jí)。資源里給的是逐點(diǎn)版本優(yōu)先保證清晰度性能不夠時(shí)再改。旋轉(zhuǎn)軸方向與渦量方向?qū)R這一步非常關(guān)鍵。eig 返回的特征向量方向在數(shù)值上是任意的正負(fù)都可能同一物理點(diǎn)在不同滿足度下可能得到完全相反的旋轉(zhuǎn)軸方向。用渦量方向做矯正能確保相鄰網(wǎng)格點(diǎn)上旋轉(zhuǎn)軸方向連續(xù)變化畫(huà)出來(lái)不會(huì)出現(xiàn)軸方向突變導(dǎo)致的渦核線斷裂。關(guān)于旋轉(zhuǎn)強(qiáng)度的符號(hào)Liutex 方法中 R 本身定義為非負(fù)方向由旋轉(zhuǎn)軸向量承擔(dān)所以這里取了虛部的絕對(duì)值。如果你習(xí)慣標(biāo)量顯示旋轉(zhuǎn)強(qiáng)度可以只取 R_mag而向量方向用 r 另外存儲(chǔ)。這兩種用法在資源里分別對(duì)應(yīng)了 Tecplot 的“云圖顯示”和“矢量顯示”兩種模式。3.3 邊界處理與相鄰點(diǎn)連續(xù)性檢查上面的循環(huán)只處理了內(nèi)部點(diǎn)邊界上的 Liutex 值沿用零或者用單側(cè)差分補(bǔ)算這一點(diǎn)要看實(shí)際需求。如果渦核不靠近邊界直接設(shè)零不影響顯示但如果關(guān)注的是管道近壁區(qū)域就必須在邊界上用單側(cè)梯度補(bǔ)算 Liutex。補(bǔ)算時(shí)需要注意邊界上特征分解結(jié)果往往很不穩(wěn)定因?yàn)閱蝹?cè)差分精度只有一階張量的實(shí)特征值會(huì)頻繁出現(xiàn)微小虛部虛部閾值在這里要適當(dāng)放松到 1e-6 級(jí)別。做完逐點(diǎn)計(jì)算我建議順手檢查一次旋轉(zhuǎn)軸方向的連續(xù)性。方法很簡(jiǎn)單取渦核中心線附近的一排點(diǎn)把 Liutex 向量的三個(gè)分量打印出來(lái)觀察其方向是否呈平滑變化。如果出現(xiàn)間隔點(diǎn)方向翻轉(zhuǎn)多半是渦量方向矯正在零渦量點(diǎn)附近失效可以再加一個(gè)約束當(dāng) |ω| 小于某個(gè)閾值比如最大渦量的 1%時(shí)沿用上一個(gè)網(wǎng)格點(diǎn)的旋轉(zhuǎn)軸方向這個(gè)修正能消除大部分渦核線上的毛刺。4. Tecplot 可視化把 Liutex 結(jié)果變成能寫(xiě)進(jìn)論文的圖4.1 從 Matlab 導(dǎo)出 Tecplot 可讀的 ASCII 數(shù)據(jù)Matlab 算出的 Liutex 三個(gè)分量要和原始速度場(chǎng)一起導(dǎo)出方便在 Tecplot 里做多變量對(duì)照。最省事的方式是寫(xiě)一個(gè)標(biāo)準(zhǔn)的 Tecplot ASCII 文件點(diǎn)格式坐標(biāo)在前速度次之Liutex 分量放在最后。代碼塊如下filename liutex_export.dat; fid fopen(filename, w); fprintf(fid, TITLE Liutex Field\n); fprintf(fid, VARIABLES X, Y, Z, U, V, W, LiutexX, LiutexY, LiutexZ\n); fprintf(fid, ZONE I%d, J%d, K%d, FPOINT\n, nx, ny, nz); for k 1:nz for j 1:ny for i 1:nx fprintf(fid, %12.6e %12.6e %12.6e %12.6e %12.6e %12.6e %12.6e %12.6e %12.6e\n, ... X(i,j,k), Y(i,j,k), Z(i,j,k), ... U(i,j,k), V(i,j,k), W(i,j,k), ... LiutexX(i,j,k), LiutexY(i,j,k), LiutexZ(i,j,k)); end end end fclose(fid);文件頭里最關(guān)鍵的兩行是 VARIABLES 和 ZONE。VARIABLES 的變量名必須和后面數(shù)據(jù)列一一對(duì)應(yīng)Tecplot 加載時(shí)按字符串匹配變量名ZONE 末尾的 FPOINT 表示數(shù)據(jù)按逐點(diǎn)方式排列而非塊狀排列選錯(cuò)會(huì)導(dǎo)致 Tecplot 把整塊數(shù)據(jù)讀成亂序。這個(gè)文件如果網(wǎng)格規(guī)模在 200 萬(wàn)點(diǎn)以上ASCII 格式會(huì)膨脹到幾百 MB讀取很慢。這時(shí)候可以改用 Tecplot 的二進(jìn)制格式或者先把數(shù)據(jù)寫(xiě)成分塊排布的 FBLOCK字段按維度排體積能小不少。4.2 云圖、等值面與切片參數(shù)怎么定數(shù)據(jù)加載進(jìn) Tecplot 后左側(cè) Zone 面板里能直接看到新增的變量。顯示 Liutex 強(qiáng)度用等值面最直觀因?yàn)樾D(zhuǎn)強(qiáng)度的數(shù)值是客觀算出來(lái)的不需要像 Q 判據(jù)那樣糾結(jié)閾值。具體取值我一般這樣定先用菜單 Compute → 計(jì)算 Liutex 向量的模長(zhǎng)放到一個(gè)新變量 R 里點(diǎn)擊 Histogram 看分布再選一個(gè)對(duì)應(yīng)渦核區(qū)域占總體積大約 5%~10% 的值作為等值面閾值。對(duì)于典型的圓柱繞流或者機(jī)翼尾跡R 的值通常在 2~20 之間低于 1 的區(qū)域基本都是剪切噪聲。插一句切片問(wèn)題直接在流場(chǎng)中間切一個(gè)平面用云圖顯示 R 的分布顏色表選用 Tecplot 默認(rèn)的“Blue-White-Red”把最小值和最大值對(duì)標(biāo)到渦核內(nèi)外這樣平面圖上能看到渦核的圓斑和 Q 判據(jù)云圖相比邊界更銳利。如果你需要顯示旋轉(zhuǎn)軸方向把 LiutexX/Y/Z 作為矢量變量疊加到切面上箭頭方向就是局部旋轉(zhuǎn)軸注意 Tecplot 矢量圖的長(zhǎng)度比例因子要調(diào)到能區(qū)分相鄰點(diǎn)方向的程度通常取網(wǎng)格間距的 0.5 倍。4.3 三維顯示里旋轉(zhuǎn)軸的呈現(xiàn)技巧三維等值面只能看強(qiáng)度分布看不出旋轉(zhuǎn)軸的連續(xù)變化。一個(gè)有效做法是畫(huà)一條過(guò)渦核中心的線沿線采樣 Liutex 方向用 Tecplot 的 Polyline 功能從渦核入口追到出口然后把沿線 Liutex 向量的三個(gè)分量畫(huà)成曲線圖。這樣可以直觀看到旋轉(zhuǎn)軸在空間里的偏轉(zhuǎn)角度這在論文里比一張彩色等值面圖更有說(shuō)服力。畫(huà)三維半透明等值面時(shí)我習(xí)慣把 R 的最小等值面閾值設(shè)為峰值 R 的 30% 到 50%讓透明層顯示外圈、不透明層顯示內(nèi)核。在 Tecplot 里設(shè)置方法是在等值面面板中添加多個(gè)等值面值然后分別調(diào)整透明度。這套顯示參數(shù)組合適合大多數(shù) FPOINT 的穩(wěn)態(tài)場(chǎng)數(shù)據(jù)。瞬態(tài)場(chǎng)數(shù)據(jù)則不建議把每一幀都導(dǎo)成 ASCII 文件最好在 Matlab 里把多幀的 Liutex 向量統(tǒng)一存成一個(gè) NETCDF 文件利用 Tecplot 對(duì) NETCDF 的原生支持做時(shí)間動(dòng)畫(huà)否則幀數(shù)一多文件管理會(huì)亂成一團(tuán)。5. Liutex 使用避坑五個(gè)容易翻車(chē)的操作與排查5.1 旋轉(zhuǎn)軸方向頻繁反轉(zhuǎn)現(xiàn)象畫(huà)渦核線矢量圖時(shí)相鄰網(wǎng)格點(diǎn)的旋轉(zhuǎn)軸方向出現(xiàn) 180 度翻轉(zhuǎn)顯示效果像是箭頭的分布隨機(jī)跳變渦核中心線上表現(xiàn)得尤其明顯。原因eig 返回的特征向量方向本來(lái)就沒(méi)有固定朝向負(fù)方向也是合法解只用渦量方向做校正在渦核中心附近渦量很小甚至為零dot(r, ω) 的符號(hào)被噪聲主導(dǎo)導(dǎo)致方向隨機(jī)翻轉(zhuǎn)。解決在方向校準(zhǔn)前先用最小渦量閾值過(guò)濾|ω| 小于最大渦量 1% 的點(diǎn)不翻轉(zhuǎn)直接用相鄰已確定的梯度方向讓渦核內(nèi)部的軸方向順著流向上游的一致性趨勢(shì)走。我實(shí)際做法是分兩遍第一遍正常算 Liutex 向量第二遍把該層每個(gè)方向的特征向量分量的第一個(gè)非零分量符號(hào)對(duì)齊到前一個(gè)點(diǎn)兩遍過(guò)后軸方向基本連續(xù)。5.2 虛部閾值設(shè)置不當(dāng)引出假旋轉(zhuǎn)現(xiàn)象在無(wú)旋的均勻剪切流動(dòng)區(qū)域計(jì)算出來(lái)的 Liutex 強(qiáng)度不為零有時(shí)候在邊界層里還能畫(huà)出大片的“偽渦”。原因單精度速度場(chǎng)的梯度噪聲達(dá)到 1e-6 量級(jí)梯度張量的特征值虛部也會(huì)出現(xiàn)微小的偽復(fù)數(shù)值如果閾值設(shè)成 1e-10這些偽復(fù)特征根就會(huì)全被當(dāng)成真實(shí)旋轉(zhuǎn)。解決檢查數(shù)據(jù)精度。單精度場(chǎng)把虛部閾值提高到 1e-6 甚至 1e-5雙精度場(chǎng)再逐步降低同時(shí)結(jié)合渦量做輔助篩選如果 Liutex 強(qiáng)度很大但該點(diǎn)渦量幾乎為零這幾乎一定是數(shù)值噪聲直接剔除。5.3 網(wǎng)格維度順序?qū)е?Tecplot 數(shù)據(jù)錯(cuò)位現(xiàn)象Tecplot 里云圖顯示時(shí)強(qiáng)度分布明明是對(duì)的但坐標(biāo)軸的方向亂了切片也切不到預(yù)期的位置。原因Matlab 里數(shù)組維度順序是 (nx, ny, nz) 對(duì)應(yīng) (x, y, z)寫(xiě)導(dǎo)出時(shí)用的是三重循環(huán)循環(huán)順序或索引順序與 ZONE 里聲明的 IJK 順序不匹配數(shù)據(jù)列當(dāng)然對(duì)不上坐標(biāo)。解決導(dǎo)出前把 ZONE 行里的 I、J、K 和循環(huán)外層對(duì)應(yīng)好。我通常統(tǒng)一用最外層 k 循環(huán)對(duì)應(yīng) Z中間 j 循環(huán)對(duì)應(yīng) Y最內(nèi)層 i 循環(huán)對(duì)應(yīng) X強(qiáng)制約束住就不會(huì)錯(cuò)。導(dǎo)出后取樣首尾若干行核對(duì) X 坐標(biāo)是否從最小漸變到最大確認(rèn)后再進(jìn) Tecplot。5.4 等值面閾值選錯(cuò)看不出渦核結(jié)構(gòu)現(xiàn)象畫(huà) R 的等值面時(shí)閾值調(diào)來(lái)調(diào)去要么看不到東西要么整個(gè)流動(dòng)區(qū)域全被包成一個(gè)大團(tuán)塊。原因R 是客觀值但沒(méi)有“歸一化”不同工況量級(jí)差很大。低速流場(chǎng) R 才 0.1高速流場(chǎng) R 超過(guò) 100拿上一步用的 5 去套新數(shù)據(jù)當(dāng)然失靈。解決每換一個(gè)案例先算 R 的最大值和直方圖把最大值的 1/3 作為第一個(gè)等值面試。另外不要直接用 Liutex 向量分量畫(huà)等值面那些分量有正負(fù)會(huì)同時(shí)出現(xiàn)正反兩組面只會(huì)把圖弄亂。5.5 非結(jié)構(gòu)化網(wǎng)格上的梯度計(jì)算偏差現(xiàn)象同一套物理流動(dòng)從非結(jié)構(gòu)網(wǎng)格插值到結(jié)構(gòu)化網(wǎng)格后Liutex 強(qiáng)度與原始求解器后處理結(jié)果差異很大渦核線形狀都對(duì)不上。原因插值過(guò)程本身會(huì)抹掉一部分速度梯度信息特別是邊界層和剪切層位置梯度大、網(wǎng)格密插值誤差容易被放大。解決先檢查原始網(wǎng)格在目標(biāo)區(qū)域的網(wǎng)格密度是否和結(jié)構(gòu)化網(wǎng)格匹配。常用做法是先把原始網(wǎng)格上速度場(chǎng)用高階插值比如 RBF 插值轉(zhuǎn)出來(lái)再用最小二乘法在目標(biāo)網(wǎng)格上重構(gòu)梯度這樣比簡(jiǎn)單線性插值保真度高很多。如果兩個(gè)網(wǎng)格差異實(shí)在太大就別強(qiáng)求 Liutex 數(shù)值完全一致重點(diǎn)比較渦核位置而不是強(qiáng)度絕對(duì)值。6. 用最小渦驗(yàn)證 Liutex 實(shí)現(xiàn)一個(gè)半小時(shí)能做完的自檢流程拿到這套資源后我強(qiáng)烈建議你先跑一遍內(nèi)置的驗(yàn)證案例而不是直接拿來(lái)算自己的流場(chǎng)。驗(yàn)證流場(chǎng)用解析的 Lamb-Oseen 渦二維足夠因此切成一個(gè)很小的三維網(wǎng)格。計(jì)算域設(shè) 60×60×3 個(gè)點(diǎn)物理尺寸 0.1m × 0.1m × 0.005m渦核半徑 rc5mm最大切向速度 Uθmax1m/s渦軸沿 z 方向。把解析速度場(chǎng)導(dǎo)入腳本輸出 Liutex 向量后檢查兩個(gè)硬指標(biāo)。第一個(gè)指標(biāo)是旋轉(zhuǎn)軸方向理論解里渦軸就是 z 方向所以 Liutex 向量的 x、y 分量應(yīng)該接近零z 分量應(yīng)該接近正值。在渦核中心那一列點(diǎn)上打印出 LiutexZ如果符號(hào)一致、數(shù)量級(jí)和理論角速度對(duì)得上說(shuō)明特征分解和方向校正部分正常。第二個(gè)指標(biāo)是旋轉(zhuǎn)強(qiáng)度峰值R 的理論峰值在渦核半徑處大約為 2 × Uθmax / rc 400 s?1實(shí)際算出來(lái)在 380~420 之間都算正常超出這個(gè)范圍就得回頭查虛部閾值和梯度差分步長(zhǎng)。驗(yàn)證通過(guò)后再跑一個(gè)三維傾斜渦的案例把渦軸旋轉(zhuǎn) 45 度檢查 Liutex 向量的三個(gè)分量是否跟隨渦軸一起偏轉(zhuǎn)。這一步能排查旋轉(zhuǎn)軸方向標(biāo)定的問(wèn)題確認(rèn)無(wú)誤后就可以放心處理真實(shí)流場(chǎng)。最后給你一個(gè)我踩過(guò)坑之后養(yǎng)成的習(xí)慣每處理一個(gè)新案例我不會(huì)直接出圖而是先把 R 的最大值、渦核中心點(diǎn)坐標(biāo)和趨向性打印出來(lái)核對(duì)速度量級(jí)和網(wǎng)格尺寸是否匹配預(yù)期。這些數(shù)據(jù)不會(huì)騙人比看圖猜結(jié)構(gòu)可靠得多。希望這次關(guān)于 Liutex 的拆解能幫你少走些彎路拿到資源后務(wù)必備好解析驗(yàn)證的數(shù)據(jù)集再動(dòng)手進(jìn)入自己的項(xiàng)目。本文還有配套的精品資源點(diǎn)擊獲取