*IF DEF,SEAICE CNVSTOP.2
C ******************************COPYRIGHT****************************** CNVSTOP.3
C (c) CROWN COPYRIGHT 1998, METEOROLOGICAL OFFICE, All Rights Reserved. CNVSTOP.4
C CNVSTOP.5
C Use, duplication or disclosure of this code is subject to the CNVSTOP.6
C restrictions as set forth in the contract. CNVSTOP.7
C CNVSTOP.8
C Meteorological Office CNVSTOP.9
C London Road CNVSTOP.10
C BRACKNELL CNVSTOP.11
C Berkshire UK CNVSTOP.12
C RG12 2SZ CNVSTOP.13
C CNVSTOP.14
C If no contract has been raised with this copy of the code, the use, CNVSTOP.15
C duplication or disclosure of it is strictly prohibited. Permission CNVSTOP.16
C to do so must first be obtained in writing from the Head of Numerical CNVSTOP.17
C Modelling at the above address. CNVSTOP.18
C ******************************COPYRIGHT****************************** CNVSTOP.19
!+Convergence stopping routine for sea-ice CNVSTOP.20
! CNVSTOP.21
! Subroutine Interface: CNVSTOP.22
subroutine cnvstop( 1,7CNVSTOP.23
*CALL ARGOINDX
CNVSTOP.24
& imt,imtm1,imtm2,jmt,jmtm1,jmtm2, CNVSTOP.25
& uice, CNVSTOP.26
& vice, CNVSTOP.27
& hice, CNVSTOP.28
& radius, CNVSTOP.29
& dphi, CNVSTOP.30
& dlambda, CNVSTOP.31
& phit, CNVSTOP.32
& hicemax,hicemid CNVSTOP.33
& ) CNVSTOP.34
IMPLICIT NONE CNVSTOP.35
! CNVSTOP.36
! Description: CNVSTOP.37
! This subroutine restricts the flow of sea-ice up a gradient in HICE CNVSTOP.38
! (GBM sea-ice depth), as pseudo-rheological resistance to convergence. CNVSTOP.39
! The flow parallel to grad(HICE) is impeded. Where HICE < HICEMID, CNVSTOP.40
! the flow is uneffected, and where HICE > HICEMAX the flow is stopped CNVSTOP.41
! completely, with ! a linear ramping of the change to the velocity CNVSTOP.42
! between HICEMID and HICEMAX. It is set up to work with the freedrift CNVSTOP.43
! scheme, being called in ICEFREEDR, but could equally well work with CNVSTOP.44
! the simple ice advection scheme, replacing ORH3F403.76 to CNVSTOP.45
! ICEDRIFT.271 (at vn4.4) of ICEDRIFT. CNVSTOP.46
! CNVSTOP.47
! Method: CNVSTOP.48
! Looping over velocity points, the maximum of surrounding values of CNVSTOP.49
! HICE is taken and labelled HICE_U. This value dictates whether, and CNVSTOP.50
! how, the ice velocity will be restricted. If HICE_U > HICEMID the CNVSTOP.51
! components of grad(HICE) and U.grad(HICE) are calculated. If CNVSTOP.52
! U.grad(HICE) is positive (i.e.the ice is deep and converging) then CNVSTOP.53
! the component of U parallel to grad(HICE) is scaled s.t. it is CNVSTOP.54
! unchanged where HICE_U = HICEMID, and zeroed where HICE_U = HICEMAX, CNVSTOP.55
! the scaling varying linearly between the two states. CNVSTOP.56
! CNVSTOP.57
! Current Code Owner: Doug Cresswell CNVSTOP.58
! CNVSTOP.59
! History: CNVSTOP.60
! Version Date Comment CNVSTOP.61
! ------- ---- ------- CNVSTOP.62
! 4.5 25.9.98 New deck. Doug Cresswell and Jonathan Gregory. CNVSTOP.63
! CNVSTOP.64
! Code Description: CNVSTOP.65
! Language: FORTRAN 77 + common extensions. CNVSTOP.66
! This code is written to UMDP3 v6 programming standards. CNVSTOP.67
C CNVSTOP.68
*CALL OTIMER
CNVSTOP.69
*CALL TYPOINDX
CNVSTOP.70
*CALL CNTLOCN
CNVSTOP.71
*CALL OARRYSIZ
CNVSTOP.72
integer CNVSTOP.73
C*CALL ARGSIZE CNVSTOP.74
C integer CNVSTOP.75
& imt ! number of tracer columns. CNVSTOP.76
&,imtm1 ! number of tracer columns - 1 CNVSTOP.77
&,imtm2 ! number of tracer columns - 2 CNVSTOP.78
&,jmt ! number of tracer rows. CNVSTOP.79
&,jmtm1 ! number of velocity rows. CNVSTOP.80
&,jmtm2 ! number of velocity rows - 1 CNVSTOP.81
real CNVSTOP.82
& uice(imt,jmtm1) ! inout zonal sea ice velocity. CNVSTOP.83
&,vice(imt,jmtm1) ! inout meridional sea ice velocity. CNVSTOP.84
&,hice(imt,jmt) ! in ice thickness. CNVSTOP.85
&,radius ! in radius of the earth in metres. CNVSTOP.86
&,dphi(jmt) ! in meridional grid spacing in radians. CNVSTOP.87
&,dlambda(imt) ! in zonal grid spacing in radians. CNVSTOP.88
&,phit(jmt) ! in latitude of mass rows in radians. CNVSTOP.89
&,hmask(imt,jmt) ! in 1.0 for ha land 0.0 for sea. CNVSTOP.90
&,umask(imt,jmtm1) ! in 1.0 for uv land 0.0 for sea. CNVSTOP.91
&,hicemax ! in max hice at which convergence is allowed. CNVSTOP.92
&,hicemid ! in min hice at which convergence is impeded. CNVSTOP.93
C Local Variables CNVSTOP.94
integer CNVSTOP.95
& i,j !Looping variables CNVSTOP.96
real CNVSTOP.97
& facex,facey,faceyp !Lengths of the E(&W),S and N sides CNVSTOP.98
!of the gridbox. CNVSTOP.99
&,gradhu(imt,jmtm1) ! U cpt of the grad h (ice depth) CNVSTOP.100
&,gradhv(imt,jmtm1) ! V cpt of the grad h (ice depth) CNVSTOP.101
&,modgradhsq(imt,jmtm1) ! Mod(Grad h), squared. CNVSTOP.102
&,udotgradh(imt,jmtm1) ! Scalar product of u and grad h CNVSTOP.103
&,hice_u(imt,jmtm1) ! h interpolated onto the U grid CNVSTOP.104
&,recip_hmax_m_hmid !1/(hicemax-hicemin) (max=100m-1) CNVSTOP.105
&,b !Scaling factor for ice between hicemid CNVSTOP.106
!and hicemax. CNVSTOP.107
C Start executable code CNVSTOP.108
if (L_OTIMER) call timer
('cnvstop',3) CNVSTOP.109
*IF DEF,MPP CNVSTOP.110
CALL SWAPBOUNDS
(UICE,IMT,JMTM1,O_EW_HALO,O_NS_HALO,1) CNVSTOP.111
CALL SWAPBOUNDS
(VICE,IMT,JMTM1,O_EW_HALO,O_NS_HALO,1) CNVSTOP.112
CALL SWAPBOUNDS
(HICE,IMT,JMT,O_EW_HALO,O_NS_HALO,1) CNVSTOP.113
*ENDIF CNVSTOP.114
recip_hmax_m_hmid=1.0/max((hicemax-hicemid),0.01) CNVSTOP.115
C Start looping over u points CNVSTOP.116
do j=J_1,J_jmt CNVSTOP.117
do i=2,imtm1 CNVSTOP.118
C Initialise arrays to zero, to be partially overwritten CNVSTOP.119
C (Note that these arrays are not made cyclic, as they are not output) CNVSTOP.120
gradhu(i,j)=0.0 CNVSTOP.121
gradhv(i,j)=0.0 CNVSTOP.122
modgradhsq(i,j)=0.0 CNVSTOP.123
udotgradh(i,j)=0.0 CNVSTOP.124
C Take hice on the u grid to be the maximum of surrounding hice's. This CNVSTOP.125
C prevents advection into deep areas up steep gradients in hice. CNVSTOP.126
hice_u(i,j)=max(hice(i,j),hice(i+1,j) CNVSTOP.127
& ,hice(i,j+1),hice(i+1,j+1)) CNVSTOP.128
C If ice is deep then procede. CNVSTOP.129
if(hice_u(i,j).gt.hicemid) then CNVSTOP.130
C Calculate gridbox side lengths CNVSTOP.131
facex = radius*dphi(j) CNVSTOP.132
facey = (radius*cos(phit(j)))*dlambda(i+1) CNVSTOP.133
faceyp = (radius*cos(phit(j+1)))*dlambda(i+1) CNVSTOP.134
C Calculate Grad h (u and v cpts) CNVSTOP.135
gradhu(i,j)=0.5*(hice(i+1,j+1)-hice(i,j+1))/faceyp CNVSTOP.136
& +0.5*(hice(i+1,j)-hice(i,j))/facey CNVSTOP.137
gradhv(i,j)=0.5*((hice(i,j+1)-hice(i,j)) CNVSTOP.138
& +(hice(i+1,j+1)-hice(i+1,j)))/facex CNVSTOP.139
C Calculate scalar product of U and Grad h CNVSTOP.140
udotgradh(i,j)=uice(i,j)*gradhu(i,j) CNVSTOP.141
& +vice(i,j)*gradhv(i,j) CNVSTOP.142
C If U is acting up the ice depth (h) gradient then CNVSTOP.143
C set the cpt of U parallel to Grad h to zero. CNVSTOP.144
if(udotgradh(i,j).gt.0.0) then CNVSTOP.145
modgradhsq(i,j)=gradhu(i,j)*gradhu(i,j) CNVSTOP.146
& +gradhv(i,j)*gradhv(i,j) CNVSTOP.147
if(hice_u(i,j).gt.hicemax)then CNVSTOP.148
uice(i,j)=uice(i,j)-(udotgradh(i,j)*gradhu(i,j)) CNVSTOP.149
& /modgradhsq(i,j) CNVSTOP.150
vice(i,j)=vice(i,j)-(udotgradh(i,j)*gradhv(i,j)) CNVSTOP.151
& /modgradhsq(i,j) CNVSTOP.152
else if(hice_u(i,j).gt.hicemid)then CNVSTOP.153
b=(hice_u(i,j)-hicemid)*recip_hmax_m_hmid CNVSTOP.154
uice(i,j)=uice(i,j)-b*(udotgradh(i,j)*gradhu(i,j)) CNVSTOP.155
& /modgradhsq(i,j) CNVSTOP.156
vice(i,j)=vice(i,j)-b*(udotgradh(i,j)*gradhv(i,j)) CNVSTOP.157
& /modgradhsq(i,j) CNVSTOP.158
endif CNVSTOP.159
endif CNVSTOP.160
endif CNVSTOP.161
C Looping - need to sort out cyclic points CNVSTOP.162
enddo CNVSTOP.163
enddo CNVSTOP.164
if(L_OCYCLIC) then CNVSTOP.165
do j=J_1,J_jmt CNVSTOP.166
uice(1,j)=uice(imtm1,j) CNVSTOP.167
uice(imt,j)=uice(2,j) CNVSTOP.168
vice(1,j)=vice(imtm1,j) CNVSTOP.169
vice(imt,j)=vice(2,j) CNVSTOP.170
enddo CNVSTOP.171
endif CNVSTOP.172
CALL SWAPBOUNDS
(UICE,IMT,JMTM1,O_EW_HALO,O_NS_HALO,1) CNVSTOP.173
CALL SWAPBOUNDS
(VICE,IMT,JMTM1,O_EW_HALO,O_NS_HALO,1) CNVSTOP.174
if (L_OTIMER) call timer
('cnvstop',4) CNVSTOP.175
end CNVSTOP.176
*ENDIF CNVSTOP.177