New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
diainsitutem.F90 in branches/UKMO/r6232_collate_bgc_diagnostics/NEMOGCM/NEMO/OPA_SRC/DIA – NEMO

source: branches/UKMO/r6232_collate_bgc_diagnostics/NEMOGCM/NEMO/OPA_SRC/DIA/diainsitutem.F90 @ 11137

Last change on this file since 11137 was 11137, checked in by jcastill, 5 years ago

Add missing files

File size: 4.7 KB
Line 
1MODULE diainsitutem
2   !!======================================================================
3   !!                       ***  MODULE  diaharm  ***
4   !! Harmonic analysis of tidal constituents
5   !!======================================================================
6   !! Initial version: D. Acreman June 2006 Original Code
7   !! History :  3.6  !  2014  (E O'Dea) for 3.6 25hour mean output
8   !!----------------------------------------------------------------------
9
10  USE dom_oce  ! ocean space and time domain
11  USE oce,            ONLY: tsn
12  USE par_oce,        ONLY: jpi, jpj, jpk, jpkm1
13  USE lbclnk,         ONLY: lbc_lnk
14  USE lib_mpp         ! MPP library
15
16  IMPLICIT NONE
17  PRIVATE
18
19  PUBLIC theta2t
20
21  !! * Accessibility
22  PUBLIC insitu_tem_alloc          ! routines called by step.F90
23
24  REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:) :: rinsitu_t
25
26   !! * Substitutions
27#  include "domzgr_substitute.h90"
28
29CONTAINS
30
31   REAL FUNCTION insitu_tem_alloc()
32      !!----------------------------------------------------------------------
33      INTEGER, DIMENSION(2) :: ierr
34      !!----------------------------------------------------------------------
35      !
36      ierr = 0
37      !
38      ALLOCATE( rinsitu_t(jpi,jpj,jpk), STAT=ierr(1) )
39         !
40      insitu_tem_alloc = MAXVAL(ierr)
41      IF( lk_mpp )   CALL mpp_sum( insitu_tem_alloc )
42      !
43  END FUNCTION insitu_tem_alloc
44
45!-----------------------------------------------------------------------------------
46!
47! Calculate the insitu temperature by integrating the adiabatic lapse rate from zero
48! to pressure at depth of tracer level. Based on UM subroutine POTTEM and function ATG
49!
50! Initial version: D. Acreman June 2006
51!
52!-----------------------------------------------------------------------------------
53
54  SUBROUTINE theta2t()
55
56    INTEGER, PARAMETER :: inum_steps=10                ! number of steps in integration
57    INTEGER            :: jstep                        ! iteration counter
58    INTEGER            :: ji, jj, jk                  ! loop indices
59    REAL(wp), DIMENSION(jpi,jpj,jpk) :: zP            ! pressure (decibars)
60    REAL(wp), DIMENSION(jpi,jpj,jpk) :: zT            ! temperature at pressure p
61    REAL(wp), DIMENSION(jpi,jpj,jpk) :: zTB           ! temperature at p-dp
62    REAL(wp), DIMENSION(jpi,jpj,jpk) :: zTA           ! temperature at p+dp
63    REAL(wp), DIMENSION(jpi,jpj,jpk) :: zDP           ! pressure step
64    REAL(wp), DIMENSION(jpi,jpj,jpk) :: zSS           ! salinity(PSU) - 35.0
65    REAL(wp), DIMENSION(jpi,jpj,jpk) :: zLAPSE        ! adiabatic lapse rate
66
67!CDIR IEXPAND (ATG)
68
69! Set the pressure interval for the integration. The integration is carried out from
70! zero (pressure at the surface) to the pressure at the depth of the tracer point. The
71! pressure in decibars is represented by the depth in metres. Pressures are "Oceanographic"
72! pressures equal to absolute pressure minus one atmosphere
73     zDP(:,:,:) = 0.0
74     DO jk = 1, jpkm1
75        DO jj = 1, jpj
76           DO ji = 1, jpi
77              ! These loops expanded for case where fsdept may be 1D
78              zDP(ji,jj,jk) = fsdept(ji,jj,jk) / real(inum_steps)
79           END DO
80        END DO
81     END DO
82
83! Salinity at each point
84     zSS(:,:,:)      = tsn(:,:,:,jp_sal) - 35.0 
85
86! Set initial values of temperature and pressure. zT is the temperature at pressure zP,
87! zTB is the temperature at pressure zP-zdP and zTA is the temperature at pressure zP+zdP
88     zT(:,:,:)  = tsn(:,:,:,jp_tem)
89     zP(:,:,:)  = 0.0         ! Pressure at surface
90     CALL ATG(zP, zT, zSS, zLAPSE)
91     zTB(:,:,:) = zT(:,:,:)  - zLAPSE(:,:,:) * zDP(:,:,:)
92
93     interation: DO jstep=1, inum_steps
94        ! Calculate lapse rate (dT/dP) and hence TA
95        CALL ATG(zP, zT, zSS, zLAPSE)
96        zTA(:,:,:) = zTB(:,:,:) + 2.0 * zLAPSE(:,:,:) * zDP(:,:,:)
97        ! Have calculated TB, T and TA for this pressure, now advance solution by dP
98        zP(:,:,:)  = zP(:,:,:)  + zDP
99        zTB(:,:,:) = zT(:,:,:)
100        zT(:,:,:)  = zTA(:,:,:)
101
102     END DO interation
103
104     rinsitu_t(:,:,:) = zT(:,:,:) * tmask(:,:,:)
105     CALL lbc_lnk( rinsitu_t,  'T', 1.0)
106
107   END SUBROUTINE theta2t
108
109   SUBROUTINE ATG(P,T,DS,zLAPSE)
110
111     REAL, INTENT(IN   ) :: P(jpi,jpj,jpk)        ! PRESSURE (DECIBARS)
112     REAL, INTENT(IN   ) :: T(jpi,jpj,jpk)        ! TEMPERATURE (DEG C)
113     REAL, INTENT(IN   ) :: DS(jpi,jpj,jpk)       ! SALINITY (PSU) -35.0
114     REAL, INTENT(  OUT) :: zLAPSE(jpi,jpj,jpk)   ! LAPSE RATE
115
116       zLAPSE = ((( -2.1687E-16*T+1.8676E-14)*T-4.6206E-13)*P &
117           + ((2.7759E-12*T-1.1351E-10)*DS+((-5.4481E-14*T    &
118           + 8.733E-12)*T-6.7795E-10)*T+1.8741E-8))*P         &
119           + (-4.2393E-8*T+1.8932E-6)*DS                      &
120           + ((6.6228E-10*T-6.836E-8)*T+8.5258E-6)*T+3.5803E-5 
121
122   END SUBROUTINE ATG
123
124END MODULE diainsitutem
Note: See TracBrowser for help on using the repository browser.