/[MITgcm]/MITgcm_contrib/SOSE/BoxAdj/code_ad/ctrl_getobcsw.F
ViewVC logotype

Diff of /MITgcm_contrib/SOSE/BoxAdj/code_ad/ctrl_getobcsw.F

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

revision 1.5 by mmazloff, Wed Jan 19 21:08:32 2011 UTC revision 1.6 by mmazloff, Fri Mar 18 18:43:05 2011 UTC
# Line 67  c     == local variables == Line 67  c     == local variables ==
67        logical obcswchanged        logical obcswchanged
68        integer obcswcount0        integer obcswcount0
69        integer obcswcount1        integer obcswcount1
       integer nk,nz  
       _RL     tmpz (nr,nsx,nsy)  
       _RL     stmp  
70    
71  cgg      _RL maskyz   (1-oly:sny+oly,nr,nsx,nsy)  cgg      _RL maskyz   (1-oly:sny+oly,nr,nsx,nsy)
72          _RL tmpfldyz (1-oly:sny+oly,nr,nsx,nsy)
73    
74        logical doglobalread        logical doglobalread
75        logical ladinit        logical ladinit
76    
77        character*(80) fnameobcsw        character*(80) fnameobcsw
78    
79  cgg(  Variables for splitting barotropic/baroclinic vels.  cmm( modes:
80        _RL ubaro        integer nk,nz
81        _RL utop        _RL     tmpz (nr,nsx,nsy)
82  cgg)        _RL     stmp
83    
84  c     == external functions ==  c     == external functions ==
85    
# Line 101  c     == end of interface == Line 99  c     == end of interface ==
99        imax = snx+olx        imax = snx+olx
100        ip1  = 1        ip1  = 1
101    
 cgg(  Initialize variables for balancing volume flux.  
       ubaro = 0.d0  
       utop = 0.d0  
 cgg)  
   
102  c--   Now, read the control vector.  c--   Now, read the control vector.
103        doglobalread = .false.        doglobalread = .false.
104        ladinit      = .false.        ladinit      = .false.
105    
106        if (optimcycle .ge. 0) then        if (optimcycle .ge. 0) then
107          ilobcsw=ilnblnk( xx_obcsw_file )         ilobcsw=ilnblnk( xx_obcsw_file )
108          write(fnameobcsw(1:80),'(2a,i10.10)')         write(fnameobcsw(1:80),'(2a,i10.10)')
109       &       xx_obcsw_file(1:ilobcsw), '.', optimcycle       &      xx_obcsw_file(1:ilobcsw), '.', optimcycle
110        endif        endif
111    
112  c--   Get the counters, flags, and the interpolation factor.  c--   Get the counters, flags, and the interpolation factor.
# Line 124  c--   Get the counters, flags, and the i Line 117  c--   Get the counters, flags, and the i
117       I                   mytime, myiter, mythid )       I                   mytime, myiter, mythid )
118    
119        do iobcs = 1,nobcs        do iobcs = 1,nobcs
120          if ( obcswfirst ) then         if ( obcswfirst ) then
121            call active_read_yz( fnameobcsw, tmpfldyz,          call active_read_yz( fnameobcsw, tmpfldyz,
122       &                         (obcswcount0-1)*nobcs+iobcs,       &                       (obcswcount0-1)*nobcs+iobcs,
123       &                         doglobalread, ladinit, optimcycle,       &                       doglobalread, ladinit, optimcycle,
124       &                         mythid, xx_obcsw_dummy )       &                       mythid, xx_obcsw_dummy )
125    
126            do bj = jtlo,jthi
127            if ( optimcycle .gt. 0) then           do bi = itlo,ithi
             if (iobcs .eq. 3) then  
 cgg         Special attention is needed for the normal velocity.  
 cgg         For the north, this is the v velocity, iobcs = 4.  
 cgg         This is done on a columnwise basis here.  
               do bj = jtlo,jthi  
                 do bi = itlo, ithi  
                   do j = jmin,jmax  
 cih    If open boundary.  
                      if ( OB_Iw(j,bi,bj) .ne. 0. ) then  
                     i = OB_Iw(J,bi,bj)  
 #ifdef ALLOW_OBCS_CONTROL_MODES  
 cih    Determine number of open vertical layers.  
                            nz = 0  
                            do k = 1,Nr  
                               nz = nz + maskW(i+ip1,j,k,bi,bj)  
                            end do  
 cih    Compute absolute velocities from the barotropic-baroclinic modes.  
 #ifdef ALLOW_CTRL_OBCS_BALANCE  
 CMM not sure if ALLOW_OBCS_CONTROL_MODES and ALLOW_CTRL_OBCS_BALANCE are  
 c     compatible - to ensure volume conservation can just set barotropic  
 c     mode amplitude to 0  
 c     however this means no inflow at every horizontal location....  
                    do k = 1,nr  
                      tmpfldyz(k,1,bi,bj)= 0.  
                    end do  
 #endif  
                            do k = 1,Nr  
                              if (k.le.nz) then  
                                 stmp = 0.  
                                 do nk = 1,nz  
                                  stmp = stmp +  
      &                            modesv(k,nk,nz)*tmpfldyz(j,nk,bi,bj)  
                                 end do  
                                 tmpz(k,bi,bj) = stmp  
                              else  
                                 tmpz(k,bi,bj) = 0.  
                              end if  
                            end do  
                            do k = 1,Nr  
                         tmpfldyz(j,k,bi,bj) = tmpz(k,bi,bj)  
      &                           *recip_hFacW(i+ip1,j,k,bi,bj)  
                            end do    
 #elif defined (ALLOW_CTRL_OBCS_BALANCE)  
 cgg         The barotropic velocity is stored in the level 1.  
                     ubaro = tmpfldyz(j,1,bi,bj)  
                     tmpfldyz(j,1,bi,bj) = 0.d0  
                     utop = 0.d0  
   
                     do k = 1,Nr  
 cgg    If cells are not full, this should be modified with hFac.  
 cgg  
 cgg    The xx field (tmpfldxz) does not contain the velocity at the  
 cgg    surface level. This velocity is not independent; it must  
 cgg    exactly balance the volume flux, since we are dealing with  
 cgg    the baroclinic velocity structure..  
                       utop = tmpfldyz(j,k,bi,bj)*  
      &                maskW(i+ip1,j,k,bi,bj) * delR(k) + utop  
 cgg    Add the barotropic velocity component.  
                       if (maskW(i+ip1,j,k,bi,bj) .ne. 0.) then  
                         tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro  
                       endif  
                     enddo  
 cgg    Compute the baroclinic velocity at level 1. Should balance flux.  
                   tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj)  
      &                                      - utop / delR(1)  
 #endif  
 cih    End open boundary.  
                         end if  
                   enddo  
                 enddo  
               enddo  
             endif  
             if (iobcs .eq. 4) then  
 cgg         Special attention is needed for the normal velocity.  
 cgg         For the north, this is the v velocity, iobcs = 4.  
 cgg         This is done on a columnwise basis here.  
               do bj = jtlo,jthi  
                 do bi = itlo, ithi  
                   do j = jmin,jmax  
 cih    If open boundary.  
                      if ( OB_Iw(j,bi,bj) .ne. 0. ) then  
                     i = OB_Iw(J,bi,bj)  
