/[MITgcm]/MITgcm/model/src/impldiff.F
ViewVC logotype

Diff of /MITgcm/model/src/impldiff.F

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

revision 1.11 by heimbach, Mon Sep 11 20:50:39 2000 UTC revision 1.12 by heimbach, Mon Nov 13 16:30:32 2000 UTC
# Line 29  C     == Routine Arguments == Line 29  C     == Routine Arguments ==
29        _RL KappaRX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL KappaRX(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
30        _RS recip_hFac(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)        _RS recip_hFac(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
31        _RL gXnm1(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)        _RL gXnm1(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
32          _RL gYnm1(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
33        INTEGER myThid        INTEGER myThid
34    
35  C     == Local variables ==  C     == Local variables ==
36        INTEGER i,j,k        INTEGER i,j,k
37        _RL a(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL a(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
38        _RL b(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL b(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
39        _RL c(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL c(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
40        _RL ckm1(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL bet(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
       _RL bet(1-Olx:sNx+Olx,1-Oly:sNy+Oly)  
41        _RL gam(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)        _RL gam(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr)
42    
43  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
44        INTEGER kkey        INTEGER kkey
45  #endif  #endif
46    
47  C--   Only need do anything if Nr>1  C--   Initialise
       IF (Nr.GT.1) THEN  
48    
49  C--    Beginning of forward sweep (top level)  C--   Old aLower
50         DO j=jMin,jMax        DO j=1-Oly,sNy+Oly
51          DO i=iMin,iMax         DO i=1-Olx,sNx+Olx
52           c(i,j)=-deltaTX*recip_hFac(i,j,1,bi,bj)*recip_drF(1)           a(i,j,1) = 0. _d 0
53       &         *KappaRX(i,j,2)*recip_drC(2)         ENDDO
54           b(i,j)=1.-c(i,j)        ENDDO
55           bet(i,j)=0.        DO k=2,Nr
56           IF (b(i,j).NE.0.) bet(i,j)=1. / b(i,j)         DO j=1-Oly,sNy+Oly
57            DO i=1-Olx,sNx+Olx
58              a(i,j,k) = -deltaTX*recip_hFac(i,j,k,bi,bj)*recip_drF(k)
59         &               *KappaRX(i,j, k )*recip_drC( k )
60          ENDDO          ENDDO
61         ENDDO         ENDDO
62          ENDDO
63    
64        ENDIF  C--   Old aUpper
65          DO k=1,Nr-1
66           DO j=1-Oly,sNy+Oly
67            DO i=1-Olx,sNx+Olx
68              c(i,j,k) = -deltaTX*recip_hFac(i,j,k,bi,bj)*recip_drF(k)
69         &               *KappaRX(i,j,k+1)*recip_drC(k+1)
70            ENDDO
71           ENDDO
72          ENDDO
73          DO j=1-Oly,sNy+Oly
74           DO i=1-Olx,sNx+Olx
75             c(i,j,Nr) = 0. _d 0
76           ENDDO
77          ENDDO
78    
79  #ifdef ALLOW_AUTODIFF_TAMC  C--   Old aCenter
80  CADJ store bet                 = comlev1_impl, key = idkey        DO k=1,Nr
81  CADJ store gXNm1(:,:,:,bi,bj)  = comlev1_impl, key = idkey         DO j=1-Oly,sNy+Oly
82  #endif          DO i=1-Olx,sNx+Olx
83              b(i,j,k) = 1. _d 0 - c(i,j,k) - a(i,j,k)
84            ENDDO
85           ENDDO
86          ENDDO
87    
88    C--   Old and new gam, bet are the same
89          DO k=1,Nr
90           DO j=1-Oly,sNy+Oly
91            DO i=1-Olx,sNx+Olx
92              bet(i,j,k) = 0. _d 0
93              gam(i,j,k) = 0. _d 0
94            ENDDO
95           ENDDO
96          ENDDO
97    
98    C--   Only need do anything if Nr>1
99          IF (Nr.GT.1) THEN
100    
101           k = 1
102    C--    Beginning of forward sweep (top level)
103         DO j=jMin,jMax         DO j=jMin,jMax
104          DO i=iMin,iMax          DO i=iMin,iMax
105           gXNm1(i,j,1,bi,bj) = gXNm1(i,j,1,bi,bj)*bet(i,j)           IF (b(i,j,1).NE.0.) bet(i,j,1) = 1. _d 0 / b(i,j,1)
106          ENDDO          ENDDO
107         ENDDO         ENDDO
108    
109          ENDIF
110    
111  C--   Middle of forward sweep  C--   Middle of forward sweep
112        IF (Nr.GT.2) THEN        IF (Nr.GT.2) THEN
113    
114         DO k=2,Nr-1  CADJ loop = sequential
115           DO k=2,Nr
116    
117  #ifdef ALLOW_AUTODIFF_TAMC  #ifdef ALLOW_AUTODIFF_TAMC
118          kkey = (idkey-1)*(Nr-2) + k-1          kkey = (idkey-1)*(Nr-2) + k-1
# Line 82  C--   Middle of forward sweep Line 120  C--   Middle of forward sweep
120    
121          DO j=jMin,jMax          DO j=jMin,jMax
122           DO i=iMin,iMax           DO i=iMin,iMax
123            ckm1(i,j)=c(i,j)            gam(i,j,k) = c(i,j,k-1)*bet(i,j,k-1)
124            a(i,j)=-deltaTX*recip_hFac(i,j,k,bi,bj)*recip_drF(k)            IF ( ( b(i,j,k) - a(i,j,k)*gam(i,j,k) ) .NE. 0.)
125       &         *KappaRX(i,j, k )*recip_drC( k )       &        bet(i,j,k) = 1. _d 0 / ( b(i,j,k) - a(i,j,k)*gam(i,j,k) )
           c(i,j)=-deltaTX*recip_hFac(i,j,k,bi,bj)*recip_drF(k)  
      &         *KappaRX(i,j,k+1)*recip_drC(k+1)  
           b(i,j)=1.-c(i,j)-a(i,j)  
          ENDDO  
         ENDDO  
   
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ store ckm1, bet = comlev1_impl_k, key = kkey  
 #endif  
   
         DO j=jMin,jMax  
          DO i=iMin,iMax  
           gam(i,j,k)=ckm1(i,j)*bet(i,j)  
          ENDDO  
         ENDDO  
   
         DO j=jMin,jMax  
          DO i=iMin,iMax  
           bet(i,j)=b(i,j)-a(i,j)*gam(i,j,k)  
           IF (bet(i,j).NE.0.) bet(i,j)=1. / bet(i,j)  
          ENDDO  
         ENDDO  
   
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ store bet                    = comlev1_impl_k, key = kkey  
 CADJ store gXNm1(:,:,k-1:k,bi,bj) = comlev1_impl_k, key = kkey  
 #endif  
   
         DO j=jMin,jMax  
          DO i=iMin,iMax  
           gXnm1(i,j,k,bi,bj)=(gXnm1(i,j,k,bi,bj)  
      &       -a(i,j)*gXnm1(i,j,k-1,bi,bj))*bet(i,j)  
126           ENDDO           ENDDO
127          ENDDO          ENDDO
128    
# Line 125  CADJ store gXNm1(:,:,k-1:k,bi,bj) = coml Line 131  CADJ store gXNm1(:,:,k-1:k,bi,bj) = coml
131        ENDIF        ENDIF
132    
133    
134        IF (Nr.GT.1) THEN        DO j=jMin,jMax
135           DO i=iMin,iMax
136  C--    End of forward sweep (bottom level)          gYNm1(i,j,1,bi,bj) = gXNm1(i,j,1,bi,bj)*bet(i,j,1)
        DO j=jMin,jMax  
         DO i=iMin,iMax  
          ckm1(i,j)=c(i,j)  
          a(i,j)=-deltaTX*recip_hFac(i,j,Nr,bi,bj)*recip_drF(Nr)  
      &        *KappaRX(i,j, Nr )*recip_drC( Nr )  
          b(i,j)=1.-a(i,j)  
         ENDDO  
        ENDDO  
   
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ store ckm1 = comlev1_impl, key = idkey  
 CADJ store a,b  = comlev1_impl, key = idkey  
 CADJ store bet  = comlev1_impl, key = idkey  
 #endif  
   
        DO j=jMin,jMax  
         DO i=iMin,iMax  
          gam(i,j,Nr)=ckm1(i,j)*bet(i,j)  
         ENDDO  
137         ENDDO         ENDDO
138          ENDDO
139          DO k=2,Nr
140         DO j=jMin,jMax         DO j=jMin,jMax
141          DO i=iMin,iMax          DO i=iMin,iMax
142           bet(i,j)=b(i,j)-a(i,j)*gam(i,j,Nr)           gYnm1(i,j,k,bi,bj) = bet(i,j,k)*
143           IF (bet(i,j).NE.0.) bet(i,j)=1. / bet(i,j)       &        (gXnm1(i,j,k,bi,bj) - a(i,j,k)*gYnm1(i,j,k-1,bi,bj))
         ENDDO  
        ENDDO  
   
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ store a,bet              = comlev1_impl, key = idkey  
 CADJ store gXnm1(:,:,:,bi,bj) = comlev1_impl, key = idkey  
 #endif  
   
        DO j=jMin,jMax  
         DO i=iMin,iMax  
          gXnm1(i,j,Nr,bi,bj)=(gXnm1(i,j,Nr,bi,bj)  
      &      -a(i,j)*gXnm1(i,j,Nr-1,bi,bj))*bet(i,j)  
144          ENDDO          ENDDO
145         ENDDO         ENDDO
146          ENDDO
147    
 #ifdef ALLOW_AUTODIFF_TAMC  
 CADJ store gam   = comlev1_impl, key = idkey  
 #endif  
148    
149  C--    Backward sweep  C--    Backward sweep
150    CADJ loop = sequential
151         DO k=Nr-1,1,-1         DO k=Nr-1,1,-1
152          DO j=jMin,jMax          DO j=jMin,jMax
153           DO i=iMin,iMax           DO i=iMin,iMax
154            gXnm1(i,j,k,bi,bj)=gXnm1(i,j,k,bi,bj)            gYnm1(i,j,k,bi,bj)=gYnm1(i,j,k,bi,bj)
155       &              -gam(i,j,k+1)*gXnm1(i,j,k+1,bi,bj)       &              -gam(i,j,k+1)*gYnm1(i,j,k+1,bi,bj)
156           ENDDO           ENDDO
157          ENDDO          ENDDO
158         ENDDO         ENDDO
159    
160        ENDIF         DO k=1,Nr
161            DO j=jMin,jMax
162             DO i=iMin,iMax
163              gXnm1(i,j,k,bi,bj)=gYnm1(i,j,k,bi,bj)
164             ENDDO
165            ENDDO
166           ENDDO
167    
168        RETURN        RETURN
169        END        END
Legend:
Removed from v.1.11  
changed lines
  Added in v.1.12
  ViewVC Help
Powered by ViewVC 1.1.22