A numerical study of tides in Titan s northern seas, Kraken and Ligeia Maria

ContentslistsavailableatScienceDirect

Icarus

journalhomepage:www.elsevier.com/locate/icarus

A numerical study of tides in Titan ^ s northern seas, Kraken and Ligeia Maria

David Vincent

^a,∗

, Özgür Karatekin

^b

, Jonathan Lambrechts

^a

, Ralph D. Lorenz

^c

, Véronique Dehant

^b,d

, Éric Deleersnijder

^e,f

a Institute of Mechanics, Materials and Civil Engineering (IMMC), Université catholique de Louvain, Avenue Georges Lemaître 4, B-1348 Louvain-la-Neuve, Belgium

b Royal Observatory of Belgium, Avenue Circulaire 3, B-1180 Bruxelles, Belgium

c Space Exploration Sector, Johns Hopkins Applied Physics Laboratory, Laurel, MD 20723, USA

d Earth and Life Institute (ELI), Université catholique de Louvain, Croix du Sud 2, B-1348 Louvain-la-Neuve, Belgium ^e

e Université catholique de Louvain, Institute of Mechanics, Materials and Civil Engineering (IMMC) & Earth and Life Institute (ELI), Avenue Georges Lemaître 4, B-1348 Louvain-la-Neuve, Belgium

f Delft Institute of Applied Mathematics (DIAM), Delft University of Technology, Mekelweg 4, 2628CD Delft, The Netherlands

a rt i c l e i n f o

Article history:

Received 13 July 2017 Revised 30 November 2017 Accepted 11 December 2017 Available online 12 December 2017 Keywords:

Titan surface seas Tides liquid body Titan oceanography

a b s t r a c t

ThetidalresponseofTitan  stwolargestnorthernseas,LigeiaMareandKrakenMare,isstudiedbymeans ofatwo-dimensional,depth-averaged,shallow watermodel,SLIM(http://www.climate.be/slim).Kraken Mareisformedoftwobasins,thenorthernonebeingassumedtobelinkedbyasinglestrait,Trevize Fretum,toLigeiaMare.Thetidalmotionstendtobeindependentofeachotherineachbasin(i.e.,Ligeia Mare,Kraken1and Kraken2)whichresultsinsharptransitionsinthe straits.Ourresultsareoverall rather similar tothoseofTokanoetal. (2014), suggestingthat a2Dmodel suchas SLIMis adequate formodellingTitan  stidesand,sinceitis(presumably)lesscomputationallydemanding,maybebetter forsensitivitystudies.Forinstance,themaximumtidalrangeinKrakenandLigeiaMariarespectivelyare 0.29mand0.14m,whichiswithintherangepredictedbyLorenzetal.(2014)althoughsmallerby0.07m andlargerby0.04mthantheestimatesofTokanoetal.(2014)(butitoccursatthesamelocation).The tidal currentsarefaster(by aboutone orderofmagnitude)inthe straitslinkingthoseMariathanin thebasinsthemselves(withamaximumof0.36m/sinthestraitlinkingKraken1andKraken2,Seldon Fretum).Adecompositionofthetidalhistoryintodifferentharmoniccomponentsiscarriedout.Exceptin specificareassuchasthestraitsandtheamphidromicpoint(s),themaintidalcomponenthasaperiod of1TitanDay. Wealsobrieflystudiedthe eigenmodesofthe northernseaswhose periodiscloseto thetidalperiod:suchmodesareverylocal.Indeed,theirmagnitudeissignificant(withrespecttothe magnitudeofthemodeselsewhereintheseas)insmallbay(s)orneartheislandsofKraken2andLigeia Mare.Theyarenotexcitedbythetidesastheydonotappearinthetidalmotion.

Asensitivityanalysistopoorlyconstrainedparameters(bottomfrictioncoefficient,depthofTrevize Fretumand attenuationfactor– thelatteris brieflydiscussedwith respecttothe values oftheLove numbersfoundintheliterature)isalsoconducted.Themodelparametersareseentohaveasignificant impactontheliquidexchangesbetweenthebasinsand,consequently,onthetidalrangeandphase,fluid velocityandlocationofamphidromicpoint(s).

© 2017ElsevierInc.Allrightsreserved.

1. Introduction

Besides Earth, Titan is the only celestial body of the solar system on which liquid-filled surface lakes and seas have been detected.Liquidhydrocarbonlakesandseashavebeenassumedto be presentonTitan^s surfacesinceVoyager1established thatthe surfaceconditionsofthisicymoonwerecloseto thetriplepoint

∗ Corresponding author.

E-mail address: david.vincent@uclouvain.be (D. Vincent).

ofmethane(Haneletal.,1981;Samuelsonetal.,1981).TheCassini spacecraft,which has been gazing at the Saturnian systemsince 2004,observed a methanecyclesimilar to the hydrologicalcycle on Earth (Atreya et al., 2006). It also detected possible surface lakes in Titan^s southern hemisphere in 2004 (only some years laterconfirmedtobeliquid)bymeansofitsImagingScienceSub- system(McEwen etal., 2005),andsurfacelakesandseas¹ inthe

1 The distinction between lake (Lacus) and sea (Mare) is based solely on the size.

https://doi.org/10.1016/j.icarus.2017.12.018 0019-1035/© 2017 Elsevier Inc. All rights reserved.

Fig. 1. Bathymetry map of Kraken and Ligeia Maria. The depth scales with the dis- tance to the nearest shore, similarly to the bathymetry of Lorenz et al. (2014) . The red letters refer respectively Seldon Fretum (A), Trevize Fretum (B), the canyons next to Vid flumina (C) and Xanthus Flumen (D), and Moray Sinus (E). (For inter- pretation of the references to color in this figure legend, the reader is referred to the web version of this article).

northernhemisphere in2006bymeansofitsRadar(Stofan etal., 2007). The distribution of these surface lakes and seas is asymmetrical with respect to Titan^s equator (Aharonson et al., 2009): there are many more lakes and seas in the northern hemisphere,wheretheyarelargeranddeeper(Hayesetal.,2008).

ThesurfacelakesandseasarenottheonlyliquidbodiesofTitan:

a global subsurface ocean has also been inferred from Cassini measurements (Iess et al., 2012). Stevenson (1992) was the first to suspect the presence of a global subsurface ocean beneath an icy shell. Since then, thisassumption has been supported by severalobservations andmeasurements such asTitan^s obliquity, itsgravityfieldanditstopography(Rappaportetal.,2008;Nimmo andBills,2010;Balandetal.,2011;2014;Hemingwayetal.,2013;

Lefevreetal.,2014)andhasbeenwidely studied(e.g.Sohletal., 1995; 2003; 2014; Beuthe, 2015). The global subsurface ocean resultsin significant deformationsof the iceshell and, hence, of thelakesandseasbottom.Suchdeformationswillreducethetidal forcingstowhichsurfaceliquidsrespond,sincetheicecrustitself partlyfollowsthechangingpotential.Thiseffectismodelledbyan attenuationfactor,

γ

2, applied to the tidal forcing (Lorenz etal., 2014;Tokano etal., 2014;Beuthe, 2015). This factor depends on theinternal structure ofTitan andis particularlysensitive tothe thickness, the rigidity and the rheological properties of the ice shell(e.g.Sohletal.,2003;Lorenzetal.,2014;Tokanoetal.,2014).

As a consequence,

γ

2 should be considered as a free parameter (Lorenzetal.,2014).

In this study, we focus on two of the northern seas, Kraken Mare and Ligeia Mare. The former is the largest sea identified on Titan, withan area of at least 400, 000km². It is centred at (68°N, 50°E) (Tokano, 2010) and it stretches from 55°N to 82°N (Lorenzetal., 2014).It isformed oftwo basins, Kraken1(in the north)andKraken2,linkedbyastraitnamedSeldonFretum(iden- tifiedby “A” in Fig.1) andaset ofsmallchannels(Lorenzetal., 2014).Itisabout17kmwideand40kmlong(Lorenzetal.,2014) whichissimilartothedimensionsoftheGibraltarstraitonEarth.

LigeiaMareis locatednorth-eastofKrakenMare:itis centredat (79°N, 112°E). Its maximal dimensionsare about 420 by 350km (Stofan et al., 2012). Recently, those seas have been found to be

linked by a strait, namely Trevize Fretum (identified by “B” in Fig. 1) (Sotin et al., 2012). The Cassini Radar altimeter also de- tected liquidfilledcanyons linked to LigeiaMare inthe northern andsouth-easternregionsofthissea(Poggialietal.,2016)(identi- fiedby “C” and“D” inFig.1).Transientfeatures called“MagicIs- lands” wereobservedinKrakenMareandLigeiaMare(Hofgartner etal.,2014;2016).InLigeiaMare,suchphenomenawereobserved in regions centered at(78°N, 123°E) and (80°N, 111°E).They are best explained by ephemeral phenomenasuch asfloating and/or suspended solids,bubbles,andwaveswhile thediurnalphenom- ena,suchasthetides,seemtoberuled out,theobservationsbe- ing all atnearly thesametrue anomaly (Hofgartner etal.,2016).

In Kraken Mare, the phenomenon was observed at (73°N, 55°E).

In thissea,the same explanations stand butthe tidescannot be ruledout.WewillstudythetidesinducedbySaturnintheabove- mentionedMariaandwe willpayspecialattentiontospecific ar- eassuchasSeldonandTrevizeFretum,theaforementionedliquid filledcanyons,andatthe“MagicIslands” location.

Variousstudies were carriedout aboutTitan^sseas andlakes:

bathymetry (e.g. Hayes et al., 2010; Ventura et al., 2012; Mas- trogiuseppe et al., 2014;Lorenz et al., 2014; Hayes, 2016), liquid composition (e.g. Brown et al., 2008; Cordier et al., 2009; Glein andShock, 2013;Tanetal., 2013;2015;Luspay-Kutietal.,2015), dissipation due to friction in Titan^s lakes and seas (e.g. Sagan andDermott,1982;DermottandSagan,1995;Sears,1995;Lorenz etal.,2014),tidalresponse(Tokano,2010;Tokanoetal.,2014;Vin- cent etal., 2016), interactions with the atmosphere(Tokanoand Lorenz,2015;2016),andphenomenasuch astheshoreline varia- tions of OntarioLacus (Turtle etal., 2011) and thetransient fea- tures inthe northernseas Hofgartner etal. (2014;2016). Among them, Tokano (2010); Lorenz et al. (2014); Tokano et al. (2014); TokanoandLorenz(2015;2016);Hofgartneretal.(2014;2016)fo- cusedonthenorthernseas.

Tokano(2010)studiedthetidalresponseofKrakenMareusing a 3D hydrostatic, baroclinic ocean circulation model based on a structuredgridwithaspatialresolutionof30km.Duetothelack of information,heused a bathymetry scaling withthe east-west distance from the closest shore (the maximum and minimum depth were respectively500m and50m) andthe shoreline was not accurately represented (for instance, the connection with Ligeia Mare and the small channels near Seldon Fretum were ignored and Mayda Insula wasconsidered to be an island). The impact ofthe iceshell deformationswasnot takeninto account.

According to his results, the tidal amplitudeis minimum in the centre of the lakeand increases withthe distance to this point.

He predicted a maximum tidal range (difference between high tideandlowtide)of4minthenorthernbay ofKraken1andhe noticed that the tide rotates clockwise. The current is generally directed fromplaces where the tide is low to places where the tide is high and its magnitude is maximum in Seldon Fretum (about 0.3m/s) and offshore. Tokano (2010) took into account the seasonal effects through the solar heating: in winter, he predicted that the lakewould be well mixedwhile it would be stratifiedin summer unlessthe evaporation cools down thelake surface.Lorenzetal.(2014)derivedananalyticalexpressionofthe bathymetry and adapted Tokano (2010) results in order to take intoaccountthesurfacedeformationsbymultiplyingthetidalam- plitudebyanattenuationfactor

γ

2.Tokanoetal.(2014)improved the previous study of Tokano (2010): the bathymetry computed by Lorenz et al. (2014) was used, an attenuation factor

γ

2=0.1 wasimplemented to representthe impact of the iceshell defor- mations and the spatial resolution was 5.2km. Ligeia Mare and Kraken Mare were linked and thesmall channels next to Seldon Fretum were considered aswidestraits separatedby islands (see Fig. 1 inTokano etal. (2014)). The largest surface displacements predicted by the model were located in the south-eastern part

of Kraken 2 and northern part of Kraken 1 with a maximum tidal rangeof0.37m. Thecrucial role ofthe(unknown) depth of Seldon Fretumwaspointedout byLorenzetal.(2014)aswell as Tokano etal. (2014). According to the latter, changing the depth of Seldon Fretum from 6m to 100m completely changes the character ofKrakenMaretides.Forinstance,itresultsin aphase lagof90° andanincreaseofthemaximumtidalrangeof17cm.

Inthiswork,we willstudythetidalresponseofTitan^snorth- ern seas by means of the state-of-the-art model SLIM, (http://

www.climate.be/slim).Ithasalreadybeenusedtopredictthetidal response ofa Titan lake,Ontario Lacus (seeVincent etal., 2016).

SLIMsolvesthe2Ddepth-averagedshallowwaterequationsonan unstructuredgrid(seeSection 2.3),whichallowsan accuraterep- resentationofthe shorelineandahigherresolution inthestraits without the computational overhead of a regular but high reso- lution.Wewill alsofocusonthetidallyinduced liquidexchanges betweenthebasins.Variousparametersarepoorlyconstrained,no- tably the roughness of the bottom and the depth of the straits linkingthebasinsoftheMaria.Theirimpact(s)andthatofTitan^s surfacedeformation willbe quantified througha sensitivityanal- ysis. We willalso briefly study the impactof the bathymetry by implementinga bathymetrysimilartothat displayedinFig.9aof Hayes(2016) inLigeiaMare andby varyingthe slopeofthebot- tominKrakenMare.Finally,theimpactoftheshapeoftheshore- lineonthetidalresponsewillbestudiedbycomparingourresults withthoseobtainedpreviouslybymeansofothermodels.Forthe referencecase, thebathymetry implemented scales withthe dis- tancetotheshore(seeSection2.3)sothatthemaximumdepthin LigeiaMare is170mandtheattenuationfactorissetto

γ

2=0.1, similarlytoTokanoetal.(2014)).

