Instability of the wake of a cylinder with STABFEM

This scripts reproduces the figures of section 3 of the paper
A Practical Review on Linear and Nonlinear Global Approaches to Flow Instabilities
by D. Fabre, V. Citro, D. Ferreira Sabino, P. Bonnefis, J. Sierra, F. Giannetti & M. Pigou.
Appl. Mech. Rev. Nov 2018, 70(6): 060802
  1. Generation of an adapted mesh
  2. Base-flow properties for Re = [2-40]
  3. Stability curves St(Re) and sigma(Re) for Re = [40-100]
  4. Determination of the instability threshold and Weakly-Nonlinear analysis
  5. Harmonic-Balance for Re = [REc-100]
  6. Self-consistent model for Re = 100

NB this script uses some depreciated syntax but is left in the present stage for legacy.

The script still runs with StabFem 3.5 (june 2020) but it is recommended for new users to look at most recent examples on the stabfem website.

Contents

CHAPTER 0 : set the global variables needed by the drivers

close all;
addpath([fileparts(fileparts(pwd)), '/SOURCES_MATLAB']);
SF_Start('verbosity',2,'ffdatadir','./WORK/');
figureformat='tif'; AspectRatio = 0.75; % aspect ration for figures; nb set to 0.53 to regenerate exactly figures from paper
system('mkdir FIGURES');
set(groot, 'defaultAxesTickLabelInterpreter','latex');
set(groot, 'defaultLegendInterpreter','latex');
tic;
NOTICE  - Initialization already done

##### CHAPTER 1 : COMPUTING THE MESH WITH ADAPTMESH PROCEDURE

disp(' STARTING ADAPTMESH PROCEDURE : ');
disp(' ');
disp(' LARGE MESH : [-40:80]x[0:40] ');
disp(' ');
%bf=SF_Init('Mesh_Cylinder.edp',[-40 80 40]);
ffmesh = SF_Mesh('Mesh_Cylinder.edp','Params',[-40 80 40],'problemtype','2D');
bf=SF_BaseFlow(ffmesh,'Re',1);
bf=SF_Adapt(bf,'Hmax',2);
bf=SF_BaseFlow(bf,'Re',10);
bf=SF_BaseFlow(bf,'Re',60);
bf=SF_Adapt(bf,'Hmax',2);
meshstrategy = 'D';
  % select 'D' or 'S'
 % 'D' will use mesh adapted on direct eigenmode (mesh M_4):
 %     this is necessary to compute correctly the structure of the mode (fig. 5a)
 %     and the energy-amplitude (fig. 7d)
 % 'S' will use mesh adapted on sensitivity (mesh M_2):
 %     figs. (5a) and (7d) will be wrong, on the other all other results
 %     will be correct and nonlinear computations will be much much faster.
if(meshstrategy=='D')
     bf=SF_BaseFlow(bf,'Re',60);
     disp('using mesh adaptated to EIGENMODE (M4) ')
     [ev,em] = SF_Stability(bf,'shift',0.04+0.76i,'nev',1,'type','D');
     bf=SF_Adapt(bf,em,'Hmax',2);
else
    disp('mesh adaptation to SENSITIVITY : ') % This is mesh M2 from the appendix
    [ev,em] = SF_Stability(bf,'shift',0.04+0.76i,'nev',1,'type','D');
    [ev,emA] = SF_Stability(bf,'shift',0.04+0.76i,'nev',1,'type','A');
    emS = SF_Sensitivity(bf,em,emA);
    bf = SF_Adapt(bf,emS,'Hmax',2);
end

    [ev,em]  = SF_Stability(bf,'shift',0.04+0.76i,'nev',1,'type','D');
    [ev,emA] = SF_Stability(bf,'shift',0.04+0.76i,'nev',1,'type','A');
    emS      = SF_Sensitivity(bf,em,emA);% compute on last mesh
 STARTING ADAPTMESH PROCEDURE : 
 
 LARGE MESH : [-40:80]x[0:40] 
 
WARNING -  Your program uses depreciated syntax 'Params' to pass parameters 
WARNING -  It is advised to switch to new method using 'Options', and modify your Freefem solvers to detect parameters using 'getARGV'  
WARNING -  Folder ./WORK/MESHBF not found : creating it !
using mesh adaptated to EIGENMODE (M4) 
WARNING -  baseflow.iter = 0 ! it seems your baseflow was projected after mesh adaptation but not recomputed !
WARNING -  baseflow.iter = 0 ! it seems your baseflow was projected after mesh adaptation but not recomputed !

CHAPTER 1b : DRAW FIGURES

plot the mesh (full size)
figure();SF_Plot(bf,'mesh','xlim',[-40 80],'ylim',[0 40]);
title('Initial mesh (full size)');
box on; pos = get(gcf,'Position'); pos(4)=pos(3)*AspectRatio;set(gcf,'Position',pos); % resize aspect ratio
set(gca,'FontSize', 18);
saveas(gca,'FIGURES/Cylinder_Mesh_Full',figureformat);

%  plot the mesh (zoom)
figure();SF_Plot(bf,'mesh','xlim',[-1.5 4.5],'ylim',[0 3]);
box on; pos = get(gcf,'Position'); pos(4)=pos(3)*AspectRatio;set(gcf,'Position',pos); % resize aspect ratio
set(gca,'FontSize', 18);
saveas(gca,'FIGURES/Cylinder_Mesh',figureformat);
plot the base flow for Re = 60
figure();
SF_Plot(bf,'p','xlim',[-1.5 4.5],'ylim',[0 3],'cbtitle','\omega_z','colormap','redblue','colorrange','centered','boundary','on','bdlabels',2,'bdcolors','k');
hold on;SF_Plot(bf,'psi','contour','only','clevels',[-.02 0 .2 1 2 5],'xlim',[-1.5 4.5],'ylim',[0 3]);
box on; pos = get(gcf,'Position'); pos(4)=pos(3)*AspectRatio;set(gcf,'Position',pos); % resize aspect ratio
set(gca,'FontSize', 18);
saveas(gca,'FIGURES/Cylinder_BaseFlowRe60',figureformat);


%  plot the eigenmode for Re = 60

figure();SF_Plot(em,'ux','xlim',[-2 8],'ylim',[0 5],'colormap','redblue','colorrange','cropcentered',...
    'boundary','on','bdlabels',2,'bdcolors','k');
box on; pos = get(gcf,'Position'); pos(4)=pos(3)*AspectRatio;set(gcf,'Position',pos); % resize aspect ratio
set(gca,'FontSize', 18);
saveas(gca,'FIGURES/Cylinder_EigenModeRe60_AdaptS',figureformat);  %