128  #ifdef ALLOW_OBCS_CONTROL_MODES  #ifdef ALLOW_OBCS_CONTROL_MODES
129              if (iobcs .gt. 2) then
130               do j = jmin,jmax
131                i = OB_Iw(j,bi,bj)
132  cih    Determine number of open vertical layers.  cih    Determine number of open vertical layers.
133                             nz = 0              nz = 0
134                             do k = 1,Nr              do k = 1,Nr
135                                nz = nz + maskS(i,j,k,bi,bj)                if (iobcs .eq. 3) then
136                             end do                  nz = nz + maskW(i+ip1,j,k,bi,bj)
137  cih    Compute absolute velocities from the barotropic-baroclinic modes.                else
138  #ifdef ALLOW_CTRL_OBCS_BALANCE                  nz = nz + maskS(i,j,k,bi,bj)
139  CMM not sure if ALLOW_OBCS_CONTROL_MODES and ALLOW_CTRL_OBCS_BALANCE are                endif
140  c     compatible - to ensure volume conservation can just set barotropic              end do
141  c     mode amplitude to 0  cih    Compute absolute velocities from the barotropic-baroclinic modes.
142  c     however this means no inflow at every horizontal location....              do k = 1,Nr
143                     do k = 1,nr               if (k.le.nz) then
144                       tmpfldyz(k,1,bi,bj)= 0.                stmp = 0.
145                     end do                do nk = 1,nz
146  #endif                 stmp = stmp +
147                             do k = 1,Nr       &         modesv(k,nk,nz)*tmpfldyz(j,nk,bi,bj)
148                               if (k.le.nz) then                end do
149                                  stmp = 0.                 tmpz(k,bi,bj) = stmp
150                                  do nk = 1,nz               else
151                                   stmp = stmp +                tmpz(k,bi,bj) = 0.
152       &                            modesv(k,nk,nz)*tmpfldyz(j,nk,bi,bj)               end if
                                 end do  
                                 tmpz(k,bi,bj) = stmp  
                              else  
                                 tmpz(k,bi,bj) = 0.  
                              end if  
                            end do    
                            do k = 1,Nr  
                         tmpfldyz(j,k,bi,bj) = tmpz(k,bi,bj)  
      &                           *recip_hFacS(i,j,k,bi,bj)  
                            end do  
 #elif defined (ALLOW_CTRL_OBCS_BALANCE)  
 cgg         The barotropic velocity is stored in the level 1.  
                     ubaro = tmpfldyz(j,1,bi,bj)  
                     tmpfldyz(j,1,bi,bj) = 0.d0  
                     utop = 0.d0  
   
                     do k = 1,Nr  
 cgg    If cells are not full, this should be modified with hFac.  
 cgg  
 cgg    The xx field (tmpfldxz) does not contain the velocity at the  
 cgg    surface level. This velocity is not independent; it must  
 cgg    exactly balance the volume flux, since we are dealing with  
 cgg    the baroclinic velocity structure..  
                       utop = tmpfldyz(j,k,bi,bj)*  
      &                maskS(i,j,k,bi,bj) * delR(k) + utop  
 cgg    Add the barotropic velocity component.  
                       if (maskS(i,j,k,bi,bj) .ne. 0.) then  
                         tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro  
                       endif  
                     enddo  
 cgg    Compute the baroclinic velocity at level 1. Should balance flux.  
                   tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj)  
      &                                      - utop / delR(1)  
 #endif  
 cih    End if iobcs = 4.  
             end if  
                   enddo  
                 enddo  
               enddo  
             endif  
           endif  
   
           do bj = jtlo,jthi  
             do bi = itlo,ithi  
               do k = 1,nr  
                 do j = jmin,jmax  
                   xx_obcsw1(j,k,bi,bj,iobcs)  = tmpfldyz (j,k,bi,bj)  
 cgg     &                                        *   maskyz (j,k,bi,bj)  
                  enddo  
               enddo  