ThegoalofthisstudyistodiscussthetidalresponseofKraken andLigeiaMaria.Thiscouldalsobe helpfulinenhancingourun- derstanding ofobserved transientevents andfor planningfuture missions. Furthermore,being able to accurately predict the tidal responseofthenorthernseascouldallowassessingsomeparame- tersfromfuturedata.Forinstance,ifwecouldmeasureaccurately the sea surfaceelevation, we could infer the value of

γ

2,which willgivefurtherinsightsintoTitan^sinternalstructure.

This articleisorganised asfollows.Section 2 briefly describes the numericalmodel used,the forcings aswell asthe mesh and bathymetry implemented. In Section 3,the tidal response of the seas is studied, with a focus on the liquid exchanges between the basins,the volumetricflowratethroughthe mainstraits,the canyons of Ligeia Mare, Moray Sinusand the magicislands phe- nomenon.Asensitivityanalysistotheattenuationfactor,thedepth of TrevizeFretum, thebathymetry scaling factor,the influenceof the artificial bathymetry, and Manning^s coefficient is conducted inSection 4.The resultsandtheapplications toandimplications forfuturemissionsarerespectivelydiscussedinSections5and6. ConclusionsaredrawninSection7.

2. Method

The modelis briefly introducedin Section 2.1(forfurther de- tails,seeVincentetal.,2016).Thetidalforcingappliedisdescribed in Section 2.2 while the computational domain, the bathymetry andthemesharedescribedinSection2.3.Thevalueoftheparam- etersofthemodelswhicharelikelytohaveasignificantimpacton theresultsarediscussedinSection2.4.

2.1. Depth-averagedmodel

The numerical model used is the same as in Vincent et al. (2016): it is the Second-generation Louvain-la- Neuve Ice-ocean Model,

SLIM

^{. This model relies on the discon-}

tinuous Galerkin finite element method (DGFEM) to solve the

2D depth-averaged shallow water equations on an unstructured mesh.An implicitRunge–Kutta schemeusing a Newton–Raphson solveris resortedto, whichallows fortime steps aslarge asfive thousandths of a Titan day (∼ 6890s) to be used. This method doesnot sufferfrom excessive numerical dissipation or spurious oscillations. Furthermore, it is highly parallelisable, local mass conservationisensuredandan efficientwetting-dryingalgorithm to deal with the tidal flats is implemented (Kärnä et al., 2011).

