Back to home page

MITgcm

 
 

    


File indexing completed on 2021-10-13 05:14:42 UTC

view on githubraw file Latest commit c59dd234 on 2021-10-12 16:00:15 UTC
01719759b8 Jean*0001 #include "CTRL_OPTIONS.h"
6b46535faa Gael*0002 
01719759b8 Jean*0003 C--  File ctrl_bound.F:
                0004 C--   Contents
                0005 C--   o ADCTRL_BOUND_3D
                0006 C--   o ADCTRL_BOUND_2D
                0007 
                0008 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0009 CBOP
6b46535faa Gael*0010 C     !ROUTINE: ADCTRL_BOUND_3D
                0011 C     !INTERFACE:
                0012       SUBROUTINE ADCTRL_BOUND_3D(
01719759b8 Jean*0013      U                fieldCur, adjFieldCur,
c59dd234b1 Jean*0014      I                mask3D, boundsVec, myThid )
01719759b8 Jean*0015 
6b46535faa Gael*0016 C     !DESCRIPTION: \bv
                0017 C     *==========================================================*
f5f1792fb5 Gael*0018 C     | started: Gael Forget gforget@mit.edu 20-Aug-2007
                0019 C     |
                0020 C     | o in forward mode: impose bounds on ctrl vector values
                0021 C     | o in adjoint mode: do nothing ... or emulate local minimum
6b46535faa Gael*0022 C     *==========================================================*
01719759b8 Jean*0023 C     \ev
6b46535faa Gael*0024 
01719759b8 Jean*0025 C     !USES:
                0026       IMPLICIT NONE
6b46535faa Gael*0027 #include "SIZE.h"
                0028 #include "EEPARAMS.h"
                0029 
