/[MITgcm]/MITgcm/pkg/streamice/streamice_vel_solve.F
ViewVC logotype

Contents of /MITgcm/pkg/streamice/streamice_vel_solve.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.10 - (show annotations) (download)
Mon Feb 16 16:46:45 2015 UTC (9 years, 3 months ago) by dgoldberg
Branch: MAIN
CVS Tags: checkpoint65j, checkpoint65k
Changes since 1.9: +4 -12 lines
require parameter streamice_use_petsc to be true even if ALLOW_PETSC cpp command is defined for streamice to use petsc for linear solveD

1 C $Header: /u/gcmpack/MITgcm/pkg/streamice/streamice_vel_solve.F,v 1.9 2014/09/09 23:09:21 jmc Exp $
2 C $Name: $
3
4 #include "STREAMICE_OPTIONS.h"
5 #ifdef ALLOW_AUTODIFF
6 # include "AUTODIFF_OPTIONS.h"
7 #endif
8
9 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
10
11 CBOP
12 SUBROUTINE STREAMICE_VEL_SOLVE( myThid, maxNLIter, maxCGiter,
13 & myiter )
14 C /============================================================\
15 C | SUBROUTINE |
16 C | o |
17 C |============================================================|
18 C | |
19 C \============================================================/
20 IMPLICIT NONE
21
22 C === Global variables ===
23 #include "SIZE.h"
24 #include "EEPARAMS.h"
25 #include "PARAMS.h"
26 #include "STREAMICE.h"
27 #include "STREAMICE_CG.h"
28
29
30 #ifdef ALLOW_AUTODIFF_TAMC
31 # include "tamc.h"
32 #endif
33
34 C !INPUT/OUTPUT ARGUMENTS
35 INTEGER myThid
36 INTEGER maxNLIter
37 INTEGER maxCGIter
38 INTEGER myIter
39
40 #ifdef ALLOW_STREAMICE
41
42 C LOCAL VARIABLES
43
44 ! real, dimension(:,:), pointer :: TAUDX, TAUDY, u_prev_iterate, v_prev_iterate, &
45 ! u_bdry_cont, v_bdry_cont, Au, Av, err_u, err_v, &
46 ! geolonq, geolatq, u_last, v_last, float_cond, H_node
47 ! type(ocean_grid_type), pointer :: G
48 ! integer :: conv_flag, i, j, k,l, iter, isym, &
49 ! isdq, iedq, jsdq, jedq, isd, ied, jsd, jed, isumstart, jsumstart, nodefloat, nsub
50 ! real :: err_max, err_tempu, err_tempv, err_init, area, max_vel, tempu, tempv, rhoi, rhow
51
52 INTEGER i, j, k, l, iter, cg_iters, bi, bj
53 INTEGER iter_numconv
54 INTEGER ikey_nl, myThidTemp
55 _RL err_max, err_tempu, err_tempv, err_init, area, err_max_fp
56 _RL max_vel, tempu, tempv, err_lastchange, cgtol
57 CHARACTER*(MAX_LEN_MBUF) msgBuf
58 LOGICAL CONVERGED
59
60 ! _RL taudx (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
61 ! _RL taudy (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
62
63 IF (STREAMICE_ppm_driving_stress) THEN
64 CALL STREAMICE_DRIVING_STRESS_PPM (myThid)
65 ELSE
66 CALL STREAMICE_DRIVING_STRESS (myThid)
67 ENDIF
68
69 #ifdef ALLOW_AUTODIFF_TAMC
70 #ifdef STREAMICE_STRESS_BOUNDARY_CONTROL
71 !$TAF STORE taudx_SI = comlev1, key=ikey_dynamics
72 !$TAF STORE taudy_SI = comlev1, key=ikey_dynamics
73 #endif
74 #endif
75
76 #ifdef STREAMICE_STRESS_BOUNDARY_CONTROL
77 _EXCH_XY_RL( taudx_SI , myThid )
78 _EXCH_XY_RL( taudy_SI , myThid )
79 CALL STREAMICE_FORCED_BUTTRESS (myThid)
80 #endif
81
82 CALL TIMER_START ('STREAMICE_VEL_SOLVE',myThid)
83
84 cgtol = streamice_cg_tol
85 CONVERGED = .false.
86 err_max = 0.
87 err_max_fp = 0.
88
89 _EXCH_XY_RL( taudx_SI , myThid )
90 _EXCH_XY_RL( taudy_SI , myThid )
91
92 DO bj = myByLo(myThid), myByHi(myThid)
93 DO bi = myBxLo(myThid), myBxHi(myThid)
94 DO j=1-OLy,sNy+OLy
95 DO i=1-OLx,sNx+OLx
96 u_old_SI (i,j,bi,bj) = U_streamice (i,j,bi,bj)
97 v_old_SI (i,j,bi,bj) = V_streamice (i,j,bi,bj)
98 ENDDO
99 ENDDO
100 ENDDO
101 ENDDO
102
103 #ifdef ALLOW_AUTODIFF_TAMC
104 !$TAF STORE U_streamice = comlev1, key=ikey_dynamics
105 !$TAF STORE V_streamice = comlev1, key=ikey_dynamics
106 #ifdef STREAMICE_HYBRID_STRESS
107 !$TAF STORE streamice_taubx = comlev1, key=ikey_dynamics
108 !$TAF STORE streamice_tauby = comlev1, key=ikey_dynamics
109 !$TAF STORE visc_streamice_full = comlev1, key=ikey_dynamics
110 #endif
111 #endif
112
113 #ifdef STREAMICE_HYBRID_STRESS
114 CALL STREAMICE_VISC_BETA_HYBRID ( myThid )
115 #else
116 CALL STREAMICE_VISC_BETA ( myThid )
117 #endif
118
119 #ifdef STREAMICE_HYBRID_STRESS
120 !$TAF STORE visc_streamice_full = comlev1, key=ikey_dynamics
121 #endif
122
123 _EXCH_XY_RL( tau_beta_eff_streamice , myThid )
124 _EXCH_XY_RL( visc_streamice , myThid )
125
126 if (STREAMICE_chkresidconvergence .or.
127 & (streamice_maxnliter_cpl.eq.0 .and. myIter.eq.0)) then
128
129 CALL STREAMICE_GET_VEL_RESID_ERR ( err_init, myThid )
130
131 #ifdef ALLOW_AUTODIFF_TAMC
132 !$TAF STORE err_init = comlev1, key=ikey_dynamics
133 #endif
134
135 WRITE(msgBuf,'(A,E15.7)') 'initial nonlinear resid (error): ',
136 & err_init
137 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
138 & SQUEEZE_RIGHT , 1)
139
140 endif !STREAMICE_chkresidconvergence
141
142 iter_numconv = 0
143 err_max = err_init
144 err_max_fp = streamice_nonlin_tol_fp * 10.
145 err_lastchange = err_init
146
147 C START NL ITER. LOOP
148 C-------------------------------------------------------------------
149
150 DO iter=1,maxNLIter
151
152 C To avoid using "exit", loop goes through all iterations
153 C but after convergence loop does nothing
154
155 #ifdef ALLOW_AUTODIFF_TAMC
156 ikey_nl = (ikey_dynamics-1)*streamice_max_nl + iter
157 #endif
158 #ifdef ALLOW_AUTODIFF_TAMC
159 !$TAF STORE cgtol = comlev1_stream_nl, key=ikey_nl
160 !$TAF STORE err_max = comlev1_stream_nl, key=ikey_nl
161 !$TAF STORE err_max_fp = comlev1_stream_nl, key=ikey_nl
162 !$TAF STORE err_tempu = comlev1_stream_nl, key=ikey_nl
163 !$TAF STORE err_tempv = comlev1_stream_nl, key=ikey_nl
164 !$TAF STORE err_lastchange = comlev1_stream_nl, key=ikey_nl
165 !$TAF STORE ru_old_si = comlev1_stream_nl, key=ikey_nl
166 !$TAF STORE rv_old_si = comlev1_stream_nl, key=ikey_nl
167 !$TAF STORE u_old_si = comlev1_stream_nl, key=ikey_nl
168 !$TAF STORE v_old_si = comlev1_stream_nl, key=ikey_nl
169 !$TAF STORE streamice_cg_a1 = comlev1_stream_nl, key=ikey_nl
170 !$TAF STORE streamice_cg_a2 = comlev1_stream_nl, key=ikey_nl
171 !$TAF STORE streamice_cg_a3 = comlev1_stream_nl, key=ikey_nl
172 !$TAF STORE streamice_cg_a4 = comlev1_stream_nl, key=ikey_nl
173 !$TAF STORE U_streamice = comlev1_stream_nl, key=ikey_nl
174 !$TAF STORE V_streamice = comlev1_stream_nl, key=ikey_nl
175 !$TAF STORE tau_beta_eff_streamice = comlev1_stream_nl, key=ikey_nl
176 !$TAF STORE visc_streamice = comlev1_stream_nl, key=ikey_nl
177 !$TAF STORE zu_old_si = comlev1_stream_nl, key=ikey_nl
178 !$TAF STORE zv_old_si = comlev1_stream_nl, key=ikey_nl
179 #endif
180
181 ! IF ((err_max .GT. streamice_nonlin_tol * err_init) .and.
182 ! & (err_max_fp .GT. streamice_nonlin_tol_fp)) THEN
183
184 #ifdef ALLOW_AUTODIFF_TAMC
185 !$TAF STORE CONVERGED = comlev1_stream_nl, key=ikey_nl
186 #endif
187
188 IF (.not.CONVERGED) THEN
189
190 iter_numconv = iter_numconv + 1
191
192 #ifdef ALLOW_AUTODIFF_TAMC
193 !$TAF STORE cgtol = comlev1_stream_nl, key=ikey_nl
194 !$TAF STORE err_max = comlev1_stream_nl, key=ikey_nl
195 !$TAF STORE err_tempu = comlev1_stream_nl, key=ikey_nl
196 !$TAF STORE err_tempv = comlev1_stream_nl, key=ikey_nl
197 !$TAF STORE err_lastchange = comlev1_stream_nl, key=ikey_nl
198 !$TAF STORE ru_old_si = comlev1_stream_nl, key=ikey_nl
199 !$TAF STORE rv_old_si = comlev1_stream_nl, key=ikey_nl
200 !$TAF STORE streamice_cg_a1 = comlev1_stream_nl, key=ikey_nl
201 !$TAF STORE streamice_cg_a2 = comlev1_stream_nl, key=ikey_nl
202 !$TAF STORE streamice_cg_a3 = comlev1_stream_nl, key=ikey_nl
203 !$TAF STORE streamice_cg_a4 = comlev1_stream_nl, key=ikey_nl
204 !$TAF STORE U_streamice = comlev1_stream_nl, key=ikey_nl
205 !$TAF STORE V_streamice = comlev1_stream_nl, key=ikey_nl
206 !$TAF STORE tau_beta_eff_streamice = comlev1_stream_nl, key=ikey_nl
207 !$TAF STORE visc_streamice = comlev1_stream_nl, key=ikey_nl
208 !$TAF STORE zu_old_si = comlev1_stream_nl, key=ikey_nl
209 !$TAF STORE zv_old_si = comlev1_stream_nl, key=ikey_nl
210 #endif
211
212 CALL STREAMICE_CG_WRAPPER(
213 & U_streamice,
214 & V_streamice,
215 & taudx_SI,
216 & taudy_SI,
217 & cgtol,
218 & cg_iters,
219 & maxCGIter,
220 & myThid )
221
222 #ifdef STREAMICE_HYBRID_STRESS
223 #ifdef ALLOW_AUTODIFF_TAMC
224 !$TAF STORE U_streamice = comlev1_stream_nl, key=ikey_nl
225 !$TAF STORE V_streamice = comlev1_stream_nl, key=ikey_nl
226 #endif
227 #endif
228
229 #ifdef STREAMICE_HYBRID_STRESS
230 CALL STREAMICE_TAUB (myThid)
231 #endif
232
233 !-----------------------------------------------------------------------------
234 WRITE(msgBuf,'(A,I5,A,I4,A)') 'streamice linear solve number',
235 & iter, " ",
236 & cg_iters,
237 & ' iterations '
238 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
239 & SQUEEZE_RIGHT , 1)
240 !-----------------------------------------------------------------------------
241
242 #ifdef ALLOW_AUTODIFF_TAMC
243 !$TAF STORE U_streamice = comlev1_stream_nl, key=ikey_nl
244 !$TAF STORE V_streamice = comlev1_stream_nl, key=ikey_nl
245 #ifdef STREAMICE_HYBRID_STRESS
246 !$TAF STORE streamice_taubx = comlev1_stream_nl, key=ikey_nl
247 !$TAF STORE streamice_tauby = comlev1_stream_nl, key=ikey_nl
248 !$TAF STORE visc_streamice_full = comlev1_stream_nl, key=ikey_nl
249 #endif
250 #endif
251
252 #ifdef STREAMICE_HYBRID_STRESS
253 CALL STREAMICE_VISC_BETA_HYBRID ( myThid )
254 #else
255 CALL STREAMICE_VISC_BETA ( myThid )
256 #endif
257
258 _EXCH_XY_RL( tau_beta_eff_streamice , myThid )
259 _EXCH_XY_RL( visc_streamice , myThid )
260
261 !----------------------CONVERGENCE TESTS-------------------------------
262
263 if (STREAMICE_chkresidconvergence .or.
264 & (streamice_maxnliter_cpl.eq.0 .and. myIter.eq.0)) then
265
266 CALL STREAMICE_GET_VEL_RESID_ERR ( err_max, myThid )
267
268 WRITE(msgBuf,'(A,E15.7)') 'err/err_init',
269 & err_max/err_init
270 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
271 & SQUEEZE_RIGHT , 1)
272
273 IF (err_max .LE. streamice_nonlin_tol * err_init) THEN
274 CONVERGED = .true.
275 ENDIF
276
277 endif
278
279 if (STREAMICE_chkfixedptconvergence .or.
280 & (streamice_maxnliter_cpl.eq.0 .and. myIter.eq.0)) then
281
282 CALL STREAMICE_GET_VEL_FP_ERR ( err_max_fp, myThid )
283
284 DO bj = myByLo(myThid), myByHi(myThid)
285 DO bi = myBxLo(myThid), myBxHi(myThid)
286 DO j=1-OLy,sNy+OLy
287 DO i=1-OLx,sNx+OLx
288 u_old_SI (i,j,bi,bj) = U_streamice (i,j,bi,bj)
289 v_old_SI (i,j,bi,bj) = V_streamice (i,j,bi,bj)
290 ENDDO
291 ENDDO
292 ENDDO
293 ENDDO
294
295 WRITE(msgBuf,'(A,1PE22.14)') 'STREAMICE_FP_ERROR =',
296 & err_max_fp
297 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
298 & SQUEEZE_RIGHT , 1)
299
300 IF (err_max_fp .LE. streamice_nonlin_tol_fp) THEN
301 CONVERGED = .true.
302 ENDIF
303
304 endif
305
306 !----------------------END CONVERGENCE TESTS-------------------------------
307
308 IF (err_max<err_lastchange*1.e-2 .and.
309 & STREAMICE_lower_cg_tol) THEN
310 cgtol = cgtol * 5.e-2
311 err_lastchange = err_max
312 WRITE(msgBuf,'(A,E15.7)') 'new cg tol: ',
313 & cgtol
314 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
315 & SQUEEZE_RIGHT , 1)
316 ENDIF
317
318 ENDIF ! (err_max .GT. streamice_nonlin_tol * err_init)
319 ENDDO
320
321
322
323 C END NL ITER. LOOP
324 C-------------------------------------------------------------------
325
326 if (iter_numconv .lt. streamice_max_nl_iter) then
327 WRITE(msgBuf,'(A,I5,A)') 'VELOCITY SOLVE CONVERGED, ',
328 & iter_numconv, ' iterations'
329 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
330 & SQUEEZE_RIGHT , 1)
331 else
332 WRITE(msgBuf,'(A,I5,A)') 'VELOCITY SOLVE NOT CONVERGED IN ',
333 & iter_numconv, ' iterations'
334 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
335 & SQUEEZE_RIGHT , 1)
336 endif
337
338 _EXCH_XY_RL(U_streamice, myThid)
339 _EXCH_XY_RL(V_streamice, myThid)
340
341 CALL TIMER_STOP ('STREAMICE_VEL_SOLVE',myThid)
342
343 #endif
344 RETURN
345 END

  ViewVC Help
Powered by ViewVC 1.1.22