Book a Demo!
CoCalc Logo Icon
StoreFeaturesDocsShareSupportNewsAboutPoliciesSign UpSign In
Download

Notebooks supporting the J. Fluid Mech. submission "Turbulent mixed convection in vertical and horizontal channels"

206 views
unlisted
ubuntu2204
Tweet
Share
Share
Kernel: Python 3 (system-wide)

import numpy as np
import sys, os, h5py
import matplotlib.pyplot as plt
import seaborn as sns
import cmocean

# Aesthetics
sns.set_theme()
sns.set_style('ticks')
sns.set_context('paper')
plt.rc('mathtext', fontset='stix')
plt.rc('font', family='serif')

cspeed = cmocean.tools.crop_by_percent(cmocean.cm.tempo, 30, which='both', N=None)

import pandas as pd

Load sheared VC data

df = pd.read_csv('../data/P-VC_global_response.csv')
df.head()

Load sheared RB data

PRB_df = pd.read_csv('../data/P-RB_global_response.csv')
PRB_df.head()

Define quantities needed (Gr=Ra/PrGr=Ra/Pr, Ri=Gr/Re2Ri=Gr/Re^2, Reδ=Wmaxδ/νRe_\delta = W_\mathrm{max} \delta/\nu, Cf=2τx/ρ0Wmax2C_f = 2\tau_x/\rho_0 W_\mathrm{max}^2)

# P-VC
df['Gr'] = df.Ra/df.Pr
df['Ri'] = df.Gr/df.Re**2

# P-RB
PRB_df['Gr'] = PRB_df.Rayleigh/PRB_df.Prandtl
PRB_df['Cf'] = 2*(PRB_df['RMS Spanwise Friction Reynolds (stafield)']/PRB_df['RMS Spanwise Reynolds (stafield)'])**2
PRB_df['Red'] = PRB_df['RMS Spanwise Reynolds (stafield)']*PRB_df['RMS Spanwise Boundary Layer (stafield)']
PRB_df['Ri'] = PRB_df.Gr/PRB_df['Bulk Reynolds']**2

Write function to solve Prandtl friction law prediction

from scipy.optimize import newton
def Cfunc(Reb):
    Cf = np.zeros(Reb.size)
    for i in range(Reb.size):
        def f(x):
            # k = 0.41
            # B = 5
            k = 0.384
            B = 4.27
            return (2/x)**0.5 - 1/k*np.log(Reb[i]/2*(x/2)**0.5) - B
        Cf[i] = newton(f, 2.5e-3)
    return Cf

cmap = cspeed
Ric = 0.25 # Richardson number cutoff value

fig, axs = plt.subplots(1,2, figsize=(5.2,2.2), layout='constrained', dpi=200)

ax = axs[0]

# Transparency of low Ri points
alph = 0.5

# P-RB
sc = ax.scatter(PRB_df.Red*(PRB_df.Ri > Ric), PRB_df.Cf*(PRB_df.Ri > Ric), c=np.log10(PRB_df.Gr), cmap=cmap, s=5, label='mixed RB')
ax.scatter(PRB_df.Red*(PRB_df.Ri <= Ric), PRB_df.Cf*(PRB_df.Ri <= Ric), c=np.log10(PRB_df.Gr), cmap=cmap, s=5, label='mixed RB', alpha=alph)

# P-VC
ax.scatter(df.Red*(df.Pr==1)*(df.Ri > Ric), df.Cfw*(df.Pr==1)*(df.Ri > Ric), c=np.log10(df.Gr), cmap=cmap, ec='k', s=25, label='mixed VC')
ax.scatter(df.Red*(df.Pr==1)*(df.Ri <=Ric), df.Cfw*(df.Pr==1)*(df.Ri <=Ric), c=np.log10(df.Gr), cmap=cmap, ec='k', s=25, alpha=alph)
ax.scatter(df.Red*(df.Pr== 4), df.Cfw*(df.Pr== 4), c=np.log10(df.Gr), cmap=cmap, ec='k', s=25, marker='^', vmin=6, vmax=8)
ax.scatter(df.Red*(df.Pr==10), df.Cfw*(df.Pr==10), c=np.log10(df.Gr), cmap=cmap, ec='k', s=25, marker='v', vmin=6, vmax=8)

fig.colorbar(sc, ax=axs, label='$\log_{10} Gr$')

