-
Notifications
You must be signed in to change notification settings - Fork 237
/
calc_eddy_stress.F
106 lines (94 loc) · 2.98 KB
/
calc_eddy_stress.F
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
#include "PACKAGES_CONFIG.h"
#include "CPP_OPTIONS.h"
#ifdef ALLOW_GMREDI
# include "GMREDI_OPTIONS.h"
#endif
SUBROUTINE CALC_EDDY_STRESS( bi, bj, myThid )
C !DESCRIPTION: \bv
C *==========================================================*
C | S/R CALC_EDDY_STRESS
C | o Calculates the eddy stress when running a residual
C | ocean model
C *==========================================================*
C | Calculates the eddy stress. Later this will be added to
C | gU the same as external sources (e.g. wind stress, bottom
C | friction, etc.
C *==========================================================*
C \ev
IMPLICIT NONE
#include "SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "GRID.h"
#include "FFIELDS.h"
#ifdef ALLOW_GMREDI
# include "GMREDI.h"
#endif
C !INPUT/OUTPUT PARAMETERS:
C == Routine arguments ==
INTEGER bi,bj
INTEGER myThid
#ifdef ALLOW_EDDYPSI
C !LOCAL VARIABLES:
C == Local variables ==
C Loop counters
INTEGER i,j,k
C Interpolated stream function and coriolis
_RL psix, psiy, coriU, coriV
C Calculate the eddy stress from the eddy induced streamfunction
#ifdef ALLOW_GMREDI
IF ( GM_InMomAsStress ) THEN
#endif
DO k=1,Nr
DO j=1-OLy,sNy+OLy-1
DO i=1-OLx+1,sNx+OLx
#ifdef ALLOW_GMREDI
psiy = op25*(GM_PsiY(i, j ,k,bi,bj)
& +GM_PsiY(i, j+1,k,bi,bj)
& +GM_PsiY(i-1,j ,k,bi,bj)
& +GM_PsiY(i-1,j+1,k,bi,bj))
#else
psiy = op25*(eddyPsiY(i, j ,k,bi,bj)
& +eddyPsiY(i, j+1,k,bi,bj)
& +eddyPsiY(i-1,j ,k,bi,bj)
& +eddyPsiY(i-1,j+1,k,bi,bj))
#endif
coriU = op5*(fcori(i-1,j,bi,bj)
& +fCori(i ,j,bi,bj))
tauxEddy(i,j,k,bi,bj) = rhoConst*coriU*psiy
ENDDO
ENDDO
DO j=1-OLy+1,sNy+OLy
DO i=1-OLx,sNx+OLx-1
#ifdef ALLOW_GMREDI
psix = op25*(GM_PsiX(i, j ,k,bi,bj)
& +GM_PsiX(i+1,j ,k,bi,bj)
& +GM_PsiX(i ,j-1,k,bi,bj)
& +GM_PsiX(i+1,j-1,k,bi,bj))
#else
psix = op25*(eddyPsiX(i, j ,k,bi,bj)
& +eddyPsiX(i+1,j ,k,bi,bj)
& +eddyPsiX(i ,j-1,k,bi,bj)
& +eddyPsiX(i+1,j-1,k,bi,bj))
#endif
coriV = op5*(fcori(i,j-1,bi,bj)
& +fCori(i,j ,bi,bj))
tauyEddy(i,j,k,bi,bj) = -rhoConst*coriV*psix
ENDDO
ENDDO
C- end k loop
ENDDO
#ifdef ALLOW_DIAGNOSTICS
IF ( useDiagnostics ) THEN
CALL DIAGNOSTICS_FILL( tauxEddy, 'TAUXEDDY',
& 0, Nr, 1, bi, bj, myThid )
CALL DIAGNOSTICS_FILL( tauyEddy, 'TAUYEDDY',
& 0, Nr, 1, bi, bj, myThid )
ENDIF
#endif
#ifdef ALLOW_GMREDI
ENDIF
#endif
#endif /* ALLOW_EDDYPSI */
RETURN
END