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.
diaar5.F90 in branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/DIA – NEMO

source: branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/DIA/diaar5.F90 @ 4400

Last change on this file since 4400 was 3211, checked in by spickles2, 13 years ago

Stephen Pickles, 11 Dec 2011

Commit to bring the rest of the DCSE NEMO development branch
in line with the latest development version. This includes
array index re-ordering of all OPA_SRC/.
  • Property svn:keywords set to Id
File size: 11.8 KB
Line 
1MODULE diaar5
2   !!======================================================================
3   !!                       ***  MODULE  diaar5  ***
4   !! AR5 diagnostics
5   !!======================================================================
6   !! History :  3.2  !  2009-11  (S. Masson)  Original code
7   !!            3.3  !  2010-10  (C. Ethe, G. Madec) reorganisation of initialisation phase + merge TRC-TRA
8   !!----------------------------------------------------------------------
9
10#if defined key_diaar5   || defined key_esopa
11   !!----------------------------------------------------------------------
12   !!   'key_diaar5'  :                           activate ar5 diagnotics
13   !!----------------------------------------------------------------------
14   !!   dia_ar5       : AR5 diagnostics
15   !!   dia_ar5_init  : initialisation of AR5 diagnostics
16   !!----------------------------------------------------------------------
17   USE oce            ! ocean dynamics and active tracers
18   USE dom_oce        ! ocean space and time domain
19   USE eosbn2         ! equation of state                (eos_bn2 routine)
20   USE lib_mpp        ! distributed memory computing library
21   USE iom            ! I/O manager library
22
23   IMPLICIT NONE
24   PRIVATE
25
26   PUBLIC   dia_ar5        ! routine called in step.F90 module
27   PUBLIC   dia_ar5_init   ! routine called in opa.F90 module
28   PUBLIC   dia_ar5_alloc  ! routine called in nemogcm.F90 module
29
30   LOGICAL, PUBLIC, PARAMETER :: lk_diaar5 = .TRUE.   ! coupled flag
31
32   REAL(wp)                         ::   vol0         ! ocean volume (interior domain)
33   REAL(wp)                         ::   area_tot     ! total ocean surface (interior domain)
34   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:  ) ::   area         ! cell surface (interior domain)
35   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:  ) ::   thick0       ! ocean thickness (interior domain)
36   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   sn0          ! initial salinity
37     
38   !! * Control permutation of array indices
39#  include "oce_ftrans.h90"
40#  include "dom_oce_ftrans.h90"
41!FTRANS sn0 :I :I :z
42
43   !! * Substitutions
44#  include "domzgr_substitute.h90"
45   !!----------------------------------------------------------------------
46   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
47   !! $Id$
48   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
49   !!----------------------------------------------------------------------
50CONTAINS
51
52   FUNCTION dia_ar5_alloc()
53      !!----------------------------------------------------------------------
54      !!                    ***  ROUTINE dia_ar5_alloc  ***
55      !!----------------------------------------------------------------------
56      INTEGER :: dia_ar5_alloc
57      !!----------------------------------------------------------------------
58      !
59      ALLOCATE( area(jpi,jpj), thick0(jpi,jpj) , sn0(jpi,jpj,jpk) , STAT=dia_ar5_alloc )
60      !
61      IF( lk_mpp             )   CALL mpp_sum ( dia_ar5_alloc )
62      IF( dia_ar5_alloc /= 0 )   CALL ctl_warn('dia_ar5_alloc: failed to allocate arrays')
63      !
64   END FUNCTION dia_ar5_alloc
65
66
67   SUBROUTINE dia_ar5( kt )
68      !!----------------------------------------------------------------------
69      !!                    ***  ROUTINE dia_ar5  ***
70      !!
71      !! ** Purpose :   compute and output some AR5 diagnostics
72      !!----------------------------------------------------------------------
73      USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released
74      USE wrk_nemo, ONLY:   zarea_ssh => wrk_2d_1 , zbotpres => wrk_2d_2   ! 2D workspace
75      USE wrk_nemo, ONLY:   zrhd      => wrk_3d_1 , zrhop    => wrk_3d_2   ! 3D      -
76      USE wrk_nemo, ONLY:   ztsn      => wrk_4d_1                          ! 4D      -
77
78      !! DCSE_NEMO: need additional directives for renamed module variables
79!FTRANS zrhd zrhop :I :I :z
80!FTRANS ztsn :I :I :z :
81
82      !
83      INTEGER, INTENT( in ) ::   kt   ! ocean time-step index
84      !
85      INTEGER  ::   ji, jj, jk                      ! dummy loop arguments
86      REAL(wp) ::   zvolssh, zvol, zssh_steric, zztmp, zarho, ztemp, zsal, zmass
87      !!--------------------------------------------------------------------
88
89      IF( wrk_in_use(2, 1,2) .OR.   &
90          wrk_in_use(3, 1,2) .OR.   &
91          wrk_in_use(4, 1)   ) THEN
92         CALL ctl_stop('dia_ar5: requested workspace arrays unavailable')   ;   RETURN
93      ENDIF
94
95      CALL iom_put( 'cellthc', fse3t(:,:,:) )
96
97      zarea_ssh(:,:) = area(:,:) * sshn(:,:)
98
99      !                                         ! total volume of liquid seawater
100      zvolssh = SUM( zarea_ssh(:,:) ) 
101      IF( lk_mpp )   CALL mpp_sum( zvolssh )
102      zvol = vol0 + zvolssh
103     
104      CALL iom_put( 'voltot', zvol               )
105      CALL iom_put( 'sshtot', zvolssh / area_tot )
106
107      !                                         ! thermosteric ssh
108      ztsn(:,:,:,jp_tem) = tn (:,:,:)
109      ztsn(:,:,:,jp_sal) = sn0(:,:,:)
110      CALL eos( ztsn, zrhd )                       ! now in situ density using initial salinity
111      !
112#if defined key_z_first
113      DO jj = 1, jpj
114         DO ji = 1, jpi
115            zbotpres(ji,jj) = 0._wp                ! no atmospheric surface pressure, levitating sea-ice
116            DO jk = 1, jpkm1
117               zbotpres(ji,jj) = zbotpres(ji,jj) + fse3t(ji,jj,jk) * zrhd(ji,jj,jk)
118            END DO
119      END DO
120#else
121      zbotpres(:,:) = 0._wp                        ! no atmospheric surface pressure, levitating sea-ice
122      DO jk = 1, jpkm1
123         zbotpres(:,:) = zbotpres(:,:) + fse3t(:,:,jk) * zrhd(:,:,jk)
124      END DO
125#endif
126      IF( .NOT.lk_vvl )   zbotpres(:,:) = zbotpres(:,:) + sshn(:,:) * zrhd(:,:,1)
127      !                                         
128      zarho = SUM( area(:,:) * zbotpres(:,:) ) 
129      IF( lk_mpp )   CALL mpp_sum( zarho )
130      zssh_steric = - zarho / area_tot
131      CALL iom_put( 'sshthster', zssh_steric )
132     
133      !                                         ! steric sea surface height
134      CALL eos( tsn, zrhd, zrhop )                 ! now in situ and potential density
135      zrhop(:,:,jpk) = 0._wp
136      CALL iom_put( 'rhop', zrhop )
137      !
138#if defined key_z_first
139      DO jj = 1, jpj
140         DO ji = 1, jpi
141            zbotpres(ji,jj) = 0._wp                ! no atmospheric surface pressure, levitating sea-ice
142            DO jk = 1, jpkm1
143               zbotpres(ji,jj) = zbotpres(ji,jj) + fse3t(ji,jj,jk) * zrhd(ji,jj,jk)
144            END DO
145      END DO
146#else
147      zbotpres(:,:) = 0._wp                        ! no atmospheric surface pressure, levitating sea-ice
148      DO jk = 1, jpkm1
149         zbotpres(:,:) = zbotpres(:,:) + fse3t(:,:,jk) * zrhd(:,:,jk)
150      END DO
151#endif
152      IF( .NOT.lk_vvl )   zbotpres(:,:) = zbotpres(:,:) + sshn(:,:) * zrhd(:,:,1)
153      !   
154      zarho = SUM( area(:,:) * zbotpres(:,:) ) 
155      IF( lk_mpp )   CALL mpp_sum( zarho )
156      zssh_steric = - zarho / area_tot
157      CALL iom_put( 'sshsteric', zssh_steric )
158     
159      !                                         ! ocean bottom pressure
160      zztmp = rau0 * grav * 1.e-4_wp               ! recover pressure from pressure anomaly and cover to dbar = 1.e4 Pa
161      zbotpres(:,:) = zztmp * ( zbotpres(:,:) + sshn(:,:) + thick0(:,:) )
162      CALL iom_put( 'botpres', zbotpres )
163
164      !                                         ! Mean density anomalie, temperature and salinity
165      ztemp = 0._wp
166      zsal  = 0._wp
167#if defined key_z_first
168      DO jj = 1, jpj
169         DO ji = 1, jpi
170            DO jk = 1, jpkm1
171#else
172      DO jk = 1, jpkm1
173         DO jj = 1, jpj
174            DO ji = 1, jpi
175#endif
176               zztmp = area(ji,jj) * fse3t(ji,jj,jk)
177               ztemp = ztemp + zztmp * tn(ji,jj,jk)
178               zsal  = zsal  + zztmp * sn(ji,jj,jk)
179            END DO
180         END DO
181      END DO
182      IF( .NOT.lk_vvl ) THEN
183         ztemp = ztemp + SUM( zarea_ssh(:,:) * tn(:,:,1) )
184         zsal  = zsal  + SUM( zarea_ssh(:,:) * sn(:,:,1) )
185      ENDIF
186      IF( lk_mpp ) THEN 
187         CALL mpp_sum( ztemp )
188         CALL mpp_sum( zsal  )
189      END IF
190      !
191      zmass = rau0 * ( zarho + zvol )                 ! total mass of liquid seawater
192      ztemp = ztemp / zvol                            ! potential temperature in liquid seawater
193      zsal  = zsal  / zvol                            ! Salinity of liquid seawater
194      !
195      CALL iom_put( 'masstot', zmass )
196      CALL iom_put( 'temptot', ztemp )
197      CALL iom_put( 'saltot' , zsal  )
198      !
199      IF( wrk_not_released(2, 1,2) .OR.   &
200          wrk_not_released(3, 1,2) .OR.   &
201          wrk_not_released(4, 1)   )   CALL ctl_stop('dia_ar5: failed to release workspace arrays')
202      !
203   END SUBROUTINE dia_ar5
204
205   !! * Reset control of array index permutation
206#  include "oce_ftrans.h90"
207#  include "dom_oce_ftrans.h90"
208!FTRANS sn0 :I :I :z
209
210   SUBROUTINE dia_ar5_init
211      !!----------------------------------------------------------------------
212      !!                  ***  ROUTINE dia_ar5_init  ***
213      !!                   
214      !! ** Purpose :   initialization for AR5 diagnostic computation
215      !!----------------------------------------------------------------------
216      USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released
217      USE wrk_nemo, ONLY:   wrk_4d_1      ! 4D workspace
218      !
219      INTEGER  ::   inum
220      INTEGER  ::   ik
221      INTEGER  ::   ji, jj, jk  ! dummy loop indices
222      REAL(wp) ::   zztmp 
223      REAL(wp), POINTER, DIMENSION(:,:,:,:) ::   zsaldta   ! Jan/Dec levitus salinity
224!FTRANS zsaldta :I :I :z :
225      !!----------------------------------------------------------------------
226      !
227      IF(wrk_in_use(4, 1) ) THEN
228         CALL ctl_stop('dia_ar5_init: requested workspace array unavailable.')   ;   RETURN
229      ENDIF
230      zsaldta => wrk_4d_1(:,:,:,1:2)
231
232      !                                      ! allocate dia_ar5 arrays
233      IF( dia_ar5_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'dia_ar5_init : unable to allocate arrays' )
234
235      area(:,:) = e1t(:,:) * e2t(:,:) * tmask_i(:,:)
236
237      area_tot = SUM( area(:,:) )   ;   IF( lk_mpp )   CALL mpp_sum( area_tot )
238
239      vol0        = 0._wp
240      thick0(:,:) = 0._wp
241      DO jk = 1, jpkm1
242         vol0        = vol0        + SUM( area (:,:) * tmask(:,:,jk) * fse3t_0(:,:,jk) )
243         thick0(:,:) = thick0(:,:) +    tmask_i(:,:) * tmask(:,:,jk) * fse3t_0(:,:,jk)
244      END DO
245      IF( lk_mpp )   CALL mpp_sum( vol0 )
246     
247      CALL iom_open ( 'data_1m_salinity_nomask', inum )
248      CALL iom_get  ( inum, jpdom_data, 'vosaline', zsaldta(:,:,:,1), 1  )
249      CALL iom_get  ( inum, jpdom_data, 'vosaline', zsaldta(:,:,:,2), 12 )
250      CALL iom_close( inum )
251      sn0(:,:,:) = 0.5_wp * ( zsaldta(:,:,:,1) + zsaldta(:,:,:,2) )       
252      sn0(:,:,:) = sn0(:,:,:) * tmask(:,:,:)
253      IF( ln_zps ) THEN               ! z-coord. partial steps
254         DO jj = 1, jpj               ! interpolation of salinity at the last ocean level (i.e. the partial step)
255            DO ji = 1, jpi
256               ik = mbkt(ji,jj)
257               IF( ik > 1 ) THEN
258                  zztmp = ( gdept_0(ik) - fsdept_0(ji,jj,ik) ) / ( gdept_0(ik) - gdept_0(ik-1) )
259                  sn0(ji,jj,ik) = ( 1._wp - zztmp ) * sn0(ji,jj,ik) + zztmp * sn0(ji,jj,ik-1)
260               ENDIF
261            END DO
262         END DO
263      ENDIF
264      !
265      IF( wrk_not_released(4, 1) )   CALL ctl_stop('dia_ar5_init: failed to release workspace array')
266      !
267   END SUBROUTINE dia_ar5_init
268
269#else
270   !!----------------------------------------------------------------------
271   !!   Default option :                                         NO diaar5
272   !!----------------------------------------------------------------------
273   LOGICAL, PUBLIC, PARAMETER :: lk_diaar5 = .FALSE.   ! coupled flag
274CONTAINS
275   SUBROUTINE dia_ar5_init    ! Dummy routine
276   END SUBROUTINE dia_ar5_init
277   SUBROUTINE dia_ar5( kt )   ! Empty routine
278      INTEGER ::   kt
279      WRITE(*,*) 'dia_ar5: You should not have seen this print! error?', kt
280   END SUBROUTINE dia_ar5
281#endif
282
283   !!======================================================================
284END MODULE diaar5
Note: See TracBrowser for help on using the repository browser.