<jupyter_start><jupyter_text>NACA0012 Tutorial![naca](http://turbmodels.larc.nasa.gov/NACA0012_grids/grid_bcs.jpg) Referenceshttp://turbmodels.larc.nasa.gov/naca0012_val.html<jupyter_code>import os
if not os.path.exists('output'):
    os.mkdir('output')
    
DATA_DIR='.'
REF_DATA_DIR='.'
from zutil import analysis
case_dict = {'n0012_257_a0p0'      : {'label'    : '0 deg SST/MG/PRECON',
                                      'case' : 1 },
             'n0012_257_a10p0'     : {'label'    : '10 deg SST/MG/PRECON',
                                      'case' : 2 },
             'n0012_257_a15p0'     : {'label'    : '15 deg SST/MG/PRECON',
                                      'case' : 3 },
             'n0012_257_a10p0_rbf' : {'label'    : '10 deg SST/MG/PRECON/RBF',
                                      'case' : 2 } }

validation_dict = {'alpha 0'  : {'case'  : 1,
                                 'alpha' : 0.0,
                                 'cpmax' : 1.2,
                                 'cpmin' : -0.6,
                                 'cfmax' : 0.01,
                                 'cfmin' : 0.0,
                                 'cpexpt': 'CP_Gregory_expdata_a0p0.csv',
                                 'cpcomp': 'n0012cp_cfl3d_sstv_a0p0.dat',
                                 'cfcomp': 'n0012cf_cfl3d_sstv_a0p0.dat',
                                 'clcomp': 0.0,
                                 'cdcomp': 0.00828},
                   'alpha 10' : {'case'  : 2,
                                 'alpha' : 10.0,
                                 'cpmax' : 2.0,
                                 'cpmin' : -6.0,
                                 'cfmax' : 0.04,
                                 'cfmin' : -0.005,
                                 'cpexpt': 'CP_Gregory_expdata_a10p0.csv',
                                 'cpcomp': 'n0012cp_cfl3d_sstv_a10p0.dat',
                                 'cfcomp': 'n0012cf_cfl3d_sstv_a10p0.dat',
                                 'clcomp': 1.0778,
                                 'cdcomp': 0.01245},
                   'alpha 15' : {'case'  : 3,
                                 'alpha' : 15.0,
                                 'cpmax' : 2.0,
                                 'cpmin' : -12.0,
                                 'cfmax' : 0.06,
                                 'cfmin' : -0.005,
                                 'cpexpt': 'CP_Gregory_expdata_a15p0.csv',
                                 'cpcomp': 'n0012cp_cfl3d_sstv_a15p0.dat',
                                 'cfcomp': 'n0012cf_cfl3d_sstv_a15p0.dat',
                                 'clcomp': 1.5060,
                                 'cdcomp': 0.02224}}<jupyter_output><empty_output><jupyter_text>Define Data Location<jupyter_code>analysis.data_init(default_data_dir=DATA_DIR,
                   default_ref_data_dir=REF_DATA_DIR)<jupyter_output><empty_output><jupyter_text>zCFD Validation and Regression<jupyter_code>validate = True

if (validate):
    tol = {}
    tol['cl'] = 0.01
    tol['cd'] = 0.0003
    valid = {}
    valid['cl'] = {}
    valid['cl']['0'] = {}
    valid['cl']['10'] = {}
    valid['cl']['15'] = {}
    valid['cl']['0']['sa'] = 0.0
    valid['cl']['10']['sa'] = 0.10778E+01
    valid['cl']['15']['sa'] = 0.15060E+01
    valid['cl']['0']['sst'] = 0.0
    valid['cl']['10']['sst'] = 0.10778E+01
    valid['cl']['15']['sst'] = 0.15060E+01 
    valid['cd'] = {}
    valid['cd']['0'] = {}
    valid['cd']['10'] = {}
    valid['cd']['15'] = {}
    valid['cd']['0']['sa'] = 0.82E-02
    valid['cd']['10']['sa'] = 0.127E-01
    valid['cd']['15']['sa'] = 0.22237E-01
    valid['cd']['0']['sst'] = 0.81E-02
    valid['cd']['10']['sst'] = 0.127E-01
    valid['cd']['15']['sst'] = 0.22237E-01<jupyter_output><empty_output><jupyter_text>Initialise Environment<jupyter_code>from zutil.plot import *
import math
# Rename colour import to prevent clash
col=cl<jupyter_output><empty_output><jupyter_text>Remote Data ConnectionThis starts paraview server on remote host and connects<jupyter_code>from zutil.post import pvserver_connect
if analysis.is_remote():
    pvserver_connect()<jupyter_output><empty_output><jupyter_text>Plot forces functions<jupyter_code>from zutil.post import cp_profile
from zutil.post import cf_profile
from zutil.post import get_fw_csv_data
from zutil.post import get_case_root, get_case_report, cp_profile_wall_from_file
from zutil.post import cf_profile_wall_from_file
from zutil.post import get_csv_data
from zutil.post import calc_lift_centre_of_action
from zutil.post import calc_drag_centre_of_action
from zutil.post import rotate_vector

def plot_theory(ax, filename, label_list):
    df = get_csv_data(filename)
    ax.scatter(df[0], df[1], color='grey', label='Gregory Re=3m free transition')
    label_list.append('Gregory Re=3m free transition')

def get_cl_alpha(filename):
    return get_fw_csv_data(filename,widths=[16,16])

current_plot = None

def cp_plot(data,pts_array,**kwargs):
    global current_plot
    ax = kwargs['axis']    
    cp_array = data.GetPointData()['cp']
    chord_array = data.GetPointData()['chord']
    current_plot, = ax.plot(chord_array, cp_array, label=kwargs['chart_label'])
    if 'data_table' in kwargs:
        data_table=kwargs['data_table']
        for ind in range(0,len(cp_array)):
            data_table.append([0,0,0,0,0,0,0,0,0])
            data_table[ind][0] = pts_array.GetPoints()[ind][0]
            data_table[ind][1] = pts_array.GetPoints()[ind][1]
            data_table[ind][2] = pts_array.GetPoints()[ind][2]
            data_table[ind][3] = chord_array[ind]
            data_table[ind][4] = cp_array[ind]
    
def cf_plot(data_array,pts,**kwargs):
    global current_plot
    ax = kwargs['axis']
    cf_array = data_array.GetPointData()['cf']
    cfmag_array = data_array.GetPointData()['cfmag']
    chord_array = data_array.GetPointData()['chord']
    current_plot, = ax.plot(chord_array, cfmag_array, label=kwargs['chart_label'])    
    if 'data_table' in kwargs:
        data_table=kwargs['data_table']
        for ind in range(0,len(cf_array)):
            data_table[ind][5] = cfmag_array[ind]
            data_table[ind][6] = cf_array[ind][0]
            data_table[ind][7] = cf_array[ind][1]
            data_table[ind][8] = cf_array[ind][2]
    
def plot_cp_profile(ax,label,file_root,alpha,beta, data_table,colour, alpha_label):
    chart_label = label + '$\mathbf{Re}$' + ' = 6m (fully turbulent)'
    force_data = cp_profile_wall_from_file(file_root,[0.0,1.0,0.0],[0.0, -0.5, 0.0],
                            axis=ax,func=cp_plot,chart_label=chart_label,
                            data_table=data_table,colour=colour)
    pforce = force_data['pressure force']
    fforce = force_data['friction force']
    pforce = rotate_vector(pforce,alpha,beta)
    fforce = rotate_vector(fforce,alpha,beta)
    total_force = [ pforce[i] + fforce[i] for i in range(0,3) ]
    pmoment = force_data['pressure moment']
    fmoment = force_data['friction moment']
    pmoment = rotate_vector(pmoment,alpha,beta)
    fmoment = rotate_vector(fmoment,alpha,beta)
    current_plot.set_label(label + ('$\mathbf{C_l=%.4f\/ C_d=%.4f}$ \n' % (total_force[2],total_force[0])))
    ax.set_title(r'NACA0012 $\mathbf{\alpha='+ ('%.1f \/' % alpha_label) 
                 + ('}$ \n'),
                 fontsize=ft.title_font_size, fontweight='normal', color =col.zeno_grey)
    return (pforce, fforce, pmoment, fmoment, force_data)

def plot_cf_profile(ax,label,file_root,alpha,beta, data_table,colour, alpha_label):
    chart_label = label + '$\mathbf{Re}$' + ' = 6m (fully turbulent)'
    force_data = cf_profile_wall_from_file(file_root,[0.0,1.0,0.0],[0.0, -0.5, 0.0],
                            axis=ax,func=cf_plot,chart_label=chart_label, 
                            data_table=data_table,colour=colour)
    pforce = force_data['pressure force']
    fforce = force_data['friction force']
    
    pforce = rotate_vector(pforce,alpha,beta)
    fforce = rotate_vector(fforce,alpha,beta)
    total_force = [ pforce[i] + fforce[i] for i in range(0,3) ]
    current_plot.set_label(label + ('$\mathbf{C_l=%.4f\/ C_d=%.4f}$ \n' % (total_force[2],total_force[0])))
    ax.set_title(r'NACA0012 $\mathbf{\alpha='+ ('%.1f \/' % alpha_label) 
                 + ('}$ \n'),
                 fontsize=ft.title_font_size, fontweight='normal', color =col.zeno_grey)
    return (pforce,fforce)

def cmy_move_moment_ref_pt(force_data):
    pforce = force_data['pressure force']
    fforce = force_data['friction force']
    pmoment = force_data['pressure moment']
    fmoment = force_data['friction moment']
    total_force = [ pforce[i] + fforce[i] for i in range(0,3) ]    
    total_moment = [ pmoment[i] + fmoment[i] for i in range(0,3) ]
    ref_point = [0.0, 0.0, 0.0]
    lift_centre,mzs0  = calc_lift_centre_of_action(total_force, total_moment, ref_point)
    drag_centre,mzs0  = calc_drag_centre_of_action(total_force, total_moment, ref_point)    
    new_ref_point = [0.25, 0.0, 0.0]
    cmy = - total_force[2] * (lift_centre[0] - 0.25) + total_force[0] * (drag_centre[2] - 0.0)
    return cmy<jupyter_output><empty_output><jupyter_text>Generate the plots<jupyter_code>from zutil.post import get_num_procs
import mpld3
from zutil.post import ProgressBar
pbar = ProgressBar()

plot = 1
angle = []
cl = {}
cd = {}
cdp = {}
cdv = {}
cmy = {}
angle = []

forces = {}
forces['cl'] = {}
forces['cd'] = {}
forces['cdp'] = {}
forces['cdv'] = {}
forces['cmy'] = {}

data_table = {}

# Analysis constants
alpha = 0.0
beta  = 0.0
reference_area = 1.0

# Create output directory if missing
if not os.path.exists('output'):
    os.makedirs('output')

for val_case in validation_dict:
    #print val_case

    comp_var = validation_dict[val_case]

    label_list = []

    fig = get_figure(plt)
    ax = fig.add_subplot(111)
    ax.grid(True)
    x_label(ax,'$\mathbf{x/c}$')
    y_label(ax,'$\mathbf{C_p}$')
    ax.axis([-0.05,1.05,comp_var['cpmax'],comp_var['cpmin']])
    set_ticks(ax)
    plot_theory(ax,os.path.join(analysis.data.ref_data_dir,'data',comp_var['cpexpt']),label_list);
    alpha = comp_var['alpha']
    data = get_cl_alpha(os.path.join(analysis.data.ref_data_dir,'data',comp_var['cpcomp']))
    ax.plot(data[0],data[1],label='CFL3D SST-V 897x257')
    label_list.append('CFL3D SST-V 897x257')

    for case in case_dict:
        case_name = case
        case_var = case_dict[case]

        if case_var['case'] == comp_var['case']:
            data_table[case_name] = []
            data_array = data_table[case_name]

            try:
                num_procs = get_num_procs(case_name)

                if num_procs != None:
                    label = str(case_dict[case]['label'])
                    alpha2 = alpha
                    
                    if label.endswith('RBF'):       
                        alpha2 = 0.0
                        
                    (pforce,fforce,pmoment,fmoment, force_data) = plot_cp_profile(ax,'zCFD '+label+' ',
                                                                                  os.path.join(analysis.data.data_dir,
                                                                                               get_case_root(case_name,str(num_procs))),
                                        alpha2,beta,data_array,'b',alpha)

                    angle.append(alpha)

                    total_force = [ pforce[i] + fforce[i] for i in range(0,3) ]    
                    total_moment = [ pmoment[i] + fmoment[i] for i in range(0,3) ]

                    cmy[case_name] = cmy_move_moment_ref_pt(force_data)

                    cl[case_name] = (alpha,total_force[2])
                    cd[case_name] = (alpha,total_force[0])
                    cdp[case_name] = (alpha,pforce[0])
                    cdv[case_name] = (alpha,fforce[0])
            except Exception as e:
                pass
    
    set_legend(ax,'best')#,label_list)
    fig.savefig('output/naca0012_case_'+str(val_case)+'_cp_profile.png',dpi=100)
    
    pbar+=5
    plot+=1

    label_list = []
    
    fig = plt.figure(figsize=(8,5),dpi=100, facecolor='w', edgecolor=col.zeno_orange)
    ax = fig.add_subplot(111)
    ax.grid(True)
    ax.set_xlabel('$\mathbf{x/c}$', fontsize=ft.axis_font_size, fontweight='bold', color=col.zeno_grey)
    ax.set_ylabel('$\mathbf{C_f}$', fontsize=ft.axis_font_size, fontweight='bold', color=col.zeno_grey)
    ax.axis([-0.05,1.05,comp_var['cfmin'],comp_var['cfmax']])
    set_ticks(ax)
    data = get_fw_csv_data(os.path.join(analysis.data.ref_data_dir,'data',comp_var['cfcomp']),widths=[27,27])
    ax.plot(data[0],data[1],label='CFL3D SST-V 897x257')
    label_list.append('CFL3D SST-V 897x257');
    
    for case in case_dict:
        case_name = case
        case_var = case_dict[case]

        if case_var['case'] == comp_var['case']:
            #print val_case, case_name
            try:
                num_procs = get_num_procs(case_name)
                
                data_array = data_table[case_name]

                if num_procs != None:
                    label = str(case_dict[case]['label'])
                    alpha2 = alpha
                    if label.endswith('RBF'):       
                        alpha2 = 0.0
                    f = plot_cf_profile(ax,'zCFD '+label+' ',
                                        os.path.join(analysis.data.data_dir,
                                                     get_case_root(case_name,str(num_procs))),
                                        alpha2,beta,data_array,'b',alpha)
            except Exception as e:
                pass
        
    set_legend(ax,'best')    
    fig.savefig(os.path.join('output', 'naca0012_case_'+str(val_case)+'_cf_profile.png'),dpi=100)
    
    pbar+=5
    plot+=1    
    
#fig.subplots_adjust(hspace=0.3) 
#fig.subplots_adjust(wspace=0.3) 
#fig.savefig("output/naca0012_cp_profile.png")
pbar.complete()
plt.show()
#from IPython.display import FileLink, display 
#display(FileLink('output/naca0012_cp_profile.png'))<jupyter_output><empty_output><jupyter_text>Compare $C_l$ versus $C_d$ curve<jupyter_code># Read Ladson data
from zutil.post import get_fw_csv_data
data = get_fw_csv_data(os.path.join(analysis.data.ref_data_dir, 'data', 'CLCD_Ladson_expdata.csv'),widths=[7,7,7])

import itertools
marker = itertools.cycle((',', '+', '.', 'o', '*','^'))
case_marker={}
for case in case_dict:
    if case not in case_marker:
        case_marker[case] = marker.next()

fig = plt.figure(figsize=(8,5),dpi=100, facecolor='w', edgecolor=col.zeno_orange)

ax = fig.add_subplot(111)
ax.set_title(r'NACA0012 $\mathbf{C_L}$ vs $\alpha$', 
             fontsize=ft.title_font_size, fontweight='normal', color=col.zeno_grey)
ax.grid(True)
ax.set_xlabel(r'$\alpha$ [degrees]', fontsize=ft.axis_font_size, fontweight='normal', color = col.zeno_grey)
ax.set_ylabel(r'$\mathbf{C_l}$', fontsize=ft.axis_font_size, fontweight='normal', color = col.zeno_grey)
ax.axis([-1.0,20,-0.2,1.8])
ax.scatter(data[0], data[1], color='green', marker="o", s=50, 
           label='Ladson ' + '$\mathbf{Re}$' + ' = 6m (fixed transition)')
val_angle=[]
val_cl_comp=[]
for val_case in validation_dict:
    case = validation_dict[val_case]
    val_angle.append(case['alpha'])
    val_cl_comp.append(case['clcomp'])
    
ax.scatter(val_angle,val_cl_comp,color='black',marker="x",s=50,
               label='CFL3D SST-V 897x257')

for case in cl:
    ax.scatter(cl[case][0],cl[case][1],color=col.zeno_orange, marker=case_marker[case], s=50,
            label='zCFD ' + str(case_dict[case]['label']))

set_legend(ax,'lower right')
set_ticks(ax)

fig.subplots_adjust(wspace=0.3) 
fig.savefig(os.path.join("output", "naca0012_clalpha.png"))
    
fig = plt.figure(figsize=(8,5),dpi=100, facecolor='w', edgecolor=col.zeno_orange)
ax = fig.add_subplot(111)
ax.set_title(r'NACA0012 $\mathbf{C_L}$ vs $\mathbf{C_D}$', 
             fontsize=ft.title_font_size, fontweight='normal', color=col.zeno_grey)
ax.grid(True)
ax.set_xlabel(r'$\mathbf{C_l}$', fontsize=ft.axis_font_size, fontweight='normal', color = col.zeno_grey)
ax.set_ylabel(r'$\mathbf{C_d}$', fontsize=ft.axis_font_size, fontweight='normal', color = col.zeno_grey)
ax.axis([-0.1,2,0,0.03])
ax.scatter(data[1], data[2], color='green', marker = "o", s=50, 
           label='Ladson ' + '$\mathbf{Re}$' + ' = 6m (fixed transition)')  
val_cl_comp=[]
val_cd_comp=[]
col_labels = [r'$\mathbf{\alpha}$', r'$\mathbf{C_L}$', r'$\mathbf{C_D}$']
table_data = []
for val_case in validation_dict:
    case = validation_dict[val_case]
    val_cd_comp.append(case['cdcomp'])
    val_cl_comp.append(case['clcomp'])
    table_data.append(['CFL3D SST-V 897x257', '%.1f' % case['alpha'], '%.4f' % case['clcomp'], '%.4f' % case['cdcomp']])
    
    
ax.scatter(val_cl_comp, val_cd_comp, marker="x", s=50, color='black', label='CFL3D SST-V 897x257')

for case in cl:        
    table_data.append([case_dict[case]['label'],'%.1f' % cl[case][0], '%.4f' % cl[case][1], '%.4f' % cd[case][1]])
    ax.scatter(cl[case][1], cd[case][1], marker=case_marker[case], s=50, color=col.zeno_orange, 
               label='zCFD ' + str(case_dict[case]['label']))
    
set_legend(ax,'upper left')
set_ticks(ax)

fig.subplots_adjust(wspace=0.3) 
fig.savefig(os.path.join("output", "naca0012_clcd.png"))

fig = plt.figure(figsize=(8,5),dpi=100, facecolor='w', edgecolor= col.zeno_orange)
ax = fig.add_subplot(111)

table_data = sorted(table_data, key=lambda x: x[1])
row_labels = [data.pop(0) for data in table_data]

ax.axis('off')
table = ax.table(cellText=table_data,
        rowLabels=row_labels,
        colLabels=col_labels,
        colWidths=[0.1,0.1,0.1],
        loc='best')


table.set_fontsize(14)
table.scale(1.5,1.5)
table_props = table.properties()
table_cells = table_props['child_artists']

for ii in xrange(len(row_labels)+1):
    for jj in xrange(-1,len(col_labels)):
        if (ii, jj) in table._cells:
            if jj == -1 or (ii == 0 and jj >= 0):
                table._cells[(ii, jj)]._text.set_fontsize(16)
                table._cells[(ii, jj)]._text.set_color(col.zeno_orange)
                table._cells[(ii, jj)].set_facecolor(col.zeno_grey)

plt.show()
#from IPython.display import FileLink, display 
#display(FileLink('output/naca0012_clcd.png'))
def save_image(file_root, label):
    renderView1 = CreateView('RenderView')
    renderView1.ViewSize = [1281, 993]
    renderView1.InteractionMode = '2D'
    renderView1.AxesGrid = 'GridAxes3DActor'
    renderView1.OrientationAxesVisibility = 0
    renderView1.CenterOfRotation = [0.6517804924152661, -0.5, 0.0651468835202671]
    renderView1.StereoType = 0
    renderView1.CameraPosition = [0.5761779123939095, -2775.4500822028735, 0.0627843028945997]
    renderView1.CameraFocalPoint = [0.5761779123939095, -0.5, 0.0627843028945997]
    renderView1.CameraViewUp = [0.0, 0.0, 1.0]
    renderView1.CameraParallelScale = 0.5865106403219323
    renderView1.Background = [0.0, 0.0, 0.0]
    n0012_pvd = PVDReader(FileName=file_root+'.pvd')
    cleantoGrid1 = CleantoGrid(Input=n0012_pvd)
    cellDatatoPointData1 = CellDatatoPointData(Input=cleantoGrid1)
    eddyLUT = GetColorTransferFunction('eddy')
    eddyLUT.RGBPoints = [0.0,    0.231, 0.298, 0.752, 
                         1241.0, 0.865, 0.865, 0.865, 
                         2484.0, 0.706, 0.016, 0.149]
    eddyLUT.ScalarRangeInitialized = 1.0
    eddyPWF = GetOpacityTransferFunction('eddy')
    eddyPWF.Points = [0.0,    0.0, 0.5, 0.0, 
                      2484.0, 1.0, 0.5, 0.0]
    eddyPWF.ScalarRangeInitialized = 1
    cellDatatoPointData1Display = Show(cellDatatoPointData1, renderView1)
    cellDatatoPointData1Display.ColorArrayName = ['POINTS', 'eddy']
    cellDatatoPointData1Display.LookupTable = eddyLUT
    cellDatatoPointData1Display.ScalarOpacityUnitDistance = 29.6
    Render()
    WriteImage(os.path.join("output", label))
    
from IPython.display import Image, display

for case in case_dict:
    case_name = case
    label = case_name+'_eddy.png'
    try:
        num_procs = get_num_procs(case_name)
        save_image(os.path.join(DATA_DIR,get_case_root(case_name,str(num_procs))),label)
        display(Image(filename=os.path.join("output", label), width=800, height=500, unconfined=True))
    except Exception as e:
        print e
        pass<jupyter_output><empty_output><jupyter_text>Cleaning up<jupyter_code>pxm = servermanager.ProxyManager()
pxm.UnRegisterProxies()
del pxm
Disconnect()

if analysis.is_remote():
    Disconnect()
    pass<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.15074568671410468, 'fraction_numerical': 0.045862169802124965, 'mean_word_length': 2.6246246246246248, 'pattern_counts': {'":': 0, '<': 37, '<?xml version=': 0, '>': 38, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 46, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '9', 'dir': 'jupyter-structured-clean-dedup', 'id': '432770', 'license': 'BSD-3-Clause', 'n_tokens_mistral': 7963, 'n_tokens_neox': 7375, 'n_words': 1092, 'path': 'naca0012/naca0012.ipynb', 'repo_name': 'zCFD/zCFD-tutorials', 'size': '20518'}