Patterns and Their Large-Scale Distortions in Marangoni Convection with Insoluble Surfactant

Abstract

Nonlinear dynamics of patterns near the threshold of long-wave monotonic Marangoni instability of conductive state in a heated thin layer of liquid covered by insoluble surfactant is considered. Pattern selection between roll and square planforms is analyzed. The dependence of pattern stability on the heat transfer from the free surface of the liquid characterized by Biot number and the gravity described by Galileo number at different surfactant concentrations is studied. Using weakly nonlinear analysis, we derive a set of amplitude equations governing the large-scale roll distortions in the presence of the surface deformation and the surfactant redistribution. These equations are used for the linear analysis of modulational instability of stationary rolls.

1. Introduction

Marangoni convection is a type of non-equilibrium process, which creates a variety of spatiotemporal periodic patterns, see [1,2], and is crucial in thin liquid layers, where the interfacial effects prevail over the bulk effects. The formation of large-scale convective patterns is a result of long-wavelength instability.

There are two kinds of longwave Marangoni instability, (i) at small Biot number with temperature as the active, slowly evolving, variable [3], and (ii) at moderate Galileo number with surface deformation as the active variable [4]. The recent research of Shklyaev et al. [5], shows that in the interval of wave numbers $k=O(B{i}^{1/2}),Bi\ll 1$ ($Bi$ is the Biot number defined later), both those variables are active. In that interval, large-scale monotonic and oscillatory instability modes exist, which produce stationary patterns and wave patterns, correspondingly [6].

It is known that when the free surface of the liquid is covered by surface-active agent (surfactant), the formation of the convective patterns significantly changes. The basic experiments describing this effect can be found in [7,8], theoretically the behavior of large-scale patterns under the influence of adsorbed insoluble surfactant is considered in [9]. When the surfactant is added, the linear stability boundaries, the types of bifurcation and selected patterns depend on the elasticity number and on the Biot number, see [10].

In the present work we investigate the large-scale Marangoni convection in a liquid layer with insoluble surfactant spread over a deformable free surface, in the interval of wavenumbers $O(B{i}^{1/2})$. We carry out a weakly nonlinear analysis of patterns on a square lattice in the Fourier space, which include rolls and squares. We consider the modulation of rolls by long-wave disturbances near the instability threshold and derive the system of the amplitude equations. The periodic patterns can be unstable with respect to a spatial modulation that breaks the spatial periodicity. The modulation of convective patterns in the neighborhood of the instability threshold is explained by Newell, Whitehead, and Segel [11,12] as interaction of disturbances with different wavenumbers close to the critical one. In the case of stationary rolls, which are chosen for the investigation in our paper, the asymptotic multiscale analysis leads to the Ginzburg-Landau equation with real coefficients. However, in the presence of slowly decaying longwave modes related to the conservation laws, that equation turns out to be insufficient. In [13] it was shown that interaction with such a mode can significantly influence the stability of patterns. In the problem under consideration, there are two conservation laws, the conservation of the liquid volume and the conservation of the total amount of the surfactant. The corresponding modification of the stability properties of rolls is considered in the present paper. The analysis of that phenomenon has not been previously carried out for Marangoni patterns in the presence of a surfactant.

2. Statement of the Problem

We consider an infinite horizontal layer with a mean thickness $d_0$, thermal diffusivity $\chi$, kinematic viscosity $\nu$, density $\rho$, dynamic viscosity $\eta =\rho \nu$, and thermal conductivity ${\mathsf{\Lambda }}_T$. The layer is confined between a rigid substrate and deformable free upper boundary. The vertical axis $z$ is directed upward. The layer is heated from below with transverse temperature gradient $-a(a>0).$ At the upper surface the insoluble surfactant is absorbed. It is convected by interfacial velocity field and diffuses over the interface but not into the bulk. ${\mathsf{\Gamma }}_0$ is the reference value of the surfactant concentration. The surface tension is a linear function of temperature and surfactant concentration, i.e., $\sigma ={\sigma }_0-{\sigma }_1(T-T_0)-{\sigma }_2(\mathsf{\Gamma }-{\mathsf{\Gamma }}_0)$ (${\sigma }_0$ is the reference value of surface tension, ${\sigma }_1=-{\partial }_T\sigma ,$ ${\sigma }_2=-{\partial }_{\mathsf{\Gamma }}\sigma$).

In [10] we derived the system of long-wave nondimensional amplitude equations for the local film thickness $H(X,Y,\tau )$, the surfactant concentration $\mathsf{\Gamma }(X,Y,\tau )$ and perturbation of the temperature disturbance $F(X,Y,\tau )$:

(1)${\partial }_{\tau }H=\nabla \cdot \left(\frac{H^3}{3}\nabla R+\frac{M{H}^2}{2}\nabla \theta +\frac{N{H}^2}{2}\nabla \mathsf{\Gamma }\right)\equiv \nabla \cdot Q_1,$

(2)$H{\partial }_{\tau }F=\nabla \cdot \left(\frac{H^4}{8}\nabla R+\frac{M{H}^3}{6}\nabla \theta +\frac{N{H}^3}{6}\nabla \mathsf{\Gamma }+H\nabla F\right)+Q_1\cdot \nabla \theta -\frac{1}{2}{\left(\nabla H\right)}^2-\beta \theta ,$

(3)${\partial }_{\tau }\mathsf{\Gamma }=\nabla \cdot \left[\mathsf{\Gamma }H\left(\frac{H}{2}\nabla R+M\nabla \theta +N\nabla \mathsf{\Gamma }\right)+L\nabla \mathsf{\Gamma }\right]\equiv \nabla \cdot Q_2.$