Solvingthe equationson unstructuredgrids isanotheradvantage of this method. Indeed, the mesh can be refined at places of particular interest (for example, in the vicinity of the shores or in the straits) without significantly increasing the computational cost.

SLIM

has already been successfully used to simulate tides invarious terrestrialdomains including the Scheldtestuary (e.g., deBryeetal.,2010),theMahakamdelta(e.g.,deBryeetal.,2011), andthe whole Great Barrier Reef (e.g., Lambrechts et al., 2008).

Ithasalso beused onTitan inorderto studyOntarioLacus(see Vincentetal.,2016).

Due to the absence of in situ measurements, models have to be used to determine the chemical composition and, hence, the density of Titan^s lakes and seas. The dominant compounds can bedeterminedfromphysicalpropertiesinferredfromCassinimea- surements such as the loss tangent and the dielectric permittiv- ity(see e.g. LeGall etal., 2016;Hayes,2016) anddifferentmod- elswere proposed (seee.g. Cordieretal.,2009; GleinandShock, 2013;Tanetal.,2013;2015;Luspay-Kutietal.,2015).Theyallpre- dicted that the seas and lakes are mainly composed of methane andethane,otherlow-molecular-masshydrocarbons,andnitrogen.

Nevertheless, depending on the model or some parameters, the polarlakesandseasappeartobe eitherethane-richormethane- rich,whichmodifiespropertiessuchastheliquiddensityandvis- cosity.Furthermore, the liquid composition varies with tempera- ture (Glein and Shock, 2013). Consequently, the liquid composi- tion can vary from one lake to another and over the year. For instance, according to Le Gall et al. (2016), Hayes et al. (2016), Ligeia Mare is methane rich and variation in backscattering be- tweenLigeiaMareandKrakenMareseemstosuggestthatKraken Mare would have moreethane than Ligeia Mare, as assumedby Lorenz(2014).ThiscouldinduceadensitygradientinTrevizeFre- tum.Nevertheless,duetothe sparseinformation aboutthecom- positionofKrakenMare,wedonottakeintoaccountsuchdensity gradientinthepresentstudy.Ontheotherhand,previousstudies haveshownthatsmallvariationsofthemeandensityhaveaneg- ligibleimpactonthetidalresponse (Vincentetal.,2016).Regard- ingtothetemporalvariationsoftheliquidcomposition,thosein- ducedbydailymodificationsofthesurfacetemperaturearenegli- gible,diurnalvariationsdisappearingatlatitudeshigherthan20°N (Cottini et al., 2012) while local variations induced by meteoro- logical eventsare beyond the scope of thisstudy, the associated timescalebeingmuchlargerthanthatofthetides.Therefore,the densityisassumedtobe spatially andtemporallyconstant inthe northernseas.Astheradarmeasurements(Hayesetal.,2016)and theindependentmicrowave radiometrydata(LeGalletal., 2016) indicatequiteclearlythatmethanedominatesthecomposition (at leastofLigeiaMare),thedensityissetto550kg/m³inagreement withTanetal.(2015).

Severalassumptions are made in order to derive the govern- ing equations solved hereinafter fromthe general massand mo- mentumconservationequations.First,densityvariationsaredisre- garded,asaforementioned.Second, the maximumdepth isabout 189m, which is much smaller than the horizontal length scale.

Such a small aspect ratiois also observed inSeldon andTrevize Fretum. The aspect ratiobeing small,the hydrostatic approxima- tionholdsvalid.Third,weneglecttheatmosphericpressuregradi- entaswe focuson tidalmotion.These assumptions allowforan

integrationoverthedepthresultinginthe2Dequationssolvedby ourmodel(seeEq.1ofVincentetal.,2016).

Eddyviscosityandbottom frictionin theshallowwaterequa- tionsare parametrised by means ofSmagorinsky^sclosure model (Smagorinsky, 1963) and an empirical Earth-based formula, the Chézy-Manning-Strickler^sformulation(see,e.g.,Lambrechtsetal., 2008).Accordingly,thebottomstressisestimatedasfollows:

τ

b=

ρ

^gn2

|

^u

|

^u

H^1/3

where

ρ

^{is the density;g}=1.352m/s² is the mean gravitational acceleration;uisthedepth-averagedvelocity;H=h+

η

^istheto-

talliquiddepthofthelakewherehisthereferenceheightofthe watercolumnand

η

^{istheseasurfaceelevation(positiveupward);}

andn=0.03sm^−1/3 is Manning^s roughnesscoefficient (it corre- spondstonaturalriverbottomonEarth).

The interactions withthe atmosphere are neglected:we take into account neither the precipitations, nor the evaporation, nor thewindstress.

2.2.Astronomicalforcings

The forcings taken into account are those induced by Saturn while those due to other moons andthe Sun are neglected. In- deed,thelatterisoneorderofmagnitudesmallerwhiletheother moonsaremuchsmallerthanSaturn,resultinginasmallerattrac- tionpotential.Theoceantidalloadingduetotheadjacentsea(s)is neglectedassuggestedbyTokanoetal.(2014).Thetidalforcingis obtainedfromthehorizontalgradientofthetidalpotential.AsTi- tanisatidallylockedmoon,thetwomaincontributionstotheas- tronomicalforcing arethatduetotheorbitaleccentricity(derived byDermottandSagan,1995)andthatduetoTitan^sobliquity(de- rivedbyTyler,2008),whoseperiodis1Titanday(TD).

Theidealtidalforcingduetothehorizontalgradientofthetidal potential would be completely effective on a totally rigid moon.

Nevertheless,Titanisnot rigid,itisevenhighlydeformablecom- paredto the Earth: the deformations of the ice shell above the global subsurface ocean are much larger than the deformations ofthe Earth crust. The solid tides of Titan andthe resulting sea bottomdisplacements are thus much moresignificant. The influ- ence of such solid tides on the tidal potential can be modelled throughtheattenuationfactor

γ

2=1+(^k2)− (^h2)^where(k₂) and(h₂)respectivelyaretherealpartoftheseconddegreetidal potentialLovenumberandtheseconddegreeradialdisplacement Love number(Sears, 1995; Sohl etal., 1995; Lorenz et al., 2014;

Tokanoetal., 2014;Beuthe, 2015). This attenuation factorrepre- sentsthefactthatTitansurfacedeformationswillreducethetidal forcingtowhichsurfaceliquidsrespond.Itisconsequentlyimple- mentedasafactormultiplyingthetidalforcingintheshallowwa- ter equations: the term representing the tidal forcing is

γ

2

∇

h

φ

where

∇

h isthehorizontaldeloperatorand

φ

^{isthetidalpoten-}

tial.Lovenumbersk₂ andh₂ respectivelycharacterizetheratioof thepotential duetotheiceshelldeformationstothetidalpoten- tialandtheratioofthesolidtideheighttotheheightofanequi- libriumocean.Suchafactorisalsousedtomodelsolid Earthand oceantides (e.g.Hendershott, 1972;Gordeevet al., 1977). Never- theless,thereareuncertaintiesabouttheseLovenumbers:theyde- pendonthesubsurfaceoceanthicknessanddensity,onTitan^sin- ternalstructure,ontherheologicalpropertiesofTitan(Sohletal., 2003),ontheiceshellrigidity(Lorenzetal.,2014;Rappaportetal., 2008), and they vary linearly with the ice shell thickness (Sohl etal., 2003;Rappaport etal., 2008).Consequently,we decidedto consider

γ

2asafreeparameter(asrecommendedbyLorenzetal., 2014). The rangeof acceptable values ofparameter

γ

2 is briefly discussedinSection2.4.Forthereferencecase,wechose

γ

2=0.1 tomatchwithTokanoetal.(2014).

2.3. Computationaldomainandbathymetry

We consider Ligeia and Kraken Maria independent from the other seas andlakes andwithoutexternal fluid input dueto hy- potheticalriver(s).TheshorelinecontourisdrawnfromtheRadar mosaicshowninLorenzetal.(2014).The resolutionoftheseim- ages vary dependingon the location: from 1 to 2km in Kraken Mare andfrom 0.3 to 1km in Ligeia Mare (Lorenz et al., 2014).