153              enddo              enddo
154                do k = 1,Nr
155                  if (iobcs .eq. 3) then
156                    tmpfldyz(j,k,bi,bj) = tmpz(k,bi,bj)
157         &            *recip_hFacW(i+ip1,j,k,bi,bj)
158                  else
159                    tmpfldyz(j,k,bi,bj) = tmpz(k,bi,bj)
160         &            *recip_hFacS(i,j,k,bi,bj)
161                  endif
162                end do
163               enddo
164              endif
165    #endif
166              do k = 1,nr
167               do j = jmin,jmax
168                xx_obcsw1(j,k,bi,bj,iobcs)  = tmpfldyz (j,k,bi,bj)
169    cgg   &                                        *   maskyz (j,k,bi,bj)
170               enddo
171            enddo            enddo
172          endif           enddo
173            enddo
174           endif
175    
176          if ( (obcswfirst) .or. (obcswchanged)) then         if ( (obcswfirst) .or. (obcswchanged)) then
177    
178            do bj = jtlo,jthi          do bj = jtlo,jthi
179             do bi = itlo,ithi           do bi = itlo,ithi
180              do k = 1,nr            do k = 1,nr
181               do j = jmin,jmax             do j = jmin,jmax
182                xx_obcsw0(j,k,bi,bj,iobcs) = xx_obcsw1(j,k,bi,bj,iobcs)              xx_obcsw0(j,k,bi,bj,iobcs) = xx_obcsw1(j,k,bi,bj,iobcs)
183                tmpfldyz (j,k,bi,bj)       = 0. _d 0              tmpfldyz (j,k,bi,bj)       = 0. _d 0
              enddo  
             enddo  
