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.
dynspg.F90 in trunk/NEMO/OPA_SRC/DYN – NEMO

source: trunk/NEMO/OPA_SRC/DYN/dynspg.F90 @ 1241

Last change on this file since 1241 was 1241, checked in by rblod, 15 years ago

Fix a stupid bug for time splitting and ensure restartability for dynspg_ts in addition, see tickets #280 and #292
  • Property svn:eol-style set to native
  • Property svn:executable set to *
  • Property svn:keywords set to Id
File size: 8.5 KB
Line 
1MODULE dynspg
2   !!======================================================================
3   !!                       ***  MODULE  dynspg  ***
4   !! Ocean dynamics:  surface pressure gradient control
5   !!======================================================================
6   !! History :  9.0  !  05-12  (C. Talandier, G. Madec)  Original code
7   !!            9.0  !  05-12  (V. Garnier)  dyn_spg_ctl: Original code
8   !!----------------------------------------------------------------------
9
10   !!----------------------------------------------------------------------
11   !!   dyn_spg     : update the dynamics trend with the lateral diffusion
12   !!   dyn_spg_ctl : initialization, namelist read, and parameters control
13   !!----------------------------------------------------------------------
14   USE oce            ! ocean dynamics and tracers variables
15   USE dom_oce        ! ocean space and time domain variables
16   USE obc_oce        ! ocean open boundary conditions
17   USE dynspg_oce     ! surface pressure gradient variables
18   USE dynspg_exp     ! surface pressure gradient     (dyn_spg_exp routine)
19   USE dynspg_ts      ! surface pressure gradient     (dyn_spg_ts  routine)
20   USE dynspg_flt     ! surface pressure gradient     (dyn_spg_flt routine)
21   USE dynspg_rl      ! surface pressure gradient     (dyn_spg_rl  routine)
22   USE trdmod         ! ocean dynamics trends
23   USE trdmod_oce     ! ocean variables trends
24   USE prtctl         ! Print control                     (prt_ctl routine)
25   USE in_out_manager ! I/O manager
26
27   IMPLICIT NONE
28   PRIVATE
29
30   PUBLIC dyn_spg         ! routine called by step module
31
32   !! * module variables
33   INTEGER ::   nspg = 0   ! type of surface pressure gradient scheme defined from lk_dynspg_...
34
35   !! * Substitutions
36#  include "domzgr_substitute.h90"
37#  include "vectopt_loop_substitute.h90"
38   !!----------------------------------------------------------------------
39   !!   OPA 9.0 , LOCEAN-IPSL (2005)
40   !! $Id$
41   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
42   !!----------------------------------------------------------------------
43
44CONTAINS
45
46   SUBROUTINE dyn_spg( kt, kindic )
47      !!----------------------------------------------------------------------
48      !!                  ***  ROUTINE dyn_spg  ***
49      !!
50      !! ** Purpose :   compute the lateral ocean dynamics physics.
51      !!----------------------------------------------------------------------
52      INTEGER, INTENT( in  ) ::   kt     ! ocean time-step index
53      INTEGER, INTENT( out ) ::   kindic ! solver flag
54      !!
55      REAL(wp) ::   z2dt                      ! temporary scalar
56      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   ztrdu, ztrdv   ! 3D workspace
57      !!----------------------------------------------------------------------
58
59      IF( kt == nit000 )   CALL dyn_spg_ctl      ! initialisation & control of options
60
61      IF( l_trddyn )   THEN                      ! temporary save of ta and sa trends
62         ztrdu(:,:,:) = ua(:,:,:)
63         ztrdv(:,:,:) = va(:,:,:)
64      ENDIF
65
66      SELECT CASE ( nspg )                       ! compute surf. pressure gradient trend and add it to the general trend
67      !                                                     
68      CASE (  0 )   ;   CALL dyn_spg_exp    ( kt )              ! explicit
69      CASE (  1 )   ;   CALL dyn_spg_ts     ( kt )              ! time-splitting
70      CASE (  2 )   ;   CALL dyn_spg_flt    ( kt, kindic )      ! filtered
71      CASE (  3 )   ;   CALL dyn_spg_rl     ( kt, kindic )      ! rigid lid
72      !                                                   
73      CASE ( -1 )                                       ! esopa: test all possibility with control print
74         ;              CALL dyn_spg_exp    ( kt )
75         ;              CALL prt_ctl( tab3d_1=ua, clinfo1=' spg0 - Ua: ', mask1=umask, &
76            &                         tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' )
77         ;              CALL dyn_spg_ts     ( kt )
78         ;              CALL prt_ctl( tab3d_1=ua, clinfo1=' spg1 - Ua: ', mask1=umask, &
79            &                         tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' )
80         ;              CALL dyn_spg_flt  ( kt, kindic )
81         ;              CALL prt_ctl( tab3d_1=ua, clinfo1=' spg2 - Ua: ', mask1=umask, &
82            &                         tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' )
83      END SELECT
84      !                   
85      IF( l_trddyn )   THEN                      ! save the horizontal diffusive trends for further diagnostics
86         SELECT CASE ( nspg )
87         CASE ( 0, 1, 3, 10, 11 )
88            ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:)
89            ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:)
90         CASE( 2, 12 )
91            z2dt = 2. * rdt
92            IF( neuler == 0 .AND. kt == nit000 ) z2dt = rdt
93            ztrdu(:,:,:) = ( ua(:,:,:) - ub(:,:,:) ) / z2dt - ztrdu(:,:,:)
94            ztrdv(:,:,:) = ( va(:,:,:) - vb(:,:,:) ) / z2dt - ztrdv(:,:,:)
95         END SELECT
96         CALL trd_mod( ztrdu, ztrdv, jpdyn_trd_spg, 'DYN', kt )
97      ENDIF
98      !                                          ! print mean trends (used for debugging)
99      IF(ln_ctl)   CALL prt_ctl( tab3d_1=ua, clinfo1=' spg  - Ua: ', mask1=umask, &
100         &                       tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' )
101      !
102   END SUBROUTINE dyn_spg
103
104
105   SUBROUTINE dyn_spg_ctl
106      !!---------------------------------------------------------------------
107      !!                  ***  ROUTINE dyn_spg_ctl  ***
108      !!               
109      !! ** Purpose :   Control the consistency between cpp options for
110      !!      surface pressure gradient schemes
111      !!----------------------------------------------------------------------
112      !! * Local declarations
113      INTEGER ::   ioptio
114      !!----------------------------------------------------------------------
115
116      ! Parameter control and print
117      ! ---------------------------
118      ! Control print
119      IF(lwp) THEN
120         WRITE(numout,*)
121         WRITE(numout,*) 'dyn_spg_ctl : choice of the surface pressure gradient scheme'
122         WRITE(numout,*) '~~~~~~~~~~~'
123         WRITE(numout,*) '     Explicit free surface                  lk_dynspg_exp = ', lk_dynspg_exp
124         WRITE(numout,*) '     Free surface with time splitting       lk_dynspg_ts  = ', lk_dynspg_ts
125         WRITE(numout,*) '     Filtered free surface cst volume       lk_dynspg_flt = ', lk_dynspg_flt
126         WRITE(numout,*) '     Rigid-lid case                         lk_dynspg_rl  = ', lk_dynspg_rl
127      ENDIF
128
129      ! Control of surface pressure gradient scheme options
130      ! ---------------------------------------------------
131      ioptio = 0
132      IF(lk_dynspg_exp)   ioptio = ioptio + 1
133      IF(lk_dynspg_ts )   ioptio = ioptio + 1
134      IF(lk_dynspg_flt)   ioptio = ioptio + 1
135      IF(lk_dynspg_rl )   ioptio = ioptio + 1
136
137      IF( ( ioptio > 1 .AND. .NOT. lk_esopa ) .OR. ioptio == 0 )   &
138           &   CALL ctl_stop( ' Choose only one surface pressure gradient scheme with a key cpp' )
139
140      IF( lk_esopa     )   nspg = -1
141      IF( lk_dynspg_exp)   nspg =  0
142      IF( lk_dynspg_ts )   nspg =  1
143      IF( lk_dynspg_flt)   nspg =  2
144      IF( lk_dynspg_rl )   nspg =  3
145      IF( nspg == 13   )   nspg =  3
146
147      IF( lk_esopa     )   nspg = -1
148
149     IF(lwp) THEN
150         WRITE(numout,*)
151         IF( nspg == -1 )   WRITE(numout,*) '     ESOPA test All scheme used except rigid-lid'
152         IF( nspg ==  0 )   WRITE(numout,*) '     explicit free surface'
153         IF( nspg ==  1 )   WRITE(numout,*) '     free surface with time splitting scheme'
154         IF( nspg ==  2 )   WRITE(numout,*) '     filtered free surface'
155         IF( nspg ==  3 )   WRITE(numout,*) '     rigid-lid'
156      ENDIF
157
158      ! Control of timestep choice
159      ! --------------------------
160      IF( lk_dynspg_ts .OR. lk_dynspg_exp ) THEN
161         IF( n_cla == 1 )   &
162           &   CALL ctl_stop( ' Crossland advection not implemented for this free surface formulation ' )
163      ENDIF
164
165#if defined key_obc
166      ! Conservation of ocean volume (key_dynspg_flt)
167      ! ---------------------------------------------
168      IF( lk_dynspg_flt ) ln_vol_cst = .true.
169
170      ! Application of Flather's algorithm at open boundaries
171      ! -----------------------------------------------------
172      IF( lk_dynspg_flt ) ln_obc_fla = .false.
173      IF( lk_dynspg_exp ) ln_obc_fla = .true.
174      IF( lk_dynspg_ts  ) ln_obc_fla = .true.
175#endif
176
177   END SUBROUTINE dyn_spg_ctl
178
179  !!======================================================================
180END MODULE dynspg
Note: See TracBrowser for help on using the repository browser.