/[MITgcm]/MITgcm_contrib/dgoldberg/streamice/streamice_vel_solve.F
ViewVC logotype

Annotation of /MITgcm_contrib/dgoldberg/streamice/streamice_vel_solve.F

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


Revision 1.14 - (hide annotations) (download)
Wed Aug 27 19:29:14 2014 UTC (10 years, 11 months ago) by dgoldberg
Branch: MAIN
CVS Tags: HEAD
Changes since 1.13: +101 -218 lines
updating contrib streamice repo with latest files, and separated out convergence checks; and parameterised maximum iteration counts and interface w shelfice for coupling

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

  ViewVC Help
Powered by ViewVC 1.1.22