Here $R=GH-S{\nabla }^{2}H$ is the pressure disturbance and $\theta =F-H$ is the temperature distribution at the free surface, $\nabla =({\partial }_X,{\partial }_Y),$$X,Y$ and $\tau$ are horizontal coordinates and time is rescaled as

$X=\epsilon x,Y=\epsilon y,\tau ={\epsilon }^{2}t,0<\epsilon \ll 1.$

This rescaling of coordinates leads to the rescaling of the wavenumber, $k=\epsilon K$. The parameter $\epsilon \ll 1$ is a ratio between the mean thickness $d_0$ and typical horizontal scale of disturbances.

The problem includes the following dimensionless parameters: $M=\raisebox{1ex}{${\sigma }_{1}a{d}_0^2$}\!\left/ \!\raisebox{-1ex}{$\eta \chi $}\right.$ is the Marangoni number, $N=\raisebox{1ex}{${\sigma }_2{d}_0{\mathsf{\Gamma }}_0$}\!\left/ \!\raisebox{-1ex}{$\eta \chi $}\right.$ is the elasticity number, $L=\raisebox{1ex}{$D_0$}\!\left/ \!\raisebox{-1ex}{$\chi $}\right.$ is the Lewis number ($D_0$ is the surfactant diffusivity), $G=\raisebox{1ex}{$g{d}_0^2$}\!\left/ \!\raisebox{-1ex}{$\nu \chi $}\right.$ is the Galileo number ($g$ is the acceleration due to gravity), $\mathsf{\Sigma }=\raisebox{1ex}{${\sigma }_0{d}_0$}\!\left/ \!\raisebox{-1ex}{$\eta \chi $}\right.$ is the inverse capillary number, $Bi=\raisebox{1ex}{$q{d}_0$}\!\left/ \!\raisebox{-1ex}{${\mathsf{\Lambda }}_T$}\right.$ is the Biot number ($q$ is the heat transfer coefficient). Here we assume $\mathsf{\Sigma }={\epsilon }^{-2}S,S=O(1)$ and $Bi={\epsilon }^2\beta .$

The linear stability analysis of the motionless state was performed in [10]. There are two instability modes: the monotonic one and the oscillatory one. The instability boundaries are determined by Equation (5) in [10] for the growth rate $\mathsf{\Lambda }(M)$ in the form:

(4)$$\begin{gathered} {\mathsf{\Lambda }}^3+c_1{\mathsf{\Lambda }}^2+c_2\mathsf{\Lambda }+c_3=0, \\ {c}_1=\left(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$3$}\right.\right)K^2\left(3+G+3L-M+3N+K^{2}S\right)+\beta , \\ \begin{array}{l}{c}_2=\left(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$144$}\right.\right)(K^4\{-24\left[2L\left(M-3\right)+3\left(M-2N\right)\right]+\\ +\left(G+K^{2}S\right)\left(48+48L-M+12N\right)\}+48K^2\left[G+3\left(L+N\right)+K^{2}S\right]\beta ),\end{array} \\ \begin{array}{l}{c}_3=\left(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$144$}\right.\right)(K^6\left\{-72LM+\left(G+K^{2}S)\left[12N-L\left(M-48\right)\right]\right)\right\}+\\ +12K^4\left(4L+N\right)\left(G+K^{2}S\right)\beta ).\end{array} \end{gathered}$$

For the monotonic mode the neutral stability curve is described by the formula

(5)$M_m(K)=\frac{\left(48+\raisebox{1ex}{$12N$}\!\left/ \!\raisebox{-1ex}{$L$}\right.\right)\left(G+K^{2}S\right)\left(\beta +K^2\right)}{K^2\left(72+G+K^{2}S\right)}.$

The critical wave number that corresponds to the minimum of the marginal stability curve of the monotonic mode does not depend on the surfactant parameters and equals the value determined in [5]:

(6)${\left(K_c^{(m)}\right)}^2=\frac{GS\beta +\sqrt{72GS\beta \left(G+72-S\beta \right)}}{S\left(72-S\beta \right)}.$

The threshold value of Marangoni number for the oscillatory mode was also found in our previous work [10]. Shklyaev et al. in [5] show that without surfactant ($N=0$) the oscillatory instability mode appears for sufficiently large values of $\beta$ in a limited range of wave numbers. The oscillatory mode is critical under conditions $\beta S>17.4$ and $G<17.2.$

In order to analyze the competition between instability modes, we use the following property of the oscillatory instability: its appearance changes the meaning of the neutral curve (5). If there is no oscillatory instability, the Equation (5) indeed determines the monotonic instability boundary: $\mathsf{\Lambda }(M,K)<0$ at $M<M_m(K)$ and $\mathsf{\Lambda }(M,K)>0$ at $M>M_m(K)$, therefore ${\left(\raisebox{1ex}{$\partial \mathsf{\Lambda }(M,K)$}\!\left/ \!\raisebox{-1ex}{$\partial M$}\right.\right)|}_{M=M_m(K)}>0$. However, if there exists an oscillatory instability with the instability boundary $M=M_{osc}(K)<M_m(K)$, Equation (5) determines the stabilization boundary for one of the two unstable monotonic modes, thus ${\left(\raisebox{1ex}{$\partial \mathsf{\Lambda }(M,K)$}\!\left/ \!\raisebox{-1ex}{$\partial M$}\right.\right)|}_{M=M_m(K)}<0$ (see Figure 5.1 (b) in [14]).

Differentiating (4) with respect to $M$, we find that at the critical value of monotonic Marangoni number

(7)${\frac{\partial \mathsf{\Lambda }}{\partial M}|}_{M=M_m}=\frac{K^4{L}^2(72+G+K^{2}S)}{{\mathrm{D}}_1}$

where

