| 7 |
|
|
| 8 |
CBOP |
CBOP |
| 9 |
SUBROUTINE STREAMICE_CG_SOLVE( |
SUBROUTINE STREAMICE_CG_SOLVE( |
| 10 |
U cg_Uin, |
U cg_Uin, ! x-velocities |
| 11 |
U cg_Vin, |
U cg_Vin, ! y-velocities |
| 12 |
I cg_Bu, |
I cg_Bu, ! force in x dir |
| 13 |
I cg_Bv, |
I cg_Bv, ! force in y dir |
| 14 |
|
I A_uu, ! section of matrix that multiplies u and projects on u |
| 15 |
|
I A_uv, ! section of matrix that multiplies v and projects on u |
| 16 |
|
I A_vu, ! section of matrix that multiplies u and projects on v |
| 17 |
|
I A_vv, ! section of matrix that multiplies v and projects on v |
| 18 |
I tolerance, |
I tolerance, |
| 19 |
O iters, |
O iters, |
| 20 |
|
I maxIter, |
| 21 |
I myThid ) |
I myThid ) |
| 22 |
C /============================================================\ |
C /============================================================\ |
| 23 |
C | SUBROUTINE | |
C | SUBROUTINE | |
| 27 |
C \============================================================/ |
C \============================================================/ |
| 28 |
IMPLICIT NONE |
IMPLICIT NONE |
| 29 |
|
|
|
C === Global variables === |
|
| 30 |
#include "SIZE.h" |
#include "SIZE.h" |
| 31 |
#include "EEPARAMS.h" |
#include "EEPARAMS.h" |
| 32 |
#include "PARAMS.h" |
#include "PARAMS.h" |
| 33 |
#include "STREAMICE.h" |
#include "STREAMICE.h" |
| 34 |
#include "STREAMICE_CG.h" |
#include "STREAMICE_CG.h" |
| 35 |
|
|
| 36 |
|
|
| 37 |
|
|
| 38 |
|
!#ifdef ALLOW_PETSC |
| 39 |
|
!#include "finclude/petsc.h" |
| 40 |
|
! UNCOMMENT IF V3.0 |
| 41 |
|
!#include "finclude/petscvec.h" |
| 42 |
|
!#include "finclude/petscmat.h" |
| 43 |
|
!#include "finclude/petscksp.h" |
| 44 |
|
!#include "finclude/petscpc.h" |
| 45 |
|
!#endif |
| 46 |
|
C === Global variables === |
| 47 |
|
|
| 48 |
|
|
| 49 |
C !INPUT/OUTPUT ARGUMENTS |
C !INPUT/OUTPUT ARGUMENTS |
| 50 |
C cg_Uin, cg_Vin - input and output velocities |
C cg_Uin, cg_Vin - input and output velocities |
| 51 |
C cg_Bu, cg_Bv - driving stress |
C cg_Bu, cg_Bv - driving stress |
| 52 |
INTEGER myThid |
INTEGER myThid |
| 53 |
INTEGER iters |
INTEGER iters |
| 54 |
|
INTEGER maxIter |
| 55 |
_RL tolerance |
_RL tolerance |
| 56 |
_RL cg_Uin (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
_RL cg_Uin (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
| 57 |
_RL cg_Vin (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
_RL cg_Vin (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
| 58 |
_RL cg_Bu (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
_RL cg_Bu (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
| 59 |
_RL cg_Bv (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
_RL cg_Bv (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
| 60 |
|
_RL |
| 61 |
|
& A_uu (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy,-1:1,-1:1), |
| 62 |
|
& A_vu (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy,-1:1,-1:1), |
| 63 |
|
& A_uv (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy,-1:1,-1:1), |
| 64 |
|
& A_vv (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy,-1:1,-1:1) |
| 65 |
|
|
| 66 |
C LOCAL VARIABLES |
C LOCAL VARIABLES |
| 67 |
INTEGER i, j, bi, bj, cg_halo, conv_flag |
INTEGER i, j, bi, bj, cg_halo, conv_flag |
| 68 |
INTEGER iter, is, js, ie, je, colx, coly, k |
INTEGER iter, is, js, ie, je, colx, coly, k |
| 69 |
_RL dot_p1, dot_p2, resid_0, alpha_k, beta_k, resid |
_RL dot_p1, dot_p2, alpha_k, beta_k, resid, resid_0 |
| 70 |
_RL dot_p1_tile (nSx,nSy) |
_RL dot_p1_tile (nSx,nSy) |
| 71 |
_RL dot_p2_tile (nSx,nSy) |
_RL dot_p2_tile (nSx,nSy) |
| 72 |
|
CHARACTER*(MAX_LEN_MBUF) msgBuf |
| 73 |
|
|
| 74 |
|
|
| 75 |
|
!#ifdef ALLOW_PETSC |
| 76 |
|
! INTEGER indices(2*(snx*nsx*sny*nsy)) |
| 77 |
|
! INTEGER n_dofs_cum_sum (0:nPx*nPy-1), idx(1) |
| 78 |
|
! _RL rhs_values(2*(snx*nsx*sny*nsy)) |
| 79 |
|
! _RL solution_values(2*(snx*nsx*sny*nsy)) |
| 80 |
|
! _RL mat_values (2*Nx*Ny,2*(snx*nsx*sny*nsy)) |
| 81 |
|
! _RL mat_values (18,1), mat_val_return(1) |
| 82 |
|
! INTEGER indices_col(18) |
| 83 |
|
! INTEGER local_dofs, global_dofs, dof_index, dof_index_col |
| 84 |
|
! INTEGER local_offset |
| 85 |
|
! Mat matrix |
| 86 |
|
! KSP ksp |
| 87 |
|
! PC pc |
| 88 |
|
! Vec rhs |
| 89 |
|
! Vec solution |
| 90 |
|
! PetscErrorCode ierr |
| 91 |
|
!#ifdef ALLOW_USE_MPI |
| 92 |
|
! integer mpiRC, mpiMyWid |
| 93 |
|
!#endif |
| 94 |
|
!#endif |
| 95 |
|
|
|
iters = streamice_max_cg_iter |
|
| 96 |
|
|
| 97 |
#ifdef ALLOW_STREAMICE |
#ifdef ALLOW_STREAMICE |
| 98 |
|
|
| 99 |
|
|
| 100 |
|
|
| 101 |
|
CALL TIMER_START ('STREAMICE_CG_SOLVE',myThid) |
| 102 |
|
#ifndef STREAMICE_SERIAL_TRISOLVE |
| 103 |
|
|
| 104 |
|
#ifdef ALLOW_PETSC |
| 105 |
|
|
| 106 |
|
CALL STREAMICE_CG_SOLVE_PETSC( |
| 107 |
|
U cg_Uin, ! x-velocities |
| 108 |
|
U cg_Vin, ! y-velocities |
| 109 |
|
I cg_Bu, ! force in x dir |
| 110 |
|
I cg_Bv, ! force in y dir |
| 111 |
|
I A_uu, ! section of matrix that multiplies u and projects on u |
| 112 |
|
I A_uv, ! section of matrix that multiplies v and projects on u |
| 113 |
|
I A_vu, ! section of matrix that multiplies u and projects on v |
| 114 |
|
I A_vv, ! section of matrix that multiplies v and projects on v |
| 115 |
|
I tolerance, |
| 116 |
|
I maxIter, |
| 117 |
|
O iters, |
| 118 |
|
I myThid ) |
| 119 |
|
|
| 120 |
|
|
| 121 |
|
#else /* ALLOW_PETSC */ |
| 122 |
|
|
| 123 |
|
|
| 124 |
|
iters = maxIter |
| 125 |
conv_flag = 0 |
conv_flag = 0 |
| 126 |
|
|
| 127 |
|
|
| 128 |
DO bj = myByLo(myThid), myByHi(myThid) |
DO bj = myByLo(myThid), myByHi(myThid) |
| 129 |
DO bi = myBxLo(myThid), myBxHi(myThid) |
DO bi = myBxLo(myThid), myBxHi(myThid) |
| 130 |
DO j=1,sNy |
DO j=1,sNy |
| 131 |
DO i=1,sNx |
DO i=1,sNx |
| 132 |
Zu_SI (i,j,bi,bj) = 0. _d 0 |
Zu_SI (i,j,bi,bj) = 0. _d 0 |
| 133 |
Zv_SI (i,j,bi,bj) = 0. _d 0 |
Zv_SI (i,j,bi,bj) = 0. _d 0 |
|
DIAGu_SI (i,j,bi,bj) = 0. _d 0 |
|
|
DIAGv_SI (i,j,bi,bj) = 0. _d 0 |
|
| 134 |
Ru_SI (i,j,bi,bj) = 0. _d 0 |
Ru_SI (i,j,bi,bj) = 0. _d 0 |
| 135 |
Rv_SI (i,j,bi,bj) = 0. _d 0 |
Rv_SI (i,j,bi,bj) = 0. _d 0 |
| 136 |
Au_SI (i,j,bi,bj) = 0. _d 0 |
Au_SI (i,j,bi,bj) = 0. _d 0 |
| 137 |
Av_SI (i,j,bi,bj) = 0. _d 0 |
Av_SI (i,j,bi,bj) = 0. _d 0 |
| 138 |
Du_SI (i,j,bi,bj) = 0. _d 0 |
Du_SI (i,j,bi,bj) = 0. _d 0 |
| 139 |
Dv_SI (i,j,bi,bj) = 0. _d 0 |
Dv_SI (i,j,bi,bj) = 0. _d 0 |
|
ubd_SI (i,j,bi,bj) = 0. _d 0 |
|
|
vbd_SI (i,j,bi,bj) = 0. _d 0 |
|
| 140 |
ENDDO |
ENDDO |
| 141 |
ENDDO |
ENDDO |
| 142 |
ENDDO |
ENDDO |
| 143 |
ENDDO |
ENDDO |
| 144 |
|
|
| 145 |
|
C FIND INITIAL RESIDUAL, and initialize r |
| 146 |
|
|
| 147 |
C PREAMBLE: get bdry contribution, find matrix diagonal, |
! #ifdef STREAMICE_CONSTRUCT_MATRIX |
|
C initialize iterates R (residual), Z, and D |
|
| 148 |
|
|
| 149 |
CALL STREAMICE_CG_BOUND_VALS( myThid, |
|
|
O ubd_SI, |
|
|
O vbd_SI) |
|
| 150 |
|
|
| 151 |
DO bj = myByLo(myThid), myByHi(myThid) |
DO bj = myByLo(myThid), myByHi(myThid) |
| 152 |
DO bi = myBxLo(myThid), myBxHi(myThid) |
DO bi = myBxLo(myThid), myBxHi(myThid) |
| 153 |
DO j=1-OLy,sNy+OLy |
DO j=1,sNy |
| 154 |
DO i=1-OLx,sNx+OLx |
DO i=1,sNx |
| 155 |
RHSu_SI (i,j,bi,bj) = cg_Bu(i,j,bi,bj) |
DO colx=-1,1 |
| 156 |
& - ubd_SI(i,j,bi,bj) |
DO coly=-1,1 |
| 157 |
RHSv_SI (i,j,bi,bj) = cg_Bv(i,j,bi,bj) |
Au_SI(i,j,bi,bj) = Au_SI(i,j,bi,bj) + |
| 158 |
& - vbd_SI(i,j,bi,bj) |
& A_uu(i,j,bi,bj,colx,coly)* |
| 159 |
ENDDO |
& cg_Uin(i+colx,j+coly,bi,bj)+ |
| 160 |
ENDDO |
& A_uv(i,j,bi,bj,colx,coly)* |
| 161 |
ENDDO |
& cg_Vin(i+colx,j+coly,bi,bj) |
|
ENDDO |
|
|
|
|
|
_EXCH_XY_RL( RHSu_SI, myThid ) |
|
|
_EXCH_XY_RL( RHSv_SI, myThid ) |
|
| 162 |
|
|
|
CALL STREAMICE_CG_ADIAG( myThid, |
|
|
O DIAGu_SI, |
|
|
O DIAGv_SI) |
|
| 163 |
|
|
| 164 |
_EXCH_XY_RL( DIAGu_SI, myThid ) |
Av_SI(i,j,bi,bj) = Av_SI(i,j,bi,bj) + |
| 165 |
_EXCH_XY_RL( DIAGv_SI, myThid ) |
& A_vu(i,j,bi,bj,colx,coly)* |
| 166 |
|
& cg_Uin(i+colx,j+coly,bi,bj)+ |
| 167 |
C FIX PROBLEM WITH PRECOND LATER |
& A_vv(i,j,bi,bj,colx,coly)* |
| 168 |
|
& cg_Vin(i+colx,j+coly,bi,bj) |
| 169 |
! DO bj = myByLo(myThid), myByHi(myThid) |
ENDDO |
| 170 |
! DO bi = myBxLo(myThid), myBxHi(myThid) |
ENDDO |
| 171 |
! DO j=1-OLy,sNy+OLy |
ENDDO |
| 172 |
! DO i=1-OLx,sNx+OLx |
ENDDO |
| 173 |
! DIAGu_SI(i,j,bi,bj)=1.0 |
ENDDO |
| 174 |
! DIAGv_SI(i,j,bi,bj)=1.0 |
ENDDO |
|
! ENDDO |
|
|
! ENDDO |
|
|
! ENDDO |
|
|
! ENDDO |
|
|
|
|
|
CALL STREAMICE_CG_ACTION( myThid, |
|
|
O Au_SI, |
|
|
O Av_SI, |
|
|
I cg_Uin, |
|
|
I cg_Vin, |
|
|
I 0, sNx+1, 0, sNy+1 ) |
|
| 175 |
|
|
| 176 |
|
|
| 177 |
_EXCH_XY_RL( Au_SI, myThid ) |
_EXCH_XY_RL( Au_SI, myThid ) |
| 178 |
_EXCH_XY_RL( Av_SI, myThid ) |
_EXCH_XY_RL( Av_SI, myThid ) |
| 179 |
|
|
| 180 |
|
|
| 181 |
DO bj = myByLo(myThid), myByHi(myThid) |
DO bj = myByLo(myThid), myByHi(myThid) |
| 182 |
DO bi = myBxLo(myThid), myBxHi(myThid) |
DO bi = myBxLo(myThid), myBxHi(myThid) |
| 183 |
DO j=1-OLy,sNy+OLy |
DO j=1-OLy,sNy+OLy |
| 184 |
DO i=1-OLx,sNx+OLx |
DO i=1-OLx,sNx+OLx |
| 185 |
Ru_SI(i,j,bi,bj)=RHSu_SI(i,j,bi,bj)- |
Ru_SI(i,j,bi,bj)=cg_Bu(i,j,bi,bj)- |
| 186 |
& Au_SI(i,j,bi,bj) |
& Au_SI(i,j,bi,bj) |
| 187 |
Rv_SI(i,j,bi,bj)=RHSv_SI(i,j,bi,bj)- |
Rv_SI(i,j,bi,bj)=cg_Bv(i,j,bi,bj)- |
| 188 |
& Av_SI(i,j,bi,bj) |
& Av_SI(i,j,bi,bj) |
| 189 |
ENDDO |
ENDDO |
| 190 |
ENDDO |
ENDDO |
| 206 |
ENDDO |
ENDDO |
| 207 |
ENDDO |
ENDDO |
| 208 |
|
|
| 209 |
|
|
| 210 |
CALL GLOBAL_SUM_TILE_RL( dot_p1_tile, dot_p1, myThid ) |
CALL GLOBAL_SUM_TILE_RL( dot_p1_tile, dot_p1, myThid ) |
| 211 |
resid_0 = sqrt(dot_p1) |
resid_0 = sqrt(dot_p1) |
| 212 |
|
|
| 213 |
DO bj = myByLo(myThid), myByHi(myThid) |
DO bj = myByLo(myThid), myByHi(myThid) |
| 214 |
DO bi = myBxLo(myThid), myBxHi(myThid) |
DO bi = myBxLo(myThid), myBxHi(myThid) |
| 215 |
|
|
| 216 |
|
WRITE(msgBuf,'(A,I1,I1,E14.7)') 'CONJ GRAD INIT RESID LOCAL, ', |
| 217 |
|
& bi,bj, dot_p1_tile(bi,bj) |
| 218 |
|
CALL PRINT_MESSAGE( msgBuf, standardMessageUnit, |
| 219 |
|
& SQUEEZE_RIGHT , 1) |
| 220 |
|
|
| 221 |
|
enddo |
| 222 |
|
enddo |
| 223 |
|
|
| 224 |
|
WRITE(msgBuf,'(A,E14.7)') 'CONJ GRAD INIT RESID, ', |
| 225 |
|
& resid_0 |
| 226 |
|
CALL PRINT_MESSAGE( msgBuf, standardMessageUnit, |
| 227 |
|
& SQUEEZE_RIGHT , 1) |
| 228 |
|
|
| 229 |
|
C CCCCCCCCCCCCCCCCCCCC |
| 230 |
|
|
| 231 |
|
DO bj = myByLo(myThid), myByHi(myThid) |
| 232 |
|
DO bi = myBxLo(myThid), myBxHi(myThid) |
| 233 |
DO j=1-OLy,sNy+OLy |
DO j=1-OLy,sNy+OLy |
| 234 |
DO i=1-OLx,sNx+OLx |
DO i=1-OLx,sNx+OLx |
| 235 |
IF (STREAMICE_umask(i,j,bi,bj).eq.1.0) |
IF (STREAMICE_umask(i,j,bi,bj).eq.1.0) |
| 241 |
ENDDO |
ENDDO |
| 242 |
ENDDO |
ENDDO |
| 243 |
|
|
|
|
|
| 244 |
cg_halo = min(OLx-1,OLy-1) |
cg_halo = min(OLx-1,OLy-1) |
| 245 |
conv_flag = 0 |
conv_flag = 0 |
| 246 |
|
|
| 263 |
c !! MAIN CG LOOP !! |
c !! MAIN CG LOOP !! |
| 264 |
c !! !! |
c !! !! |
| 265 |
c !!!!!!!!!!!!!!!!!! |
c !!!!!!!!!!!!!!!!!! |
| 266 |
|
|
| 267 |
|
|
| 268 |
|
|
| 269 |
|
|
| 270 |
c ! initially, b-grid data is valid up to 3 halo nodes out -- right? (check for MITgcm!!) |
c ! initially, b-grid data is valid up to 3 halo nodes out -- right? (check for MITgcm!!) |
| 271 |
|
|
| 272 |
print *, "BEGINNING MAIN CG LOOP" |
WRITE(msgBuf,'(A)') 'BEGINNING MAIN CG LOOP' |
| 273 |
|
CALL PRINT_MESSAGE( msgBuf, standardMessageUnit, |
| 274 |
|
& SQUEEZE_RIGHT , 1) |
| 275 |
|
|
| 276 |
! IF(STREAMICE_construct_matrix) CALL STREAMICE_CG_MAKE_A(myThid) |
! IF(STREAMICE_construct_matrix) CALL STREAMICE_CG_MAKE_A(myThid) |
|
#ifdef STREAMICE_CONSTRUCT_MATRIX |
|
|
|
|
|
CALL STREAMICE_CG_MAKE_A(myThid) |
|
| 277 |
|
|
|
#endif |
|
| 278 |
|
|
| 279 |
do iter = 1, streamice_max_cg_iter |
do iter = 1, maxIter |
| 280 |
if (resid .gt. tolerance*resid_0) then |
if (resid .gt. tolerance*resid_0) then |
| 281 |
|
|
| 282 |
c to avoid using "exit" |
c to avoid using "exit" |
| 293 |
DO i=1-OLx,sNx+OLx |
DO i=1-OLx,sNx+OLx |
| 294 |
Au_SI(i,j,bi,bj) = 0. _d 0 |
Au_SI(i,j,bi,bj) = 0. _d 0 |
| 295 |
Av_SI(i,j,bi,bj) = 0. _d 0 |
Av_SI(i,j,bi,bj) = 0. _d 0 |
|
! Du_SI(i,j,bi,bj) = Real(i) |
|
|
! Dv_SI(i,j,bi,bj) = 0.0 |
|
| 296 |
ENDDO |
ENDDO |
| 297 |
ENDDO |
ENDDO |
| 298 |
ENDDO |
ENDDO |
| 300 |
|
|
| 301 |
! IF (STREAMICE_construct_matrix) THEN |
! IF (STREAMICE_construct_matrix) THEN |
| 302 |
|
|
| 303 |
#ifdef STREAMICE_CONSTRUCT_MATRIX |
! #ifdef STREAMICE_CONSTRUCT_MATRIX |
| 304 |
|
|
| 305 |
DO bj = myByLo(myThid), myByHi(myThid) |
DO bj = myByLo(myThid), myByHi(myThid) |
| 306 |
DO bi = myBxLo(myThid), myBxHi(myThid) |
DO bi = myBxLo(myThid), myBxHi(myThid) |
| 309 |
DO colx=-1,1 |
DO colx=-1,1 |
| 310 |
DO coly=-1,1 |
DO coly=-1,1 |
| 311 |
Au_SI(i,j,bi,bj) = Au_SI(i,j,bi,bj) + |
Au_SI(i,j,bi,bj) = Au_SI(i,j,bi,bj) + |
| 312 |
& streamice_cg_A1(i,j,bi,bj,colx,coly)* |
& A_uu(i,j,bi,bj,colx,coly)* |
| 313 |
& Du_SI(i+colx,j+coly,bi,bj)+ |
& Du_SI(i+colx,j+coly,bi,bj)+ |
| 314 |
& streamice_cg_A2(i,j,bi,bj,colx,coly)* |
& A_uv(i,j,bi,bj,colx,coly)* |
| 315 |
& Dv_SI(i+colx,j+coly,bi,bj) |
& Dv_SI(i+colx,j+coly,bi,bj) |
| 316 |
Av_SI(i,j,bi,bj) = Av_SI(i,j,bi,bj) + |
Av_SI(i,j,bi,bj) = Av_SI(i,j,bi,bj) + |
| 317 |
& streamice_cg_A3(i,j,bi,bj,colx,coly)* |
& A_vu(i,j,bi,bj,colx,coly)* |
| 318 |
& Du_SI(i+colx,j+coly,bi,bj)+ |
& Du_SI(i+colx,j+coly,bi,bj)+ |
| 319 |
& streamice_cg_A4(i,j,bi,bj,colx,coly)* |
& A_vv(i,j,bi,bj,colx,coly)* |
| 320 |
& Dv_SI(i+colx,j+coly,bi,bj) |
& Dv_SI(i+colx,j+coly,bi,bj) |
| 321 |
ENDDO |
ENDDO |
| 322 |
ENDDO |
ENDDO |
| 326 |
ENDDO |
ENDDO |
| 327 |
|
|
| 328 |
! else |
! else |
| 329 |
#else |
! #else |
| 330 |
|
! |
| 331 |
CALL STREAMICE_CG_ACTION( myThid, |
! CALL STREAMICE_CG_ACTION( myThid, |
| 332 |
O Au_SI, |
! O Au_SI, |
| 333 |
O Av_SI, |
! O Av_SI, |
| 334 |
I Du_SI, |
! I Du_SI, |
| 335 |
I Dv_SI, |
! I Dv_SI, |
| 336 |
I is,ie,js,je) |
! I is,ie,js,je) |
| 337 |
|
! |
| 338 |
! ENDIF |
! ! ENDIF |
| 339 |
|
! |
| 340 |
#endif |
! #endif |
| 341 |
|
|
|
! ! if (iter.eq.1) then |
|
|
! ! CALL WRITE_FLD_XY_RL ("Au2","",Au_SI,0,myThid) |
|
|
! ! CALL WRITE_FLD_XY_RL ("Av2","",Av_SI,0,myThid) |
|
|
! ! CALL WRITE_FLD_XY_RL ("Du","",Du_SI,0,myThid) |
|
|
! ! |
|
|
! ! CALL WRITE_FLD_XY_RL("Dv","",Au_SI,0,myThid) |
|
|
! ! |
|
|
! ! CALL WRITE_FLD_XY_RL("DiagU1","",Diagu_SI,0,myThid) |
|
|
! ! CALL WRITE_FLD_XY_RL("DiagV1","",Diagv_SI,0,myThid) |
|
|
! ! CALL WRITE_FLD_XY_RL("DiagU2","", |
|
|
! ! & streamice_cg_A1(i,j,bi,bj,0,0),0,myThid) |
|
|
! ! CALL WRITE_FLD_XY_RL("DiagV2","", |
|
|
! ! & streamice_cg_A4(i,j,bi,bj,0,0),0,myThid) |
|
|
! ! CALL WRITE_FLD_XY_RL("DiagV2","",Diagv_SI,0,myThid) |
|
|
! ! |
|
|
! ! endif |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
! ! if (iter.eq.1) then |
|
|
! ! CALL WRITE_FLD_XY_RL ("Au1","",Au_SI,0,myThid) |
|
|
! ! CALL WRITE_FLD_XY_RL ("Av1","",Av_SI,0,myThid) |
|
|
! ! |
|
|
! ! endif |
|
|
|
|
|
|
|
| 342 |
|
|
| 343 |
DO bj = myByLo(myThid), myByHi(myThid) |
DO bj = myByLo(myThid), myByHi(myThid) |
| 344 |
DO bi = myBxLo(myThid), myBxHi(myThid) |
DO bi = myBxLo(myThid), myBxHi(myThid) |
| 493 |
_EXCH_XY_RL( cg_Vin, myThid ) |
_EXCH_XY_RL( cg_Vin, myThid ) |
| 494 |
endif |
endif |
| 495 |
|
|
| 496 |
|
|
| 497 |
endif |
endif |
| 498 |
enddo ! end of CG loop |
enddo ! end of CG loop |
| 499 |
|
|
| 500 |
c to avoid using "exit" |
c to avoid using "exit" |
| 501 |
c if iters has reached max_iters there is no convergence |
c if iters has reached max_iters there is no convergence |
| 502 |
|
|
| 503 |
IF (iters .lt. streamice_max_cg_iter) THEN |
IF (iters .lt. maxIter) THEN |
| 504 |
conv_flag = 1 |
conv_flag = 1 |
| 505 |
ENDIF |
ENDIF |
| 506 |
|
|
| 507 |
DO bj = myByLo(myThid), myByHi(myThid) |
! DO bj = myByLo(myThid), myByHi(myThid) |
| 508 |
DO bi = myBxLo(myThid), myBxHi(myThid) |
! DO bi = myBxLo(myThid), myBxHi(myThid) |
| 509 |
DO j=1-OLy,sNy+OLy |
! DO j=1-OLy,sNy+OLy |
| 510 |
DO i=1-OLy,sNx+OLy |
! DO i=1-OLy,sNx+OLy |
| 511 |
IF (STREAMICE_umask(i,j,bi,bj).eq.3.0) |
! IF (STREAMICE_umask(i,j,bi,bj).eq.3.0) |
| 512 |
& cg_Uin(i,j,bi,bj)=u_bdry_values_SI(i,j,bi,bj) |
! & cg_Uin(i,j,bi,bj)=u_bdry_values_SI(i,j,bi,bj) |
| 513 |
IF (STREAMICE_vmask(i,j,bi,bj).eq.3.0) |
! IF (STREAMICE_vmask(i,j,bi,bj).eq.3.0) |
| 514 |
& cg_Vin(i,j,bi,bj)=v_bdry_values_SI(i,j,bi,bj) |
! & cg_Vin(i,j,bi,bj)=v_bdry_values_SI(i,j,bi,bj) |
| 515 |
ENDDO |
! ENDDO |
| 516 |
ENDDO |
! ENDDO |
| 517 |
ENDDO |
! ENDDO |
| 518 |
ENDDO |
! ENDDO |
| 519 |
|
! |
| 520 |
|
! _EXCH_XY_RL( cg_Uin, myThid ) |
| 521 |
|
! _EXCH_XY_RL( cg_Vin, myThid ) |
| 522 |
|
|
| 523 |
_EXCH_XY_RL( cg_Uin, myThid ) |
#endif /* ifndef ALLOW_PETSC */ |
|
_EXCH_XY_RL( cg_Vin, myThid ) |
|
| 524 |
|
|
| 525 |
#endif |
#else /* STREAMICE_SERIAL_TRISOLVE */ |
|
RETURN |
|
|
END |
|
| 526 |
|
|
| 527 |
|
iters = 0 |
| 528 |
|
|
| 529 |
|
CALL STREAMICE_TRIDIAG_SOLVE( |
| 530 |
|
U cg_Uin, ! x-velocities |
| 531 |
|
U cg_Vin, |
| 532 |
|
U cg_Bu, ! force in x dir |
| 533 |
|
I A_uu, ! section of matrix that multiplies u and projects on u |
| 534 |
|
I STREAMICE_umask, |
| 535 |
|
I myThid ) |
| 536 |
|
|
| 537 |
|
#endif |
| 538 |
|
|
| 539 |
|
CALL TIMER_STOP ('STREAMICE_CG_SOLVE',myThid) |
| 540 |
|
|
| 541 |
|
|
| 542 |
|
#endif |
| 543 |
|
RETURN |
| 544 |
|
END |