Back to home page

MITgcm

 
 

    


File indexing completed on 2018-03-02 18:43:13 UTC

view on githubraw file Latest commit add29e06 on 2018-01-31 20:35:05 UTC
6060ec2938 Dimi*0001 #include "SBO_OPTIONS.h"
85d26bffb0 Jean*0002 #ifdef ALLOW_SEAICE
                0003 # include "SEAICE_OPTIONS.h"
                0004 #endif
6060ec2938 Dimi*0005 
ef080e1d37 Dimi*0006 CBOP
                0007 C !ROUTINE: SBO_CALC
                0008 
                0009 C !INTERFACE: ==========================================================
94c8eb5701 Jean*0010       SUBROUTINE SBO_CALC( myTime, myIter, myThid )
ef080e1d37 Dimi*0011 
                0012 C !DESCRIPTION: \bv
85d26bffb0 Jean*0013 C     *==========================================================*
6060ec2938 Dimi*0014 C     | SUBROUTINE SBO_CALC                                      |
                0015 C     | o Do SBO diagnostic output.                              |
85d26bffb0 Jean*0016 C     *==========================================================*
6060ec2938 Dimi*0017 C     | NOTE: The following subtleties are ignored for time      |
                0018 C     | being but may need revisiting at some point in time.     |
                0019 C     | 1) The model is volume-preserving and Boussinesq so      |
                0020 C     |    quantities like oceanic mass need to be interpreted   |
85d26bffb0 Jean*0021 C     |    with some care.  We remove spurious mass variations   |
7e75f24202 Jean*0022 C     |    using the Greatbatch correction. Real freshwater      |
                0023 C     |    fluxes retained in mass load.                         |
                0024 C     | 2) The sea surface height variable etaN might lag other  |
                0025 C     |    prognostic variables by half a time step.  This lag   |
                0026 C     |    is ignored in SBO computations.                       |
                0027 C     | 3) OAM due to currents assumes constant density          |
85d26bffb0 Jean*0028 C     |    (=rhoConst), rms differences using variable density   |
                0029 C     |    is less than 1%, assuming rhoConst is a good measure  |
                0030 C     |    of the actual mean density                            |
7e75f24202 Jean*0031 C     | 4) Seaice motion added to OAMC.  Seaice mass is in OAMP  |
                0032 C     |    and COM.   Net freshwater flux is between atmosphere  |
85d26bffb0 Jean*0033 C     |    and liquid ocean plus seaice.  I.e. changes in seaice |
                0034 C     |    mass due to melt/freeze with liquid ocean do not      |
                0035 C     |    change net freshwater flux.                           |
                0036 C     *==========================================================*
6060ec2938 Dimi*0037 
222b416016 Jean*0038 C=======================================================================
                0039 C
85d26bffb0 Jean*0040 C     Based on ftp://euler.jpl.nasa.gov/sbo/software/calc_sbo2.f
222b416016 Jean*0041 C     Written  by Richard Gross (Richard.Gross@jpl.nasa.gov)
85d26bffb0 Jean*0042 C     Reference
222b416016 Jean*0043 C           Gross, R. S., F. O. Bryan, Y. Chao, J. O. Dickey, S. L. Marcus,
                0044 C           R. M. Ponte, and R. Tokmakian, The IERS Special Bureau for the
                0045 C           Oceans, in IERS Technical Note on the IERS Global Geophysical
                0046 C           Fluids Center, edited by B. Chao, in press, Observatoire de Paris,
85d26bffb0 Jean*0047 C           Paris, France, 2002.
222b416016 Jean*0048 C
85d26bffb0 Jean*0049 C     June 10, 2001: Modified for online computations in MIT GCM UV
                0050 C              by Dimitris Menemenlis (Menemenlis@jpl.nasa.gov)
                0051 C     Jan 7, 2014: Modified for real freshwater flux and coordinates other
                0052 C              than spherical polar by Katy Quinn (kquinn@aer.com)
222b416016 Jean*0053 C
85d26bffb0 Jean*0054 C     Purpose
                0055 C           calc_sbo calculates the core products of the IERS Special Bureau
                0056 C           for the Oceans including oceanic mass, center-of-mass, and angular
                0057 C           momentum.
222b416016 Jean*0058 C
                0059 C=======================================================================
                0060 C \ev