(8)$\begin{array}{l}{\mathrm{D}}_1=12\{K^6\left(4L^2-N+LN\right)S^2-(36\left[5G-24\right]L^2+G\left[G+72\right]N+36\left[G-24\right]LN)\beta \\ +K^{4}S\left(4\left[2G+27\right]L^2+2\left[G+18\right]LN-N\left[72+2G+S\beta \right]\right)\\ +K^2(4L^2\left[216+27G+G^2-45S\beta \right]+LN\left[864+36G+G^2-36S\beta \right]\\ -N\left[72G+G^2+72S\beta +2GS\beta \right])\}.\end{array}$

One can show that even if the oscillatory instability is absent at $N=0$, it can appear for $N>N_m$, where the value of the elasticity number $N_m=N_m(G,\beta ,K)$ is determined by the relation ${\mathrm{D}}_1=0$.

Without loss of generality we can fix the value of the inverse capillary number equal to unity, i.e., $S=1$, which corresponds to the choice $\epsilon \equiv \raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$\sqrt{\mathsf{\Sigma }}$}\right.$. While the elasticity number, $N$, can vary in a relatively large interval, the value of the surfactant Lewis number, which is the ratio of mass diffusivity to thermal diffusivity, has significant limitation. As an example of insoluble surfactant we consider $C_{12}E{O}_8$ produced by Nikko Chemical Co. (Tokyo, Japan), for its properties see [15,16]. For our calculation in this paper the Lewis number is taken equal to $L=0.003$. The value of the wave number is taken equal to $K_c^{(m)}$. Then we obtain $N_m$ as a function on parameters $G$ and $\beta$, see Figure 1. For $N<N_m$ the critical Marangoni number $M_m$ is a threshold of the monotonic instability of the equilibrium state.

3. Weakly Nonlinear Analysis of Monotonic Mode

We study the nonlinear dynamics of perturbations in the neighborhood of $M_m$, $M=M_m+{\delta }^2{M}_2$ ($\delta$ is a small parameter of supercriticality). Here we consider spatially periodic patterns corresponding to a square lattice in the Fourier space.

We expand functions $H,F,\mathsf{\Gamma }$ in powers of the small parameter of supercriticality $\delta$,

(9)$H=1+\delta h_1+{\delta }^2{h}_2+{\delta }^3{h}_3+\dots ,$

(10)$F=1+\delta f_1+{\delta }^2{f}_2+{\delta }^3{f}_3+\dots ,$

(11)$\mathsf{\Gamma }=1+\delta {\gamma }_1+{\delta }^2{\gamma }_2+{\delta }^3{\gamma }_3+\dots ,$

and rescale the time variable

(12)${\tau }_2={\delta }^2\tau .$

Substituting these expansions into (1)–(3) and collecting terms with the same order of parameter $\delta$, we obtain the following system of equations at the leading order:

(13)$\left(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$3$}\right.\right){\nabla }^2{R}_1+\left(\raisebox{1ex}{$M_m$}\!\left/ \!\raisebox{-1ex}{$2$}\right.\right){\nabla }^2(f_1-h_1)+\left(\raisebox{1ex}{$N$}\!\left/ \!\raisebox{-1ex}{$2$}\right.\right){\nabla }^2{\gamma }_1=0,$

(14)$\left(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$8$}\right.\right){\nabla }^2{R}_1+\left(\raisebox{1ex}{$M_m$}\!\left/ \!\raisebox{-1ex}{$6$}\right.\right){\nabla }^2(f_1-h_1)+\left(\raisebox{1ex}{$N$}\!\left/ \!\raisebox{-1ex}{$6$}\right.\right){\nabla }^2{\gamma }_1+{\nabla }^2{f}_1-\beta \left(f_1-h_1\right)=0,$

(15)$\left(\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$2$}\right.\right){\nabla }^2{R}_1+M_m{\nabla }^2\left(f_1-h_1\right)+(N+L){\nabla }^2{\gamma }_1=0.$

Here the Marangoni number $M_m$ is taken from (5), $R_1=G{h}_1-S{\nabla }^2{h}_1.$ The solution of this system is presented in the form of squares and rolls:

(16)$h_1=A_1\left({\tau }_2\right)e^{iKX}+A_2\left({\tau }_2\right)e^{iKY}+c.c.,$

(17)$f_1={\alpha }_1\left(A_1\left({\tau }_2\right)e^{iKX}+A_2\left({\tau }_2\right)e^{iKY}\right)+c.c.,$

(18)${\gamma }_1={\alpha }_2\left(A_1\left({\tau }_2\right)e^{iKX}+A_2\left({\tau }_2\right)e^{iKY}\right)+c.c.,$

where ${\alpha }_1=-\raisebox{1ex}{$\left(G{K}^2+K^{4}S-72\beta \right)$}\!\left/ \!\raisebox{-1ex}{$72\left(K^2+\beta \right)$}\right.$ and ${\alpha }_2=\raisebox{1ex}{$\left(G+K^{2}S\right)$}\!\left/ \!\raisebox{-1ex}{$6L$}\right.$.

The solution at the second order of the expansion is presented in the form:

(19)$h_2=C_1\left({\tau }_2\right)e^{2iKX}+C_2\left({\tau }_2\right)e^{2iKY}+C_3\left({\tau }_2\right)e^{iK(X-Y)}+C_4\left(\tau \right)e^{iK(X+Y)}+c.c.,$

(20)$f_2=F_0+F_1\left({\tau }_2\right)e^{2iKX}+F_2\left({\tau }_2\right)e^{2iKY}+F_3\left({\tau }_2\right)e^{iK(X-Y)}+F_4\left(\tau \right)e^{iK(X+Y)}+c.c.,$

(21)${\gamma }_2=D_1\left({\tau }_2\right)e^{2iKX}+D_2\left({\tau }_2\right)e^{2iKY}+D_3\left({\tau }_2\right)e^{iK(X-Y)}+D_4\left(\tau \right)e^{iK(X+Y)}+c.c..$

