version 1.0.0.1 (78.3 MB) by
Oliver Schmidt

Implements the frequency domain form of proper orthogonal decomposition (POD)

SPOD() is a Matlab implementation of the frequency domain form of proper orthogonal decomposition (POD, also known as principle component analysis or Karhunen-Loève decomposition) called spectral proper orthogonal decomposition (SPOD). SPOD is derived from a space-time POD problem for stationary flows [1,2] and leads to modes that each oscillate at a single frequency. SPOD modes represent dynamic structures that optimally account for the statistical variability of stationary random processes.

The large-eddy simulation data provided along with this example is a subset of the database of a Mach 0.9 turbulent jet described in [3] and was calculated using the unstructured flow solver Charles developed at Cascade Technologies. If you are using the database in your research or teaching, please include explicit mention of Brès et al. [3]. The test database consists of 5000 snapshots of the symmetric component (m=0) of a round turbulent jet.

spod.m is a stand-alone Matlab function with no toolbox dependencies. All other Matlab files contained in this repository are related to the six examples that demonstrate the functionality of the code (see file descriptions below). A physical interpretation of the results obtained from the examples can be found in [4].

Download

------------

Repository zip file with examples (81.5 MB): https://github.com/SpectralPOD/spod_matlab/archive/master.zip

Matlab function only (15 KB): https://raw.githubusercontent.com/SpectralPOD/spod_matlab/master/spod.m

File Description

------------------

spod.m - Spectral proper orthogonal decomposition in Matlab

example_1.m - Inspect data and plot SPOD spectrum

example_2.m - Plot SPOD spectrum and inspect SPOD modes

example_3.m - Specify spectral estimation parameters and use weighted inner product

example_4.m - Calculate the SPOD of large data and save results on hard drive

example_5.m - Calculate full SPOD spectrum of large data

example_6.m - Calculate and plot confidence intervals for SPOD eigenvalues

jet_data/getjet.m - Interfaces external data source with SPOD() (examples 4-5)

utils/trapzWeightsPolar.m - Integration weight matrix for cylindrical coordinates (examples 3-6)

utils/jetLES.mat - Mach 0.9 turbulent jet test database

LICENSE.txt - License

References

-------------

[1] Towne, A., Schmidt, O. T., Colonius, T., Spectral proper orthogonal decomposition and its relationship to dynamic mode decomposition and resolvent analysis, J. of Fluid Mech. 847, 821–867, 2018

[2] Lumley, J. L., Stochastic tools in turbulence, Academic Press, 1970

[3] G. A. Brès, P. Jordan, M. Le Rallic, V. Jaunet, A. V. G. Cavalieri, A. Towne, S. K. Lele, T. Colonius, O. T. Schmidt, Importance of the nozzle-exit boundary-layer state in subsonic turbulent jets, J. of Fluid Mech. 851, 83-124, 2018

[4] Schmidt, O. T., Towne, A., Rigas, G., Colonius, T., Bres, G. A., Spectral analysis of jet turbulence, J. of Fluid Mech. 855, 953–982, 2018

Towne, A., Schmidt, O. T., Colonius, T., Spectral proper orthogonal decomposition and its relationship to dynamic mode decomposition and resolvent analysis, J. of Fluid Mech. 847, 821–867, 2018

Created with
R2017b

Compatible with any release

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!Create scripts with code, output, and formatted text in a single executable document.

Lican WangChi Young MoonChristian Amor RodríguezOliver SchmidtHi Christian,

Regarding your first question: increasing the frequency bins by reducing the number of snapshots per block is one option, as you mentioned (you're sacrificing frequency resolution doing that though). A better option is to use more data, if available/computable, to get more realizations of the Fourier transform. Regarding your second question: assuming that the sampling frequency is fixed, your only option here is to increase the number of snapshots per block. It really sounds like you'll need more data to get nicely converged modes...

Christian Amor RodríguezHi Oliver. Do you know any form to reduce the noise that could appear on the figures? I'm trying to calculate the first SPOD mode from some frequencies and I can obtain a cleaner figure from low frequencies instead from highest. I reduced the number of snapshots per block, reaching some clean figures but only for certain frequencies.

Also I'm trying to replicate some results that I have previously obtained using a DMD on my data, so I want to calculate the SPOD modes using an array with some frequencies that I founded relevant after calculate the DMD modes. It is possible to introduce manually some frequencies instead of those that the program calculate itself?

Thanks, Christian.

Oliver SchmidtOmid,

Have you tried out the examples? Each of the 6 examples builds on the previous one and introduces new features. The examples also show you how to visualize the results.

Oliver

omid nematollahi maleHi Oliver, Do you have any tutorial to work with this program? I appreciate you if you can prepare a tutorial on youtube or anywhere.

Best,

Omid

Oliver SchmidtYonatan,

Below this message is a code snippet that you can append at the end of example_2.m to animate the modes on screen.

Cheers, Oliver

%% Animate the same modes.

% Note how all wavepackets travel at approximately the same phase

% speed c_ph. The reason is that their streamwise wavenumber k_x changes

% with frequency such that c_ph = omega/k_x is approximately constant.

figure

nt = 30;

T = 1/f(10); % period of the 10th frequency

time = linspace(0,T,nt); % animate over one period

count = 1;

for ti = 1:nt

for fi = [10 15 20]

for mi = [1 2]

subplot(3,2,count)

pcolor(x,r,real(squeeze(P(fi,:,:,mi)*exp(2i*pi*f(fi)*time(ti))))), shading interp, axis equal tight, caxis(max(abs(caxis))*[-1 1])

xlabel('x'), ylabel('r'), title(['$f=' num2str(f(fi),'%.2f$') ', mode ' num2str(mi) ', $\lambda=' num2str(L(fi,mi),'%.2g$')])

xlim([0 10]); ylim([0 2])

count = count + 1;

hold on

end

end

drawnow

hold off

count = 1;

end

Yonatan CadavidI wonder if is it possible to see the time evolution or oscilation of one coherent structure or several coherent structures with SPOD.