fl ow into the North Aegean Sea Modeling the buoyancy-driven Black Sea Waterout ScienceDirect

ORIGINAL RESEARCH ARTICLE

Modeling the buoyancy-driven Black Sea Water out flow into the North Aegean Sea

Nikolaos Kokkos, Georgios Sylaios *

DepartmentofEnvironmentalEngineering,DemocritusUniversityofThrace,Xanthi,Greece

Received30July2015;accepted29December2015 Availableonline20January2016

KEYWORDS

Hydrodynamicmodel;

Remotesensingdata;

Modelvalidation;

Potentialenergy anomaly;

NorthAegeanSea

Summary Athree-dimensionalnumericalmodelwasappliedtosimulatetheBlackSeaWater (BSW) outfluxand spreading over theNorth Aegean Sea, and itsimpact oncirculation and stratification—mixing dynamics. Model results were validated against satellite-derived sea surfacetemperatureand in-situtemperatureandsalinityprofiles.Further,themodelresults werepost-processedintermsofthepotentialenergyanomaly,f,analyzingthefactorscontrib- utingtoitschange.ItoccursthatBSWcontributessignificantlyontheThracianSeawatercolumn stratification,butitssignalreducesintherestoftheNorthAegeanSea.TheBSWbuoyancyflux contributedto thechangeoff in theThracianSeaby 1.23
10^3Wm^3 in thewinterand 7.9
10^4Wm^3 in the summer, significantly higher than the corresponding solarheat flux contribution(1.41
10^5Wm^3and7.4
10^5Wm^3,respectively).Quantificationofthef- advectivetermcrossingthenorth-westernBSWbranch(tothenorthofLemnosIsland),depicteda strongnon-linearrelationtotherelativevorticityofSamothrakiAnticyclone.Similaranalysisfor thesouth-westernbranchillustratedarelationshipbetweenthef-advectivetermsignandthe relativevorticity inthe Sporadessystem. Thef-mixingterm increasesitssignificance under strongwinds(>15ms^1),tendingtodestroysurfacemeso-scaleeddies.

#2016InstituteofOceanologyofthePolishAcademyofSciences.Productionandhostingby Elsevier Sp. z o.o. This is an open access article under the CC BY-NC-ND license (http://

creativecommons.org/licenses/by-nc-nd/4.0/).

PeerreviewundertheresponsibilityofInstituteofOceanologyofthePolishAcademyofSciences.

* Correspondingauthorat:Vas.Sofias12,67100Xanthi,Greece.Tel.:+302541079398;fax:+302541079398.

E-mailaddress:gsylaios@env.duth.gr(G.Sylaios).

Availableonlineatwww.sciencedirect.com

ScienceDirect

jo u rn al ho m e p age :w w w. el s ev i er. co m / lo c a te /o c e an o

http://dx.doi.org/10.1016/j.oceano.2015.12.003

0078-3234/#2016InstituteofOceanologyofthePolishAcademyofSciences.ProductionandhostingbyElsevierSp.zo.o.Thisisanopen accessarticleundertheCCBY-NC-NDlicense(http://creativecommons.org/licenses/by-nc-nd/4.0/).

1. Introduction

TheNorthAegeanSea(NAS)(Fig.1)actsasadilutionbasin, directlyaffectedbythebuoyancyoutflowofthelow-salinity, nutrient-rich Black SeaWater (BSW)exported throughthe Dardanelles Strait, at the north-eastern part of the basin (AndroulidakisandKourafalou,2011).AlongtheDardanelles, astrongly-stratifiedtwo-layersystemwithopposingflowsis formed,consistingoftheexitingtotheNASsurfaceBSWlayer and the outflowing into the Sea of Marmara bottom layer (Ünlüataetal.,1990).Thesignificantverticalshearstresses developedatthesystem'sinterfacefavortheupwardentrain- ment,thusreturningpartsofNASwaterbacktotheAegean (KarnaskaandMaderich,2008).Recentestimatesoftheannual BSWdischargereportvaluesrangingbetween38,820m³s^1 (Karnaska and Maderich, 2008) and 42,790m³s^1 (Tuğrul etal., 2002), with approximately 67% of its volume being directlytransferredfromtheMarmaraSeaand33%attributed toNASwaterentrainment.Afteritsoutflow,BSWoccupiesthe first50mofthewatercolumn,overlyingon thenorthward flowingwarmandhighlysalineLevantineIntermediateWater (LIW),between adepth of50m and 400m,andtheNorth Aegean Deep Water (NADW) from 400m to the bottom (Zervakis et al., 2000). Its impact on surface dynamics of NASexhibitsstrongseasonality,especiallyduringearlyspring, when the mean monthly outflux through the Dardanelles Straitsreaches450km³perseason(57,740m³s^1),asaresult oftheincreasedriverrunoffandprecipitationovertheBlack SeaandtheraisedlocalentrainmentfluxesattheDardanelles (Tuğruletal.,2002).

TheexitofBSWfromtheDardanellesproducesacyclonic flow,bifurcatingatLemnosIsland,withthesouthernbranch beingstrongerduringthesummer,undertheinfluenceofthe annual northerly Etesian winds, and the northern branch coveringtheentireThracianSeacontinentalshelf(Vlasenko etal.,1996;ZervakisandGeorgopoulos,2002).InwinterBSW coreconcentratesat thenortherncoast of LemnosIsland, whereit bifurcates primarily to the north-westand occa- sionallytothesouthwest,undertheinfluenceofnorth-east- erly(bora-type)gales.Themostprominentsurfacepatterns in NAS include the thermohaline BSW-LIW frontal zone,

affectedstrongly by seasonal meteorology(Sylaios, 2011), a semi-permanent anticyclone of variable strength and dimensions developed in the Thracian Sea, around Samo- thrakiandpossiblyImvrosIslands(Cordero,1999;Theocharis andGeorgopoulos,1993;ZervakisandGeorgopoulos,2002) andarapidlychangingcyclone—anticyclonesystem,covering theupper200mof theSporades Basin.Thislatter system appearssuppliedeitherbythehighersalinitywatersmoving fromthesouthernAegeanSeaorbythefresherwaterofBSW (Kontoyiannisetal.,2003).Othereddyfeaturesofvariable strengthandsizeappeardependenton BSWdischargeand ThermaikosGulffreshwateroutflows(Olsonetal.,2007).

Severalpreviousmodelingeffortsexist,usingmostlythe PrincetonOceanModel(POM),appliedundervariablespatial resolution grids and different boundary conditions for the BSWoutflux.Zodiatisetal.(1996)usedatwo-dimensional vertically-integrated model to simulate the synoptic-scale flowpatternsoftheNASandutilizedNOAA-AVHRRthermal imagesformodelvalidation.Korresetal.(2002)coupleda high resolutionatmospheric model with POMand assessed model'sskillinforecastingthesatellite-recordedseasurface temperature(SST)field.Kourafalou andBarbopoulos(2003) usedthemodifiedsigma-coordinatedPOMinahighresolution mode(2.5km
2.5km),nestedatitssouthernboundarytoa coarsegridMediterraneanmodelandaregional,intermediate resolution model. BSW discharge boundary condition was determined by the net monthly inflow-outflow budget at the Dardanelles, ranging from 5000m³s^1 in Decemberto 15,000m³s^1inJune.KourafalouandTsiaras(2007)applied POMina1.6km
1.6kmgridimposingatwo-layer(inflow- outflow) transportrate at the Dardanelles boundary. Tzali etal.(2010)appliedPOMinaseriesofBSWoutflowscenarios, toassesstheimpactofBSWinflowrateontheestablishedNAS circulation. Finally, Androulidakis and Kourafalou (2011) applied theHybrid CoordinateOcean Model (HYCOM) on a veryhighresolutiongrid(1/508
1/508)andintroducedanew mathematicalschemeforBSWoutflowrates.Modelvalidation involvedthecomparisonofmodelresultstosatelliteandinsitu data.

