/[MITgcm]/MITgcm/pkg/generic_advdiff/gad_dst2u1_adv_y.F
ViewVC logotype

Annotation of /MITgcm/pkg/generic_advdiff/gad_dst2u1_adv_y.F

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


Revision 1.3 - (hide annotations) (download)
Mon Jun 12 16:06:28 2006 UTC (17 years, 11 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint58j_post, checkpoint58i_post
Changes since 1.2: +17 -8 lines
change formulation to be less sensitive to truncation errors.

1 jmc 1.3 C $Header: /u/gcmpack/MITgcm/pkg/generic_advdiff/gad_dst2u1_adv_y.F,v 1.2 2006/06/07 01:55:14 heimbach Exp $
2 jmc 1.1 C $Name: $
3    
4     #include "GAD_OPTIONS.h"
5    
6     CBOP
7     C !ROUTINE: GAD_DST2U1_ADV_Y
8    
9     C !INTERFACE: ==========================================================
10 jmc 1.3 SUBROUTINE GAD_DST2U1_ADV_Y(
11 jmc 1.1 I bi,bj,k, advectionScheme, deltaTloc,
12     I vTrans, vVel,
13     I tracer,
14     O vT,
15     I myThid )
16    
17     C !DESCRIPTION:
18     C Calculates the area integrated meridional flux due to advection
19 jmc 1.3 C of a tracer using second-order Direct Space and Time (DST-2)
20 jmc 1.1 C interpolation (=Lax-Wendroff) or simple 1rst order upwind scheme.
21    
22     C !USES: ===============================================================
23     IMPLICIT NONE
24     #include "SIZE.h"
25     #include "GRID.h"
26     #include "GAD.h"
27    
28     C !INPUT PARAMETERS: ===================================================
29     C bi,bj :: tile indices
30     C k :: vertical level
31     C advectionScheme :: advection scheme to use: either 2nd Order DST
32     C or 1rst Order Upwind
33     C vTrans :: meridional volume transport
34     C vVel :: meridional flow
35     C tracer :: tracer field
36     C myThid :: thread number
37     INTEGER bi,bj, k, advectionScheme
38     _RL deltaTloc
39     _RL vTrans(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
40     _RL vVel (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
41     _RL tracer(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
42     INTEGER myThid
43    
44     C !OUTPUT PARAMETERS: ==================================================
45     C vT :: meridional advective flux
46     _RL vT (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
47    
48     C !LOCAL VARIABLES: ====================================================
49     C i,j :: loop indices
50     C yLimit :: centered (vs upwind) fraction
51     C vFld :: velocity [m/s], meridional component
52     C vCFL :: Courant-Friedrich-Levy number
53     INTEGER i,j
54 jmc 1.3 _RL vFld, vCFL, yLimit, vAbs
55 jmc 1.1 CEOP
56    
57     yLimit = 0. _d 0
58     IF ( advectionScheme.EQ.ENUM_DST2 ) yLimit = 1. _d 0
59    
60     DO i=1-Olx,sNx+Olx
61     vT(i,1-Oly)=0.
62     ENDDO
63     DO j=1-Oly+1,sNy+Oly
64     DO i=1-Olx,sNx+Olx
65    
66     c vFld = vVel(i,j,k,bi,bj)
67     vFld = vTrans(i,j)*recip_dxG(i,j,bi,bj)
68 heimbach 1.2 & *recip_drF(k)*_recip_hFacS(i,j,k,bi,bj)
69 jmc 1.1 vCFL = ABS(vFld*deltaTloc*recip_dyC(i,j,bi,bj))
70    
71 jmc 1.3 c vT(i,j) =
72     c & vTrans(i,j)*(tracer(i,j-1)+tracer(i,j))*0.5 _d 0
73     c & + ( 1. _d 0 - yLimit*(1. _d 0 - vCFL) )*ABS(vTrans(i,j))
74     c & *(tracer(i,j-1)-tracer(i,j))*0.5 _d 0
75     C-- above formulation produces large truncation error when:
76     C 1rst.O upWind and v > 0 & |tracer(i,j-1)| << |tracer(i,j)|
77     C or v < 0 & |tracer(i,j-1)| >> |tracer(i,j)|
78     C-- change to a more robust expression:
79     vAbs = ABS(vTrans(i,j))
80     & *( 1. _d 0 - yLimit*(1. _d 0 - vCFL) )
81     vT(i,j) = ( vTrans(i,j)+vAbs )* 0.5 _d 0 * tracer(i,j-1)
82     & + ( vTrans(i,j)-vAbs )* 0.5 _d 0 * tracer(i,j)
83    
84 jmc 1.1 ENDDO
85     ENDDO
86    
87     RETURN
88     END

  ViewVC Help
Powered by ViewVC 1.1.22