pkg/ecco/cost_obcsn.F

8f7d13d0c9 Jean*0001 #include "ECCO_OPTIONS.h"
6805a315c1 Gael*0002 #ifdef ALLOW_CTRL
                0003 # include "CTRL_OPTIONS.h"
                0004 #endif
5001c65f45 Patr*0005
15f6a15ab5 Mart*0006 CBOP
                0007 C     !ROUTINE: COST_OBCSN
                0008 C     !INTERFACE:
5001c65f45 Patr*0009       subroutine cost_obcsn(
720a211d38 Ou W*0010      I                       startrec, endrec,
                0011      I                       myTime, myIter, myThid )
5001c65f45 Patr*0012
15f6a15ab5 Mart*0013 C     !DESCRIPTION: \bv
5001c65f45 Patr*0014 c     ==================================================================
                0015 c     SUBROUTINE cost_obcsn
                0016 c     ==================================================================
                0017 c
                0018 c     o cost function contribution obc
                0019 c
                0020 c     ==================================================================
                0021 c     SUBROUTINE cost_obcsn
                0022 c     ==================================================================
15f6a15ab5 Mart*0023 C     \ev
                0024
                0025 C     !USES:
5001c65f45 Patr*0026
                0027       implicit none
                0028
                0029 c     == global variables ==
                0030 #include "EEPARAMS.h"
                0031 #include "SIZE.h"
                0032 #include "PARAMS.h"
                0033 #include "GRID.h"
f75c953e99 Jean*0034 c#ifdef ALLOW_OBCS
                0035 c# include "OBCS_GRID.h"
                0036 c#endif
5001c65f45 Patr*0037
c509d7e04a Gael*0038 #ifdef ALLOW_CAL
                0039 # include "cal.h"
                0040 #endif
                0041 #ifdef ALLOW_CTRL
                0042 # include "CTRL_SIZE.h"
4d72283393 Mart*0043 # include "CTRL.h"
edcd27be69 Mart*0044 # include "CTRL_DUMMY.h"
65754df434 Mart*0045 # include "OPTIMCYCLE.h"
c509d7e04a Gael*0046 # include "CTRL_OBCS.h"
                0047 #endif
5001c65f45 Patr*0048
15f6a15ab5 Mart*0049 C     !INPUT/OUTPUT PARAMETERS:
720a211d38 Ou W*0050       integer startrec, endrec
                0051       _RL     myTime
                0052       integer myIter
                0053       integer myThid
5001c65f45 Patr*0054
c509d7e04a Gael*0055 #if (defined (ALLOW_CTRL) && defined (ALLOW_OBCS))
                0056
f75c953e99 Jean*0057 #ifdef ALLOW_OBCSN_COST_CONTRIBUTION
720a211d38 Ou W*0058 C     !FUNCTIONS:
f75c953e99 Jean*0059       integer  ilnblnk
                0060       external ilnblnk
                0061
15f6a15ab5 Mart*0062 C     !LOCAL VARIABLES:
5001c65f45 Patr*0063       integer bi,bj
f75c953e99 Jean*0064       integer i,k
5001c65f45 Patr*0065       integer imin,imax
                0066       integer irec
                0067       integer iobcs
15f6a15ab5 Mart*0068       integer nrec
                0069       integer ilfld
                0070       integer igg
5001c65f45 Patr*0071       _RL fctile
                0072       _RL fcthread
                0073       _RL dummy
15f6a15ab5 Mart*0074       _RL gg
                0075       _RL tmpx
9f5240b52a Jean*0076       _RL tmpfield (1-OLx:sNx+OLx,Nr,nSx,nSy)
                0077       _RL maskxz   (1-OLx:sNx+OLx,Nr,nSx,nSy)
de57a2ec4b Mart*0078       character*(MAX_LEN_FNAM) fnamefld
5001c65f45 Patr*0079       logical doglobalread
                0080       logical ladinit
                0081 #ifdef ECCO_VERBOSE
                0082       character*(MAX_LEN_MBUF) msgbuf
                0083 #endif
15f6a15ab5 Mart*0084 CEOP
5001c65f45 Patr*0085
                0086       imin = 1
9f5240b52a Jean*0087       imax = sNx
5001c65f45 Patr*0088
                0089 c--   Read tiled data.
                0090       doglobalread = .false.
                0091       ladinit      = .false.
                0092
15f6a15ab5 Mart*0093 c     Number of records to be used.
                0094       nrec = endrec-startrec+1
                0095
f75c953e99 Jean*0096 c     jp1 = 0
5001c65f45 Patr*0097       fcthread = 0. _d 0
                0098
15f6a15ab5 Mart*0099 #ifdef ECCO_VERBOSE
720a211d38 Ou W*0100       _BEGIN_MASTER( myThid )
15f6a15ab5 Mart*0101       write(msgbuf,'(a)') ' '
720a211d38 Ou W*0102       call print_message( msgbuf, standardMessageUnit,
                0103      &                    SQUEEZE_RIGHT, myThid )
15f6a15ab5 Mart*0104       write(msgbuf,'(a)') ' '
720a211d38 Ou W*0105       call print_message( msgbuf, standardMessageUnit,
                0106      &                    SQUEEZE_RIGHT, myThid )
15f6a15ab5 Mart*0107       write(msgbuf,'(a,i9.8)')
                0108      &  ' cost_obcsn: number of records to process: ',nrec
720a211d38 Ou W*0109       call print_message( msgbuf, standardMessageUnit,
                0110      &                    SQUEEZE_RIGHT, myThid )