6060ec2938 Dimi*0061 
ef080e1d37 Dimi*0062 C !USES: ===============================================================
85d26bffb0 Jean*0063       IMPLICIT NONE
6060ec2938 Dimi*0064 C     === Global variables ===
                0065 #include "SIZE.h"
                0066 #include "EEPARAMS.h"
                0067 #include "PARAMS.h"
                0068 #include "GRID.h"
                0069 #include "DYNVARS.h"
fac3e5bea0 Patr*0070 #include "FFIELDS.h"
6060ec2938 Dimi*0071 #include "SBO.h"
85d26bffb0 Jean*0072 #ifdef ALLOW_SEAICE
                0073 # include "SEAICE_SIZE.h"
                0074 # include "SEAICE.h"
                0075 #endif
6060ec2938 Dimi*0076 
ef080e1d37 Dimi*0077 C !INPUT PARAMETERS: ===================================================
6060ec2938 Dimi*0078 C     == Routine arguments ==
94c8eb5701 Jean*0079 C     myTime  :: Current time of simulation ( s )
                0080 C     myIter  :: Iteration number
                0081 C     myThid  :: Number of this instance of SBO_CALC
                0082       _RL     myTime
6060ec2938 Dimi*0083       INTEGER myIter, myThid
                0084 
                0085 #ifdef ALLOW_SBO
                0086 
ef080e1d37 Dimi*0087 C !LOCAL VARIABLES: ====================================================
222b416016 Jean*0088 C     bi, bj    :: array indices
                0089 C     i         :: index over longitude grid points
                0090 C     j         :: index over latitude  grid points
                0091 C     k         :: index over layers
                0092 C     lat       :: latitude  of grid point (radians)
                0093 C     lon       :: longitude of grid point (radians)
85d26bffb0 Jean*0094 C     darea     :: element of surface area (m**2)
222b416016 Jean*0095 C     dvolume   :: element of volume (m**3)
85d26bffb0 Jean*0096 C     ae        :: mean radius of Earth (m) (PREM value)
                0097 C     sbo_omega :: mean angular velocity of Earth (rad/s)
                0098 C     UE,VN     :: geographic (east,north) ocean velocities at cell centers (m/s)
                0099 C     UEice,VNice  :: geographic (east,north) seaice velocities at cell centers (m/s)
                0100 C     Mload     :: total mass load (kg/m**2)
                0101 C     GCload    :: mass load for Greatbatch correction (kg/m**2)
                0102 C     FWLoad    :: real freshwater flux mass load (kg/m**2)
                0103       integer bi, bj, i, j, k
                0104       _RL lat, lon, darea, dvolume
                0105       _RL ae, sbo_omega
9ef95b2e9f Ed H*0106       PARAMETER ( ae        = 6.3710 _d 6    )
                0107       PARAMETER ( sbo_omega = 7.292115 _d -5 )