Coefficients in (19)–(21) can be determined from the system of equations at the second order. The elimination of the secular terms in the system of the third order using the solvability condition yields a set of the Landau equations that govern the evolution of the complex amplitudes $A_1$ and $A_2$

(22)$\frac{\partial A_1}{\partial {\tau }_2}={\kappa }_0{A}_1+{\kappa }_1{\left|A_1\right|}^2{A}_1+{\kappa }_2{\left|A_2\right|}^2{A}_1,$

(23)$\frac{\partial A_2}{\partial {\tau }_2}={\kappa }_0{A}_2+{\kappa }_1{\left|A_2\right|}^2{A}_2+{\kappa }_2{\left|A_1\right|}^2{A}_2$

Here

(24)${\kappa }_0=\frac{K^4{L}^2{\left(72+G+K^2{S}^2\right)}^2{M}_2}{12\left(K^2{\mathrm{D}}_2-\beta {\mathrm{D}}_3\right)},$

where

$\begin{array}{l}{\mathrm{D}}_2=-N\left(G+K^{2}S\right)\left(72+G+K^{2}S\right)+LN\left[864+36G+G^2+2\left(18+G\right)K^{2}S+K^4{S}^2\right]\\ +4L^2\left[216+27G+G^2+\left(27+2G\right)K^{2}S+K^4{S}^2\right],\end{array}$

${\mathrm{D}}_3=36LN\left(-24+G+K^{2}S\right)+N\left(G+K^{2}S\right)\left(72+G+K^{2}S\right)+36L^2\left(-24+5G+5K^{2}S\right).$

In the case without surfactant $\left(N=0\right)$ the expression for ${\kappa }_0$ is equal to that determined by Equation (43) in Ref. [17]. The self-interaction coefficient, ${\kappa }_1$, as well as the cross-interaction one, ${\kappa }_2$, are real.

We have two kinds of steady state solutions in addition to the trivial solution ($A_1=A_2=0$) corresponding to the motionless state. They describe the rolls (if one of $A_1$ or $A_2$ is zero) and squares (both $A_1$ and $A_2$ are non-zero). As usual, the pattern selection is determined by the signs of ${\kappa }_1,{\kappa }_1-{\kappa }_2,{\kappa }_1+{\kappa }_2$. The line ${\kappa }_1=0$ divides the region of parameters into regions of subcritical (${\kappa }_1>0$) and supercritical (${\kappa }_1<0$) bifurcations of the roll patterns. The structures can be stable in the supercritical region if ${\kappa }_1$ and ${\kappa }_1+{\kappa }_2$ are negative. Rolls are selected if ${\kappa }_1-{\kappa }_2>0$, in the opposite case the squares are stable.

We consider the influence of the concentration of the surfactant on the pattern selection at very low values of the elasticity number. As shown in Figure 2, even low concentration of the surfactant remarkably changes the pattern selection of the structures on the square lattice in the Fourier space. There are two regions of stationary pattern stability, one shown in Figure 2a and another one, for small values of $\beta$, in Figure 2b. Here the set #1 (red colored) belongs to the case without surfactant. The set #2 (blue colored) represents the case with $N={10}^{-6}$ and the set #3 (black colored) corresponds to the case $N={10}^{-5}$. The marks on the Figure show the supercritical stability regions, and they include the number corresponding to the set number and the letter, “R” for stable rolls and “S” for stable squares. Here the solid lines are boundaries between supercritical and subcritical bifurcations for rolls (${\kappa }_1=0$). The dashed lines are similar boundaries for squares (${\kappa }_1+{\kappa }_2=0$), and the dot-dashed lines, which correspond to the relation ${\kappa }_1={\kappa }_2$, separate the regions of the stable rolls and stable squares. Figure 2b shows domains at small $\beta$. Here the rolls are mainly selected, for greater values of $\beta$ the squares are selected. Even a low concentration of the surfactant diminishes the stability domains. This effect is enhanced with the growth of the elasticity number.

Therefore, it was found that the most typical form of patterns is the roll pattern. Note that in a rectangular container with one size much longer than another one or in a narrow annular container of large radius, that pattern is the only possible kind of patterns.

In the next section, we investigate the instability of roll patterns with respect to one-dimensional disturbances breaking the spatial periodicity of patterns (Eckhaus instability).

4. Modulation of Roll Patterns

4.1. Amplitude Equations

We consider the roll patterns near the instability threshold,

(25)$M=M_m+{\delta }^2{M}_2,M_m=M\left(K_c\right).$

In that region, the motionless state is unstable with respect to disturbances in the interval of wavenumbers of the width $O(\delta )$ around $K=K_c$, and their growth rate is $O\left({\delta }^2\right)$.

Following the approach of Newell, Whitehead, and Segel [11,12], we introduce two spatial scales: (i) the scale $X_0=X$ corresponding to the basic wavelength of patterns, and (ii) the scale $X_1=\delta X$ corresponding to the spatial modulation of patterns due to the superposition of Fourier components within and around the instability region. Also, we rescale the time variable: ${\tau }_2={\delta }^2\tau .$

The differential operators are transformed in the following way

${\partial }_X={\partial }_{X_0}+\delta {\partial }_{X_1},{\partial }_{XX}={\partial }_{X_0{X}_0}+2\delta {\partial }_{X_0{X}_1}+{\delta }^2{\partial }_{X_1{X}_1},{\partial }_{\tau }={\delta }^2{\partial }_{{\tau }_2}.$

Therefore, the solution at the leading order is

(26)$h_1=A\left(X_1,{\tau }_2\right)e^{i{K}_c{X}_0}+A^{\ast }\left(X_1,{\tau }_2\right)e^{-i{K}_c{X}_0},f_1={\alpha }_1{h}_1,{\gamma }_1={\alpha }_2{h}_1$