figure();SF_Plot(emA,'ux','xlim',[-2 8],'ylim',[0 5],'colormap','redblue','colorrange','cropcentered','boundary','on','bdlabels',2,'bdcolors','k');
box on; pos = get(gcf,'Position'); pos(4)=pos(3)*AspectRatio;set(gcf,'Position',pos); % resize aspect ratio
set(gca,'FontSize', 18);
saveas(gca,'FIGURES/Cylinder_EigenModeAdjRe60',figureformat);

figure();SF_Plot(emS,'S','xlim',[-2 4],'ylim',[0 3],'colormap','ice','boundary','on','bdlabels',2,'bdcolors','k');
hold on;
SF_Plot(bf,'psi','contour','only','clevels', [0 0],'xlim',[-2 4],'ylim',[0 3],...
'colormap','ice','colorrange',[min(real(emS.S)), max(real(emS.S))]);

box on; pos = get(gcf,'Position'); pos(4)=pos(3)*AspectRatio;set(gcf,'Position',pos); % resize aspect ratio
set(gca,'FontSize', 18);
saveas(gca,'FIGURES/Cylinder_SensitivityRe60',figureformat);
WARNING -  Error in SF_Plot : Field psi does not exist
WARNING -  Error in SF_Plot : Field psi does not exist

Demonstration of "spectrum explorator" mode

plotoptions = {'ux','xlim',[-2 8],'ylim',[0 5],'colormap','redblue','colorrange','cropcentered',...
                'boundary','on','bdlabels',2,'bdcolors','k'};
ev = SF_Stability(bf,'nev',20,'shift',0.74i,'PlotSpectrum',true,'PlotModeOptions',plotoptions);
figure(100); box on; pos = get(gcf,'Position'); pos(4)=pos(3)*AspectRatio;set(gcf,'Position',pos); % resize aspect ratio
set(gca,'FontSize', 14);
saveas(gca,'FIGURES/SpectrumExplorator',figureformat);


%figure();SF_Plot(emE,'endogeneity.re','contour','endogeneity.im','xlim',[-2 4],'ylim',[0 3],...
%    'colormap','redblue','cstyle','patchdashedneg','boundary','on','bdlabels',2,'bdcolors','k');
%hold on;
%SF_Plot(bf,'psi','contour','on','clevels', [0 0],'CColor','b','CStyle','monochrome',...
%    'ColorBar','off','xlim',[-2 4],'ylim',[0 3],'colormap','redblue',...
%    'colorrange',[min(real(emE.endogeneity)), max(real(emE.endogeneity))]);
%box on; pos = get(gcf,'Position'); pos(4)=pos(3)*AspectRatio;set(gcf,'Position',pos); % resize aspect ratio
%set(gca,'FontSize', 18);
%saveas(gca,'FIGURES/Cylinder_EndogeneityRe60',figureformat);
%pause(0.1);
WARNING -  baseflow.iter = 0 ! it seems your baseflow was projected after mesh adaptation but not recomputed !

CHAPTER 2 : DESCRIPTION OF BASE FLOW PROPERTIES (range 2-50)

Re_BF = [2 : 2: 50];
Fx_BF = []; Lx_BF = [];
    for Re = Re_BF
        bf = SF_BaseFlow(bf,'Re',Re);
        Fx_BF = [Fx_BF,bf.Fx];
        Lx_BF = [Lx_BF,bf.Lx];
    end


% chapter 2B : figures


figure(22);hold off;
LiteratureData=csvread('./literature_data/fig3a_Lw_giannetti2004.csv'); %read literature data
plot(LiteratureData(:,1),LiteratureData(:,2)+0.5,'r+-','LineWidth',2);
plot(Re_BF,Fx_BF,'b+-','LineWidth',2);
xlabel('Re');ylabel('Fx');
box on; pos = get(gcf,'Position'); pos(4)=pos(3)*AspectRatio;set(gcf,'Position',pos); % resize aspect ratio
set(gca,'FontSize', 18);
saveas(gca,'FIGURES/Cylinder_Fx_baseflow',figureformat);

figure(23);hold off;
plot(Re_BF,Lx_BF,'b+-','LineWidth',2);
xlabel('Re');ylabel('Lx');
box on; pos = get(gcf,'Position'); pos(4)=pos(3)*AspectRatio;set(gcf,'Position',pos); % resize aspect ratio
set(gca,'FontSize', 18);
saveas(gca,'FIGURES/Cylinder_Lx_baseflow',figureformat);
pause(0.1);

CHAPTER 3 : COMPUTING STABILITY BRANCH

    disp('COMPUTING STABILITY BRANCH')

% LOOP OVER RE FOR BASEFLOW + EIGENMODE
Re_LIN = [40 : 2: 100];
bf=SF_BaseFlow(bf,'Re',40);
[ev,em] = SF_Stability(bf,'shift',-.03+.72i,'nev',1,'type','D');

Fx_LIN = []; Lx_LIN = [];lambda_LIN=[];
    for Re = Re_LIN
        bf = SF_BaseFlow(bf,'Re',Re);
        Fx_LIN = [Fx_LIN,bf.Fx];
        Lx_LIN = [Lx_LIN,bf.Lx];
        [ev,em] = SF_Stability(bf,'nev',1,'shift','cont');
        lambda_LIN = [lambda_LIN ev];
    end
completed_lambda = 1;
COMPUTING STABILITY BRANCH

CHAPTER 3b : figures

figure(20);
LiteratureData=csvread('./literature_data/fig4a_sigma_giannetti2004.csv'); %read literature data
plot(LiteratureData(:,1),LiteratureData(:,2),'r+-','LineWidth',2);
plot(Re_LIN,real(lambda_LIN),'b+-');
xlabel('Re');ylabel('$\sigma$','Interpreter','latex');
box on; pos = get(gcf,'Position'); pos(4)=pos(3)*AspectRatio;set(gcf,'Position',pos); % resize aspect ratio
set(gca,'FontSize', 18);
saveas(gca,'FIGURES/Cylinder_Sigma_Re',figureformat);

figure(21);hold off;
LiteratureData=csvread('./literature_data/fig4b_st_giannetti2004.csv'); %read literature data
plot(LiteratureData(:,1),LiteratureData(:,2),'r+-','LineWidth',2);
plot(Re_LIN,imag(lambda_LIN)/(2*pi),'b+-');
xlabel('Re');ylabel('St');
box on; pos = get(gcf,'Position'); pos(4)=pos(3)*AspectRatio;set(gcf,'Position',pos); % resize aspect ratio
set(gca,'FontSize', 18);
saveas(gca,'FIGURES/Cylinder_Strouhal_Re',figureformat);
pause(0.1);


