*IF DEF,S40_1A UMSLAB1A.2
C ******************************COPYRIGHT****************************** GTS2F400.10783
C (c) CROWN COPYRIGHT 1995, METEOROLOGICAL OFFICE, All Rights Reserved. GTS2F400.10784
C GTS2F400.10785
C Use, duplication or disclosure of this code is subject to the GTS2F400.10786
C restrictions as set forth in the contract. GTS2F400.10787
C GTS2F400.10788
C Meteorological Office GTS2F400.10789
C London Road GTS2F400.10790
C BRACKNELL GTS2F400.10791
C Berkshire UK GTS2F400.10792
C RG12 2SZ GTS2F400.10793
C GTS2F400.10794
C If no contract has been raised with this copy of the code, the use, GTS2F400.10795
C duplication or disclosure of it is strictly prohibited. Permission GTS2F400.10796
C to do so must first be obtained in writing from the Head of Numerical GTS2F400.10797
C Modelling at the above address. GTS2F400.10798
C ******************************COPYRIGHT****************************** GTS2F400.10799
C GTS2F400.10800
C*LL UMSLAB1A.3
CLL SUBROUTINE UMSLAB UMSLAB1A.4
CLL ----------------- UMSLAB1A.5
CLL UMSLAB1A.6
CLL THIS ROUTINE IS FOR USE WITH THE 'SLAB' OCEAN MODEL ONLY. UMSLAB1A.7
CLL UMSLAB1A.8
CLL THE TEMPERATURE OF THE MIXED LAYER OCEAN IS UPDATED UMSLAB1A.9
CLL ACCORDING TO THE EQUATION UMSLAB1A.10
CLL UMSLAB1A.11
CLL CW * RHOW * DTM NET HEAT FLUX IN + HEAT CONVERGENCE UMSLAB1A.12
CLL --- = ------------------------------------ UMSLAB1A.13
CLL DT H UMSLAB1A.14
CLL UMSLAB1A.15
CLL THIS ROUTINE FORMS PART OF SYSTEM COMPONENT P40. UMSLAB1A.16
CLL IT CAN BE COMPILED BY CFT77, BUT DOES NOT CONFORM TO ANSI UMSLAB1A.17
CLL FORTRAN77 STANDARDS, BECAUSE OF THE INLINE COMMENTS. UMSLAB1A.18
CLL IT ADHERES TO THE STANDARDS OF DOCUMENTATION PAPER 3, VERSION 5. UMSLAB1A.19
CLL UMSLAB1A.20
CLL ALL QUANTITIES IN THIS ROUTINE ARE IN S.I. UNITS UNLESS UMSLAB1A.21
CLL OTHERWISE STATED. UMSLAB1A.22
CLL UMSLAB1A.23
CLL DIFFUSION OF SLAB TEMPERATURES ADDED SJT1F304.823
CLL SJT1F304.824
CLL CALLED BY: SLABCNTL UMSLAB1A.24
CLL UMSLAB1A.25
CLL WRITTEN BY C.A.SENIOR (14/1/91) UMSLAB1A.26
CLL MODIFIED BY A.B.KEEN (02/02/93) UMSLAB1A.27
CLL MODIFIED BY C.A.SENIOR(22/03/93) SJT1F304.825
CLL MODIFIED BY C.A.SENIOR(14/09/93) SJT1F304.826
CLL MODIFIED BY J.F.THOMSON (23/05/94) Thermodynamics altered to SJT1F304.827
CLL allow addition of dynamics. SJT1F304.828
CLL Version Description of change SJC1F400.97
CLL 4.0 Vertical SST advection included (R.Carnell) and SJC1F400.98
CLL change order of heat convergence calculation for SJC1F400.99
CLL ice dynamics. J.F.Crossley SJC1F400.100
!LL 4.4 04/08/97 Add missing ARGOINDX to various argument lists. SDR1F404.97
!LL D. Robinson. SDR1F404.98
!LL SDR1F404.99
CLL VERSION NUMBER 1.1 UMSLAB1A.28
CLL REVIEWER: W.INGRAM (01/03/93) UMSLAB1A.29
CLL UMSLAB1A.30
CLLEND--------------------------------------------------------------- UMSLAB1A.31
C*L UMSLAB1A.32
SUBROUTINE UMSLAB( 1,2SDR1F404.100
*CALL ARGOINDX 
SDR1F404.101
+ ICY, SDR1F404.102
+ ATMSFLUX, SJT1F304.830
+ ADJHCONV, SJT1F304.831
+ SNOWRATE, SJT1F304.832
+ CARYHEAT, SJT1F304.833
+ SLABTEMP, SJT1F304.834
+ wtsfc, SJC1F400.101
+ wtbase, SJC1F400.102
+ LAND, SJT1F304.835
+ AICE, SJT1F304.836
+ NEWICE, SJT1F304.837
+ OIFLUX, SJT1F304.838
+ DTADV,DTDIFF, SJC1F400.103
+ L1,L2,DT,DZ1,u_field, SJC1F400.104
+ JROWS,JROWSM1,ICOLS, SJT1F304.840
+ TCLIMC, SJT1F304.841
+ HCLIM, SJT1F304.842
+ CALIB, SJT1F304.843
+ L_THERM,L_IDYN,L_IDRIF, SJT1F304.844
+ L_SLBADV,LGLOBAL, SJT1F304.845
+ EDDYDIFF, SJT1F304.846
+ DELTA_LONG,DELTA_LAT,BASE_LAT, SJT1F304.847
+ AH, SJT1F304.848
+ OPENSEA, SJT1F304.849
+ UCURRENT, SJT1F304.850
+ VCURRENT, SJT1F304.851
+ COS_P_LATITUDE,COS_U_LATITUDE, SJT1F304.852
+ SEC_P_LATITUDE,SIN_U_LATITUDE) SJT1F304.853
C UMSLAB1A.36
IMPLICIT NONE UMSLAB1A.37
C UMSLAB1A.38
*CALL TYPOINDX SDR1F404.103
INTEGER L1, ! IN SIZE OF DATA VECTORS UMSLAB1A.39
+ L2 ! IN AMOUNT OF DATA TO BE PROCESSED UMSLAB1A.40
+,JROWS ! IN NO OF ROWS N-S SJT1F304.854
+,JROWSM1 ! IN NO OF ROWS N-S - 1 SJT1F304.855
+,ICOLS ! IN NO OF COLUMNS E-W SJT1F304.856
+,u_field ! In points in u field SJC1F400.105
C UMSLAB1A.41
REAL ATMSFLUX(L1) ! IN NET DOWNWARD SURFACE HEAT FLUX (W M-2) SJT1F304.857
+,ADJHCONV(L1) ! IN ADJUSTED HEAT CONVERGENCE IN W M-2. SJT1F304.858
+,SNOWRATE(L1) ! IN RATE OF SNOWFALL, IN KG M-2 S-1. SJT1F304.859
+,TCLIMC(L1) ! IN CLIMATOLOGICAL SEA SURFACE TEMPS C SJT1F304.860
+,HCLIM(L1) ! IN CLIMATOLOGICAL SEA-ICE EXTENTS M SJT1F304.861
C UMSLAB1A.48
REAL SJT1F304.862
+ COS_P_LATITUDE(L1)! IN COS LATITUDE ON P GRID SJT1F304.863
+,COS_U_LATITUDE(L1)! IN COS LATITUDE ON UV GRID SJT1F304.864
+,SEC_P_LATITUDE(L1)! 1.0/COS LATITUDE ON P GRID SJT1F304.865
+,SIN_U_LATITUDE(L1)! IN SIN LATITUDE ON UV GRID SJT1F304.866
C UMSLAB1A.52
REAL DT ! IN TIMESTEP FOR UPDATING THE SLAB OCEAN IN S. SJT1F304.867
+,DZ1 ! IN THICKNESS OF THE SLAB OCEAN IN METRES. SJT1F304.868
+,AH ! DIFFUSION COEFFICENT SJT1F304.869
+,EDDYDIFF ! IN EDDY DIFFUSION COEFF FOR OIFLUX CALCULATION SJT1F304.870
+,UCURRENT(ICOLS,JROWSM1) ! IN ZONAL SURFACE CURRENT (M/S) SJT1F304.871
+,VCURRENT(ICOLS,JROWSM1) ! IN MERIDIONAL SURFACE CURRENT (M/S) SJT1F304.872
C SJT1F304.873
REAL SJT1F304.874
& DELTA_LONG ! IN EW grid spacing (degrees) SJT1F304.875
&,DELTA_LAT ! IN NS grid spacing (degrees) SJT1F304.876
&,BASE_LAT ! IN latitude of first row (degrees) SJT1F304.877
C SJT1F304.878
LOGICAL ICY(L1) ! IN TRUE IF BOX CONTAINS SEA-ICE. SJT1F304.879
+,CALIB ! IN TRUE IF CALIBRATION EXPT SJT1F304.880
+,OPENSEA(L1) ! IN TRUE IF BOX CONTAINS OPEN SEA (NO ICE) SJT1F304.881
+,LAND(L1) ! IN TRUE AT LAND POINTS SJT1F304.882
+,NEWICE(L1) ! OUT TRUE IF ICE IS FORMING SJT1F304.883
+,L_THERM ! IN TRUE FOR COUPLED LIKE ICE THERMODYNAMICS SJT1F304.884
+,L_IDYN ! IN TRUE FOR CAV FLUID ICE DYNAMICS SJT1F304.885
+,L_IDRIF ! IN TRUE FOR ICE DEPTH ADVECTION SJT1F304.886
+,L_SLBADV ! IN TRUE FOR SLAB TEMP ADVECTION SJT1F304.887
+,LGLOBAL ! IN TRUE IF GLOBAL MODEL SJT1F304.888
C SJT1F304.889
REAL SLABTEMP(L1) ! INOUT TEMPERATURE OF THE SLAB OCEAN. SJT1F304.890
+,wtsfc(l1) ! OUT W x sfc slab temp SJC1F400.106
+,wtbase(l1) ! OUT w x base slab temp SJC1F400.107
+,AICE(L1) ! INOUT ICE FRACTION SJT1F304.891
+,OIFLUX(L1) ! OUT OCEAN TO ICE HEAT FLUX SJT1F304.892
+,CARYHEAT(L1) ! OUT ZERO EXCEPT AT POINTS WHERE ICE IS SJT1F304.893
+ ! ABOUT TO FORM, WHERE IT SHOULD BE THE UMSLAB1A.55
+ ! NEGATIVE HEAT FLUX UMSLAB1A.56
+,DTADV(L1) ! OUT HEATING RATE DUE TO ADVECTION K/S SJC1F400.108
+,DTDIFF(L1) ! OUT HEATING RATE DUE TO DIFFUSION K/S SJC1F400.109
C* UMSLAB1A.57
C Include COMDECKS UMSLAB1A.58
C UMSLAB1A.59
*CALL C_SLAB UMSLAB1A.60
*CALL C_LHEAT UMSLAB1A.61
*CALL C_0_DG_C UMSLAB1A.62
C UMSLAB1A.63
C UMSLAB1A.64
C VARIABLES LOCAL TO THIS ROUTINE ARE NOW DEFINED. UMSLAB1A.65
C UMSLAB1A.66
INTEGER J ! LOOP COUNTERS UMSLAB1A.67
C UMSLAB1A.68
REAL TFREEZE ! FREEZING POINT OF SEAWATER IN C SJT1F304.894
+,SLABD(L1) ! TEMPORARY STORE FOR SLABTEMP-TCLIMC SJT1F304.895
+,TWORK(L1) ! TEMPORARY STORE FOR SLABTEMP SJC1F400.110
+,DTTOFLUX ! CONVERSION FACTOR FROM K TO WM-2 SJT1F304.896
+,FLUXTODT ! CONVERSION FACTOR FROM WM-2 TO K SJT1F304.897
+,DEPTHTOK ! CONVERSION FACTOR FROM M ICE TO K SJT1F304.898
+,AHDT ! DIFFUSION COEFFICENT * TIMESTEP SJT1F304.900
C UMSLAB1A.75
PARAMETER(TFREEZE=TFS-ZERODEGC) UMSLAB1A.76
C UMSLAB1A.77
C SJT1F304.901
C SET VARIOUS CONSTANTS FOR USE IN CALCULATING THE HEAT SJT1F304.902
C CONVERGENCE SJT1F304.903
C SJT1F304.904
DTTOFLUX = RHOCP * DZ1 / DT SJT1F304.905
FLUXTODT = 1.0 / DTTOFLUX SJT1F304.906
DEPTHTOK = RHOICE * LF / ( RHOCP * DZ1 ) SJT1F304.907
C SJT1F304.908
DO J = 1,L2 UMSLAB1A.78
C UMSLAB1A.79
C 1. SET TEMPERATURE OF ICE POINTS TO TFREEZE (if coupled model SJT1F304.909
C ice thermodynamics is not required) SJT1F304.910
C AND INITIALISE CARYHEAT AND ADJHCONV (IN CALIBRATION) SJT1F304.911
C and initialise newice for coupled model ice thermo. SJT1F304.912
C UMSLAB1A.81
IF ( ICY(J) .and. .not. l_therm) THEN SJT1F304.913
SLABTEMP(J) = TFREEZE UMSLAB1A.94
END IF UMSLAB1A.95
C UMSLAB1A.96
IF ( CALIB ) THEN SJT1F304.914
ADJHCONV(J) = 0.0 SJT1F304.915
END IF SJT1F304.916
C SJT1F304.917
CARYHEAT(J) = 0.0 SJT1F304.918
newice(j) = .false. SJT1F304.919
C SJT1F304.920
C 2. CALCULATE THE CHANGE IN SLAB TEMPERATURE SJT1F304.921
C and ocean to ice heat flux SJT1F304.922
C SJT1F304.923
C PERFORM OPERATIONS OVER SEA POINTS ONLY UMSLAB1A.97
C (for coupled ice thermodynamics, sea points means sea and sea ice) SJT1F304.924
C UMSLAB1A.98
C First coupled thermodynamics SJT1F304.925
IF (L_THERM) THEN SJT1F304.926
IF ( .not. LAND(J) ) THEN SJT1F304.927
oiflux(j) = rhocp * eddydiff / (0.5*dz1) SJT1F304.928
& * ( slabtemp(j) -tfreeze ) SJT1F304.929
if (opensea(j)) oiflux(j) = 0.0 SJT1F304.930
SLABTEMP(J) = SLABTEMP(J) + ( ATMSFLUX(J) SJT1F304.931
& - LF*SNOWRATE(J)*(1.0-aice(j)) SJT1F304.932
& - oiflux(j) ) * fluxtodt SJT1F304.933
C UMSLAB1A.100
ENDIF SJT1F304.934
C Then original thermodynamic treatment SJT1F304.935
ELSE SJT1F304.936
IF ( OPENSEA(J) ) THEN SJT1F304.937
SLABTEMP(J) = SLABTEMP(J) + ( ATMSFLUX(J) SJT1F304.938
+ - LF*SNOWRATE(J) ) * FLUXTODT SJT1F304.939
ENDIF SJT1F304.940
ENDIF SJT1F304.941
C End J loop. SJT1F304.942
END DO SJT1F304.943
C UMSLAB1A.102
C 3. CALCULATE ADVECTION OF SLAB OCEAN TEMPERATURE SJT1F304.944
C (and diagnose heating rate) SJC1F400.111
IF (L_SLBADV) THEN SJT1F304.945
do j=1,l1 SJC1F400.112
twork(j) = slabtemp(j) SJC1F400.113
end do SJC1F400.114
CALL SLAB_T_UV( SJT1F304.946
*CALL ARGOINDX SDR1F404.104
& L1,L2,ICOLS,JROWS,JROWSM1,LAND SJT1F304.947
& ,LGLOBAL,u_field,dz1 SJC1F400.115
& ,DELTA_LAT,DELTA_LONG,BASE_LAT,DT SJT1F304.949
& ,COS_U_LATITUDE,COS_P_LATITUDE SJT1F304.950
& ,SEC_P_LATITUDE,SIN_U_LATITUDE SJT1F304.951
& ,UCURRENT,VCURRENT SJT1F304.952
& ,SLABTEMP SJT1F304.953
& ,OPENSEA SJT1F304.954
& ,wtsfc SJC1F400.116
& ,wtbase SJC1F400.117
& ) SJT1F304.955
do j=1,l1 SJC1F400.118
dtadv(j) = ( slabtemp(j)-twork(j) ) / dt SJC1F400.119
end do SJC1F400.120
ENDIF SJT1F304.956
C UMSLAB1A.106
C UMSLAB1A.110
C USE TEMPORARY SPACE TO HOLD SLABTEMP-TCLIMC FOR USE SJT1F304.957
C IN DIFFUSION SUBROUTINE SJT1F304.958
C SJT1F304.959
DO J=1,L2 SJT1F304.960
IF ( .NOT. CALIB ) THEN SJT1F304.961
SLABD(J) = SLABTEMP(J) - TCLIMC(J) SJT1F304.962
ENDIF SJT1F304.963
END DO SJT1F304.964
C SJT1F304.965
C SJT1F304.966
C 4. CALCULATE DEL2 DIFFUSION OF SLAB TEMPERATURES SJT1F304.967
C (DIFFUSE SLABTEMP-TCLIMC - SO DONT CALL IN SJT1F304.968
C CALIBRATION EXPERIMENT) SJT1F304.969
C SJT1F304.970
IF ( .NOT. CALIB) THEN SJT1F304.971
IF ( AH .GT. 0.0) THEN SJT1F304.972
AHDT= AH * DT SJT1F304.973
CALL SLABDIFF(SLABD, SJT1F304.974
+ OPENSEA, SJT1F304.975
+ L1,L2, SJT1F304.976
+ JROWS,ICOLS, SJT1F304.977
+ AHDT, SJT1F304.978
+ DELTA_LONG,DELTA_LAT,BASE_LAT, SJT1F304.979
+ COS_P_LATITUDE,COS_U_LATITUDE,SEC_P_LATITUDE) SJT1F304.980
ENDIF SJT1F304.981
ENDIF SJT1F304.982
C SJT1F304.983
C RECALCULATE SLABTEMP IN NON-CALIBRATION EXPERIMENT SJT1F304.984
C PUT DIFFUSION INCREMENT INTO NEW STORE TO ADD TO HEAT CONV SJT1F304.985
C SJT1F304.986
C 5. IF MIXED LAYER TEMPERATURE IS BELOW FREEZING THEN ICE FORMATION SJT1F304.987
C IS JUST STARTING. CALCULATE THE HEAT DEFICIT THAT ARISES UMSLAB1A.112
C FROM THE TEMPERATURE CHANGE AND PASS TO SLABICE THROUGH UMSLAB1A.113
C CARYHEAT(J) FOR ICE FORMATION TO OCCUR. SET SLABTEMP TO UMSLAB1A.114
C TFREEZE UMSLAB1A.115
C First coupled thermodynamics SJT1F304.988
C SJT1F304.989
IF (l_therm) then SJT1F304.990
DO J = 1,L2 SJT1F304.991
IF ( .NOT. CALIB) THEN SJT1F304.992
DTDIFF(J) = SLABD(J) - SLABTEMP(J) + TCLIMC(J) SJC1F400.121
ADJHCONV(J) = ADJHCONV(J) + DTDIFF(J) * DTTOFLUX SJC1F400.122
dtdiff(j) = dtdiff(j) / dt SJC1F400.123
ENDIF SJT1F304.995
C SJT1F304.996
IF ( .NOT.LAND(J) ) THEN SJT1F304.997
C SJT1F304.998
IF ( SLABTEMP(J) .LT. TFREEZE ) THEN SJC1F400.124
CARYHEAT(J) = SLABTEMP(J) - TFREEZE SJC1F400.125
CARYHEAT(J) = (CARYHEAT(J) * RHOCP * DZ1) / DT SJC1F400.126
SLABTEMP(J) = TFREEZE SJC1F400.127
NEWICE(J) = (.NOT. ICY(J)) SJC1F400.128
ENDIF SJC1F400.129
IF ( CALIB ) THEN SJT1F304.999
IF ( l_idyn .or. l_idrif ) THEN SJT1F304.1000
ADJHCONV(J) = ( TCLIMC(J) - SJT1F304.1001
+ SLABTEMP(J) ) * DTTOFLUX SJT1F304.1002
ELSE SJT1F304.1003
ADJHCONV(J) = ( ( TCLIMC(J) - DEPTHTOK * HCLIM(J) ) - SJT1F304.1004
+ SLABTEMP(J) ) * DTTOFLUX SJT1F304.1005
END IF SJT1F304.1006
END IF SJT1F304.1007
C SJT1F304.1008
SLABTEMP(J) = SLABTEMP (J) + ADJHCONV(J) * FLUXTODT SJT1F304.1009
C UMSLAB1A.116
ENDIF UMSLAB1A.122
C End J loop. SJT1F304.1011
END DO SJT1F304.1012
C Then original thermodynamic treatment SJT1F304.1013
ELSE SJT1F304.1014
DO J = 1,L2 SJT1F304.1015
IF ( .NOT. CALIB) THEN SJT1F304.1016
DTDIFF(J) = SLABD(J) - SLABTEMP(J) + TCLIMC(J) SJC1F400.130
ADJHCONV(J) = ADJHCONV(J) + DTDIFF(J) * DTTOFLUX SJC1F400.131
dtdiff(j) = dtdiff(j) / dt SJC1F400.132
ENDIF SJT1F304.1019
C UMSLAB1A.123
IF ( OPENSEA(J) ) THEN SJT1F304.1020
C SJT1F304.1021
IF ( CALIB ) THEN SJT1F304.1022
ADJHCONV(J) = ( ( TCLIMC(J) - DEPTHTOK * HCLIM(J) ) - SJT1F304.1023
+ SLABTEMP(J) ) * DTTOFLUX SJT1F304.1024
END IF SJT1F304.1025
C SJT1F304.1026
SLABTEMP(J) = SLABTEMP (J) + ADJHCONV(J) * FLUXTODT SJT1F304.1027
C SJT1F304.1028
IF ( SLABTEMP(J) .LT. TFREEZE ) THEN SJT1F304.1029
CARYHEAT(J) = SLABTEMP(J) - TFREEZE SJT1F304.1030
CARYHEAT(J) = (CARYHEAT(J) * RHOCP * DZ1) / DT SJT1F304.1031
SLABTEMP(J) = TFREEZE SJT1F304.1032
NEWICE(J) = (.NOT. ICY(J)) SJT1F304.1033
ENDIF SJT1F304.1034
ENDIF SJT1F304.1035
C End J loop. SJT1F304.1036
END DO UMSLAB1A.124
C End IF block for coupled thermodynamics. SJT1F304.1037
ENDIF SJT1F304.1038
RETURN UMSLAB1A.125
END UMSLAB1A.126
*ENDIF UMSLAB1A.127