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

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

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


Revision 1.1 - (show annotations) (download)
Thu Mar 29 15:59:21 2012 UTC (13 years, 3 months ago) by heimbach
Branch: MAIN
Initial check-in of Dan Goldberg's streamice package

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
25 C !INPUT/OUTPUT ARGUMENTS
26 INTEGER myThid
27
28 #ifdef ALLOW_STREAMICE
29
30 C LOCAL VARIABLES
31
32 ! real, dimension(:,:), pointer :: TAUDX, TAUDY, u_prev_iterate, v_prev_iterate, &
33 ! u_bdry_cont, v_bdry_cont, Au, Av, err_u, err_v, &
34 ! geolonq, geolatq, u_last, v_last, float_cond, H_node
35 ! type(ocean_grid_type), pointer :: G
36 ! integer :: conv_flag, i, j, k,l, iter, isym, &
37 ! isdq, iedq, jsdq, jedq, isd, ied, jsd, jed, isumstart, jsumstart, nodefloat, nsub
38 ! real :: err_max, err_tempu, err_tempv, err_init, area, max_vel, tempu, tempv, rhoi, rhow
39
40 INTEGER conv_flag, i, j, k, l, iter, cg_iters, bi, bj
41 INTEGER iter_numconv
42 _RL err_max, err_tempu, err_tempv, err_init, area
43 _RL max_vel, tempu, tempv, err_lastchange, cgtol
44 CHARACTER*(MAX_LEN_MBUF) msgBuf
45 ! _RL taudx (1-OLx:sNx+OLx,1-OLy:sNy+Oly,nSx,nSy)
46 ! _RL taudy (1-OLx:sNx+OLx,1-OLy:sNy+Oly,nSx,nSy)
47
48 CALL STREAMICE_DRIVING_STRESS (myThid)
49
50 cgtol = streamice_cg_tol
51
52 ! CALL WRITE_FULLARRAY_RL ("taudy_SI",taudy_SI,1,0,0,1,0,myThid)
53
54 _EXCH_XY_RL ( taudx_SI , myThid )
55 _EXCH_XY_RL ( taudy_SI , myThid )
56
57 ! CALL WRITE_FULLARRAY_RL ("taudy_SI_2",taudy_SI,1,0,0,1,0,myThid)
58
59 DO bj = myByLo(myThid), myByHi(myThid)
60 DO bi = myBxLo(myThid), myBxHi(myThid)
61 DO j=1-OLy,sNy+OLy
62 DO i=1-OLx,sNx+OLx
63 u_old_SI (i,j,bi,bj) = U_streamice (i,j,bi,bj)
64 v_old_SI (i,j,bi,bj) = V_streamice (i,j,bi,bj)
65 ENDDO
66 ENDDO
67 ENDDO
68 ENDDO
69
70 CALL STREAMICE_VISC_BETA ( myThid )
71
72 _EXCH_XY_RL ( tau_beta_eff_streamice , myThid )
73 _EXCH_XY_RL ( visc_streamice , myThid )
74
75 DO bj = myByLo(myThid), myByHi(myThid)
76 DO bi = myBxLo(myThid), myBxHi(myThid)
77 DO j=1,sNy
78 DO i=1,sNx
79 Au_SI (i,j,bi,bj) = 0. _d 0
80 Av_SI (i,j,bi,bj) = 0. _d 0
81 ubd_SI (i,j,bi,bj) = 0. _d 0
82 vbd_SI (i,j,bi,bj) = 0. _d 0
83 ENDDO
84 ENDDO
85 ENDDO
86 ENDDO
87
88 CALL STREAMICE_CG_BOUND_VALS( myThid,
89 O ubd_SI,
90 O vbd_SI)
91
92 CALL STREAMICE_CG_ACTION( myThid,
93 O Au_SI,
94 O Av_SI,
95 I U_streamice,
96 I V_streamice,
97 I 0, sNx+1, 0, sNy+1 )
98
99 err_init = 0. _d 0
100 err_tempu = 0. _d 0
101 err_tempv = 0. _d 0
102
103 DO bj = myByLo(myThid), myByHi(myThid)
104 DO bi = myBxLo(myThid), myBxHi(myThid)
105 DO j=1,sNy
106 DO i=1,sNx
107 IF (STREAMICE_umask(i,j,bi,bj).eq.1) THEN
108 err_tempu =
109 & ABS (Au_SI(i,j,bi,bj)+ubd_SI(i,j,bi,bj) -
110 & taudx_SI(i,j,bi,bj))
111 ENDIF
112 IF (STREAMICE_vmask(i,j,bi,bj).eq.1) THEN
113 err_tempv = MAX( err_tempu,
114 & ABS (Av_SI(i,j,bi,bj)+vbd_SI(i,j,bi,bj) -
115 & taudy_SI(i,j,bi,bj)))
116 ENDIF
117 IF (err_tempv .ge. err_init) err_init = err_tempv
118 ENDDO
119 ENDDO
120 ENDDO
121 ENDDO
122
123 CALL GLOBAL_MAX_R8 (err_init, myThid)
124
125 iter_numconv = 0
126 err_max = err_init
127 err_lastchange = err_init
128
129 DO iter=1,streamice_max_nl_iter
130
131 C To avoid using "exit", loop goes through all iterations
132 C but after convergence loop does nothing
133
134 IF (err_max .GT. streamice_nonlin_tol * err_init) THEN
135
136 iter_numconv = iter_numconv + 1
137
138 CALL STREAMICE_CG_SOLVE(
139 & U_streamice,
140 & V_streamice,
141 & taudx_SI,
142 & taudy_SI,
143 & myThid , cgtol , cg_iters)
144
145
146
147 WRITE(msgBuf,'(A,I5,A,I4,A)') 'streamice linear solve number',
148 & iter, " ",
149 & cg_iters,
150 & ' iterations '
151 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
152 & SQUEEZE_RIGHT , 1)
153
154 CALL STREAMICE_VISC_BETA ( myThid )
155
156 _EXCH_XY_RL ( tau_beta_eff_streamice , myThid )
157 _EXCH_XY_RL ( visc_streamice , myThid )
158
159 DO bj = myByLo(myThid), myByHi(myThid)
160 DO bi = myBxLo(myThid), myBxHi(myThid)
161 DO j=1,sNy
162 DO i=1,sNx
163 Au_SI (i,j,bi,bj) = 0. _d 0
164 Av_SI (i,j,bi,bj) = 0. _d 0
165 ubd_SI (i,j,bi,bj) = 0. _d 0
166 vbd_SI (i,j,bi,bj) = 0. _d 0
167 ENDDO
168 ENDDO
169 ENDDO
170 ENDDO
171
172 CALL STREAMICE_CG_BOUND_VALS( myThid,
173 O ubd_SI,
174 O vbd_SI)
175
176 CALL STREAMICE_CG_ACTION( myThid,
177 O Au_SI,
178 O Av_SI,
179 I U_streamice,
180 I V_streamice,
181 I 0, sNx+1, 0, sNy+1 )
182
183 err_max = 0. _d 0
184 err_tempu = 0. _d 0
185 err_tempv = 0. _d 0
186
187 DO bj = myByLo(myThid), myByHi(myThid)
188 DO bi = myBxLo(myThid), myBxHi(myThid)
189 DO j=1,sNy
190 DO i=1,sNx
191 IF (STREAMICE_umask(i,j,bi,bj).eq.1) THEN
192 err_tempu =
193 & ABS (Au_SI(i,j,bi,bj)+ubd_SI(i,j,bi,bj) -
194 & taudx_SI(i,j,bi,bj))
195 ENDIF
196 IF (STREAMICE_vmask(i,j,bi,bj).eq.1) THEN
197 err_tempv = MAX( err_tempu,
198 & ABS (Av_SI(i,j,bi,bj)+vbd_SI(i,j,bi,bj) -
199 & taudy_SI(i,j,bi,bj)))
200 ENDIF
201 IF (err_tempv .ge. err_max) err_max = err_tempv
202 ENDDO
203 ENDDO
204 ENDDO
205 ENDDO
206
207 CALL GLOBAL_MAX_R8 (err_max, myThid)
208
209 WRITE(msgBuf,'(A,F11.7)') 'err/err_init',
210 & err_max/err_init
211 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
212 & SQUEEZE_RIGHT , 1)
213
214 IF (err_max<err_lastchange*1.e-2 .and.
215 & STREAMICE_lower_cg_tol) THEN
216 cgtol = cgtol * 5.e-2
217 err_lastchange = err_max
218 WRITE(msgBuf,'(A,F11.7)') 'new cg tol: ',
219 & cgtol
220 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
221 & SQUEEZE_RIGHT , 1)
222 ENDIF
223
224 ENDIF ! (err_max .GT. streamice_nonlin_tol * err_init)
225 ENDDO
226
227 if (iter_numconv .lt. streamice_max_nl_iter) then
228 PRINT *, "VELOCITY SOLVE CONVERGED, ", iter_numconv,
229 & " iterations"
230 else
231 PRINT *, "VELOCITY SOLVE DID NOT CONVERGE IN ",
232 & iter, " iterations"
233 endif
234
235 _EXCH_XY_RL (U_streamice, myThid)
236 _EXCH_XY_RL (V_streamice, myThid)
237
238 #endif
239 RETURN
240 END
241

  ViewVC Help
Powered by ViewVC 1.1.22