Particular attentionis paidtoSeldon andTrevizeFretum (respec- tively point A and B in Fig. 1). The latter is elongated (200km, five timeslonger than Seldon Fretum) and narrow (about15km, 2km lessthan Seldon Fretum). As recently discovered, there are liquid filled canyons in VidFlumina (Point C inFig. 1) and Xan- thus Flumen (Point DinFig. 1),near LigeiaMare (Poggiali etal., 2016). Nevertheless,astheydonotimpactthetidesinthenorth- ernseas,theyare notincludedinourdomainexceptforstudying thetidalmotiontherein(seeSection3.3).AtVidFluminalocation, thewidthoftheimplementedcanyonsislessthan1300m,which isslightlylargerthantheobservationofPoggialietal.(2016)(they predicted a maximum width of 1000m). Nonetheless, this value is only reached locally, at intersection, most of the canyons be- ing narrower than 1000m and the mean width is close to the value predicted by Poggiali et al. (2016) (700m). Near Xanthus Flumen, there is a strait and a small widening before the start of the canyon. The canyon maximum widthis about 6700m (at themouthofthecanyon),whichmatcheswiththeobservationsof Poggialietal.(2016).

Theunstructuredmeshesaregeneratedbymeansof

GMSH

(see Geuzaine andRemacle (2009)), which is one of themost widely usedopensourceunstructuredmeshgenerator.Thespatialresolu- tionin theseasranges from1kmto12km.We usethreerefine- mentcriteria:

• The local element size is proportional to the celerity of the longsurfacegravitywaves,c=



gh,assuggestedbyHenryand Walters(1993),Legrandetal.(2006).

• Themeshisrefinedneartheshores.

• ThespatialresolutionisincreasedinTrevizeFretumandSeldon Fretum(theminimumgridsizeisstill1kmbutthemaximum gridsizeis3kmintheseareas).

A bathymetry profile along the T91 flyby (in Ligeia Mare) was derived from Cassini measurements by Mastrogiuseppe et al. (2014) and a bathymetry map of Ligeia Mare isdisplayedin Fig.9aofHayes(2016).Nevertheless,dueto theopacityofKrakenMaretotheradarwavesnobathymetrypro- fileisavailableexceptinMoraySinus.Consequently,nomeasured bathymetry isavailable for KrakenMare. Lorenzet al.(2014) de- rived a hypothetical bathymetry with a simple bottom shape: a constant slope was assumed away from the shore. Nevertheless, due to variations in the shoreline contour, we could not use this bathymetry. Consequently, we build our own bathymetry similarly to the method used by Lorenz et al. (2014) (see Fig.

2 in Tokano et al., 2014): the depth scales with the distance to the nearest shore (i.e., depth=

α

× distance), the scaling factor,

α

^{, being computed in such a way that the maximum depth in}

LigeiaMareisthesameasinLorenzetal.(2014)(i.e.170 m)(see Fig. 1). It results insmall local variationsof thebathymetry due to smallmodifications of theshoreline orthe fact that an island is taken into account or ignored. For instance, the maximum depth of the implemented bathymetry is 189m and is located in Kraken 2 instead ofKraken 1 according to the bathymetry of Lorenz et al. (2014). We will briefly discuss to what extent an artificialbathymetryisjustifiableinSection4.4bycomparingour resultswith thoseobtained fora bathymetry ofLigeia Mare that mimicsFig.9aofHayes(2016).

2.4. Modelparameters

Several parameters are likely to have a significant impact on themodelresults.Asensitivityanalysiswillbeconductedwithre- specttofourofthem:theattenuationfactor

γ

2,thebathymetryof KrakenMare,thedepthofTrevizeFretumH,andManning^srough- nesscoefficientn.Thelatterisstudiedwithinarangefrom0.01to 0.06sm^−1/3,whichrespectivelycorrespondstosmoothman-made channel andtonaturalchannelswithstonesonEarth.Theresults arepresentedinSection4.3.

As explained in Section 2.3, the bathymetry is not derived from measurements but is assumed to scale with the distance to the shore, the scaling factor,

α

^{, being set to reach the same}

maximum depth in Ligeia Mare as in the bathymetry map of Lorenzetal.(2014).However,thereisnoevidencetosuggestthat the scaling factor is the samein Kraken Mare. Consequently, we consider two additional cases: in case1, we multiplied thescal- ing factorby 2(i.e.,depth=2

α

× distance) and, incase2,we di- vided it by 2 (i.e., depth=0.5

α

× distance). The former and lat- ter caseare referred to asbathy A and bathyB. Forbathy B,the maximum depth is170mand liesin LigeiaMare; themaximum depth inKraken 1andKraken 2respectively are 92mand 95m.

ForbathyA,themaximumdepthinKraken1andKraken2respec- tivelyare365mand379m.Itresultsinsignificantmodificationsof thevolume ratiobetweenKrakenandLigeiaMare.The variations ofthetidalresponseinducedbysuchmodificationsaredetailedin Section 4.5. Furthermore,aspointed out by Tokanoet al.(2014), thedepthofSeldonFretum(whichisunderconstrained)hasasig- nificantimpactonthetidalcurrentinsidethestrait.Itislikelythat the depth of Trevize Fretum, which is also not well constrained, hasasimilar impact.Consequently,weconsiderthe depthofthis strait as a free parameter. We carried out simulations for four cases: the referencecase andthree constant depths. Forthe ref- erence case, the depth increases with the distance to the shore, which resultsina maximumdepth of15m anda bottomprofile whichhastheshapeofa“∨”.Thethreeothercasescorrespondto a bottomprofilewhichhastheshapeofa“ࣶ” and whosedepths respectivelyare10m,20m,and50m.Theimpactofthisparame- terisdetailedinSection4.2.

TheattenuationfactorvalueisafunctionoftheLovenumbers k₂ andh₂ which,inturn,dependsontheinternalstructureofTi- tan. Unfortunately, thelatter is poorly constrained,which results in uncertaintiesaboutthese Love numbers.Iess etal. (2012) de- rived some valuesofk₂ fromCassini measurements by meansof three analysismodels whileSohl etal.(2003),Sohl etal.(2014), Baland et al. (2014), Lefevre et al. (2014), Beuthe (2015), among others,usedmodels toreconstructTitan^s internalstructure from observationsandmeasurementssuchasthegravityfield,thesec- onddegreetidalpotential Lovenumberk₂ andthemomentofin- ertia (whichwasrecentlyquestionedby Hemingwayetal.(2013), Balandetal.(2014)andLefevreetal.(2014)).Asaresult,different valuescanbefounddependingontheworkinghypothesesandthe modelused(seeTable1).

Beuthe(2015)predictedvaluesofk₂matchingwiththosecom- puted by Iess et al.(2012) from observations. By settingthe ice shell relative thickness in order to correspond with the value of k₂ obtained by Iess et al. (2012), the corresponding atten- uation factors range from 0 to 0.21 depending on the density of the crust and subsurface ocean (See Fig. 5 in Beuthe, 2015).

Tokanoetal.(2014) computed

γ

2 fordifferentinternal structures andonly retainedthose whichsatisfy themoment ofinertia and the Lovenumberk₂ predictedbyIessetal.(2012),whichresults in

γ

2∈[0.1,0.2]dependingontheiceshellthickness.However,the value retained fromthemoment ofinertia wasobtainedby hav- ing recoursetothe hydrostaticassumption,which isquestionable dueto thesignificant degree-three signal observedin thegravity

field (Hemingway etal., 2013; Lefevre et al., 2014;Baland etal., 2014). All these models assume a thin and homogeneous crust.

However, the crust could be non-homogeneous: there could be someclathratesoutgassingmethaneintheatmospherefromtime to time (Tobie et al., 2006). Such clathrates at the base of the iceshellmighthaveasmallershear moduluswhichwouldcause largerdeformationsoftheshellandtheLovenumberk₂wouldbe larger(Rappaportetal.,2008).

In the light of the values found in the literature, we decided tostudythetidalresponseforthreeadditionalvaluesof

γ

2:0.05, 0.2,and0.3(seeSection4.1).

3. Results

