pkg/frazil/frazil_calc_rhs.F

e0b3e1bdd8 Dimi*0001 #include "FRAZIL_OPTIONS.h"
                0002
                0003 CBOP
                0004 C !ROUTINE: FRAZIL_CALC_RHS
                0005
                0006 C !INTERFACE: ==========================================================
                0007       SUBROUTINE FRAZIL_CALC_RHS(
                0008      I                     myTime, myIter, myThid )
                0009
                0010 C !DESCRIPTION:
                0011 C Check water temperature and if colder than freezing
                0012 C point bring excess negative heat to the surface.
                0013
                0014 C !USES: ===============================================================
                0015       IMPLICIT NONE
                0016 #include "EEPARAMS.h"
33e5e9881a Dimi*0017 #include "SIZE.h"
e0b3e1bdd8 Dimi*0018 #include "PARAMS.h"
                0019 #include "GRID.h"
33e5e9881a Dimi*0020 #include "DYNVARS.h"
e0b3e1bdd8 Dimi*0021 #include "FFIELDS.h"
33e5e9881a Dimi*0022 #ifdef ALLOW_FRAZIL
                0023 # include "FRAZIL.h"
                0024 #endif
e0b3e1bdd8 Dimi*0025
                0026 C     !INPUT/OUTPUT PARAMETERS:
33e5e9881a Dimi*0027 C     myTime :: current time in simulation
                0028 C     myIter :: current iteration number in simulation
                0029 C     myThid :: my Thread Id number
                0030       _RL     myTime
e0b3e1bdd8 Dimi*0031       INTEGER myIter
                0032       INTEGER myThid
33e5e9881a Dimi*0033 CEOP
e0b3e1bdd8 Dimi*0034
                0035 #ifdef ALLOW_FRAZIL
33e5e9881a Dimi*0036
                0037 C     !LOCAL VARIABLES:
                0038 C     Tfreezing :: freezing threshold temperature
e0b3e1bdd8 Dimi*0039       INTEGER bi,bj,i,j,k,kTop
                0040       _RL Tfreezing, Tresid, pLoc, sLoc, tLoc
                0041       _RL a0, a1, a2, b
                0042       PARAMETER( a0 = -0.0575   _d  0 )
                0043       PARAMETER( a1 = 1.710523  _d -3 )
                0044       PARAMETER( a2 = -2.154996 _d -4 )
                0045       PARAMETER( b  = -7.53     _d -4 )
                0046
                0047       _RL SW_TEMP
                0048       EXTERNAL SW_TEMP
                0049
                0050       DO bj=myByLo(myThid),myByHi(myThid)
                0051        DO bi=myBxLo(myThid),myBxHi(myThid)
33e5e9881a Dimi*0052
                0053 C     Initialize FrazilForcingT to zero
                0054         DO k=1,Nr
                0055          DO j=1-Oly,sNy+OLy
                0056           DO i=1-Olx,sNx+Olx
                0057            FrazilForcingT(i,j,k,bi,bj) = 0. _d 0
                0058           ENDDO
                0059          ENDDO
                0060         ENDDO
                0061
                0062 C     Check for water below freezing point.
                0063         DO k = 2, Nr
                0064          DO j=1-OLy,sNy+OLy
                0065           DO i=1-OLx,sNx+OLx
e0b3e1bdd8 Dimi*0066            IF ( maskC(i,j,k-1,bi,bj) .NE. 0. _d 0 .AND.
                0067      &          maskC(i,j,k,  bi,bj) .NE. 0. _d 0 ) THEN
                0068
                0069             pLoc = ABS(RC(k))
                0070             sLoc = MAX(salt(i,j,k,bi,bj), 0. _d 0)
                0071             tLoc = SW_TEMP(sLoc,theta(i,j,k,bi,bj),pLoc,0. _d 0)
                0072
                0073 C Freezing point of seawater
                0074 C   REFERENCE: UNESCO TECH. PAPERS IN THE MARINE SCIENCE NO. 28. 1978
                0075 C   EIGHTH REPORT JPOTS
                0076 C   ANNEX 6 FREEZING POINT OF SEAWATER F.J. MILLERO PP.29-35.
                0077 C
                0078 C  UNITS:
                0079 C         PRESSURE      P          DECIBARS
                0080 C         SALINITY      S          PSS-78
                0081 C         TEMPERATURE   TF         DEGREES CELSIUS
                0082 C         FREEZING PT.
                0083 C************************************************************
                0084 C  CHECKVALUE: TF= -2.588567 DEG. C FOR S=40.0, P=500. DECIBARS
                0085             Tfreezing = (a0 + a1*sqrt(sLoc) + a2*sLoc) * sLoc + b*pLoc
                0086
                0087             IF (tLoc .LT. Tfreezing) THEN
                0088 C     Move the negative heat to surface level.
33e5e9881a Dimi*0089              kTop = kSurfC(i,j,bi,bj)
                0090              Tresid = ( Tfreezing - tloc )
610e3cd50d Dimi*0091      &            * HeatCapacity_Cp * rUnit2mass
                0092      &            * drF(k) * _hFacC(i,j,k,bi,bj)
33e5e9881a Dimi*0093              FrazilForcingT(i,j,k,bi,bj) = Tresid / dTtracerLev(k)
610e3cd50d Dimi*0094              FrazilForcingT(i,j,kTop,bi,bj) =
33e5e9881a Dimi*0095      &            FrazilForcingT(i,j,kTop,bi,bj)
610e3cd50d Dimi*0096      &            - Tresid / dTtracerLev(kTop)
e0b3e1bdd8 Dimi*0097             ENDIF
                0098            ENDIF
                0099           ENDDO
                0100          ENDDO
                0101         ENDDO
                0102        ENDDO
                0103       ENDDO
                0104
bb86b07795 Dimi*0105 # ifdef ALLOW_DIAGNOSTICS
                0106       IF ( useDiagnostics ) THEN
                0107          CALL DIAGNOSTICS_FILL( FrazilForcingT, 'FrzForcT',
                0108      &                          0,Nr, 0, 1, 1, myThid )
                0109       ENDIF
                0110 # endif /* ALLOW_DIAGNOSTICS */
                0111
e0b3e1bdd8 Dimi*0112 #endif /* ALLOW_FRAZIL */
                0113
                0114       RETURN
                0115       END