/[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.2 - (hide annotations) (download)
Wed May 2 02:36:01 2012 UTC (13 years, 2 months ago) by heimbach
Branch: MAIN
Changes since 1.1: +102 -5 lines
Various updates to streamice code

1 heimbach 1.1 C $Header: /u/gcmpack/MITgcm/pkg/streamice/streamice_init_varia.F,v 1.6 2011/06/29 16:24:10 dng Exp $
2     C $Name: $
3    
4     #include "STREAMICE_OPTIONS.h"
5    
6     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
7    
8     CBOP
9     SUBROUTINE STREAMICE_VEL_SOLVE( myThid )
10     C /============================================================\
11     C | SUBROUTINE |
12     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 heimbach 1.2 #ifdef ALLOW_AUTODIFF_TAMC
25     # include "tamc.h"
26     #endif
27 heimbach 1.1
28     C !INPUT/OUTPUT ARGUMENTS
29     INTEGER myThid
30    
31     #ifdef ALLOW_STREAMICE
32    
33     C LOCAL VARIABLES
34    
35     ! real, dimension(:,:), pointer :: TAUDX, TAUDY, u_prev_iterate, v_prev_iterate, &
36     ! u_bdry_cont, v_bdry_cont, Au, Av, err_u, err_v, &
37     ! geolonq, geolatq, u_last, v_last, float_cond, H_node
38     ! type(ocean_grid_type), pointer :: G
39     ! integer :: conv_flag, i, j, k,l, iter, isym, &
40     ! isdq, iedq, jsdq, jedq, isd, ied, jsd, jed, isumstart, jsumstart, nodefloat, nsub
41     ! real :: err_max, err_tempu, err_tempv, err_init, area, max_vel, tempu, tempv, rhoi, rhow
42    
43     INTEGER conv_flag, i, j, k, l, iter, cg_iters, bi, bj
44     INTEGER iter_numconv
45 heimbach 1.2 INTEGER ikey_nl
46 heimbach 1.1 _RL err_max, err_tempu, err_tempv, err_init, area
47     _RL max_vel, tempu, tempv, err_lastchange, cgtol
48     CHARACTER*(MAX_LEN_MBUF) msgBuf
49     ! _RL taudx (1-OLx:sNx+OLx,1-OLy:sNy+Oly,nSx,nSy)
50     ! _RL taudy (1-OLx:sNx+OLx,1-OLy:sNy+Oly,nSx,nSy)
51    
52     CALL STREAMICE_DRIVING_STRESS (myThid)
53    
54     cgtol = streamice_cg_tol
55    
56     ! CALL WRITE_FULLARRAY_RL ("taudy_SI",taudy_SI,1,0,0,1,0,myThid)
57    
58     _EXCH_XY_RL ( taudx_SI , myThid )
59     _EXCH_XY_RL ( taudy_SI , myThid )
60    
61     ! CALL WRITE_FULLARRAY_RL ("taudy_SI_2",taudy_SI,1,0,0,1,0,myThid)
62    
63     DO bj = myByLo(myThid), myByHi(myThid)
64     DO bi = myBxLo(myThid), myBxHi(myThid)
65     DO j=1-OLy,sNy+OLy
66     DO i=1-OLx,sNx+OLx
67     u_old_SI (i,j,bi,bj) = U_streamice (i,j,bi,bj)
68     v_old_SI (i,j,bi,bj) = V_streamice (i,j,bi,bj)
69     ENDDO
70     ENDDO
71     ENDDO
72     ENDDO
73    
74 heimbach 1.2 #ifdef ALLOW_AUTODIFF_TAMC
75     !$TAF STORE U_streamice = comlev1, key=ikey_dynamics
76     !$TAF STORE V_streamice = comlev1, key=ikey_dynamics
77     #endif
78    
79 heimbach 1.1 CALL STREAMICE_VISC_BETA ( myThid )
80    
81     _EXCH_XY_RL ( tau_beta_eff_streamice , myThid )
82     _EXCH_XY_RL ( visc_streamice , myThid )
83    
84     DO bj = myByLo(myThid), myByHi(myThid)
85     DO bi = myBxLo(myThid), myBxHi(myThid)
86     DO j=1,sNy
87     DO i=1,sNx
88     Au_SI (i,j,bi,bj) = 0. _d 0
89     Av_SI (i,j,bi,bj) = 0. _d 0
90     ubd_SI (i,j,bi,bj) = 0. _d 0
91     vbd_SI (i,j,bi,bj) = 0. _d 0
92     ENDDO
93     ENDDO
94     ENDDO
95     ENDDO
96    
97     CALL STREAMICE_CG_BOUND_VALS( myThid,
98     O ubd_SI,
99     O vbd_SI)
100    
101 heimbach 1.2 #ifdef ALLOW_AUTODIFF_TAMC
102     !$TAF STORE U_streamice = comlev1, key=ikey_dynamics
103     !$TAF STORE V_streamice = comlev1, key=ikey_dynamics
104     #endif
105    
106 heimbach 1.1 CALL STREAMICE_CG_ACTION( myThid,
107     O Au_SI,
108     O Av_SI,
109     I U_streamice,
110     I V_streamice,
111     I 0, sNx+1, 0, sNy+1 )
112    
113     err_init = 0. _d 0
114    
115     DO bj = myByLo(myThid), myByHi(myThid)
116     DO bi = myBxLo(myThid), myBxHi(myThid)
117 heimbach 1.2 #ifdef ALLOW_AUTODIFF_TAMC
118     !$TAF STORE err_init = comlev1, key=ikey_dynamics
119     #endif
120 heimbach 1.1 DO j=1,sNy
121     DO i=1,sNx
122 heimbach 1.2 err_tempu = 0. _d 0
123     err_tempv = 0. _d 0
124 heimbach 1.1 IF (STREAMICE_umask(i,j,bi,bj).eq.1) THEN
125     err_tempu =
126     & ABS (Au_SI(i,j,bi,bj)+ubd_SI(i,j,bi,bj) -
127     & taudx_SI(i,j,bi,bj))
128     ENDIF
129     IF (STREAMICE_vmask(i,j,bi,bj).eq.1) THEN
130     err_tempv = MAX( err_tempu,
131     & ABS (Av_SI(i,j,bi,bj)+vbd_SI(i,j,bi,bj) -
132     & taudy_SI(i,j,bi,bj)))
133     ENDIF
134     IF (err_tempv .ge. err_init) err_init = err_tempv
135     ENDDO
136     ENDDO
137     ENDDO
138     ENDDO
139 heimbach 1.2 #ifdef ALLOW_AUTODIFF_TAMC
140     !$TAF STORE err_init = comlev1, key=ikey_dynamics
141     #endif
142 heimbach 1.1
143     CALL GLOBAL_MAX_R8 (err_init, myThid)
144    
145     iter_numconv = 0
146     err_max = err_init
147     err_lastchange = err_init
148    
149 heimbach 1.2 C START NL ITER. LOOP
150     C-------------------------------------------------------------------
151    
152 heimbach 1.1 DO iter=1,streamice_max_nl_iter
153    
154     C To avoid using "exit", loop goes through all iterations
155     C but after convergence loop does nothing
156    
157 heimbach 1.2 #ifdef ALLOW_AUTODIFF_TAMC
158     ikey_nl = (ikey_dynamics-1)*streamice_max_nl + iter
159     #endif
160     #ifdef ALLOW_AUTODIFF_TAMC
161     !$TAF STORE cgtol = comlev1_stream_nl, key=ikey_nl
162     !$TAF STORE err_max = comlev1_stream_nl, key=ikey_nl
163     !$TAF STORE err_tempu = comlev1_stream_nl, key=ikey_nl
164     !$TAF STORE err_tempv = comlev1_stream_nl, key=ikey_nl
165     !$TAF STORE err_lastchange = comlev1_stream_nl, key=ikey_nl
166     !$TAF STORE ru_old_si = comlev1_stream_nl, key=ikey_nl
167     !$TAF STORE rv_old_si = comlev1_stream_nl, key=ikey_nl
168     !$TAF STORE streamice_cg_a1 = comlev1_stream_nl, key=ikey_nl
169     !$TAF STORE streamice_cg_a2 = comlev1_stream_nl, key=ikey_nl
170     !$TAF STORE streamice_cg_a3 = comlev1_stream_nl, key=ikey_nl
171     !$TAF STORE streamice_cg_a4 = comlev1_stream_nl, key=ikey_nl
172     !$TAF STORE U_streamice = comlev1_stream_nl, key=ikey_nl
173     !$TAF STORE V_streamice = comlev1_stream_nl, key=ikey_nl
174     !$TAF STORE tau_beta_eff_streamice = comlev1_stream_nl, key=ikey_nl
175     !$TAF STORE visc_streamice = comlev1_stream_nl, key=ikey_nl
176     !$TAF STORE zu_old_si = comlev1_stream_nl, key=ikey_nl
177     !$TAF STORE zv_old_si = comlev1_stream_nl, key=ikey_nl
178     #endif
179    
180 heimbach 1.1 IF (err_max .GT. streamice_nonlin_tol * err_init) THEN
181    
182     iter_numconv = iter_numconv + 1
183    
184 heimbach 1.2 #ifdef ALLOW_AUTODIFF_TAMC
185     !$TAF STORE cgtol = comlev1_stream_nl, key=ikey_nl
186     !$TAF STORE err_max = comlev1_stream_nl, key=ikey_nl
187     !$TAF STORE err_tempu = comlev1_stream_nl, key=ikey_nl
188     !$TAF STORE err_tempv = comlev1_stream_nl, key=ikey_nl
189     !$TAF STORE err_lastchange = comlev1_stream_nl, key=ikey_nl
190     !$TAF STORE ru_old_si = comlev1_stream_nl, key=ikey_nl
191     !$TAF STORE rv_old_si = comlev1_stream_nl, key=ikey_nl
192     !$TAF STORE streamice_cg_a1 = comlev1_stream_nl, key=ikey_nl
193     !$TAF STORE streamice_cg_a2 = comlev1_stream_nl, key=ikey_nl
194     !$TAF STORE streamice_cg_a3 = comlev1_stream_nl, key=ikey_nl
195     !$TAF STORE streamice_cg_a4 = comlev1_stream_nl, key=ikey_nl
196     !$TAF STORE U_streamice = comlev1_stream_nl, key=ikey_nl
197     !$TAF STORE V_streamice = comlev1_stream_nl, key=ikey_nl
198     !$TAF STORE tau_beta_eff_streamice = comlev1_stream_nl, key=ikey_nl
199     !$TAF STORE visc_streamice = comlev1_stream_nl, key=ikey_nl
200     !$TAF STORE zu_old_si = comlev1_stream_nl, key=ikey_nl
201     !$TAF STORE zv_old_si = comlev1_stream_nl, key=ikey_nl
202     #endif
203    
204 heimbach 1.1 CALL STREAMICE_CG_SOLVE(
205     & U_streamice,
206     & V_streamice,
207     & taudx_SI,
208     & taudy_SI,
209 heimbach 1.2 & cgtol,
210     & cg_iters,
211     & myThid )
212 heimbach 1.1
213    
214    
215     WRITE(msgBuf,'(A,I5,A,I4,A)') 'streamice linear solve number',
216     & iter, " ",
217     & cg_iters,
218     & ' iterations '
219     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
220     & SQUEEZE_RIGHT , 1)
221    
222 heimbach 1.2 #ifdef ALLOW_AUTODIFF_TAMC
223     !$TAF STORE U_streamice = comlev1_stream_nl, key=ikey_nl
224     !$TAF STORE V_streamice = comlev1_stream_nl, key=ikey_nl
225     #endif
226    
227 heimbach 1.1 CALL STREAMICE_VISC_BETA ( myThid )
228    
229 heimbach 1.2 #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    
234 heimbach 1.1 _EXCH_XY_RL ( tau_beta_eff_streamice , myThid )
235     _EXCH_XY_RL ( visc_streamice , myThid )
236    
237     DO bj = myByLo(myThid), myByHi(myThid)
238     DO bi = myBxLo(myThid), myBxHi(myThid)
239     DO j=1,sNy
240     DO i=1,sNx
241     Au_SI (i,j,bi,bj) = 0. _d 0
242     Av_SI (i,j,bi,bj) = 0. _d 0
243     ubd_SI (i,j,bi,bj) = 0. _d 0
244     vbd_SI (i,j,bi,bj) = 0. _d 0
245     ENDDO
246     ENDDO
247     ENDDO
248     ENDDO
249    
250     CALL STREAMICE_CG_BOUND_VALS( myThid,
251     O ubd_SI,
252     O vbd_SI)
253    
254 heimbach 1.2 #ifdef ALLOW_AUTODIFF_TAMC
255     !$TAF STORE U_streamice = comlev1_stream_nl, key=ikey_nl
256     !$TAF STORE V_streamice = comlev1_stream_nl, key=ikey_nl
257     #endif
258    
259 heimbach 1.1 CALL STREAMICE_CG_ACTION( myThid,
260     O Au_SI,
261     O Av_SI,
262     I U_streamice,
263     I V_streamice,
264     I 0, sNx+1, 0, sNy+1 )
265    
266 heimbach 1.2 #ifdef ALLOW_AUTODIFF_TAMC
267     !$TAF STORE U_streamice = comlev1_stream_nl, key=ikey_nl
268     !$TAF STORE V_streamice = comlev1_stream_nl, key=ikey_nl
269     #endif
270    
271 heimbach 1.1 err_max = 0. _d 0
272    
273 heimbach 1.2 #ifdef ALLOW_AUTODIFF_TAMC
274     !$TAF STORE err_max = comlev1_stream_nl, key=ikey_nl
275     #endif
276 heimbach 1.1 DO bj = myByLo(myThid), myByHi(myThid)
277     DO bi = myBxLo(myThid), myBxHi(myThid)
278 heimbach 1.2 #ifdef ALLOW_AUTODIFF_TAMC
279     !$TAF STORE err_max = comlev1_stream_nl, key=ikey_nl
280     #endif
281 heimbach 1.1 DO j=1,sNy
282     DO i=1,sNx
283 heimbach 1.2 err_tempu = 0. _d 0
284     err_tempv = 0. _d 0
285 heimbach 1.1 IF (STREAMICE_umask(i,j,bi,bj).eq.1) THEN
286     err_tempu =
287     & ABS (Au_SI(i,j,bi,bj)+ubd_SI(i,j,bi,bj) -
288     & taudx_SI(i,j,bi,bj))
289     ENDIF
290     IF (STREAMICE_vmask(i,j,bi,bj).eq.1) THEN
291     err_tempv = MAX( err_tempu,
292     & ABS (Av_SI(i,j,bi,bj)+vbd_SI(i,j,bi,bj) -
293     & taudy_SI(i,j,bi,bj)))
294     ENDIF
295     IF (err_tempv .ge. err_max) err_max = err_tempv
296     ENDDO
297     ENDDO
298     ENDDO
299     ENDDO
300    
301     CALL GLOBAL_MAX_R8 (err_max, myThid)
302    
303     WRITE(msgBuf,'(A,F11.7)') 'err/err_init',
304     & err_max/err_init
305     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
306     & SQUEEZE_RIGHT , 1)
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,F11.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 heimbach 1.2 C END NL ITER. LOOP
322     C-------------------------------------------------------------------
323    
324 heimbach 1.1 if (iter_numconv .lt. streamice_max_nl_iter) then
325     PRINT *, "VELOCITY SOLVE CONVERGED, ", iter_numconv,
326     & " iterations"
327     else
328     PRINT *, "VELOCITY SOLVE DID NOT CONVERGE IN ",
329     & iter, " iterations"
330     endif
331    
332     _EXCH_XY_RL (U_streamice, myThid)
333     _EXCH_XY_RL (V_streamice, myThid)
334    
335     #endif
336     RETURN
337     END
338    

  ViewVC Help
Powered by ViewVC 1.1.22