disp(' ');
disp('       cpu time for Linear calculations : ');
tlin = toc;
disp([ '   ' num2str(tlin) ' seconds']);
tic;
 
       cpu time for Linear calculations : 
   1985.0862 seconds

CHAPTER 4 : computation of weakly nonlinear expansion

disp(' ');
disp('######     ENTERING NONLINEAR PART       ####### ');
disp(' ');


% plot the eigenmode for Re = 60
%em.xlim = [-2 8]; em.ylim=[0,5];
%figure();SF_Plot(em,'ux1','colorrange','cropcentered','xlim',[-2 8],'ylim',[0 5]);
%title('Eigenmode for Re=60');
%box on; pos = get(gcf,'Position'); pos(4)=pos(3)*AspectRatio;set(gcf,'Position',pos); % resize aspect ratio
%set(gca,'FontSize', 18);
%saveas(gca,'FIGURES/Cylinder_EigenModeRe60_AdaptD',figureformat);  %




% DETERMINATION OF THE INSTABILITY THRESHOLD
disp('COMPUTING INSTABILITY THRESHOLD');
bf=SF_BaseFlow(bf,'Re',50);
[ev,em] = SF_Stability(bf,'shift',+.75i,'nev',1,'type','D');
[bf,em]=SF_FindThreshold(bf,em);
Rec = bf.Re;  Fxc = bf.Fx;
Lxc=bf.Lx;    Omegac=imag(em.lambda);


[ev,em] = SF_Stability(bf,'shift',1i*Omegac,'nev',1,'type','D');
[ev,emA] = SF_Stability(bf,'shift',1i*Omegac,'nev',1,'type','A');
[wnl,meanflow,mode] = SF_WNL(bf,em,'Adjoint',emA,'Retest',47.); % Here to generate a starting point for the next chapter
 
######     ENTERING NONLINEAR PART       ####### 
 
COMPUTING INSTABILITY THRESHOLD
WARNING -  option "AdjointType" is no longer operational
WARNING -  In SFcore_MoveDataFiles : File ./WORK/WNL_results.txt does not exit ! 

PLOTS of WNL predictions

epsilon2_WNL = -0.003:.0001:.005; % will trace results for Re = 40-55 approx.
Re_WNL = 1./(1/Rec-epsilon2_WNL);
A_WNL = wnl.Aeps.*real(sqrt(epsilon2_WNL));
Fy_WNL = wnl.Fyeps.*real(sqrt(epsilon2_WNL))*2; % factor 2 because of complex conjugate
omega_WNL =Omegac + epsilon2_WNL.*imag(wnl.Lambda) ...
                  - epsilon2_WNL.*(epsilon2_WNL>0)*real(wnl.Lambda)*imag(wnl.nu0+wnl.nu2)/real(wnl.nu0+wnl.nu2)  ;
Fx_WNL = wnl.Fx0 + wnl.Fxeps2.*epsilon2_WNL  ...
                 + wnl.FxA20.*real(wnl.Lambda)/real(wnl.nu0+wnl.nu2)*epsilon2_WNL.*(epsilon2_WNL>0) ;

figure(20);hold on;
plot(Re_WNL,real(wnl.Lambda)*epsilon2_WNL,'g--','LineWidth',2);hold on;

figure(21);hold on;
plot(Re_WNL,omega_WNL/(2*pi),'g--','LineWidth',2);hold on;
xlabel('Re');ylabel('St');

figure(22);hold on;
plot(Re_WNL,Fx_WNL,'g--','LineWidth',2);hold on;
xlabel('Re');ylabel('Fx');

figure(24); hold on;
plot(Re_WNL,abs(Fy_WNL),'g--','LineWidth',2);
xlabel('Re');ylabel('Fy')

figure(25);hold on;
plot(Re_WNL,A_WNL,'g--','LineWidth',2);
xlabel('Re');ylabel('AE')

pause(0.1);

CHAPTER 5 : SELF CONSISTENT / HB1

disp('SC quasilinear model on the range [Rec , 100]');
Re_HB = [Rec 47 47.5 48 49 50 52.5 55 60 65 70 75 80 85 90 95 100];

% THE STARTING POINT HAS BEEN GENERATED ABOVE, WHEN PERFORMING THE WNL ANALYSIS
Res = 47. ;

 Lx_HB = [Lxc]; Fx_HB = [Fxc]; omega_HB = [Omegac]; Aenergy_HB  = [0]; Fy_HB = [0];

[meanflow,mode] = SF_HB1(meanflow,mode,'sigma',0.,'Re',Res);

for Re = Re_HB(2:end)
    [meanflow,mode] = SF_HB1(meanflow,mode,'Re',Re);
    Lx_HB = [Lx_HB meanflow.Lx];
    Fx_HB = [Fx_HB meanflow.Fx];
    omega_HB = [omega_HB imag(mode.lambda)];
    Aenergy_HB  = [Aenergy_HB mode.AEnergy];
    Fy_HB = [Fy_HB mode.Fy];

    if(Re==60)
       figure();
       SF_Plot(meanflow,'vort','xlim',[-1.5 4.5],'ylim',[0 3],'cbtitle','\omega_z','colormap','redblue','colorrange','centered','boundary','on','bdlabels',2,'bdcolors','k');
       hold on;SF_Plot(meanflow,'psi','contour','only','clevels',[-.02 0 .2 1 2 5],'xlim',[-1.5 4.5],'ylim',[0 3]);
       box on; pos = get(gcf,'Position'); pos(4)=pos(3)*AspectRatio;set(gcf,'Position',pos); % resize aspect ratio
       set(gca,'FontSize', 18);
       saveas(gca,'FIGURES/Cylinder_MeanFlowRe60',figureformat);
    end
end
SC quasilinear model on the range [Rec , 100]
Note: To improve runtime build MEX function fftri2gridfast() from fftri2gridfast.c

chapter 5b : figures

%load('Cylinder_AllFigures.mat');

figure(21);hold off;
plot(Re_LIN,imag(lambda_LIN)/(2*pi),'b+-');
hold on;
plot(Re_WNL,omega_WNL/(2*pi),'g--','LineWidth',2);hold on;
plot(Re_HB,omega_HB/(2*pi),'r-','LineWidth',2);
plot(Rec,Omegac/2/pi,'ko');
LiteratureData=csvread('./literature_data/fig7a_st_Re_experience.csv'); %read literature data
plot(LiteratureData(:,1),LiteratureData(:,2),'ko','LineWidth',2);xlabel('Re');ylabel('St');
box on; pos = get(gcf,'Position'); pos(4)=pos(3)*AspectRatio;set(gcf,'Position',pos); % resize aspect ratio
set(gca,'FontSize', 18);
legend('Linear','WNL','HB1','Ref. [23]','Location','northwest');
saveas(gca,'FIGURES/Cylinder_Strouhal_Re_HB',figureformat);

