ESI for Suppression of Self-Stratification in Colloidal Mixtures with High Peclet Numbers by Schulz et al in Soft Matter

Model and comparison of model with experimental findings for stratification

by Malin Schulz, Richard Brinkhuis, Carol Crean, Richard P. Sear and Joseph L. Keddie

This Jupyter notebook has the Python code for a model that predicts when stratification will and will not occur, in the parameter space of the volume fraction of small colloid, $\phi_S$ and Peclet number for the small colloid, Pe$_S$. It also briefly summarises the equations, see main text for more details.

import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
from scipy import stats

Fit power law to measurements for viscosity at pH=12 and shear rate 0.1/s

We fit a power law

$$ \eta=a[\mbox{ASE %}]^b $$

to a set of experimental data at pH=12 and at the lowest shear rate we could reach, $0.1/$s. Here ASE = alkali-swellable emulsion - our thickener, and % ASE $=100\phi_{CP}$. Note that as viscosity measurements are done without colloids present, there $\phi_T=\phi_{CP}$.

We do the fit to logs of the data, i.e.,

$$ \log_{10}\eta=\log_{10}a+b\log_{10}[\mbox{ASE %}] $$

data is of form 4 measurements of viscosity in Pa s at pH=12 at % ASE $\equiv 100\phi_T$.

ASE_percent=np.zeros(4)
ASE_percent[:]=[1.5,2.5,3,4]
log10ASE_percent=np.log10(ASE_percent)
log10Viscosity=np.zeros(4)
log10Viscosity[:]=[-1.20615,-0.17145,0.41039,0.679547]
Viscosity=np.zeros(4)
Viscosity=10.0**log10Viscosity
print(' % ASE         ',ASE_percent)
print(' Viscosity ',np.round(Viscosity,4),' Pa s')
print('')
gradient, intercept, r_value, p_value, std_err = stats.linregress(log10ASE_percent,log10Viscosity)
print('intercept, gradient ',intercept,gradient)
print('exponent b = ',round(gradient,4))
a=10.0**intercept
print('     and a = ',round(a,4))

 % ASE          [1.5 2.5 3.  4. ]
 Viscosity  [0.0622 0.6738 2.5727 4.7813]  Pa s

intercept, gradient  -1.97591338472721 4.606782537301664
exponent b =  4.6068
     and a =  0.0106

plot data and fit to check that it looks OK

n_pts=100
ASE_percent_plot = np.linspace(1, 10, n_pts)
log10ASE_percent_plot=np.log10(ASE_percent_plot)
Viscosity_plot=10.0**(intercept+gradient*log10ASE_percent_plot)
plt.semilogy(ASE_percent_plot/100.0,Viscosity_plot,lw=4)
plt.semilogy(ASE_percent/100.0,Viscosity,lw=0,marker='o',markersize=20)

plt.ylabel(r'viscosity  Pa s', fontsize=24)
plt.xlabel(r'$\phi_T$', fontsize=24)
plt.xticks(fontsize=20)
plt.yticks(fontsize=20)
plt.xlim([0,0.06])
#plt.ylim([0,c_t[90]])
plt.show()
#
testASE=2.5#/0.7
testViscosity=10.0**(intercept+gradient*np.log10(testASE))
print('predicted viscosity at',round(testASE,4),' % ASE = ',testViscosity)

predicted viscosity at 2.5  % ASE =  0.7199622045495464

Now plot figures for functions of time during evaporation, for a particular experiment

Evaporation_speed= 140.0e-9 #[m/s]
print('Evpaoration speed E ',np.round(Evaporation_speed*1e9),' nm/s')
initial_filmthickness=300.0e-6 #[m]
print('initial film thickness ',np.round(initial_filmthickness*1e6),' nm/s')

t_evap=initial_filmthickness/Evaporation_speed
print('(hypothetical) time for film thickness to go to zero',t_evap,' s')
initial_ASE=1.0
initial_effASE=initial_ASE/0.7
print('initial ASE % ',initial_ASE,' and effective % ASE ',initial_effASE,' = ASE %/0.7 as 30% volume occupied by colloid')
print('initial phi_T ',initial_ASE/100.0,' and effective phi_T ',initial_effASE/100.0,' = phi_T/0.7 as 30% volume occupied by colloid')

