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.
agrif_lim3_update.F90 in branches/2017/dev_merge_2017/NEMOGCM/NEMO/NST_SRC – NEMO

source: branches/2017/dev_merge_2017/NEMOGCM/NEMO/NST_SRC/agrif_lim3_update.F90 @ 9454

Last change on this file since 9454 was 9454, checked in by clem, 6 years ago

make sure that agrif is working with sea ice in any conditions

File size: 10.7 KB
Line 
1#define TWO_WAY
2
3MODULE agrif_lim3_update
4   !!=====================================================================================
5   !!                       ***  MODULE agrif_lim3_update ***
6   !! Nesting module :  update surface ocean boundary condition over ice from a child grid
7   !! Sea-Ice model  :  LIM 3.6 Sea ice model time-stepping
8   !!=====================================================================================
9   !! History :  2.0   !  04-2008  (F. Dupont)  initial version
10   !!            3.4   !  08-2012  (R. Benshila, C. Herbaut) update and EVP
11   !!            3.6   !  05-2016  (C. Rousset)  Add LIM3 compatibility
12   !!----------------------------------------------------------------------
13#if defined key_agrif && defined key_lim3
14   !!----------------------------------------------------------------------
15   !!   'key_lim3'  :                                 LIM 3.6 sea-ice model
16   !!   'key_agrif' :                                 AGRIF library
17   !!----------------------------------------------------------------------
18   !!   agrif_update_lim3  : update sea-ice on boundaries or total
19   !!                        child domain for velocities and ice properties
20   !!   update_tra_ice     : sea-ice properties
21   !!   update_u_ice       : zonal      ice velocity
22   !!   update_v_ice       : meridional ice velocity
23   !!----------------------------------------------------------------------
24   USE dom_oce
25   USE sbc_oce
26   USE agrif_oce
27   USE ice
28   USE agrif_ice 
29   USE phycst , ONLY: rt0
30
31   IMPLICIT NONE
32   PRIVATE
33
34   PUBLIC   agrif_update_lim3   ! called by agrif_user.F90
35
36   !!----------------------------------------------------------------------
37   !! NEMO/NST 4.0 , LOCEAN-IPSL (2017)
38   !! $Id: agrif_lim3_update.F90 6204 2016-01-04 13:47:06Z cetlod $
39   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
40   !!----------------------------------------------------------------------
41CONTAINS
42
43   SUBROUTINE agrif_update_lim3( kt )
44      !!----------------------------------------------------------------------
45      !!                     *** ROUTINE agrif_update_lim3 ***
46      !! ** Method  :   Call the hydrostaticupdate pressure at the boundary or the entire domain
47      !!
48      !! ** Action : - Update (u_ice,v_ice) and ice tracers
49      !!----------------------------------------------------------------------
50      INTEGER, INTENT(in) :: kt
51      !!----------------------------------------------------------------------
52      !
53      !! clem: I think the update should take place each time the ocean sees the surface forcings
54      !!       (but maybe I am wrong and we should update every rhot time steps)
55      IF( ( MOD( (kt-nit000)/nn_fsbc + 1, Agrif_irhot() * Agrif_Parent(nn_fsbc) / nn_fsbc ) /=0 ) .AND. (kt /= 0) ) RETURN ! do not update if nb of child time steps differ from time refinement
56                                                                                                                           ! i.e. update only at the parent time step
57      IF( nn_ice == 0 ) RETURN   ! do not update if child domain does not have ice
58      !
59      Agrif_SpecialValueFineGrid    = -9999.
60      Agrif_UseSpecialValueInUpdate = .TRUE.
61# if defined TWO_WAY
62      CALL Agrif_Update_Variable( tra_ice_id , procname = update_tra_ice  )
63      CALL Agrif_Update_Variable( u_ice_id   , procname = update_u_ice    )
64      CALL Agrif_Update_Variable( v_ice_id   , procname = update_v_ice    )
65
66!      CALL Agrif_Update_Variable( tra_ice_id , locupdate=(/0,2/), procname = update_tra_ice  )
67!      CALL Agrif_Update_Variable( u_ice_id   , locupdate=(/0,1/), procname = update_u_ice    )
68!      CALL Agrif_Update_Variable( v_ice_id   , locupdate=(/0,1/), procname = update_v_ice    )
69# endif
70      Agrif_SpecialValueFineGrid    = 0.
71      Agrif_UseSpecialValueInUpdate = .FALSE.
72      !
73   END SUBROUTINE agrif_update_lim3
74
75
76   SUBROUTINE update_tra_ice( ptab, i1, i2, j1, j2, k1, k2, before )
77      !!-----------------------------------------------------------------------
78      !!                        *** ROUTINE update_tra_ice ***
79      !! ** Method  : Compute the mass properties on the fine grid and recover
80      !!              the properties per mass on the coarse grid
81      !!-----------------------------------------------------------------------
82      INTEGER                               , INTENT(in   ) ::   i1, i2, j1, j2, k1, k2
83      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) ::   ptab
84      LOGICAL                               , INTENT(in   ) ::   before
85      !!
86      INTEGER  :: ji, jj, jk, jl, jm
87      !!-----------------------------------------------------------------------
88      ! it is ok not to multiply by e1*e2 since we conserve tracers here (same as in the ocean).
89      IF( before ) THEN
90         jm = 1
91         DO jl = 1, jpl
92            ptab(i1:i2,j1:j2,jm  ) = a_i (i1:i2,j1:j2,jl)
93            ptab(i1:i2,j1:j2,jm+1) = v_i (i1:i2,j1:j2,jl)
94            ptab(i1:i2,j1:j2,jm+2) = v_s (i1:i2,j1:j2,jl)
95            ptab(i1:i2,j1:j2,jm+3) = sv_i(i1:i2,j1:j2,jl)
96            ptab(i1:i2,j1:j2,jm+4) = oa_i(i1:i2,j1:j2,jl)
97            ptab(i1:i2,j1:j2,jm+5) = a_ip(i1:i2,j1:j2,jl)
98            ptab(i1:i2,j1:j2,jm+6) = v_ip(i1:i2,j1:j2,jl)
99            ptab(i1:i2,j1:j2,jm+7) = t_su(i1:i2,j1:j2,jl)
100            jm = jm + 8
101            DO jk = 1, nlay_s
102               ptab(i1:i2,j1:j2,jm) = e_s(i1:i2,j1:j2,jk,jl)   ;   jm = jm + 1
103            END DO
104            DO jk = 1, nlay_i
105               ptab(i1:i2,j1:j2,jm) = e_i(i1:i2,j1:j2,jk,jl)   ;   jm = jm + 1
106            END DO
107         END DO
108         !
109         DO jk = k1, k2
110            WHERE( tmask(i1:i2,j1:j2,1) == 0. )   ptab(i1:i2,j1:j2,jk) = Agrif_SpecialValueFineGrid 
111         END DO
112         !
113      ELSE
114         !
115         jm = 1
116         DO jl = 1, jpl
117            !
118            DO jj = j1, j2
119               DO ji = i1, i2
120                  IF( ptab(ji,jj,jm) /= Agrif_SpecialValueFineGrid ) THEN
121                     a_i (ji,jj,jl) = ptab(ji,jj,jm  ) * tmask(ji,jj,1)
122                     v_i (ji,jj,jl) = ptab(ji,jj,jm+1) * tmask(ji,jj,1)
123                     v_s (ji,jj,jl) = ptab(ji,jj,jm+2) * tmask(ji,jj,1)
124                     sv_i(ji,jj,jl) = ptab(ji,jj,jm+3) * tmask(ji,jj,1)
125                     oa_i(ji,jj,jl) = ptab(ji,jj,jm+4) * tmask(ji,jj,1)
126                     a_ip(ji,jj,jl) = ptab(ji,jj,jm+5) * tmask(ji,jj,1)
127                     v_ip(ji,jj,jl) = ptab(ji,jj,jm+6) * tmask(ji,jj,1)
128                     t_su(ji,jj,jl) = ptab(ji,jj,jm+7) * tmask(ji,jj,1)
129                  ENDIF
130               END DO
131            END DO
132            jm = jm + 8
133            !
134            DO jk = 1, nlay_s
135               WHERE( ptab(i1:i2,j1:j2,jm) /= Agrif_SpecialValueFineGrid )
136                  e_s(i1:i2,j1:j2,jk,jl) = ptab(i1:i2,j1:j2,jm) * tmask(i1:i2,j1:j2,1)
137               ENDWHERE
138               jm = jm + 1
139            END DO
140            !
141            DO jk = 1, nlay_i
142               WHERE( ptab(i1:i2,j1:j2,jm) /= Agrif_SpecialValueFineGrid )
143                  e_i(i1:i2,j1:j2,jk,jl) = ptab(i1:i2,j1:j2,jm) * tmask(i1:i2,j1:j2,1)
144               ENDWHERE
145               jm = jm + 1
146            END DO
147            !
148         END DO
149         !
150         ! integrated values
151         vt_i(i1:i2,j1:j2) = SUM(      v_i(i1:i2,j1:j2,:)           , dim=3 )
152         vt_s(i1:i2,j1:j2) = SUM(      v_s(i1:i2,j1:j2,:)           , dim=3 )
153         at_i(i1:i2,j1:j2) = SUM(      a_i(i1:i2,j1:j2,:)           , dim=3 )
154         et_s(i1:i2,j1:j2) = SUM( SUM( e_s(i1:i2,j1:j2,:,:), dim=4 ), dim=3 )
155         et_i(i1:i2,j1:j2) = SUM( SUM( e_i(i1:i2,j1:j2,:,:), dim=4 ), dim=3 )
156
157         at_ip(i1:i2,j1:j2) = SUM( a_ip(i1:i2,j1:j2,:), dim=3 ) ! melt ponds
158         vt_ip(i1:i2,j1:j2) = SUM( v_ip(i1:i2,j1:j2,:), dim=3 )
159         !
160         ato_i(i1:i2,j1:j2) = 1._wp - at_i(i1:i2,j1:j2)         ! open water fraction 
161
162         DO jl = 1, jpl
163            WHERE( tmask(i1:i2,j1:j2,1) == 0._wp )   t_su(i1:i2,j1:j2,jl) = rt0   ! to avoid a division by 0 in sbcblk.F90
164         END DO
165         
166      ENDIF
167      !
168   END SUBROUTINE update_tra_ice
169
170
171   SUBROUTINE update_u_ice( ptab, i1, i2, j1, j2, before )
172      !!-----------------------------------------------------------------------
173      !!                        *** ROUTINE update_u_ice ***
174      !! ** Method  : Update the fluxes and recover the properties (C-grid)
175      !!-----------------------------------------------------------------------
176      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
177      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
178      LOGICAL                         , INTENT(in   ) ::   before
179      !!
180      REAL(wp) ::   zrhoy   ! local scalar
181      !!-----------------------------------------------------------------------
182      !
183      IF( before ) THEN
184         zrhoy = Agrif_Rhoy()
185         ptab(:,:) = e2u(i1:i2,j1:j2) * u_ice(i1:i2,j1:j2) * zrhoy
186         WHERE( umask(i1:i2,j1:j2,1) == 0._wp )   ptab(:,:) = Agrif_SpecialValueFineGrid
187      ELSE
188         WHERE( ptab(i1:i2,j1:j2) /= Agrif_SpecialValueFineGrid )
189            u_ice(i1:i2,j1:j2) = ptab(i1:i2,j1:j2) / e2u(i1:i2,j1:j2) * umask(i1:i2,j1:j2,1)
190         ENDWHERE
191      ENDIF
192      !
193   END SUBROUTINE update_u_ice
194
195
196   SUBROUTINE update_v_ice( ptab, i1, i2, j1, j2, before )
197      !!-----------------------------------------------------------------------
198      !!                    *** ROUTINE update_v_ice ***
199      !! ** Method  : Update the fluxes and recover the properties (C-grid)
200      !!-----------------------------------------------------------------------
201      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
202      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
203      LOGICAL                         , INTENT(in   ) ::   before
204      !!
205      REAL(wp) ::   zrhox   ! local scalar
206      !!-----------------------------------------------------------------------
207      !
208      IF( before ) THEN
209         zrhox = Agrif_Rhox()
210         ptab(:,:) = e1v(i1:i2,j1:j2) * v_ice(i1:i2,j1:j2) * zrhox
211         WHERE( vmask(i1:i2,j1:j2,1) == 0._wp )   ptab(:,:) = Agrif_SpecialValueFineGrid
212      ELSE
213         WHERE( ptab(i1:i2,j1:j2) /= Agrif_SpecialValueFineGrid )
214            v_ice(i1:i2,j1:j2) = ptab(i1:i2,j1:j2) / e1v(i1:i2,j1:j2) * vmask(i1:i2,j1:j2,1)
215         ENDWHERE
216      ENDIF
217      !
218   END SUBROUTINE update_v_ice
219
220#else
221   !!----------------------------------------------------------------------
222   !!   Empty module                                             no sea-ice
223   !!----------------------------------------------------------------------
224CONTAINS
225   SUBROUTINE agrif_lim3_update_empty
226      WRITE(*,*)  'agrif_lim3_update : You should not have seen this print! error?'
227   END SUBROUTINE agrif_lim3_update_empty
228#endif
229
230   !!======================================================================
231END MODULE agrif_lim3_update
Note: See TracBrowser for help on using the repository browser.