(${\alpha }_1$ and ${\alpha }_2$ are the same as in the previous section). At the second order the solution is

(27)$h_2=h_0\left(X_1,{\tau }_2\right)+B_1\left(X_1,{\tau }_2\right)e^{2i{K}_c{X}_0}+B_1^{\ast }\left(X_1,{\tau }_2\right)e^{-2i{K}_c{X}_0},$

(28)$f_2=f_0\left(X_1,{\tau }_2\right)+B_2\left(X_1,{\tau }_2\right)e^{2i{K}_c{X}_0}+B_2^{\ast }\left(X_1,{\tau }_2\right)e^{-2i{K}_c{X}_0},$

(29)${\gamma }_2={\gamma }_0\left(X_1,{\tau }_2\right)+B_3\left(X_1,{\tau }_2\right)e^{2i{K}_c{X}_0}+B_3^{\ast }\left(X_1,{\tau }_2\right)e^{-2i{K}_c{X}_0}.$

The equation for the heat transfer, (2), at the second order gives us the relation between $h_0$ and $f_0$:

(30)$f_0=h_0-\frac{K_c{}^2}{\beta }{\left|A\right|}^2.$

At the third order the solvability condition yields the evolution equation for the amplitude $A$ similar to that obtained by Golovin et al. [13] in the problem of the interaction of the short-wave and long-wave instability modes:

(31)${\partial }_{{\tau }_2}A={\kappa }_{0}A+{\kappa }_1{\left|A\right|}^{2}A+{\mu }_1{\partial }_{X_1{X}_1}A+{\mu }_2{h}_{0}A+{\mu }_3{\gamma }_{0}A.$

The coefficients ${\mu }_i$ $\left(i=1,2,3\right)$ are cumbersome and are not presented here. The equation for the surface deformation $h_0$ which can be obtained from the kinematic Equation (1), as well as the equation for perturbations of surfactant concentration ${\gamma }_0$, are obtained at the fourth order as (here we eliminated $f_0$ using (31)):

(32)${\partial }_{{\tau }_2}{h}_0=a_1{\partial }_{X_1{X}_1}{\left|A\right|}^2+a_2{\partial }_{X_1{X}_1}{h}_0+a_3{\partial }_{X_1{X}_1}{\gamma }_0,$

(33)${\partial }_{{\tau }_2}{\gamma }_0=b_1{\partial }_{X_1{X}_1}{\left|A\right|}^2+b_2{\partial }_{X_1{X}_1}{h}_0+b_3{\partial }_{X_1{X}_1}{\gamma }_0,$

where

(34)$$\begin{gathered} a_1=\frac{1}{3}\left(G+7K^{2}S\right)-\frac{M_m{K}_c{}^2}{2\beta },a_2=\frac{G}{3},a_3=\frac{N}{2}, \\ {b}_1=\frac{1}{36L}\left[G\left(12L-G\right)+4K_c{}^2\left(G+21L\right)S+5K_c{}^4{S}^2\right]-\frac{M_m{K}_c{}^2}{\beta }, \\ {b}_2=\frac{G}{2},b_3=L+N. \end{gathered}$$

The obtained system of Equations (31)–(33) describes the evolution of modulated rolls near the instability threshold. Recall that $A$ is the amplitude of the rolls (see (26)), while $h_0$ and ${\gamma }_0$ are long-wave deformation of the free surface and concentration of the surfactant, respectively (see (27),(29)). Therefore, additional terms ${\mu }_2{h}_{0}A$ and ${\mu }_3{\gamma }_{0}A$ reflect the influence of large-scale surface deformation and inhomogeneity of surfactant concentration on the roll-like Marangoni convection. In Equations (32) and (33) the term ${\partial }_{X_1{X}_1}{\left|A\right|}^2$ is due to a large-scale flow and a surfactant flux generated by the spatial inhomogeneity of the convection roll amplitude.

System (32)–(33) is solved under global conditions on $h_0$ and ${\gamma }_0$:

(35)$⟨h_0⟩=0,⟨{\gamma }_0⟩=0.$

Here $⟨\dots ⟩$ is averaging over $X_1$; the appearance of convection cannot change the total volume of the liquid and the amount of the surfactant.

4.2. Equation for the Growth Rate

Now we consider the linear stability analysis of the spatially periodic roll patterns, in the framework of the system of Equations (31)–(34).

Let us consider the family of stationary roll patterns

(36)$A_0=\sqrt{\frac{q^2{\mu }_1-{\kappa }_0}{{\kappa }_1}}{e}^{iq{X}_1},h_0={\gamma }_0=0\left({\kappa }_1<0\right)$

The parameter $q$ is the rescaled deviation of the wave number from the critical value $K_c$; $q^2<q_0^2=\raisebox{1ex}{${\kappa }_0$}\!\left/ \!\raisebox{-1ex}{${\mu }_1$}\right.$. Solutions (36), (26) describe the stationary spatially periodic supercritical convection with constant amplitude, without interfacial deformation and disturbance of the surfactant concentration independent of $X_0$.

To describe the roll stability, we add small perturbations to the roll amplitude and phase in the form $A=A_0\left(1+r\right)e^{i\left(q{X}_1+\phi \right)},$ where $r$ and $\phi$ are real functions of ${\tau }_2$ and $X_1$.

Linearizing Equations (32)–(34), we obtain the following evolution equations for $r,$ $\phi ,$ $h_0,$ and ${\gamma }_0:$

(37)${\partial }_{{\tau }_2}r=2\left(q^2{\mu }_1-{\kappa }_0\right)r+{\mu }_1{\partial }_{X_1{X}_1}r-2{\mu }_{1}q{\partial }_{X_1}\phi +{\mu }_2{h}_0+{\mu }_3{\gamma }_0,$

