Contents

Equilibrium shape and oscillation modes of a liquid bridge resulting from the coalescence of two droplers

This script reproduces the main results from the paper by Chireux et al. (Phys. fluids, 2015)

The script reproduces figures 5, 11, 12 of the paper plus a collection of other results

close all;
addpath([fileparts(fileparts(pwd)),'/SOURCES_MATLAB']);
SF_Start('verbosity',2,'storagemode',1);
SF_core_setopt('eigensolver','SLEPC'); % this set of programs uses arpack ; to be rewritten someday with Slepc
figureformat = 'png';
mkdir('FIGURES');
NOTICE  - Initialization already done

CHAPTER 1 : creation of initial mesh for cylindrical bridge, L=4

L = 4;
density=20;
 %creation of an initial mesh (with volume corresponding to coalescence of two spherical caps)
freesurf = SF_Init('MeshInit_Bridge.edp','Params',[0 L density],'problemtype','3DFreeSurfaceStatic','cleanworkdir','yes');

V = pi*L/2*(1+L^2/12);
gamma = 1;
rhog = 0;
freesurf = SF_Deform(freesurf,'V',V,'typestart','pined','typeend','pined','gamma',gamma,'rhog',rhog);

%Vol0 = freesurf.Vol;

CHAPTER 2 : Eigenvalue computation for m=0 and m=1 for L=4

[evm0,emm0] =  SF_Stability(freesurf,'nev',10,'m',0,'sort','SIA');
[evm1,emm1] =  SF_Stability(freesurf,'nev',10,'m',1,'sort','SIA');

% PLOT RESULTS
figure(1);
E=0.15;
SF_Plot(emm0(5),'phi.im','title',{'Mode m=0,s',['\omega_r = ',num2str(-imag(evm0(5)))] } );hold on;
SF_Plot_ETA(emm0(5),'Amp',E,'style','r');hold off;

figure(4);
E=0.15;
subplot(1,4,1);
SF_Plot(emm0(3),'phi.im','title',{'Mode m=0,a',['\omega_r = ',num2str(-imag(evm0(3)))] } );
hold on;SF_Plot_ETA(emm0(3),'Amp',E,'style','r');hold off;
subplot(1,4,2);
SF_Plot(emm0(5),'phi.im','title',{'Mode m=0,s',['\omega_r = ',num2str(-imag(evm0(5)))] } );hold on;
SF_Plot_ETA(emm0(5),'Amp',E,'style','r');hold off;
subplot(1,4,3);
SF_Plot(emm1(1),'phi.im','title',{'Mode m=1,s',['\omega_r = ',num2str(-imag(evm1(1)))] } );hold on;
SF_Plot_ETA(emm1(1),'Amp',E,'style','r');hold off;
subplot(1,4,4);
SF_Plot(emm1(3),'phi.im','title',{'Mode m=1,a',['\omega_r = ',num2str(-imag(evm1(3)))] } );hold on;
SF_Plot_ETA(emm1(3),'Amp',E,'style','r');hold off;

pos = get(gcf,'Position'); pos(3)=pos(4)*2.6;set(gcf,'Position',pos); % resize aspect ratio
saveas(gcf,'FIGURES/Bridges_NV_Eigenmodes_phi_cyl_L3_5',figureformat);


figure(2);
plot(imag(evm0),real(evm0),'ro',imag(evm1),real(evm1),'bo');
title('Cylindrical bridge, L= 4 : spectra for m=0 (red) and m=1 (blue)');
xlabel('\omega_r');ylabel('\omega_i');

figure(3); hold off;
plot(freesurf.mesh.S0,real(emm0(3).eta));hold on;
plot(freesurf.mesh.S0,real(emm0(5).eta));
% note that for m=0 the two first modes are spurious, so we take modes 3 and 5
plot(freesurf.mesh.S0,real(emm1(1).eta));
plot(freesurf.mesh.S0,real(emm1(3).eta));
% on the other hand for m=1 the first modes are regular
title('Cylindrical bridge, L= 4 : structure of the four simplest modes eta(s) ');
xlabel('s');ylabel('\eta(s)');
legend('m=0,a','m=0,s','m=1,s','m=1,a');

pause(1);

CHAPTER 3 : Loop over L/R

% CHAPTER 3a : First loop in the interval [3.5,7] (increasing values)
L = 3.5;
freesurf = SF_Init('MeshInit_Bridge.edp','Params',[0 L density],'problemtype','3DFreeSurfaceStatic'); %% creation of an initial mesh (cylindrical liquid bridge)
tabL = [2:.1:3.5 ,3.5 :.25 :7];
for i=17:31
    L = tabL(i);Lans = tabL(i-1);
    freesurf = SF_MeshStretch(freesurf,'Yratio',L/Lans);
    V = pi*L/2*(1+L^2/12);
    freesurf = SF_Deform(freesurf,'V',V,'rhog',0);
    tabV(i) = freesurf.Vol; tabP(i) = freesurf.P0;
    tabR(i) = freesurf.mesh.rsurf(round(freesurf.mesh.nsurf+1)/2);
    figure(30);
    plot([freesurf.mesh.rsurf; NaN; -freesurf.mesh.rsurf],[freesurf.mesh.zsurf; NaN; freesurf.mesh.zsurf]); hold on;
    pause(0.1);
    evm0 =  SF_Stability(freesurf,'nev',10,'m',0,'sort','SIA');
    evm1 =  SF_Stability(freesurf,'nev',10,'m',1,'sort','SIA');
    tabEVm0(:,i) = evm0;
    tabEVm1(:,i) = evm1;