184             enddo             enddo
185            enddo            enddo
186             enddo
187            enddo
188    
189            call active_read_yz( fnameobcsw, tmpfldyz,          call active_read_yz( fnameobcsw, tmpfldyz,
190       &                         (obcswcount1-1)*nobcs+iobcs,       &                       (obcswcount1-1)*nobcs+iobcs,
191       &                         doglobalread, ladinit, optimcycle,       &                       doglobalread, ladinit, optimcycle,
192       &                         mythid, xx_obcsw_dummy )       &                       mythid, xx_obcsw_dummy )
193    
194            if ( optimcycle .gt. 0) then          do bj = jtlo,jthi
195              if (iobcs .eq. 3) then           do bi = itlo,ithi
 cgg         Special attention is needed for the normal velocity.  
 cgg         For the north, this is the v velocity, iobcs = 4.  
 cgg         This is done on a columnwise basis here.  
               do bj = jtlo,jthi  
                 do bi = itlo, ithi  
                   do j = jmin,jmax  
 cih    If open boundary.  
                      if ( OB_Iw(j,bi,bj) .ne. 0. ) then  
                     i = OB_Iw(J,bi,bj)  
 #ifdef ALLOW_OBCS_CONTROL_MODES  
 cih    Determine number of open vertical layers.  
                            nz = 0  
                            do k = 1,Nr  
                               nz = nz + maskW(i+ip1,j,k,bi,bj)  
                            end do  
 cih    Compute absolute velocities from the barotropic-baroclinic modes.  
 #ifdef ALLOW_CTRL_OBCS_BALANCE  
 CMM not sure if ALLOW_OBCS_CONTROL_MODES and ALLOW_CTRL_OBCS_BALANCE are  
 c     compatible - to ensure volume conservation can just set barotropic  
 c     mode amplitude to 0  
 c     however this means no inflow at every horizontal location....  
                    do k = 1,nr  
                      tmpfldyz(k,1,bi,bj)= 0.  
                    end do  
 #endif  
                            do k = 1,Nr  
                              if (k.le.nz) then  
                                 stmp = 0.  
                                 do nk = 1,nz  
                                  stmp = stmp +  
      &                            modesv(k,nk,nz)*tmpfldyz(j,nk,bi,bj)  
                                 end do  
                                 tmpz(k,bi,bj) = stmp  
                              else  
                                 tmpz(k,bi,bj) = 0.  
                              end if  
                            end do  
                            do k = 1,Nr  
                         tmpfldyz(j,k,bi,bj) = tmpz(k,bi,bj)  
      &                           *recip_hFacW(i+ip1,j,k,bi,bj)  
                            end do    
 #elif defined (ALLOW_CTRL_OBCS_BALANCE)  
 cgg         The barotropic velocity is stored in the level 1.  
                     ubaro = tmpfldyz(j,1,bi,bj)  
                     tmpfldyz(j,1,bi,bj) = 0.d0  
                     utop = 0.d0  
   
                     do k = 1,Nr  
 cgg    If cells are not full, this should be modified with hFac.  
 cgg  
 cgg    The xx field (tmpfldxz) does not contain the velocity at the  
 cgg    surface level. This velocity is not independent; it must  
 cgg    exactly balance the volume flux, since we are dealing with  
 cgg    the baroclinic velocity structure..  
                       utop = tmpfldyz(j,k,bi,bj)*  
      &                maskW(i+ip1,j,k,bi,bj) * delR(k) + utop  
 cgg    Add the barotropic velocity component.  
                       if (maskW(i+ip1,j,k,bi,bj) .ne. 0.) then  
                         tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro  
                       endif  
                     enddo  
 cgg    Compute the baroclinic velocity at level 1. Should balance flux.  
                   tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj)  
      &                                      - utop / delR(1)  
 #endif  
 cih    End open boundary.  
                         end if  
                   enddo  
                 enddo  
               enddo  
             endif  
             if (iobcs .eq. 4) then  
 cgg         Special attention is needed for the normal velocity.  
 cgg         For the north, this is the v velocity, iobcs = 4.  
 cgg         This is done on a columnwise basis here.  
               do bj = jtlo,jthi  
                 do bi = itlo, ithi  
                   do j = jmin,jmax  
 cih    If open boundary.  
                      if ( OB_Iw(j,bi,bj) .ne. 0. ) then  
                     i = OB_Iw(J,bi,bj)  