Inthissection,wefirstdescribethetidalresponseofLigeiaand KrakenMariapredictedby ourmodel.Forthispurpose,westudy the decomposition in different harmonic components as well as some of thenorthern seas eigenmodes. Then, we focuson some specificareassuchasTrevizeandSeldon Fretum,thecanyonsob- servedbyPoggialietal.(2016),MoraySinusandtheMagicislands locations.The bathymetryisthatshowninFig.1,theattenuation factoris

γ

2=0.1andManning^scoefficientissetto0.03sm^−1/3. 3.1. TidalresponseofKrakenandLigeiaMaria

Thetidalforcing(seeSection2.2)appliedrotatesanticlockwise withan exactperiodof1TD (notshown). Itis notunidirectional overtheseas.Indeed,therecanbeasignificantdifferenceinforc- ingorientationbetweenLigeiaMareandKrakenMareorevenbe- tween thesouthern and northernregion ofKraken Mare. Forin- stance, 0.25TD after perikron/apokron (point of the orbit where Titan isnearest/farthest to Saturn), theangle betweenthe direc- tionoftheforcinginthenorth-easternandsouth-westernregions can be larger than 90°. The forcing in t^∗ and t^∗+T/2 has the samemagnitudebutitisorientedinoppositedirections.InLigeia Mare, the maximum is about 3.59× 10⁻⁷m/s² and occurs 0.1TD before perikron/apokron in the eastern part of the sea while, in Kraken Mare, it is about 5.46× 10⁻⁷m/s² andit occurs 0.145TD after perikron/apokron inthe south-eastern partof the southern basin.

Themaximum/minimumseasurfaceelevation(

η

^)is ± 0.145m.

Itoccurs0.31TD afterperikron/apokron inthesouth-eastern part ofKraken2,0.165TDafterthemaximumforcingmagnitudewhich is directed south-westward (see the video inthe additional con- tent,availableatNorthernSeas_sse_reference_case.mp4).Themax- imum/minimumseasurfaceelevationinKraken1is ± 0.115mand occurs0.18TD before perikron/apokron atthe north-western end ofthebasin. In LigeiaMare, themaximum/minimum is ± 0.07m andoccurs0.05TDbeforeperikron/apokronattheeasternshoreof LigeiaMare.TheseasurfaceelevationinKraken2andLigeiaMare can be positive/negative over the whole basin, which suggests strong liquid exchanges through the straits (see Section 3.2). It lastslongerinLigeiaMare(about0.2TD)thaninKraken2(where theseasurfaceelevationispositive/negativeoverthewholebasin during respectively 0.095TD and 0.06TD) and does not occur in Kraken1dueto thepresenceofthe amphidromicpoint. Thepe- riodoftimeduringwhichtheelevationispositiveoverthewhole basin is not the same as that duringwhich it is negative which suggeststhat the volumetric flow rategoing through Seldon and TrevizeFretumisnotthesameatt^∗thanatt^∗+T/2.Thefactthat thetidalmotionisindependentineachbasinresultsinstrongsea surfaceelevationgradientinthe straits.Thisgradient islarger in the small channels near Seldon Fretum (the maximum is about 2.22× 10⁻⁵) than inSeldon andTrevize Fretum whereit can re- spectively reach about 1× 10⁻⁵ and 7× 10⁻⁶. The gradient in-

Table 1

Real part of the Love numbers (  ( k 2 ) and  ( h 2 )) and attenuation factor ( γ2 ) corre- sponding to the presence of a subsurface ocean (as established by Sohl et al. (2003) , Sohl et al. (2014) , Baland et al. (2014) , Lefevre et al. (2014) , Beuthe (2015) from ob- servations and measurements) found in the literature. k 2 and h 2 respectively are the second degree tidal potential Love number and the second degree radial displacement Love number.

Authors  ( k 2 )  ( h 2 ) γ2

Sohl et al. (1995) 0.36 1.19 0.17

Sohl et al. (2003) ^a [0.39, 0.32] [1.25, 1.05] [0.13, 0.27]

[0.39, 0.35] [1.28, 1.15] [0.12, 0.2]

Nimmo and Bills (2010) – 1.28 –

Iess et al. (2012) ^b 0.589 ± 0.075 – –

0.67 ± 0.09 – –

0.637 ± 0.112 – –

Sohl et al. (2014) 0.437 1.29 0.147

Beuthe (2015) ^c 0 . 7 − 0 . 5 1 . 7 − 1 . 3 0 − 0 . 2 0 . 7 − 0 . 41 1 . 7 − 1 . 2 0 − 0 . 21

a Various ice thickness and two ammonia concentrations.

b Different analysis models.

c Two crust densities and various relative thickness of the ice shell.

Table 2

Maximum amplitude ² of the first tidal component and relation between the amplitude of the other tidal components and that of the first tidal component. In the basins, these values are computed at three points: (80 °N, 110 °E) in Ligeia Mare, (70 °N, 50 °E) in Kraken 1, and (62.5 °N, 40 °E) in Kraken 2. In the straits, they are computed where the relation ^comp2_comp1 is maximum.

Period [ TD ] Ligeia Mare ^b Kraken 1 ^c Kraken 2 ^d Seldon Fretum Trevize Fretum Small channels

Comp 1 1 0.027 m 0.023 m 0.056 m 0.018 m 0.008 m 0.001 m

Comp 2 ^a 0.5 1% 1% 1% 96% 62% 378%

Comp 3 ^a 1/3 2% 4% 3% 3% 18% 58%

Comp 4 ^a 0.25 0% 0% 0% 21% 10% 118%

Comp 5 ^a 0.2 0% 1% 1% 1% 1% 48%

Comp 6 ^a 1/6 0% 0% 0% 4% 1% 54%

Comp 7 ^a 1/7 0% 0% 0% 1% 1% 34%

a Amplitude of the tidal component expressed as a percentage of the amplitude of the first tidal component.

b at (80 °N, 110 °E).

c at (70 °N, 50 °E).

d at (62.5 °N, 40 °E).

creasesanddecreasestwiceonatidalcycleandcanbeverysmall duringsome shortperiodsoftime (see,forexample,theseasur- faceelevationinthevideointheadditionalcontent).

While studyingtheseasurfaceelevationofseas,itiscommon to studythe tidal components.Indeed, the sea surface elevation ata point can be viewed asa sumof severalcomponents writ- tenas

η

(^t,x,y)=E(^x,y)^cos(

ω

^t+

φ

(^x,y))^{whereE(x,y)istheam-} plitude²,

ω

=2

π

/T isthe angularperiod,and

φ

^{(x,y) isthe tidal}

phase.On Titan, the first,second, third,and fourthtidal compo- nentsrespectivelyhaveaperiodof1,0.5,0.33,and0.25TD.Fourier analysis shows that the first tidal component (whose period is 1TD) isresponsibleformorethan 90%oftheseasurfacemotion, exceptforsomesmallareasinthestraitsandattheamphidromic point of the first tidal component (see Fig. 2). The second and fourthcomponentsare significant inthe straits(see Table 2) but canbe neglected elsewherewhile thethirdcomponent issignifi- cantinTrevizeFretum andwherethetidalrangeduetothefirst tidal component is weak (i.e. near the amphidromic point and in the south-western part of Kraken 2 (see Fig. 2)). The other tidalcomponentsareinsignificantincomparisontotheaforemen- tionedcomponents. In the small straits nearSeldon Fretum, the tidalrange ofthe second andfourthtidal componentsare larger thanthatofthefirst one.Nevertheless,itremains small(atmost 0.017m forthe second). The fact that the amplitude ofthe oth- erscomponentsissignificantwithrespecttothefirstoneinthese straits is due to the small amplitude of the latter at this loca-

2 The amplitude of a tidal component is half of the tidal range which is the dif- ference between high tide and low tide.

Fig. 2. Tidal range (in m ) and tidal phase (white lines, with a 15 ° spacing) of the first tidal component (whose period is 1 TD ) in Kraken and Ligeia Maria. The white arrows show the direction of the tide. The tidal range is the surface elevation dif- ference between high tide and low tide. There is an amphidromic point in Kraken 1 (at 65.213 °N, 66.459 °E) but not in Kraken 2 or Ligeia Mare. The tidal phase at each end of the Seldon and Trevize Fretum are not the same, which results in strong transition in the straits. This is also observed for the tidal range. Such phenomenon suggests independent tidal motion in each basin. The tidal range is maximum in the south-eastern part of Kraken 2.

