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

source: trunk/NEMO/OPA_SRC/SOL/solfet.F90 @ 3

Last change on this file since 3 was 3, checked in by opalod, 20 years ago

Initial revision

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 6.8 KB
Line 
1MODULE solfet
2   !!======================================================================
3   !!                     ***  MODULE  solfet
4   !! Ocean solver :  Finite Elements Tearing & Interconnecting solver
5   !!=====================================================================
6#if defined key_feti
7   !!----------------------------------------------------------------------
8   !!   'key_feti' :                                            FETI solver
9   !!----------------------------------------------------------------------
10   !!   sol_fet     : FETI solver
11   !!----------------------------------------------------------------------
12   !! * Modules used
13   USE oce             ! ocean dynamics and tracers variables
14   USE dom_oce         ! ocean space and time domain variables
15   USE sol_oce         ! ocean solver
16   USE lib_mpp         ! distribued memory computing
17   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
18
19   IMPLICIT NONE
20   PRIVATE
21
22   !! * Routine accessibility
23   PUBLIC sol_fet              ! ???
24   !!----------------------------------------------------------------------
25
26CONTAINS
27
28   SUBROUTINE sol_fet( kindic )
29      !!---------------------------------------------------------------------
30      !!                  ***  ROUTINE sol_fet  ***
31      !!
32      !! ** Purpose :
33      !!     Solve the ellipic equation for the barotropic stream function
34      !!     system (default option) or the transport divergence system
35      !!     ("key_dynspg_fsc") using a Finite Elements Tearing &
36      !!      Interconnecting (FETI) approach.
37      !!     In the former case, the barotropic stream function trend has a
38      !!     zero boundary condition along all coastlines (i.e. continent
39      !!     as well as islands) while in the latter the boundary condition
40      !!     specification is not required.
41      !!
42      !! ** Method :
43      !!      Resolution of the elliptic equation by a Dual formulation of
44      !!      the Schur Complement Method or Finite Elements Tearing &
45      !!      Interconnecting (FETI) approach
46      !!
47      !! ** Action :
48      !!
49      !! ** References :
50      !!      Guyon, M, Roux, F-X, Chartier, M and Fraunie, P, 1994 :
51      !!      A domain decomposition solver to compute the barotropic
52      !!      component of an OGCM in the parallel processing field.
53      !!      Ocean Modelling, issue 105, december 94.
54      !!
55      !! History :
56      !!        !  97-02  (M. Guyon)  original code
57      !!   8.5  !  02-08  (G. Madec)  F90: Free form
58      !!----------------------------------------------------------------------
59      !! * Arguments
60      INTEGER, INTENT( inout ) ::   kindic   ! solver indicator, < 0 if the conver-
61      !                                      ! gence is not reached: the model is
62      !                                      ! stopped in step
63
64      !! * Local declarations
65      INTEGER ::   ji, jj                    ! dummy loop indices
66      INTEGER ::   iit0, inn, iju 
67      REAL(wp) ::   zgcad, zgwgt
68      !!----------------------------------------------------------------------
69      !!  OPA 8.5, LODYC-IPSL (2002)
70      !!----------------------------------------------------------------------
71
72
73      ! Norme(gcr)   = (gcb, gcb)
74      ! gcdes = gsr
75     
76      ! bmask field is filtering the differents contribution on the
77      ! non-overlapping interface
78     
79      gcb(:,:) = bmask(:,:) * gcb(:,:)
80     
81      ! Mpp: sum over all the global domain
82     
83      ! copy the right hand member
84     
85      CALL feti_vmov(noeuds,gcb(1,1),wfeti(may))
86     
87      ! conservation of descent direction if ntest = 0
88     
89      IF(ntest /= 0) mjj0=0
90      iit0 = mjj0
91      !                      --->
92      !    resolution of the Grad(PS) equation by a Dual formulation of the
93      !    Schur Complement Method or Finite Elements Tearing & Interconnecting
94      !    (FETI) approach
95      !    interface problem (Lagrange  multiplier) : PCPG algorithm
96      !    local problem (trend of the 2D potential field) : LU factorization
97      !    preconditioner : lumped
98      !    optimisation : Krylov initialisation + Krylov correction
99
100      CALL feti_dualschur(noeuds,nifmat+1,njfmat+1,wfeti(maan),   &
101          npblo,wfeti(mablo),ninterf,   &
102          ninterfc,nni,nnic,mfet(mandvois),mfet(mandvoisc),   &
103          mfet(maplistin),mfet(malistin),   &
104          wfeti(mapoids),wfeti(miax),   &
105          wfeti(maz),wfeti(may),   &
106          wfeti(mabitw),wfeti(mautilu),   &
107          wfeti(malambda),wfeti(mag),   &
108          wfeti(mapg),wfeti(mamg),nitmax,nmaxd,mjj0,   &
109          wfeti(mawj),   &
110          wfeti(madwj),wfeti(madwwj),wfeti(magamm),   &
111          wfeti(mawork),   &
112          wfeti(mabufin),wfeti(mabufout),narea,epsilo,   &
113          ndlblo,mfet(malisblo),ndkerep,   &
114          wfeti(maxnul),wfeti(maynul),numit0ete,nitmaxete,   &
115          wfeti(maeteg),wfeti(maeteag),wfeti(maeted),   &
116          wfeti(maetead),wfeti(maeteadd),wfeti(maetegamm),   &
117          wfeti(mansp),   &
118          wfeti(maetev),wfeti(maetew),nnih,mfet(maplistih),   &
119          wfeti(magh),   &
120          wfeti(maw),wfeti(madw),   &
121          res,kindic,inn)
122
123      ! number of iteration of the pcg to solve the interface pb
124
125      inn =  mjj0 - iit0
126
127      ! test of convergence
128      IF( res < epsilo .OR. inn == nmax ) THEN
129          niter = inn
130          ncut  = 999
131      ENDIF
132
133      ! indicator of non-convergence or explosion
134      IF( inn == nmax .OR. rr > 1.e+20 ) kindic = -2
135      IF( ncut == 999 ) GOTO 999
136
137999   CONTINUE
138
139      !  2. Output in gcx
140      !  -----------------
141     
142      CALL feti_vmov( noeuds, wfeti(miax), gcx )
143
144      ! boundary conditions   !!bug ???  check  arguments...
145#   if defined key_dynspg_fsc
146#      if defined key_mpp
147      !   Mpp: export boundary values to neighbouring processors
148      CALL lbc_lnk( gcx, 'S', 1. )
149#      else
150      !   mono- or macro-tasking: W-point, >0, 2D array, no slab
151      IF( nperio /= 0 ) THEN
152         CALL lbc_lnk( gcx, 'T', 1. )
153      ENDIF
154#      endif
155#   else
156#      if defined key_mpp
157      !   Mpp: export boundary values to neighbouring processors
158      CALL lbc_lnk( gcx, 'G', 1. )
159#      else
160      !   mono- or macro-tasking: W-point, >0, 2D array, no slab
161      IF( nperio /= 0 ) THEN
162         CALL lbc_lnk( gcx, 'F', 1. )
163      ENDIF
164#      endif
165#   endif
166
167   END SUBROUTINE sol_fet
168
169#else
170   !!----------------------------------------------------------------------
171   !!   Default option :                                       Empty module
172   !!----------------------------------------------------------------------
173CONTAINS
174   SUBROUTINE sol_fet( kindic )        ! Empty routine
175      INTEGER, INTENT( inout ) ::   kindic  ! solver problem
176      kindic = -100
177   END SUBROUTINE sol_fet
178#endif
179
180   !!=====================================================================
181END MODULE solfet
Note: See TracBrowser for help on using the repository browser.