JournalofMarineSystems65(2007)376–399 www.elsevier.com/locate/jmarsys An inverse model for calculation of global volume transport from wind and hydrographic data ⁎ Peter C. Chu , Chenwu Fan NavalOceanAnalysisandPredictionLaboratory,DepartmentofOceanography,NavalPostgraduateSchool,Monterey,CA93943,USA Received24September2004;accepted9June2005 Availableonline28November2006 Abstract TheP-vectorinversemethodhasbeensuccessfullyusedtoinverttheabsolutevelocityfromhydrographicdatafortheextra- equatorial hemispheres, but not for the equatorial region since it is based on the geostrophic balance. A smooth interpolation schemeacrosstheequatorisdevelopedinthisstudytobringtogetherthetwoalreadyknownsolutions(P-vectors)fortheextra- equatorialhemispheres.Thismodelcontainsfourmajorcomponents:(a)theP-vectorinversemodeltoobtainthesolutionsforthe two extra-equatorial hemispheres, (b) the objective method to determine the Ψ-values at individual islands, (c) the Poisson equation-solvertoobtaintheΠ-valuesovertheequatorialregionfromthevolumetransportvorticityequation,and(d)thePoisson equation-solver to obtain the Ψ and depth-integrated velocity field (U, V) over the globe from the Poisson Ψ-equation. The Poissonequation-solverissimilartotheboxmodeldevelopedbyWunsch.Thus,thismethodcombinesthestrengthfrombothbox andP-vector models.Thecalculated depth-integrated velocityand Ψ-field agreewell with earlier studies. Published byElsevier B.V. Keywords:Inversemodel;Volumetransportstreamfunction;Volumetransportvorticity;Wind-drivencirculation;Density-drivencirculation;Depth- integratedvelocity;Poissonequation;P-vector 1. Introduction Sverdrup relation in Antarctic regions. Godfrey (1989) used the Sverdrup model with climatological annual Winds, density, and friction determine volume trans- winds(HellermanandRosenstein1983)tocalculatethe port in oceans. The wind-driven volume transport has meandepth-integratedstreamfunctionfortheworldocean beenestimatedusingtheSverdrup(1947)relationforthe under two assumptions: (1) the ocean is stagnant below interiorocean,andusingtheStommel(1948)andMunk some depth, and (2) all major undersea topographic (1950) linear frictional ocean model for the intensive featuressuchasmid-oceanridgesliebelowthatdepth. westernboundarycurrents.LeetmaaandBunker(1978), The density-driven volume transport has been calcu- Meyers (1980), andGodfrey and Golding (1981) latedbyseveralauthors.Thedensityfielddirectlydeter- calculated similarly for the North Atlantic, tropical minesthegeostrophicvelocityrelativetothebottomflow. Pacific, and Indian Oceans. Baker (1982) examined the Thebottomvelocity(u−H,v−H)isusuallycalculatedusing the β-spiral (Stommeland Schott,1977),Box (Wunsch, ⁎ 1978),andP-vector(Chu,1995;Chuetal.,1998a,b;Chu, Correspondingauthor. E-mailaddress:[email protected](P.C.Chu). 2000;Chuetal.,2001a,b;Chu,2006)models. 0924-7963/$-seefrontmatter.PublishedbyElsevierB.V. doi:10.1016/j.jmarsys.2005.06.010 Report Documentation Page Form Approved OMB No. 0704-0188 Public reporting burden for the collection of information is estimated to average 1 hour per response, including the time for reviewing instructions, searching existing data sources, gathering and maintaining the data needed, and completing and reviewing the collection of information. Send comments regarding this burden estimate or any other aspect of this collection of information, including suggestions for reducing this burden, to Washington Headquarters Services, Directorate for Information Operations and Reports, 1215 Jefferson Davis Highway, Suite 1204, Arlington VA 22202-4302. Respondents should be aware that notwithstanding any other provision of law, no person shall be subject to a penalty for failing to comply with a collection of information if it does not display a currently valid OMB control number. 1. REPORT DATE 3. DATES COVERED SEP 2004 2. REPORT TYPE 00-00-2004 to 00-00-2004 4. TITLE AND SUBTITLE 5a. CONTRACT NUMBER An inverse model for calculation of global volume transport from wind 5b. GRANT NUMBER and hydrographic data 5c. PROGRAM ELEMENT NUMBER 6. AUTHOR(S) 5d. PROJECT NUMBER 5e. TASK NUMBER 5f. WORK UNIT NUMBER 7. PERFORMING ORGANIZATION NAME(S) AND ADDRESS(ES) 8. PERFORMING ORGANIZATION Naval Postgraduate School,Department of REPORT NUMBER Oceanography,Monterey,CA,93943 9. SPONSORING/MONITORING AGENCY NAME(S) AND ADDRESS(ES) 10. SPONSOR/MONITOR’S ACRONYM(S) 11. SPONSOR/MONITOR’S REPORT NUMBER(S) 12. DISTRIBUTION/AVAILABILITY STATEMENT Approved for public release; distribution unlimited 13. SUPPLEMENTARY NOTES 14. ABSTRACT The P-vector inverse method has been successfully used to invert the absolute velocity from hydrographic data for the extraequatorial hemispheres, but not for the equatorial region since it is based on the geostrophic balance. A smooth interpolation scheme across the equator is developed in this study to bring together the two already known solutions (P-vectors) for the extraequatorial hemispheres. This model contains four major components: (a) the P-vector inverse model to obtain the solutions for the two extra-equatorial hemispheres, (b) the objective method to determine the Ψ-values at individual islands, (c) the Poisson equation-solver to obtain the Π-values over the equatorial region from the volume transport vorticity equation, and (d) the Poisson equation-solver to obtain the Ψ and depth-integrated velocity field (U, V) over the globe from the Poisson Ψ-equation. The Poisson equation-solver is similar to the box model developed by Wunsch. Thus, this method combines the strength from both box and P-vector models. The calculated depth-integrated velocity and Ψ-field agree well with earlier studies. 15. SUBJECT TERMS 16. SECURITY CLASSIFICATION OF: 17. LIMITATION OF 18. NUMBER 19a. NAME OF ABSTRACT OF PAGES RESPONSIBLE PERSON a. REPORT b. ABSTRACT c. THIS PAGE Same as 24 unclassified unclassified unclassified Report (SAR) Standard Form 298 (Rev. 8-98) Prescribed by ANSI Std Z39-18 P.C.Chu,C.Fan/JournalofMarineSystems65(2007)376–399 377 Since the geostrophic balance is used, these models The horizontal diffusivity A can be estimated by h provide the solutions for the extra-equatorial hemi- Smargrinsky parameterization, spheres,butnotfortheequatorialregion.Animproved D inverse method is developed in this study to bring A ¼ DxDyjjVþðjVÞTj; ð6Þ h 2 togetherthetwoknownsolutionsfortheextra-equatorial hemispheresacrosstheequatorandtoestablishaglobal where velocity dataset. The rest of the paper is outlined as " # (cid:1) (cid:3) (cid:1) (cid:3) (cid:1) (cid:3) 1=2 follows. Section 2 describes the basic theory of the Au 2 1 Av Au 2 Av 2 model.Sections3–5describedepth-integratedvelocity, jjVþðjVÞTju Ax þ2 AxþAy þ Ay : volume transport streamfunction, and volume transport vorticity. Section 6 depicts theΨ-Poisson equation and Here, the nondimensional parameter D varies from itssolver.Section7depictsthemodelsensitivity.Section 0.1 to 0.2 (Mellor, 2003); we set D=0.15. The 8providestheglobalcirculationcharacteristics.Section horizontal grid in this study is 1°×1°, i.e., (Δx,Δy) 9presentstheconclusions. ∼100 km. Let the spatial variability of the velocity be scaledby 0.1 ms−1, we have 2. Dynamics 0:1 m s−1 2.1. Basic equations jjVþðjVÞTjf2(cid:2) ¼2(cid:2)10−6s−1: ð7Þ 105m Let (x, y, z) be the coordinates with x-axis in the Substitutionof (7)into (6) leads to zonal direction (eastward positive), y-axis in the latitudinal direction (northward positive), and z-axis in A ¼1:5(cid:2)103m2s−1: ð8Þ thevertical(upwardpositive).Theunitvectorsalongthe h three axes are represented by (i, j, k). For the extra- equatorialregion,thelinearsteadystatesystemwiththe 2.2. Ekman number Boussinesq approximationisgiven by The Ekman number can identify the relative A2u −fðv−v Þ¼A þA j2u; ð1Þ importance of the horizontal gradient of the Reynolds g zAz2 h stress versus the Coriolis force, fðu−ugÞ¼AzAA2zv2þAhj2v; ð2Þ E ¼OðOjAðjhfjV2jVÞjÞ¼jfAjhL2; ð9Þ Ap ¼−qg; ð3Þ where L is the characteristic horizontal length scale. In Az this study, the motion with L larger than 200 km is considered. For extra-equatorial regions (north of 8°N Au Av Aw and south of 8°S), þ þ ¼0; ð4Þ Ax Ay Az jfjN0:2(cid:2)10−4s−1; whereρisthein-situdensity;f=2Ωsinφ,istheCoriolis parameter,ΩtheEarthrotationrate,andφthelatitude. and theEkman number isestimatedby V=(u, v), is the horizontal velocity; w is the vertical velocity; ∇=i∂/∂x+j∂/∂y, is the horizontal gradient Eb 1:5(cid:2)103m2s−1 ¼1:875(cid:2)10−3; operator; V =(u , v ), is the geostrophic velocity ð0:2(cid:2)10−4s−1Þ(cid:2)ð2(cid:2)105mÞ2 g g g representing thehorizontal pressure(p) gradients whichshowsthatthehorizontalgradientoftheReynolds 1 Ap 1 Ap stresscanbeneglectedagainsttheCoriolisforce,i.e., u ¼− ; V ¼ ; ð5Þ g f q0Ay g f q0Ax Ah ¼0; ð10aÞ whereρ0isthecharacteristicvalue(1025kg/m3)ofthe forextra-equatorialregions. sea water density. The two coefficients (Az, Ah) are the Fortheequatorialregionsespeciallyneartheequator, vertical and horizontal eddy diffusivities. |f| is very small. The Ekman number is not a small 378 P.C.Chu,C.Fan/JournalofMarineSystems65(2007)376–399 Fig.1.Bathymetryofworldoceans. parameter.Thehorizontalviscousforce,(A ∇2u,A ∇2v), where z=−H(x, y) represents the ocean bottom, and h h cannot be neglected against the Coriolis force in the z=0referstotheoceansurface.Depth-integrationof(1) equatorialregion,thatis, and(2)fromtheoceanbottomtotheoceansurfaceleads to A p0; ð10bÞ h Au Au −fðV−VgÞ¼AzAzjz¼g−(cid:3)AzAzjz¼−H þAhj2U ð13Þ fortheequatorialregions. −2Ahju−H djH−Ahu−Hj2H; 3. Depth-integrated velocity Av Av Let (U, V) and (U , V ) be the depth-integrated fðU−UgÞ¼AzAzjz¼g−(cid:3)AzAzjz¼−H þAhj2V ð14Þ horizontal velocity g g −2Ahjv−H djH−Ahv−Hj2H; Z 0 where (u−H, v−H) are velocity components at the ocean ðU;VÞ¼ ðu;vÞdz; ð11Þ bottom. −H Theturbulentmomentumfluxattheoceansurfaceis calculated by and geostrophic velocity, Z (cid:1) (cid:3) 0 Au Av ðs ;s Þ ðUg;VgÞ¼ ðug;vgÞdz; ð12Þ Az Az;Az jz¼g ¼ xq y ; ð15Þ −H 0 P.C.Chu,C.Fan/JournalofMarineSystems65(2007)376–399 379 Fig.2.BoundaryconditionsofΨfortheglobalocean. where (τ ,τ ) are the surface wind stress components. is the density-driven transport. Re-arranging (13) and x y The turbulent momentum flux at the ocean bottom is (14), we have parameterized by (cid:1) (cid:3) s Az AAuz;qAAvzffiffiffiffiffijffizffi¼ffiffiffi−ffiffiHffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Ahj2U þfV ¼fVdenþfVb−qx0þAhQ1; ð21Þ ¼CD u2−H þv2−Hðu−H;v−HÞ; ð16Þ −Ahj2V þfU ¼fUdenþfUbþqsy (cid:4)AhQ2; ð22Þ 0 where C =0.0025 isthedrag coefficient. D The thermal wind relation can be obtained from where vertical integration of the hydrostatic balance Eq. (3) (cid:3) fromthebottom(−H)toanydepth(z)andthentheuse of thegeostrophic Eq. (5) Q1uð2ju−H djH þu−Hj2HÞ; Z (cid:3) g z Aq ug ¼u−H þfq0 −H AydzV; ð17Þ Q2uð2jv−H djH þv−Hj2HÞ; Z g z Aq and vg ¼v−H−fq0 −H AxdzV: ð18Þ (cid:1) qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi(cid:3) C Substitutionof(17)and(18)intothesecondequation Ub ¼ H− fD u2−H þv2−H u−H; ð23Þ of (12) leads to (cid:1) qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi(cid:3) ðUg;VgÞ¼ðUdenþHu−H;VdenþHv−HÞ; ð19Þ Vb ¼ H þCfD u2−H þv2−H v−H; where ðU ;V Þ are the transport due to the bottom currents, or simply den de(cid:1)nZ Z Z Z (cid:3) g 0 z Aq 0 z Aq called the bottom transport. With the known bottom ¼fq0 −H −H AydzVdz; − −H −H AxdzVdz ; (vUel,oVci)tycavnecbtoerd(uet−eHrm,vi−nHe)d,tfhroemdetphteh-winitnedg,radteednsviteyl,ocaintyd ð20Þ topographic data. 380 P.C.Chu,C.Fan/JournalofMarineSystems65(2007)376–399 Fig.3.ComputedΨ-valuesforeachcontinent/island:(a)annualmean,(b)January,and(c)July. P.C.Chu,C.Fan/JournalofMarineSystems65(2007)376–399 381 For the extra-equatorial regions, the horizontal dif- is the volume transport vorticity. Eq. (28) is called the fusion can be neglected [see (8)]. Eqs. (22) and (21) PoissonΨ-equation. become 5. Volume transport vorticity U⁎ ¼U þU þ sy ; ð24Þ den b fq 0 5.1. Volume transport vorticity equation V⁎ ¼V þV − sx : ð25Þ den b fq Summationofthedifferentiationof(21)withrespect 0 to y and the differentiation of (22) with respect to x Withtheknown(u−H,v−H),thedepth-integratedflow gives the volume transport vorticityequation, (U⁎,V⁎)maybedirectlycalculatedusing(24)and(25). (cid:1) (cid:3) ⁎ ⁎ b 1 As As However, the computed (U , V ) field using (24) and j2C ¼ ðV−V −V Þ− y− x (25) is quite noisy and cannot not be the final product. Ah(cid:1) den b(cid:3) Ahq0 Ax Ay Thus, the subscript ‘⁎’ is used to represent the interim þ AQ2−AQ1 ; ð30Þ depth-integratedvelocitycalculatedusing(24)and(25). Ax Ay 4. Volume transport streamfunction where β=df/dy, and (28) is used. Integrationofthecontinuityequationwithrespectto z from thebottomto thesurfaceyields, 5.2. Extra-equatorial region AU AH AV AH Ax þu−H Ax þ Ay þv−H Ay −w−H ¼0: ð26Þ For the extra-equatorial region, the horizontal diffusion can be neglected [see (8)]. Substitution of A h With the assumption that the water flows following into 0 leads to thebottom topography, (cid:5) (cid:6) AH AH Cu1 AðfVdenÞ−AðfUdenÞ w−H ¼u−H Ax þv−H Ay ; f (cid:5) Ax Ay (cid:6) (cid:5) (cid:1) (cid:3) (cid:1) (cid:3)(cid:6) 1 AðfV Þ AðfU Þ 1 A s A s Eq. (26) becomes þ b − b − x þ y : f Ax Ay f Ax q Ax q AU AV 0 0 þ ¼0; ð31Þ Ax Ay which leads to the definition of the volume transport Similarly, (30) becomes streamfunctionΨ, (cid:1) (cid:3) 1 As As U ¼−AAWy ; V ¼AAWx : ð27Þ bðV−Vden−VbÞ¼q0 Axy− Ayx ; ð32Þ which isthe Sverdrup relation. Subtractionofthedifferentiationof(22)withrespect In (31), (U , V ) depend on ρ only; (τ ,τ ) are to y from the differentiation of (21) with respect to x den den x y windstresscomponents;and(U ,V )aredeterminedby gives b b the bottom current velocity (u−H, v−H). The P-vector j2W¼P; ð28Þ inverse method (Chu, 1995; Chu et al., 1998a,b; Chu, where 2000)isusedtodetermine(u−H,v−H)fromhydrographic (cid:5) (cid:6) (cid:5) (cid:6) data(seeAppendixA).Inthisstudy,theclimatological 1 AðfV Þ AðfU Þ 1 AðfV Þ AðfU Þ Cu den − den þ b − b : hydrographic data (Levitus et al., 1994) are used to f Ax Ay f Ax Ay compute (U , V ) [see (20)]. The climatological (cid:5) (cid:1) (cid:3) (cid:1) (cid:3)(cid:6) den den 1 A s A s surface wind stress (τ ,τ ) data are obtained from the − x þ y x y f Ax q Ax q ComprehensiveOcean-AtmosphereDataSet(COADS, (cid:1) 0 (cid:3) 0 A AQ AQ da Silva et al., 1994). The bottom topography is þ h 1þ 2 ; obtained from the Navy's Digital Bathymetry Data f Ax Ay Base 5-min (DBDB5) (Fig. 1). The volume transport ð29Þ vorticityΠ isquite noisy. 382 P.C.Chu,C.Fan/JournalofMarineSystems65(2007)376–399 Fig.4.Globalvolumetransportstreamfunction(Ψ)computedfromtheinversemodel:(a)annualmean,(b)January,and(c)July. P.C.Chu,C.Fan/JournalofMarineSystems65(2007)376–399 383 5.3. Equatorial region (between 8°S and 8°N) boundary values of the vorticity Eq. (30). Here, the forcing terms [righthand-side of (30)] are calculated Let the volume transport vorticity Π calculated with the assumptions that (1) f=f(8°N) north of using(31)at8°Nand8°Sasthenorthernandsouthern the equator, and f=f(8°S) south of the equator, and Fig.5.Globaldepth-integratedvelocity(U,V)vectorscomputedfromtheinversemodel:(a)annualmean,(b)January,and(c)July. 384 P.C.Chu,C.Fan/JournalofMarineSystems65(2007)376–399 (2) (U, V) are calculated by (24) and (25). With the subjectivelysetupinsomeearlierstudies.Forexample, given forcing terms and the northern and southern in calculating the geostrophic transport in the Pacific boundary conditions and the cyclic eastern and Ocean,Reid(1997)setupΨ-valuetobe0forAntarctic, western boundary conditions, the volume transport 135Sv(1Sv=106m3s−1)forAustralia,and130Svfor vorticity Eq. (30) can be solved in the equatorial America.Incalculatingthegeostrophictransportinthe regionbetween8°Nand8°S.ThecomputedΠ-fieldis SouthAtlanticOcean,Reid(1989)setupΨ-valuetobe quite smooth. 0 for Antarctic, 132 Sv for Africa, and 130 Sv for America. Such a treatment subjectively prescribes 6. Poisson Ψ-equation 130Sv throughtheDrakePassage and 132Svthrough area between Africa and Antarctica. 6.1. Boundary conditions AnobjectivemethoddepictedinAppendix-Bisused to determine Ψ-values at islands. Fig. 3 shows the Withthecomputedglobalvolumetransportvorticity distribution of Ψ-value for each continent/island (Π), the Poisson Ψ-equation (28) can be solved when computed from the annual, January, and July mean the boundary conditions are given. No flow over the hydrographicandwinddata.Takingtheannualmeanas Antarctic Continent leads to the southern boundary an example, we have: 0 for the American Continent, condition 157.30 Sv for Antarctica, −21.74 Sv for Australia, −27.17 Sv for Madagascar, and −21.74 Sv for New W¼C ; at the southern boundary y¼y : ð33Þ Guinea. 1 s 7. Model sensitivity Nohorizontalconvergenceofthe2-dimensionalflow (U, V) at the North Pole (treated as an island) leads to Withthegiven valuesattheboundariesandislands, thenorthern boundarycondition we solve the Poisson Ψ-equation (28) with climatolog- icalannualandmonthlyΠ-fieldsandobtainannualand W¼C ; at the northern boundary y¼y ; ð34Þ 2 n where C and C are the constants to be determined. 1 2 Thecyclicboundaryconditionisappliedtothewestern and the eastern boundaries (Fig. 2). We integrate ∂Ψ/ ∂y=−U⁎withrespecttoyalongthewestern(oreastern) boundary from the southern boundary (Ψ=0) to the northernboundary toobtain Z y ⁎ Wj ðyÞ¼− U ðx ;yVÞdyV: ð35Þ west west ys The cyclic boundary condition leads to Wj ðyÞ¼Wj ðyÞ: ð36Þ east west Thus, thenorthern boundaryconditionis given by Z yn ⁎ W¼− U ðx ;yÞdy¼C : ð37Þ west 2 ys 6.2. Ψ-Values at islands Before solving the Poisson Ψ-equation (28) with theboundaryconditions(33)–(35)and(37),weneedto Fig.6.MonthlyvolumetransportthroughtheDrakePassagewitha know the Ψ-values at all islands. These values were smallseasonalvariation.