tion(seeTable2).Intheremainderofthispaper,unlessotherwise specified, the expressions “tidalrange” and “tidal phase” refer to therangeandphaseofthefirsttidalcomponent.Anamphidromic point ofthe firsttidal componentislocated inthe south-eastern partofKraken1,at(65.21°N,66.47°E)(seeFig.2).Atthislocation, themostsignificant componentisthe thirdonebuttheresulting tidalrangeremainsweak:atmost2.5× 10⁻³m.Asaconsequence, thetidesinthisbasinrotatesaroundthispointandthetidalphase rangesfrom−180^◦ to180°.There isnosuch pointinLigeiaMare orKraken2.InLigeiaMare,thetidalphaseislargerthan−86^◦and smallerthan19° while,inKraken2itrangesfrom−208^◦to−52^◦. Inbothbasins, thecotidallines(linesofconstant tidalphase)are not rectilinear,whichisduetothefrictionwhichishigherinthe vicinityoftheshorelineandislandsastheseaisshalloweratthese locations.

Fig. 2 showsthat, in Kraken2, the tidalrange increases with the distance to a point on the south-western shore and ismax- imum (0.288m) at the south-eastern end. In Kraken 1, the tidal range increaseswiththe distanceto theamphidromic point and, hence,ishigheronthenorthern shore(maximum: 0.225m) than in the middle and southern part of the sea. The co-range lines (linesofsametidalrange)haveanovoidshapewiththelongitudi- nalaxisparalleltothemaximumwidthdirectionandwhosecen- terislocatedattheamphidromicpoint.InLigeiaMare,itincreases withthedistancetoapointlocatedon thesouthernshoreofthe sea, east ofthe Ligeia outlet ofTrevize Fretum and ismaximum (0.14m)inthesouth-easternpartofthesea.Thetidalmotionsre- main independent of each other in each basin. Indeed,a transi- tiontakesplaceineachstraitlinkingthebasins.Thistransitionis sharperinthestraitslinkingKraken1andKraken2thaninTrevize Fretum.Intheformer,therearebufferareaswherethetidalrange issmallerthanelsewhereinthestrait,whichresultsinsharptran- sitioninitsphase.Inthelatter,thetransitionissmootherforboth tidalrangeandphasealthoughitdoesnottakeplaceatthesame location(thetransitionoftidalphaseoccursmuchclosertoLigeia Marethanthatofthetidalrange).Thissuggeststhat,toalargeex- tent, thetidal motioninKraken1 andLigeiaMare canbe repre- sentedbyasloshingwavemode.Indeed,inthesebasins,thetidal phase isnearlyuniformoverlarge geographicalareaswhilethere are strong variations(sometimes by about 180°) over the straits connectingthosebasinswithKraken1(seeFig.2).

ThenormofthevelocityismaximuminthestraitssuchasTre- vize Fretum, Seldon Fretum, the small channels next to the lat- ter, and the straitnorth of Mayda Insula (see Fig. 3). Except for these specific locations, the largest speed is observed nearshore but the fluid velocity is one order of magnitudesmaller than in thestraits(about0.027m/s).InSeldonFretum,thefluidspeedhas two maxima located in the middle of the straits, 0.364m/s and 0.355m/s,whichoccurrespectively0.06TD afterperikron/apokron while the speed is the lowest (0.08m/s, which remains larger than everywhere else in the seas, excluding the straits) 0.15TD before perikron/apokron. The currentis mainly unidirectional to- wards Kraken 1 or 2 depending on the period of time (see Fig. 4(a–f)). The transition between two unidirectional flows oc- curs when the fluid speed at this location is the weakest. Once the flow is unidirectional, it accelerates to reach its maximum speed and then decelerates until it starts turning again. In Tre- vize Fretum, the maximum speeds, 0.235m/s and 0.234m/s, oc- cur respectively 0.3TD after perikron/apokron in the southern part of Trevize Fretum while the speed is the lowest (about 0.065m/s) 0.07TD after perikron/apokron. As in Seldon Fretum, the flow is mostly unidirectional (see Fig. 4(g–l)) and the flow turns when its speed is the weakest. The period of time dur- ing which the flow is unidirectional is longer in Trevize Fretum than in Seldon Fretum. These periods of time are concomitant

with each other as long as the flow is unidirectional in Seldon Fretum.

ThefactthatthecurrentisalmostunidirectionalintheSeldon andTrevizeFretum is also shownin Fig.5: the tidalellipses³ of themaintidalcomponentarenearlyrectilinearattheselocations, whichmeansthatthereisafasttransitionbetweenthetwomax- imumspeedsata givenpoint. Elongatedtidalellipses alsooccur elsewhereintheMaria butthevelocityismuch smaller.At loca- tionswheretheflowspeedissmallerthan0.001m/s(forinstance, inthewesternpartofKraken2,thesouthwesternpartofKraken 1,andthenorthernpartofLigeia(seeFig.3)),theellipsesarenear circularwhichindicatesthattheflowvelocityrotatesover1TDin- stead of being nil during the transition between two directions.

The orientation of the flow and, hence, of the tidal ellipses de- pendsonthelocation. Inthestraits,theellipsesarealignedwith themainaxisofthestraitwhile,intheirvicinity,theypointatthe nearest mouth of the strait. Offshore, where the speed is larger than 0.001m/s, the main axisof the ellipses is alignedwith the streamlinesatperikron/apokron(seeFig.3).

The tidalforcing studied in thispaperdoes not seemto gen- erate resonance in the northern seas (at least during the 15TD studied). Nevertheless, it does not mean that eigenmodes of the northernseascould notbe excitedby otherphenomena. Inorder to study the eigenmodes, we performed a discrete modal analy- sis of the Discontinuous Galerkin solution of the linear shallow water equations following the method of Bernard et al. (2008). The eigenvalues andeigenvectors are computedby means of the

linalg.eig

function of

scipy

packageofpythonwhichsolves the generalized eigenvalue problem M^dU_dt =AU where M is the massmatrix,A is thediscontinuousgalerkindiscretisation of the linear shallow water space operators for the domain, and U= (

η

,u,

v

)^{isthesolution.Theperiodsoftheeigenmodesrangefrom} 2.55× 10⁻³TDto15171TD.Nevertheless,mostoftheeigenmodes (90%)haveaperiodsmallerthan25TD(seeFig.6(a)).Someofthe eigenmodeshaveaperiodcloseto1TDbuttheseeigenmodescon- sistinlocalizedmodificationsoftheseasurfaceelevationinsome ofthesmallbays,alongthenorthernshoreofLigeiaMareornear theislandsofKraken2andLigeiaMare.Theselocalizedphenom- enadonotplayasignificantroleintheglobaltidesofthenorthern seasandtheseeigenmodesarenotnoticeablewhilelookingatthe tidalmotion. Some ofthe eigenmodes correspondto phenomena such as coastal trapped wavesor funnel shaped waves buttheir resonancefrequency is not excited by the tides.A complete dis- cussionofthenorthernseaseigenmodesisofinterestbutbeyond thescopeofthisarticle.

In Kraken 1, 2, and Ligeia Mare, the Rossby radius of defor- mation (R= √_gh

f where g=1.352m/s² is the mean gravitational acceleration atTitan^s surface, h is the depth, and fis the Cori- olis parameter) ranges are respectively [91.3, 1776.3]km, [98.5, 1938.2]km,and[93.5,1655.3]km.Thesevaluesarelargerthanthe characteristiclength scaleof thebasinsexceptforareas nearthe shorewheretheseasareshallower.Thismeansthatrotatingshal- low water waves such asKelvin waves could appear near shore, Rbeingsmallerthanthecharacteristiclength scaleattheseloca- tions.Other coastal trapped wavessuch astheedge wavescould occurinthenorthernseaswhilecontinentalshelfwavesareruled out because there is no shelf. Some of the eigenmodes of the northernseas correspondto suchcoastal trappedwavesbutthey are not excited by the tides as their periods are much smaller than the tidal period. Nevertheless, coastal trapped waves could appearduetoother forcingssuchasthewindbutthisisbeyond thescopeofthisarticle.Furthermore,theevolutionoftheseasur-

3 A tidal ellipse is the locus of the end of the velocity vector associated with one tidal component over a tidal period at one point.