Thepresentpaperpresentstheresultsofathree-dimen- sionalhydrodynamicmodelforthethermohalineNAScircu- lation, placing particular emphasis on (a) the model validation procedure, by comparing surface model results withextensivesatellitedatasetsandCTDprofilesfromfield campaigns,(b)thestratification—mixingconditionsofBSW alongitsrouteoverNAS,(c)thequantificationofindividual termsofthepotentialenergyanomalygeneralequation,and (d) thedescription of thepermanent andsemi-permanent surface meso-scale patterns, developed in NAS and their inter-linktoBSWbuoyancyspreading.

2. Material and methods

The hydrodynamic Estuary and Lake Computer Model (ELCOM), a three-dimensional numerical model developed bytheCentreforWaterResearchattheUniversityofWestern Australia(HodgesandDallimore,2001)wasusedformodeling theNAS.Thismodelhasmostlybeen usedto simulatethe hydrodynamicsoflakesandreservoirs,butithasalso been Figure1 TheNorthAegeanSeaanditsphysiographicbasins.

appliedinthemodelingofbuoyancyfluxesandgeophysical flows at large-scale enclosed basins, estuaries and coastal seas,asthePersianGulf(Alosairietal.,2011),theNorthern AdriaticSea(Spillmanetal., 2007)andtheRedSea(Barry etal.,2009).Kopasakisetal.(2012)modeledsuccessfullyNAS watercirculationandpollutantsaccumulation,whileKamidis etal.(2011)validatedthemodeltosimulatethetransportand diffusion processes in Nestos River plume. Further, Sylaios etal. (2013) utilizedthe abovemodel-validated results to assessthealong-andcross-shorecirculationandthestratifi- cation-destratificationprocessesintheshallowandelongated ThassosPassage(ThracianSea).

ELCOMsolvestheReynolds-AveragedNavierStokesequa- tions,usingboththehydrostaticandtheBoussinesqapprox- imations. The numerical method applies a semi-implicit formulationonafinite-volumeframeworkusingrectangular CartesiancellsinanArakawa-Cgridstencil(Hodges,2000).

Horizontalgridspacingisfixed,whereastheverticalspacing mayvary as a functionofdepth, butremains horizontally uniformandfixedovertime.Free-surfaceevolutionineach gridcell'scolumnissolvedbytheverticalintegrationofthe conservation of mass equation for incompressible flow appliedtothekinematicboundarycondition.

ThenumericalschemefollowstheadaptedTRIMapproach (CasulliandCheng,1992)withmodificationsforscalarcon- servation, numerical diffusion, and implementation of a mixed-layerturbulenceclosure(Hodges,2000).Convective terms are calculated using a third order Euler—Lagrange scheme,whiletheULTIMATE-QUICKESTschemeisintroduced for the advection of scalars. The model produces the dynamics ofstratified waterbodies with externalenviron- mentalforcing,suchastidalforcing,windstresses,surface thermalforcingaswellasinflowsandoutflows(Hodgesand Dallimore,2001)andthusseemsappropriateforapplication in NAS. Heat exchange through the water's surface is governed by standard bulk transfer models corrected for non-neutralatmosphericstabilityeffects(ImbergerandPat- terson,1990).Themodeldoesnotfollowthecommonly-used eddy diffusivityapproach, but adoptsa unique 1D mixed- layermodelforcomputingtheverticalmixingofmomentum and scalars (Laval et al., 2003), an approach particularly properforstratifiedwaterbodiesasNAS.

2.2. NASmodelset-up

2.2.1. Bathymetryandcomputationaldomain

Modeldomainextendsfrom38.68Nto41.08Nandfrom22.58E to 27.08E,thus coveringbothshelfand deepareasof NAS (Fig. 1). The bottom bathymetry was digitized in a Geo- graphic Information System (MapInfo v10.0), to develop a Digital Elevation Model (DEM) from the Hellenic Navy 1:200,000 bathymetric chart. Model's computational grid wasdevelopedusingaspatialmappingtool(VerticalMapper forMapInfo)andtheapplicationoflineartriangulationinter- polation.Themodeledarea wasdiscretizedintoauniform high resolution horizontal grid consisting of 1km
1km orthogonal cells (0.0098
0.0098), allowing the accurate representation ofarea'scomplex bathymetry andtopogra- phy.Thewatercolumnateach horizontalcellwasdivided into20exponentiallystretchedlayers,consistingofthinner surface layers with gradually increasing thickness towards

the bottom. All variables were located at fixed z-levels (z-coordinate model), thus producing a computational domaincomprisedby3,109,280rectangularfinite-volumes.

2.2.2. Initialandboundaryconditions

Aconstantinitial conditionoverthegrid domain(T=108C andS=38.0psu)wasdeterminedfortheambientwaterof NAS. Model boundary conditions involved hydrologic, meteorological and tidal forcing. Hydrologic forcing was determinedonameanmonthlybasisbyapplyingriverdis- chargeandwatertemperatureforallmainriversofGreece and Turkey, as reported by Skoulikidis et al. (1998). BSW buoyancy flux was considered seasonally variable, with a mean annual value of 42,222m³s^1 or 1331km³yr^1, as computed by Tuğrul et al. (2002). BSW daily temperature wasextractedfromMODISSSTsatellitedatafortheSeaof Marmara.BSW salinityat the Dardanellesboundary varied seasonally,withameanannualvalueof28.03psu,following Türkoğlu(2010).

