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.
icbclv.F90 in NEMO/branches/UKMO/NEMO_4.0.1_icesheet_and_river_coupling/src/OCE/ICB – NEMO

source: NEMO/branches/UKMO/NEMO_4.0.1_icesheet_and_river_coupling/src/OCE/ICB/icbclv.F90 @ 12576

Last change on this file since 12576 was 12576, checked in by dancopsey, 2 years ago

Merge in iceberg calving stuff from dev_r5518_coupling_GSI7_GSI8_landice from the start of the branch to revision 6023.

File size: 11.7 KB
Line 
1MODULE icbclv
2   !!======================================================================
3   !!                       ***  MODULE  icbclv  ***
4   !! Icebergs:  calving routines for iceberg calving
5   !!======================================================================
6   !! History : 3.3.1 !  2010-01  (Martin&Adcroft) Original code
7   !!            -    !  2011-03  (Madec)          Part conversion to NEMO form
8   !!            -    !                            Removal of mapping from another grid
9   !!            -    !  2011-04  (Alderson)       Split into separate modules
10   !!            -    !  2011-05  (Alderson)       budgets into separate module
11   !!----------------------------------------------------------------------
12
13   !!----------------------------------------------------------------------
14   !!   icb_clv_flx   : transfer input flux of ice into iceberg classes
15   !!   icb_clv       : calve icebergs from stored ice
16   !!----------------------------------------------------------------------
17   USE par_oce        ! NEMO parameters
18   USE dom_oce        ! NEMO ocean domain
19   USE phycst         ! NEMO physical constants
20   USE lib_mpp        ! NEMO MPI library, lk_mpp in particular
21   USE lbclnk         ! NEMO boundary exchanges for gridded data
22
23   USE icb_oce        ! iceberg variables
24   USE icbdia         ! iceberg diagnostics
25   USE icbutl         ! iceberg utility routines
26   USE icb_oce        ! iceberg parameters
27
28   USE sbc_oce        ! for icesheet freshwater input variables
29   USE in_out_manager 
30   USE iom 
31
32   IMPLICIT NONE
33   PRIVATE
34
35   PUBLIC   icb_clv_flx  ! routine called in icbstp.F90 module
36   PUBLIC   icb_clv      ! routine called in icbstp.F90 module
37
38   !!----------------------------------------------------------------------
39   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
40   !! $Id$
41   !! Software governed by the CeCILL license (see ./LICENSE)
42   !!----------------------------------------------------------------------
43CONTAINS
44
45   SUBROUTINE icb_clv_flx( kt )
46      !!----------------------------------------------------------------------
47      !!                 ***  ROUTINE icb_clv_flx  ***
48      !!
49      !! ** Purpose :   accumulate ice available for calving into class arrays
50      !!
51      !!----------------------------------------------------------------------
52      INTEGER, INTENT(in) ::   kt
53      !
54      REAL(wp)      ::   zcalving_used, zdist, zfact
55      REAL(wp)      ::   zgreenland_calving_sum, zantarctica_calving_sum 
56      INTEGER       ::   jn, ji, jj                    ! loop counters
57      INTEGER       ::   imx                           ! temporary integer for max berg class
58      LOGICAL, SAVE ::   ll_first_call = .TRUE.
59      !!----------------------------------------------------------------------
60      !
61      ! Adapt calving flux and calving heat flux from coupler for use here
62      ! Use interior mask: so no bergs in overlap areas and convert from km^3/year to kg/s
63      ! this assumes that input is given as equivalent water flux so that pure water density is appropriate
64
65      zfact = ( (1000._wp)**3 / ( NINT(rday) * nyear_len(1) ) ) * rn_rho_bergs
66      berg_grid%calving(:,:) = src_calving(:,:) * zfact * tmask_i(:,:) * tmask(:,:,1)
67
68      ! Heat in units of W/m2, and mask (just in case)
69      berg_grid%calving_hflx(:,:) = src_calving_hflx(:,:) * tmask_i(:,:) * tmask(:,:,1)
70
71      IF( lk_oasis) THEN
72        ! nn_coupled_iceshelf_fluxes uninitialised unless lk_oasis=true
73        IF( nn_coupled_iceshelf_fluxes .gt. 0 ) THEN
74          ll_write = ((MOD( kt, sn_cfctl%ptimincr ) == 0) .OR. ( kt == nitend )) .AND. lwp
75          ! Adjust total calving rates so that sum of iceberg calving and iceshelf melting in the northern
76          ! and southern hemispheres equals rate of increase of mass of greenland and antarctic ice sheets
77          ! to preserve total freshwater conservation in coupled models without an active ice sheet model.
78
79           zgreenland_calving_sum = SUM( berg_grid%calving(:,:) * greenland_icesheet_mask(:,:) )
80           IF( lk_mpp ) CALL mpp_sum( zgreenland_calving_sum )
81           WHERE( greenland_icesheet_mask(:,:) == 1.0 )                                                                                 &
82          &    berg_grid%calving(:,:) = berg_grid%calving(:,:) * greenland_icesheet_mass_rate_of_change * rn_greenland_calving_fraction &
83          &                                     / ( zgreenland_calving_sum + 1.0e-10_wp )
84
85           ! check
86           IF(ll_write) WRITE(numout, *) 'Greenland iceberg calving climatology (kg/s) : ',zgreenland_calving_sum
87           zgreenland_calving_sum = SUM( berg_grid%calving(:,:) * greenland_icesheet_mask(:,:) )
88           IF( lk_mpp ) CALL mpp_sum( zgreenland_calving_sum )
89           IF(ll_write) WRITE(numout, *) 'Greenland iceberg calving adjusted value (kg/s) : ',zgreenland_calving_sum
90
91           zantarctica_calving_sum = SUM( berg_grid%calving(:,:) * antarctica_icesheet_mask(:,:) )
92           IF( lk_mpp ) CALL mpp_sum( zantarctica_calving_sum )
93           WHERE( antarctica_icesheet_mask(:,:) == 1.0 )                                                                              &
94           berg_grid%calving(:,:) = berg_grid%calving(:,:) * antarctica_icesheet_mass_rate_of_change * rn_antarctica_calving_fraction &
95          &                           / ( zantarctica_calving_sum + 1.0e-10_wp )
96
97           ! check
98           IF(ll_write) WRITE(numout, *) 'Antarctica iceberg calving climatology (kg/s) : ',zantarctica_calving_sum
99           zantarctica_calving_sum = SUM( berg_grid%calving(:,:) * antarctica_icesheet_mask(:,:) )
100           IF( lk_mpp ) CALL mpp_sum( zantarctica_calving_sum )
101           IF(ll_write) WRITE(numout, *) 'Antarctica iceberg calving adjusted value (kg/s) : ',zantarctica_calving_sum
102
103        ENDIF
104      ENDIF
105   
106      CALL iom_put( 'berg_calve', berg_grid%calving(:,:) )
107
108
109      IF( ll_first_call .AND. .NOT. l_restarted_bergs ) THEN      ! This is a hack to simplify initialization
110         ll_first_call = .FALSE.
111         !do jn=1, nclasses
112         !  where (berg_grid%calving==0.) berg_grid%stored_ice(:,:,jn)=0.
113         !end do
114         DO jj = 2, jpjm1
115            DO ji = 2, jpim1
116               IF( berg_grid%calving(ji,jj) /= 0._wp )                                          &    ! Need units of J
117                  berg_grid%stored_heat(ji,jj) = SUM( berg_grid%stored_ice(ji,jj,:) ) *         &    ! initial stored ice in kg
118                     &                   berg_grid%calving_hflx(ji,jj) * e1e2t(ji,jj) / berg_grid%calving(ji,jj)   ! J/s/m2 x m^2
119                     !                                                                                             ! = J/s/calving in kg/s
120            END DO
121         END DO
122      ENDIF
123
124      ! assume that all calving flux must be distributed even if distribution array does not sum
125      ! to one - this may not be what is intended, but it's what you've got
126      DO jj = 1, jpj
127         DO ji = 1, jpi
128            imx = berg_grid%maxclass(ji,jj)
129            zdist = SUM( rn_distribution(1:nclasses) ) / SUM( rn_distribution(1:imx) )
130            DO jn = 1, imx
131               berg_grid%stored_ice(ji,jj,jn) = berg_grid%stored_ice(ji,jj,jn)     &
132                  &                           + berg_dt * berg_grid%calving(ji,jj) * rn_distribution(jn) * zdist
133            END DO
134         END DO
135      END DO
136
137      ! before changing the calving, save the amount we're about to use and do budget
138      zcalving_used = SUM( berg_grid%calving(:,:) )
139      berg_grid%tmp(:,:) = berg_dt * berg_grid%calving_hflx(:,:) * e1e2t(:,:) * tmask_i(:,:)
140      berg_grid%stored_heat (:,:) = berg_grid%stored_heat (:,:) + berg_grid%tmp(:,:)
141      CALL icb_dia_income( kt,  zcalving_used, berg_grid%tmp )
142      !
143   END SUBROUTINE icb_clv_flx
144
145
146   SUBROUTINE icb_clv( kt )
147      !!----------------------------------------------------------------------
148      !!                 ***  ROUTINE icb_clv  ***
149      !!
150      !! ** Purpose :   This routine takes a stored ice field and calves to the ocean,
151      !!                so the gridded array stored_ice has only non-zero entries at selected
152      !!                wet points adjacent to known land based calving points
153      !!
154      !! ** method  : - Look at each grid point and see if there's enough for each size class to calve
155      !!                If there is, a new iceberg is calved.  This happens in the order determined by
156      !!                the class definition arrays (which in the default case is smallest first)
157      !!                Note that only the non-overlapping part of the processor where icebergs are allowed
158      !!                is considered
159      !!----------------------------------------------------------------------
160      INTEGER, INTENT(in) ::   kt
161      INTEGER       ::   ji, jj, jn   ! dummy loop indices
162      INTEGER       ::   icnt, icntmax
163      TYPE(iceberg) ::   newberg
164      TYPE(point)   ::   newpt
165      REAL(wp)      ::   zday, zcalved_to_berg, zheat_to_berg
166      !!----------------------------------------------------------------------
167      !
168      icntmax = 0
169      zday    = REAL(nday_year,wp) + REAL(nsec_day,wp)/86400.0_wp
170      !
171      DO jn = 1, nclasses
172         DO jj = nicbdj, nicbej
173            DO ji = nicbdi, nicbei
174               !
175               icnt = 0
176               !
177               DO WHILE (berg_grid%stored_ice(ji,jj,jn) >= rn_initial_mass(jn) * rn_mass_scaling(jn) )
178                  !
179                  newpt%lon = glamt(ji,jj)         ! at t-point (centre of the cell)
180                  newpt%lat = gphit(ji,jj)
181                  newpt%xi  = REAL( mig(ji), wp )
182                  newpt%yj  = REAL( mjg(jj), wp )
183                  !
184                  newpt%uvel = 0._wp               ! initially at rest
185                  newpt%vvel = 0._wp
186                  !                                ! set berg characteristics
187                  newpt%mass           = rn_initial_mass     (jn)
188                  newpt%thickness      = rn_initial_thickness(jn)
189                  newpt%width          = first_width         (jn)
190                  newpt%length         = first_length        (jn)
191                  newberg%mass_scaling = rn_mass_scaling     (jn)
192                  newpt%mass_of_bits   = 0._wp                          ! no bergy
193                  !
194                  newpt%year   = nyear
195                  newpt%day    = zday
196                  newpt%heat_density = berg_grid%stored_heat(ji,jj) / berg_grid%stored_ice(ji,jj,jn)   ! This is in J/kg
197                  !
198                  CALL icb_utl_incr()
199                  newberg%number(:) = num_bergs(:)
200                  !
201                  CALL icb_utl_add( newberg, newpt )
202                  !
203                  zcalved_to_berg = rn_initial_mass(jn) * rn_mass_scaling(jn)           ! Units of kg
204                  !                                ! Heat content
205                  zheat_to_berg           = zcalved_to_berg * newpt%heat_density             ! Units of J
206                  berg_grid%stored_heat(ji,jj) = berg_grid%stored_heat(ji,jj) - zheat_to_berg
207                  !                                ! Stored mass
208                  berg_grid%stored_ice(ji,jj,jn) = berg_grid%stored_ice(ji,jj,jn) - zcalved_to_berg
209                  !
210                  icnt = icnt + 1
211                  !
212                  CALL icb_dia_calve(ji, jj, jn,  zcalved_to_berg, zheat_to_berg )
213               END DO
214               icntmax = MAX( icntmax, icnt )
215            END DO
216         END DO
217      END DO
218      !
219      DO jn = 1, nclasses
220         CALL lbc_lnk( 'icbclv', berg_grid%stored_ice(:,:,jn), 'T', 1._wp )
221      END DO
222      CALL lbc_lnk( 'icbclv', berg_grid%stored_heat, 'T', 1._wp )
223      !
224      IF( nn_verbose_level > 0 .AND. icntmax > 1 )   WRITE(numicb,*) 'icb_clv: icnt=', icnt,' on', narea
225      !
226   END SUBROUTINE  icb_clv
227
228   !!======================================================================
229END MODULE icbclv
Note: See TracBrowser for help on using the repository browser.