Fig. 3. Tidal current streamlines at Perikron (a) and Apokron (b). The empty areas correspond to areas where the predicted velocity field magnitude is smaller than 0.001 m/s.

At these periods, the main current pattern is a flux going from the south-eastern part of Kraken 2 towards the north of Kraken 1 and from Moray Sinus towards Ligeia Mare (Panel (b)) or in the opposite direction (Panel (a)).

face elevationover one tidal period does not show evidence for anycoastaltrappedwaves(seethevideoofthetidalmotionover 1TD in additional content). Such waves are, hence,unlikely, un- lesstheyhaveaperiodsmallerthan0.05TD(∼ 6.9× 10⁴s)(which isthenumericaltimestep) orlargerthan0.5TD(∼ 6.9× 10⁵s)or they occur on a really small spatialscale. Topographic planetary wavescouldoccurif|

∇

hH|≥ HR0(LeBlondandMysak,1978)where Histhedepth, R₀=Rtan

φ

0 (RisTitan radiusand

φ

0 isthelati- tude), and |

∇

hH|∼ 0.002 is the bathymetry gradient. This condi- tionbeingnever encountered withthe bathymetry implemented, such waves are not expected. Offshore,the shallow water waves phasespeedtendstoc=



ghandtobenondispersive.InTrevize andSeldonFretum,capillarywavescouldoccur–thustidalcurrent couldgeneratesurfaceroughnessdetectableinradarobservations.

Indeed,Hayesetal.(2013)calculated aspeed forcapillarywaves (3cmwavelength)of0.11m/s,whichislowerthantheflowveloc- ityinthosestraits.Otherphenomena,whichcouldnotbedetected byspacecraftinstrument,suchassurfacegravitywavesandturbu- lence,could alsooccurintheseseasandmightplayarole inthe transferofenergyandmomentumtothesewaves.

It is noteworthy that the system of two basins connected by a channel is a classic Helmholtz oscillator. A terrestrial example ofsuch behaviour is the identification ofa 48–72hr periodicity in currents through the Straits of Mackinac, between the Great Lakes,MichiganandHuron(AndersonandSchwab,2013).Thetwo major Kraken basins and Seldon Fretum form an analogoussys- tem.Thenaturalperiodcanbe calculated(Eq.8ofAndersonand Schwab (2013)) from the channel cross-section (A_c∼ 0.2km², for

∼ 10m depth), the area of the basins (A_K1∼ 2.4× 10⁵km² and A_K2∼ 1.3× 10⁵km²), thechannel lengthl=40kmandgravityg= 1.352m/s² as T=2

π

^{ l}

Acg

 1 AK1+_A¹

K2

. With these values, the pre-

dictedperiodisabout8.8Earthdays.Asthisperiodisnotan in- tegerdivisorofthetidal period,itwouldtake severaltidalcycles toobserve anyresonance phenomenawiththe tides. The Ligeia- Trevize-Krakensystemhasanevenlongerperiod(about18.7Earth days) and would be more heavily damped by the long tortuous channel,sotheHelmholtzmodewillbe insignificantthere.Other

possibleresonancephenomenon couldbe thequarterwavelength resonance. The latter can be observed on shelves in the Earth ocean. The speed of long waves is givenby c∼



gh.For a tidal periodof1TD(about15.95Earthdays),theresultingquarterwave- lengthwidthrangesfrom283kmnearshoreupto5500kmatthe maximumdepth.Thisincreaseisnotlinear.Forinstance,thequar- ter wavelengthwidth is higher than 1400km at 1.7km offshore.

Consequently, such resonance phenomenon is unlikely in Titan^s northernseas.

3.2. Fluidexchangesbetweenthebasins

In this sub-section, we focus on the tidally induced fluid ex- changesbetweenLigeiaMare andKrakenMareandbetweenboth basinsofKrakenMare.Therearemuchmoreliquidexchangesbe- tweenbothbasinsofKrakenMarethanbetweenKrakenMareand Ligeia Mare: the maximum volumetric flow rate through Seldon Fretum is about three times larger than in Trevize Fretum (see Fig. 7) andthetotal volumetric flow ratethrough Seldon Fretum (inonewayoranother)overoneTitandayis32.7km³/TDwhileit is10.5km³/TDinTrevizeFretum.Theformermatchesthesemiana- lyticestimateofLorenzetal.(2014);theyconsideredawedge-like volume of liquid of 10–30km³ in each half period. Fig. 7 shows thatthereisanasymmetryintheflowratethroughbothstraits:it doesnotbehaveasasinusoid.Indeed,thebehaviourofthevolu- metricflowrateisnotthesameatitsmaximumandatits mini- mum.Thisisduetothefactthat theflow rateinthestraitsdoes notonlydependonthevelocitybutalsoontheseasurfaceeleva- tionatthislocation.Forinstance,inTrevizeFretum, theflowrate increases faster after its minimum than after its maximum (the slopeof thegraphissteeper).Thiscan be linkedto thefact that the sea surfaceelevation is positive in TrevizeFretum when the volumetric flow rate is maximum. Thislarger total depth of liq- uidresultsinaslowerdecreaseofthevolumetricflowratewhen the velocity decreases. Despitethisasymmetry, thedaily average isnil.Thevolumetricflowratesexplainthesmallervolumevaria- tion inLigeiaMare. Indeed,the latterisaboutthree timeslarger in both basins of Kraken Mare than in Ligeia Mare (see Fig. 8).

Fig. 4. Streamlines (a, b, c, d, g, h, i, j) and orientation (e, f, k, l) of the tidal current (in m/s) in Seldon and Trevize Fretum at four times of day: Perikron (a, e, h), T /4 after Perikron (b, i), Apokron (c, f, j), and T /4 after Apokron (d, k). The fluid speed is maximum (0.364 m/s) in Seldon Fretum 0.06 TD after perikron. The other areas where the speed is high are, in first instance, the other straits and, in second instance, the nearshore regions corresponding to specific shapes in the shoreline such as cape, peninsula and bay. The velocity in the strait is unidirectional for most of a Titan day (see Panels (e, f, k, l)).

Consequently, there isa strong correlation betweenvolume vari- ation in Kraken 1 and Kraken 2: the volume variation is maxi- mum in one basin soon after/before the minimum in the other andviceversa(see Fig.8).Thiscanbe explainedbythe factthat the tidal forcing magnitude is more significant in the surround- ingsofSeldon FretumthannearTrevizeFretumandbytheshape of the straits: TrevizeFretum is much more elongated than Sel- don Fretum andthe small channelsnext to Seldon Fretum allow foradditionalliquidexchanges.Thevolumevariationofeachbasin behaves asa sinusoid whoseperiod is 1TD.Although being sig- nificant, the volume variation is quite weak with respect to the

volume of the basins. Indeed, forthis bathymetry, the total vol- umeofLigeiaMare,Kraken1andKraken2arerespectivelyabout 4500km³,13000km³ and7200km³ while themaximum volume variationsoveraTitandaypredictedbyourmodelwiththesepa- rametersareabout2.5km³,8.9km³,and9.2km³respectively.

We also studied the fluid exchanges between the basins and Mariaover aperiod oftime of150TD (about 6.55Earthyearsor 0.22 Titan years) by means of a passive tracer: the tracer con- centration is set to 1 in Ligeia Mare and Kraken 2 and to 0 in Kraken1.Thismeans thatwherever thetracerislarger than0in Kraken1 (or smaller than1 in the othersbasins), there issome

Fig. 5. Major and minor axes of the tidal ellipses in the northern seas (Panel (a)), Trevize Fretum (Panel (b)) and Seldon Fretum (Panel (c)). The tidal ellipses represent the orientation of the first tidal component of the current over a whole period. Red and blue ellipses respectively indicate an anti-clockwise/clockwise rotation while the color scales with the magnitude of the maximum velocity over 1 TD . The smaller the minor axis, the more elongated the ellipses and the more unidirectional the velocity. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

Fig. 6. Periods of the northern seas eigenmodes (Panel (a), in TD ) with a zoom from 0 to 1 TD (Panel (b)) and from 0.99 to 1.01 TD (Panel (c)). The tidal period and its first fifth integer divisors are in red. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

liquidoriginatingfrom theadjacent basin(s). Such volumesgrow overtime. Theygrow faster inKrakenMare than in LigeiaMare.