Meteorologicaltime-seriesforafour-yearperiod(2005—

2008)witha3hrstime-stepwereacquiredfromtheNOAA database (http://ready.arl.noaa.gov/READYamet.php) at 18
18spatialresolution,todescribetheatmosphere-ocean boundary conditions. Datasets included wind speed and direction,atmosphericrelativehumidity,cloudcover,baro- metric pressure, atmospheric temperature, precipitation andsolarradiationintensity.

Boundarytemperature,salinity,currentsandtidaleleva- tionsatthe southernmodelboundary wereobtainedfrom dailyMyOcean products(www.myocean.eu)usingtheMed- iterranean Sea Physics Re-analysis database (1987—2012).

Water temperature, salinity and velocity profiles were acquired at 72 layers with a daily time-step at a spatial resolutionof 0.06258 along thesouthern boundary line of themodel'scomputationalgrid.Initially,theseprofileswere verticallyspline-interpolated to match the layeredconfig- uration of ELCOM model. At a second stage, a horizontal interpolationalgorithm,appliedtoeachlayer,filledthegaps ofhorizontalmodel'sdiscretization.

TheNASmodelwasrunforfourconsecutiveyearswitha time-stepof3min.Yearone(2005)wasusedtodevelopand stabilizetheflowfield, especiallydueto BSWentry,while years2006—2008toobtain,compareandpost-processmodel results.Theseconsistedofthethree-dimensionalflowfield, the temperature and salinity fields and the sea-surface heightevolution.

2.3. Modelvalidation procedureandcriteria Modelvalidationwasachievedintwomanners:(a)bycom- paring sea-surface temperature model results at 44 regu- larly-spacedmodelgridcellswiththeSSTdataobtainedby the International Group for High Resolution SST (GHRSST, www.ghrsst.org), a platform combining daily updates of near-real-timelevel2data(SSTsatobservedpixels),level 3data(griddedinspace, butwithnogap-filling) andlevel 4 data (gap-free objective analysis data), with a 0.258
0.258 spatial resolution (Martin et al., 2012); and (b)bycomparingmodelproducedtemperatureandsalinity profiles with the 2006 CTD casts of summer NAS cruises (Sylaios,2011;Sylaios,unpublisheddata).

Toevaluatethevalidityofmodelresultswithsatelliteand in-situ data, a set of criteria was established (Jedrasik, 2005):

(a)thecorrelationcoefficient,r,definedas:

r¼covðoi;m_iÞ

sosm ¼o_im_iom

sosm ; (1)

where oiandmirepresenttheobservedand modeled pairedvalues(withi=1,2,...,N),respectively,ando,m andsoandsmtherespectivemeansandstandarddevia- tionsoftheobservedandmodeleddatasets.Thisisthe non-dimensionalmeasureofco-variationoftheobserved andmodeledvalues,reachingitsmaximumatunity.Such criterionisinsensitivetothetotalmodelbias.

(b)themeansquarederror,MSE,definedas:

MSE¼s²_o ð1r²Þþ s_m sor

2

þðmoÞ² s²_o

#

;

"

(2)

wherethe firsttermisthe squaredcorrelationcoeffi- cient;thesecondtermdescribestheconditionalmodel biasC,expressingthecorrelationbetweenmodelerror andthevalues simulated bythemodel;andthethird termtheunconditionalmodelbiasB,definedastheratio ofabsolutebiastothesquaredstandarddeviationofthe observations.ForincreasedmodelreliabilityMSEshould beminimum.

(c)the Nash-Sutcliffe effectiveness coefficientE, defined as:

E¼1MSE

s²_o : (3)

It is considered as a dimensionless transformation of MSE, also called a quadratic model score. Its value increasesuptounitywithincreasingmodelquality.

(d)thespecialcorrelationcoefficient:

Rs¼

ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 MSE

s²_oþm²: s

(4)

A perfect fit is reached when Rs=1, meaning that MSE=0.Rsbetterrepresentsthefittingofmodelsimu- lations to the observed ones, compared to the total squareerror(Jedrasik,2005),and

(e)themeanmodifiedmodelbias,definedastheratioof themeanmodeledtothemeanobservedcounterpart, indicatingtheover-orunder-estimationtendencyofthe simulations.

2.4. Modelresultspost-processing

To assessthe impactof BSW buoyancy flow pathover the near-surface water column (0—200m depth) of NAS, the potentialenergyanomalyf(PEA,inJm^3)wasconsidered throughmodelresults.PEA-valuesexpressthestabilityofthe water column, defined as the amount of energy per unit volumeneededtoverticallyhomogenizethewatercolumn, as(Simpson,1981):

f¼1 h

Z 0

hgzðrrÞdz; (5)

wheregisthe gravitationalacceleration[ms^2], zis the verticalcoordinatepointingupwardattheseasurface[m], r is the depth-averaged density [kgm^3], r is the local density [kgm^3] at depth z, and h is the total water depth [m]. By definition f is depth independent, but it varieshorizontallyandintime.Itoccursthatfiszerofora fully-mixedwater column,andpositiveforstablestratifi- cation.

AlthoughPEAexplainstheinstantstateofwatercolumnin termsof mixingandstratification,thetemporalchange of PEAmayexplaintheinteractedprocessesrelatedtomechan- icalmixing(windandtidalstirring)andstratifyingmechan- ismsasthesolarheatfluxandthefreshwaterbuoyancyflux (Simpsonetal.,1991).Therefore,thechangeofPEAintime maybeexpressedas:

@f

@t

 

TOTAL

¼dksr_ajW³j

h ekbr_wju³_bj h þagQ_S

2Cp

þ 1 320

g²h⁴ NZr_w

@r_w

@x

2

; (6)

wherethefirsttermrepresentsthewindstirringinfluence andthesecondtermrepresentsthetidalstirringeffect,both tending to decrease f over time. d and e are empirically determined mixing coefficients (d=0.039 and e=0.0038) expressingtheefficiencyintheconversionofwind andthe tidally-generated turbulent kinetic energy into potential energy(Lund-Hansenetal.,1996),ksandkbarethesurface and bottom drag coefficients (ks=6.4
10^5 and kb=2.5
10^3), ra is the air density (1.2kgm^3), rw is the water density [kgm^3], h is the water column depth [m] and W and ub represent the wind speed [ms^1] and bottomcurrentvelocity[ms^1].Thethirdandfourthterms expressthestratifyingsolarheatfluxandfreshwaterbuoy- ancy flux, respectively,tending to increase thef-value in time. a and Cp are thermal expansion coefficients (a=1.6
10^48C^1 at 98C and Cp=4.0
10³Jkg^18C^1), gisthegravityacceleration[ms^2],QSistheincidentsolar heatfluxattheseasurface[Wm^2]andNZisthevertical eddyviscositycoefficient[m²s^1],determinedasNZ=ghu, whereg=3.3
10^3anduthedepth-averagedtidalspeed [ms^1].

Inorderto derivetherelativeimpactofeachindividual termonthechangeofthepotentialenergyanomalyofNAS watercolumn,thecombinedEq.(6)wasintegratedintime, for a period of a month, using the mean daily values of respectiveparameters,thus:

f_TOTAL¼dksr_a 1 h

Z t

0jWj³dtþ ag 2Cp

Z _t

0

Q_Sdtþ 1 320

g²h⁴ NZr_w

Z t 0

@r_w

@x

2

dtekbr_w1 h

Z t 0

jubj³dt; (7) whereterms1to4refertotherelativeimpactofthewind, theincidentsolarradiation,theBSWbuoyancyfluxandthe tide on the total potential energy anomaly of the water column.Theterm(@rw/@x)referstothe model-computed depth-averaged (0—200m) density gradient between BSW exit point at the Dardanelles and each selected studiedsitethroughoutNAS.Theimpactofthetidalterm was considered as negligible forthe micro-tidalenviron- ment of NAS.

2.5. RelativecontributionofPEAequationterms alongameridionaltransect

To comprehendtheimpactof advectionofBSWpulseson NAS, the relative magnitude of eachterm in the general dynamic f-equation was quantified, utilizing the ELCOM modelresultsforthe upper130mlayer,correspondingto theupper5layers,alongameridionaltransect.Thegeneral f-equation derived by Burchard and Hofmeister (2008) reads:

@f

@t¼r|fflfflfflfflfflffl{zfflfflfflfflfflffl}hðufÞ

A

þg hrhr:

Z ₀

h~uzdz

|fflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflffl}

B

g h

Z ₀

h hh 2z



~urh~rdz

|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}

C

g h

Z ₀

h hh 2z



~w@z~rdz

|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}

D

þr₀ h

Z ₀

hPbdz

|fflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflffl}

E

r₀

2ðP_b^sþP_b^bÞ

|fflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflffl}

F

þg h

Z 0

h hh 2z



Qdz

|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}

G

þg h

Z 0

h hh 2z



rhðKhrhrÞdz

|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}

H

;

