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

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

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

revision 1.17 by jmc, Sat Aug 9 01:02:28 2008 UTC revision 1.18 by mlosch, Tue Jul 19 12:53:24 2011 UTC
# Line 53  C     == Local variables == Line 53  C     == Local variables ==
53        INTEGER i,j        INTEGER i,j
54        _RL refTemp,refSalt,tP,sP        _RL refTemp,refSalt,tP,sP
55        _RL t1, t2, t3, s1, s3o2, p1, p2, sp5, p1t1        _RL t1, t2, t3, s1, s3o2, p1, p2, sp5, p1t1
56          _RL ct, sa, sqrtsa, p
57        _RL drhoP0dtheta, drhoP0dthetaFresh, drhoP0dthetaSalt        _RL drhoP0dtheta, drhoP0dthetaFresh, drhoP0dthetaSalt
58        _RL dKdtheta, dKdthetaFresh, dKdthetaSalt, dKdthetaPres        _RL dKdtheta, dKdthetaFresh, dKdthetaSalt, dKdthetaPres
59        _RL locPres(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL locPres(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
# Line 267  C     of sea water at p Line 268  C     of sea water at p
268           ENDDO           ENDDO
269        ENDDO        ENDDO
270    
271          ELSEIF ( equationOfState.EQ.'TEOS10' ) THEN
272    
273           CALL PRESSURE_FOR_EOS(
274         I      bi, bj, iMin, iMax, jMin, jMax,  kRef,
275         O      locPres,
276         I      myThid )
277    
278           CALL FIND_RHOTEOS(
279         I      iMin, iMax, jMin, jMax, locPres,
280         I      theta(1-OLx,1-OLy,k,bi,bj), salt(1-OLx,1-OLy,k,bi,bj),
281         O      rhoLoc, rhoDen,
282         I      myThid )
283    
284           DO j=jMin,jMax
285            DO i=iMin,iMax
286             ct      = theta(i,j,k,bi,bj)
287             sa      = salt(i,j,k,bi,bj)
288             IF ( sa .GT. 0. _d 0 ) THEN
289              sqrtsa = SQRT(sa)
290             ELSE
291              sa     = 0. _d 0
292              sqrtsa = 0. _d 0
293             ENDIF
294             p       = locPres(i,j)*SItodBar
295    
296             dnum_dtheta = teos(02)
297         &   + ct*(2.*teos(03) + 3.*teos(04)*ct)  
298         &   + sa*(teos(06) + 2.*teos(07)*ct
299         &   + sqrtsa*(teos(09) + ct*(2.*teos(10) + 3.*teos(11)*ct)))
300         &   + p*(     teos(13) + 2.*teos(14)*ct  + sa*2.*teos(16)
301         &        + p*(teos(18) + 2.*teos(19)*ct))
302            
303             dden_dtheta = teos(22)
304         &   + ct*(2.*teos(23) + ct*(3.*teos(24) + 4.*teos(25)*ct))
305         &   + sa*(teos(27)
306         &   + ct*(2.*teos(28) + ct*(3.*teos(29) + 4.*teos(30)*ct))
307         &   + sqrtsa*(teos(32)
308         &   + ct*(2.*teos(33) + ct*(3.*teos(34) + 4.*teos(35)*ct))))  
309         &   + p*(teos(38) + ct*(2.*teos(39) + 3.*teos(40)*ct)
310         &   + teos(42)
311         &   + p*(teos(44) + 2.*teos(45)*ct + teos(46)*sa
312         &   + p*teos(48) ))
313    
314             alphaLoc(i,j)    = rhoDen(i,j)*(dnum_dtheta
315         &        - (rhoLoc(i,j)*rhoDen(i,j))*dden_dtheta)
316            
317            ENDDO
318           ENDDO
319    
320        ELSE        ELSE
321           WRITE(*,*) 'FIND_ALPHA: equationOfState = ',equationOfState           WRITE(*,*) 'FIND_ALPHA: equationOfState = ',equationOfState
322           STOP 'FIND_ALPHA: "equationOfState" has illegal value'           STOP 'FIND_ALPHA: "equationOfState" has illegal value'
# Line 315  C     == Local variables == Line 365  C     == Local variables ==
365        INTEGER i,j        INTEGER i,j
366        _RL refTemp,refSalt,tP,sP        _RL refTemp,refSalt,tP,sP
367        _RL t1, t2, t3, s1, s3o2, p1, sp5, p1t1        _RL t1, t2, t3, s1, s3o2, p1, sp5, p1t1
368          _RL ct, sa, sqrtsa, p
369        _RL drhoP0dS        _RL drhoP0dS
370        _RL dKdS, dKdSSalt, dKdSPres        _RL dKdS, dKdSSalt, dKdSPres
371        _RL locPres(1-Olx:sNx+Olx,1-Oly:sNy+Oly)        _RL locPres(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
# Line 506  C     of sea water at p Line 557  C     of sea water at p
557              ENDDO              ENDDO
558           ENDDO           ENDDO
559    
560          ELSEIF ( equationOfState.EQ.'TEOS10' ) THEN
561    
562           CALL PRESSURE_FOR_EOS(
563         I      bi, bj, iMin, iMax, jMin, jMax,  kRef,
564         O      locPres,
565         I      myThid )
566    
567           CALL FIND_RHOTEOS(
568         I      iMin, iMax, jMin, jMax, locPres,
569         I      theta(1-OLx,1-OLy,k,bi,bj), salt(1-OLx,1-OLy,k,bi,bj),
570         O      rhoLoc, rhoDen,
571         I      myThid )
572          
573           DO j=jMin,jMax
574            DO i=iMin,iMax
575             ct      = theta(i,j,k,bi,bj)
576             sa      = salt(i,j,k,bi,bj)
577             IF ( sa .GT. 0. _d 0 ) THEN
578              sqrtsa = SQRT(sa)
579             ELSE
580              sa     = 0. _d 0
581              sqrtsa = 0. _d 0
582             ENDIF
583             p       = locPres(i,j)*SItodBar
584    
585             dnum_dsalt =  teos(05) + ct*(teos(06) + teos(07)*ct)
586         &   + 1.5*sqrtsa*(teos(08)
587         &           + ct*(teos(09) + ct*(teos(10) + teos(11)*ct)))
588         &   + p*(teos(15) + teos(16)*ct + p*teos(20))
589    
590             dden_dsalt = teos(26)
591         &   + ct*(teos(27) + ct*(teos(28) + ct*(teos(29) + teos(30)*ct)))
592         &   + 2.*teos(36)*sa
593         &   + 1.5*sqrtsa*(teos(31) + ct*(teos(32) + ct*(teos(33)
594         &                          + ct*(teos(34) + teos(35)*ct))))
595         &   + p*(teos(41) + teos(42)*ct + p*teos(46))
596    
597             betaLoc(i,j) = rhoDen(i,j)*( dnum_dsalt
598         &        - (rhoLoc(i,j)*rhoDen(i,j))*dden_dsalt )
599    
600            ENDDO
601           ENDDO
602    
603        ELSE        ELSE
604           WRITE(*,*) 'FIND_BETA: equationOfState = ',equationOfState           WRITE(*,*) 'FIND_BETA: equationOfState = ',equationOfState
605           STOP 'FIND_BETA: "equationOfState" has illegal value'           STOP 'FIND_BETA: "equationOfState" has illegal value'
Legend:
Removed from v.1.17  
changed lines
  Added in v.1.18
  ViewVC Help
Powered by ViewVC 1.1.22