text: '''
  airfoil_module library + geo_module
  '''
  
  from steapy import *
  import os
  import numpy
  from numpy import *
  import math
  from scipy import integrate, linalg
  from matplotlib import pyplot
  from pylab import *
  
  
  '''
  Defining the characteristics for NACA airfoil
  For symmetrical airfoil leave m, p 
  '''
  
  x=1
  m=0.03
  p=0.4
  t=0.12
  c=1
  
  
  '''
  Choose the type of the airfoil
  SYMMETRICAL
      OR
  CAMBERED
  '''
  #coordinates_symmetrical(x, t, c)
  coordinate = coordinates_cambered(x, m, p, t, c)
  
  '''
  Comment or Incomment the file path for the specific type of airfoil
  
  '''
  # load geometry from data file
  #naca_filepath = os.path.abspath('NACA'+'_'+'00'+str(int(t*100))+'.txt') #SYMMETRICAL
  naca_filepath = os.path.abspath('NACA'+'_'+str(int(m*100))+str(int(p*10))+str(int(t*100))+'.txt') #CAMBERED
  with open(naca_filepath, 'r') as infile:
      x, y = numpy.loadtxt(infile, dtype=float, unpack=True)
      
  # plot geometry
  width = 10
  pyplot.figure(figsize=(width, width))
  pyplot.grid()
  pyplot.xlabel('x', fontsize=16)
  pyplot.ylabel('y', fontsize=16)
  pyplot.plot(x, y, color='k', linestyle='-', linewidth=2)
  pyplot.axis('scaled', adjustable='box')
  pyplot.xlim(-0.1, 1.1)
  pyplot.ylim(-0.5, 0.5);
  savefig('Geo_NACA'+'_'+str(int(m*100))+str(int(p*10))+str(int(t*100))+'.pdf')
  show()
      
    
  # discretize geoemetry into panels
  panels = define_panels(x, y, N=40)
          
  # plot discretized geometry
  width = 10
  pyplot.figure(figsize=(width, width))
  pyplot.grid()
  pyplot.xlabel('x', fontsize=16)
  pyplot.ylabel('y', fontsize=16)
  pyplot.plot(x, y, color='k', linestyle='-', linewidth=2)
  pyplot.plot(numpy.append([panel.xa for panel in panels], panels[0].xa),
              numpy.append([panel.ya for panel in panels], panels[0].ya),
              linestyle='-', linewidth=1, marker='o', markersize=6, color='green')
  pyplot.axis('scaled', adjustable='box')
  pyplot.xlim(-0.1, 1.1)
  pyplot.ylim(-0.5, 0.5);
  savefig('Disc_Geo_NACA'+'_'+str(int(m*100))+str(int(p*10))+str(int(t*100))+'.pdf')
  show()
  
  
  # define freestream conditions
  freestream = Freestream(u_inf=1.0, alpha=5.0)
  
  #solve the source contribution
  A_source = source_contribution_normal(panels)
  
  #solve the vortex contribution
  B_vortex = vortex_contribution_normal(panels)
  
  #build LHS and RHS
  A = build_singularity_matrix(A_source, B_vortex)
  b = build_freestream_rhs(panels, freestream)
  
  # solve for singularity strengths
  strengths = numpy.linalg.solve(A, b)
  
  # store source strength on each panel
  for i , panel in enumerate(panels):
      panel.sigma = strengths[i]
      
  # store circulation density
  gamma = strengths[-1]
  
  # tangential velocity at each panel center.
  compute_tangential_velocity(panels, freestream, gamma, A_source, B_vortex)
  
  # surface pressure coefficient
  compute_pressure_coefficient(panels, freestream)
  
  # plot surface pressure coefficient
  pyplot.figure(figsize=(10, 6))
  pyplot.grid()
  pyplot.xlabel('$x$', fontsize=16)
  pyplot.ylabel('$C_p$', fontsize=16)
  pyplot.plot([panel.xc for panel in panels if panel.loc == 'upper'],
              [panel.cp for panel in panels if panel.loc == 'upper'],
              label='upper surface',
              color='r', linestyle='-', linewidth=2, marker='o', markersize=6)
  pyplot.plot([panel.xc for panel in panels if panel.loc == 'lower'],
              [panel.cp for panel in panels if panel.loc == 'lower'],
              label= 'lower surface',
              color='b', linestyle='-', linewidth=1, marker='o', markersize=6)
  pyplot.legend(loc='best', prop={'size':16})
  pyplot.xlim(-0.1, 1.1)
  pyplot.ylim(1.0, -2.0)
  pyplot.title('Number of panels: {}'.format(panels.size), fontsize=16);
  savefig('Press_coeff_NACA'+'_'+str(int(m*100))+str(int(p*10))+str(int(t*100))+'.pdf')
  show()
  
  
  # calculate the accuracy
  accuracy = sum([panel.sigma*panel.length for panel in panels])
  print('sum of singularity strengths: {:0.6f}'.format(accuracy))
  
  # compute the chord and lift coefficient
  c = abs(max(panel.xa for panel in panels)
          - min(panel.xa for panel in panels))
  cl = (gamma*sum(panel.length for panel in panels)
        / (0.5*freestream.u_inf*c))
  print('lift coefficient: CL = {:0.3f}'.format(cl))
      
          
  
  # define a mesh grid
  nx, ny = 20, 20  # number of points in the x and y directions
  x_start, x_end = -1.0, 2.0
  y_start, y_end = -0.5, 0.5
  X, Y = numpy.meshgrid(numpy.linspace(x_start, x_end, nx),
                        numpy.linspace(y_start, y_end, ny))
  
  # compute the velocity field on the mesh grid
  u, v = get_velocity_field(panels, freestream, X, Y)
  
  # plot the velocity field
  width = 10
  pyplot.figure(figsize=(width, width))
  pyplot.xlabel('x', fontsize=16)
  pyplot.ylabel('y', fontsize=16)
  pyplot.streamplot(X, Y, u, v,
                    density=2, linewidth=1, arrowsize=1, arrowstyle='->')
  pyplot.fill([panel.xc for panel in panels],
              [panel.yc for panel in panels],
              color='k', linestyle='solid', linewidth=2, zorder=2)
  pyplot.axis('scaled', adjustable='box')
  pyplot.xlim(x_start, x_end)
  pyplot.ylim(y_start, y_end)
  #pyplot.title('Streamlines around a NACA 0012 airfoil (AoA = ${}^o$)'.format(alpha), fontsize=16);
  savefig('Velo_NACA'+'_'+str(int(m*100))+str(int(p*10))+str(int(t*100))+'.pdf')
  show()
                                                        
                                              
  # compute the pressure field
  cp = 1.0 - (u**2+v**2)/freestream.u_inf**2
  
  # plot the pressure field
  width = 10
  pyplot.figure(figsize=(width, width))
  pyplot.xlabel('x', fontsize=16)
  pyplot.ylabel('y', fontsize=16)
  contf = pyplot.contourf(X, Y, cp,
                          levels=numpy.linspace(-2.0, 1.0, 5000), extend='both')
  cbar = pyplot.colorbar(contf,
                         orientation='horizontal',
                         shrink=0.5, pad = 0.1,
                         ticks=[-2.0, -1.0, 0.0, 1.0])
  cbar.set_label('$C_p$', fontsize=16)
  pyplot.fill([panel.xc for panel in panels],
              [panel.yc for panel in panels],
              color='k', linestyle='solid', linewidth=2, zorder=2)
  pyplot.axis('scaled', adjustable='box')
  pyplot.xlim(x_start, x_end)
  pyplot.ylim(y_start, y_end)
  pyplot.title('Contour of pressure field', fontsize=16);
  savefig('Press_Field_NACA'+'_'+str(int(m*100))+str(int(p*10))+str(int(t*100))+'.pdf')
  show()                                                      
  