(8) whereuand~uarethedepth-averagedandthedeviationfrom the depth-mean horizontal velocity vectors, respectively [ms^1],~wisthedeviationfromdepth-meanverticalvelocity [ms^1],histheelevationofseasurfaceabovethemean[m], rand~rarethedepth-averagedandthedeviationfromthe depth-meandensities,respectively[kgm^3],r0istherefer- encedensity[kgm^3],Qisthesourcetermfordensitydueto heating [Wm^2], Kh is the horizontal eddy diffusivity [m²s^1], Pb is the vertical buoyancy flux (=(g/r0)KZ(@r/

@z)),whileP^S_b andP_b^barethesurfaceandbottombuoyancy fluxes[m²s^3],respectively,and5histhehorizontalgradi- entoperator.

In the above equation the described processes are:

f-advection (term A), representing the advection of potentialenergyanomalybythedepth-averagedhorizontal velocityvector;thedepth-meanstraining(termB),repre- sentingthestrainingofthedepth-averagedhorizontalden- sity gradient due to the deviation from the depth-mean velocity vector; the non-mean straining (term C), repre- senting thestraining producedby the deviationfromthe depth-averaged horizontal density gradient; the vertical advection (term D), produced by the deviation fromthe linear vertical velocity proceeding from the kinematic boundaryconditionsimposedonthewatercolumnsurface andbottom;theverticalmixing(termE),expressedbythe integratedverticalbuoyancyflux;thesurfaceandbottom densityfluxes(termF),bothincreasingf;theheatingdue to short-wave radiation (term G), representing inner sources or sinks of potentialdensity; and the divergence ofhorizontalturbulentdensityfluxes(term H),increasing fattheupperhalfofthewatercolumn.

For our analysis, only terms A, B, C, E, F and G were accountedforasthemostimportantfactorscontrollingthe stratification—mixingdynamics of thearea. FollowingHoi- tink et al. (2011), contributions from the temporal and spatialvariationsofverticalvelocityinthedeviationofmean density(termD)andtheimpact onstratificationfrom the horizontaldivergenceofhorizontalturbulentdensityfluxes (termH)wereconsideredasnegligible.

TermE(verticalmixing)wassimplifiedinto:

@f

@t



Term-E

¼CdGr₀ju³j

h ; (9)

whereCdisthedragcoefficient(=0.0025)andGisthemixing efficiencycoefficient(=0.04).

TermsFandGinvolving innersources ofdensitydueto absorption of solar radiation and the surface and bottom buoyancyfluxesappeardominatedbysurfaceheating.These terms could be parameterized according to Wiles et al.

(2006), involving meteorological data(derived from NOAA database) as the incoming short wave radiation, relative humidity,wind speedandseaalbedofor thederivationof theupwardheatfluxduetoevaporation,thelong-waveback radiationandthesensibleheattransfer,asfollowing:

@f

@t



TermsF&G

¼ ag

2CpðQ⁰STQuÞ; (10) whereCpistheheatcapacity(=4.0
10³Jkg^18C^1),aisthe thermalexpansioncoefficient(=1.51
10^4K^1),Q⁰_Sisthe sub-surfaceincidentsolarradiation,Tistheheatreleased whenradiationreachestheseabottom,andQuisthesumof upwardheatfluxduetoevaporation,long-waveback-radia- tionandsensibleheatflux(Qu=QE+QB+QC).Q⁰_S isdeter- minedbyQ⁰_S¼Q_Sð1AÞ,whereAistheseasurfacealbedo (=0.15)andQSistheincidentsolarradiation.Tiscomputed asa functionof theeffectivediffusive attenuation coeffi- cientforshortwaveradiation(k=0.3)andtotalwaterdepth, as:

T¼12ð1a1Þ

kh ð1e^khÞ; (11)

wherea1=0.55,representingthefractionofheatabsorbed atthethinsurfacelayer.Theupwardheatfluxtermswere determinedusing relativehumidity andwindvelocity data fromNOAAdatabasecombinedwithseasurfacetemperature datafromtheELCOMmodelandapplyingthestandardheat fluxformulas(Gill,1982).Allcalculationswereperformedin MATLAB10.

3. Results

3.1. Modelresultsvalidation

Atthefirstvalidationstage,daily-averagedSSTmodelresults werecomparedtosatelliterecordedSSTdailydataobtained fromtheGHRSSTdataset(n=365peryear).Foreachstation, theproducedscatterpointsfollowedcloselythediagonalin thesimulated—observeddataplots(Fig.2aandb),providing quite successful statistical measures. The spatially-mean validationcriteria,presentedinTable1,indicatethatmod- el's performance improves over time (years 2006—2008).

Meanconditionalandunconditionalmodelbiasesofthemean squarederrorforalltheseyearswerecomputedat0.040and 0.099, respectively. From years 2006 to 2008 the mean modifiedmodelbiasturnsfromover-toslightunder-estima- tionofthesatellite-derivedSST(Table1).Thisperformance varies spatially and temporally. In theThracian Sea slight winterandautumnover-estimationexists,changingtofair under-estimationin spring andsummer (Fig. 2c). Summer under-estimationappearsmoreprominentin theSporades

Basin(Fig.2d). Suchunderestimationcould alsobeattrib- utedtothedaily-averagingofmodeleddata,theimpactof relativelypoormeteorologicalconditionsresolutionandthe fixedhorizontal eddyviscosity and diffusivitycoefficients.

Fig.3illustratesthespatialvariabilityinthemodelvalidation statisticalmeasuresforyear2008.Highercorrelationcoeffi- cientswere obtainedin the Thracian Seaand Thermaikos Gulf, where the model slightly under-estimated observa- tions. Lower correlation but improved modified bias was achievedinthecentralpartsofNASandintheChiosBasin.

At the second validation stage, the model's ability to reproducetheverticaldistributionoftemperatureandsali- nityfieldswastested.Year2006CTDcastsinNAS(103sam- pling sites) were compared to temperature and salinity modelresults,upto200mdepth.In-situdatawereinitially layer-averaged, to match model discretization, and then directly compared to model outputs (n=457). Validation statisticsshowninTable1indicatearathergoodagreement and a fair under-estimation of in-situ observed variables.

Indicative comparative profiles of water temperature and salinityforsitesinThracianSeaandSporadesBasinareshown inFig.4.

3.2. Modelresultsdescription

Thesimulatedsurfaceflowexhibitssimilarpatternstothose produced by other recent NAS models (Androulidakis and Kourafalou, 2011; Androulidakis et al., 2012; Kopasakis et al., 2012) and in-situ observations (Eronat and Sayin, 2014;Sayinetal.,2011;Sylaios,2011).Inwinter(09February 2008), theDardanelles low salinity flow bifurcatesaround LemnosIsland,withitsnorthernbranchcrossingtheLemnos- ImvrosPassage,movingeastwardstowardsChalkidikiPenin- sula (Fig. 5). The flow around Athos Peninsula and into SingitikosBayagreestotheLagrangianobservationsofOlson etal.(2007).UnderEkman transport,thesouthernDarda- nelles branchflows around the southern coastsof Lemnos Island,turningwestwardstofeedaweakanticycloneinthe SporadesBasin.Thispatternisconsistenttotheobservations madebyOlsonetal.(2007)usingdrifterswhenreportedon the enhanced westward Ekman drift associated with the strong northerly winds. A strong anti-cyclonicflow with a largediameter(70km)isformedinthevicinityoftheBSW- LIWfront,tothenorth-westofLesvosIsland.

Table 1 Spatially-mean statistical measures for ELCOM modelvalidationintermsofSST.

Parameter 2006 2007 2008

Satellite-derivedSST

Correlationcoefficient,r 0.917 0.924 0.941 Meansquared-error,MSE 4.985 5.060 3.838 Nash-Sutcliffeeffectiveness

coefficient,E

0.774 0.779 0.818 Specialcorrelation