(38)${\partial }_{{\tau }_2}\phi =2{\mu }_{1}q{\partial }_{X_1}r+{\mu }_2{\partial }_{X_1{X}_1}\phi ,$

(39)${\partial }_{{\tau }_2}{h}_0=2a_1\left(\frac{{\kappa }_0-q^2{\mu }_1}{{\kappa }_1}\right){\partial }_{X_1{X}_1}r+a_2{\partial }_{X_1{X}_1}{h}_0+a_3{\partial }_{X_1{X}_1}{\gamma }_0,$

(40)${\partial }_{{\tau }_2}{\gamma }_0=2b_1\left(\frac{{\kappa }_0-q^2{\mu }_1}{{\kappa }_1}\right){\partial }_{X_1{X}_1}r+b_2{\partial }_{X_1{X}_1}{h}_0+b_3{\partial }_{X_1{X}_1}{\gamma }_0.$

Though the disturbances are real, one can perform the stability analysis using normal modes in the form

$\left(r,\phi ,h_0,{\gamma }_0\right)=\left(\tilde{r},\tilde{\phi },\widetilde{h_0},\widetilde{{\gamma }_0}\right)e^{\lambda {\tau }_2+ik{X}_1},$

where the eigenvector $\left(\tilde{r},\tilde{\phi },\widetilde{h_0},\widetilde{{\gamma }_0}\right)$ and the eigenvalue $\lambda (k)$ can be complex.

Substituting these modes into (37)–(40) we obtain the following equation for the growth rate of the side-band instability:

(41)${\lambda }^4+c_{0,1}{\lambda }^3+c_{0,2}{\lambda }^2+c_{0,3}\lambda +c_{0,4}=0,$

where

$c_{0,1}=c_{1,0}+c_{1,2}{k}^2,c_{0,2}=c_{2,2}{k}^2+c_{2,4}{k}^4,c_{0,3}=c_{3,4}{k}^4+c_{3,6}{k}^6,c_{0,4}=c_{4,6}{k}^6+c_{4,8}{k}^8\mathbf{}$

and

$$\begin{gathered} c_{1,0}=2\left({\kappa }_0-q^2{\mu }_1\right),c_{1,2}=a_2+b_3+2{\mu }_1, \\ {c}_{2,2}=2\left(a_2+b_3\right)\left({\kappa }_0-q^2{\mu }_1\right)+2{\mu }_1\left({\kappa }_0-3q^2{\mu }_1\right)-\frac{2}{{\kappa }_1}\left(a_1{\mu }_2+b_1{\mu }_3\right)\left({\kappa }_0-q^2{\mu }_1\right), \\ {c}_{2,4}={\mu }_1^2+2\left(a_2+b_3\right){\mu }_1+a_2{b}_3-a_3{b}_2, \\ \begin{array}{l}{c}_{3,4}=\frac{2}{{\kappa }_1}\{\left({\kappa }_0-q^2{\mu }_1\right)[{\kappa }_1\left(a_2{b}_3-a_3{b}_2\right)-{\mu }_1\left(a_1{\mu }_2+b_1{\mu }_3\right)+{\mu }_2\left(a_3{b}_1-a_1{b}_3\right)+\\ +{\mu }_3\left(a_1{b}_2-a_2{b}_1\right)]+\left({\kappa }_0-3q^2{\mu }_1\right)\left(a_2+b_3\right){\kappa }_1{\mu }_1\},\end{array} \\ {c}_{3,6}=\left(a_2+b_3\right){\mu }_1^2+2{\mu }_1\left(a_2{b}_3-a_3{b}_2\right), \\ {c}_{4,6}=2{\mu }_1\left(a_2{b}_3-a_3{b}_2\right)\left({\kappa }_0-3q^2{\mu }_1\right)+\frac{2{\mu }_1}{{\kappa }_1}\left({\kappa }_0-q^2{\mu }_1\right)\left[\left(a_3{b}_1-a_1{b}_3\right){\mu }_2+\left(a_1{b}_2-a_2{b}_1\right){\mu }_3\right], \\ {c}_{4,8}=\left(a_2{b}_3-a_3{b}_2\right){\mu }_1^2. \end{gathered}$$

It is known, see [11], that the limit $k^2\ll 1$ is crucial for the long-wave modulations of the patterns. Assuming $\lambda ={\lambda }_0+\mathsf{\Lambda }{k}^2+o(k^2)$, we obtain at the leading order,

${\lambda }_0^4+c_{1,0}{\lambda }_0^3=0.$

We find that one root,

${\lambda }_0=-c_{1,0}=2\left(q^2{\mu }_1-{\kappa }_0\right),$

corresponding to amplitude modulation of the roll pattern, is negative, while other three roots tend to zero when $k\to 0$. For those roots, we obtain at the leading order

(42)$c_{1,0}{\mathsf{\Lambda }}_{}^3+c_{2,2}{\mathsf{\Lambda }}_{}^2+c_{3,4}\mathsf{\Lambda }+c_{4,6}=0.$

This equation describes the interaction of three “soft” (Goldstone) modes corresponding to definite symmetries of the problem, phase disturbances (translational invariance of the envelope function), surface deformation (conservation of the volume) and surfactant concentration disturbances (conservation of the amount of surfactant).

4.3. Stability in the Absence of Surfactant

First, let us consider the roll modulation in the absence of surfactant $\left(N=0\right)$. In that case, coefficients ${\mu }_3$ and $a_3$ in (32), (33) are equal to zero, hence those equations do not contain ${\gamma }_0$. Equations (37)–(39) form a closed system, while Equation (40), which has no physical meaning, has the obvious solution

