fig11.ipynb

Kernel: Python 3 (system-wide)

import numpy as np
import matplotlib.pyplot as plt
import cmocean 
from matplotlib import cm
from matplotlib.colors import ListedColormap, LinearSegmentedColormap
import scipy as sp
import scipy.special 
import seaborn as sns
import scipy.ndimage
from matplotlib.gridspec import GridSpec
%matplotlib inline
from matplotlib import rc
rc('font',**{'family':'sans-serif','sans-serif':['Helvetica'], 'size':20})
rc('text', usetex=True)
rc('text.latex', preamble=r'\usepackage{mathrsfs}')
from matplotlib.gridspec import GridSpec

RC = {'figure.figsize':(10,5),
      'axes.facecolor':'white',
      'axes.edgecolor':' .15',
      'axes.linewidth':' 1.25',
      'axes.grid' : True,
      'grid.color': '.8',
      'font.size' : 15,
     "xtick.major.size": 4, "ytick.major.size": 4}

plt.rcParams.update(RC)
sns.set_palette(sns.color_palette('dark'), n_colors=None, desat=None, color_codes=False)

#Define non-dimensional parameters for the simulation 
Re = 2.478883e+03 
Fr = 1.105914e+00
Ri = 1/((Fr/(2*np.pi)))**2
Pr = 7
ch = 50 #input column height

#Define dimensionless functions f and g from equations (2.14) and (2.15) in the manuscript
#These can be easily modified by the user
a = 1
b = 0.8
c = 0.9
d = 0.9 
def f(x):
    return 19/8 + 11/8*np.tanh(a*np.log(x)-b)

def g(x):
    return 2 + np.tanh(c*np.log(x)-d)

#Model name
model = 'shearr_mixed_all_ch50'
data_loc = 'processed_data/test/Pr7_decaying/'
output_loc = 'PCNN_chi/'

t = ['T=0.5','T=1', 'T=2', 'T=4', 'T=6', 'T=7.7']
#Inputs are normalised during training for faster convergence, this must be taken into account when taking derivatives 
#with respect to the raw inputs
sigma_X = 0.32755166052145035 
sigma_Y = 0.26732916889081537
xdim, zdim = 500, 500
sns.set_palette("tab10")
gs=GridSpec(1,2, width_ratios=[1,1], wspace=0.1, hspace=0.12)
fig=plt.figure(figsize=(12,3.6))
ax1=fig.add_subplot(gs[0,0]) 
ax2=fig.add_subplot(gs[0,1]) 
rc('font',**{'family':'sans-serif','sans-serif':['Helvetica'], 'size':20})
x = np.linspace(-25,25,50)
for i in range(6):
    pcnn_mean_chi = np.load(output_loc+model+'_means_'+t[i]+'.npy').reshape((xdim,zdim))
    gradX = np.load(output_loc+model+'_grads_'+t[i]+'.npy')[:,:,1].sum(axis=1)
    gradY = np.load(output_loc+model+'_grads_'+t[i]+'.npy')[:,:,0].sum(axis=1)
    D_X = np.log(10)*np.reshape(gradX,(xdim,zdim))*10**(pcnn_mean_chi[:,:])/sigma_X
    D_Y = np.log(10)*np.reshape(gradY,(xdim,zdim))*10**(pcnn_mean_chi[:,:])/sigma_Y
    histX = np.histogram(D_X.flatten(), bins=50, range=(-5,11), density=True)
    histY = np.histogram(D_Y.flatten(), bins=50, range=(-5,5), density=True)
    ax1.fill_between(histX[1][:-1],histX[0], alpha=0.4+i/10, zorder=1-i, label=t[i])
    ax2.fill_between(histY[1][:-1],histY[0], alpha=0.4+i/10, zorder=1-i, label=t[i])
plt.legend(fontsize=10)
rc('font',**{'family':'sans-serif','sans-serif':['Helvetica'], 'size':20})
ax2.set_xlim(-1,2)
ax1.set_xlim(-5,10)
ax1.set_ylim(0,0.75)
ax2.set_ylim(0,3.1)
ax1.axvline(3, color='k', linestyle='--')
ax1.axvline(1, color='k', linestyle='--')
ax1.set_ylabel('$\\mathrm{p.d.f.}$')
ax2.set_xlabel('$D_Y$')
ax1.set_xlabel('$D_X$')
ax2.set_yticklabels([])
ax1.annotate('$a)$', (-6,0.77), annotation_clip=False)
ax1.annotate('$b)$', (10.5,0.77), annotation_clip=False)
fig.show()