source: branches/2014/dev_r4650_UKMO14.12_STAND_ALONE_OBSOPER/NEMOGCM/NEMO/OPA_SRC/SOL/solver.F90 @ 5600

Last change on this file since 5600 was 5600, checked in by andrewryan, 5 years ago

merged in latest version of trunk alongside changes to SAO_SRC to be compatible with latest OBS

  • Property svn:keywords set to Id
File size: 6.9 KB
Line 
1MODULE solver
2   !!======================================================================
3   !!                     ***  MODULE  solver  ***
4   !! Ocean solver :  initialization of ocean solver
5   !!=====================================================================
6   !! History :  OPA  ! 1990-10  (G. Madec)  Original code           
7   !!                 ! 1993-02  (O. Marti)                         
8   !!                 ! 1997-02  (G. Madec)  local depth inverse computation
9   !!                 ! 1998-10  (G. Roullet, G. Madec)  free surface
10   !!   NEMO     1.0  ! 2003-07  (G. Madec)  free form, F90
11   !!            3.2  ! 2009-07  (R. Benshila) suppression of rigid-lid & FETI solver
12   !!----------------------------------------------------------------------
13#if defined key_dynspg_flt   ||   defined key_esopa 
14   !!----------------------------------------------------------------------
15   !!   'key_dynspg_flt'                              filtered free surface
16   !!----------------------------------------------------------------------
17   !!   solver_init: solver initialization
18   !!----------------------------------------------------------------------
19   USE oce             ! ocean dynamics and tracers variables
20   USE dom_oce         ! ocean space and time domain variables
21   USE zdf_oce         ! ocean vertical physics variables
22   USE sol_oce         ! solver variables
23   USE dynspg_oce      ! choice/control of key cpp for surface pressure gradient
24   USE solmat          ! matrix of the solver
25   USE in_out_manager  ! I/O manager
26   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
27   USE lib_mpp         ! MPP library
28   USE timing          ! timing
29
30   IMPLICIT NONE
31
32   !!----------------------------------------------------------------------
33   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
34   !! $Id$
35   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
36   !!----------------------------------------------------------------------
37CONTAINS
38
39   SUBROUTINE solver_init( kt )
40      !!----------------------------------------------------------------------
41      !!                  ***  ROUTINE solver_init  ***
42      !!                   
43      !! ** Purpose :   Initialization of the elliptic solver
44      !!     
45      !! ** Method  :   a solver is required when using the filtered free
46      !!              surface.
47      !!
48      !! ** Action  : - c_solver_pt : nature of the gridpoint at which the solver is applied
49      !!
50      !! References : Jensen, 1986: Adv. Phys. Oceanogr. Num. Mod.,Ed. O Brien,87-110.
51      !!----------------------------------------------------------------------
52      INTEGER, INTENT(in) :: kt
53      INTEGER             ::   ios       ! Local integer output status for namelist read
54      !!
55      NAMELIST/namsol/ nn_solv, nn_sol_arp, nn_nmin, nn_nmax, nn_nmod, rn_eps, rn_resmax, rn_sor
56      !!----------------------------------------------------------------------
57      !
58      IF( nn_timing == 1 )  CALL timing_start('solver_init')
59      !
60
61      IF(lwp) THEN                  !* open elliptic solver statistics file (only on the printing processors)
62         CALL ctl_opn( numsol, 'solver.stat', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea )
63      ENDIF
64
65      REWIND( numnam_ref )              ! Namelist namsol in reference namelist : Elliptic solver / free surface
66      READ  ( numnam_ref, namsol, IOSTAT = ios, ERR = 901)
67901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsol in reference namelist', lwp )
68
69      REWIND( numnam_cfg )              ! Namelist namsol in configuration namelist : Elliptic solver / free surface
70      READ  ( numnam_cfg, namsol, IOSTAT = ios, ERR = 902 )
71902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsol in configuration namelist', lwp )
72      IF(lwm) WRITE ( numond, namsol )
73
74      IF(lwp) THEN                  !* Control print
75         WRITE(numout,*)
76         WRITE(numout,*) 'solver_init : solver to compute the surface pressure gradient'
77         WRITE(numout,*) '~~~~~~~~~~~'
78         WRITE(numout,*) '   Namelist namsol : set solver parameters'
79         WRITE(numout,*) '      type of elliptic solver            nn_solv    = ', nn_solv
80         WRITE(numout,*) '      absolute/relative (0/1) precision  nn_sol_arp = ', nn_sol_arp
81         WRITE(numout,*) '      minimum iterations for solver      nn_nmin    = ', nn_nmin
82         WRITE(numout,*) '      maximum iterations for solver      nn_nmax    = ', nn_nmax
83         WRITE(numout,*) '      frequency for test                 nn_nmod    = ', nn_nmod
84         WRITE(numout,*) '      absolute precision of solver       rn_eps     = ', rn_eps
85         WRITE(numout,*) '      absolute precision for SOR solver  rn_resmax  = ', rn_resmax
86         WRITE(numout,*) '      optimal coefficient of sor         rn_sor     = ', rn_sor
87         WRITE(numout,*)
88      ENDIF
89      eps = rn_eps
90
91      !                              ! allocate solver arrays
92      IF( .NOT. lk_agrif .OR. .NOT. ln_rstart) THEN
93         IF( sol_oce_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'solver_init : unable to allocate sol_oce arrays' )
94         gcx (:,:) = 0.e0
95         gcxb(:,:) = 0.e0
96      ENDIF
97
98      SELECT CASE( nn_solv )          !* parameter check
99      !
100      CASE ( 1 )                          ! preconditioned conjugate gradient solver
101         IF(lwp) WRITE(numout,*) '   a preconditioned conjugate gradient solver is used'
102         IF( jpr2di /= 0 .AND. jpr2dj /= 0 )   CALL ctl_stop( ' jpr2di and jpr2dj should be equal to zero' )
103         !
104      CASE ( 2 )                          ! successive-over-relaxation solver
105         IF(lwp) WRITE(numout,*) '   a successive-over-relaxation solver with extra outer halo is used'
106         IF(lwp) WRITE(numout,*) '   with jpr2di =', jpr2di, ' and  jpr2dj =', jpr2dj
107         IF( .NOT. lk_mpp .AND. jpr2di /= 0 .AND. jpr2dj /= 0 ) THEN
108             CALL ctl_stop( 'jpr2di and jpr2dj are not equal to zero',   &
109             &              'In this case the algorithm should be used only with the key_mpp_... option' )
110         ELSE
111            IF( ( ( jperio == 1 .OR. jperio == 4 .OR. jperio == 6 ) .OR. ( jpni /= 1 ) ) &
112              &  .AND. ( jpr2di /= jpr2dj ) )   CALL ctl_stop( 'jpr2di should be equal to jpr2dj' )
113         ENDIF
114         !
115      CASE DEFAULT                        ! error in parameter
116         WRITE(ctmp1,*) '          bad flag value for nn_solv = ', nn_solv
117         CALL ctl_stop( ctmp1 )
118      END SELECT
119      !
120
121      !                             !* Grid-point at which the solver is applied
122!!gm  c_solver_pt should be removed: nomore bsf, only T-point is used
123      c_solver_pt = 'T'                   ! always T-point (ssh solver only, not anymore bsf)
124
125      CALL sol_mat( kt )            !* Construction of the elliptic system matrix
126      !
127      IF( nn_timing == 1 )  CALL timing_stop('solver_init')
128      !
129   END SUBROUTINE solver_init
130#endif
131
132   !!======================================================================
133END MODULE solver
Note: See TracBrowser for help on using the repository browser.