Evpaoration speed E  140.0  nm/s
initial film thickness  300.0  nm/s
(hypothetical) time for film thickness to go to zero 2142.8571428571427  s
initial ASE %  1.0  and effective % ASE  1.4285714285714286  = ASE %/0.7 as 30% volume occupied by colloid
initial phi_T  0.01  and effective phi_T  0.014285714285714285  = phi_T/0.7 as 30% volume occupied by colloid

$$ t_{evap}=\frac{H}{v_{evap}} $$

is a nominal evaporation time, in the sense that as we start with 30% (involatile) solids then all the solvent has evaporated in 70% of $t_{evap}$. We are assuming a constant rate of descent of the water/air interface. So the fraction of the total initial volume of the film occupied by the solvent or continuous phase, after $t$ seconds of evaporation at assumed constant speed $v_{evap}$ is

$$ f_{CP}(t)=0.7-t/t_{evap} $$

which starts at $0.7$ because the continuous phase only occupies 70% of the initial volume, with the rest occupied by solids. After a time $t$ the effective concentration of ASE is

$$ \mbox{effective ASE concentration}=\phi_T/f_{CP}(t) $$

t = np.linspace(0, 0.6999*t_evap, n_pts)#[s]



fraction_CP_t=0.7-t/t_evap

effASE_of_t=initial_ASE/fraction_CP_t





plt.plot( t/60.0,effASE_of_t/100.0,lw=4)#,'k')
plt.ylabel(r' $\phi_{CP}$', fontsize=24)
plt.xlabel(r'$t$    [ min ]', fontsize=24)
plt.xticks(fontsize=20)
plt.yticks([0,0.05,0.1,0.15,0.2],fontsize=20)
plt.xlim([0,25])
plt.ylim([0,0.2])
plt.tight_layout()
plt.savefig('figureS3a.png',dpi=900)

plt.show()

########################

def Viscosity_fit(ASE):
    Viscosity=10.0**(intercept+gradient*np.log10(ASE))
#    print(intercept,gradient,Viscosity)
    return Viscosity

log10effASE=np.log10(effASE_of_t)
Viscosity_of_t=Viscosity_fit(effASE_of_t)




plt.semilogy(t/60,Viscosity_of_t,lw=4)

plt.ylabel(r'viscosity  [ Pa s ]', fontsize=24)
plt.xlabel(r'$t$    [ min ]', fontsize=24)
plt.xticks(fontsize=20)
plt.yticks(fontsize=22)
plt.ylim([0.01,1e4])
plt.xlim([0,25])
plt.tight_layout()
plt.savefig('figure3a.png',dpi=900)
plt.show()

RSmall=14.0e-9



def StokesEinstein(T,effectiveASE):
    if(effectiveASE==0):
# in Pa s at T=60 C for water
        viscosity=4.658e-4
    else:
# with ASE so use fit
        viscosity=Viscosity_fit(effectiveASE)
#
    D=1.381e-23*T/(6.0*np.pi*viscosity*RSmall)
    return D
# T in K
T=303.15
D_of_t=np.zeros(n_pts)
for i in range(0,n_pts):
    D_of_t[i]=StokesEinstein(T,effASE_of_t[i])
#
print('initial value of D         ',D_of_t[0],' m^2/s')

plt.semilogy( t/60.0,D_of_t*1e12,lw=4)#,'k')
plt.ylabel(r'$D$   $\mu$m$^2$/s', fontsize=24)
plt.xlabel(r'$t$    [ min ]', fontsize=24)
plt.xticks(fontsize=20)
plt.yticks(fontsize=20)
plt.xlim([0,25])
plt.ylim([1e-6,1])
plt.tight_layout()
plt.savefig('figureS3b.png',dpi=900)
plt.show()
#
phijam=0.64
phi0=0.1
vjam=Evaporation_speed/(1.0-phi0/phijam)
print('speed of descent of jammed layer v_jam ',np.round(vjam*1e9),' nm/s')
lambdaG_of_t=D_of_t/vjam
#
plt.semilogy( t/60.0,lambdaG_of_t*1e9,lw=4)#,'k')
plt.ylabel(r'$\lambda_G$   [ nm ] ', fontsize=24)

plt.xlabel(r'$t$    [ min ]', fontsize=24)

plt.xticks(fontsize=20)
plt.yticks(fontsize=20)
plt.xlim([0,25])
plt.ylim([10,3000])
plt.tight_layout()
plt.savefig('figure3b.png',dpi=900)
plt.show()

initial value of D          2.902256769542028e-13  m^2/s

speed of descent of jammed layer v_jam  166.0  nm/s