196  #ifdef ALLOW_OBCS_CONTROL_MODES  #ifdef ALLOW_OBCS_CONTROL_MODES
197              if (iobcs .gt. 2) then
198               do j = jmin,jmax
199                i = OB_Iw(j,bi,bj)
200  cih    Determine number of open vertical layers.  cih    Determine number of open vertical layers.
201                             nz = 0              nz = 0
202                             do k = 1,Nr              do k = 1,Nr
203                                nz = nz + maskS(i,j,k,bi,bj)                if (iobcs .eq. 3) then
204                             end do                  nz = nz + maskW(i+ip1,j,k,bi,bj)
205  cih    Compute absolute velocities from the barotropic-baroclinic modes.                else
206  #ifdef ALLOW_CTRL_OBCS_BALANCE                  nz = nz + maskS(i,j,k,bi,bj)
207  CMM not sure if ALLOW_OBCS_CONTROL_MODES and ALLOW_CTRL_OBCS_BALANCE are                endif
208  c     compatible - to ensure volume conservation can just set barotropic              end do
209  c     mode amplitude to 0  cih    Compute absolute velocities from the barotropic-baroclinic modes.
210  c     however this means no inflow at every horizontal location....              do k = 1,Nr
211                     do k = 1,nr               if (k.le.nz) then
212                       tmpfldyz(k,1,bi,bj)= 0.                stmp = 0.
213                     end do                do nk = 1,nz
214  #endif                 stmp = stmp +
215                             do k = 1,Nr       &         modesv(k,nk,nz)*tmpfldyz(j,nk,bi,bj)
216                               if (k.le.nz) then                end do
217                                  stmp = 0.                 tmpz(k,bi,bj) = stmp
218                                  do nk = 1,nz               else
219                                   stmp = stmp +                 tmpz(k,bi,bj) = 0.
220       &                            modesv(k,nk,nz)*tmpfldyz(j,nk,bi,bj)               end if
                                 end do  
                                 tmpz(k,bi,bj) = stmp  
                              else  
                                 tmpz(k,bi,bj) = 0.  
                              end if  
                            end do    
                            do k = 1,Nr  
                         tmpfldyz(j,k,bi,bj) = tmpz(k,bi,bj)  
      &                           *recip_hFacS(i,j,k,bi,bj)  
                            end do  
 #elif defined (ALLOW_CTRL_OBCS_BALANCE)  
 cgg         The barotropic velocity is stored in the level 1.  
                     ubaro = tmpfldyz(j,1,bi,bj)  
                     tmpfldyz(j,1,bi,bj) = 0.d0  
                     utop = 0.d0  
   
                     do k = 1,Nr  
 cgg    If cells are not full, this should be modified with hFac.  
 cgg  
 cgg    The xx field (tmpfldxz) does not contain the velocity at the  
 cgg    surface level. This velocity is not independent; it must  
 cgg    exactly balance the volume flux, since we are dealing with  
 cgg    the baroclinic velocity structure..  
                       utop = tmpfldyz(j,k,bi,bj)*  
      &                maskS(i,j,k,bi,bj) * delR(k) + utop  
 cgg    Add the barotropic velocity component.  
                       if (maskS(i,j,k,bi,bj) .ne. 0.) then  
                         tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro  
                       endif  
                     enddo  
 cgg    Compute the baroclinic velocity at level 1. Should balance flux.  
                   tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj)  
      &                                      - utop / delR(1)  
 #endif  
 cih    End if iobcs = 4.  
             end if  
                   enddo  
                 enddo  
               enddo  
             endif  
           endif  
   
           do bj = jtlo,jthi  
             do bi = itlo,ithi  
               do k = 1,nr  
                 do j = jmin,jmax  
                   xx_obcsw1 (j,k,bi,bj,iobcs) = tmpfldyz (j,k,bi,bj)  
 cgg     &                                        *   maskyz (j,k,bi,bj)  
                  enddo  
               enddo  