01719759b8 Jean*0030 C     !INPUT/OUTPUT PARAMETERS:
c59dd234b1 Jean*0031       _RL fieldCur   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
                0032       _RL adjFieldCur(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
                0033       _RS mask3D     (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
f5f1792fb5 Gael*0034       _RL boundsVec(5)
01719759b8 Jean*0035       INTEGER myThid
54efc7cf00 Gael*0036 
                0037 #ifdef ALLOW_ADCTRLBOUND
01719759b8 Jean*0038 C     !LOCAL VARIABLES:
                0039       INTEGER bi,bj,i,j,k
f5f1792fb5 Gael*0040       _RL x0,x0p5,l0p5,x1,x2,x2p5,l2p5,x3
1f293db3d3 Jean*0041       _RL tmpCur,xCur,adxCur
01719759b8 Jean*0042 CEOP
f5f1792fb5 Gael*0043 
                0044       x0=boundsVec(1)
                0045       x1=boundsVec(2)
4b9076472b Gael*0046         x0p5=(x0+x1)/2.0
                0047         l0p5=(x1-x0)/2.0
f5f1792fb5 Gael*0048       x2=boundsVec(3)
                0049       x3=boundsVec(4)
4b9076472b Gael*0050         x2p5=(x2+x3)/2.0
                0051         l2p5=(x3-x2)/2.0
f5f1792fb5 Gael*0052 
fae6796590 Jean*0053 C  x0<x1<x2<x3  => ctrl_bound and adctrl_bound   act on xx/adxx
01719759b8 Jean*0054 C  x0=x3        => ctrl_bound and adctrl_bound   DO nothing
fae6796590 Jean*0055 C  otherwise    => error
f5f1792fb5 Gael*0056 
01719759b8 Jean*0057       IF ( x0.LT.x3 ) THEN
                0058        IF ( (x0.LT.x1).AND.(x1.LT.x2).AND.(x2.LT.x3) ) THEN
                0059 
                0060         DO bj=myByLo(myThid), myByHi(myThid)
                0061          DO bi=myBxLo(myThid), myBxHi(myThid)
                0062 
                0063           DO k = 1,Nr
                0064            DO j = 1,sNy
                0065             DO i = 1,sNx
c59dd234b1 Jean*0066               IF (mask3D(i,j,k,bi,bj).NE.0.) THEN
01719759b8 Jean*0067                 xCur=fieldCur(i,j,k,bi,bj)
                0068                 adxCur=adjFieldCur(i,j,k,bi,bj)
                0069                 IF ( (xCur.GT.x2).AND.(adxCur.LT.0.) ) THEN
                0070                   tmpCur=1.0
                0071                   adjFieldCur(i,j,k,bi,bj)=abs(adxCur)*
                0072      &            min((xCur-x2p5)/l2p5,tmpCur)
                0073                 ENDIF
                0074                 IF ( (xCur.LT.x1).AND.(adxCur.GT.0.) ) THEN
                0075                   tmpCur=-1.0
                0076                   adjFieldCur(i,j,k,bi,bj)=abs(adxCur)*
                0077      &            max((xCur-x0p5)/l0p5,tmpCur)
                0078                 ENDIF
                0079               ENDIF
                0080             ENDDO
                0081            ENDDO
                0082           ENDDO
                0083 
                0084          ENDDO
                0085         ENDDO
                0086 
                0087        ELSE
                0088           PRINT*,"boundsVec is not self-consistent"
                0089           STOP 'ABNORMAL END: S/R ADCTRL_BOUND_3D'
                0090        ENDIF
f5f1792fb5 Gael*0091       ENDIF
                0092 
01719759b8 Jean*0093 #endif /* ALLOW_ADCTRLBOUND */
6b46535faa Gael*0094 
3c0a47830d Jean*0095       RETURN
                0096       END
6b46535faa Gael*0097 
3c0a47830d Jean*0098 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
01719759b8 Jean*0099 CBOP
6b46535faa Gael*0100 C     !ROUTINE: ADCTRL_BOUND_2D
                0101 C     !INTERFACE:
                0102       SUBROUTINE ADCTRL_BOUND_2D(
01719759b8 Jean*0103      U                fieldCur, adjFieldCur,
7b8b86ab99 Timo*0104      I                mask2D, boundsVec, myThid )
01719759b8 Jean*0105 
6b46535faa Gael*0106 C     !DESCRIPTION: \bv
                0107 C     *==========================================================*
f5f1792fb5 Gael*0108 C     | started: Gael Forget gforget@mit.edu 20-Aug-2007
                0109 C     |
                0110 C     | o in forward mode: impose bounds on ctrl vector values
01719759b8 Jean*0111 C     | o in adjoint mode: DO nothing ... or emulate local minimum
6b46535faa Gael*0112 C     *==========================================================*
01719759b8 Jean*0113 C     \ev
6b46535faa Gael*0114 
01719759b8 Jean*0115 C     !USES:
                0116       IMPLICIT NONE
6b46535faa Gael*0117 #include "SIZE.h"
                0118 #include "EEPARAMS.h"
                0119 
01719759b8 Jean*0120 C     !INPUT/OUTPUT PARAMETERS:
c59dd234b1 Jean*0121       _RL fieldCur   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0122       _RL adjFieldCur(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0123       _RS mask2D     (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
f5f1792fb5 Gael*0124       _RL boundsVec(5)
01719759b8 Jean*0125       INTEGER myThid
f5f1792fb5 Gael*0126 
                0127 #ifdef ALLOW_ADCTRLBOUND
01719759b8 Jean*0128 C     !LOCAL VARIABLES:
                0129       INTEGER bi,bj,i,j
                0130       _RL x0,x0p5,l0p5,x1,x2,x2p5,l2p5,x3
                0131       _RL tmpCur,xCur,adxCur
                0132 CEOP
f5f1792fb5 Gael*0133 
                0134       x0=boundsVec(1)
                0135       x1=boundsVec(2)
4b9076472b Gael*0136         x0p5=(x0+x1)/2.0
                0137         l0p5=(x1-x0)/2.0
f5f1792fb5 Gael*0138       x2=boundsVec(3)
                0139       x3=boundsVec(4)
4b9076472b Gael*0140         x2p5=(x2+x3)/2.0
                0141         l2p5=(x3-x2)/2.0
f5f1792fb5 Gael*0142 
fae6796590 Jean*0143 C  x0<x1<x2<x3  => ctrl_bound and adctrl_bound   act on xx/adxx
01719759b8 Jean*0144 C  x0=x3        => ctrl_bound and adctrl_bound   DO nothing
fae6796590 Jean*0145 C  otherwise    => error
f5f1792fb5 Gael*0146 
01719759b8 Jean*0147       IF ( x0.LT.x3 ) THEN
                0148        IF ( (x0.LT.x1).AND.(x1.LT.x2).AND.(x2.LT.x3) ) THEN
                0149 
                0150         DO bj=myByLo(myThid), myByHi(myThid)
                0151          DO bi=myBxLo(myThid), myBxHi(myThid)
                0152 
                0153            DO j = 1,sNy
                0154             DO i = 1,sNx
7b8b86ab99 Timo*0155               IF (mask2D(i,j,bi,bj).NE.0.) THEN
01719759b8 Jean*0156                 xCur=fieldCur(i,j,bi,bj)
                0157                 adxCur=adjFieldCur(i,j,bi,bj)
                0158                 IF ( (xCur.GT.x2).AND.(adxCur.LT.0.) ) THEN
                0159                   tmpCur=1.0
                0160                   adjFieldCur(i,j,bi,bj)=abs(adxCur)*
                0161      &            min((xCur-x2p5)/l2p5,tmpCur)
                0162                 ENDIF
                0163                 IF ( (xCur.LT.x1).AND.(adxCur.GT.0.) ) THEN
                0164                   tmpCur=-1.0
                0165                   adjFieldCur(i,j,bi,bj)=abs(adxCur)*
                0166      &            max((xCur-x0p5)/l0p5,tmpCur)
                0167                 ENDIF
                0168               ENDIF
                0169             ENDDO
                0170            ENDDO
                0171 
                0172          ENDDO
                0173         ENDDO
                0174 
                0175        ELSE
                0176           PRINT*,"boundsVec is not self-consistent"
                0177           STOP 'ABNORMAL END: S/R ADCTRL_BOUND_2D'
                0178        ENDIF
f5f1792fb5 Gael*0179       ENDIF
01719759b8 Jean*0180 
                0181 #endif /* ALLOW_ADCTRLBOUND */
6b46535faa Gael*0182 
3c0a47830d Jean*0183       RETURN
                0184       END