end

% CHAPTER 3b : Second loop in the interval [3.5,2] (decreasing values)
L = 3.5;
freesurf = SF_Init('MeshInit_Bridge.edp','Params',[0 L density],'problemtype','3DFreeSurfaceStatic'); %% creation of an initial mesh (cylindrical liquid bridge)
for i=16:-1:1
    L = tabL(i);Lans = tabL(i+1);
    freesurf = SF_MeshStretch(freesurf,'Yratio',L/Lans);
    V = pi*L/2*(1+L^2/12);
    freesurf = SF_Deform(freesurf,'V',V,'rhog',0);
    tabV(i) = freesurf.Vol; tabP(i) = freesurf.P0;
    tabR(i) = freesurf.mesh.rsurf(round(freesurf.mesh.nsurf+1)/2);
    figure(30);
    plot([freesurf.mesh.rsurf; NaN; -freesurf.mesh.rsurf],[freesurf.mesh.zsurf; NaN; freesurf.mesh.zsurf]); hold on;
    pause(0.1);
    evm0 =  SF_Stability(freesurf,'nev',10,'m',0,'sort','SIA');
    evm1 =  SF_Stability(freesurf,'nev',10,'m',1,'sort','SIA');
    tabEVm0(:,i) = evm0;
    tabEVm1(:,i) = evm1;
end

Chapter 3c : PLOTS for equilibrium shape

figure(30);
title('A few equilibrium shapes for liquid bridges');
axis equal;xlabel('r');ylabel('z');

% PLOTS
figure(31);
hold on;
subplot(1,2,1); plot(tabL,tabV./(pi*tabL),'bo-');
xlabel('L^*');ylabel('V^*'); title('Volume vs. Length');
subplot(1,2,2); plot(tabL,tabR,'bo-');
xlabel('L^*');ylabel('R_m/R'); title('Bridge radius vs. Length')

Chapter 3d : PLOTS for oscillation frequencies

figure(32);hold on;
for num=1:6
    plot(tabL,imag(tabEVm0(num,:)),'ro-',tabL,imag(tabEVm1(num,:)),'bo-');
end
title('frequencies of m=0 (red) and m=1 (blue) modes vs. L');
xlabel('L');ylabel('\omega_r');
ylim([0 6]);
box on; pos = get(gcf,'Position'); pos(4)=pos(3)*1;set(gcf,'Position',pos); % resize aspect ratio
set(gca,'FontSize', 14);
saveas(gca,'FIGURES/Bridges_NV_coal_omega',figureformat);

figure(33);hold on;
for num=1:6
    plot(tabL,imag(tabEVm0(num,:)).*tabL.^1.5,'ro-',tabL,imag(tabEVm1(num,:)).*tabL.^1.5,'bo-');
end
title('rescaled frequencies of m=0 (red) and m=1 (blue) modes vs. L');
xlabel('L');ylabel('\omega_r/\omega_L');
ylim([0 22]);
box on; pos = get(gcf,'Position'); pos(4)=pos(3)*1;set(gcf,'Position',pos); % resize aspect ratio
set(gca,'FontSize', 14);
saveas(gca,'FIGURES/Bridges_NV_L3_5_coal_omegaL',figureformat);

un petit dernier...

[evm1,emm1] =  SF_Stability(freesurf,'nev',10,'m',1,'sort','SIA')
evm1 =

   0.0000 - 1.8572i
   0.0000 + 1.8572i
  -0.0000 - 3.9412i
   0.0000 + 3.9412i
  -0.0000 - 7.3157i
   0.0000 + 7.3157i
   0.0000 -11.9115i
  -0.0000 +11.9115i
   0.0000 -17.3896i
  -0.0000 +17.3896i


emm1 = 

  1x10 struct array with fields:

    mesh
    filename
    DataDescription
    datatype
    datastoragemode
    datadescriptors
    meshfilename
    baseflowfilename
    lambda
    INDEXING
    m
    ur
    uz
    ut
    eta
    DetaDs
    phi
    p
    type

figure(1);
E=0.15;
SF_Plot(emm1(2),'phi.im','title',{'Mode m=0,s',['\omega_r = ',num2str(-imag(evm1(2)))] } );hold on;
SF_Plot_ETA(emm0(5),'Amp',E,'style','r');hold off;


% [[PUBLISH]] (this tag is to enable automatic publication as html ; don't touch this line)
% [[FIGURES]] (this tag is to generate figures for the manual)
%
fclose(fopen('SCRIPT_LiquidBridges.success','w+'));