85d26bffb0 Jean*0108       _RL UE(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
                0109       _RL VN(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
                0110       _RL UEice(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0111       _RL VNice(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0112       _RL Mload(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
7e75f24202 Jean*0113       _RL GCload, FWload
85d26bffb0 Jean*0114 C Tiled global sums
7e75f24202 Jean*0115       _RL tile_FWload(nSx,nSy)
85d26bffb0 Jean*0116       _RL tile_sboarea(nSx,nSy)
                0117       _RL tile_GCload(nSx,nSy)
                0118       _RL tile_mass(nSx,nSy)
                0119       _RL tile_xcom(nSx,nSy)
                0120       _RL tile_ycom(nSx,nSy)
                0121       _RL tile_zcom(nSx,nSy)
                0122       _RL tile_xoamc(nSx,nSy)
                0123       _RL tile_yoamc(nSx,nSy)
                0124       _RL tile_zoamc(nSx,nSy)
                0125       _RL tile_xoamp(nSx,nSy)
                0126       _RL tile_yoamp(nSx,nSy)
                0127       _RL tile_zoamp(nSx,nSy)
                0128       _RL tile_xoamc_si(nSx,nSy)
                0129       _RL tile_yoamc_si(nSx,nSy)
                0130       _RL tile_zoamc_si(nSx,nSy)
                0131       _RL tile_mass_si(nSx,nSy)
                0132       _RL tile_mass_fw(nSx,nSy)
                0133       _RL tile_xcom_fw(nSx,nSy)
                0134       _RL tile_ycom_fw(nSx,nSy)
                0135       _RL tile_zcom_fw(nSx,nSy)
                0136       _RL tile_xoamp_fw(nSx,nSy)
                0137       _RL tile_yoamp_fw(nSx,nSy)
                0138       _RL tile_zoamp_fw(nSx,nSy)
                0139       _RL tile_mass_gc (nSx,nSy)
                0140 C Pre-computed cos(lat), sin(lat), cos(lon), sin(lon)
                0141       _RL COSlat(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0142       _RL SINlat(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0143       _RL COSlon(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0144       _RL SINlon(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
ef080e1d37 Dimi*0145 CEOP
6060ec2938 Dimi*0146 
85d26bffb0 Jean*0147 C Initialize variables to be computed---------------------------------
55a4640a63 Jean*0148 C-    note: only done once (by master thread) for var in common block
                0149       _BEGIN_MASTER(myThid)
85d26bffb0 Jean*0150 
6060ec2938 Dimi*0151       xoamc = 0.0
                0152       yoamc = 0.0
                0153       zoamc = 0.0
                0154       xoamp = 0.0
                0155       yoamp = 0.0
                0156       zoamp = 0.0
                0157       mass  = 0.0
                0158       xcom  = 0.0
                0159       ycom  = 0.0
                0160       zcom  = 0.0
55a4640a63 Jean*0161       sboarea = 0.0
94c8eb5701 Jean*0162 
85d26bffb0 Jean*0163       xoamc_si = 0.0
                0164       yoamc_si = 0.0
                0165       zoamc_si = 0.0
                0166       mass_si  = 0.0
                0167 
                0168       xoamp_fw = 0.0
                0169       yoamp_fw = 0.0
                0170       zoamp_fw = 0.0
                0171       mass_fw  = 0.0
                0172       xcom_fw  = 0.0
                0173       ycom_fw  = 0.0
                0174       zcom_fw  = 0.0
                0175 
                0176       mass_gc  = 0.0
                0177 
55a4640a63 Jean*0178       _END_MASTER(myThid)
                0179 
85d26bffb0 Jean*0180 C Get geographic (East,North) velocities------------------------------
                0181 
                0182       CALL ROTATE_UV2EN_RL(
                0183      U          uVel, vVel,
                0184      U          UE, VN,
                0185      I          .TRUE., .TRUE., .FALSE., Nr, mythid )
                0186 
                0187 #ifdef ALLOW_SEAICE
                0188       IF ( useSEAICE ) THEN
                0189         CALL ROTATE_UV2EN_RL(
                0190      U          UICE, VICE,
                0191      U          UEice, VNice,
                0192      I          .TRUE., .TRUE., .FALSE., 1, mythid )
                0193       ELSE
                0194 #else /* ALLOW_SEAICE */
                0195       IF ( .TRUE. ) THEN
                0196 #endif /* ALLOW_SEAICE */
                0197         DO bj = myByLo(myThid), myByHi(myThid)
                0198          DO bi = myBxLo(myThid), myBxHi(myThid)
                0199           DO j=1-OLy,sNy+OLy
                0200            DO i=1-OLx,sNx+OLx
                0201              UEice(i,j,bi,bj) = 0.
                0202              VNice(i,j,bi,bj) = 0.
                0203            ENDDO
                0204           ENDDO
                0205          ENDDO
                0206         ENDDO
                0207       ENDIF
                0208 
                0209 C Calculate mass load-------------------------------------------------
                0210 C     Calculate mass load (Mload), Greatbatch correction for
                0211 C     spurious mass but spatial mean freshwater flux retained.
                0212 C     Mload *needs* to be total mass (for center of mass), so add
                0213 C     back missing time invariant term: -R_low*rhoConst
                0214 
                0215 c Calculate freshwater load
                0216 C Calculate Greatbatch correction load over global ocean volume
                0217 c Note: no halo regions in i,j loops, do not want to double book sums
7e75f24202 Jean*0218       FWload = 0.0
85d26bffb0 Jean*0219       GCload = 0.0
6060ec2938 Dimi*0220       DO bj = myByLo(myThid), myByHi(myThid)
222b416016 Jean*0221        DO bi = myBxLo(myThid), myBxHi(myThid)
7e75f24202 Jean*0222         tile_FWload(bi,bj) = 0.0
85d26bffb0 Jean*0223         tile_GCload(bi,bj) = 0.0
                0224         tile_sboarea(bi,bj) = 0.0
                0225         DO j = 1, sNy
                0226          DO i = 1, sNx
                0227             darea = rA(i,j,bi,bj)*maskC(i,j,1,bi,bj)
                0228             tile_sboarea(bi,bj) = tile_sboarea(bi,bj) + darea
7e75f24202 Jean*0229             tile_FWload(bi,bj) =  tile_FWload(bi,bj) +
                0230      &           rhoConst*etaN(i,j,bi,bj)*darea +
                0231      &           sIceLoad(i,j,bi,bj)*darea
85d26bffb0 Jean*0232             DO k = 1, Nr
7e75f24202 Jean*0233                dvolume = rA(i,j,bi,bj)*drF(k)*hFacC(i,j,k,bi,bj)
85d26bffb0 Jean*0234                tile_GCload(bi,bj) = tile_GCload(bi,bj) +
                0235      &              rhoInSitu(i,j,k,bi,bj) * dvolume
                0236             ENDDO
6060ec2938 Dimi*0237          ENDDO
85d26bffb0 Jean*0238         ENDDO
                0239        ENDDO
                0240       ENDDO
7e75f24202 Jean*0241       CALL GLOBAL_SUM_TILE_RL( tile_FWload  , FWload  , myThid )
85d26bffb0 Jean*0242       CALL GLOBAL_SUM_TILE_RL( tile_sboarea , sboarea , myThid )
                0243       CALL GLOBAL_SUM_TILE_RL( tile_GCload  , GCload  , myThid )
7e75f24202 Jean*0244       FWload = FWload/sboarea
85d26bffb0 Jean*0245       GCload = -1.0 * GCload/sboarea
                0246 
                0247 c Total mass load with freshwater flux
                0248       DO bj = myByLo(myThid), myByHi(myThid)
                0249        DO bi = myBxLo(myThid), myBxHi(myThid)
                0250         DO j = 1-OLy, sNy+OLy
                0251          DO i = 1-OLx, sNx+OLx
7e75f24202 Jean*0252             Mload(i,j,bi,bj) =  rhoConst*etaN(i,j,bi,bj) +
                0253      &           sIceLoad(i,j,bi,bj) +
                0254      &           GCload - R_low(i,j,bi,bj)*rhoConst
                0255             DO k = 1, Nr
                0256                Mload(i,j,bi,bj) = Mload(i,j,bi,bj) +
                0257      &           rhoInSitu(i,j,k,bi,bj)*drF(k)*hFacC(i,j,k,bi,bj)
                0258             ENDDO
85d26bffb0 Jean*0259          ENDDO
                0260         ENDDO
222b416016 Jean*0261        ENDDO
6060ec2938 Dimi*0262       ENDDO
                0263 
85d26bffb0 Jean*0264 c Pre-compute cos(lat), sin(lat), cos(lon), sin(lon)
222b416016 Jean*0265       DO bj = myByLo(myThid), myByHi(myThid)
                0266        DO bi = myBxLo(myThid), myBxHi(myThid)
85d26bffb0 Jean*0267         DO j = 1-OLy, sNy+OLy
                0268          DO i = 1-OLx, sNx+OLx
                0269             lat = yC(i,j,bi,bj) * deg2rad
                0270             lon = xC(i,j,bi,bj) * deg2rad
                0271             COSlat(i,j,bi,bj) = COS(lat)
                0272             SINlat(i,j,bi,bj) = SIN(lat)
                0273             COSlon(i,j,bi,bj) = COS(lon)
                0274             SINlon(i,j,bi,bj) = SIN(lon)
                0275          ENDDO
                0276         ENDDO
                0277        ENDDO
                0278       ENDDO
222b416016 Jean*0279 
85d26bffb0 Jean*0280 C Main loops----------------------------------------------------------
                0281 C     loop over all grid points, accumulating mass, com, oam
                0282 C     Note: no halo regions in i,j loops, do not want to double book sums
222b416016 Jean*0283 
85d26bffb0 Jean*0284       DO bj = myByLo(myThid), myByHi(myThid)
                0285        DO bi = myBxLo(myThid), myBxHi(myThid)
                0286 
                0287 C Initialize tile sums---------------------------------
                0288         tile_xoamc(bi,bj) = 0.0
                0289         tile_yoamc(bi,bj) = 0.0
                0290         tile_zoamc(bi,bj) = 0.0
                0291         tile_xoamp(bi,bj) = 0.0
                0292         tile_yoamp(bi,bj) = 0.0
                0293         tile_zoamp(bi,bj) = 0.0
                0294         tile_mass(bi,bj) = 0.0
                0295         tile_xcom(bi,bj) = 0.0
                0296         tile_ycom(bi,bj) = 0.0
                0297         tile_zcom(bi,bj) = 0.0
                0298 
                0299         tile_xoamc_si(bi,bj) = 0.0
                0300         tile_yoamc_si(bi,bj) = 0.0
                0301         tile_zoamc_si(bi,bj) = 0.0
                0302         tile_mass_si(bi,bj) = 0.0
                0303 
                0304         tile_xoamp_fw(bi,bj) = 0.0
                0305         tile_yoamp_fw(bi,bj) = 0.0
                0306         tile_zoamp_fw(bi,bj) = 0.0
                0307         tile_mass_fw(bi,bj) = 0.0
                0308         tile_xcom_fw(bi,bj) = 0.0
                0309         tile_ycom_fw(bi,bj) = 0.0
                0310         tile_zcom_fw(bi,bj) = 0.0
                0311 
                0312         tile_mass_gc(bi,bj) = 0.0
222b416016 Jean*0313 
                0314         DO j = 1, sNy
85d26bffb0 Jean*0315          DO i = 1, sNx
                0316 
                0317            IF ( maskC(i,j,1,bi,bj) .NE. 0. ) THEN
                0318 
                0319 C     horizontal area
                0320               darea = rA(i,j,bi,bj)*maskC(i,j,1,bi,bj)
                0321 
                0322 C     accumulate mass of oceans, Greatbatch correction, seaice
                0323               tile_mass(bi,bj) = tile_mass(bi,bj) +
                0324      &          Mload(i,j,bi,bj)*darea
                0325               tile_mass_gc(bi,bj) = tile_mass_gc(bi,bj) +
                0326      &          GCload*darea
                0327               tile_mass_si(bi,bj) = tile_mass_si(bi,bj) +
                0328      &          sIceLoad(i,j,bi,bj)*darea
                0329 
                0330 C     accumulate center-of-mass of oceans (need to divide by total mass at end)
                0331               tile_xcom(bi,bj) = tile_xcom(bi,bj) +
                0332      &          Mload(i,j,bi,bj)*COSlat(i,j,bi,bj)*COSlon(i,j,bi,bj)
                0333      &          * ae * darea
                0334               tile_ycom(bi,bj) = tile_ycom(bi,bj) +
                0335      &          Mload(i,j,bi,bj)*COSlat(i,j,bi,bj)*SINlon(i,j,bi,bj)
                0336      &          * ae * darea
                0337               tile_zcom(bi,bj) = tile_zcom(bi,bj) +
                0338      &          Mload(i,j,bi,bj)*SINlat(i,j,bi,bj)
                0339      &          * ae * darea
                0340 
                0341 C     accumulate oceanic angular momentum due to currents (need depth integral too)
                0342 C     Note: depth integral goes from k=1,Nr.  hFacC takes care of R_low and etaN (as per JMC)
                0343               DO k = 1, Nr
                0344                  dvolume = rA(i,j,bi,bj)*drF(k)
                0345      &             * maskC(i,j,k,bi,bj)*hFacC(i,j,k,bi,bj)
                0346                  tile_xoamc(bi,bj) = tile_xoamc(bi,bj) +
                0347      &             ( VN(i,j,k,bi,bj)*SINlon(i,j,bi,bj) -
                0348      &             UE(i,j,k,bi,bj)*
                0349      &             SINlat(i,j,bi,bj)*COSlon(i,j,bi,bj) )
                0350      &             * rhoConst * ae * dvolume
                0351                  tile_yoamc(bi,bj) = tile_yoamc(bi,bj) +
                0352      &             (-VN(i,j,k,bi,bj)*COSlon(i,j,bi,bj) -
                0353      &             UE(i,j,k,bi,bj)*
                0354      &             SINlat(i,j,bi,bj)*SINlon(i,j,bi,bj) )
                0355      &             * rhoConst * ae * dvolume
                0356                  tile_zoamc(bi,bj) = tile_zoamc(bi,bj) +
                0357      &             UE(i,j,k,bi,bj)*COSlat(i,j,bi,bj)
                0358      &             * rhoConst * ae * dvolume
                0359               ENDDO
                0360 
                0361 C     accumulate sea angular momentum due to motion (one layer, so no depth integral needed)
                0362               tile_xoamc_si(bi,bj) = tile_xoamc_si(bi,bj) +
                0363      &          ( VNice(i,j,bi,bj)*SINlon(i,j,bi,bj) -
                0364      &          UEice(i,j,bi,bj)*
                0365      &          SINlat(i,j,bi,bj)*COSlon(i,j,bi,bj) )
                0366      &          * sIceLoad(i,j,bi,bj) * ae * darea
                0367               tile_yoamc_si(bi,bj) = tile_yoamc_si(bi,bj) +
                0368      &          (-VNice(i,j,bi,bj)*COSlon(i,j,bi,bj) -
                0369      &          UEice(i,j,bi,bj)*
                0370      &          SINlat(i,j,bi,bj)*SINlon(i,j,bi,bj) )
                0371      &          * sIceLoad(i,j,bi,bj) * ae * darea
                0372               tile_zoamc_si(bi,bj) = tile_zoamc_si(bi,bj) +
                0373      &          UEice(i,j,bi,bj)*COSlat(i,j,bi,bj)
                0374      &          * sIceLoad(i,j,bi,bj) * ae * darea
                0375 
                0376 C     accumulate oceanic angular momentum due to pressure
                0377               tile_xoamp(bi,bj) = tile_xoamp(bi,bj) -
                0378      &          SINlat(i,j,bi,bj)*COSlat(i,j,bi,bj)*COSlon(i,j,bi,bj)
                0379      &          * sbo_omega * Mload(i,j,bi,bj) * ae*ae * darea
                0380               tile_yoamp(bi,bj) = tile_yoamp(bi,bj) -
                0381      &          SINlat(i,j,bi,bj)*COSlat(i,j,bi,bj)*SINlon(i,j,bi,bj)
                0382      &          * sbo_omega * Mload(i,j,bi,bj) * ae*ae * darea
                0383               tile_zoamp(bi,bj) = tile_zoamp(bi,bj) +
                0384      &          COSlat(i,j,bi,bj) * COSlat(i,j,bi,bj)
                0385      &          * sbo_omega * Mload(i,j,bi,bj) * ae*ae * darea
                0386 
                0387 C     accumulate mass of real freshwater flux
                0388               tile_mass_fw(bi,bj) = tile_mass_fw(bi,bj) +
                0389      &          FWload * darea
                0390 
                0391 C     accumulate center-of-mass of real freshwater flux (need to divide by total FW mass at end)
                0392               tile_xcom_fw(bi,bj) = tile_xcom_fw(bi,bj) +
                0393      &          FWload * COSlat(i,j,bi,bj) * COSlon(i,j,bi,bj)
                0394      &          * ae * darea
                0395               tile_ycom_fw(bi,bj) = tile_ycom_fw(bi,bj) +
                0396      &          FWload * COSlat(i,j,bi,bj) * SINlon(i,j,bi,bj)
                0397      &          * ae * darea
                0398               tile_zcom_fw(bi,bj) = tile_zcom_fw(bi,bj) +
                0399      &          FWload * SINlat(i,j,bi,bj)
                0400      &          * ae * darea
                0401 
                0402 C     accumulate oceanic angular momentum due to real freshwater flux
                0403               tile_xoamp_fw(bi,bj) = tile_xoamp_fw(bi,bj) -
                0404      &          SINlat(i,j,bi,bj)*COSlat(i,j,bi,bj)*COSlon(i,j,bi,bj)
                0405      &          * sbo_omega * FWload * ae*ae * darea
                0406               tile_yoamp_fw(bi,bj) = tile_yoamp_fw(bi,bj) -
                0407      &          SINlat(i,j,bi,bj)*COSlat(i,j,bi,bj)*SINlon(i,j,bi,bj)
                0408      &          * sbo_omega * FWload * ae*ae * darea
                0409               tile_zoamp_fw(bi,bj) = tile_zoamp_fw(bi,bj) +
                0410      &          COSlat(i,j,bi,bj) * COSlat(i,j,bi,bj)
                0411      &          * sbo_omega * FWload * ae*ae * darea
                0412 
                0413 C     end if over ocean
                0414            ENDIF
                0415 C     end loop over i,j
                0416          ENDDO
222b416016 Jean*0417         ENDDO
                0418 
                0419 C     end loop over bi,bj
                0420        ENDDO
                0421       ENDDO
                0422 
85d26bffb0 Jean*0423 C     sum all global values across model tiles
                0424       CALL GLOBAL_SUM_TILE_RL( tile_mass  , mass  , myThid )
                0425       CALL GLOBAL_SUM_TILE_RL( tile_xcom  , xcom  , myThid )
                0426       CALL GLOBAL_SUM_TILE_RL( tile_ycom  , ycom  , myThid )
                0427       CALL GLOBAL_SUM_TILE_RL( tile_zcom  , zcom  , myThid )
                0428       CALL GLOBAL_SUM_TILE_RL( tile_xoamc , xoamc , myThid )
                0429       CALL GLOBAL_SUM_TILE_RL( tile_yoamc , yoamc , myThid )
                0430       CALL GLOBAL_SUM_TILE_RL( tile_zoamc , zoamc , myThid )
                0431       CALL GLOBAL_SUM_TILE_RL( tile_xoamp , xoamp , myThid )
                0432       CALL GLOBAL_SUM_TILE_RL( tile_yoamp , yoamp , myThid )
                0433       CALL GLOBAL_SUM_TILE_RL( tile_zoamp , zoamp , myThid )
                0434 
                0435       CALL GLOBAL_SUM_TILE_RL( tile_xoamc_si , xoamc_si , myThid )
                0436       CALL GLOBAL_SUM_TILE_RL( tile_yoamc_si , yoamc_si , myThid )
                0437       CALL GLOBAL_SUM_TILE_RL( tile_zoamc_si , zoamc_si , myThid )
                0438       CALL GLOBAL_SUM_TILE_RL( tile_mass_si , mass_si , myThid )
                0439 
                0440       CALL GLOBAL_SUM_TILE_RL( tile_mass_fw  , mass_fw  , myThid )
                0441       CALL GLOBAL_SUM_TILE_RL( tile_xcom_fw  , xcom_fw  , myThid )
                0442       CALL GLOBAL_SUM_TILE_RL( tile_ycom_fw  , ycom_fw  , myThid )
                0443       CALL GLOBAL_SUM_TILE_RL( tile_zcom_fw  , zcom_fw  , myThid )
                0444       CALL GLOBAL_SUM_TILE_RL( tile_xoamp_fw , xoamp_fw , myThid )
                0445       CALL GLOBAL_SUM_TILE_RL( tile_yoamp_fw , yoamp_fw , myThid )
                0446       CALL GLOBAL_SUM_TILE_RL( tile_zoamp_fw , zoamp_fw , myThid )
                0447       CALL GLOBAL_SUM_TILE_RL( tile_mass_gc  , mass_gc  , myThid )
6060ec2938 Dimi*0448 
222b416016 Jean*0449 C     finish calculating center-of-mass of oceans
55a4640a63 Jean*0450 C-    note: only master thread updates/modifies var in common block
                0451       _BEGIN_MASTER(myThid)
                0452 
85d26bffb0 Jean*0453       IF ( mass.NE.zeroRL ) THEN
                0454         xcom = xcom / mass
                0455         ycom = ycom / mass
                0456         zcom = zcom / mass
                0457       ENDIF
                0458 
                0459       IF ( mass_fw.NE.zeroRL ) THEN
                0460         xcom_fw = xcom_fw / mass_fw
                0461         ycom_fw = ycom_fw / mass_fw
                0462         zcom_fw = zcom_fw / mass_fw
                0463       ENDIF
                0464 
                0465 C     Add seaice OAMC to total OAMC
                0466       xoamc = xoamc + xoamc_si
                0467       yoamc = yoamc + yoamc_si
                0468       zoamc = zoamc + zoamc_si
6060ec2938 Dimi*0469 
55a4640a63 Jean*0470       _END_MASTER(myThid)
                0471 
f6bfe3bad8 Alis*0472 #endif /* ALLOW_SBO */
6060ec2938 Dimi*0473 
222b416016 Jean*0474       RETURN
                0475       END