Having attracted much attention over the past decades, electrospray is a method for producing continuous stream of monodisperse droplets in ambient air. By exerting a potential difference between a nozzle ejecting a continuous stream of droplets and a substrate positioned directly below the nozzle, various EHD modes can be observed. Deionized (DI) water EHD modes were classified in several researches, thus, the main focus of current study is on electrospray of non-Newtonian fluids. Viscoelastic fluid depicting both viscous and elastic behaviour is chosen as the non-Newtonian fluid in this research in order to investigate the role of viscoelasticity on EHD modes. First, related researches are reviewed and governing equations including log-conformation reformulation, two-phase flow equations, leaky dielectric model and Navier-Stokes equations are examined. Viscoelastic solutions used in electrospray tests comprised of three different concentrations of aqueous PAA solutions whose physical properties were thoroughly measured. Experimental tests were conducted utilizing a syringe pump and a high-voltage source, subsequently, various EHD modes were photographed by a high-speed camera and obtained images were processed by MATLAB software to investigate the effect of flowrate, applied voltage and solution concentration on EHD modes, jet profile, droplet diameter and its corresponding distribution. By classifying various EHD modes into dimensionless diagrams, it can be concluded that if solution concentration or flowrate increases, dimensionless regions where dripping and bead mode are observed plummet while stable jet region grows. Furthermore, contrary to DI water, when applied voltage surges, instead of multi-jet mode stick-jet mode is observed where jet sticks to outer surface of nozzle and jet’s asymptotic thickness diminishes. Distinguishable from dripping mode in DI water, dripping mode in viscoelastic solutions is identifiable by the presence of constantly thinning viscoelastic filament. Additionally, it can be asserted that by rising electric field strength, droplet diameter reduces and droplet frequency escalates while droplet diameter distribution becomes wider. Also obtained by MATLAB image processing codes, curve fitted jet profiles reveal that by increasing viscoelasticity or exerted voltage, jet approaches annular region of nozzle and it becomes thinner in jet’s asymptotic region, nevertheless, rising flowrate thickens jet’s asymptotic radius. Thereafter, several numerical simulations were conducted with COMSOL multiphysics software where viscoelastic flow past a sphere, viscoelastic droplet, DI water jet and viscoelastic jet were modeled and solved. Used for verifying the authenticity of viscoelastic equations and their implimentation in software, simulation of viscoelastic flow past a sphere benchmark was also stabilized beyond previously suggested treshold. Next, circulatory flow inside droplet, pressure jump on droplet interface and viscoelastic filament thinning were all successfully modeled in viscoelastic droplet simulation and obtained results were compared with experimental data where good agreement was seen between two data sets. Also simulated numerically, circulatory flow inside cone, pressure contours, space charge density contours and surface plot of log-conformation terms are presented and discussed in subsection allocated to DI water and viscoelastic solution jet simulations. Afterwards, retrieved jet profiles from numerical simulations and experimental tests are plotted in a single diagram where little discrepancy was seen between two data sets and the overall effect of viscoelasticity on jet profile was similar between numerical and experimental results. Finally, yielded results are discussed and suggestions for future work are made. Keywords Electrospray, Electrohydrodynamic, Viscoelastic Fluid, Multiphase Flow, Leaky Dielectric Model, Numerical Simulation, Rheology