Back to home page

MITgcm

 
 

    


File indexing completed on 2023-09-21 05:10:46 UTC

view on githubraw file Latest commit 96b00645 on 2023-09-20 15:15:14 UTC
5ca83cd8f7 Dani*0001 #include "STREAMICE_OPTIONS.h"
                0002 
                0003 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0004 
                0005 CBOP
f991e74a3b Jean*0006       SUBROUTINE ADSTREAMICE_CG_SOLVE(
c8fca1659b Jean*0007      U                               U_state,    ! velocities - need to be recalc ed
5ca83cd8f7 Dani*0008      I                               cg_Bu,      ! adjoint of vel (input)
c8fca1659b Jean*0009      U                               V_state,    ! velocities - need to be recalc ed
5ca83cd8f7 Dani*0010      I                               cg_Bv,      ! adjoint of vel (input)
                0011      I                               Bu_state,   ! to recalc velocities
                0012      U                               cg_Uin,     ! adjoint of RHS (output)
                0013      I                               Bv_state,   ! to recalc velocities
                0014      U                               cg_Vin,     ! adjoint of RHS (output)
                0015      I                               A_uu,       ! section of matrix that multiplies u and projects on u
                0016      U                               adA_uu,     ! adjoint of matrix coeffs (output)
                0017      I                               A_uv,       ! section of matrix that multiplies v and projects on u
                0018      U                               adA_uv,     ! adjoint of matrix coeffs (output)
                0019      I                               A_vu,       ! section of matrix that multiplies u and projects on v
                0020      U                               adA_vu,     ! adjoint of matrix coeffs (output)
                0021      I                               A_vv,       ! section of matrix that multiplies v and projects on v
                0022      U                               adA_vv,     ! adjoint of matrix coeffs (output)
f991e74a3b Jean*0023      I                               tolerance,
d2cdb9260d Dani*0024      I                               maxiters,
5ca83cd8f7 Dani*0025      I                               myThid )
f991e74a3b Jean*0026 C     *============================================================*
                0027 C     | SUBROUTINE                                                 |
5ca83cd8f7 Dani*0028 C     | o                                                          |
f991e74a3b Jean*0029 C     *============================================================*
                0030 
                0031 C     !USES:
5ca83cd8f7 Dani*0032       IMPLICIT NONE
                0033 
                0034 C     === Global variables ===
                0035 #include "SIZE.h"
                0036 #include "EEPARAMS.h"
                0037 #include "PARAMS.h"
                0038 #include "STREAMICE.h"
                0039 #include "STREAMICE_CG.h"
                0040 
                0041 C     !INPUT/OUTPUT ARGUMENTS
                0042 C     cg_Uin, cg_Vin - input and output velocities
                0043 C     cg_Bu, cg_Bv - driving stress
                0044       _RL cg_Uin (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0045       _RL cg_Vin (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0046       _RL cg_Bu (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0047       _RL cg_Bv (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0048       _RL U_state (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0049       _RL V_state (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0050       _RL Bu_state (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0051       _RL Bv_state (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
f991e74a3b Jean*0052       _RL
5ca83cd8f7 Dani*0053      & A_uu (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy,-1:1,-1:1),
                0054      & A_vu (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy,-1:1,-1:1),
                0055      & A_uv (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy,-1:1,-1:1),
                0056      & A_vv (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy,-1:1,-1:1),
                0057      & adA_uu (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy,-1:1,-1:1),
                0058      & adA_vu (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy,-1:1,-1:1),
                0059      & adA_uv (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy,-1:1,-1:1),
                0060      & adA_vv (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy,-1:1,-1:1)
                0061       _RL tolerance
d2cdb9260d Dani*0062       INTEGER maxiters
5ca83cd8f7 Dani*0063       INTEGER myThid
                0064 
f991e74a3b Jean*0065 C     !LOCAL VARIABLES
96b006450c dngo*0066       INTEGER i, j, bi, bj, conv_flag, tmpiter
                0067       INTEGER colx, coly
5ca83cd8f7 Dani*0068       _RL Utemp (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0069       _RL Vtemp (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0070       _RL UtempSt (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0071       _RL VtempSt (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0072       _RL ad_tolerance
                0073       CHARACTER*(MAX_LEN_MBUF) msgBuf
f991e74a3b Jean*0074 CEOP
5ca83cd8f7 Dani*0075 
96b006450c dngo*0076 c       iters = streamice_max_cg_iter
5ca83cd8f7 Dani*0077 
                0078 #ifdef ALLOW_STREAMICE
                0079 
                0080       WRITE(msgBuf,'(A)') 'CALLING MANUAL CG ADJOINT'
                0081        CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
                0082      &                     SQUEEZE_RIGHT , 1)
                0083 
96b006450c dngo*0084 c      print *, "GOT HERE myThid=", myThid, tolerance
d2cdb9260d Dani*0085 
5ca83cd8f7 Dani*0086       conv_flag = 0
                0087       ad_tolerance = 1.e-14
                0088 
                0089       DO bj = myByLo(myThid), myByHi(myThid)
                0090        DO bi = myBxLo(myThid), myBxHi(myThid)
f991e74a3b Jean*0091         DO j=1-OLy,sNy+OLy
                0092          DO i=1-OLx,sNx+OLx
5ca83cd8f7 Dani*0093           Utemp (i,j,bi,bj) =
                0094      &     cg_Uin (i,j,bi,bj)
                0095           Vtemp (i,j,bi,bj) =
                0096      &     cg_Vin (i,j,bi,bj)
                0097           UtempSt (i,j,bi,bj) =
                0098      &     U_state (i,j,bi,bj)
                0099           VtempSt (i,j,bi,bj) =
                0100      &     V_state (i,j,bi,bj)
                0101          ENDDO
                0102         ENDDO
                0103        ENDDO
                0104       ENDDO
                0105 
96b006450c dngo*0106 c      print *, "GOT HERE 2 myThid=", myThid, tolerance
d2cdb9260d Dani*0107 
59f7526048 Dani*0108 #ifndef ALLOW_OPENAD
f991e74a3b Jean*0109       CALL STREAMICE_CG_SOLVE(
5ca83cd8f7 Dani*0110      &  U_state,
                0111      &  V_state,
                0112      &  Bu_state,
                0113      &  Bv_state,
                0114      &  A_uu,
                0115      &  A_uv,
                0116      &  A_vu,
f991e74a3b Jean*0117      &  A_vv,
                0118      &  tolerance,
5ca83cd8f7 Dani*0119      &  tmpiter,
d2cdb9260d Dani*0120      &  maxiters,
5ca83cd8f7 Dani*0121      &  myThid )
59f7526048 Dani*0122 #endif
5ca83cd8f7 Dani*0123 
96b006450c dngo*0124 c      print *, "GOT HERE 3 myThid=", myThid, tolerance
d2cdb9260d Dani*0125 
5ca83cd8f7 Dani*0126       tmpiter = 0
                0127 
                0128       _EXCH_XY_RL( cg_Bu, myThid )
                0129       _EXCH_XY_RL( cg_Bv, myThid )
                0130 
f991e74a3b Jean*0131       CALL STREAMICE_CG_SOLVE(
5ca83cd8f7 Dani*0132      &  cg_Uin,
                0133      &  cg_Vin,
                0134      &  cg_Bu,
                0135      &  cg_Bv,
                0136      &  A_uu,
                0137      &  A_uv,
                0138      &  A_vu,
f991e74a3b Jean*0139      &  A_vv,
                0140      &  ad_tolerance,
5ca83cd8f7 Dani*0141      &  tmpiter,
d2cdb9260d Dani*0142      &  maxiters,
5ca83cd8f7 Dani*0143      &  myThid )
                0144 
96b006450c dngo*0145 c      print *, "GOT HERE 4 myThid=", myThid, tolerance
d2cdb9260d Dani*0146 
5ca83cd8f7 Dani*0147       _EXCH_XY_RL( cg_Uin, myThid )
f991e74a3b Jean*0148       _EXCH_XY_RL( cg_Vin, myThid )
                0149       _EXCH_XY_RL( U_state, myThid )
                0150       _EXCH_XY_RL( V_state, myThid )
5ca83cd8f7 Dani*0151 
                0152       DO bj = myByLo(myThid), myByHi(myThid)
                0153        DO bi = myBxLo(myThid), myBxHi(myThid)
5721e29731 Dani*0154         DO j=1,sNy
                0155          DO i=1,sNx
5ca83cd8f7 Dani*0156           DO colx=-1,1
                0157            DO coly=-1,1
                0158 
                0159             if (STREAMICE_umask(i,j,bi,bj).eq.1) then
                0160              if (STREAMICE_umask(i+colx,j+coly,bi,bj).eq.1) then
f991e74a3b Jean*0161                 adA_uu(i,j,bi,bj,colx,coly) =
                0162      &           adA_uu(i,j,bi,bj,colx,coly) -
5ca83cd8f7 Dani*0163      &           cg_Uin(i,j,bi,bj) *
                0164      &           U_state(i+colx,j+coly,bi,bj)
                0165 
                0166              endif
                0167              if (STREAMICE_vmask(i+colx,j+coly,bi,bj).eq.1) then
                0168                 adA_uv(i,j,bi,bj,colx,coly) =
f991e74a3b Jean*0169      &           adA_uv(i,j,bi,bj,colx,coly) -
5ca83cd8f7 Dani*0170      &           cg_Uin(i,j,bi,bj) *
                0171      &           V_state(i+colx,j+coly,bi,bj)
                0172              endif
                0173             endif
                0174 
                0175             if (STREAMICE_vmask(i,j,bi,bj).eq.1) then
                0176              if (STREAMICE_umask(i+colx,j+coly,bi,bj).eq.1) then
                0177                 adA_vu(i,j,bi,bj,colx,coly) =
f991e74a3b Jean*0178      &           adA_vu(i,j,bi,bj,colx,coly) -
5ca83cd8f7 Dani*0179      &           cg_Vin(i,j,bi,bj) *
                0180      &           U_state(i+colx,j+coly,bi,bj)
                0181              endif
                0182              if (STREAMICE_vmask(i+colx,j+coly,bi,bj).eq.1) then
                0183                 adA_vv(i,j,bi,bj,colx,coly) =
f991e74a3b Jean*0184      &           adA_vv(i,j,bi,bj,colx,coly) -
5ca83cd8f7 Dani*0185      &           cg_Vin(i,j,bi,bj) *
                0186      &           V_state(i+colx,j+coly,bi,bj)
                0187              endif
                0188             endif
                0189 
                0190            enddo
                0191           enddo
                0192          enddo
                0193         enddo
                0194        enddo
                0195       enddo
                0196 
                0197       DO bj = myByLo(myThid), myByHi(myThid)
                0198        DO bi = myBxLo(myThid), myBxHi(myThid)
f991e74a3b Jean*0199         DO j=1-OLy,sNy+OLy
                0200          DO i=1-OLx,sNx+OLx
5ca83cd8f7 Dani*0201           if (i.lt.1.or.i.gt.sNx.or.
                0202      &        j.lt.1.or.j.gt.sNy) then
                0203            cg_Uin (i,j,bi,bj) = 0.0
                0204            cg_Vin (i,j,bi,bj) = 0.0
f991e74a3b Jean*0205 
5ca83cd8f7 Dani*0206            DO colx=-1,1
f991e74a3b Jean*0207             DO coly=-1,1
5ca83cd8f7 Dani*0208              ada_uu(i,j,bi,bj,colx,coly)=0.0
                0209              ada_uv(i,j,bi,bj,colx,coly)=0.0
                0210              ada_vu(i,j,bi,bj,colx,coly)=0.0
                0211              ada_vv(i,j,bi,bj,colx,coly)=0.0
                0212             enddo
                0213            enddo
                0214 
                0215           endif
                0216           cg_Uin (i,j,bi,bj) =
f991e74a3b Jean*0217      &     cg_Uin (i,j,bi,bj) +
5ca83cd8f7 Dani*0218      &     Utemp (i,j,bi,bj)
                0219           cg_Vin (i,j,bi,bj) =
f991e74a3b Jean*0220      &     cg_Vin (i,j,bi,bj) +
5ca83cd8f7 Dani*0221      &     Vtemp (i,j,bi,bj)
                0222           cg_bu (i,j,bi,bj) = 0.
                0223           cg_bv (i,j,bi,bj) = 0.
                0224           U_state (i,j,bi,bj) =
                0225      &     UtempSt (i,j,bi,bj)
                0226           V_state (i,j,bi,bj) =
                0227      &     VtempSt (i,j,bi,bj)
                0228          ENDDO
                0229         ENDDO
                0230        ENDDO
                0231       ENDDO
                0232 
                0233       WRITE(msgBuf,'(A,I5,A)') 'DONE WITH MANUAL CG ADJOINT:',tmpiter,
                0234      & 'ITERS'
                0235        CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
                0236      &                     SQUEEZE_RIGHT , 1)
                0237 
                0238 #endif
                0239       RETURN
                0240       END