${\gamma }_0=\widetilde{{\gamma }_0}{e}^{\lambda {\tau }_2+ik{X}_1},\lambda =-(L+N)k^2,$

hence

$\mathsf{\Lambda }=-(L+N).$

Therefore, the dispersion relation (42) can be presented as

$\left({\tilde{c}}_{1,0}{\mathsf{\Lambda }}^2+{\tilde{c}}_{2,2}\mathsf{\Lambda }+{\tilde{c}}_{3,4}\right)\left(\mathsf{\Lambda }+L+N\right)=0.$

The analysis of the dispersion relation is reduced to the consideration of a quadratic equation,

(43)${\tilde{c}}_{1,0}{\mathsf{\Lambda }}^2+{\tilde{c}}_{2,2}\mathsf{\Lambda }+{\tilde{c}}_{3,4}=0,$

where ${\tilde{c}}_{1,0}=2({\kappa }_0-q^2{\mu }_1),$${\tilde{c}}_{2,2}=2\left(a_2-\frac{a_1{\mu }_2}{{\kappa }_1}\right)\left({\kappa }_0-q^2{\mu }_1\right)+2{\mu }_1\left({\kappa }_0-3q^2{\mu }_1\right),$ and ${\tilde{c}}_{3,4}=\frac{2{\mu }_1}{{\kappa }_{{}_1}}\left[a_2{\kappa }_1\left({\kappa }_0-3q^2{\mu }_1\right)-a_1{\mu }_2\left({\kappa }_0-q^2{\mu }_1\right)\right].$

Equation (43) for the growth rate has two roots. Let us plot the growth rates ${\mathsf{\Lambda }}_{1,2}$ as functions of $q^2$ for several values of $\beta$ at fixed parameter $G=1$ (parameters $S=1,M_2=1$).

Figure 3a shows the growth rate at $\beta =5.$ For $q^2<q_m^2\approx 0.193,$ one root is positive and another one negative; for $q^2>q_m^2$, both roots are positive. We find that in this case the rolls with arbitrary values of $q$ are monotonically unstable. Figure 3b presents case of $\beta =10,$ where one growth rate is negative for all $q$, but the second one is positive in whole interval $q^2<q_0^2\approx$0.112, where the roll solution (37) exists. Again, all the rolls are monotonically unstable. At $\beta =10$(Figure 3c), in the range $q^2<q_m^2\approx 0.0245$ the rolls are stable (both ${\mathsf{\Lambda }}_{1,2}$ are negative); for $q^2$ between $q_m^2$ and $q_0^2=0.103,$ one of the growth rates is positive, therefore we have unstable rolls. Starting from $\beta =16.7,$ the instability is oscillatory on the boundary of the stability interval of rolls. At Figure 3d the dashed line describes the imaginary part of the growth rate. The rolls are stable when $q^2<q_{osc}^2$ and unstable (oscillatory or monotonically) when $q_{osc}^2<q^2<q_0^2.$

The stability map of the stationary rolls with respect to long-wave modulation is shown in Figure 4. The black solid line here indicates the boundary between supercritical and subcritical regions for the rolls. We consider only the region of supercritical rolls (case ${\kappa }_1<0$). The previous investigation of Shklyaev et al. [5] shows the existence of long-wave oscillatory instability of the conductive state above the red solid line on Figure 4. Therefore, the analysis of the modulational instability of rolls has to be done only between the black solid line and the red solid line, where there are three regions:

• (i). region 1, where all the rolls are unstable with respect to modulations;

• (ii). region 2, where the rolls are stable within the interval $q^2<q_m^2$ and monotonically unstable for $q_m^2<q^2<q_0^2$;

• (iii). region 3, where the rolls are stable within the interval $q^2<q_{osc}^2,$ oscillatory unstable for $q^2$ slightly above $q_{osc}^2$, and oscillatory or monotonically unstable for any $q_{osc}^2<q^2<q_0^2.$

All regions are shown in Figure 4a. The values of $q_m^2$ and $q_{osc}^2$ can be found analytically.

The boundary of the monotonic Eckhaus instability (${\mathsf{\Lambda }}_0=0$) is given by the relation ${\tilde{c}}_{3,4}\left(q_m^2\right)=0$, hence

$q_m^2=\frac{{\kappa }_0}{{\mu }_1}\left(\frac{a_2{\kappa }_1-a_1{\mu }_2}{3a_2{\kappa }_1-a_1{\mu }_2}\right).$

The quantity $q_{osc}^2$ is determined by the relation ${\tilde{c}}_{2,2}\left(q_{osc}^2\right)=0,$ hence

$q_{osc}^2=\frac{{\kappa }_0\left(a_2{\kappa }_1+{\kappa }_1{\mu }_1-a_1{\mu }_2\right)}{{\mu }_1\left(a_2{\kappa }_1+3{\kappa }_1{\mu }_1-a_1{\mu }_2\right)}.$

Figure 4b shows an additional region of stable rolls in the neighborhood of the origin (small values of $\beta$ and $G$).

4.4. Influence of Insoluble Surfactant

As shown in our previous research [10], even rather small quantity of insoluble surfactant on the free surface of liquid can significantly change the instability thresholds and creates new oscillatory regimes. That instability is determined by both the surface property of the surfactant (${\sigma }_2$) and its concentration (${\mathsf{\Gamma }}_0$) through the elasticity number $N$ that includes both factors.

In the present paper, we consider the influence of a very small amount of insoluble surfactant on the modulation stability of stationary rolls. For this goal we use dispersion relation (42). From this relation the monotonic Eckhaus instability is determined by the relation $c_{4,6}=0$ or

(44)$q_m^2=q_0^2\frac{1+ℂ}{3+ℂ},$

where $q_0$ is the width of the existence interval of the periodic solutions ($0\le q^2\le q_0^2$), and

