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