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

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

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

revision 1.4 by dgoldberg, Thu Jul 26 02:14:25 2012 UTC revision 1.11 by dgoldberg, Wed Jan 9 21:56:18 2013 UTC
# Line 43  C     LOCAL VARIABLES Line 43  C     LOCAL VARIABLES
43        INTEGER conv_flag, i, j, k, l, iter, cg_iters, bi, bj        INTEGER conv_flag, i, j, k, l, iter, cg_iters, bi, bj
44        INTEGER iter_numconv        INTEGER iter_numconv
45        INTEGER ikey_nl        INTEGER ikey_nl
46        _RL err_max, err_tempu, err_tempv, err_init, area        _RL err_max, err_tempu, err_tempv, err_init, area, err_max_fp
47        _RL max_vel, tempu, tempv, err_lastchange, cgtol        _RL max_vel, tempu, tempv, err_lastchange, cgtol
48        CHARACTER*(MAX_LEN_MBUF) msgBuf        CHARACTER*(MAX_LEN_MBUF) msgBuf
49  !       _RL taudx (1-OLx:sNx+OLx,1-OLy:sNy+Oly,nSx,nSy)  !       _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)  !       _RL taudy (1-OLx:sNx+OLx,1-OLy:sNy+Oly,nSx,nSy)
51    
52        CALL STREAMICE_DRIVING_STRESS (myThid)        IF (STREAMICE_ppm_driving_stress) THEN
53           CALL STREAMICE_DRIVING_STRESS_PPM (myThid)
54          ELSE
55           CALL STREAMICE_DRIVING_STRESS (myThid)
56          ENDIF
57    
58        cgtol = streamice_cg_tol        cgtol = streamice_cg_tol
59    
# Line 74  C     LOCAL VARIABLES Line 78  C     LOCAL VARIABLES
78  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
79  !$TAF STORE U_streamice = comlev1, key=ikey_dynamics  !$TAF STORE U_streamice = comlev1, key=ikey_dynamics
80  !$TAF STORE V_streamice = comlev1, key=ikey_dynamics  !$TAF STORE V_streamice = comlev1, key=ikey_dynamics
81    #ifdef STREAMICE_HYBRID_STRESS
82    !$TAF STORE streamice_taubx = comlev1, key=ikey_dynamics
83    !$TAF STORE streamice_tauby = comlev1, key=ikey_dynamics
84    !$TAF STORE visc_streamice_full = comlev1, key=ikey_dynamics
85    #endif
86  #endif  #endif
87    
88        CALL STREAMICE_VISC_BETA ( myThid )  #ifdef STREAMICE_HYBRID_STRESS
89           CALL STREAMICE_VISC_BETA_HYBRID ( myThid )
90    #else
91           CALL STREAMICE_VISC_BETA ( myThid )
92    #endif
93    
94    #ifdef STREAMICE_HYBRID_STRESS
95    !$TAF STORE visc_streamice_full = comlev1, key=ikey_dynamics
96    #endif
97    
98        _EXCH_XY_RL ( tau_beta_eff_streamice , myThid )        _EXCH_XY_RL ( tau_beta_eff_streamice , myThid )
99        _EXCH_XY_RL ( visc_streamice , myThid )        _EXCH_XY_RL ( visc_streamice , myThid )
# Line 99  C     LOCAL VARIABLES Line 116  C     LOCAL VARIABLES
116       O    ubd_SI,       O    ubd_SI,
117       O    vbd_SI)       O    vbd_SI)
118    
119    !      CALL WRITE_FLD_XY_RL("u_bound_cont","",ubd_SI,0,myThid)
120    !      CALL WRITE_FLD_XY_RL("v_bound_cont","",vbd_SI,0,myThid)
121    !      CALL WRITE_FLD_XY_RL("taudx_u","",taudx_SI,0,myThid)
122    !      CALL WRITE_FLD_XY_RL("taudx_v","",taudy_SI,0,myThid)
123    
124  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
125  !$TAF STORE U_streamice = comlev1, key=ikey_dynamics  !$TAF STORE U_streamice = comlev1, key=ikey_dynamics
126  !$TAF STORE V_streamice = comlev1, key=ikey_dynamics  !$TAF STORE V_streamice = comlev1, key=ikey_dynamics
# Line 111  C     LOCAL VARIABLES Line 133  C     LOCAL VARIABLES
133       I    V_streamice,       I    V_streamice,
134       I    0, sNx+1, 0, sNy+1 )       I    0, sNx+1, 0, sNy+1 )
135    
136    
137    
138        err_init = 0. _d 0        err_init = 0. _d 0
139    
140        DO bj = myByLo(myThid), myByHi(myThid)        DO bj = myByLo(myThid), myByHi(myThid)
# Line 124  C     LOCAL VARIABLES Line 148  C     LOCAL VARIABLES
148            err_tempv = 0. _d 0            err_tempv = 0. _d 0
149            IF (STREAMICE_umask(i,j,bi,bj).eq.1) THEN            IF (STREAMICE_umask(i,j,bi,bj).eq.1) THEN
150             err_tempu =             err_tempu =
151       &      ABS (Au_SI(i,j,bi,bj)+ubd_SI(i,j,bi,bj) -       &      ABS (Au_SI(i,j,bi,bj)+0*ubd_SI(i,j,bi,bj) -
152       &           taudx_SI(i,j,bi,bj))       &           taudx_SI(i,j,bi,bj))
153    !            print *, "err_temp_u", err_tempu
154            ENDIF            ENDIF
155            IF (STREAMICE_vmask(i,j,bi,bj).eq.1) THEN            IF (STREAMICE_vmask(i,j,bi,bj).eq.1) THEN
156             err_tempv = MAX( err_tempu,             err_tempv = MAX( err_tempu,
157       &      ABS (Av_SI(i,j,bi,bj)+vbd_SI(i,j,bi,bj) -       &      ABS (Av_SI(i,j,bi,bj)+0*vbd_SI(i,j,bi,bj) -
158       &           taudy_SI(i,j,bi,bj)))       &           taudy_SI(i,j,bi,bj)))
159            ENDIF            ENDIF
160            IF (err_tempv .ge. err_init) err_init = err_tempv            IF (err_tempv .ge. err_init) THEN
161                err_init = err_tempv
162              ENDIF
163           ENDDO           ENDDO
164          ENDDO          ENDDO
165         ENDDO         ENDDO
# Line 151  C     LOCAL VARIABLES Line 178  C     LOCAL VARIABLES
178    
179        iter_numconv = 0        iter_numconv = 0
180        err_max = err_init        err_max = err_init
181          err_max_fp = streamice_nonlin_tol_fp * 10.
182        err_lastchange = err_init        err_lastchange = err_init
183    
184  C START NL ITER. LOOP  C START NL ITER. LOOP
# Line 167  C       but after convergence loop does Line 195  C       but after convergence loop does
195  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
196  !$TAF STORE cgtol = comlev1_stream_nl, key=ikey_nl  !$TAF STORE cgtol = comlev1_stream_nl, key=ikey_nl
197  !$TAF STORE err_max = comlev1_stream_nl, key=ikey_nl  !$TAF STORE err_max = comlev1_stream_nl, key=ikey_nl
198    !$TAF STORE err_max_fp = comlev1_stream_nl, key=ikey_nl
199  !$TAF STORE err_tempu = comlev1_stream_nl, key=ikey_nl  !$TAF STORE err_tempu = comlev1_stream_nl, key=ikey_nl
200  !$TAF STORE err_tempv = comlev1_stream_nl, key=ikey_nl  !$TAF STORE err_tempv = comlev1_stream_nl, key=ikey_nl
201  !$TAF STORE err_lastchange = comlev1_stream_nl, key=ikey_nl  !$TAF STORE err_lastchange = comlev1_stream_nl, key=ikey_nl
202  !$TAF STORE ru_old_si = comlev1_stream_nl, key=ikey_nl  !$TAF STORE ru_old_si = comlev1_stream_nl, key=ikey_nl
203  !$TAF STORE rv_old_si = comlev1_stream_nl, key=ikey_nl  !$TAF STORE rv_old_si = comlev1_stream_nl, key=ikey_nl
204    !$TAF STORE u_old_si = comlev1_stream_nl, key=ikey_nl
205    !$TAF STORE v_old_si = comlev1_stream_nl, key=ikey_nl
206  !$TAF STORE streamice_cg_a1 = comlev1_stream_nl, key=ikey_nl  !$TAF STORE streamice_cg_a1 = comlev1_stream_nl, key=ikey_nl
207  !$TAF STORE streamice_cg_a2 = comlev1_stream_nl, key=ikey_nl  !$TAF STORE streamice_cg_a2 = comlev1_stream_nl, key=ikey_nl
208  !$TAF STORE streamice_cg_a3 = comlev1_stream_nl, key=ikey_nl  !$TAF STORE streamice_cg_a3 = comlev1_stream_nl, key=ikey_nl
# Line 184  C       but after convergence loop does Line 215  C       but after convergence loop does
215  !$TAF STORE zv_old_si = comlev1_stream_nl, key=ikey_nl  !$TAF STORE zv_old_si = comlev1_stream_nl, key=ikey_nl
216  #endif  #endif
217    
218         IF (err_max .GT. streamice_nonlin_tol * err_init) THEN         IF ((err_max .GT. streamice_nonlin_tol * err_init) .and.
219         &     (err_max_fp .GT. streamice_nonlin_tol_fp)) THEN
220    
221         iter_numconv = iter_numconv + 1         iter_numconv = iter_numconv + 1
222    
# Line 230  C       but after convergence loop does Line 262  C       but after convergence loop does
262       &       cg_iters,       &       cg_iters,
263       &       myThid )       &       myThid )
264    
265          #ifdef STREAMICE_HYBRID_STRESS
266    #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    #endif
271    
272    #ifdef STREAMICE_HYBRID_STRESS
273            CALL STREAMICE_TAUB (myThid)
274    #endif
275    
276         WRITE(msgBuf,'(A,I5,A,I4,A)') 'streamice linear solve number',         WRITE(msgBuf,'(A,I5,A,I4,A)') 'streamice linear solve number',
277       &                       iter, " ",       &                       iter, " ",
# Line 242  C       but after convergence loop does Line 283  C       but after convergence loop does
283  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
284  !$TAF STORE U_streamice = comlev1_stream_nl, key=ikey_nl  !$TAF STORE U_streamice = comlev1_stream_nl, key=ikey_nl
285  !$TAF STORE V_streamice = comlev1_stream_nl, key=ikey_nl  !$TAF STORE V_streamice = comlev1_stream_nl, key=ikey_nl
286    #ifdef STREAMICE_HYBRID_STRESS
287    !$TAF STORE streamice_taubx = comlev1_stream_nl, key=ikey_nl
288    !$TAF STORE streamice_tauby = comlev1_stream_nl, key=ikey_nl
289    !$TAF STORE visc_streamice_full = comlev1_stream_nl, key=ikey_nl
290    #endif
291  #endif  #endif
292    
293    #ifdef STREAMICE_HYBRID_STRESS
294           CALL STREAMICE_VISC_BETA_HYBRID ( myThid )
295    #else
296         CALL STREAMICE_VISC_BETA ( myThid )         CALL STREAMICE_VISC_BETA ( myThid )
297    #endif
298    
299    
300  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
# Line 284  C       but after convergence loop does Line 334  C       but after convergence loop does
334       I    V_streamice,       I    V_streamice,
335       I    0, sNx+1, 0, sNy+1 )       I    0, sNx+1, 0, sNy+1 )
336    
337          if (iter .eq. streamice_max_nl_iter) then
338          CALL WRITE_FLD_XY_RL("u_bound_cont_A","",ubd_SI,0,myThid)
339          CALL WRITE_FLD_XY_RL("v_bound_cont_A","",vbd_SI,0,myThid)
340          CALL WRITE_FLD_XY_RL("u_bound_cont_B","",Au_SI,0,myThid)
341          CALL WRITE_FLD_XY_RL("v_bound_cont_B","",Av_SI,0,myThid)
342          endif
343    
344  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
345  !$TAF STORE U_streamice = comlev1_stream_nl, key=ikey_nl  !$TAF STORE U_streamice = comlev1_stream_nl, key=ikey_nl
346  !$TAF STORE V_streamice = comlev1_stream_nl, key=ikey_nl  !$TAF STORE V_streamice = comlev1_stream_nl, key=ikey_nl
347  #endif  #endif
348    
349         err_max = 0. _d 0         err_max = 0. _d 0
350           err_max_fp = 0. _d 0
351    
352  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
353  !$TAF STORE err_max = comlev1_stream_nl, key=ikey_nl  !$TAF STORE err_max = comlev1_stream_nl, key=ikey_nl
# Line 305  C       but after convergence loop does Line 363  C       but after convergence loop does
363             err_tempv = 0. _d 0             err_tempv = 0. _d 0
364             IF (STREAMICE_umask(i,j,bi,bj).eq.1) THEN             IF (STREAMICE_umask(i,j,bi,bj).eq.1) THEN
365              err_tempu =              err_tempu =
366       &       ABS (Au_SI(i,j,bi,bj)+ubd_SI(i,j,bi,bj) -       &       ABS (Au_SI(i,j,bi,bj)+0*ubd_SI(i,j,bi,bj) -
367       &            taudx_SI(i,j,bi,bj))       &            taudx_SI(i,j,bi,bj))
368             ENDIF             ENDIF
369             IF (STREAMICE_vmask(i,j,bi,bj).eq.1) THEN             IF (STREAMICE_vmask(i,j,bi,bj).eq.1) THEN
370              err_tempv = MAX( err_tempu,              err_tempv = MAX( err_tempu,
371       &       ABS (Av_SI(i,j,bi,bj)+vbd_SI(i,j,bi,bj) -       &       ABS (Av_SI(i,j,bi,bj)+0*vbd_SI(i,j,bi,bj) -
372       &            taudy_SI(i,j,bi,bj)))       &            taudy_SI(i,j,bi,bj)))
373             ENDIF             ENDIF
374             IF (err_tempv .ge. err_max) err_max = err_tempv  !           if (err_tempu.ge.1.e2.or.err_tempv.ge.1.e2) THEN
375    !            print *, "FOUND MAX ", i,j,err_tempu,err_tempv,
376    !     &      ubd_SI(i,j,bi,bj),vbd_SI(i,j,bi,bj)
377    !           endif
378               IF (err_tempv .ge. err_max) THEN
379                err_max = err_tempv
380               ENDIF
381            ENDDO            ENDDO
382           ENDDO           ENDDO
383          ENDDO          ENDDO
# Line 321  C       but after convergence loop does Line 385  C       but after convergence loop does
385    
386         CALL GLOBAL_MAX_R8 (err_max, myThid)         CALL GLOBAL_MAX_R8 (err_max, myThid)
387    
388         WRITE(msgBuf,'(A,F11.7)') 'err/err_init',         DO bj = myByLo(myThid), myByHi(myThid)
389            DO bi = myBxLo(myThid), myBxHi(myThid)
390             DO j=1,sNy
391              DO i=1,sNx
392               err_tempu = 0. _d 0
393               err_tempv = 0. _d 0
394               IF (STREAMICE_umask(i,j,bi,bj).eq.1) THEN
395                err_tempu =
396         &       ABS (U_streamice(i,j,bi,bj)-u_old_SI(i,j,bi,bj))
397               ENDIF
398               IF (STREAMICE_vmask(i,j,bi,bj).eq.1) THEN
399                err_tempv = MAX( err_tempu,
400         &       ABS (V_streamice(i,j,bi,bj)-v_old_SI(i,j,bi,bj)))
401               ENDIF
402               IF (err_tempv .ge. err_max_fp) err_max_fp = err_tempv
403              ENDDO
404             ENDDO
405            ENDDO
406           ENDDO
407    
408           CALL GLOBAL_MAX_R8 (err_max_fp, myThid)
409           WRITE(msgBuf,'(A,E15.7)') '||x_i-x_{i-1}||_inf',
410         &                       err_max_fp
411           CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
412         &                    SQUEEZE_RIGHT , 1)
413    
414           DO bj = myByLo(myThid), myByHi(myThid)
415            DO bi = myBxLo(myThid), myBxHi(myThid)
416             DO j=1-OLy,sNy+OLy
417              DO i=1-OLx,sNx+OLx
418               u_old_SI (i,j,bi,bj) = U_streamice (i,j,bi,bj)
419               v_old_SI (i,j,bi,bj) = V_streamice (i,j,bi,bj)
420              ENDDO
421             ENDDO
422            ENDDO
423           ENDDO
424    
425           WRITE(msgBuf,'(A,E15.7)') 'err/err_init',
426       &                       err_max/err_init       &                       err_max/err_init
427         CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,         CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
428       &                    SQUEEZE_RIGHT , 1)       &                    SQUEEZE_RIGHT , 1)
# Line 344  C END NL ITER. LOOP Line 445  C END NL ITER. LOOP
445  C-------------------------------------------------------------------  C-------------------------------------------------------------------
446    
447        if (iter_numconv .lt. streamice_max_nl_iter) then        if (iter_numconv .lt. streamice_max_nl_iter) then
448         PRINT *, "VELOCITY SOLVE CONVERGED, ", iter_numconv,         WRITE(msgBuf,'(A,I5,A)') 'VELOCITY SOLVE CONVERGED, ',
449       & " iterations"       &         iter_numconv, ' iterations'
450           CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
451         &                     SQUEEZE_RIGHT , 1)
452        else        else
453         PRINT *, "VELOCITY SOLVE DID NOT CONVERGE IN ",         WRITE(msgBuf,'(A,I5,A)') 'VELOCITY SOLVE NOT CONVERGED IN ',
454       & iter, " iterations"       &         iter_numconv, ' iterations'
455           CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
456         &                     SQUEEZE_RIGHT , 1)
457        endif        endif
458    
 !        DO bj = myByLo(myThid), myByHi(myThid)  
 !         DO bi = myBxLo(myThid), myBxHi(myThid)  
 !          DO j=1,sNy  
 !           DO i=1,sNx  
 !            U_streamice (i,j,bi,bj) = 0. _d 0  
 !            V_streamice (i,j,bi,bj) = 0. _d 0  
 !           ENDDO  
 !          ENDDO  
 !         ENDDO  
 !        ENDDO  
 !  
 !        CALL STREAMICE_CG_WRAPPER(  
 !      &       U_streamice,  
 !      &       V_streamice,  
 !      &       taudx_SI,  
 !      &       taudy_SI,  
 !      &       cgtol,  
 !      &       cg_iters,  
 !      &       myThid )  
   
   
459        _EXCH_XY_RL (U_streamice, myThid)        _EXCH_XY_RL (U_streamice, myThid)
460        _EXCH_XY_RL (V_streamice, myThid)        _EXCH_XY_RL (V_streamice, myThid)
461    
Legend:
Removed from v.1.4  
changed lines
  Added in v.1.11
  ViewVC Help
Powered by ViewVC 1.1.22