coefficient,Rs

0.993 0.993 0.995 Modifiedmodelbias,MMB 1.063 1.024 0.966 Temperature Salinity Year2006CTDcasts

Correlationcoefficient,r 0.859 0.794 Meansquared-error,MSE 6.410 1.383 Nash-Sutcliffeeffectiveness

coefficient,E

0.402 0.138

Specialcorrelation coefficient,Rs

0.990 0.999

Modifiedmodelbias,MMB 0.921 0.975

Figure2 ScatterdiagramsandtemporalvariabilitybetweenthemodeledandtheGHRSSTsea-surfacetemperaturedataduring 2008.Subplots(a)and(c)comparedataintheThracianSea(circle),and(b)and(d)intheSporadesBasin(diamond).

Inspring(08May2008)thewindinfluencelessensandas BSWoutflowremainsatsubstantiallevels(50,000m³s^1), its northwestern branch achieves speeds up to 1 ms^1, stronglyaffectingtheThracianSea(Fig.6).Thislowsalinity water (33.7—34.5) forms the Samothraki Anticyclone, a distinctive feature in Thracian Sea. Part of the BSW flow reachesthesouthwesternpartofThassosIsland,andbifur- cates into a branchmoving towards Strymonikos Gulf and anothertowardsAthosPeninsula.Similarfeatureswerealso producedbythemodelofAndroulidakisetal.(2012).

Acirculationsnapshotrepresentingthesummerconditions ofmoderatesouth-easternwindsanddiminishedBSWoutflow isshowninFig.7(05July2008).BSWflowfollowsawestward pathway between Lemnos and Imvros Islands, with speeds ranging from 0.7to 0.95ms^1. As flow moves northwards

towardsThassosIsland,itseparatesdevelopingatwincyclo- nic—anti-cyclonicsystemtothewestandnorth-westofSamo- thrakiIsland.Thistwincyclone—anticyclonesystemwasalso simulatedbyAndroulidakisetal.(2012)duringthesummer period(August2003).TheBSW-LIWfrontalzonewasfoundto the south and east of Lemnos Island, allowing the south- westernwatertransfertowardsSkirosBasin.Again,Androu- lidakisetal.(2012)pointsoutthestrongsouth-westernpro- pagation of surface waters and the presence of Sporades Anticycloneduringtheentiresummer.

3.3. PEA distributionandchange

ThePEA quantifies the deficit in thepotential energy due to stratification over a water column of 200m depth, as Figure3 Spatialvariabilityof(a)thecorrelationcoefficientand(b)themodifiedmodelbias,betweenthemodeledandtheGHRSST sea-surfacetemperaturedatasets.

Figure4 Modeledandobservedprofiles(surfaceto200mdepth)intheThracianSea(circle)for(a)watertemperatureand(b) salinityandintheSporadesBasin(diamond)for(c)watertemperatureand(d)salinity,insummer2006.

comparedtothatofthefullymixedwatercolumn.ThePEA distributioninNASincreasesduetotheBSWbuoyancyoutflux andthesolarradiationeffectanddiminishesunderthewind mixinginfluence.Fig.8presentsthePEAtemporalvariability atselectedsitesintheThracianSea,inLemnosBasin,tothe south of Lemnos Island and in Sporades Basin for year 2008.Results illustrate thePEA characteristic bell-shaped curveresultingfromtheincidentheatfluxseasonality.Inthe ThracianSeaandLemnosBasinsites,thesuddenPEAchanges are attributed to increased BSW bulges reaching the site (positivechange)andtowindimpactandtheLIWoriginated masses entrapped through eddies at the site (negative change). For example, the significant PEA reduction pro- ducedintheThracianSeaandLemnosBasinbetween26June 2008 and 04 July 2008 is attributed to the intrusion and entrapment of saltier Levantine origin water, while wind mixing(10ms^1)isresponsibleforasimilarsuddenchange on12October2008.MovingtothesouthofLemnosIslandand towards Sporades Basin, the impact of BSW on the water

column stratification reduces, and PEA values peak at 650Jm^3and580Jm^3,respectively.

Tounderstandthemechanismsresponsibleforthestrati- fication—mixingdynamicsoverNAS,theindividualtermsof Eq. (6) were explored and their relative contribution on f-change was accounted for. In the Thracian Sea, during January2008,meanwatercolumndensityof1028.11kgm^3, mean density gradient of 1.41
10^5kgm^4, mean wind speed of 5.8ms^1 and mean solar radiation of 92Wm^2 were considered. The wind impact on f-change may be estimatedat3.62
10^6Wm^3,theinfluenceofsolarradia- tionis an order of magnitude higher(1.41
10^5Wm^3), whilethebuoyancyfluxeffectinducedbyBSWcontributes themost(1.23
10^3Wm^3).DuringJuly2008,solarheat flux reached 330Wm^2 thus depth-averaged density reducedto 1026.80kgm^3,whilethedepth-averagedden- sity gradient decreased slightly at 1.18
10^3kgm^4. As windspeedsdecreased,thewindeffectonf-changereached 2.1
10^6Wm^3,theinfluenceofsolarheatfluxincreased approximatelybyafactoroffour(7.4
10^5Wm^3),while the BSW impact almost halved compared to winter (7.9
10^4Wm^3), butremained the dominant stratifica- tionmechanism.

Integratingtheaboveresultsoveramonthlyperiodduring 2006—2008,the relative contribution of eachterm on the total PEA may be examined. Table 2 presents the mean monthlyf-valuesforThracianSea,asobtainedfromcalcu- latedverticaldensitydistributionandtheuseofEq.(5),and the integrated monthly f-values as obtained through Eq.(7).Theresultsobtainedbybothmethodsappearrather comparative. Moreover, the relative contribution of each mechanism (solar radiation, buoyancy effect and wind impact)onthemean-monthlyPEAvalueisalsoshown.The solarradiationcontributiononthewatercolumnPEAdepicts an average value of 130.79Jm^3, ranging between 55.28Jm^3inJanuaryand55.86Jm^3inDecember upto 208.02Jm^3 in July. These findings appear in agreement withtheNorthAegeansea-atmosphereheatfluxesanalysis developedby Poulosetal.(1997)andAndroulidakis etal.

(2012).Inthesamemonth,theBSWbuoyancycontribution Figure7 SurfaceflowandsalinitypatternsintheNorthAege- anSeaasproducedbyELCOMmodelon05July2008.

Figure6 SurfaceflowandsalinitypatternsintheNorthAe- geanSeaasproducedbyELCOMmodelon08May2008.

Figure5 SurfaceflowandsalinitypatternsintheNorthAege- anSeaasproducedbyELCOMmodelon09February2008.

alsoreceivesits maximumvalue(611.98Jm^3),leadingto the most stratified water column havinga f-stratification valueof811.27Jm^3.Thewindimpactreachesitsmaximum f-contributioninJanuary(22.19Jm^3),producingarather well-mixed water column in the Thracian Sea, since the fTOTALvalueislimitedto94.94Jm^3.

3.4. QuantificationofPEA equationterms

Fig.9presents thetemporalvariationin themagnitude of the examined potential energy equation terms crossing theNorthernLemnos—the ThracianSeameridionaltransect during2008(A:40.4758N,25.2418E;B:40.0228N,25.2418E).

The df/dt term varies between 1.014
10^4Wm^3 and +1.389
10^4Wm^3(Fig.9a),whilethef-advectiveterm,

