*IF DEF,OCEAN @DYALLOC.4660
C ******************************COPYRIGHT****************************** GTS2F400.8785
C (c) CROWN COPYRIGHT 1995, METEOROLOGICAL OFFICE, All Rights Reserved. GTS2F400.8786
C GTS2F400.8787
C Use, duplication or disclosure of this code is subject to the GTS2F400.8788
C restrictions as set forth in the contract. GTS2F400.8789
C GTS2F400.8790
C Meteorological Office GTS2F400.8791
C London Road GTS2F400.8792
C BRACKNELL GTS2F400.8793
C Berkshire UK GTS2F400.8794
C RG12 2SZ GTS2F400.8795
C GTS2F400.8796
C If no contract has been raised with this copy of the code, the use, GTS2F400.8797
C duplication or disclosure of it is strictly prohibited. Permission GTS2F400.8798
C to do so must first be obtained in writing from the Head of Numerical GTS2F400.8799
C Modelling at the above address. GTS2F400.8800
C ******************************COPYRIGHT****************************** GTS2F400.8801
C GTS2F400.8802
C*LL SFCADD.2
CLL Subroutine SFCADD SFCADD.3
CLL Can run on any FORTRAN 77 compiler with long lower case variables SFCADD.4
CLL SFCADD.5
CLL The code must be pre-compiled by the UPDOC system. SFCADD.6
CLL Option A indicates that the Unified Model version is to be used. SFCADD.7
CLL The default is the COX standard code. SFCADD.8
CLL SFCADD.9
CLL Author: S J Foreman SFCADD.10
CLL Date: 17 September 1990 SFCADD.11
CLL Reviewer: D L Roberts SFCADD.12
CLL Review date: 17 September 1990 SFCADD.13
CLL Version 1.01 date 17 September 1990 SFCADD.14
CLL Version 3.4 4/8/94 Remove ice depth flux correction and separate OJT0F304.66
CLL ice depth haney forcing from SST/SSS. (JFT) OJT0F304.67
CLL Version 4.0 7.7.95 Take heat from layer 2 if flux correction OJG0F400.31
CLL causes layer 1 to freeze (J. Murphy, coded by J. Gregory) OJG0F400.32
! 3.5 16.01.95 Remove *IF dependency. R.Hill ORH1F305.400
CLL 4.1 22.5.96 J.M.Gregory If heat flux correction cannot be OJG0F401.16
CLL applied because it would freeze layer 1, record it in HEATSINK OJG0F401.17
CLL SFCADD.15
CLL Programming standards use Cox naming convention for Cox variables SFCADD.16
CLL with the addition that lower case variables are local to the SFCADD.17
CLL routine. SFCADD.18
CLL Otherwise follows UM doc paper 4 version 1. SFCADD.19
CLL SFCADD.20
CLL This forms part of UM brick P4. SFCADD.21
CLL SFCADD.22
CLL This routine calculates the change in temperature and salinity SFCADD.23
CLL due to non-penetrating surface fluxes. SFCADD.24
CLL SFCADD.25
CLL Programming: Brick P4, paper 3, version number 1. SFCADD.26
CLL SFCADD.27
CLL Subroutine dependencies: NONE SFCADD.28
CLL SFCADD.29
CLLEND --------------------------------------------------------------- SFCADD.30
C* SFCADD.31
C*L -------------------------- Arguments ---------------------------- SFCADD.32
C SFCADD.33
SUBROUTINE SFCADD (TA, QFLUX, SFLUX, SREF, 2SFCADD.34
& IMT, KM, NT, ORH1F305.401
+ DZ, C2DTTS, SPECIFIC_HEAT,RHO_WATER SFCADD.36
&, fluxcorh, fluxcorw ORH1F305.402
&,TFREEZE,FM OJG0F400.33
&, anom_ice_heat, icy ORH1F305.403
&,VTCO2_FLUX,VALK_FLUX,TB OJP0F404.949
&,heatsink OJG0F401.13
+ ,salref) OJL1F405.1
C SFCADD.45
IMPLICIT NONE SFCADD.46
C SFCADD.47
*CALL CNTLOCN 
ORH1F305.404
*CALL OARRYSIZ ORH1F305.405
*CALL OTRACPNT ! Required for carbon model tracer pointers OJP0F404.950
C Define constants for array sizes SFCADD.48
C SFCADD.49
INTEGER SFCADD.50
+ IMT ! IN Number of points in horizontal SFCADD.51
+, KM ! IN Number of layers in model SFCADD.52
+, NT ! IN Number of tracers SFCADD.53
C SFCADD.54
C Physical arguments SFCADD.55
C SFCADD.56
LOGICAL ORH1F305.406
& icy(IMT_IHY) ! IN TRUE at seaice points ORH1F305.407
ORH1F305.408
REAL SFCADD.65
+ TA (IMT, KM, NT) ! INOUT Tracer values SFCADD.66
&, TB (IMT, KM, NT) ! IN Tracer values OJP0F404.951
+, QFLUX (IMT) ! IN Non-penetrating heat flux (SI) SFCADD.67
+, SFLUX (IMT) ! IN Precip less evap (SI) SFCADD.68
+, SREF (IMT) ! IN Reference salinity SFCADD.69
+, DZ (KM) ! IN Layer thicknesses (cm) SFCADD.70
+, C2DTTS ! IN Timestep (S) SFCADD.71
+, SPECIFIC_HEAT ! IN Specific heat capacity (SI) SFCADD.72
+, RHO_WATER ! IN density of sea water (SI) SFCADD.73
&, fluxcorh (IMT_FLX) ! IN Heat flux correction -sea (SI) ORH1F305.409
&, fluxcorw (IMT_FLX) ! IN Salt flux correction (SI) ORH1F305.410
&,anom_ice_heat(IMT_IHY)! IN Anomalous seaice heat flux (SI) ORH1F305.411
+, FM(IMT,KM) ! IN Masking array for tracer points OJG0F400.34
+, TFREEZE ! IN Freezing point of sea water, in Celsius OJG0F400.35
&, VTCO2_FLUX(IMT),VALK_FLUX(IMT) ! OUT Virtual Air-sea fluxes of OJP0F404.952
& ! TCO2 and ALK driven by PME. OJP0F404.953
+,HEATSINK(IMT) ! OUT Heat 'lost' when T set to TFREEZE at bottom OJG0F401.14
+ ! model level OJG0F401.15
C* SFCADD.82
real salref OJL1F405.2
C SFCADD.83
C Locally defined variables SFCADD.84
C SFCADD.85
INTEGER SFCADD.86
+ i ! Horizontal loop index SFCADD.87
+, inc ! Indicator for applying flux correction OJG0F400.36
C SFCADD.88
REAL SFCADD.89
+ con_temp ! Conversion W/m2 to temp SFCADD.90
+, con_salt ! Conversion P-E (kg/m2/s) to salt SFCADD.91
&, fluxtoth (IMT) ! Total heat flux correction ORH1F305.412
+, deltat ! IN SST increment arising from flux corr OJG0F400.37
+, deltat1 ! IN T increment to be transferred to level 2 OJG0F400.38
+, tval ! IN Temporary store for value of TA(i,1,1) OJG0F400.39
C SFCADD.98
CL 1.1 Define conversion factors SFCADD.99
C SFCADD.100
C Heat SI to temperature change SFCADD.101
C SFCADD.102
con_temp = C2DTTS/DZ(1)*(100.0/SPECIFIC_HEAT)/RHO_WATER SFCADD.103
con_salt = -100.*(C2DTTS/DZ(1))/RHO_WATER SFCADD.104
SFCADD.105
IF (L_FLUXCORR) THEN ORH1F305.413
C SFCADD.107
C At non-ice points, add sea-ice heat flux correction to Haney SFCADD.108
C term to obtain total heat flux correction. SFCADD.109
C SFCADD.110
DO 2000 I=1,IMT SFCADD.111
if (L_SALFLUXFIX) then OJL1F405.13
fluxtoth(I)=0.0 OJL1F405.14
else OJL1F405.15
fluxtoth(I)=fluxcorh(I) OJL1F405.16
endif OJL1F405.17
2000 CONTINUE SFCADD.117
ENDIF ORH1F305.414
IF (L_SEAICE.AND.L_IHANEY) THEN ORH1F305.415
C SFCADD.120
C At non-ice points, add anomalous seaice heat flux in the manner SFCADD.121
C of a flux-correction term to surface layer of ocean SFCADD.122
C SFCADD.123
DO I = 1, IMT ORH1F305.416
IF (.NOT.icy(I)) THEN SFCADD.125
fluxtoth(I)=anom_ice_heat(I) SFCADD.126
ELSE SFCADD.127
fluxtoth(I)=0.0 SFCADD.128
ENDIF SFCADD.129
ENDDO ! over I ORH1F305.417
ENDIF ORH1F305.418
C SFCADD.132
CL 2.1 Add increments SFCADD.133
DO I = 1, IMT ORH1F305.419
TA(i,1,1) = TA(i,1,1) + con_temp*QFLUX(i) ORH1F305.420
if (L_REFSAL) then OJL1F405.3
TA(i,1,2)=TA(i,1,2)+con_salt*SFLUX(i)*salref OJL1F405.4
else OJL1F405.5
TA(i,1,2) = TA(i,1,2)+con_salt*SFLUX(i)*(0.035+SREF(i)) OJL1F405.6
endif OJL1F405.7
IF (L_OCARBON) THEN OJP0F404.954
VTCO2_FLUX(i)= OJP0F404.955
& con_salt*SFLUX(i)*TB(i,1,TCO2_TRACER)/C2DTTS OJP0F404.956
TA(i,1,TCO2_TRACER)= OJP0F404.957
& TA(i,1,TCO2_TRACER)+VTCO2_FLUX(i)*C2DTTS OJP0F404.958
ENDIF OJP0F404.959
IF (L_OCARBON) THEN OJP0F404.960
VALK_FLUX(i)= OJP0F404.961
& con_salt*SFLUX(i)*TB(i,1,ALK_TRACER)/C2DTTS OJP0F404.962
TA(i,1,ALK_TRACER)= OJP0F404.963
& TA(i,1,ALK_TRACER)+VALK_FLUX(i)*C2DTTS OJP0F404.964
ENDIF OJP0F404.965
ENDDO ORH1F305.422
ORH1F305.423
IF (L_FLUXCORR) THEN ORH1F305.424
! Add in heat and water flux corrections ORH1F305.425
DO I = 1, IMT ORH1F305.426
C Apply heat flux correction to layer 1 OJG0F400.40
C If application of flux correction term causes SST to fall below OJG0F400.41
C freezing point of ice, reset SST to TFREEZE by taking heat from OJG0F400.42
C the layer below OJG0F400.43
deltat=con_temp*fluxtoth(I) OJG0F400.44
tval=TA(i,1,1) OJG0F400.45
IF(deltat.ge.0.0) THEN OJG0F400.46
inc=1 OJG0F400.47
ELSEIF ((tval+deltat).ge.TFREEZE) THEN OJG0F400.48
inc=1 OJG0F400.49
ELSEIF (tval.le.TFREEZE) THEN OJG0F400.50
inc=2 OJG0F400.51
ELSE OJG0F400.52
inc=3 OJG0F400.53
ENDIF OJG0F400.54
C Case 1: flux correction leaves layer 1 above freezing. This must be OJG0F400.55
C the case if the flux correction is positive, because upon entry layer OJG0F400.56
C 1 should never be below freezing. OJG0F400.57
IF(inc.eq.1) THEN OJG0F400.58
TA(i,1,1)=tval+deltat OJG0F400.59
C Case 2: Layer 1 was freezing and flux correction would leave it below. OJG0F400.60
C In this case, we apply the flux correction to layer 2 instead. OJG0F400.61
ELSEIF(inc.eq.2) THEN OJG0F400.62
if (fm(i,2).eq.0.) then OJG0F401.18
heatsink(i)=heatsink(i)+fluxtoth(i) OJG0F401.19
else OJG0F401.20
TA(i,2,1)=TA(i,2,1)+(DZ(1)/DZ(2))*deltat OJG0F401.21
endif OJG0F401.22
C Case 3: Layer 1 was above freezing but flux correction freezes it. OJG0F400.64
C We leave layer 1 at freezing and take the rest of the flux correction OJG0F400.65
C from layer 2. OJG0F400.66
ELSE OJG0F400.67
deltat1=tval-TFREEZE+deltat OJG0F400.68
TA(i,1,1)=TFREEZE OJG0F400.69
if (fm(i,2).eq.0.) then OJG0F401.23
heatsink(i)=heatsink(i)+deltat1/con_temp OJG0F401.24
else OJG0F401.25
TA(i,2,1)=TA(i,2,1)+(DZ(1)/DZ(2))*deltat1 OJG0F401.26
endif OJG0F401.27
ENDIF OJG0F400.71
if (L_REFSAL) then OJL1F405.8
TA(i,1,2)=TA(i,1,2)+con_salt*fluxcorw(I)*salref OJL1F405.9
else OJL1F405.10
TA(i,1,2)=TA(i,1,2)+con_salt*fluxcorw(I)*(0.035+SREF(I)) OJL1F405.11
endif OJL1F405.12
ENDDO ORH1F305.429
ENDIF ORH1F305.430
ORH1F305.431
IF (L_IHANEY) THEN ORH1F305.432
! Add anomalous seaice heat flux where ice is absent ORH1F305.433
DO I = 1, IMT ORH1F305.434
TA(i,1,1) = TA(i,1,1) + con_temp*fluxtoth(I) ORH1F305.435
ENDDO ORH1F305.436
ENDIF ORH1F305.437
C SFCADD.154
CL Return from SFCADD SFCADD.155
C SFCADD.156
RETURN SFCADD.157
END SFCADD.158
*ENDIF @DYALLOC.4661