Forinstance, after 150TD, at least 1% of the liquid at one loca- tionisoriginatingfromanotherbasininmostofKrakenMare(see Fig.9(a)).Theareaswherethisisnotthecasearealmostenclosed bayssuchasthatwestofMaydaInsula.Onthecontrary,inLigeia,

Fig. 7. Volumetric flow rate through Seldon (black line) and Trevize Fretum (blue line) over three Titan days. In Seldon and Trevize Fretum, it is respectively positive towards Kraken 1 and Ligeia Mare. The volumetric flow rate through the former is about three times larger than that through the latter. The net difference over 1 TD is 0 km ³ . (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article).

therearesometracesofliquidoriginatingfromKrakenMareonly inthesouthernpartofthesea,nearTrevizeFretum.Ifwelookat higherthresholds, we can seethat the liquidexchangesbetween LigeiaMare andKrakenMare remainsconfinedto TrevizeFretum andthe small channelsat Ligeia outlet. On the other side, there isa wide area (abouthalf ofthe basin)inboth basinsofKraken Marewhereatleast5% oftheliquidisoriginatingfromtheother

Fig. 8. Volume variation in the basins of Kraken and Ligeia Maria over three Titan days. The variation is larger in both basins of Kraken Mare than in Ligeia Mare and the total variation of each basin over a Titan day is nil.

basin,whichisduetothemuchmoresignificantliquidexchanges throughSeldonFretum(seeFig.9(b)).Theareawhereatleast10%

of the liquidis originating fromthe other basin remains located near Seldon Fretum in Kraken 1 whileit spreads over the west- ernpartofKraken2(seeFig.9(c)).Thiscanbeduetothesmaller volume ofliquidand smallerarea ofKraken2 (about54% ofthe volume of Kraken1). Unless a precipitation event(which is un- likely inKraken2andin thesouthern region ofKraken1) oran evaporation rate peak (which seems to occur during the spring) occurs, our model gives a good approximation of the tidally in- duced liquid exchanges between the basins of Kraken Mare. In- deed, the density driven circulation predicted in Seldon Fretum byTokanoandLorenz(2016)areeitherinsignificant(forthether- mallyforcedcirculation)oratleastoneorderofmagnitudesmaller (seeFig.14aofTokanoandLorenz(2016))thanthetidallyinduced flow(forprecipitation/evaporationdrivenflows).Furthermore,the formeroccurs onamuch largertime scale thanthe diurnaltime scaleofthetides.Ontheotherhand,therecouldbeamoresignif- icant difference incomposition betweenLigeia andKrakenMaria (Lorenz,2014).Astheshape ofthe straitalmost reducesthetidal flowtoabackanforthmotion,asignificantdensitygradientinthe straitcouldmodifytheflowinTrevizeFretumand,hence,theliq- uidexchanges.Butstudyingthisisbeyondthescopeofthepresent article.Furthermore,other factors could modifytheseresults.For example,inthecaseofapuretidalflowinandoutofLigeiaMare throughTrevizeFretum, thenthesameparcelofliquidiswashed back and forth and the net mixing is small: a parcel of liquid comes throughthechannel andsitsatthemouth,thenissucked back through. However, if there is, for instance, a strong wind- generatedgyreinLigeiaMare,thenwhateverKrakenMareliquidis hosedoutofTrevizeFretumintoLigeiaMareisquicklysweptaway andmixedintoLigeiaMare,andwhat issuckedback intoTrevize FretumisnowmostlyLigeiaMareliquid,notthepreviouslyKraken parcel. Consequently,wecan onlyassessthatthevariationofliq- uidcompositionobservedinMoraySinus(Hayesetal.,2011)with respecttoKrakenMarecannotbeexplainedbytidallyinducedliq- uidexchangesoveratimescalesmallerthan150TD.Thisvariation couldberelatedtoliquidexchangeswithLigeiaMareoveramuch longerperiodortoliquidexchangesinteractingwithwindgyre(s) (orsimilarphenomenon)atone(atleast)ofTrevizeFretumoutlets.

3.3. ThecanyonsofLigeiaMare

In this section, we investigate the tidal current and sea sur- faceelevationinthecanyons studiedbyPoggialietal.(2016).We show that the tidesdo not play a significant role infilling these

canyons.As no informationisavailable aboutthe liquiddepth in thesecanyons,thebathymetryisderivedsimilarlytothatofLigeia Mare.Near VidFlumina,itresults inamaximum depthof2.3m, the mean depth being about 1m while the maximum depth is about6.5mnearXhantusFlumen. Thisbathymetrydoesnottake intoaccounttheslopeofthebottomofthecanyons,thelatterbe- ing unknown:Poggialiet al.(2016)predictedthat the liquidsur- faceisatthesamelevelasinLigeiaMare withan errorofabout 0.7m.

The maximum sea surface elevation in the canyons is about 0.034m, which does not significantly modify the liquid depth at thislocationandcannotbedetectedbyanyCassiniinstrument.At thetime ofthe flybyanalysed by Poggiali etal.(2016) (T91,at a trueanomalyof

ν

=69^◦),thetidewashighinthecanyons:thesea surfaceelevationwasbetween0.009mand0.02minthecanyons near Vid Flumina while it wasbetween 0.027m and 0.031m in thoselocated next to XanthusFlumen. Unfortunately, such varia- tionsare toosmallto be observedby meansofthe Radarwhose precision along the vertical axis is about 0.7m (Poggiali et al., 2016).There isnodiscontinuitybetweenthetidesinthecanyons andinLigeiaMare:theliquidmotioninthecanyonsisinducedby thetidesofLigeiaMare.Inbothcanyons,theonlysignificanttidal componentisthefirstharmonic(seeFig.10).

The maximum velocity in the canyons near Vid Flumina is 0.064m/sandoccursatthemouthofthedendriticnetwork,where the canyon is narrow. Except for this location andat the mouth ofthe southern canyon, thespeed islower than 0.04m/s. In the canyonnearXanthusFlumen,thespeedislowerthan0.012m/s.

Inthelight ofsuch predictionsandaswedonottakeintoac- counttheslopeinthecanyon,itishighlyunlikely thattheliquid fillingthesecanyonsisduetothetides.Nonetheless,thoseconclu- sionsmaybeputintoquestioniftherearesignificantdiscrepancies betweenthebathymetryimplementedinthenumericalmodeland therealbathymetry.

3.4.Moraysinusandthe“MagicIslands” phenomenon

Inthissection,wewilldetailthetidalresponseofthenorthern seas inregions wheretransient features were observed and ina baywestofTrevizeFretum,MoraySinus.The“MagicIslands” phe- nomenonwas observed both in Kraken andLigeia Maria. Never- theless,Hofgartneretal.(2016)suggestedthat, inLigeia, thephe- nomenonisindependentfromthediurnaltide.Consequently,itis notofaparticularinteresttostudyaccuratelytheseareas.

In Kraken Mare, the “Magic Island” was observed at (73°N, 55°E)betweenT91andT104flybys.ThetrueanomalyofTitandur- ingtheseflybys was69° and246° respectively,which correspond to a variation of 177°, almost half a Titan day. Consequently, di- urnalphenomenon such asthetides cannot beruled out. Inthis area,tidalcomponentsother thanthediurnaloneare notsignifi- cantand,hence,canbedisregarded.Theseasurfaceelevationvari- ationbetweenbothflybysisabout0.1m,whichistoosmalltoex- plainthephenomenon.The tidalcurrentsinthisarea areweaker than0.029m/s(themaximumoccursneartheshoresofthenearby island),whichissmallerthanthespeedforcapillarywaves(3cm wavelength)calculatedbyHayesetal.(2013)(0.11m/s).

InMoray Sinus,thebathymetryisamphitheater-shapedwitha maximumdepth(45m) atthe inlet.Themaintidalcomponentis diurnalandthe tidalrangerangesfrom0.17to 0.225m.The flow is reallyweak exceptin the narrow part atthe northern end of MoraySinuswhereitcanreach0.058m/s.Inbothareas,nocapil- larywavesareexpected,thetidalflowbeingtoosmall.Neverthe- less,withdifferentbathymetryassumptions,theflowspeedcould behigher.Inthisconnection,wenotethat thespeedforcapillary waves(3cm wavelength) calculated by Hayes etal. (2013) could bereached–thustidalcurrentscould conceivablygeneratesurface