figure(22);hold off;
plot(Re_LIN,Fx_LIN,'b+-');
hold on;
plot(Re_WNL,Fx_WNL,'g--','LineWidth',2);hold on;
plot(Re_HB,Fx_HB,'r+-','LineWidth',2);
plot(Rec,Fxc,'ro')
xlabel('Re');ylabel('Fx');
box on; pos = get(gcf,'Position'); pos(4)=pos(3)*AspectRatio;set(gcf,'Position',pos); % resize aspect ratio
set(gca,'FontSize', 18);
legend('BF','WNL','HB1','Location','south');
saveas(gca,'FIGURES/Cylinder_Cx_Re_HB',figureformat);

figure(23);hold off;
plot(Re_LIN,Lx_LIN,'b+-');
hold on;
plot(Re_HB,Lx_HB,'r+-','LineWidth',2);
LiteratureData=csvread('./literature_data/fig7e_Lx_mean.csv'); %read literature data
plot(LiteratureData(:,1),LiteratureData(:,2),'ko','LineWidth',2);
plot(Rec,Lxc,'ro','LineWidth',2);
xlabel('Re');ylabel('Lx');
box on; pos = get(gcf,'Position'); pos(4)=pos(3)*AspectRatio;set(gcf,'Position',pos); % resize aspect ratio
set(gca,'FontSize', 18);
legend('BF','HB1','Ref. [5]','Location','northwest');
saveas(gca,'FIGURES/Cylinder_Lx_Re_HB',figureformat);

figure(24);hold off;
plot(Re_WNL,abs(Fy_WNL),'g--','LineWidth',2);
hold on;
plot(Re_HB,real(Fy_HB),'r+-','LineWidth',2);
%title('Harmonic Balance results');
xlabel('Re');  ylabel('Fy')
box on;  pos = get(gcf,'Position');  pos(4)=pos(3)*AspectRatio;  set(gcf,'Position',pos); % resize aspect ratio
set(gca,'FontSize', 18);
legend('WNL','HB1','Location','south');
saveas(gca,'FIGURES/Cylinder_Cy_Re_SC',figureformat);

figure(25);hold off;
plot(Re_WNL,A_WNL,'g--','LineWidth',2);
hold on;
plot(Re_HB,Aenergy_HB,'r+-','LineWidth',2);
LiteratureData=csvread('./literature_data/fig7d_energy_amplitude.csv'); %read literature data
plot(LiteratureData(:,1),LiteratureData(:,2),'ko','LineWidth',2);
%title('Harmonic Balance results');
xlabel('Re');ylabel('A_E')
box on; pos = get(gcf,'Position'); pos(4)=pos(3)*AspectRatio; set(gcf,'Position',pos); % resize aspect ratio
set(gca,'FontSize', 18);
legend('WNL','HB1','Ref. [5]','Location','south');
if(meshstrategy=='D')
    filename = 'FIGURES/Cylinder_Energy_Re_SC_AdaptD';
else
    filename = 'FIGURES/Cylinder_Energy_Re_SC_AdaptS';
end
saveas(gca,filename,figureformat);


tnolin = toc;
disp(' ');
disp('       cpu time for Nonlinear calculations : ');
disp([ '   ' num2str(tnolin) ' seconds']);

disp(' ');
disp('Total cpu time for the linear & nonlinear calculations and generation of all figures : ');
disp([ '   ' num2str(tlin+tnolin) ' seconds']);
 
       cpu time for Nonlinear calculations : 
   5838.8994 seconds
 
Total cpu time for the linear & nonlinear calculations and generation of all figures : 
   7823.9855 seconds

Summary

SF_Status

