從DNS數據中提取出了壓力數據,我想對其進行傅里葉變換,看看在頻域里能否找到一些關於壓力的有用的信息。我分別用了Matlab和Tecplot進行傅里葉變換,並對結果進行了比較,二者得到了一致的結果。現將過程記錄一下。
Matlab:
Matlab提供了現成的傅里葉變換工具:fft,在開始變換前,需要將數據的直流分量(DC),(對我的case來說)也就是壓力的平均值減去,否則FFT之后會在0Hz頻率處有一個很強的頻率分量,有可能會導致其他頻率分量被淹沒。
圖1. 壓力隨時間的變化
圖1 顯示了將要進行傅里葉變換的壓力數據。
首先用Matlab進行傅里葉變換。
dt = 7.7e-5 % DNS中的時間步長
freq_save = 100 % 數據保存頻率為100,即100個時間步(dt)保存一次
由此可以計算出,我的采樣周期為:
T = dt * freq_save = 7.7e-5 * 100 = 0.0077
那么采樣頻率為:
Fs = 1 / T = 130
代碼:
1 %% load pressure data and plot it 2 clc; 3 clear all; 4 close all; 5 6 delt_t=5; % = c/u_infinity = 1/0.2 = 5, reference time 7 8 load('time_pressure_1-7.mat') 9 load('time_pressure_8-15.mat') 10 P=[time_pressuer_1_7(:,2:8) time_pressuer_8_15(:,2:9)]; 11 time=(time_pressuer_1_7(:,1)-time_pressuer_1_7(1,1))*delt_t; % real time 12 13 tap_pressure = P(:,9); 14 15 figure 16 plot(time, tap_pressure) 17 18 Fs = 130; % 1 / (7.7*e-5*100) 19 T = 1/Fs; 20 L = length(tap_pressure); 21 22 P_mean = mean(tap_pressure); % the mean value will appear in the 0Hz of fft with a peak value 23 % need to be removed from the original data 24 25 fft_p = fft(tap_pressure-P_mean); % detrent the data, remove the mean value 26 p2 = abs(fft_p/L); % amplitude 27 p1 = p2(1:L/2+1); % single side amplitude 28 p1(2:end-1) = 2*p1(2:end-1); 29 30 f = Fs*(0:(L/2))/L; 31 figure 32 plot(f(1:150),p1(1:150)) % here I just plot the low frequency components, thee is nearly no high frequence components 33 title('single-sided amplitude spectrum of y') 34 xlabel('f (Hz)') 35 ylabel('|p1(f)|') 36 37 [m,n] = max(p1); 38 peak_freq = f(n); % find the peak frequence
代碼中有相應的解釋。
上圖:
(a) fft for original pressure
(b) fft after remove the mean pressure
圖2 頻域中的壓力
圖2(a)是沒有減去壓力平均值的fft,可以看到 0Hz 頻率處有一個峰值,這個峰值其實對應的就是壓力的平均值(數據的直流分量),因為平均值很大(相對),導致其他頻率分量完全被淹沒,觀察不到。(b) 就是去除平均值之后的fft,0Hz 頻率處不再有峰值, 其他頻率分量也都顯示了出來。橫坐標是頻率,縱坐標表示各頻率分量的幅值。
Tecplot:
Tecplot中做傅里葉變換(fft)要相對簡單一些,在做變換之前同樣需要減去壓力的平均值。在Tecplot中,可以用積分的方法計算某個量的平均值。
利用 Analysis -- Perform Integration 選項,type選擇average,求平均值。再利用自定義公式求出減去平均值后的脈動壓力。
然后,選擇Data -- Fourier Transform,進行傅里葉變換,注意選擇正確的自變量和因變量。
上圖:
可以看到Tecplot得到了和Matlab相同的結果。