accountingfortheeffectofBSWbuoyancyoutfluxcrossingthe transectintheeast-to-westdirection(positivevalues),ranges from0.914
10^4Wm^3to+3.926
10^4Wm^3.Positive values prevail throughout the year, indicating that BSW is mostly transferred westward, towards the Thracian Sea (Fig.9b).Thiswestwardtransferredf-advectivefluxappears relatedtothepresenceofSamothrakiAnticyclone,aneddy systemwith 50—60km diameter spread over the Thracian Sea. The depth-mean straining term, explained as the effect of tidal shear on the vertically constant horizontal density gradients appears of lower magnitude than the f-advective term, ranging between 0.671
10^5Wm^3 and+0.440
10^5Wm^3(Fig.9c).Similarly, limitedis the variability of the non-mean straining (term C), while the vertical mixing term (term E, Fig. 9e), expressed as

Table2 TheThracianSeameanmonthlyPEA,ascomputed bymodelresultsandmonthly-integratedPEAand relativePEA- contributionofindividualtermsasderivedfromEq.(7)during2006—2008.

Months Meanmonthly

PEA[Jm^3]

f_TOTAL[Jm^3] f_SOLAR[Jm^3] f_BOUYANCY[Jm^3] f_WIND[Jm^3]

January 103.83 94.94 55.28 62.27 22.19

February 104.04 115.97 70.55 54.89 9.17

March 224.71 206.22 120.10 95.87 9.57

April 308.49 317.18 156.73 163.05 2.57

May 613.32 601.80 206.76 398.15 3.08

June 786.51 673.00 200.03 475.02 2.02

July 755.01 811.27 208.02 611.98 8.70

August 774.92 761.55 195.60 567.40 1.42

September 662.56 633.02 141.20 497.36 5.51

October 431.28 455.27 96.64 368.29 9.64

November 249.55 292.79 62.79 235.60 5.56

December 163.08 172.46 55.86 121.66 5.03

Figure8 Temporalchangeofthepotentialenergyanomalyat(a)theThracianSea(solidblackline),(b)theLemnosBasin(red dashedline),(c)thesouthofLemnosIsland(bluedottedline)and(d)theSporadesBasin(greendashed-dottedline)during2008.The mapshowstheexactlocationofstationsattheThracianSea(circle),theLemnosBasin(triangle),thesouthofLemnosIsland(star)and theSporadesBasin(diamond).(Forinterpretationofthereferencestocolorinthisfigurelegend,thereaderisreferredtotheweb versionofthisarticle.)

theintegratedverticalbuoyancyflux,fluctuatesbetweennear zeroand0.012
10^3Wm^3,obtainingincreasednegative valuesundertheinfluenceofstrongsouthernandnorth-east- ernwinds(15—20ms^1).Theimpactofsolarheating,dom- inatingtermsG&F,increasesduringthesummerobtaininga maximumvalueof0.155
10^4Wm^3andreducestonear zerovaluesinthewinter(Fig.9f).

Thevariabilityofthef-equationtermscrossinganindica- tivemeridionaltransectfortheBSWsouth-westernbranch(C:

39.5738N, 25.2418E; D: 39.1238N, 25.2418E), is shown in Fig.10.Thepositivesigninthef-advectiveterm(Fig.10b)

is related to the strengthening of the SW branch, with increasedflowupto0.6ms^1movingfromthesouth-eastern capeofLemnostowardsAgiosEfstratiosIsland.TheSporades anticyclone appears mostly supported by high positive f-advective values, as on 09 March 2008 when USW- branch=0.72ms^1, the fadvective term 0.60
10^4Wm^3 and an anticyclonic flow is formed in Sporades(Fig. 11a).

Similarconditions prevailon20March 2008(fadvectiveterm

0.83
10^4Wm^3),asshowninFig.11b.Underanegative sign,theSWbranchofBSWisweakorcompletelyabsent,and the flowat the southof Lemnos Island isnarrowed to the Figure9 Temporalvariabilityof(a)thedf/dt-term(
10^4Wm^3),(b)thef-advectiveterm(
10^4Wm^3),(c)thef-depth-mean strainingterm(
10^5Wm^3),(d)thef-non-meanstrainingterm(
10^5Wm^3),(e)thef-verticalmixingterm(
10^3Wm^3)and(f) thesurfaceheatingbuoyancyintroductionterms(
10^3Wm^3),crossingthenorth-westBSWbranch.Alltermsweresmoothedusinga low-pass-filterof0.2Hzcutofffrequency.Blacklinesinsubplot(b)representperiodsofSamothrakiAnticyclonepresenceinThracian Sea.

Figure10 Temporalvariabilityof(a)thedf/dt-term(
10^4Wm^3),(b)thef-advectiveterm(
10^4Wm^3),(c)thef-depth- meanstrainingterm(
10^5Wm^3),(d)thef-non-meanstrainingterm(
10^5Wm^3),(e)thef-verticalmixingterm(
10^3Wm^3) and(f)thesurfaceheatingbuoyancyintroductionterms(
10^3Wm^3),crossingthesouth-westBSWbranch.Alltermsweresmoothed usingalow-pass-filterof0.2Hzcutofffrequency.

coastal zone, with a westward direction. Strong negative fadvectivetermvaluesappear linkedto theoccurrenceof a cyclonic flow in NE Sporades Basin, as on 28 May 2008 (fadvective0.80
10^4Wm^3; Fig. 11c) and 02 June 2008 (fadvective1.30
10^4Wm^3; Fig. 11d). Sudden changes in the f-advective term sign occur mostly during August and September (Fig. 10b), leading to consequent changesonthevorticitysignofSporadesBasinflow.

4. Discussion

Inthispaper,ELCOMmodelrathersuccessfullysimulatedthe variabilityofsurfacecirculationandwatermassesdistribu- tionintheNorthAegeanSea,producingtherapidlychanging surfacemeso-scalepatternsandthewatercolumndynamics undertheBSW,meteorologicalandheatexchangeinfluence.

Asseenfrommodelresults,theBSWbifurcationatthenorth andsouthofLemnosIslandandtherelativestrengthofflow betweenthesebranches, determinesthevariability ofthe abovedescribedmeso-scalepatternsandthewatercolumn stratification-mixing processes, expressed bythe potential energy anomalyof the studiedsystem. The north-western BSWflowaffectsmostlythecirculationintheThracianSea, feedingsub-basinscale gyresandflowsalongtheThracian coastline(the CoastalCurrent) andbetweenLemnosBasin andChalkidikiPeninsula(theRimCurrent).Thestrongsouth- westernBSWflowenhancesthehorizontaldensitygradients acrosstheBSW-LIWfrontalzoneandaffectstheSkirosand Sporades Basins inducingcyclonic—anti-cyclonic flows (the SporadesEddy).Theabovedescribedresultsappearinagree- menttothecirculationpatternsdescribedexplicitlybyTzali et al. (2010), Androulidakis and Kourafalou (2011) and Androulidakisetal.(2012).

Figure 11 Surface water circulation in North Aegean Sea (a) on 09 March 2008, when the f_advective term of SW-branch

0.60
10^4Wm^3,relatedtoan anticyclonein SporadesBasin, (b)on20 March2008,whenthefadvectivetermof SW-branch

0.83
10^4Wm^3,relatedtoananticycloneintheSporadesBasin,(c)on28May2008,whenthef_advectivetermofSW-branch

0.80
10^4Wm^3,relatedtoacyclonicflowintheSporadesBasin,and(d)on2June2008,whenthefadvectivetermofSW-branch

1.30
10^4Wm^3,relatedtoacyclonicflowintheSporadesBasin.