# Add theoretical comparison
R = 10**np.linspace(0.5,3.5,101)
ax.loglog(R, 8*R**-1, 'k--')
R = 10**np.linspace(2.5,3.5,101)
ax.loglog(R, Cfunc(R), 'k-')

# Labels and aesthetics
ax.grid(True)
ax.set(
    xlabel=r'$Re_\delta = W_\mathrm{max}\delta/\nu$',
    ylabel='$C_f$',
    xlim = [8e0,1.2e3],
    ylim=[1e-2,1]
)
ax.set_ylabel('$C_f$', labelpad=0)
ax.annotate('$(a)$', (-0.05, 1.05), xycoords='axes fraction', ha='right', va='bottom')

# Add inset axes for example vertical velocity profile
axin = ax.inset_axes([0.65,0.7,0.35,0.3])

# Load example profile for Gr=10^7, Pr=1, Re=10^{3.5}
with h5py.File('../data/profile_record.h5','r') as fp:
    grp = 'Gr7_Sc1_Re3.50'
    xm = fp[grp+'/xm'][:]
    vybar = fp[grp+'/vybar'][:]

# Use antisymmetry to average profile
n = xm.size
vys = 0.5*(vybar[:n//2] - vybar[-1:-1-n//2:-1])
axin.plot(xm[:n//2], vys)

# Calculate peak velocity and location, adding labels there
vmx = np.max(vys)
xmx = xm[np.argmax(vys)]
axin.set_xticks([0,xmx], labels=['0','$\delta$'])
axin.set_yticks([0,vmx], labels=['0','$W_\mathrm{max}$'], )
axin.set_xlim([0,0.2])
axin.set_ylim(bottom=0, top=1.4*vmx)
axin.set_title('$\overline{w}(y)$')
axin.grid(True)


ax = axs[1]

# P-RB
ax.scatter(PRB_df.Red/PRB_df['Red0']*(PRB_df.Ri > Ric), PRB_df.Cf/PRB_df['Cf0']*(PRB_df.Ri > Ric), c=np.log10(PRB_df.Gr), cmap=cmap, s=5)
ax.scatter(PRB_df.Red/PRB_df['Red0']*(PRB_df.Ri <=Ric), PRB_df.Cf/PRB_df['Cf0']*(PRB_df.Ri <=Ric), c=np.log10(PRB_df.Gr), cmap=cmap, s=5, alpha=alph)

# P-VC
ax.scatter(df.Red/df.Red0*(df.Pr==1)*(df.Ri > Ric), df.Cfw/df.Cfw0*(df.Pr==1)*(df.Ri > Ric), c=np.log10(df.Gr), cmap=cmap, ec='k', s=25, label='$Ri > \\frac{1}{4}$')
ax.scatter(df.Red/df.Red0*(df.Pr==1)*(df.Ri <=Ric), df.Cfw/df.Cfw0*(df.Pr==1)*(df.Ri <=Ric), c=np.log10(df.Gr), cmap=cmap, ec='k', s=25, alpha=alph, label='$Ri < \\frac{1}{4}$')
ax.scatter(df.Red/df.Red0*(df.Pr==4), df.Cfw/df.Cfw0*(df.Pr==4), c=np.log10(df.Gr), cmap=cmap, ec='k', s=25, marker='^')
ax.scatter(df.Red/df.Red0*(df.Pr==10), df.Cfw/df.Cfw0*(df.Pr==10), c=np.log10(df.Gr), cmap=cmap, ec='k', s=25, marker='v')

# Add theoretical comparison
R = 10**np.linspace(-1,1,101)
ax.loglog(R, R**-1, 'k--')

# Labels and aesthetics
ax.grid(True)
ax.legend()
ax.set(
    xlabel=r'$Re_\delta/ Re_{\delta,0}$',
    ylabel='$C_f/ C_0$',
    xlim=[5e-1, 10],
    ylim=[1e-1, 2e0]
)
ax.set_ylabel('$C_f/ C_0$', labelpad=0)
ax.annotate('$(Re_\\delta/Re_{\\delta,0})^{-1}$', (4, 1.8e-1), ha='center', va='center', rotation=-42)
ax.annotate('$(b)$', (-0.05, 1.05), xycoords='axes fraction', ha='right', va='bottom')

# fig.savefig('Convective_friction_coefficients.pdf')

plt.show()
Image in a Jupyter notebook