/[MITgcm]/MITgcm/pkg/ctrl/ctrl_getobcsw.F
ViewVC logotype

Diff of /MITgcm/pkg/ctrl/ctrl_getobcsw.F

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

revision 1.6 by mlosch, Fri Dec 3 00:48:57 2004 UTC revision 1.10 by mlosch, Mon Mar 7 09:24:10 2011 UTC
# Line 1  Line 1 
1    C $Header$
2    C $Name$
3    
4  #include "CTRL_CPPOPTIONS.h"  #include "CTRL_CPPOPTIONS.h"
5  #ifdef ALLOW_OBCS  #ifdef ALLOW_OBCS
# Line 67  c     == local variables == Line 69  c     == local variables ==
69        integer obcswcount1        integer obcswcount1
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
# Line 107  c--   Now, read the control vector. Line 110  c--   Now, read the control vector.
110    
111        if (optimcycle .ge. 0) then        if (optimcycle .ge. 0) then
112          ilobcsw=ilnblnk( xx_obcsw_file )          ilobcsw=ilnblnk( xx_obcsw_file )
113          write(fnameobcsw(1:80),'(2a,i10.10)')          write(fnameobcsw(1:80),'(2a,i10.10)')
114       &       xx_obcsw_file(1:ilobcsw), '.', optimcycle       &       xx_obcsw_file(1:ilobcsw), '.', optimcycle
115        endif        endif
116    
# Line 118  c--   Get the counters, flags, and the i Line 121  c--   Get the counters, flags, and the i
121       O                   obcswcount0,obcswcount1,       O                   obcswcount0,obcswcount1,
122       I                   mytime, myiter, mythid )       I                   mytime, myiter, mythid )
123    
124    CML      print *,'ml-getobcs: ',myIter,obcswfirst,obcswchanged,
125    CML     &     obcswcount0,obcswcount1,obcswfac
126        do iobcs = 1,nobcs        do iobcs = 1,nobcs
127          if ( obcswfirst ) then          if ( obcswfirst ) then
128            call active_read_yz_loc( fnameobcsw, tmpfldyz,            call active_read_yz( fnameobcsw, tmpfldyz,
129       &                         (obcswcount0-1)*nobcs+iobcs,       &                         (obcswcount0-1)*nobcs+iobcs,
130       &                         doglobalread, ladinit, optimcycle,       &                         doglobalread, ladinit, optimcycle,
131       &                         mythid, xx_obcsw_dummy )       &                         mythid, xx_obcsw_dummy )
132    
133  #ifdef ALLOW_CTRL_OBCS_BALANCE  #ifdef ALLOW_CTRL_OBCS_BALANCE
134    
135            if ( optimcycle .gt. 0) then                      if ( optimcycle .gt. 0) then
136              if (iobcs .eq. 3) then              if (iobcs .eq. 3) then
137  cgg         Special attention is needed for the normal velocity.  cgg         Special attention is needed for the normal velocity.
138  cgg         For the north, this is the v velocity, iobcs = 4.  cgg         For the north, this is the v velocity, iobcs = 4.
# Line 144  cgg         The barotropic velocity is s Line 149  cgg         The barotropic velocity is s
149    
150                      do k = 1,Nr                      do k = 1,Nr
151  cgg    If cells are not full, this should be modified with hFac.  cgg    If cells are not full, this should be modified with hFac.
152  cgg      cgg
153  cgg    The xx field (tmpfldxz) does not contain the velocity at the  cgg    The xx field (tmpfldyz) does not contain the velocity at the
154  cgg    surface level. This velocity is not independent; it must  cgg    surface level. This velocity is not independent; it must
155  cgg    exactly balance the volume flux, since we are dealing with  cgg    exactly balance the volume flux, since we are dealing with
156  cgg    the baroclinic velocity structure..  cgg    the baroclinic velocity structure..
157                        utop = tmpfldyz(j,k,bi,bj)*                        utop = tmpfldyz(j,k,bi,bj)*
158       &                maskW(i+ip1,j,k,bi,bj) * delR(k) + utop       &                maskW(i+ip1,j,k,bi,bj) * delR(k) + utop
159  cgg    Add the barotropic velocity component.  cgg    Add the barotropic velocity component.
160                        if (maskW(i+ip1,j,k,bi,bj) .ne. 0.) then                        if (maskW(i+ip1,j,k,bi,bj) .ne. 0.) then
161                          tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro                          tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro
162                        endif                        endif
163                      enddo                      enddo
164  cgg    Compute the baroclinic velocity at level 1. Should balance flux.  cgg    Compute the baroclinic velocity at level 1. Should balance flux.
165                    tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj)                    tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj)
166       &                                      - utop / delR(1)       &                                      - utop / delR(1)
167                    enddo                    enddo
168                  enddo                  enddo
# Line 179  cgg         The barotropic velocity is s Line 184  cgg         The barotropic velocity is s
184    
185                      do k = 1,Nr                      do k = 1,Nr
186  cgg    If cells are not full, this should be modified with hFac.  cgg    If cells are not full, this should be modified with hFac.
187  cgg      cgg
188  cgg    The xx field (tmpfldxz) does not contain the velocity at the  cgg    The xx field (tmpfldyz) does not contain the velocity at the
189  cgg    surface level. This velocity is not independent; it must  cgg    surface level. This velocity is not independent; it must
190  cgg    exactly balance the volume flux, since we are dealing with  cgg    exactly balance the volume flux, since we are dealing with
191  cgg    the baroclinic velocity structure..  cgg    the baroclinic velocity structure..
192                        utop = tmpfldyz(j,k,bi,bj)*                        utop = tmpfldyz(j,k,bi,bj)*
193       &                maskS(i,j,k,bi,bj) * delR(k) + utop       &                maskS(i,j,k,bi,bj) * delR(k) + utop
194  cgg    Add the barotropic velocity component.  cgg    Add the barotropic velocity component.
195                        if (maskS(i,j,k,bi,bj) .ne. 0.) then                        if (maskS(i,j,k,bi,bj) .ne. 0.) then
196                          tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro                          tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro
197                        endif                        endif
198                      enddo                      enddo
199  cgg    Compute the baroclinic velocity at level 1. Should balance flux.  cgg    Compute the baroclinic velocity at level 1. Should balance flux.
200                    tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj)                    tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj)
201       &                                      - utop / delR(1)       &                                      - utop / delR(1)
202                    enddo                    enddo
203                  enddo                  enddo
# Line 215  cgg     & Line 220  cgg     &
220          endif          endif
221    
222          if ( (obcswfirst) .or. (obcswchanged)) then          if ( (obcswfirst) .or. (obcswchanged)) then
             
 cgg(    This is a terribly long way to do it. However, the dimensions do not exactly  
 cgg     match up. I will blame Fortran for the ugliness.  
223    
224            do bj = jtlo,jthi            do bj = jtlo,jthi
225              do bi = itlo,ithi             do bi = itlo,ithi
226                do k = 1,nr              do k = 1,nr
227                  do j = jmin,jmax               do j = jmin,jmax
228                    tmpfldyz(j,k,bi,bj) = xx_obcsw1(j,k,bi,bj,iobcs)                xx_obcsw0(j,k,bi,bj,iobcs) = xx_obcsw1(j,k,bi,bj,iobcs)
229                  enddo                tmpfldyz (j,k,bi,bj)       = 0. _d 0
230                enddo               enddo
             enddo  
           enddo  
   
           call exf_swapffields_yz( tmpfldyz2, tmpfldyz, mythid)  
   
           do bj = jtlo,jthi  
             do bi = itlo,ithi  
               do k = 1,nr  
                 do j = jmin,jmax  
                   xx_obcsw0(j,k,bi,bj,iobcs) = tmpfldyz2(j,k,bi,bj)  
                 enddo  
               enddo  
231              enddo              enddo
232               enddo
233            enddo            enddo
234    
235            call active_read_yz_loc( fnameobcsw, tmpfldyz,            call active_read_yz( fnameobcsw, tmpfldyz,
236       &                         (obcswcount1-1)*nobcs+iobcs,       &                         (obcswcount1-1)*nobcs+iobcs,
237       &                         doglobalread, ladinit, optimcycle,       &                         doglobalread, ladinit, optimcycle,
238       &                         mythid, xx_obcsw_dummy )       &                         mythid, xx_obcsw_dummy )
239    
240  #ifdef ALLOW_CTRL_OBCS_BALANCE  #ifdef ALLOW_CTRL_OBCS_BALANCE
241    
242            if ( optimcycle .gt. 0) then                      if ( optimcycle .gt. 0) then
243              if (iobcs .eq. 3) then              if (iobcs .eq. 3) then
244  cgg         Special attention is needed for the normal velocity.  cgg         Special attention is needed for the normal velocity.
245  cgg         For the north, this is the v velocity, iobcs = 4.  cgg         For the north, this is the v velocity, iobcs = 4.
# Line 265  cgg         The barotropic velocity is s Line 256  cgg         The barotropic velocity is s
256    
257                      do k = 1,Nr                      do k = 1,Nr
258  cgg    If cells are not full, this should be modified with hFac.  cgg    If cells are not full, this should be modified with hFac.
259  cgg      cgg
260  cgg    The xx field (tmpfldxz) does not contain the velocity at the  cgg    The xx field (tmpfldyz) does not contain the velocity at the
261  cgg    surface level. This velocity is not independent; it must  cgg    surface level. This velocity is not independent; it must
262  cgg    exactly balance the volume flux, since we are dealing with  cgg    exactly balance the volume flux, since we are dealing with
263  cgg    the baroclinic velocity structure..  cgg    the baroclinic velocity structure..
264                        utop = tmpfldyz(j,k,bi,bj)*                        utop = tmpfldyz(j,k,bi,bj)*
265       &                maskW(i+ip1,j,k,bi,bj) * delR(k) + utop       &                maskW(i+ip1,j,k,bi,bj) * delR(k) + utop
266  cgg    Add the barotropic velocity component.  cgg    Add the barotropic velocity component.
267                        if (maskW(i+ip1,j,k,bi,bj) .ne. 0.) then                        if (maskW(i+ip1,j,k,bi,bj) .ne. 0.) then
268                          tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro                          tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro
269                        endif                        endif
270                      enddo                      enddo
271  cgg    Compute the baroclinic velocity at level 1. Should balance flux.  cgg    Compute the baroclinic velocity at level 1. Should balance flux.
272                      tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj)                      tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj)
273       &                                      - utop / delR(1)       &                                      - utop / delR(1)
274                    enddo                    enddo
275                  enddo                  enddo
# Line 300  cgg         The barotropic velocity is s Line 291  cgg         The barotropic velocity is s
291    
292                      do k = 1,Nr                      do k = 1,Nr
293  cgg    If cells are not full, this should be modified with hFac.  cgg    If cells are not full, this should be modified with hFac.
294  cgg      cgg
295  cgg    The xx field (tmpfldxz) does not contain the velocity at the  cgg    The xx field (tmpfldyz) does not contain the velocity at the
296  cgg    surface level. This velocity is not independent; it must  cgg    surface level. This velocity is not independent; it must
297  cgg    exactly balance the volume flux, since we are dealing with  cgg    exactly balance the volume flux, since we are dealing with
298  cgg    the baroclinic velocity structure..  cgg    the baroclinic velocity structure..
299                        utop = tmpfldyz(j,k,bi,bj)*                        utop = tmpfldyz(j,k,bi,bj)*
300       &                maskS(i,j,k,bi,bj) * delR(k) + utop       &                maskS(i,j,k,bi,bj) * delR(k) + utop
301  cgg    Add the barotropic velocity component.  cgg    Add the barotropic velocity component.
302                        if (maskS(i,j,k,bi,bj) .ne. 0.) then                        if (maskS(i,j,k,bi,bj) .ne. 0.) then
303                          tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro                          tmpfldyz(j,k,bi,bj) = tmpfldyz(j,k,bi,bj)+ ubaro
304                        endif                        endif
305                      enddo                      enddo
306  cgg    Compute the baroclinic velocity at level 1. Should balance flux.  cgg    Compute the baroclinic velocity at level 1. Should balance flux.
307                      tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj)                      tmpfldyz(j,1,bi,bj) = tmpfldyz(j,1,bi,bj)
308       &                                      - utop / delR(1)       &                                      - utop / delR(1)
309                    enddo                    enddo
310                  enddo                  enddo
# Line 335  cgg     & Line 326  cgg     &
326            enddo            enddo
327          endif          endif
328    
329  c--     Add control to model variable.  c--   Add control to model variable.
330          do bj = jtlo, jthi          do bj = jtlo, jthi
331             do bi = itlo, ithi           do bi = itlo, ithi
332  c--        Calculate mask for tracer cells (0 => land, 1 => water).  c--   Calculate mask for tracer cells (0 => land, 1 => water).
333                do k = 1,nr            do k = 1,nr
334                   do j = 1,sny             do j = 1,sny
335                      i = OB_Iw(j,bi,bj)              i = OB_Iw(j,bi,bj)
336                      if (iobcs .EQ. 1) then              if (iobcs .EQ. 1) then
337                         OBWt(j,k,bi,bj) = OBWt (j,k,bi,bj)               OBWt(j,k,bi,bj) = OBWt (j,k,bi,bj)
338       &                 + obcswfac            *xx_obcsw0(j,k,bi,bj,iobcs)       &            + obcswfac            *xx_obcsw0(j,k,bi,bj,iobcs)
339       &                 + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)       &            + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
340                         OBWt(j,k,bi,bj) = OBWt(j,k,bi,bj)               OBWt(j,k,bi,bj) = OBWt(j,k,bi,bj)
341       &                      *maskW(i+ip1,j,k,bi,bj)       &            *maskW(i+ip1,j,k,bi,bj)
342                      else if (iobcs .EQ. 2) then              else if (iobcs .EQ. 2) then
343                         OBWs(j,k,bi,bj) = OBWs (j,k,bi,bj)               OBWs(j,k,bi,bj) = OBWs (j,k,bi,bj)
344       &                 + obcswfac            *xx_obcsw0(j,k,bi,bj,iobcs)       &            + obcswfac            *xx_obcsw0(j,k,bi,bj,iobcs)
345       &                 + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)       &            + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
346                         OBWs(j,k,bi,bj) = OBWs(j,k,bi,bj)               OBWs(j,k,bi,bj) = OBWs(j,k,bi,bj)
347       &                      *maskW(i+ip1,j,k,bi,bj)       &            *maskW(i+ip1,j,k,bi,bj)
348                      else if (iobcs .EQ. 3) then              else if (iobcs .EQ. 3) then
349                         OBWu(j,k,bi,bj) = OBWu (j,k,bi,bj)               OBWu(j,k,bi,bj) = OBWu (j,k,bi,bj)
350       &                 + obcswfac            *xx_obcsw0(j,k,bi,bj,iobcs)       &            + obcswfac            *xx_obcsw0(j,k,bi,bj,iobcs)
351       &                 + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)       &            + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
352                         OBWu(j,k,bi,bj) = OBWu(j,k,bi,bj)               OBWu(j,k,bi,bj) = OBWu(j,k,bi,bj)
353       &                      *maskW(i+ip1,j,k,bi,bj)       &            *maskW(i+ip1,j,k,bi,bj)
354                      else if (iobcs .EQ. 4) then              else if (iobcs .EQ. 4) then
355                         OBWv(j,k,bi,bj) = OBWv (j,k,bi,bj)               OBWv(j,k,bi,bj) = OBWv (j,k,bi,bj)
356       &                 + obcswfac            *xx_obcsw0(j,k,bi,bj,iobcs)       &            + obcswfac            *xx_obcsw0(j,k,bi,bj,iobcs)
357       &                 + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)       &            + (1. _d 0 - obcswfac)*xx_obcsw1(j,k,bi,bj,iobcs)
358                         OBWv(j,k,bi,bj) = OBWv(j,k,bi,bj)               OBWv(j,k,bi,bj) = OBWv(j,k,bi,bj)
359       &                      *maskS(i,j,k,bi,bj)       &            *maskS(i,j,k,bi,bj)
360                      endif              endif
                  enddo  
               enddo  
361             enddo             enddo
362              enddo
363             enddo
364          enddo          enddo
365    
366  C--   End over iobcs loop  C--   End over iobcs loop
Legend:
Removed from v.1.6  
changed lines
  Added in v.1.10
  ViewVC Help
Powered by ViewVC 1.1.22