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.
asmbkg.F90 in branches/UKMO/AMM15_v3_6_STABLE_package_collate_BGC_DA/NEMOGCM/NEMO/OPA_SRC/ASM – NEMO

source: branches/UKMO/AMM15_v3_6_STABLE_package_collate_BGC_DA/NEMOGCM/NEMO/OPA_SRC/ASM/asmbkg.F90 @ 10622

Last change on this file since 10622 was 10622, checked in by dford, 5 years ago

Implement biogeochemistry assimilation for FABM-ERSEM.
File size: 13.2 KB
Line 
1MODULE asmbkg
2   !!======================================================================
3   !!                       ***  MODULE asmtrj -> asmbkg  ***
4   !! Assimilation trajectory interface: Write to file the background state and the model state trajectory
5   !!======================================================================
6   !! History :       ! 2007-03 (M. Martin)  Met. Office version
7   !!                 ! 2007-04 (A. Weaver)  asm_trj_wri, original code
8   !!                 ! 2007-03 (K. Mogensen)  Adapt to NEMOVAR and use IOM instead of IOIPSL
9   !!                 ! 2007-04 (A. Weaver)  Name change (formally asmbkg.F90). Distinguish
10   !!                                        background states in Jb term and at analysis time.
11   !!                                        Include state trajectory routine (currently empty)
12   !!                 ! 2007-07 (A. Weaver)  Add tke_rst and flt_rst for case nitbkg=0
13   !!                 ! 2009-03 (F. Vigilant)  Add hmlp (zdfmxl) for no tracer nmldp=2
14   !!                 ! 2009-06 (F. Vigilant) asm_trj_wri: special case when kt=nit000-1
15   !!                 ! 2009-07 (F. Vigilant) asm_trj_wri: add computation of eiv at restart
16   !!                 ! 2010-01 (A. Vidard) split asm_trj_wri into tam_trj_wri and asm_bkg_wri
17   !!----------------------------------------------------------------------
18
19   !!----------------------------------------------------------------------
20   !!   'key_asminc' : Switch on the assimilation increment interface
21   !!----------------------------------------------------------------------
22   !!   asm_bkg_wri  : Write out the background state
23   !!   asm_trj_wri  : Write out the model state trajectory (used with 4D-Var)
24   !!----------------------------------------------------------------------
25   USE oce                ! Dynamics and active tracers defined in memory
26   USE sbc_oce            ! Ocean surface boundary conditions
27   USE zdf_oce            ! Vertical mixing variables
28   USE zdfddm             ! Double diffusion mixing parameterization
29   USE ldftra_oce         ! Lateral tracer mixing coefficient defined in memory
30   USE ldfslp             ! Slopes of neutral surfaces
31   USE tradmp             ! Tracer damping
32#if defined key_zdftke
33   USE zdftke             ! TKE vertical physics
34#endif
35   USE eosbn2             ! Equation of state (eos_bn2 routine)
36   USE zdfmxl             ! Mixed layer depth
37   USE dom_oce, ONLY :   ndastp
38   USE sol_oce, ONLY :   gcx   ! Solver variables defined in memory
39   USE in_out_manager     ! I/O manager
40   USE iom                ! I/O module
41   USE asmpar             ! Parameters for the assmilation interface
42   USE zdfmxl             ! mixed layer depth
43#if defined key_traldf_c2d
44   USE ldfeiv             ! eddy induced velocity coef.      (ldf_eiv routine)
45#endif
46#if defined key_lim2
47   USE ice_2
48#endif
49#if defined key_lim3
50   USE ice
51#endif
52   USE asminc, ONLY: ln_avgbkg
53#if defined key_top
54   USE asmbgc, ONLY: asm_bgc_bkg_alloc, &
55      &              asm_bgc_bkg_tavg,  &
56      &              asm_bgc_bkg_wri
57#endif
58   IMPLICIT NONE
59   PRIVATE
60   
61   PUBLIC   asm_bkg_wri   !: Write out the background state
62
63  !! * variables for calculating time means
64   REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   tn_tavg  , sn_tavg 
65   REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   un_tavg  , vn_tavg
66   REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   avt_tavg
67#if defined key_zdfgls || key_zdftke
68   REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   en_tavg
69#endif
70   REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:)   ::   sshn_tavg
71   REAL(wp),SAVE :: numtimes_tavg     ! No of times to average over
72
73   !!----------------------------------------------------------------------
74   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
75   !! $Id$
76   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
77   !!----------------------------------------------------------------------
78CONTAINS
79
80   SUBROUTINE asm_bkg_wri( kt )
81      !!-----------------------------------------------------------------------
82      !!                  ***  ROUTINE asm_bkg_wri ***
83      !!
84      !! ** Purpose : Write to file the background state for later use in the
85      !!              inner loop of data assimilation or for direct initialization
86      !!              in the outer loop.
87      !!
88      !! ** Method  : Write out the background state for use in the Jb term
89      !!              in the cost function and for use with direct initialization
90      !!              at analysis time.
91      !!-----------------------------------------------------------------------
92      INTEGER, INTENT( IN ) :: kt               ! Current time-step
93      !
94      CHARACTER (LEN=50) :: cl_asmbkg
95      CHARACTER (LEN=50) :: cl_asmdin
96      LOGICAL :: llok          ! Check if file exists
97      INTEGER :: inum          ! File unit number
98      REAL(wp) :: zdate        ! Date
99      INTEGER :: ierror
100      !!-----------------------------------------------------------------------
101
102      ! If creating an averaged assim bkg, initialise on first timestep
103      IF ( ln_avgbkg .AND. kt == ( nn_it000 - 1) ) THEN
104         ! Allocate memory
105         ALLOCATE( tn_tavg(jpi,jpj,jpk), STAT=ierror )
106         IF( ierror > 0 ) THEN
107            CALL ctl_stop( 'asm_wri_bkg: unable to allocate tn_tavg' )   ;   RETURN
108         ENDIF
109         tn_tavg(:,:,:)=0
110         ALLOCATE( sn_tavg(jpi,jpj,jpk), STAT=ierror )
111         IF( ierror > 0 ) THEN
112            CALL ctl_stop( 'asm_wri_bkg: unable to allocate sn_tavg' )   ;   RETURN
113         ENDIF
114         sn_tavg(:,:,:)=0
115         ALLOCATE( un_tavg(jpi,jpj,jpk), STAT=ierror )
116         IF( ierror > 0 ) THEN
117            CALL ctl_stop( 'asm_wri_bkg: unable to allocate un_tavg' )   ;   RETURN
118         ENDIF
119         un_tavg(:,:,:)=0
120         ALLOCATE( vn_tavg(jpi,jpj,jpk), STAT=ierror )
121         IF( ierror > 0 ) THEN
122            CALL ctl_stop( 'asm_wri_bkg: unable to allocate vn_tavg' )   ;   RETURN
123         ENDIF
124         vn_tavg(:,:,:)=0
125         ALLOCATE( sshn_tavg(jpi,jpj), STAT=ierror )
126         IF( ierror > 0 ) THEN
127            CALL ctl_stop( 'asm_wri_bkg: unable to allocate sshn_tavg' )   ;   RETURN
128         ENDIF
129         sshn_tavg(:,:)=0
130#if defined key_zdftke
131         ALLOCATE( en_tavg(jpi,jpj,jpk), STAT=ierror )
132         IF( ierror > 0 ) THEN
133            CALL ctl_stop( 'asm_wri_bkg: unable to allocate en_tavg' )   ;   RETURN
134         ENDIF
135         en_tavg(:,:,:)=0
136#endif
137         ALLOCATE( avt_tavg(jpi,jpj,jpk), STAT=ierror )
138         IF( ierror > 0 ) THEN
139            CALL ctl_stop( 'asm_wri_bkg: unable to allocate avt_tavg' )   ;   RETURN
140         ENDIF
141         avt_tavg(:,:,:)=0
142         
143         numtimes_tavg = REAL ( nitavgbkg_r -  nn_it000 + 1 )
144      ENDIF
145     
146#if defined key_top
147      ! Allocate BGC average arrays whatever, to save code repetition later
148      IF ( kt == ( nn_it000 - 1) ) THEN
149         CALL asm_bgc_bkg_alloc
150      ENDIF
151#endif
152
153      ! If creating an averaged assim bkg, sum the contribution every timestep
154      IF ( ln_avgbkg ) THEN
155         IF (lwp) THEN
156              WRITE(numout,*) 'asm_wri_bkg : Summing assim bkg fields at timestep ',kt
157              WRITE(numout,*) '~~~~~~~~~~~~ '
158         ENDIF
159
160         tn_tavg(:,:,:)        = tn_tavg(:,:,:) + tsn(:,:,:,jp_tem) / numtimes_tavg
161         sn_tavg(:,:,:)        = sn_tavg(:,:,:) + tsn(:,:,:,jp_sal) / numtimes_tavg
162         sshn_tavg(:,:)        = sshn_tavg(:,:) + sshn (:,:) / numtimes_tavg
163         un_tavg(:,:,:)        = un_tavg(:,:,:) + un(:,:,:) / numtimes_tavg
164         vn_tavg(:,:,:)        = vn_tavg(:,:,:) + vn(:,:,:) / numtimes_tavg
165         avt_tavg(:,:,:)        = avt_tavg(:,:,:) + avt(:,:,:) / numtimes_tavg
166#if defined key_zdftke
167         en_tavg(:,:,:)       = en_tavg(:,:,:) + en(:,:,:) / numtimes_tavg
168#endif
169#if defined key_top
170         CALL asm_bgc_bkg_tavg( kt, numtimes_tavg )
171#endif
172      ENDIF
173     
174
175      ! Write out background at time step nitbkg_r or nitavgbkg_r
176      IF ( ( .NOT. ln_avgbkg .AND. (kt == nitbkg_r) ) .OR. &
177      &          ( ln_avgbkg .AND. (kt == nitavgbkg_r) ) ) THEN
178         !
179         WRITE(cl_asmbkg, FMT='(A,".nc")' ) TRIM( c_asmbkg )
180         cl_asmbkg = TRIM( cl_asmbkg )
181         INQUIRE( FILE = cl_asmbkg, EXIST = llok )
182         !
183         IF( .NOT. llok ) THEN
184            IF(lwp) WRITE(numout,*) ' Setting up assimilation background file '// TRIM( c_asmbkg )
185            !
186            !                                      ! Define the output file       
187            CALL iom_open( c_asmbkg, inum, ldwrt = .TRUE., kiolib = jprstlib)
188            !
189            !
190            ! Write the information
191            IF ( ln_avgbkg ) THEN
192               IF( nitavgbkg_r == nit000 - 1 ) THEN      ! Treat special case when nitavgbkg = 0
193                  zdate = REAL( ndastp )
194#if defined key_zdftke
195                  ! lk_zdftke=T :   Read turbulent kinetic energy ( en )
196                  IF(lwp) WRITE(numout,*) ' Reading TKE (en) from restart...'
197                  CALL tke_rst( nit000, 'READ' )               ! lk_zdftke=T :   Read turbulent kinetic energy ( en )
198
199#endif
200               ELSE
201                  zdate = REAL( ndastp )
202               ENDIF
203               CALL iom_rstput( kt, nitavgbkg_r, inum, 'rdastp' , zdate   )
204               CALL iom_rstput( kt, nitavgbkg_r, inum, 'un'     , un_tavg )
205               CALL iom_rstput( kt, nitavgbkg_r, inum, 'vn'     , vn_tavg )
206               CALL iom_rstput( kt, nitavgbkg_r, inum, 'tn'     , tn_tavg )
207               CALL iom_rstput( kt, nitavgbkg_r, inum, 'sn'     , sn_tavg )
208               CALL iom_rstput( kt, nitavgbkg_r, inum, 'sshn'   , sshn_tavg)
209#if defined key_zdftke
210               CALL iom_rstput( kt, nitavgbkg_r, inum, 'en'     , en_tavg )
211#endif
212               CALL iom_rstput( kt, nitavgbkg_r, inum, 'avt'    , avt_tavg)
213               !
214            ELSE
215               IF( nitbkg_r == nit000 - 1 ) THEN      ! Treat special case when nitbkg = 0
216                  zdate = REAL( ndastp )
217#if defined key_zdftke
218                  ! lk_zdftke=T :   Read turbulent kinetic energy ( en )
219                  IF(lwp) WRITE(numout,*) ' Reading TKE (en) from restart...'
220                  CALL tke_rst( nit000, 'READ' )               ! lk_zdftke=T :   Read turbulent kinetic energy ( en )
221
222#endif
223               ELSE
224                  zdate = REAL( ndastp )
225               ENDIF
226               CALL iom_rstput( kt, nitbkg_r, inum, 'rdastp' , zdate   )
227               CALL iom_rstput( kt, nitbkg_r, inum, 'un'     , un                )
228               CALL iom_rstput( kt, nitbkg_r, inum, 'vn'     , vn                )
229               CALL iom_rstput( kt, nitbkg_r, inum, 'tn'     , tsn(:,:,:,jp_tem) )
230               CALL iom_rstput( kt, nitbkg_r, inum, 'sn'     , tsn(:,:,:,jp_sal) )
231               CALL iom_rstput( kt, nitbkg_r, inum, 'sshn'   , sshn              )
232#if defined key_zdftke
233               CALL iom_rstput( kt, nitbkg_r, inum, 'en'     , en                )
234#endif
235               CALL iom_rstput( kt, nitbkg_r, inum, 'avt'    , avt               )
236               !
237            ENDIF
238           
239#if defined key_top
240            CALL asm_bgc_bkg_wri( kt, inum, ln_avgbkg )
241#endif
242            CALL iom_close( inum )
243
244         ENDIF
245         !
246      ENDIF
247
248      !                                !-------------------------------------------
249      IF( kt == nitdin_r ) THEN        ! Write out background at time step nitdin_r
250         !                             !-----------------------------------========
251         !
252         WRITE(cl_asmdin, FMT='(A,".nc")' ) TRIM( c_asmdin )
253         cl_asmdin = TRIM( cl_asmdin )
254         INQUIRE( FILE = cl_asmdin, EXIST = llok )
255         !
256         IF( .NOT. llok ) THEN
257            IF(lwp) WRITE(numout,*) ' Setting up assimilation background file '// TRIM( c_asmdin )
258            !
259            !                                      ! Define the output file       
260            CALL iom_open( c_asmdin, inum, ldwrt = .TRUE., kiolib = jprstlib)
261            !
262            IF( nitdin_r == nit000 - 1 ) THEN      ! Treat special case when nitbkg = 0
263
264               zdate = REAL( ndastp )
265            ELSE
266               zdate = REAL( ndastp )
267            ENDIF
268            !
269            !                                      ! Write the information
270            CALL iom_rstput( kt, nitdin_r, inum, 'rdastp' , zdate             )
271            CALL iom_rstput( kt, nitdin_r, inum, 'un'     , un                )
272            CALL iom_rstput( kt, nitdin_r, inum, 'vn'     , vn                )
273            CALL iom_rstput( kt, nitdin_r, inum, 'tn'     , tsn(:,:,:,jp_tem) )
274            CALL iom_rstput( kt, nitdin_r, inum, 'sn'     , tsn(:,:,:,jp_sal) )
275            CALL iom_rstput( kt, nitdin_r, inum, 'sshn'   , sshn              )
276#if defined key_lim2 || defined key_lim3
277            IF(( nn_ice == 2 ) .OR. ( nn_ice == 3 )) THEN
278          IF(ALLOCATED(frld)) THEN
279                  CALL iom_rstput( kt, nitdin_r, inum, 'iceconc', 1.0 - frld(:,:)   )
280               ELSE
281        CALL ctl_warn('Ice concentration not written to background as ice variable frld not allocated on this timestep')
282          ENDIF
283            ENDIF
284#endif
285            !
286            CALL iom_close( inum )
287         ENDIF
288         !
289      ENDIF
290      !                   
291   END SUBROUTINE asm_bkg_wri
292
293   !!======================================================================
294END MODULE asmbkg
Note: See TracBrowser for help on using the repository browser.