(45)$ℂ=\frac{\left(a_3{b}_1-a_1{b}_3\right){\mu }_2-\left(a_2{b}_1-a_1{b}_2\right){\mu }_3}{\left(a_2{b}_3-a_3{b}_2\right){\kappa }_1}.$

If the interaction of the roll convection with the deformations and disturbances of the surfactant concentration is switched off, i.e., ${\mu }_2={\mu }_3=0$ or $a_1=a_2=0$, then $ℂ=0$, and the classical result

$q_m^2=\frac{1}{3}{q}_0^2$

is recovered. The expression for the oscillatory Eckhaus instability boundary $q_{osc}^2$ is cumbersome, and it is not presented here.

For plotting stability maps we take two values of the elasticity parameter: $N={10}^{-8}$ and $N={10}^{-7}.$ In both cases the Lewis number was fixed at $L=0.003.$ The stability maps are shown in Figure 5. Here Figure 5a corresponds to the case $N={10}^{-8}$ and Figure 5b to $N={10}^{-7}$.

Comparing the stability maps with surfactant and without surfactant we conclude that the region of stable rolls shrinks with addition of the insoluble surfactant while the subcriticality region expands. The behavior of the growth rates in each of these regions (“1”, “2”, and “3”) are similar to those describing in the case without surfactant. The shrinking of the region of stable supercritical rolls means that with addition of an insoluble surfactant the appearance of a subcritical instability becomes more typical. Therefore, the linear stability analysis become inefficient, and the nonlinear analysis is unavoidable.

5. Conclusions

In the present paper we have studied a stationary Marangoni convection without and with insoluble surfactant spread over a deformable free surface. The analysis of nonlinear dynamics of perturbations near the instability threshold shows the existence of two kinds of steady supercritical periodic structures–rolls and squares. We investigated the dependence of their stability regions on the heat transfer at the liquid surface, characterized by Biot number, deformability of the free surface of the liquid characterized by Galileo number, and the action of the surfactant on the surface tension characterized by the elasticity number. The presented theory can be applied to different active agents on different surfaces, independently of the physico-chemical mechanism of the surface-tension decrease, on the condition of low solubility of the surfactant (the influence of the surfactant on the properties of the bulk liquid is assumed to be negligible).

It is shown how the regions of pattern selection change under the influence of the surfactant.

In addition, we carried out the linear analysis of modulational instability of stationary rolls near the bifurcation point. Using the weakly nonlinear analysis, the set of amplitude equations that describe large-scale one-dimensional (longitudinal) distortions of periodic roll patterns is derived. The stability maps of the stationary rolls with and without surfactant have been presented. The results indicate the existence of three regions of the supercritical rolls: (i) region where the rolls are unstable with respect to modulations; (ii) region where the rolls are stable within the interval $q^2<q_m^2$ and monotonically unstable for $q_m^2<q^2<q_0^2$; (iii) region where the rolls are stable within the interval $q^2<q_{osc}^2,$ oscillatory unstable for $q^2$ slightly above $q_{osc}^2$, and oscillatory or monotonically unstable for any $q_{osc}^2<q^2<q_0^2.$

To check the results of described theory we propose performing an experiment on Marangoni convection in a rectangular container with one size much longer than another one or in a narrow annular container of large radius.

Author Contributions

A.B.M.: software, formal analysis, investigation, writing–original draft, writing–review & editing and visualization; A.A.N.: conceptualization, methodology, validation, writing–original draft, review & editing. Both authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Acknowledgments

Authors acknowledge the support by the Israel Science Foundation (Grant No. 843/18).

Conflicts of Interest

The authors declare no conflict of interest.

Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Figures

Figure 1. Critical elasticity number Nm vs. parameter β for different values of parameter G at the critical wave number Kc(m). Curve #1 is for G=2, curve #2 is for G=1 and curve #3 is for G=0.5. Other parameters are: S=1 and L=0.003.

Figure 2. Pattern selection domains for stability on the square lattice. The solid lines are for κ1=0, the dashed lines are for κ1+κ2=0 and dot-dashed lines are for κ2−κ1=0. Here the red color set (#1) is for the case without surfactant (N=0 ), the blue colored set (#2) is for N=10−6 and the black colored set (#3) is for N=10−5. The stability domains are marked, “R” is for rolls, “S” is for squares. Other parameters: S=1 and L=0.003.

Figure 3. Growth rates Λ1,2 of Equation (44) versus q2 for different values of parameter β. Panel (a) β=5, (b) β=10, (c) β=15, and (d) β=20. Solid lines are for real part of Λ1,2, dashed lines for the imaginary one. Other parameters: G=1,S=1,M2=1.

Figure 4. Stability maps of the stationary rolls without insoluble surfactant. Panel (a) Solid black line is a boundary between supercritical and subcritical roll regions (“subcr. R”). Region “1” is a region of unstable rolls; region “2”–rolls are stable within q2<qm2 and monotonically unstable for qm2<q2<q02; region “3”-rolls are stable within the interval q2<qosc2, oscillatory unstable for q2 slightly above qosc2, and oscillatory or monotonically unstable for any qosc2<q2<q02. Panel (b) the zoomed-in region of panel (a) at small values of β and G (“stab. R” is a region of stable rolls).

Figure 5. Stability maps of the modulated stationary rolls with insoluble surfactant. Solid black line is a boundary between supercritical and subcritical roll regions (“subcr. R”). Region “1” is a region of unstable rolls; region “2”–rolls are stable within q2<qm2 and monotonically unstable for qm2<q2<q02; region “3”-rolls are stable within the interval q2<qosc2, oscillatory unstable for q2 slightly above qosc2, and oscillatory or monotonically unstable for any qosc2<q2<q02. Panel (a) N=10−8, panel (b) N=10−7. Other parameters: L=0.003,S=1.