The potential energyanomaly variabilityillustratedthe importanceofBSWoutfluxinthestratificationconditionsof thewatercolumn.Ourresultssuggestthatduringthewinter theimpactofBSW-inducedbuoyancyonwatercolumnstra- tificationishigherbyalmosttwoordersofmagnitudethan thatofsolarradiation.Duringspringandsummer,solarheat fluxgraduallyincreases, allowingBSWbuoyancyimpact as thedominantf-changeterm,determiningthewatercolumn stratificationintheThracianSea,butonlybeinganorderof magnitudehigherthanthesolarradiationeffect.Atthearea of Sporades, solar radiation and buoyancy effects appear comparableduringthesummerperiod,whileatChiosBasin f-temporalvariabilityshowslimitedmean-monthlyfluctua- tion and solar radiation dominates the water column dynamics.

TocomprehendtheimpactofBSWpulsesadvectiononthe North Aegean Sea, the general dynamic f-equation, as derivedbyBurchardandHofmeister(2008),wassolvedalong a meridional transect in the Thracian Sea, utilizing the ELCOMmodelresultsfortheupper130m.Thef-advective termcrossingtheNWbranchappearsmostlypositive,imply- ingthe westward transport of the buoyant jet. The term fluctuates strongly in time andreceives its highest values duringspring andsummer, related to theoccurrence of a meso-scale eddy spread between Samothraki and Thassos Islands.Basedonthemodeledflowfield,therelativevorti- cityof thissystem, knownas SamothrakiAnticyclone, was calculated,as:

z¼@v

@x@u

@y (12)

andthenassociatedtothef-advectiveterm.Relativevor- ticity in the Thracian Sea exhibits negative values (mean z=1.1
10^5s^1), indicating the anticyclonic nature of thecirculation.Afifth-orderpolynomialregressionbetween relative vorticity and cross-transect f-advection is shown (Fig. 12), indicating that Samothraki Anticyclone was fed bythe north-western branchof the BSW plume. Asimilar analysiswasalsoperformedbySoosaaretal.(2014)onthe anticycloniceddyofthesouthernGulfofRigaanditsrela- tionshiptothewindandhorizontaldensitygradients.Results depictedthatthisanticyclonewasmostlyfedbythebuoy- ancy field, being enhanced or reversed by the dominant

winds.Sylaios(2011)explainedthatundernortherlywinds, theanticyclonewaspushedtowardsthenorth-westofLem- nos Island, while under the influence of south to south- westerlywinds, itmoved to the north-westof Samothraki Island.Asshownbythetemporalvariationofthef—vertical mixingterm(Fig.9e),intheNorthAegeanSea,strongwinds (>15ms^1)tend to destroythe anticyclonicpattern,pro- motingverticalwatercolumnmixing.Indeed,thenegative peaksinthetemporalvariabilityofterm-E(20February08;

21April08;30May08;18June08;11July08;20July08;

14 September 08) correspond to short-to-medium scale stormswithwindspeedsbetween15and20ms^1.

Basedon modelresultsfor theThracian Sea,thebaro- clinic Rossby radius of deformation for Samothraki Antic- yclonewascomputed,as:

R1¼c1

jfj; (13)

wherec1isthegravitybaroclinicwavespeed,estimatedby:

c1¼1 p

Z 0

HNðzÞdz; (14)

with N(z) the Brunt-Väisälä frequency. The produced R1

values during the occurrence of increased positive values inthecross-transectf-advectivetermexhibitsameanvalue of7.50.3kmthroughout2008.Thisresultisinagreement with the findings of Androulidakis and Kourafalou (2011), duringtheirnumericalexperimental tests.

5. Conclusions

In this paper the three-dimensional model ELCOM was adapted, implemented and validated, aiming to derive thehydrodynamicfieldintheNorthAegeanSea.TheBlack SeaWateroutflowandspreadinggovernssurfacehydrody- namics,withtwodistinctbranches(NWandSW)atLemnos Island. Intensive post-processing on the validated model resultswascarriedout,aimingtolinkBSWbuoyancytransfer to NAS water column dynamics and surface meso-scale patterns.

The potential energy anomaly illustrated a significant spatial variability in NAS, due to the variable impact of BSW throughout the system's surface. PEA forcing analysis indicatedthatBSWinducesbuoyancycomparabletothesolar heatingimpactonthewinterwatercolumnstratificationof theThracianSea.However,inspringandsummer,theBSW influence in the Thracian Sea appears up to three times higherthanthecorrespondingsolarheatingeffect.InSpor- ades, solar radiation and buoyancy impact seem of equal importanceforsummerwatercolumnstratification.Analysis andquantificationoftheindividualPEA termsexhibitsthe impactofBSW-inducedbuoyancyonthesurfacemeso-scale patternsofNAS.Indeed,thef-advectivetermcrossingthe NW branch exhibits strong relation to the occurrence of Samothraki Anticyclone. A non-linear regression between eddy's relativevorticity andtheNWf-advective term was developed, explaining the impact of BSW NW branch on SamothrakiAnticyclone.

Similarly, a SW branch enhancement, indicated by the highly positivef-advective values, appearsrelated to the intensification of Sporades Anticyclone. On the contrary, Figure 12 Non-linear regression between the f-advective

term crossing the north-west BSW branch and the relative vorticityofSamothrakiAnticyclone.

negativevaluesinthef-advectivetermappearrelatedtothe occurrenceofacycloneinSporadesBasin.Thisvorticitysign changeoccurredmostlyduringAugustandSeptember2008, associated to thesuddenchanges in thef-advection term crossingtheSWBSWbranch.Asshownbythevariabilityof thef—verticalmixingterm,strongwinds(>15ms^1)tendto destroytheabovemeso-scaleeddysystems,thuspromoting verticalwatercolumnmixinginNAS.

References

Alosairi,Y.,Imberger,J.,Falconer,R.A.,2011.Mixingandflushingin the PersianGulf(ArabianGulf).J.Geophys.Res.116,C03029, 1—14,http://dx.doi.org/10.1029/2010JC006769.

Androulidakis,Y.S.,Kourafalou,V.H.,2011.Evolutionofabuoyant outflowinthepresenceofcomplextopography:theDardanelles plume(NorthAegeanSea).J.Geophys.Res.116,C04019,http://

dx.doi.org/10.1029/2010JC006316.

Androulidakis,Y.S.,Krestenitis,Y.N.,Kourafalou,V.H.,2012.Con- nectivityof NorthAegean circulationto the Black Sea water budget.Cont.ShelfRes.48,8—26,http://dx.doi.org/10.1016/

j.csr.2012.08.019.

Barry,M.E.,Houdré,F.M.,McAlister,A.B.,Botelho,D.A.,2009.Three dimensional hydraulicandwaterqualitymodellingoftheRed Sea:challengesandlearnings. In:Proceedingsofthe18thWorld IMACSCongressandMODSIM09InternationalCongressonModel- lingandSimulation,“InterfacingModellingandSimulationwith MathematicalandComputationalSciences”,4142—4148.

Burchard, H., Hofmeister, R., 2008. A dynamic equation for the potentialenergyanomalyforanalysingmixingandstratification inestuariesandcoastalseas.Estuar.Coast.ShelfSci.77(4),679— 687,http://dx.doi.org/10.1016/j.csr.2012.08.01910.1016/j.ecss.

2007.10.025.

Casulli,V.,Cheng,R.T.,1992.Semi-implicitfinitedifferencemeth- odsforthreedimensionalshallowwaterflow.Int.J.Numer.Meth.

Fluids15,629—648.

