Air quality forecast systems need reliable and accurate representations of the planetary boundary layer (PBL) to perform well. An important question is how accurately numerical weather prediction models can reproduce conditions in diverse synoptic flow types. Here, observations from the summer 2014 HygrA-CD (Hygroscopic Aerosols to Cloud Droplets) experimental campaign are used to validate simulations from the Weather Research and Forecasting (WRF) model over the complex, urban terrain of the Greater Athens Area. Three typical atmospheric flow types were identified during the 39-day campaign based upon two-day backward trajectories: Continental, Etesians, and Saharan. It is shown the numerical model simulations differ dramatically depending on the PBL scheme, atmospheric dynamics, and meteorological parameter (e.g., 2-m air temperature). Eight PBL schemes from WRF version 3.4 are tested with daily simulations on an inner domain at 1-km grid spacing. Near-surface observations of 2-m air temperature and relative humidity, and 10-m wind speed are collected from multiple meteorological stations. Estimates of the PBL height come from measurements using a multiwavelength Raman lidar, with an adaptive extended Kalman filter technique. Vertical profiles of atmospheric variables are obtained from radiosonde launches, along with PBL heights calculated using bulk Richardson number. Daytime maximum PBL heights ranged from 2.57km during Etesian flows, to as low as 0.37km during Saharan flows. The largest differences between model and observations are found with simulated PBL height during Saharan synoptic flows. During the daytime, campaign-averaged near-surface variables show WRF tended to have a cool, moist bias with higher simulated wind speeds than the observations, especially near the coast. It is determined that non-local PBL schemes give the most agreeable solutions when compared with observations.