Newer
Older
1001
1002
1003
1004
1005
1006
1007
1008
1009
1010
1011
1012
1013
1014
1015
1016
1017
1018
1019
1020
1021
1022
1023
1024
1025
1026
1027
1028
1029
1030
1031
1032
1033
1034
1035
1036
1037
1038
1039
1040
1041
1042
1043
1044
1045
1046
1047
1048
1049
1050
1051
1052
1053
1054
1055
1056
1057
1058
1059
1060
1061
1062
1063
1064
1065
1066
1067
1068
1069
1070
1071
1072
1073
1074
1075
1076
1077
1078
1079
1080
1081
1082
1083
1084
1085
1086
1087
1088
1089
1090
1091
1092
1093
1094
1095
1096
1097
1098
1099
1100
1101
1102
1103
1104
1105
1106
1107
1108
1109
1110
1111
1112
1113
1114
1115
1116
1117
1118
1119
1120
1121
1122
1123
1124
1125
1126
1127
1128
1129
1130
1131
1132
1133
1134
1135
1136
1137
1138
1139
1140
1141
1142
1143
1144
1145
1146
1147
1148
1149
1150
1151
1152
1153
1154
1155
1156
1157
1158
1159
1160
1161
pv_rhs(ji,jj,jk) = pv_rhs(ji,jj,jk) + zva
END_2D
! ! ===============
END DO ! End of slab
! ! ===============
END SUBROUTINE vor_eeT
SUBROUTINE dyn_vor_init
!!---------------------------------------------------------------------
!! *** ROUTINE dyn_vor_init ***
!!
!! ** Purpose : Control the consistency between cpp options for
!! tracer advection schemes
!!----------------------------------------------------------------------
INTEGER :: ji, jj, jk ! dummy loop indices
INTEGER :: ioptio, ios ! local integer
REAL(wp) :: zmsk ! local scalars
!!
NAMELIST/namdyn_vor/ ln_dynvor_ens, ln_dynvor_ene, ln_dynvor_enT, ln_dynvor_eeT, &
& ln_dynvor_een, nn_e3f_typ , ln_dynvor_mix, ln_dynvor_msk
!!----------------------------------------------------------------------
!
IF(lwp) THEN
WRITE(numout,*)
WRITE(numout,*) 'dyn_vor_init : vorticity term : read namelist and control the consistency'
WRITE(numout,*) '~~~~~~~~~~~~'
ENDIF
!
READ ( numnam_ref, namdyn_vor, IOSTAT = ios, ERR = 901)
901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdyn_vor in reference namelist' )
READ ( numnam_cfg, namdyn_vor, IOSTAT = ios, ERR = 902 )
902 IF( ios > 0 ) CALL ctl_nam ( ios , 'namdyn_vor in configuration namelist' )
IF(lwm) WRITE ( numond, namdyn_vor )
!
IF(lwp) THEN ! Namelist print
WRITE(numout,*) ' Namelist namdyn_vor : choice of the vorticity term scheme'
WRITE(numout,*) ' enstrophy conserving scheme ln_dynvor_ens = ', ln_dynvor_ens
WRITE(numout,*) ' f-point energy conserving scheme ln_dynvor_ene = ', ln_dynvor_ene
WRITE(numout,*) ' t-point energy conserving scheme ln_dynvor_enT = ', ln_dynvor_enT
WRITE(numout,*) ' energy conserving scheme (een using e3t) ln_dynvor_eeT = ', ln_dynvor_eeT
WRITE(numout,*) ' enstrophy and energy conserving scheme ln_dynvor_een = ', ln_dynvor_een
WRITE(numout,*) ' e3f = averaging /4 (=0) or /sum(tmask) (=1) nn_e3f_typ = ', nn_e3f_typ
WRITE(numout,*) ' mixed enstrophy/energy conserving scheme ln_dynvor_mix = ', ln_dynvor_mix
WRITE(numout,*) ' masked (=T) or unmasked(=F) vorticity ln_dynvor_msk = ', ln_dynvor_msk
ENDIF
!!gm this should be removed when choosing a unique strategy for fmask at the coast
! If energy, enstrophy or mixed advection of momentum in vector form change the value for masks
! at angles with three ocean points and one land point
IF(lwp) WRITE(numout,*)
IF(lwp) WRITE(numout,*) ' change fmask value in the angles (T) ln_vorlat = ', ln_vorlat
IF( ln_vorlat .AND. ( ln_dynvor_ene .OR. ln_dynvor_ens .OR. ln_dynvor_mix ) ) THEN
DO_3D( 1, 0, 1, 0, 1, jpk )
IF( tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) &
& + tmask(ji,jj ,jk) + tmask(ji+1,jj ,jk) == 3._wp ) fmask(ji,jj,jk) = 1._wp
END_3D
!
CALL lbc_lnk( 'dynvor', fmask, 'F', 1._wp ) ! Lateral boundary conditions on fmask
!
ENDIF
!!gm end
ioptio = 0 ! type of scheme for vorticity (set nvor_scheme)
IF( ln_dynvor_ens ) THEN ; ioptio = ioptio + 1 ; nvor_scheme = np_ENS ; ENDIF
IF( ln_dynvor_ene ) THEN ; ioptio = ioptio + 1 ; nvor_scheme = np_ENE ; ENDIF
IF( ln_dynvor_enT ) THEN ; ioptio = ioptio + 1 ; nvor_scheme = np_ENT ; ENDIF
IF( ln_dynvor_eeT ) THEN ; ioptio = ioptio + 1 ; nvor_scheme = np_EET ; ENDIF
IF( ln_dynvor_een ) THEN ; ioptio = ioptio + 1 ; nvor_scheme = np_EEN ; ENDIF
IF( ln_dynvor_mix ) THEN ; ioptio = ioptio + 1 ; nvor_scheme = np_MIX ; ENDIF
!
IF( ioptio /= 1 ) CALL ctl_stop( ' use ONE and ONLY one vorticity scheme' )
!
IF(lwp) WRITE(numout,*) ! type of calculated vorticity (set ncor, nrvm, ntot)
ncor = np_COR ! planetary vorticity
SELECT CASE( n_dynadv )
CASE( np_LIN_dyn )
IF(lwp) WRITE(numout,*) ' ==>>> linear dynamics : total vorticity = Coriolis'
nrvm = np_COR ! planetary vorticity
ntot = np_COR ! - -
CASE( np_VEC_c2 )
IF(lwp) WRITE(numout,*) ' ==>>> vector form dynamics : total vorticity = Coriolis + relative vorticity'
nrvm = np_RVO ! relative vorticity
ntot = np_CRV ! relative + planetary vorticity
CASE( np_FLX_c2 , np_FLX_ubs )
IF(lwp) WRITE(numout,*) ' ==>>> flux form dynamics : total vorticity = Coriolis + metric term'
nrvm = np_MET ! metric term
ntot = np_CME ! Coriolis + metric term
!
SELECT CASE( nvor_scheme ) ! pre-computed gradients for the metric term:
CASE( np_ENT ) !* T-point metric term : pre-compute di(e2u)/2 and dj(e1v)/2
ALLOCATE( di_e2u_2(jpi,jpj), dj_e1v_2(jpi,jpj) )
DO_2D( 0, 0, 0, 0 )
di_e2u_2(ji,jj) = ( e2u(ji,jj) - e2u(ji-1,jj ) ) * 0.5_wp
dj_e1v_2(ji,jj) = ( e1v(ji,jj) - e1v(ji ,jj-1) ) * 0.5_wp
END_2D
CALL lbc_lnk( 'dynvor', di_e2u_2, 'T', -1.0_wp , dj_e1v_2, 'T', -1.0_wp ) ! Lateral boundary conditions
!
CASE DEFAULT !* F-point metric term : pre-compute di(e2u)/(2*e1e2f) and dj(e1v)/(2*e1e2f)
ALLOCATE( di_e2v_2e1e2f(jpi,jpj), dj_e1u_2e1e2f(jpi,jpj) )
DO_2D( 1, 0, 1, 0 )
di_e2v_2e1e2f(ji,jj) = ( e2v(ji+1,jj ) - e2v(ji,jj) ) * 0.5 * r1_e1e2f(ji,jj)
dj_e1u_2e1e2f(ji,jj) = ( e1u(ji ,jj+1) - e1u(ji,jj) ) * 0.5 * r1_e1e2f(ji,jj)
END_2D
CALL lbc_lnk( 'dynvor', di_e2v_2e1e2f, 'F', -1.0_wp , dj_e1u_2e1e2f, 'F', -1.0_wp ) ! Lateral boundary conditions
END SELECT
!
END SELECT
#if defined key_qco || defined key_linssh
SELECT CASE( nvor_scheme ) ! qco or linssh cases : pre-computed a specific e3f_0 for some vorticity schemes
CASE( np_ENS , np_ENE , np_EEN , np_MIX )
!
ALLOCATE( e3f_0vor(jpi,jpj,jpk) )
!
SELECT CASE( nn_e3f_typ )
CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4)
DO_3D( 0, 0, 0, 0, 1, jpk )
e3f_0vor(ji,jj,jk) = ( e3t_0(ji ,jj+1,jk)*tmask(ji ,jj+1,jk) &
& + e3t_0(ji+1,jj+1,jk)*tmask(ji+1,jj+1,jk) &
& + e3t_0(ji ,jj ,jk)*tmask(ji ,jj ,jk) &
& + e3t_0(ji+1,jj ,jk)*tmask(ji+1,jj ,jk) ) * 0.25_wp
END_3D
CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask)
DO_3D( 0, 0, 0, 0, 1, jpk )
zmsk = (tmask(ji,jj+1,jk) +tmask(ji+1,jj+1,jk) &
& + tmask(ji,jj ,jk) +tmask(ji+1,jj ,jk) )
!
IF( zmsk /= 0._wp ) THEN
e3f_0vor(ji,jj,jk) = ( e3t_0(ji ,jj+1,jk)*tmask(ji ,jj+1,jk) &
& + e3t_0(ji+1,jj+1,jk)*tmask(ji+1,jj+1,jk) &
& + e3t_0(ji ,jj ,jk)*tmask(ji ,jj ,jk) &
& + e3t_0(ji+1,jj ,jk)*tmask(ji+1,jj ,jk) ) / zmsk
ELSE ; e3f_0vor(ji,jj,jk) = 0._wp
ENDIF
END_3D
END SELECT
!
CALL lbc_lnk( 'dynvor', e3f_0vor, 'F', 1._wp )
! ! insure e3f_0vor /= 0
WHERE( e3f_0vor(:,:,:) == 0._wp ) e3f_0vor(:,:,:) = e3f_0(:,:,:)
!
END SELECT
!
#endif
IF(lwp) THEN ! Print the choice
WRITE(numout,*)
SELECT CASE( nvor_scheme )
CASE( np_ENS ) ; WRITE(numout,*) ' ==>>> enstrophy conserving scheme (ENS)'
CASE( np_ENE ) ; WRITE(numout,*) ' ==>>> energy conserving scheme (Coriolis at F-points) (ENE)'
CASE( np_ENT ) ; WRITE(numout,*) ' ==>>> energy conserving scheme (Coriolis at T-points) (ENT)'
IF( ln_dynadv_vec ) CALL ctl_warn('dyn_vor_init: ENT scheme may not work in vector form')
CASE( np_EET ) ; WRITE(numout,*) ' ==>>> energy conserving scheme (EEN scheme using e3t) (EET)'
CASE( np_EEN ) ; WRITE(numout,*) ' ==>>> energy and enstrophy conserving scheme (EEN)'
CASE( np_MIX ) ; WRITE(numout,*) ' ==>>> mixed enstrophy/energy conserving scheme (MIX)'
END SELECT
ENDIF
!
END SUBROUTINE dyn_vor_init
!!==============================================================================
END MODULE dynvor