Appendix to : A hydrodynamic instability is used to create aesthetically appealing patterns in painting

Appendix to: A hydrodynamic instability is used to create aesthetically appealing patterns in painting Sandra Zetina, Francisco A. Godı́nez, Roberto Zenit3,∗ 1 Instituto de Investigaciones Estéticas, Universidad Nacional Autónoma de México, México D.F. 04510, México 2 Instituto de Ingenierı́a, Universidad Nacional Autónoma de México, México D.F. 04510, México 3 Instituto de Investigaciones en Materiales, Universidad Nacional Autónoma de México, México D.F. 04510, México ∗ E-mail: zenit@unam.mx

where D= d/dz, w is the vertical velocity component, µ and ρ are the viscosity and density of the fluid, respectively. n and k = k 2 x + k 2 y are frequency and wave number, respectively, of the normal disturbance solutions on x, y and t, given by exp(ik x x + ik y y + nt) If the fluid is confined, the following conditions must be satisfied at the boundaries Additionally, w, Dw and µ(D 2 + k 2 )w must be continuos across the interface between the two fluids. If surface tension is to be considered, the following relation must also be satisfied: where σ is the surface tension. ∆ s (x) indicates the difference between the values of the property x on each side of the interface. The subscript s indicates the value of the property at the interface. This relation arises from the continuity of stresses across the interface.
A thin viscous layer composed of two miscible fluids on top of each other For the flow that occurs in Siqueiros' accidental painting technique, we consider a thin dual layer of two viscous miscible fluids of constant, but different, densities and viscosities, as depicted in Fig.1. The fluids, separated at z = 0, have densities, ρ 1 and ρ 2 , and viscosities, µ 1 and µ 2 . The subscripts 1 and 2 refer to the fluids at the bottom and top, respectively. In this case, σ = 0; therefore, this part is dropped form the rest of the analysis. Assuming that the layer is confined on the bottom and top, the following boundary conditions need to be satisfied: Note that in the accidental painting technique flow the top layer is, in fact, not confined. To model the free surface, the properties of the interface of the top layer with the ambient air would also have to be considered (a 3-fluid layer). This analysis is also considered below. As will be shown later, despite the much longer mathematical description, considering the slip-condition on the top layer does not change the nature of the instability significantly. Therefore, we will retain the confined case as a simpler solution for which more solutions can be readily obtained.  In each of the regions, the general instability equation (Eqn. 1) reduces to: Since both density and viscosity are constant, on each side of the interface, the expression above can be further simplified to: The general solution of Eqn.(8) is: where w 1 and w 2 are the solutions for the bottom (−H 1 < z < 0) and top (0 < z < H 2 ) layers, respectively. We also use q 2 1 = k 2 + n/ν 1 and q 2 2 = k 2 + n/ν 2 , where ν = µ/ρ. Therefore, at z = 0, the interfacial conditions are: where w 0 and (Dw) 0 are the common values of w and Dw at z = 0. From the boundary conditions (Eqns. 5 and 6), considering the w-solutions (Eqns. 9 and 10), we have: After some algebra, we can obtain where are constants. In this manner the number of unknown coefficients is reduced to only four (A 1 , B 1 , A 2 and B 2 ).
From the interfacial conditions (Eqns. 11,12,13, and 14), we have: Equations (15) through (34) can be written as: Equation (38) can be written as: where where At = α 2 − α 1 is the Atwood number and H = H 1 + H 2 . Now, written in matrix form, the equation to solve is where C is: The entire expression for R appears on Solving Eqn.46 the functional relation between n and k can be obtained, resulting in the so-called dispersion relation. From it, the parametric range of system instability can be readily obtained. We solved the equation using Mathematica [2]. We first obtained the solution for the physical properties for the experiment shown in Fig. 5 of the main text; then we varied the physical properties artificially to observe changes in the predicted dispersion relationship. Fig. 8 (of the main text) shows the calculated dispersion relationships. The thick continuous line shows the result obtained for the black-white unstable layer, Table ??.

A triple layer: exact experimental conditions
To model the system completely and more accurately, a free surface condition has to be considered at the upper edge of the top layer. To do that, we can consider the same scheme described in the previous section but also model the motion of the fluid (air on top of the upper layer). Such an approach is briefly described below. In that case, we have two dispersion relations, that have to be solved simultaneously.
Let us consider again the scheme shown in Fig. 1. For this case, we also consider that for z > H 2 there is another fluid layer with properties ρ 3 and µ 3 ; also at z = H 2 , there is a non zero value of the surface tension, σ= 40 mPa m. If we consider the properties of air (ρ 3 =1 kg/m 3 and µ 3 = 2 × 10 −5 Pa s.), we have to consider: where, now, w 3 is the solution for the upper most layer (H 2 < z < H 3 ). Clearly, q 2 3 = k 2 + n/ν 3 . For this case, we consider that H 3 → ∞, therefore A 3 = C 3 = 0 in Eqn. (49).
Hence, in addition to the boundary conditions for the first two layers (Eqns. 11-14), at z = H 2 we need to consider: Table 1. S1 Table. Matrix form of R, from Eqn. 44.
Therefore, we have a system of 10 unknowns (A 1 , B 1 , C 1 , D 1 , A 2 , B 2 , C 2 , D 2 , B 3 , D 3 ) with 10 equations. In the same way as with the case with two layers, we can write the system in a matrix form: The vector C 3 contains the unknown coefficients and R 3 contains all the conditions for the solution to exists. By making det(R 3 ) = 0 we obtain the desired dispersion relation.
Evidently, the solution of this equation is much more elaborated. First of all, the numerical solution is more time-intensive, but also we found that the solution is also numerically unstable; the initial conditions have to be chosen carefully to avoid spurious results. The most important result of this analysis is that we found that the nature of the instability remains largely unchanged. By making the upper layer not to be confined, the size and speed of the most unstable modes increases slightly, but remain in the same order of magnitude. We have calculated only one case, shown in Fig. 9 (of the main text), which corresponds to the physical properties listed on Table 1 (of the main document). As can be seen in the figure, the agreement between the prediction and the measurement is even closer.