Cordero, S.G.,1999.The useofthermalsatellitedataindense water formationstudiesin theMediterraneanSea.J.Marine Syst. 20 (1—4), 175—186, http://dx.doi.org/10.1016/S0924- 7963(98)00081-5.

Eronat,C.,Sayin,E.,2014.Temporalevolutionofthewaterchar- acteristicsinthebaysalongtheeasterncoastoftheAegeanSea:

Saros,İzmir,andGökovabays.TurkishJ.EarthSci.23,53—66, http://dx.doi.org/10.3906/yer-1307-4.

Gill,A.E.,1982.Atmosphere—OceanDynamics.Int.Geophys.Ser., vol.30.AcademicPress,London,662pp.

Hodges, B.R., 2000. Numerical Techniques in CWR-ELCOM (Code Release v.1). CWR Manuscript WP1422 BH. Cent. Water Res., Univ.WesternAustralia,Perth,37pp.

Hodges,B.,Dallimore,C.,2001.EstuaryandLakeComputerModel:

ELCOMScience Manual CodeVersion 2.0.0.Cent.Water Res., Univ.WesternAustralia,Perth.

Hoitink, A.J.F., vanMaren, D.S., Hoekstra, P., 2011. Mixingand stratificationinatropicaltidalembaymentsubjecttoadistrib- utedfreshwatersource.J.MarineSyst.88(1),34—44,http://dx.

doi.org/10.1016/j.jmarsys.2011.02.015.

Imberger,J.,Patterson,J.C.,1990.Physicallimnology. In:Wu,T.

(Ed.),AdvancesinAppliedMechanics,vol.27,302—475.

Jedrasik,J.,2005.Validationofthehydrodynamicpartoftheecohy- drodynamicmodelforthesouthernBaltic. Oceanologia47 (4), 543—566.

Kamidis, N., Sylaios, G., Tsihrintzis, V.A., 2011. Modelling NestosRiverplumedynamicsusingELCOM.Desalin.WaterTreat.

33(1—3),22—35,http://dx.doi.org/10.5004/dwt.2011.2627.

Karnaska,Y.,Maderich,V., 2008.Modellingofseasonalexchange flowsthroughtheDardanellesStrait.Estuar.Coast.ShelfSci.79 (3),449—458,http://dx.doi.org/10.1016/j.ecss.2008.04.019.

Kontoyiannis,H.,Kourafalou,V.H.,Papadopoulos,V.,2003.Seasonal characteristicsofthehydrologyandcirculationintheNorthwest AegeanSea(easternMediterranean):observationsandmodeling.

J. Geophys. Res. 108 (C9), 3302, http://dx.doi.org/10.1029/

2001JC001132.

Kopasakis,K.I.,Georgoulas,A.N.,Angelidis,P.B.,Kostovinos,N.E., 2012.Numericalmodelingofthelong-termtransport,dispersion, andaccumulationofBlackSeapollutantsintotheNorthAegean coastalwaters.Estuar.Coast.35(6),1530—1550,http://dx.doi.

org/10.1007/s12237-012-9540-9.

Korres,G.,Lascaratos,A.,Hatziapostolou,E.,Katsafados,P.,2002.

Towardsan ocean forecasting system for the Aegean Sea.J.

Atmos.OceanSci.8(2—3),173—200.

Kourafalou,V.H.,Barbopoulos,K.,2003.Highresolutionsimulations ontheNorthAegeanSeaseasonalcirculation.Ann.Geophys.21 (1),251—265,http://dx.doi.org/10.5194/angeo-21-251-2003.

Kourafalou,V.,Tsiaras,K.,2007.Anestedcirculationmodelforthe NorthAegeanSea.OceanSci.3(1), 1—16,http://dx.doi.org/

10.5194/os-3-1-2007.

Laval,B.,Imberger,J.,Hodges,B.R.,Stocker,R.,2003.Modeling circulation in lakes: spatial and temporal variations. Limnol.

Oceanogr.48 (3),983—994,http://dx.doi.org/10.4319/lo.2003.

48.3.0983.

Lund-Hansen, L.C., Skyum, P., Christiansen, C., 1996. Modes of stratification ina semi-enclosed bay at theNorth Sea—Baltic Seatransition.Estuar.Coast.ShelfSci.42(1),45—54,http://

dx.doi.org/10.1006/ecss.1996.0004.

Martin,M.,Dash,P.,Ignatov,A.,Banzon,V.,Beggs,H.,Brasnett,B., Cayula,J.-F.,Cummings,J.,Donlon,C.,Gentemann,C.,Grum- bine,R.,Ishizaki,S.,Maturi,E.,Reynolds,R.W.,Roberts-Jones, J., 2012. Groupfor HighResolutionSea Surfacetemperature (GHRSST)analysis fields inter-comparisons. Part 1:A GHRSST multi-product ensemble (GMPE). Deep-Sea Res. Pt. II 77—80, 21—30,http://dx.doi.org/10.1016/j.dsr2.2012.04.013.

Olson,D.B.,Kourafalou,V.H.,Johns,W.E.,Samuels,G.,Veneziani, M., 2007. Aegeansurfacecirculation from a satellite-tracked drifterarray.J.Phys.Oceanogr.37(7),1898—1917,http://dx.

doi.org/10.1175/JPO3028.1.

Poulos,S.E.,Drakopoulos,P.G.,Collins,M.B.,1997.Seasonalvari- abilityinseasurfaceoceanographicconditions intheAegean Sea(EasternMediterranean):anoverview.J. MarineSyst.13 (1—4), 225—244, http://dx.doi.org/10.1016/S0924-7963(96) 00113-3.

Sayin,E.,Eronat,C.,Uc¸kac¸,S¸.,Beşiktepe,S¸.T.,2011.Hydrographyof theeasternpartoftheAegeanSeaduringtheEastMediterranean Transient(EMT).J.MarineSyst.88(4),502—515,http://dx.doi.

org/10.1016/j.jmarsys.2011.06.005.

Simpson, J.H., 1981. The shelf sea fronts: implications of their existenceandbehavior.Philos.T.R.Soc.A302,531—546.

Simpson, J.H., Sharples, J.,Rippeth, T.P., 1991.A prescriptive modelof stratificationinducedbyfreshwaterrunoff.Estuar.

Coast.ShelfSci.33,23—35,http://dx.doi.org/10.1016/0272- 7714(91)90068-M.

Skoulikidis,N.T.,Bertahas,I.,Koussouris,T.,1998.Theenvironmen- talstateoffreshwater resourcesinGreece(riversandlakes).

Environ.Geol.36(1—2),1—17.

Soosaar,E.,Maljutenko,I.,Raudsepp,U.,Elken,J.,2014.Aninves- tigationofanticycloniccirculationinthesouthernGulfofRiga duringthespringperiod.Cont.ShelfRes.78,75—84,http://dx.

doi.org/10.1016/j.csr.2014.02.009.

Spillman,C.M.,Imberger,J.,Hamilton,D.P.,Hipsey,M.R.,Romero, J.R.,2007.ModellingtheeffectsofPoRiverdischarge,internal nutrientcyclingandhydrodynamicsonbiogeochemistryofthe NorthernAdriaticSea.J.MarineSyst.68(1—2),167—200,http://

dx.doi.org/10.1016/j.jmarsys.2006.11.006.

Sylaios,G.,2011. Meteorologicalinfluenceonthesurfacehydro- graphicpatternsoftheNorthAegeanSea.Oceanologia53(1), 1—24,http://dx.doi.org/10.5697/oc.53-1.057.