15f6a15ab5 Mart*0111       write(msgbuf,'(a)') ' '
720a211d38 Ou W*0112       call print_message( msgbuf, standardMessageUnit,
                0113      &                    SQUEEZE_RIGHT, myThid )
                0114       _END_MASTER( myThid )
15f6a15ab5 Mart*0115 #endif
5001c65f45 Patr*0116
15f6a15ab5 Mart*0117       if (optimcycle .ge. 0) then
                0118         ilfld=ilnblnk( xx_obcsn_file )
de57a2ec4b Mart*0119         write(fnamefld,'(2a,i10.10)')
15f6a15ab5 Mart*0120      &       xx_obcsn_file(1:ilfld), '.', optimcycle
                0121       endif
5001c65f45 Patr*0122
15f6a15ab5 Mart*0123 c--   Loop over records.
                0124       do irec = 1,nrec
                0125
720a211d38 Ou W*0126 #ifdef ALLOW_AUTODIFF
                0127         call active_read_xz( fnamefld, tmpfield, irec,
                0128      &                       doglobalread, ladinit, optimcycle,
                0129      &                       myThid, xx_obcsn_dummy )
                0130 #else
                0131         CALL READ_REC_XZ_RL( fnamefld, ctrlprec, Nr,
                0132      &                       tmpfield, irec, 1, myThid )
                0133 #endif
15f6a15ab5 Mart*0134
720a211d38 Ou W*0135 cgg     Need to solve for iobcs would have been.
                0136         gg    = (irec-1)/nobcs
                0137         igg   = int(gg)
                0138         iobcs = irec - igg*nobcs
                0139
                0140 #ifdef ALLOW_AUTODIFF
                0141         call active_read_xz( 'maskobcsn', maskxz, iobcs,
                0142      &                       doglobalread, ladinit, 0,
                0143      &                       myThid, dummy )
                0144 #else
                0145         CALL READ_REC_XZ_RL( 'maskobcsn', ctrlprec, Nr,
                0146      &                       maskxz, iobcs, 1, myThid )
                0147 #endif
15f6a15ab5 Mart*0148
80451941d6 Jean*0149 c--     Loop over this thread s tiles.
9f5240b52a Jean*0150         do bj = myByLo(myThid), myByHi(myThid)
                0151           do bi = myBxLo(myThid), myBxHi(myThid)
5001c65f45 Patr*0152
15f6a15ab5 Mart*0153 c--         Determine the weights to be used.
5001c65f45 Patr*0154             fctile = 0. _d 0
15f6a15ab5 Mart*0155
                0156             do k = 1, Nr
                0157               do i = imin,imax
f75c953e99 Jean*0158 c                j = OB_Jn(I,bi,bj)
15f6a15ab5 Mart*0159 cgg                 if (maskS(i,j+jp1,k,bi,bj) .ne. 0.) then
                0160                   tmpx = tmpfield(i,k,bi,bj)
                0161 CMM                  fctile = fctile + wobcsn2(i,k,bi,bj,iobcs)
                0162                   fctile = fctile + wobcsn(k,iobcs)
                0163      &                          *tmpx*tmpx*maskxz(i,k,bi,bj)
                0164 cgg                endif
                0165 CMM                 if (wobcsn2(i,k,bi,bj,iobcs)*maskxz(i,k,bi,bj).ne.0.)
                0166                  if (wobcsn(k,iobcs)*maskxz(i,k,bi,bj).ne.0.)
                0167      &                    num_obcsn(bi,bj) = num_obcsn(bi,bj) + 1. _d 0
                0168 cgg                print*,'S fctile',fctile
                0169               enddo
5001c65f45 Patr*0170             enddo
15f6a15ab5 Mart*0171
5001c65f45 Patr*0172             objf_obcsn(bi,bj) = objf_obcsn(bi,bj) + fctile
15f6a15ab5 Mart*0173             fcthread         = fcthread + fctile
5001c65f45 Patr*0174           enddo
                0175         enddo
                0176
                0177 #ifdef ECCO_VERBOSE
                0178 c--     Print cost function for all tiles.
6637358eea Jean*0179         _GLOBAL_SUM_RL( fcthread , myThid )
5001c65f45 Patr*0180         write(msgbuf,'(a)') ' '
720a211d38 Ou W*0181         call print_message( msgbuf, standardMessageUnit,
                0182      &                      SQUEEZE_RIGHT, myThid )
5001c65f45 Patr*0183         write(msgbuf,'(a,i8.8)')
                0184      &    ' cost_obcsn: irec = ',irec
720a211d38 Ou W*0185         call print_message( msgbuf, standardMessageUnit,
                0186      &                      SQUEEZE_RIGHT, myThid )
5001c65f45 Patr*0187         write(msgbuf,'(a,a,d22.15)')
                0188      &    ' global cost function value',
                0189      &    ' (obcsn) = ',fcthread
720a211d38 Ou W*0190         call print_message( msgbuf, standardMessageUnit,
                0191      &                      SQUEEZE_RIGHT, myThid )
5001c65f45 Patr*0192         write(msgbuf,'(a)') ' '
720a211d38 Ou W*0193         call print_message( msgbuf, standardMessageUnit,
                0194      &                      SQUEEZE_RIGHT, myThid )
5001c65f45 Patr*0195 #endif
                0196
                0197       enddo
                0198 c--   End of loop over records.
                0199
f75c953e99 Jean*0200 #endif /* ALLOW_OBCSN_COST_CONTRIBUTION */
5001c65f45 Patr*0201
c509d7e04a Gael*0202 #endif /* ALLOW_CTRL and ALLOW_OBCS */
                0203
5001c65f45 Patr*0204       return
                0205       end