forked from altMITgcm/MITgcm
-
Notifications
You must be signed in to change notification settings - Fork 0
/
gad_osc_mul_x.F
81 lines (57 loc) · 2.21 KB
/
gad_osc_mul_x.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
C $Header: /u/gcmpack/MITgcm/pkg/generic_advdiff/gad_osc_mul_x.F,v 1.1 2016/03/13 01:44:02 jmc Exp $
C $Name: $
# include "GAD_OPTIONS.h"
SUBROUTINE GAD_OSC_MUL_X(ix,hh,mask,ohat,scal)
C |================================================================|
C | OSC_MUL_X: evaluate WENO oscillation weights in X. |
C |================================================================|
implicit none
C =============================================== global variables
# include "SIZE.h"
integer ix,hh
_RL mask(1-OLx:sNx+OLx)
_RL ohat(1:2,
& 1-OLx:sNx+OLx)
_RL scal(1:2)
integer ii
_RL dels,dfs1,dfs2
_RL osum,zero,mval
_RL oval,omin,omax
C =============================== calc. WENO oscillation weighting
c omin = +huge(+1. _d 0)
c omax = -huge(+1. _d 0)
omin = +1. _d 99
omax = -1. _d 99
zero = 1. _d -20
mval = 1. _d + 0
do ii = ix-hh, ix+hh
C =============================== calc. derivatives centred on II.
dels = (ii - ix) * 2. _d 0
dfs1 = ohat(1,ii)
dfs2 = ohat(2,ii)
dfs1 = dfs1 + dfs2 * dels
C =============================== oscl. = NORM(H^N * D^N/DX^N(F)).
oval = (2. _d 0 * dfs1)**2
& + (4. _d 0 * dfs2)**2
if (oval.lt.omin) omin = oval
if (oval.gt.omax) omax = oval
C =============================== any mask across oscil. stencil
mval = mval * mask(ii)
end do
if (mval .gt. 0. _d 0) then
C =============================== calc. WENO-style profile weights
scal(1) = 1. _d 5
& / (omax + zero)**3
scal(2) = 1. _d 0
& / (omin + zero)**3
osum = scal(1) + scal(2)
scal(1) = scal(1) / osum
scal(2) = scal(2) / osum
else
C =============================== default to MONO. profile weights
scal(1) = 0. _d 0
scal(2) = 1. _d 0
end if
return
c end subroutine GAD_OSC_MUL_X
end