New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
isfhdiv.F90 in NEMO/trunk/src/OCE/ISF – NEMO

source: NEMO/trunk/src/OCE/ISF/isfhdiv.F90 @ 13295

Last change on this file since 13295 was 13295, checked in by acc, 4 years ago

Replace do-loop macros in the trunk with alternative forms with greater flexibility for extra halo applications. This alters a lot of routines but does not change any behaviour or results. do_loop_substitute.h90 is greatly simplified by this change. SETTE results are identical to those with the previous revision
File size: 6.5 KB
Line 
1MODULE isfhdiv
2   !!======================================================================
3   !!                       ***  MODULE  isfhdiv  ***
4   !! ice shelf horizontal divergence module :  update the horizontal divergence
5   !!                   with the ice shelf melt and coupling correction
6   !!======================================================================
7   !! History :  4.0  !  2019-09  (P. Mathiot) Original code
8   !!----------------------------------------------------------------------
9
10   !!----------------------------------------------------------------------
11   !!   isf_hdiv    : update the horizontal divergence with the ice shelf
12   !!                 melt and coupling correction
13   !!----------------------------------------------------------------------
14
15   USE isf_oce                ! ice shelf
16
17   USE dom_oce                ! time and space domain
18   USE phycst , ONLY: r1_rho0 ! physical constant
19   USE in_out_manager         !
20
21   IMPLICIT NONE
22
23   PRIVATE
24
25   PUBLIC isf_hdiv
26   !! * Substitutions
27#  include "do_loop_substitute.h90"
28#  include "domzgr_substitute.h90"
29
30CONTAINS
31
32   SUBROUTINE isf_hdiv( kt, Kmm, phdiv )
33      !!----------------------------------------------------------------------
34      !!                  ***  SUBROUTINE isf_hdiv  ***
35      !!       
36      !! ** Purpose :   update the horizontal divergence with the ice shelf contribution
37      !!                (parametrisation, explicit, ice sheet coupling conservation
38      !!                 increment)
39      !!
40      !!----------------------------------------------------------------------
41      REAL(wp), DIMENSION(:,:,:), INTENT( inout ) ::   phdiv   ! horizontal divergence
42      !!----------------------------------------------------------------------
43      INTEGER, INTENT(in) :: kt
44      INTEGER, INTENT(in) :: Kmm      !  ocean time level index
45      !
46      IF ( ln_isf ) THEN
47         !
48         ! ice shelf cavity contribution
49         IF ( ln_isfcav_mlt ) CALL isf_hdiv_mlt(misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav, fwfisf_cav, fwfisf_cav_b, phdiv)
50         !
51         ! ice shelf parametrisation contribution
52         IF ( ln_isfpar_mlt ) CALL isf_hdiv_mlt(misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par, fwfisf_par, fwfisf_par_b, phdiv)
53         !
54         ! ice sheet coupling contribution
55         IF ( ln_isfcpl .AND. kt /= 0 ) THEN
56            !
57            ! Dynamical stability at start up after change in under ice shelf cavity geometry is achieve by correcting the divergence.
58            ! This is achieved by applying a volume flux in order to keep the horizontal divergence after remapping
59            ! the same as at the end of the latest time step. So correction need to be apply at nit000 (euler time step) and
60            ! half of it at nit000+1 (leap frog time step).
61            IF ( kt == nit000   ) CALL isf_hdiv_cpl(Kmm, risfcpl_vol       , phdiv)
62            IF ( kt == nit000+1 ) CALL isf_hdiv_cpl(Kmm, risfcpl_vol*0.5_wp, phdiv)
63            !
64            ! correct divergence every time step to remove any trend due to coupling
65            ! conservation option
66            IF ( ln_isfcpl_cons ) CALL isf_hdiv_cpl(Kmm, risfcpl_cons_vol, phdiv)
67            !
68         END IF
69         !
70      END IF
71      !
72   END SUBROUTINE isf_hdiv
73
74   SUBROUTINE isf_hdiv_mlt(ktop, kbot, phtbl, pfrac, pfwf, pfwf_b, phdiv)
75      !!----------------------------------------------------------------------
76      !!                  ***  SUBROUTINE sbc_isf_div  ***
77      !!       
78      !! ** Purpose :   update the horizontal divergence with the ice shelf inflow
79      !!
80      !! ** Method  :   pfwf is positive (outflow) and expressed as kg/m2/s
81      !!                increase the divergence
82      !!
83      !! ** Action  :   phdivn   increased by the ice shelf outflow
84      !!----------------------------------------------------------------------
85      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: phdiv
86      !!----------------------------------------------------------------------
87      INTEGER , DIMENSION(jpi,jpj), INTENT(in   ) :: ktop , kbot
88      REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) :: pfrac, phtbl
89      REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) :: pfwf , pfwf_b
90      !!----------------------------------------------------------------------
91      INTEGER  ::   ji, jj, jk   ! dummy loop indices
92      INTEGER  ::   ikt, ikb 
93      REAL(wp), DIMENSION(jpi,jpj) :: zhdiv
94      !!----------------------------------------------------------------------
95      !
96      !==   fwf distributed over several levels   ==!
97      !
98      ! compute integrated divergence correction
99      zhdiv(:,:) = 0.5_wp * ( pfwf(:,:) + pfwf_b(:,:) ) * r1_rho0 / phtbl(:,:)
100      !
101      ! update divergence at each level affected by ice shelf top boundary layer
102      DO_2D( 1, 1, 1, 1 )
103         ikt = ktop(ji,jj)
104         ikb = kbot(ji,jj)
105         ! level fully include in the ice shelf boundary layer
106         DO jk = ikt, ikb - 1
107            phdiv(ji,jj,jk) = phdiv(ji,jj,jk) + zhdiv(ji,jj)
108         END DO
109         ! level partially include in ice shelf boundary layer
110         phdiv(ji,jj,ikb) = phdiv(ji,jj,ikb) + zhdiv(ji,jj) * pfrac(ji,jj)
111      END_2D
112      !
113   END SUBROUTINE isf_hdiv_mlt
114
115   SUBROUTINE isf_hdiv_cpl(Kmm, pqvol, phdiv)
116      !!----------------------------------------------------------------------
117      !!                  ***  SUBROUTINE isf_hdiv_cpl  ***
118      !!       
119      !! ** Purpose :   update the horizontal divergence with the ice shelf
120      !!                coupling conservation increment
121      !!
122      !! ** Method  :   pqvol is positive (outflow) and expressed as m3/s
123      !!                increase the divergence
124      !!
125      !! ** Action  :   phdivn   increased by the ice shelf outflow
126      !!
127      !!----------------------------------------------------------------------
128      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: phdiv
129      !!----------------------------------------------------------------------
130      INTEGER,                          INTENT(in)    :: Kmm     ! ocean time level index
131      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in   ) :: pqvol
132      !!----------------------------------------------------------------------
133      INTEGER :: jk
134      !!----------------------------------------------------------------------
135      !
136      DO jk=1,jpk 
137         phdiv(:,:,jk) =  phdiv(:,:,jk) + pqvol(:,:,jk) * r1_e1e2t(:,:)   &
138            &                             / e3t(:,:,jk,Kmm)
139      END DO
140      !
141   END SUBROUTINE isf_hdiv_cpl
142
143END MODULE isfhdiv
Note: See TracBrowser for help on using the repository browser.