% [[PUBLISH]] (this tag is to enable automatic publication as html ; don't touch it !)
% [[FIGURES]] (this tag is to enable automatic generation of figures for manual or article ; don't touch it !)

%
fclose(fopen('SCRIPT_CYLINDER_ALLFIGURES.success','w+'));
################################################################################
...  SUMMARY OF YOUR DATABASE FOLDER :    ./WORK/
#################################################################################
 
.... CONTENT OF DIRECTORY ./WORK/MESHES :
     (list of meshes previously created/adapted ; couples of .msh/.ff2m files )
 
Index | Name              | generation mode | Date                 | Nv      
1     | FFMESH_000001.msh | initial         | 22-mars-2023 18:00:48 | 2192    
2     | FFMESH_000002.msh | adapted         | 22-mars-2023 18:00:59 | 1867    
3     | FFMESH_000003.msh | adapted         | 22-mars-2023 18:01:17 | 2691    
4     | FFMESH_000004.msh | adapted         | 22-mars-2023 18:01:35 | 7886    
#################################################################################
 
     (list of base flows associated to newly computed meshes ; couples of .txt/.ff2m files )
     NB : these base flows are most likely simply projected and not recomputed, it is not recommended to use them
 
 Index | Name              | Type         | Date                 | Mesh file         | Re         | Fx        
 1     | FFDATA_000003.txt | BaseFlow     | 22-mars-2023 18:01:01 | FFMESH_000002.msh | 1          | 5.5507    
 2     | FFDATA_000004.txt | BaseFlow     | 22-mars-2023 18:01:19 | FFMESH_000003.msh | 60         | 0.63959   
 3     | FFDATA_000005.txt | BaseFlow     | 22-mars-2023 18:01:41 | FFMESH_000004.msh | 60         | 0.64333   
 
#################################################################################
#################################################################################
 
.... CONTENT OF DIRECTORY ./WORK/BASEFLOWS
     (couples of .txt/.ff2m files )
 
 Index | Name              | Type         | Date                 | Mesh file         | Re         | Fx        
 1     | FFDATA_000001.txt | BaseFlow     | 22-mars-2023 18:00:57 | FFMESH_000001.msh | 1          | 5.5507    
 2     | FFDATA_000002.txt | BaseFlow     | 22-mars-2023 18:01:08 | FFMESH_000002.msh | 10         | 1.3892    
 3     | FFDATA_000003.txt | BaseFlow     | 22-mars-2023 18:01:15 | FFMESH_000002.msh | 60         | 0.63959   
 4     | FFDATA_000004.txt | BaseFlow     | 22-mars-2023 18:01:27 | FFMESH_000003.msh | 60         | 0.64333   
 5     | FFDATA_000005.txt | BaseFlow     | 22-mars-2023 18:04:31 | FFMESH_000004.msh | 2          | 3.3986    
 6     | FFDATA_000006.txt | BaseFlow     | 22-mars-2023 18:04:56 | FFMESH_000004.msh | 4          | 2.2532    
 7     | FFDATA_000007.txt | BaseFlow     | 22-mars-2023 18:05:21 | FFMESH_000004.msh | 6          | 1.8046    
 8     | FFDATA_000008.txt | BaseFlow     | 22-mars-2023 18:05:45 | FFMESH_000004.msh | 8          | 1.5539    
 9     | FFDATA_000009.txt | BaseFlow     | 22-mars-2023 18:06:10 | FFMESH_000004.msh | 10         | 1.39      
 10    | FFDATA_000010.txt | BaseFlow     | 22-mars-2023 18:06:35 | FFMESH_000004.msh | 12         | 1.2725    
 11    | FFDATA_000011.txt | BaseFlow     | 22-mars-2023 18:06:59 | FFMESH_000004.msh | 14         | 1.1832    
 12    | FFDATA_000012.txt | BaseFlow     | 22-mars-2023 18:07:24 | FFMESH_000004.msh | 16         | 1.1126    
 13    | FFDATA_000013.txt | BaseFlow     | 22-mars-2023 18:07:49 | FFMESH_000004.msh | 18         | 1.0549    
 14    | FFDATA_000014.txt | BaseFlow     | 22-mars-2023 18:08:13 | FFMESH_000004.msh | 20         | 1.0066    
 15    | FFDATA_000015.txt | BaseFlow     | 22-mars-2023 18:08:35 | FFMESH_000004.msh | 22         | 0.96551   
 16    | FFDATA_000016.txt | BaseFlow     | 22-mars-2023 18:09:00 | FFMESH_000004.msh | 24         | 0.92995   
 17    | FFDATA_000017.txt | BaseFlow     | 22-mars-2023 18:09:22 | FFMESH_000004.msh | 26         | 0.89879   
 18    | FFDATA_000018.txt | BaseFlow     | 22-mars-2023 18:09:45 | FFMESH_000004.msh | 28         | 0.87119   
 19    | FFDATA_000019.txt | BaseFlow     | 22-mars-2023 18:10:09 | FFMESH_000004.msh | 30         | 0.84652   
 20    | FFDATA_000020.txt | BaseFlow     | 22-mars-2023 18:10:34 | FFMESH_000004.msh | 32         | 0.82429   
 21    | FFDATA_000021.txt | BaseFlow     | 22-mars-2023 18:10:59 | FFMESH_000004.msh | 34         | 0.80411   
 22    | FFDATA_000022.txt | BaseFlow     | 22-mars-2023 18:11:23 | FFMESH_000004.msh | 36         | 0.7857    
 23    | FFDATA_000023.txt | BaseFlow     | 22-mars-2023 18:11:47 | FFMESH_000004.msh | 38         | 0.76879   
 24    | FFDATA_000024.txt | BaseFlow     | 22-mars-2023 18:12:12 | FFMESH_000004.msh | 40         | 0.7532    
 25    | FFDATA_000025.txt | BaseFlow     | 22-mars-2023 18:12:37 | FFMESH_000004.msh | 42         | 0.73876   
 26    | FFDATA_000026.txt | BaseFlow     | 22-mars-2023 18:13:01 | FFMESH_000004.msh | 44         | 0.72533   
 27    | FFDATA_000027.txt | BaseFlow     | 22-mars-2023 18:13:24 | FFMESH_000004.msh | 46         | 0.71279   
 28    | FFDATA_000028.txt | BaseFlow     | 22-mars-2023 18:13:49 | FFMESH_000004.msh | 48         | 0.70105   
 29    | FFDATA_000029.txt | BaseFlow     | 22-mars-2023 18:14:11 | FFMESH_000004.msh | 50         | 0.69002   
 30    | FFDATA_000030.txt | BaseFlow     | 22-mars-2023 18:14:35 | FFMESH_000004.msh | 40         | 0.7532    
 31    | FFDATA_000031.txt | BaseFlow     | 22-mars-2023 18:15:02 | FFMESH_000004.msh | 40         | 0.7532    
 32    | FFDATA_000032.txt | BaseFlow     | 22-mars-2023 18:15:39 | FFMESH_000004.msh | 42         | 0.73876   
 33    | FFDATA_000033.txt | BaseFlow     | 22-mars-2023 18:16:18 | FFMESH_000004.msh | 44         | 0.72533   
 34    | FFDATA_000034.txt | BaseFlow     | 22-mars-2023 18:16:57 | FFMESH_000004.msh | 46         | 0.71279   
 35    | FFDATA_000035.txt | BaseFlow     | 22-mars-2023 18:17:36 | FFMESH_000004.msh | 48         | 0.70105   
 36    | FFDATA_000036.txt | BaseFlow     | 22-mars-2023 18:18:10 | FFMESH_000004.msh | 50         | 0.69002   
 37    | FFDATA_000037.txt | BaseFlow     | 22-mars-2023 18:18:45 | FFMESH_000004.msh | 52         | 0.67964   
 38    | FFDATA_000038.txt | BaseFlow     | 22-mars-2023 18:19:22 | FFMESH_000004.msh | 54         | 0.66984   
 39    | FFDATA_000039.txt | BaseFlow     | 22-mars-2023 18:19:58 | FFMESH_000004.msh | 56         | 0.66056   
 40    | FFDATA_000040.txt | BaseFlow     | 22-mars-2023 18:20:36 | FFMESH_000004.msh | 58         | 0.65177   
 41    | FFDATA_000041.txt | BaseFlow     | 22-mars-2023 18:21:13 | FFMESH_000004.msh | 60         | 0.64341   
 42    | FFDATA_000042.txt | BaseFlow     | 22-mars-2023 18:21:53 | FFMESH_000004.msh | 62         | 0.63546   
 43    | FFDATA_000043.txt | BaseFlow     | 22-mars-2023 18:22:32 | FFMESH_000004.msh | 64         | 0.62787   
 44    | FFDATA_000044.txt | BaseFlow     | 22-mars-2023 18:23:06 | FFMESH_000004.msh | 66         | 0.62063   
 45    | FFDATA_000045.txt | BaseFlow     | 22-mars-2023 18:23:42 | FFMESH_000004.msh | 68         | 0.61371   
 46    | FFDATA_000046.txt | BaseFlow     | 22-mars-2023 18:24:18 | FFMESH_000004.msh | 70         | 0.60707   
 47    | FFDATA_000047.txt | BaseFlow     | 22-mars-2023 18:24:55 | FFMESH_000004.msh | 72         | 0.60071   
 48    | FFDATA_000048.txt | BaseFlow     | 22-mars-2023 18:25:32 | FFMESH_000004.msh | 74         | 0.59461   
 49    | FFDATA_000049.txt | BaseFlow     | 22-mars-2023 18:26:10 | FFMESH_000004.msh | 76         | 0.58874   
 50    | FFDATA_000050.txt | BaseFlow     | 22-mars-2023 18:26:47 | FFMESH_000004.msh | 78         | 0.58309   
 51    | FFDATA_000051.txt | BaseFlow     | 22-mars-2023 18:27:24 | FFMESH_000004.msh | 80         | 0.57764   
 52    | FFDATA_000052.txt | BaseFlow     | 22-mars-2023 18:27:58 | FFMESH_000004.msh | 82         | 0.5724    
 53    | FFDATA_000053.txt | BaseFlow     | 22-mars-2023 18:28:33 | FFMESH_000004.msh | 84         | 0.56734   
 54    | FFDATA_000054.txt | BaseFlow     | 22-mars-2023 18:29:10 | FFMESH_000004.msh | 86         | 0.56244   
 55    | FFDATA_000055.txt | BaseFlow     | 22-mars-2023 18:29:49 | FFMESH_000004.msh | 88         | 0.55772   
 56    | FFDATA_000056.txt | BaseFlow     | 22-mars-2023 18:30:27 | FFMESH_000004.msh | 90         | 0.55314   
 57    | FFDATA_000057.txt | BaseFlow     | 22-mars-2023 18:31:06 | FFMESH_000004.msh | 92         | 0.54871   
 58    | FFDATA_000058.txt | BaseFlow     | 22-mars-2023 18:31:44 | FFMESH_000004.msh | 94         | 0.54442   
 59    | FFDATA_000059.txt | BaseFlow     | 22-mars-2023 18:32:24 | FFMESH_000004.msh | 96         | 0.54026   
 60    | FFDATA_000060.txt | BaseFlow     | 22-mars-2023 18:33:00 | FFMESH_000004.msh | 98         | 0.53622   
 61    | FFDATA_000061.txt | BaseFlow     | 22-mars-2023 18:33:35 | FFMESH_000004.msh | 100        | 0.53231   
 62    | FFDATA_000062.txt | BaseFlow     | 22-mars-2023 18:34:22 | FFMESH_000004.msh | 50         | 0.69002   
 63    | FFDATA_000063.txt | BaseFlow     | 22-mars-2023 19:03:54 | FFMESH_000004.msh | 46.6244    | 0.70906   
 
#################################################################################
 
.... CONTENT OF DIRECTORY ./WORK/EIGENMODES
     (couples of .txt/.ff2m files )
 
 Index | Name              | Type         | Date                 | Mesh file         | Re         | lambda          
 1     | FFDATA_000001.txt | EigenMode    | 22-mars-2023 18:01:33 | FFMESH_000003.msh | 60         | 0.047117+0.74492i
 2     | FFDATA_000002.txt | EigenMode    | 22-mars-2023 18:01:53 | FFMESH_000004.msh | 60         | 0.047043+0.74504i
 3     | FFDATA_000003.txt | EigenModeA   | 22-mars-2023 18:02:05 | FFMESH_000004.msh | 60         | 0.047043-0.74504i
 4     | FFDATA_000004.txt | EigenMode    | 22-mars-2023 18:03:50 | FFMESH_000004.msh | 60         | 0.047043+0.74504i
 5     | FFDATA_000005.txt | EigenMode    | 22-mars-2023 18:03:50 | FFMESH_000004.msh | 60         | -0.093808+0.74678i
 6     | FFDATA_000006.txt | EigenMode    | 22-mars-2023 18:03:50 | FFMESH_000004.msh | 60         | -0.092115+0.71194i
 7     | FFDATA_000007.txt | EigenMode    | 22-mars-2023 18:03:51 | FFMESH_000004.msh | 60         | -0.095964+0.7874i
 8     | FFDATA_000008.txt | EigenMode    | 22-mars-2023 18:03:51 | FFMESH_000004.msh | 60         | -0.08862+0.66856i
 9     | FFDATA_000009.txt | EigenMode    | 22-mars-2023 18:03:51 | FFMESH_000004.msh | 60         | -0.12749+0.71795i
 10    | FFDATA_000010.txt | EigenMode    | 22-mars-2023 18:03:52 | FFMESH_000004.msh | 60         | -0.12892+0.75994i
 11    | FFDATA_000011.txt | EigenMode    | 22-mars-2023 18:03:52 | FFMESH_000004.msh | 60         | -0.10515+0.83064i
 12    | FFDATA_000012.txt | EigenMode    | 22-mars-2023 18:03:52 | FFMESH_000004.msh | 60         | -0.12771+0.67634i
 13    | FFDATA_000013.txt | EigenMode    | 22-mars-2023 18:03:53 | FFMESH_000004.msh | 60         | -0.12915+0.80118i
 14    | FFDATA_000014.txt | EigenMode    | 22-mars-2023 18:03:53 | FFMESH_000004.msh | 60         | -0.14684+0.7365i
 15    | FFDATA_000015.txt | EigenMode    | 22-mars-2023 18:03:53 | FFMESH_000004.msh | 60         | -0.091068+0.62139i
 16    | FFDATA_000016.txt | EigenMode    | 22-mars-2023 18:03:54 | FFMESH_000004.msh | 60         | -0.14418+0.70003i
 17    | FFDATA_000017.txt | EigenMode    | 22-mars-2023 18:03:54 | FFMESH_000004.msh | 60         | -0.14912+0.75633i
 18    | FFDATA_000018.txt | EigenMode    | 22-mars-2023 18:03:54 | FFMESH_000004.msh | 60         | -0.14624+0.78064i
 19    | FFDATA_000019.txt | EigenMode    | 22-mars-2023 18:03:54 | FFMESH_000004.msh | 60         | -0.15824+0.71006i
 20    | FFDATA_000020.txt | EigenMode    | 22-mars-2023 18:03:55 | FFMESH_000004.msh | 60         | -0.14874+0.80288i
 21    | FFDATA_000021.txt | EigenMode    | 22-mars-2023 18:03:55 | FFMESH_000004.msh | 60         | -0.12642+0.63639i
 22    | FFDATA_000022.txt | EigenMode    | 22-mars-2023 18:03:55 | FFMESH_000004.msh | 60         | -0.12934+0.84348i
 23    | FFDATA_000023.txt | EigenMode    | 22-mars-2023 18:03:56 | FFMESH_000004.msh | 60         | -0.16617+0.74143i
 24    | FFDATA_000024.txt | EigenMode    | 22-mars-2023 18:14:49 | FFMESH_000004.msh | 40         | -0.030477+0.72027i
 25    | FFDATA_000025.txt | EigenMode    | 22-mars-2023 18:15:16 | FFMESH_000004.msh | 40         | -0.030477+0.72027i
 26    | FFDATA_000026.txt | EigenMode    | 22-mars-2023 18:15:52 | FFMESH_000004.msh | 42         | -0.020665+0.72472i
 27    | FFDATA_000027.txt | EigenMode    | 22-mars-2023 18:16:33 | FFMESH_000004.msh | 44         | -0.011401+0.72861i
 28    | FFDATA_000028.txt | EigenMode    | 22-mars-2023 18:17:11 | FFMESH_000004.msh | 46         | -0.0026448+0.732i
 29    | FFDATA_000029.txt | EigenMode    | 22-mars-2023 18:17:49 | FFMESH_000004.msh | 48         | 0.0056396+0.73493i
 30    | FFDATA_000030.txt | EigenMode    | 22-mars-2023 18:18:23 | FFMESH_000004.msh | 50         | 0.013485+0.73745i
 31    | FFDATA_000031.txt | EigenMode    | 22-mars-2023 18:18:58 | FFMESH_000004.msh | 52         | 0.020921+0.73959i
 32    | FFDATA_000032.txt | EigenMode    | 22-mars-2023 18:19:35 | FFMESH_000004.msh | 54         | 0.027975+0.74138i
 33    | FFDATA_000033.txt | EigenMode    | 22-mars-2023 18:20:12 | FFMESH_000004.msh | 56         | 0.034672+0.74286i
 34    | FFDATA_000034.txt | EigenMode    | 22-mars-2023 18:20:49 | FFMESH_000004.msh | 58         | 0.041033+0.74404i
 35    | FFDATA_000035.txt | EigenMode    | 22-mars-2023 18:21:26 | FFMESH_000004.msh | 60         | 0.047079+0.74495i
 36    | FFDATA_000036.txt | EigenMode    | 22-mars-2023 18:22:08 | FFMESH_000004.msh | 62         | 0.052829+0.7456i
 37    | FFDATA_000037.txt | EigenMode    | 22-mars-2023 18:22:43 | FFMESH_000004.msh | 64         | 0.058299+0.74602i
 38    | FFDATA_000038.txt | EigenMode    | 22-mars-2023 18:23:20 | FFMESH_000004.msh | 66         | 0.063506+0.74621i
 39    | FFDATA_000039.txt | EigenMode    | 22-mars-2023 18:23:55 | FFMESH_000004.msh | 68         | 0.068463+0.74619i
 40    | FFDATA_000040.txt | EigenMode    | 22-mars-2023 18:24:32 | FFMESH_000004.msh | 70         | 0.073184+0.74598i
 41    | FFDATA_000041.txt | EigenMode    | 22-mars-2023 18:25:09 | FFMESH_000004.msh | 72         | 0.077683+0.74558i
 42    | FFDATA_000042.txt | EigenMode    | 22-mars-2023 18:25:46 | FFMESH_000004.msh | 74         | 0.081969+0.74501i
 43    | FFDATA_000043.txt | EigenMode    | 22-mars-2023 18:26:24 | FFMESH_000004.msh | 76         | 0.086055+0.74428i
 44    | FFDATA_000044.txt | EigenMode    | 22-mars-2023 18:27:00 | FFMESH_000004.msh | 78         | 0.089951+0.74339i
 45    | FFDATA_000045.txt | EigenMode    | 22-mars-2023 18:27:35 | FFMESH_000004.msh | 80         | 0.093666+0.74235i
 46    | FFDATA_000046.txt | EigenMode    | 22-mars-2023 18:28:12 | FFMESH_000004.msh | 82         | 0.097211+0.74118i
 47    | FFDATA_000047.txt | EigenMode    | 22-mars-2023 18:28:46 | FFMESH_000004.msh | 84         | 0.10059+0.73987i
 48    | FFDATA_000048.txt | EigenMode    | 22-mars-2023 18:29:25 | FFMESH_000004.msh | 86         | 0.10382+0.73844i
 49    | FFDATA_000049.txt | EigenMode    | 22-mars-2023 18:30:03 | FFMESH_000004.msh | 88         | 0.1069+0.7369i  
 50    | FFDATA_000050.txt | EigenMode    | 22-mars-2023 18:30:41 | FFMESH_000004.msh | 90         | 0.10984+0.73524i
 51    | FFDATA_000051.txt | EigenMode    | 22-mars-2023 18:31:19 | FFMESH_000004.msh | 92         | 0.11265+0.73347i
 52    | FFDATA_000052.txt | EigenMode    | 22-mars-2023 18:31:58 | FFMESH_000004.msh | 94         | 0.11534+0.73161i
 53    | FFDATA_000053.txt | EigenMode    | 22-mars-2023 18:32:37 | FFMESH_000004.msh | 96         | 0.11791+0.72966i
 54    | FFDATA_000054.txt | EigenMode    | 22-mars-2023 18:33:13 | FFMESH_000004.msh | 98         | 0.12036+0.72761i
 55    | FFDATA_000055.txt | EigenMode    | 22-mars-2023 18:33:49 | FFMESH_000004.msh | 100        | 0.12271+0.72548i
 56    | FFDATA_000056.txt | EigenMode    | 22-mars-2023 18:34:37 | FFMESH_000004.msh | 50         | 0.013485+0.73745i
 57    | FFDATA_000057.txt | EigenModeD   | 22-mars-2023 19:03:56 | UNKNOWN.msh       | 46.6244    | 0+0.73295i      
 58    | FFDATA_000058.txt | EigenMode    | 22-mars-2023 19:04:11 | FFMESH_000004.msh | 46.6244    | -8.3191e-08+0.73295i
 59    | FFDATA_000059.txt | EigenModeA   | 22-mars-2023 19:04:29 | FFMESH_000004.msh | 46.6244    | -8.3191e-08-0.73295i
 60    | FFDATA_000060.txt | HarmonicMode | 22-mars-2023 19:05:05 | FFMESH_000004.msh | 47         | 0+0.73923i      
 61    | FFDATA_000061.txt | SecondHarmonicMode | 22-mars-2023 19:05:06 | FFMESH_000004.msh | 47         | 0+1.4785i       
 
 
 You have previously computed 58 eigenvalues in this session
 The full list will be available as field  'EIGENVALUES' of this function's result
 
#################################################################################
 
.... CONTENT OF DIRECTORY ./WORK/MEANFLOWS
     (couples of .txt/.ff2m files )
 
 Index | Name              | Type         | Date                 | Mesh file         | Re         | Fx        
 1     | FFDATA_000001.txt | MeanFlow     | 22-mars-2023 19:05:03 | FFMESH_000004.msh | 47         | 0.71189   
 2     | MeanFlow_Re47_Omegax0.txt | MeanFlow     | 22-mars-2023 19:12:45 | FFMESH_000004.msh | 47         | 0.71001   
 3     | HBMode1_Re47_Omegax0.txt | HBmode       | 22-mars-2023 19:12:47 | FFMESH_000004.msh | 47         | ???       
 4     | MeanFlow_Re47.5_Omegax0.txt | MeanFlow     | 22-mars-2023 19:16:36 | FFMESH_000004.msh | 47.5       | 0.70756   
 5     | HBMode1_Re47.5_Omegax0.txt | HBmode       | 22-mars-2023 19:16:37 | FFMESH_000004.msh | 47.5       | ???       
 6     | MeanFlow_Re48_Omegax0.txt | MeanFlow     | 22-mars-2023 19:20:25 | FFMESH_000004.msh | 48         | 0.70515   
 7     | HBMode1_Re48_Omegax0.txt | HBmode       | 22-mars-2023 19:20:26 | FFMESH_000004.msh | 48         | ???       
 8     | MeanFlow_Re49_Omegax0.txt | MeanFlow     | 22-mars-2023 19:24:09 | FFMESH_000004.msh | 49         | 0.70049   
 9     | HBMode1_Re49_Omegax0.txt | HBmode       | 22-mars-2023 19:24:11 | FFMESH_000004.msh | 49         | ???       
 10    | MeanFlow_Re50_Omegax0.txt | MeanFlow     | 22-mars-2023 19:28:10 | FFMESH_000004.msh | 50         | 0.69602   
 11    | HBMode1_Re50_Omegax0.txt | HBmode       | 22-mars-2023 19:28:11 | FFMESH_000004.msh | 50         | ???       
 12    | MeanFlow_Re52.5_Omegax0.txt | MeanFlow     | 22-mars-2023 19:32:16 | FFMESH_000004.msh | 52.5       | 0.68558   
 13    | HBMode1_Re52.5_Omegax0.txt | HBmode       | 22-mars-2023 19:32:17 | FFMESH_000004.msh | 52.5       | ???       
 14    | MeanFlow_Re55_Omegax0.txt | MeanFlow     | 22-mars-2023 19:36:08 | FFMESH_000004.msh | 55         | 0.67609   
 15    | HBMode1_Re55_Omegax0.txt | HBmode       | 22-mars-2023 19:36:10 | FFMESH_000004.msh | 55         | ???       
 16    | MeanFlow_Re60_Omegax0.txt | MeanFlow     | 22-mars-2023 19:39:45 | FFMESH_000004.msh | 60         | 0.65949   
 17    | HBMode1_Re60_Omegax0.txt | HBmode       | 22-mars-2023 19:39:47 | FFMESH_000004.msh | 60         | ???       
 18    | MeanFlow_Re65_Omegax0.txt | MeanFlow     | 22-mars-2023 19:43:43 | FFMESH_000004.msh | 65         | 0.64543   
 19    | HBMode1_Re65_Omegax0.txt | HBmode       | 22-mars-2023 19:43:45 | FFMESH_000004.msh | 65         | ???       
 20    | MeanFlow_Re70_Omegax0.txt | MeanFlow     | 22-mars-2023 19:47:43 | FFMESH_000004.msh | 70         | 0.63339   
 21    | HBMode1_Re70_Omegax0.txt | HBmode       | 22-mars-2023 19:47:44 | FFMESH_000004.msh | 70         | ???       
 22    | MeanFlow_Re75_Omegax0.txt | MeanFlow     | 22-mars-2023 19:51:23 | FFMESH_000004.msh | 75         | 0.62295   
 23    | HBMode1_Re75_Omegax0.txt | HBmode       | 22-mars-2023 19:51:24 | FFMESH_000004.msh | 75         | ???       
 24    | MeanFlow_Re80_Omegax0.txt | MeanFlow     | 22-mars-2023 19:55:21 | FFMESH_000004.msh | 80         | 0.61382   
 25    | HBMode1_Re80_Omegax0.txt | HBmode       | 22-mars-2023 19:55:22 | FFMESH_000004.msh | 80         | ???       
 26    | MeanFlow_Re85_Omegax0.txt | MeanFlow     | 22-mars-2023 19:59:10 | FFMESH_000004.msh | 85         | 0.60576   
 27    | HBMode1_Re85_Omegax0.txt | HBmode       | 22-mars-2023 19:59:11 | FFMESH_000004.msh | 85         | ???       
 28    | MeanFlow_Re90_Omegax0.txt | MeanFlow     | 22-mars-2023 20:03:10 | FFMESH_000004.msh | 90         | 0.59859   
 29    | HBMode1_Re90_Omegax0.txt | HBmode       | 22-mars-2023 20:03:12 | FFMESH_000004.msh | 90         | ???       
 30    | MeanFlow_Re95_Omegax0.txt | MeanFlow     | 22-mars-2023 20:06:59 | FFMESH_000004.msh | 95         | 0.59218   
 31    | HBMode1_Re95_Omegax0.txt | HBmode       | 22-mars-2023 20:07:00 | FFMESH_000004.msh | 95         | ???       
 32    | MeanFlow_Re100_Omegax0.txt | MeanFlow     | 22-mars-2023 20:10:56 | FFMESH_000004.msh | 100        | 0.58642   
 33    | HBMode1_Re100_Omegax0.txt | HBmode       | 22-mars-2023 20:10:58 | FFMESH_000004.msh | 100        | ???       
 
#################################################################################
 
.... CONTENT OF DIRECTORY ./WORK/MESHBF
     (couples of .txt/.ff2m files )
 
 Index | Name              | Type         | Date                 | Mesh file         | Re         | Fx        
 1     | FFDATA_000001.txt | BaseFlow     | 22-mars-2023 18:00:57 | FFMESH_000001.msh | 1          | 5.5507    
 
#################################################################################
 
.... CONTENT OF DIRECTORY ./WORK/MISC
     (couples of .txt/.ff2m files )
   NB this directory may contain various secondary files produced by StabFem,
      such as flowfields projected after adaptation but not recomputed, adaptation masks, etc... 
 
 Index | Name              | Type         | Date                 | (no bf or mesh indication)
 1     | FFDATA_000001.txt | Sensitivity  | 22-mars-2023 18:02:12 | (unknown)       
 2     | FFDATA_000002.txt | WNL          | 22-mars-2023 19:05:06 | (unknown)       
 
 
#################################################################################

ans = 

  struct with fields:

         MESHES: [1x4 struct]
      BASEFLOWS: [1x63 struct]
     EIGENMODES: [1x61 struct]
    EIGENVALUES: [1x1 struct]
      MEANFLOWS: [1x33 struct]
         MESHBF: [1x1 struct]
           MISC: [1x2 struct]