221              enddo              enddo
222            enddo              do k = 1,Nr
223          endif                if (iobcs .eq. 3) then
224                    tmpfldyz(j,k,bi,bj) = tmpz(k,bi,bj)
225  c--     Add control to model variable.       &            *recip_hFacW(i+ip1,j,k,bi,bj)
226          do bj = jtlo, jthi                else
227             do bi = itlo, ithi                  tmpfldyz(j,k,bi,bj) = tmpz(k,bi,bj)
228  c--        Calculate mask for tracer cells (0 => land, 1 => water).       &            *recip_hFacS(i,j,k,bi,bj)
229                do k = 1,nr                endif
230                   do j = 1,sny              end do
231                      i = OB_Iw(j,bi,bj)             enddo
232                      if (iobcs .EQ. 1) then            endif
233                         OBWt(j,k,bi,bj) = OBWt (j,k,bi,bj)  #endif
234       &                 + obcswfac            *xx_obcsw0(j,k,bi,bj,iobcs)            do k = 1,nr
235       &                 + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)             do j = jmin,jmax
236                         OBWt(j,k,bi,bj) = OBWt(j,k,bi,bj)              xx_obcsw1 (j,k,bi,bj,iobcs) = tmpfldyz (j,k,bi,bj)
237       &                      *maskW(i+ip1,j,k,bi,bj)  cgg   &                                        *   maskyz (j,k,bi,bj)
                     else if (iobcs .EQ. 2) then  
                        OBWs(j,k,bi,bj) = OBWs (j,k,bi,bj)  
      &                 + obcswfac            *xx_obcsw0(j,k,bi,bj,iobcs)  
      &                 + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)  
                        OBWs(j,k,bi,bj) = OBWs(j,k,bi,bj)  
      &                      *maskW(i+ip1,j,k,bi,bj)  
                     else if (iobcs .EQ. 3) then  
                        OBWu(j,k,bi,bj) = OBWu (j,k,bi,bj)  
      &                 + obcswfac            *xx_obcsw0(j,k,bi,bj,iobcs)  
      &                 + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)  
                        OBWu(j,k,bi,bj) = OBWu(j,k,bi,bj)  
      &                      *maskW(i+ip1,j,k,bi,bj)  
                     else if (iobcs .EQ. 4) then  
                        OBWv(j,k,bi,bj) = OBWv (j,k,bi,bj)  
      &                 + obcswfac            *xx_obcsw0(j,k,bi,bj,iobcs)  
      &                 + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)  
                        OBWv(j,k,bi,bj) = OBWv(j,k,bi,bj)  
      &                      *maskS(i,j,k,bi,bj)  
                     endif  
                  enddo  
               enddo  
238             enddo             enddo
239              enddo
240             enddo
241          enddo          enddo
242           endif
243    
244    c--   Add control to model variable.
245           do bj = jtlo, jthi
246            do bi = itlo, ithi
247    c--   Calculate mask for tracer cells (0 => land, 1 => water).
248             do k = 1,nr
249              do j = 1,sny
250               i = OB_Iw(j,bi,bj)
251               if (iobcs .EQ. 1) then
252                OBWt(j,k,bi,bj) = OBWt (j,k,bi,bj)
253         &           + obcswfac            *xx_obcsw0(j,k,bi,bj,iobcs)
254         &           + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
255                OBWt(j,k,bi,bj) = OBWt(j,k,bi,bj)
256         &           *maskW(i+ip1,j,k,bi,bj)
257               else if (iobcs .EQ. 2) then
258                OBWs(j,k,bi,bj) = OBWs (j,k,bi,bj)
259         &           + obcswfac            *xx_obcsw0(j,k,bi,bj,iobcs)
260         &           + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
261                OBWs(j,k,bi,bj) = OBWs(j,k,bi,bj)
262         &           *maskW(i+ip1,j,k,bi,bj)
263               else if (iobcs .EQ. 3) then
264                OBWu(j,k,bi,bj) = OBWu (j,k,bi,bj)
265         &           + obcswfac            *xx_obcsw0(j,k,bi,bj,iobcs)
266         &           + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
267                OBWu(j,k,bi,bj) = OBWu(j,k,bi,bj)
268         &           *maskW(i+ip1,j,k,bi,bj)
269               else if (iobcs .EQ. 4) then
270                OBWv(j,k,bi,bj) = OBWv (j,k,bi,bj)
271         &           + obcswfac            *xx_obcsw0(j,k,bi,bj,iobcs)
272         &           + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
273                OBWv(j,k,bi,bj) = OBWv(j,k,bi,bj)
274         &           *maskS(i,j,k,bi,bj)
275               endif
276              enddo
277             enddo
278            enddo
279           enddo
280          
281  C--   End over iobcs loop  C--   End over iobcs loop
282        enddo        enddo
283    
Legend:
Removed from v.1.5  
changed lines
  Added in v.1.6
  ViewVC Help
Powered by ViewVC 1.1.22