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

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

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

revision 1.16 by jmc, Thu Jan 24 17:22:06 2008 UTC revision 1.17 by jmc, Mon Aug 11 01:19:45 2008 UTC
# Line 2  C $Header$ Line 2  C $Header$
2  C $Name$  C $Name$
3    
4  #include "CPP_OPTIONS.h"  #include "CPP_OPTIONS.h"
5  #define EXCLUDE_EOS_CHECK  #undef INCLUDE_EOS_CHECK
6    
7    C--  File ini_eos.F: Routines to initialise Equation of State
8    C--   Contents
9    C--   o INI_EOS
10    C--   o EOS_CHECK
11    
12    C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
13  CBOP  CBOP
14  C     !ROUTINE: INI_EOS  C     !ROUTINE: INI_EOS
15  C     !INTERFACE:  C     !INTERFACE:
# Line 23  C     == Global variables == Line 29  C     == Global variables ==
29  #include "EEPARAMS.h"  #include "EEPARAMS.h"
30  #include "PARAMS.h"  #include "PARAMS.h"
31  #include "EOS.h"  #include "EOS.h"
 #include "GRID.h"  
 #include "DYNVARS.h"  
32    
33  C     !INPUT/OUTPUT PARAMETERS:  C     !INPUT/OUTPUT PARAMETERS:
34  C     == Routine arguments ==  C     == Routine arguments ==
# Line 249  C     5. secant bulk modulus K of sea wa Line 253  C     5. secant bulk modulus K of sea wa
253    
254        ELSE        ELSE
255    
256           WRITE(msgbuf,'(3a)') ' INI_EOS: equationOfState = "',           WRITE(msgbuf,'(3A)') ' INI_EOS: equationOfState = "',
257       &        equationOfState,'"'       &        equationOfState,'"'
258           CALL PRINT_ERROR( msgbuf, myThid )           CALL PRINT_ERROR( msgbuf, myThid )
259           STOP 'ABNORMAL END: S/R INI_EOS'           STOP 'ABNORMAL END: S/R INI_EOS'
# Line 258  C     5. secant bulk modulus K of sea wa Line 262  C     5. secant bulk modulus K of sea wa
262    
263  C--   Check EOS initialisation:  C--   Check EOS initialisation:
264    
265        CALL CHECK_EOS( myThid )        CALL EOS_CHECK( myThid )
266    
267        _END_MASTER( myThid )        _END_MASTER( myThid )
268        _BARRIER        _BARRIER
# Line 267  C--   Check EOS initialisation: Line 271  C--   Check EOS initialisation:
271        END        END
272    
273  CBOP  CBOP
274  C     !ROUTINE: CHECK_EOS  C     !ROUTINE: EOS_CHECK
275  C     !INTERFACE:  C     !INTERFACE:
276        SUBROUTINE CHECK_EOS( myThid )        SUBROUTINE EOS_CHECK( myThid )
277  C     !DESCRIPTION: \bv  C     !DESCRIPTION: \bv
278  C     *==========================================================*  C     *==========================================================*
279  C     | SUBROUTINE CHECK_EOS  C     | SUBROUTINE EOS_CHECK
280  C     | o check the equation of state.  C     | o check the equation of state.
281  C     *==========================================================*  C     *==========================================================*
282  C     \ev  C     \ev
# Line 284  C     !USES: Line 288  C     !USES:
288  #include "EEPARAMS.h"  #include "EEPARAMS.h"
289  #include "PARAMS.h"  #include "PARAMS.h"
290  #include "EOS.h"  #include "EOS.h"
291    #include "GRID.h"
292    #include "DYNVARS.h"
293    
294  C     !INPUT/OUTPUT PARAMETERS:  C     !INPUT/OUTPUT PARAMETERS:
295  C     == Routine arguments ==  C     == Routine arguments ==
296  C     myThid -  Number of this instance of CHECK_EOS  C     myThid ::  Number of this instance of EOS_CHECK
297        INTEGER myThid        INTEGER myThid
298    
299  #ifndef EXCLUDE_EOS_CHECK  #ifdef INCLUDE_EOS_CHECK
300  C     !LOCAL VARIABLES:  C     !LOCAL VARIABLES:
301  C     == Local variables ==  C     == Local variables ==
302  C     bi,bj  - Loop counters  C     bi,bj  - Loop counters
303  C     I,J,K  C     i,j,k
304        INTEGER bi, bj        INTEGER bi, bj
305        INTEGER imin, imax, jmin, jmax        INTEGER iMin, iMax, jMin, jMax
306        INTEGER  i, j, k        INTEGER  i, j, k
307        _RL tFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)        _RL tFld   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
308        _RL sFld(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)        _RL sFld   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
309          _RL pFld   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
310        _RL rhoLoc (1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL rhoLoc (1-Olx:sNx+Olx,1-Oly:sNy+Oly)
311        _RL bulkMod(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL bulkMod(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
       _RL psave  
312    
313        INTEGER ncheck, kcheck        INTEGER ncheck, kcheck
314        PARAMETER ( ncheck = 13 )        PARAMETER ( ncheck = 13 )
315        _RL tloc(ncheck), ptloc(ncheck), sloc(ncheck), ploc(ncheck)        _RL tLoc(ncheck), ptLoc(ncheck), sLoc(ncheck), pLoc(ncheck)
316        _RL rloc(ncheck), bloc(ncheck)        _RL rLoc(ncheck), bLoc(ncheck)
317          _RS mskSave
318          _RS rC_Save
319    
320        CHARACTER*(MAX_LEN_MBUF) msgBuf        CHARACTER*(MAX_LEN_MBUF) msgBuf
321    
322        DATA tloc        DATA tLoc
323       &     /3.25905152915860 _d 0, 20.38687090048638 _d 0,       &     /3.25905152915860 _d 0, 20.38687090048638 _d 0,
324       &     25.44820830309568 _d 0, 20.17368557065936 _d 0,       &     25.44820830309568 _d 0, 20.17368557065936 _d 0,
325       &     13.43397459640398 _d 0,       &     13.43397459640398 _d 0,
# Line 319  C     I,J,K Line 327  C     I,J,K
327       &      5.               _d 0, 25.               _d 0,       &      5.               _d 0, 25.               _d 0,
328       &      5.               _d 0, 25.               _d 0,       &      5.               _d 0, 25.               _d 0,
329       &      5.               _d 0, 25.               _d 0/,       &      5.               _d 0, 25.               _d 0/,
330       &     ptloc       &     ptLoc
331       &     /3.               _d 0, 20.               _d 0,       &     /3.               _d 0, 20.               _d 0,
332       &     25.               _d 0, 20.               _d 0,       &     25.               _d 0, 20.               _d 0,
333       &     12.               _d 0,       &     12.               _d 0,
# Line 327  C     I,J,K Line 335  C     I,J,K
335       &      5.               _d 0, 25.               _d 0,       &      5.               _d 0, 25.               _d 0,
336       &      4.03692566635316 _d 0, 22.84661726775120 _d 0,       &      4.03692566635316 _d 0, 22.84661726775120 _d 0,
337       &      3.62720389416752 _d 0, 22.62420229124846 _d 0/       &      3.62720389416752 _d 0, 22.62420229124846 _d 0/
338       &     sloc       &     sLoc
339       &     /35.5 _d 0, 35. _d 0,       &     /35.5 _d 0, 35. _d 0,
340       &      35.0 _d 0, 20. _d 0,       &      35.0 _d 0, 20. _d 0,
341       &      40.0 _d 0,       &      40.0 _d 0,
# Line 335  C     I,J,K Line 343  C     I,J,K
343       &      35.  _d 0, 35. _d 0,       &      35.  _d 0, 35. _d 0,
344       &       0.  _d 0,  0. _d 0,       &       0.  _d 0,  0. _d 0,
345       &      35.  _d 0, 35. _d 0/       &      35.  _d 0, 35. _d 0/
346       &     ploc       &     pLoc
347       &     /300. _d 5,  200. _d 5,       &     /300. _d 5,  200. _d 5,
348       &      200. _d 5,  100. _d 5,       &      200. _d 5,  100. _d 5,
349       &      800. _d 5,       &      800. _d 5,
# Line 343  C     I,J,K Line 351  C     I,J,K
351       &        0. _d 0,    0. _d 0,       &        0. _d 0,    0. _d 0,
352       &     1000. _d 5, 1000. _d 5,       &     1000. _d 5, 1000. _d 5,
353       &     1000. _d 5, 1000. _d 5/       &     1000. _d 5, 1000. _d 5/
354        DATA rloc        DATA rLoc
355       &     /1041.83267  _d 0, 1033.213387 _d 0,       &     /1041.83267  _d 0, 1033.213387 _d 0,
356       &      1031.654229 _d 0, 1017.726743 _d 0,       &      1031.654229 _d 0, 1017.726743 _d 0,
357       &      1062.928258 _d 0,       &      1062.928258 _d 0,
# Line 351  C     I,J,K Line 359  C     I,J,K
359       &      1027.67547  _d 0, 1023.34306  _d 0,       &      1027.67547  _d 0, 1023.34306  _d 0,
360       &      1044.12802  _d 0, 1037.90204  _d 0,       &      1044.12802  _d 0, 1037.90204  _d 0,
361       &      1069.48914  _d 0, 1062.53817  _d 0/       &      1069.48914  _d 0, 1062.53817  _d 0/
362       &     bloc       &     bLoc
363       &     /   -1.00000 _d 0,    -1.00000 _d 0,       &     /   -1.00000 _d 0,    -1.00000 _d 0,
364       &         -1.00000 _d 0,    -1.00000 _d 0,       &         -1.00000 _d 0,    -1.00000 _d 0,
365       &         -1.00000 _d 0,       &         -1.00000 _d 0,
# Line 364  C     I,J,K Line 372  C     I,J,K
372        bi   = 1        bi   = 1
373        bj   = 1        bj   = 1
374        k    = 1        k    = 1
375        imin = 1        iMin = 1
376        imax = 1        iMax = 1
377        jmin = 1        jMin = 1
378        jmax = 1        jMax = 1
379        i    = 1        i    = 1
380        j    = 1        j    = 1
381        IF ( equationOfState.NE.'LINEAR'        IF ( equationOfState.NE.'LINEAR'
382       &     .AND. equationOfState.NE.'POLY3' ) THEN       &     .AND. equationOfState.NE.'POLY3'
383         &     .AND. equationOfState.NE.'IDEALG' ) THEN
384  C     check nonlinear EOS  C     check nonlinear EOS
385           WRITE(msgBuf,'(a,a)')          WRITE(msgBuf,'(A,A)')
386       &        'check_eos: Check the equation of state: Type ',       &        'EOS_CHECK: Check the equation of state: Type ',
387       &        equationOfState       &        equationOfState
388           CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,          CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
389       &        SQUEEZE_RIGHT , myThid )       &       SQUEEZE_RIGHT , myThid )
390           psave = pressure(i,j,k,bi,bj)  C-    to get the right presure out of S/R PRESSURE_FOR_EOS:
391           DO kcheck = 1,ncheck          mskSave = maskC(i,j,k,bi,bj)
392              pressure(i,j,k,bi,bj) = ploc(kcheck)          rC_Save = rC(k)
393              IF ( equationOfState.NE.'UNESCO' ) THEN          maskC(i,j,k,bi,bj) = 1.
394                 tFld(i,j,k,bi,bj) = ptloc(kcheck)          totPhiHyd(i,j,k,bi,bj) = 0.
395              ELSE  C-    + set rC accordingly
396                 tFld(i,j,k,bi,bj) = tloc(kcheck)          DO kcheck = 1,ncheck
397              ENDIF            IF ( usingZCoords ) THEN
398              sFld(i,j,k,bi,bj)    = sloc(kcheck)              rC(k) = -pLoc(kcheck)*recip_rhoConst*recip_gravity
399              rholoc(i,j)          =  0. _d 0            ELSE
400              bulkMod(i,j)         = -1. _d 0              rC(k) =  pLoc(kcheck)
401              ENDIF
402              CALL FIND_RHO(            IF ( equationOfState.NE.'UNESCO' ) THEN
403       &           bi, bj, iMin, iMax, jMin, jMax,  k, k,               tFld(i,j) = ptLoc(kcheck)
404       &           tFld, sFld, rholoc, myThid )            ELSE
405                 tFld(i,j) = tLoc(kcheck)
406              ENDIF
407              sFld(i,j)    = sLoc(kcheck)
408              pFld(i,j)    = pLoc(kcheck)
409              rhoLoc(i,j)  =  0. _d 0
410    
411              CALL FIND_RHO_2D(
412         I           iMin, iMax, jMin, jMax, k,
413         I           tFld, sFld,
414         O           rhoLoc,
415         I           k, bi, bj, myThid )
416    
417              IF ( equationOfState.EQ.'MDJWF' ) THEN
418                bulkMod(i,j) = -1. _d 0
419              ELSE
420              CALL FIND_BULKMOD(              CALL FIND_BULKMOD(
421       &           bi, bj, imin, imax, jmin, jmax, k, k,       I           iMin, iMax, jMin, jMax,
422       &           tFld, sFld, bulkMod, myThid )       I           pFld, tFld, sFld,
423         O           bulkMod,
424              WRITE(msgBuf,       I           myThid )
425       &           '(a4,f4.1,a5,f4.1,a6,f5.0,a5,a3,f10.5,1x,f11.5)')            ENDIF
426       &           'rho(', sFld(i,j,k,bi,bj), ' PSU,',  
427       &           tFld(i,j,k,bi,bj), ' degC,',            WRITE(msgBuf,'(2(A,F4.1),A,F5.0,2A,F10.5,1X,F11.5)')
428  c    &           pressure(i,j,k,bi,bj)*SItoBar, ' bar)',' = ',       &         'rho(', sFld(i,j), ' PSU,',
429       &           rloc(kcheck), bloc(kcheck)       &         tFld(i,j), ' degC,',
430              CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,       &         pLoc(kcheck)*SItoBar, ' bar)',
431       &           SQUEEZE_RIGHT , myThid )       &         ' = ', rLoc(kcheck), bLoc(kcheck)
432              WRITE(msgBuf,'(a14,a22,f10.5,1x,f11.5)')            CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
433       &           'rho(find_rho) ',       &         SQUEEZE_RIGHT , myThid )
434       &           ' = ', rholoc(i,j)+rhoConst, bulkMod(i,j)            WRITE(msgBuf,'(A,F10.5,1X,F11.5)')
435              CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,       &         ' rho(find_rho_2d)     = ',
436       &           SQUEEZE_RIGHT , myThid )       &         rhoLoc(i,j)+rhoConst, bulkMod(i,j)
437              CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
438              CALL FIND_RHO_SCALAR( tFld(i,j,k,bi,bj), sLoc(kcheck),       &         SQUEEZE_RIGHT , myThid )
439       &           pLoc(kcheck), rhoLoc(i,j), myThid )  
440              bulkMod(i,j) = 0. _d 0            CALL FIND_RHO_SCALAR( tFld(i,j), sLoc(kcheck),
441              WRITE(msgBuf,'(a21,a15,f10.5,1x,f11.5)')       &         pLoc(kcheck), rhoLoc(i,j), myThid )
442       &           'rho(find_rho_scalar) ',            WRITE(msgBuf,'(A,F10.5,1X,F11.5)')
443       &           ' = ', rholoc(i,j)+rhoConst, bulkMod(i,j)       &         ' rho(find_rho_scalar) = ',
444              CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,       &         rhoLoc(i,j)+rhoConst
445       &           SQUEEZE_RIGHT , myThid )  c    &         rhoLoc(i,j)+rhoConst, bulkMod(i,j)
446              CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
447         &         SQUEEZE_RIGHT , myThid )
448    
449           ENDDO          ENDDO
450  C     end check nonlinear EOS  C     end check nonlinear EOS
451  c        pressure(i,j,k,bi,bj) = psave          maskC(i,j,k,bi,bj) = mskSave
452            rC(k) = rC_Save
453    
454           WRITE(msgBuf,'(A)') 'end check the equation of state'          WRITE(msgBuf,'(A)') 'EOS_CHECK: Done'
455           CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,          CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
456       &        SQUEEZE_RIGHT , myThid )       &       SQUEEZE_RIGHT , myThid )
457    
458        ENDIF        ENDIF
459  #endif /* EXCLUDE_EOS_CHECK */  #endif /* INCLUDE_EOS_CHECK */
460    
461        RETURN        RETURN
462        END        END
   

Legend:
Removed from v.1.16  
changed lines
  Added in v.1.17

  ViewVC Help
Powered by ViewVC 1.1.22