/[MITgcm]/MITgcm/pkg/seaice/seaice_ocean_stress.F
ViewVC logotype

Annotation of /MITgcm/pkg/seaice/seaice_ocean_stress.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.24 - (hide annotations) (download)
Fri May 29 14:51:21 2009 UTC (15 years ago) by mlosch
Branch: MAIN
CVS Tags: checkpoint61p, checkpoint61q
Changes since 1.23: +15 -103 lines
  - clean up computation of Hibler+Bryan (1987) stress coupling for the case
    of LSR (this change is expected to change the results slightly because now
    it uses slightly difference moduli from the second last LSR solution, but
    that is more consistent with the stress computations; this part of the code
    is not tested in the verification experiments)

1 mlosch 1.24 C $Header: /u/gcmpack/MITgcm/pkg/seaice/seaice_ocean_stress.F,v 1.23 2009/05/29 10:18:03 mlosch Exp $
2 mlosch 1.1 C $Name: $
3    
4     #include "SEAICE_OPTIONS.h"
5    
6     CStartOfInterface
7 jmc 1.19 SUBROUTINE SEAICE_OCEAN_STRESS(
8 mlosch 1.1 I myTime, myIter, myThid )
9     C /==========================================================\
10     C | SUBROUTINE SEAICE_OCEAN_STRESS |
11     C | o Calculate ocean surface stresses |
12     C | - C-grid version |
13     C |==========================================================|
14     C \==========================================================/
15     IMPLICIT NONE
16    
17     C === Global variables ===
18     #include "SIZE.h"
19     #include "EEPARAMS.h"
20     #include "PARAMS.h"
21 dimitri 1.21 #include "DYNVARS.h"
22 mlosch 1.5 #include "GRID.h"
23 mlosch 1.1 #include "FFIELDS.h"
24     #include "SEAICE.h"
25     #include "SEAICE_PARAMS.h"
26    
27     C === Routine arguments ===
28     C myTime - Simulation time
29     C myIter - Simulation timestep number
30     C myThid - Thread no. that called this routine.
31     _RL myTime
32     INTEGER myIter
33     INTEGER myThid
34     CEndOfInterface
35    
36 jmc 1.19 #ifdef SEAICE_CGRID
37 mlosch 1.1 C === Local variables ===
38     C i,j,bi,bj - Loop counters
39    
40     INTEGER i, j, bi, bj
41 mlosch 1.5 _RL SINWAT, COSWAT, SINWIN, COSWIN
42 mlosch 1.24 _RL fuIceLoc, fvIceLoc
43 mlosch 1.4 _RL areaW, areaS
44 mlosch 1.1
45     c introduce turning angle (default is zero)
46     SINWAT=SIN(SEAICE_waterTurnAngle*deg2rad)
47     COSWAT=COS(SEAICE_waterTurnAngle*deg2rad)
48 mlosch 1.5 SINWIN=SIN(SEAICE_airTurnAngle*deg2rad)
49     COSWIN=COS(SEAICE_airTurnAngle*deg2rad)
50 mlosch 1.1
51 mlosch 1.5 IF ( useHB87StressCoupling ) THEN
52     C
53 jmc 1.19 C use an intergral over ice and ocean surface layer to define
54 mlosch 1.5 C surface stresses on ocean following Hibler and Bryan (1987, JPO)
55 jmc 1.19 C
56 mlosch 1.5 DO bj=myByLo(myThid),myByHi(myThid)
57     DO bi=myBxLo(myThid),myBxHi(myThid)
58 mlosch 1.24 DO J=1,sNy
59     DO I=1,sNx
60 jmc 1.19 C average wind stress over ice and ocean and apply averaged wind
61 mlosch 1.14 C stress and internal ice stresses to surface layer of ocean
62 mlosch 1.24 areaW = 0.5 * (AREA(I,J,1,bi,bj) + AREA(I-1,J,1,bi,bj))
63     & * SEAICEstressFactor
64     areaS = 0.5 * (AREA(I,J,1,bi,bj) + AREA(I,J-1,1,bi,bj))
65     & * SEAICEstressFactor
66     fu(I,J,bi,bj)=(ONE-areaW)*fu(I,J,bi,bj)
67     & + areaW*taux(I,J,bi,bj)
68     & + stressDivergenceX(I,J,bi,bj) * SEAICEstressFactor
69     fv(I,J,bi,bj)=(ONE-areaS)*fv(I,J,bi,bj)
70     & + areaS*tauy(I,J,bi,bj)
71     & + stressDivergenceY(I,J,bi,bj) * SEAICEstressFactor
72 mlosch 1.14 ENDDO
73 mlosch 1.24 ENDDO
74 mlosch 1.5 ENDDO
75     ENDDO
76 jmc 1.17
77 mlosch 1.5 ELSE
78 jmc 1.17 C else: useHB87StressCoupling=F
79 mlosch 1.5
80 jmc 1.19 C-- Compute ice-affected wind stress (interpolate to U/V-points)
81     C by averaging wind stress and ice-ocean stress according to
82 mlosch 1.5 C ice cover
83 mlosch 1.1 DO bj=myByLo(myThid),myByHi(myThid)
84     DO bi=myBxLo(myThid),myBxHi(myThid)
85     DO j=1,sNy
86     DO i=1,sNx
87 mlosch 1.18 fuIceLoc=HALF*( DWATN(I,J,bi,bj)+DWATN(I-1,J,bi,bj) )*
88 jmc 1.19 & COSWAT *
89 dimitri 1.21 & ( UICE(I,J,1,bi,bj)-uVel(I,J,1,bi,bj) )
90 jmc 1.19 & - SIGN(SINWAT, _fCori(I,J,bi,bj)) * 0.5 _d 0 *
91 mlosch 1.6 & ( DWATN(I ,J,bi,bj) *
92 dimitri 1.21 & 0.5 _d 0*(vIce(I ,J ,1,bi,bj)-vVel(I ,J ,1,bi,bj)
93     & +vIce(I ,J+1,1,bi,bj)-vVel(I ,J+1,1,bi,bj))
94 mlosch 1.6 & + DWATN(I-1,J,bi,bj) *
95 dimitri 1.21 & 0.5 _d 0*(vIce(I-1,J ,1,bi,bj)-vVel(I-1,J ,1,bi,bj)
96     & +vIce(I-1,J+1,1,bi,bj)-vVel(I-1,J+1,1,bi,bj))
97 mlosch 1.1 & )
98 mlosch 1.18 fvIceLoc=HALF*( DWATN(I,J,bi,bj)+DWATN(I,J-1,bi,bj) )*
99 mlosch 1.6 & COSWAT *
100 dimitri 1.21 & ( VICE(I,J,1,bi,bj)-vVel(I,J,1,bi,bj) )
101 mlosch 1.6 & + SIGN(SINWAT, _fCori(I,J,bi,bj)) * 0.5 _d 0 *
102     & ( DWATN(I,J ,bi,bj) *
103 dimitri 1.21 & 0.5 _d 0*(uIce(I ,J ,1,bi,bj)-uVel(I ,J ,1,bi,bj)
104     & +uIce(I+1,J ,1,bi,bj)-uVel(I+1,J ,1,bi,bj))
105 mlosch 1.6 & + DWATN(I,J-1,bi,bj) *
106 dimitri 1.21 & 0.5 _d 0*(uIce(I ,J-1,1,bi,bj)-uVel(I ,J-1,1,bi,bj)
107     & +uIce(I+1,J-1,1,bi,bj)-uVel(I+1,J-1,1,bi,bj))
108 mlosch 1.1 & )
109 mlosch 1.4 areaW = 0.5 _d 0 * (AREA(I,J,1,bi,bj) + AREA(I-1,J,1,bi,bj))
110 mlosch 1.9 & * SEAICEstressFactor
111 mlosch 1.4 areaS = 0.5 _d 0 * (AREA(I,J,1,bi,bj) + AREA(I,J-1,1,bi,bj))
112 mlosch 1.9 & * SEAICEstressFactor
113 mlosch 1.11 fu(I,J,bi,bj)=(ONE-areaW)*fu(I,J,bi,bj)+areaW*fuIceLoc
114     fv(I,J,bi,bj)=(ONE-areaS)*fv(I,J,bi,bj)+areaS*fvIceLoc
115 mlosch 1.1 ENDDO
116     ENDDO
117     ENDDO
118     ENDDO
119 mlosch 1.5 ENDIF
120 mlosch 1.1 CALL EXCH_UV_XY_RS(fu, fv, .TRUE., myThid)
121 mlosch 1.3
122 jmc 1.17 #endif /* SEAICE_CGRID */
123 mlosch 1.1
124     RETURN
125     END

  ViewVC Help
Powered by ViewVC 1.1.22