ROMS
Loading...
Searching...
No Matches
ini_adjust_mod Module Reference

Functions/Subroutines

subroutine, public ini_adjust (ng, tile, linp, lout)
 
subroutine ini_adjust_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, linp, lout, rmask, umask, vmask, tl_t, tl_u, tl_v, tl_ubar, tl_vbar, tl_zeta, t, u, v, ubar, vbar, zeta)
 
subroutine, public rp_ini_adjust (ng, tile, linp, lout)
 
subroutine rp_ini_adjust_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, linp, lout, rmask, umask, vmask, ad_t, ad_u, ad_v, ad_zeta, t, u, v, zeta, tl_t, tl_u, tl_v, tl_zeta)
 
subroutine, public load_adtotl (ng, tile, linp, lout, add)
 
subroutine load_adtotl_tile (ng, tile, lbi, ubi, lbj, ubj, lbij, ubij, imins, imaxs, jmins, jmaxs, linp, lout, add, rmask, umask, vmask, ad_t_obc, ad_u_obc, ad_v_obc, ad_ubar_obc, ad_vbar_obc, ad_zeta_obc, ad_ustr, ad_vstr, ad_tflux, ad_t, ad_u, ad_v, ad_ubar, ad_vbar, ad_zeta, tl_t_obc, tl_u_obc, tl_v_obc, tl_ubar_obc, tl_vbar_obc, tl_zeta_obc, tl_ustr, tl_vstr, tl_tflux, tl_t, tl_u, tl_v, tl_ubar, tl_vbar, tl_zeta)
 
subroutine, public load_tltoad (ng, tile, linp, lout, add)
 
subroutine load_tltoad_tile (ng, tile, lbi, ubi, lbj, ubj, lbij, ubij, imins, imaxs, jmins, jmaxs, linp, lout, add, rmask, umask, vmask, tl_t_obc, tl_u_obc, tl_v_obc, tl_ubar_obc, tl_vbar_obc, tl_zeta_obc, tl_ustr, tl_vstr, tl_tflux, tl_t, tl_u, tl_v, tl_ubar, tl_vbar, tl_zeta, ad_t_obc, ad_u_obc, ad_v_obc, ad_ubar_obc, ad_vbar_obc, ad_zeta_obc, ad_ustr, ad_vstr, ad_tflux, ad_t, ad_u, ad_v, ad_ubar, ad_vbar, ad_zeta)
 
subroutine, public ini_perturb (ng, tile, linp, lout)
 
subroutine ini_perturb_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nstp, nnew, linp, lout, rmask, umask, vmask, bed_thick, hz, h, om_v, on_u, zice, z_r, z_w, zt_avg1, ad_t, ad_u, ad_v, t, u, v, ad_ubar, ad_vbar, ad_zeta, ubar, vbar, zeta)
 
subroutine, public ad_ini_perturb (ng, tile, kinp, kout, ninp, nout)
 
subroutine ad_ini_perturb_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, kinp, kout, ninp, nout, rmask, umask, vmask, e_t, e_u, e_v, e_zeta, hz, tl_zeta, tl_t, tl_u, tl_v, ad_t, ad_u, ad_v, ad_zeta)
 
subroutine, public tl_ini_perturb (ng, tile, linp, lout)
 
subroutine tl_ini_perturb_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nstp, nnew, linp, lout, rmask, umask, vmask, tl_bed_thick, tl_hz, h, tl_h, om_v, on_u, zice, tl_z_r, tl_z_w, zt_avg1, tl_zt_avg1, ubar, vbar, zeta, ad_t, ad_u, ad_v, tl_t, tl_u, tl_v, ad_ubar, ad_vbar, ad_zeta, tl_ubar, tl_vbar, tl_zeta)
 

Function/Subroutine Documentation

◆ ad_ini_perturb()

subroutine, public ini_adjust_mod::ad_ini_perturb ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) kinp,
integer, intent(in) kout,
integer, intent(in) ninp,
integer, intent(in) nout )

Definition at line 2974 of file ini_adjust.F.

2975!
2976!=======================================================================
2977! !
2978! Initialize adjoint state variables with tangent linear state scaled !
2979! by the energy norm, as required by the Generalized stability Theory !
2980! propagators. !
2981! !
2982!=======================================================================
2983!
2984 USE mod_param
2985# ifdef SOLVE3D
2986 USE mod_coupling
2987# endif
2988 USE mod_grid
2989 USE mod_ocean
2990!
2991! Imported variable declarations.
2992!
2993 integer, intent(in) :: ng, tile, Kinp, Kout, Ninp, Nout
2994!
2995! Local variable declarations.
2996!
2997 character (len=*), parameter :: MyFile = &
2998 & __FILE__//", ad_ini_perturb"
2999!
3000# include "tile.h"
3001!
3002# ifdef PROFILE
3003 CALL wclock_on (ng, iadm, 7, __line__, myfile)
3004# endif
3005 CALL ad_ini_perturb_tile (ng, tile, &
3006 & lbi, ubi, lbj, ubj, &
3007 & imins, imaxs, jmins, jmaxs, &
3008 & kinp, kout, ninp, nout, &
3009# ifdef MASKING
3010 & grid(ng) % rmask, &
3011 & grid(ng) % umask, &
3012 & grid(ng) % vmask, &
3013# endif
3014# ifdef BNORM
3015# ifdef SOLVE3D
3016 & ocean(ng) % e_t, &
3017 & ocean(ng) % e_u, &
3018 & ocean(ng) % e_v, &
3019# else
3020 & ocean(ng) % e_ubar, &
3021 & ocean(ng) % e_vbar, &
3022# endif
3023 & ocean(ng) % e_zeta, &
3024# endif
3025# ifdef SOLVE3D
3026 & grid(ng) % Hz, &
3027# else
3028 & grid(ng) % h, &
3029 & ocean(ng) % tl_ubar, &
3030 & ocean(ng) % tl_vbar, &
3031# endif
3032 & ocean(ng) % tl_zeta, &
3033# ifdef SOLVE3D
3034 & ocean(ng) % tl_t, &
3035 & ocean(ng) % tl_u, &
3036 & ocean(ng) % tl_v, &
3037 & ocean(ng) % ad_t, &
3038 & ocean(ng) % ad_u, &
3039 & ocean(ng) % ad_v, &
3040# else
3041 & ocean(ng) % ad_ubar, &
3042 & ocean(ng) % ad_vbar, &
3043# endif
3044 & ocean(ng) % ad_zeta)
3045# ifdef PROFILE
3046 CALL wclock_off (ng, iadm, 7, __line__, myfile)
3047# endif
3048!
3049 RETURN
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer, parameter iadm
Definition mod_param.F:665
recursive subroutine wclock_off(ng, model, region, line, routine)
Definition timers.F:148
recursive subroutine wclock_on(ng, model, region, line, routine)
Definition timers.F:3

References ad_ini_perturb_tile(), mod_grid::grid, mod_param::iadm, mod_ocean::ocean, wclock_off(), and wclock_on().

Referenced by propagator_mod::propagator_fsv(), propagator_mod::propagator_hop(), propagator_mod::propagator_hso(), propagator_mod::propagator_op(), and propagator_mod::propagator_so().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ ad_ini_perturb_tile()

subroutine ini_adjust_mod::ad_ini_perturb_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) kinp,
integer, intent(in) kout,
integer, intent(in) ninp,
integer, intent(in) nout,
real(r8), dimension(lbi:,lbj:), intent(in) rmask,
real(r8), dimension(lbi:,lbj:), intent(in) umask,
real(r8), dimension(lbi:,lbj:), intent(in) vmask,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(in) e_t,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) e_u,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) e_v,
real(r8), dimension(lbi:,lbj:,:), intent(in) e_zeta,
real(r8), dimension(lbi:,lbj:,:), intent(in) hz,
real(r8), dimension(lbi:,lbj:,:), intent(in) tl_zeta,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(in) tl_t,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) tl_u,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) tl_v,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) ad_t,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) ad_u,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) ad_v,
real(r8), dimension(lbi:,lbj:,:), intent(inout) ad_zeta )
private

Definition at line 3053 of file ini_adjust.F.

3082!***********************************************************************
3083!
3084 USE mod_param
3085 USE mod_scalars
3086!
3087! Imported variable declarations.
3088!
3089 integer, intent(in) :: ng, tile
3090 integer, intent(in) :: LBi, UBi, LBj, UBj
3091 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
3092 integer, intent(in) :: Kinp, Kout, Ninp, Nout
3093!
3094# ifdef ASSUMED_SHAPE
3095# ifdef MASKING
3096 real(r8), intent(in) :: rmask(LBi:,LBj:)
3097 real(r8), intent(in) :: umask(LBi:,LBj:)
3098 real(r8), intent(in) :: vmask(LBi:,LBj:)
3099# endif
3100# ifdef SOLVE3D
3101 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
3102# else
3103 real(r8), intent(in) :: h(LBi:,LBj:)
3104 real(r8), intent(in) :: tl_ubar(LBi:,LBj:,:)
3105 real(r8), intent(in) :: tl_vbar(LBi:,LBj:,:)
3106# endif
3107 real(r8), intent(in) :: tl_zeta(LBi:,LBj:,:)
3108# ifdef SOLVE3D
3109 real(r8), intent(in) :: tl_t(LBi:,LBj:,:,:,:)
3110 real(r8), intent(in) :: tl_u(LBi:,LBj:,:,:)
3111 real(r8), intent(in) :: tl_v(LBi:,LBj:,:,:)
3112# else
3113 real(r8), intent(inout) :: ad_ubar(LBi:,LBj:,:)
3114 real(r8), intent(inout) :: ad_vbar(LBi:,LBj:,:)
3115# endif
3116 real(r8), intent(inout) :: ad_zeta(LBi:,LBj:,:)
3117# ifdef SOLVE3D
3118 real(r8), intent(inout) :: ad_t(LBi:,LBj:,:,:,:)
3119 real(r8), intent(inout) :: ad_u(LBi:,LBj:,:,:)
3120 real(r8), intent(inout) :: ad_v(LBi:,LBj:,:,:)
3121# endif
3122# ifdef BNORM
3123 real(r8), intent(in) :: e_zeta(LBi:,LBj:,:)
3124# ifdef SOLVE3D
3125 real(r8), intent(in) :: e_t(LBi:,LBj:,:,:,:)
3126 real(r8), intent(in) :: e_u(LBi:,LBj:,:,:)
3127 real(r8), intent(in) :: e_v(LBi:,LBj:,:,:)
3128# else
3129 real(r8), intent(in) :: e_ubar(LBi:,LBj:,:)
3130 real(r8), intent(in) :: e_vbar(LBi:,LBj:,:)
3131# endif
3132# endif
3133# else
3134# ifdef MASKING
3135 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
3136 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
3137 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
3138# endif
3139# ifdef SOLVE3D
3140 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
3141# else
3142 real(r8), intent(in) :: h(LBi:UBi,LBj:UBj)
3143 real(r8), intent(in) :: tl_ubar(LBi:UBi,LBj:UBj,:)
3144 real(r8), intent(in) :: tl_vbar(LBi:UBi,LBj:UBj,:)
3145# endif
3146 real(r8), intent(in) :: tl_zeta(LBi:UBi,LBj:UBj,:)
3147# ifdef SOLVE3D
3148 real(r8), intent(in) :: tl_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
3149 real(r8), intent(in) :: tl_u(LBi:UBi,LBj:UBj,N(ng),2)
3150 real(r8), intent(in) :: tl_v(LBi:UBi,LBj:UBj,N(ng),2)
3151# else
3152 real(r8), intent(inout) :: ad_ubar(LBi:UBi,LBj:UBj,:)
3153 real(r8), intent(inout) :: ad_vbar(LBi:UBi,LBj:UBj,:)
3154# endif
3155 real(r8), intent(inout) :: ad_zeta(LBi:UBi,LBj:UBj,:)
3156# ifdef SOLVE3D
3157 real(r8), intent(inout) :: ad_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
3158 real(r8), intent(inout) :: ad_u(LBi:UBi,LBj:UBj,N(ng),2)
3159 real(r8), intent(inout) :: ad_v(LBi:UBi,LBj:UBj,N(ng),2)
3160# endif
3161# ifdef BNORM
3162 real(r8), intent(in) :: e_zeta(LBi:UBi,LBj:UBj,NSA)
3163# ifdef SOLVE3D
3164 real(r8), intent(in) :: e_t(LBi:UBi,LBj:UBj,N(ng),NSA,NT(ng))
3165 real(r8), intent(in) :: e_u(LBi:UBi,LBj:UBj,N(ng),NSA)
3166 real(r8), intent(in) :: e_v(LBi:UBi,LBj:UBj,N(ng),NSA)
3167# else
3168 real(r8), intent(in) :: e_ubar(LBi:UBi,LBj:UBj,NSA)
3169 real(r8), intent(in) :: e_vbar(LBi:UBi,LBj:UBj,NSA)
3170# endif
3171# endif
3172# endif
3173!
3174! Local variable declarations.
3175!
3176 integer :: i, j
3177# ifdef SOLVE3D
3178 integer :: itrc, k
3179# endif
3180 real(r8) :: cff, scale
3181!
3182# include "set_bounds.h"
3183!
3184!-----------------------------------------------------------------------
3185! Initialize adjoint state variables with tangent linear state scaled
3186! by the energy norm.
3187!-----------------------------------------------------------------------
3188
3189# ifdef FULL_GRID
3190# define I_R_RANGE IstrT,IendT
3191# define I_U_RANGE IstrP,IendT
3192# define J_R_RANGE JstrT,JendT
3193# define J_V_RANGE JstrP,JendT
3194# else
3195# define I_R_RANGE Istr,Iend
3196# define I_U_RANGE IstrU,Iend
3197# define J_R_RANGE Jstr,Jend
3198# define J_V_RANGE JstrV,Jend
3199# endif
3200!
3201! Free-surface.
3202!
3203# ifndef BNORM
3204 scale=0.5_r8*g*rho0
3205# endif
3206 DO j=j_r_range
3207 DO i=i_r_range
3208# ifdef BNORM
3209 IF (e_zeta(i,j,1).gt.0.0_r8) THEN
3210 scale=1.0_r8/(e_zeta(i,j,1)*e_zeta(i,j,1))
3211 ELSE
3212 scale=1.0_r8
3213 END IF
3214# endif
3215 ad_zeta(i,j,kout)=scale*tl_zeta(i,j,kinp)
3216# ifdef MASKING
3217 ad_zeta(i,j,kout)=ad_zeta(i,j,kout)*rmask(i,j)
3218# endif
3219 END DO
3220 END DO
3221
3222# ifndef SOLVE3D
3223!
3224! 2D u-momentum component.
3225!
3226 cff=0.25_r8*rho0
3227 DO j=j_r_range
3228 DO i=i_u_range
3229# ifdef BNORM
3230 IF (e_ubar(i,j,1).gt.0.0_r8) THEN
3231 scale=1.0_r8/(e_ubar(i,j,1)*e_ubar(i,j,1))
3232 ELSE
3233 scale=1.0_r8
3234 ENDIF
3235# else
3236 scale=cff*(h(i-1,j)+h(i,j))
3237# endif
3238 ad_ubar(i,j,kout)=scale*tl_ubar(i,j,kinp)
3239# ifdef MASKING
3240 ad_ubar(i,j,kout)=ad_ubar(i,j,kout)*umask(i,j)
3241# endif
3242 END DO
3243 END DO
3244!
3245! 2D v-momentum component.
3246!
3247 cff=0.25_r8*rho0
3248 DO j=j_v_range
3249 DO i=i_r_range
3250# ifdef BNORM
3251 IF (e_vbar(i,j,1).gt.0.0_r8) THEN
3252 scale=1.0_r8/(e_vbar(i,j,1)*e_vbar(i,j,1))
3253 ELSE
3254 scale=1.0_r8
3255 ENDIF
3256# else
3257 scale=cff*(h(i,j-1)+h(i,j))
3258# endif
3259 ad_vbar(i,j,kout)=scale*tl_vbar(i,j,kinp)
3260# ifdef MASKING
3261 ad_vbar(i,j,kout)=ad_vbar(i,j,kout)*vmask(i,j)
3262# endif
3263 END DO
3264 END DO
3265# else
3266!
3267! 3D u-momentum component.
3268!
3269 cff=0.25_r8*rho0
3270 DO k=1,n(ng)
3271 DO j=j_r_range
3272 DO i=i_u_range
3273# ifdef BNORM
3274 IF (e_u(i,j,k,1).gt.0.0_r8) THEN
3275 scale=1.0_r8/(e_u(i,j,k,1)*e_u(i,j,k,1))
3276 ELSE
3277 scale=1.0_r8
3278 END IF
3279# else
3280 scale=cff*(hz(i-1,j,k)+hz(i,j,k))
3281# endif
3282 ad_u(i,j,k,nout)=scale*tl_u(i,j,k,ninp)
3283# ifdef MASKING
3284 ad_u(i,j,k,nout)=ad_u(i,j,k,nout)*umask(i,j)
3285# endif
3286 END DO
3287 END DO
3288 END DO
3289!
3290! 3D v-momentum component.
3291!
3292 cff=0.25_r8*rho0
3293 DO k=1,n(ng)
3294 DO j=j_v_range
3295 DO i=i_r_range
3296# ifdef BNORM
3297 IF (e_v(i,j,k,1).gt.0.0_r8) THEN
3298 scale=1.0_r8/(e_v(i,j,k,1)*e_v(i,j,k,1))
3299 ELSE
3300 scale=1.0_r8
3301 END IF
3302# else
3303 scale=cff*(hz(i,j-1,k)+hz(i,j,k))
3304# endif
3305 ad_v(i,j,k,nout)=scale*tl_v(i,j,k,ninp)
3306# ifdef MASKING
3307 ad_v(i,j,k,nout)=ad_v(i,j,k,nout)*vmask(i,j)
3308# endif
3309 END DO
3310 END DO
3311 END DO
3312!
3313! Tracers. For now, use salinity scale for passive tracers.
3314!
3315 DO itrc=1,nt(ng)
3316 IF (itrc.eq.itemp) THEN
3317 cff=0.5_r8*rho0*tcoef(ng)*tcoef(ng)*g*g/bvf_bak
3318 ELSE IF (itrc.eq.isalt) THEN
3319 cff=0.5_r8*rho0*scoef(ng)*scoef(ng)*g*g/bvf_bak
3320 ELSE
3321 cff=0.5_r8*rho0*scoef(ng)*scoef(ng)*g*g/bvf_bak
3322 END IF
3323 DO k=1,n(ng)
3324 DO j=j_r_range
3325 DO i=i_r_range
3326# ifdef BNORM
3327 IF (e_t(i,j,k,1,itrc).gt.0.0_r8) THEN
3328 scale=1.0_r8/(e_t(i,j,k,1,itrc)*e_t(i,j,k,1,itrc))
3329 ELSE
3330 scale=1.0_r8
3331 END IF
3332# else
3333 scale=cff*hz(i,j,k)
3334# endif
3335 ad_t(i,j,k,nout,itrc)=scale*tl_t(i,j,k,ninp,itrc)
3336# ifdef MASKING
3337 ad_t(i,j,k,nout,itrc)=ad_t(i,j,k,nout,itrc)*rmask(i,j)
3338# endif
3339 END DO
3340 END DO
3341 END DO
3342 END DO
3343# endif
3344
3345# undef I_R_RANGE
3346# undef I_U_RANGE
3347# undef J_R_RANGE
3348# undef J_V_RANGE
3349!
3350 RETURN
integer, dimension(:), allocatable n
Definition mod_param.F:479
integer, dimension(:), allocatable nt
Definition mod_param.F:489
real(dp) bvf_bak
real(r8), dimension(:), allocatable tcoef
integer isalt
integer itemp
real(dp) g
real(dp) rho0
real(r8), dimension(:), allocatable scoef

References mod_scalars::bvf_bak, mod_scalars::g, mod_scalars::isalt, mod_scalars::itemp, mod_scalars::rho0, mod_scalars::scoef, and mod_scalars::tcoef.

Referenced by ad_ini_perturb().

Here is the caller graph for this function:

◆ ini_adjust()

subroutine public ini_adjust_mod::ini_adjust ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) linp,
integer, intent(in) lout )

Definition at line 68 of file ini_adjust.F.

69!
70!=======================================================================
71! !
72! This routine adds 4D-Var inner loop tangent linear increments to !
73! nonlinear state initial conditions. The boundary condition and !
74! barotropic/baroclinic coupling, if any, are processed latter in !
75! routine "ini_fields" before time-stepping. !
76! !
77! On Input: !
78! !
79! ng Nested grid number. !
80! Linp Tangent linear state time index to add. !
81! Lout Nonlinear state time index to update. !
82! tile Domain partition. !
83! !
84!=======================================================================
85!
86 USE mod_param
87 USE mod_grid
88 USE mod_ocean
89!
90! Imported variable declarations.
91!
92 integer, intent(in) :: ng, tile, Linp, Lout
93!
94! Local variable declarations.
95!
96 character (len=*), parameter :: MyFile = &
97 & __FILE__
98!
99# include "tile.h"
100!
101# ifdef PROFILE
102 CALL wclock_on (ng, inlm, 7, __line__, myfile)
103# endif
104 CALL ini_adjust_tile (ng, tile, &
105 & lbi, ubi, lbj, ubj, &
106 & imins, imaxs, jmins, jmaxs, &
107 & linp, lout, &
108# ifdef MASKING
109 & grid(ng) % rmask, &
110 & grid(ng) % umask, &
111 & grid(ng) % vmask, &
112# endif
113# ifdef SOLVE3D
114 & ocean(ng) % tl_t, &
115 & ocean(ng) % tl_u, &
116 & ocean(ng) % tl_v, &
117# endif
118 & ocean(ng) % tl_ubar, &
119 & ocean(ng) % tl_vbar, &
120 & ocean(ng) % tl_zeta, &
121# ifdef SOLVE3D
122 & ocean(ng) % t, &
123 & ocean(ng) % u, &
124 & ocean(ng) % v, &
125# endif
126 & ocean(ng) % ubar, &
127 & ocean(ng) % vbar, &
128 & ocean(ng) % zeta)
129# ifdef PROFILE
130 CALL wclock_off (ng, inlm, 7, __line__, myfile)
131# endif
132!
133 RETURN
integer, parameter inlm
Definition mod_param.F:662

References mod_grid::grid, ini_adjust_tile(), mod_param::inlm, mod_ocean::ocean, wclock_off(), and wclock_on().

Referenced by i4dvar_mod::analysis(), convolve_mod::error_covariance(), and rbl4dvar_mod::increment().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ ini_adjust_tile()

subroutine ini_adjust_mod::ini_adjust_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) linp,
integer, intent(in) lout,
real(r8), dimension(lbi:,lbj:), intent(in) rmask,
real(r8), dimension(lbi:,lbj:), intent(in) umask,
real(r8), dimension(lbi:,lbj:), intent(in) vmask,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(in) tl_t,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) tl_u,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) tl_v,
real(r8), dimension(lbi:,lbj:,:), intent(in) tl_ubar,
real(r8), dimension(lbi:,lbj:,:), intent(in) tl_vbar,
real(r8), dimension(lbi:,lbj:,:), intent(in) tl_zeta,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) t,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) u,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) v,
real(r8), dimension(lbi:,lbj:,:), intent(inout) ubar,
real(r8), dimension(lbi:,lbj:,:), intent(inout) vbar,
real(r8), dimension(lbi:,lbj:,:), intent(inout) zeta )
private

Definition at line 137 of file ini_adjust.F.

152!***********************************************************************
153!
154 USE mod_param
155!
156! Imported variable declarations.
157!
158 integer, intent(in) :: ng, tile
159 integer, intent(in) :: LBi, UBi, LBj, UBj
160 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
161 integer, intent(in) :: Linp, Lout
162!
163# ifdef ASSUMED_SHAPE
164# ifdef MASKING
165 real(r8), intent(in) :: rmask(LBi:,LBj:)
166 real(r8), intent(in) :: umask(LBi:,LBj:)
167 real(r8), intent(in) :: vmask(LBi:,LBj:)
168# endif
169# ifdef SOLVE3D
170 real(r8), intent(in) :: tl_t(LBi:,LBj:,:,:,:)
171 real(r8), intent(in) :: tl_u(LBi:,LBj:,:,:)
172 real(r8), intent(in) :: tl_v(LBi:,LBj:,:,:)
173# endif
174 real(r8), intent(in) :: tl_ubar(LBi:,LBj:,:)
175 real(r8), intent(in) :: tl_vbar(LBi:,LBj:,:)
176 real(r8), intent(in) :: tl_zeta(LBi:,LBj:,:)
177# ifdef SOLVE3D
178 real(r8), intent(inout) :: t(LBi:,LBj:,:,:,:)
179 real(r8), intent(inout) :: u(LBi:,LBj:,:,:)
180 real(r8), intent(inout) :: v(LBi:,LBj:,:,:)
181# endif
182 real(r8), intent(inout) :: ubar(LBi:,LBj:,:)
183 real(r8), intent(inout) :: vbar(LBi:,LBj:,:)
184 real(r8), intent(inout) :: zeta(LBi:,LBj:,:)
185# else
186# ifdef MASKING
187 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
188 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
189 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
190# endif
191# ifdef SOLVE3D
192 real(r8), intent(in) :: tl_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
193 real(r8), intent(in) :: tl_u(LBi:UBi,LBj:UBj,N(ng),2)
194 real(r8), intent(in) :: tl_v(LBi:UBi,LBj:UBj,N(ng),2)
195# endif
196 real(r8), intent(in) :: tl_ubar(LBi:UBi,LBj:UBj,:)
197 real(r8), intent(in) :: tl_vbar(LBi:UBi,LBj:UBj,:)
198 real(r8), intent(in) :: tl_zeta(LBi:UBi,LBj:UBj,:)
199# ifdef SOLVE3D
200 real(r8), intent(inout) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
201 real(r8), intent(inout) :: u(LBi:UBi,LBj:UBj,N(ng),2)
202 real(r8), intent(inout) :: v(LBi:UBi,LBj:UBj,N(ng),2)
203# endif
204 real(r8), intent(inout) :: ubar(LBi:UBi,LBj:UBj,:)
205 real(r8), intent(inout) :: vbar(LBi:UBi,LBj:UBj,:)
206 real(r8), intent(inout) :: zeta(LBi:UBi,LBj:UBj,:)
207# endif
208!
209! Local variable declarations.
210!
211 integer :: i, j
212# ifdef SOLVE3D
213 integer :: itrc, k
214# endif
215
216# include "set_bounds.h"
217!
218!-----------------------------------------------------------------------
219! Adjust initial conditions for 2D state variables by adding tangent
220! linear increments from data assimilation.
221!-----------------------------------------------------------------------
222!
223# ifndef SOLVE3D
224 DO j=jstrt,jendt
225 DO i=istrp,iendt
226 ubar(i,j,lout)=ubar(i,j,lout)+tl_ubar(i,j,linp)
227# ifdef MASKING
228 ubar(i,j,lout)=ubar(i,j,lout)*umask(i,j)
229# endif
230 END DO
231 END DO
232
233 DO j=jstrp,jendt
234 DO i=istrt,iendt
235 vbar(i,j,lout)=vbar(i,j,lout)+tl_vbar(i,j,linp)
236# ifdef MASKING
237 vbar(i,j,lout)=vbar(i,j,lout)*vmask(i,j)
238# endif
239 END DO
240 END DO
241# endif
242
243 DO j=jstrt,jendt
244 DO i=istrt,iendt
245 zeta(i,j,lout)=zeta(i,j,lout)+tl_zeta(i,j,linp)
246# ifdef MASKING
247 zeta(i,j,lout)=zeta(i,j,lout)*rmask(i,j)
248# endif
249 END DO
250 END DO
251
252# ifdef SOLVE3D
253!
254!-----------------------------------------------------------------------
255! Adjust initial conditions for 3D state variables by adding tangent
256! linear increments from data assimilation.
257!-----------------------------------------------------------------------
258!
259 DO k=1,n(ng)
260 DO j=jstrt,jendt
261 DO i=istrp,iendt
262 u(i,j,k,lout)=u(i,j,k,lout)+tl_u(i,j,k,linp)
263# ifdef MASKING
264 u(i,j,k,lout)=u(i,j,k,lout)*umask(i,j)
265# endif
266 END DO
267 END DO
268
269 DO j=jstrp,jendt
270 DO i=istrt,iendt
271 v(i,j,k,lout)=v(i,j,k,lout)+tl_v(i,j,k,linp)
272# ifdef MASKING
273 v(i,j,k,lout)=v(i,j,k,lout)*vmask(i,j)
274# endif
275 END DO
276 END DO
277 END DO
278
279 DO itrc=1,nt(ng)
280 DO k=1,n(ng)
281 DO j=jstrt,jendt
282 DO i=istrt,iendt
283 t(i,j,k,lout,itrc)=t(i,j,k,lout,itrc)+ &
284 & tl_t(i,j,k,linp,itrc)
285# ifdef MASKING
286 t(i,j,k,lout,itrc)=t(i,j,k,lout,itrc)*rmask(i,j)
287# endif
288 END DO
289 END DO
290 END DO
291 END DO
292# endif
293!
294 RETURN

References mod_grid::grid, ini_adjust_tile(), mod_param::inlm, mod_ocean::ocean, wclock_off(), and wclock_on().

Referenced by ini_adjust(), and ini_adjust_tile().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ ini_perturb()

subroutine, public ini_adjust_mod::ini_perturb ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) linp,
integer, intent(in) lout )

Definition at line 2450 of file ini_adjust.F.

2451!
2452!=======================================================================
2453! !
2454! Add a perturbation to nonlinear state variables according to the !
2455! outer and inner loop iterations. The added term is a function of !
2456! the steepest descent direction (grad(J)) times the perturbation !
2457! amplitude "p" (controlled by inner loop). !
2458! !
2459!=======================================================================
2460!
2461 USE mod_param
2462# ifdef SOLVE3D
2463 USE mod_coupling
2464# endif
2465 USE mod_grid
2466 USE mod_ocean
2467# if defined SEDIMENT && defined SED_MORPH && defined SOLVE3D
2468 USE mod_sedbed
2469# endif
2470 USE mod_stepping
2471!
2472! Imported variable declarations.
2473!
2474 integer, intent(in) :: ng, Linp, Lout, tile
2475!
2476! Local variable declarations.
2477!
2478 character (len=*), parameter :: MyFile = &
2479 & __FILE__//", ini_perturb"
2480!
2481# include "tile.h"
2482!
2483# ifdef PROFILE
2484 CALL wclock_on (ng, inlm, 7, __line__, myfile)
2485# endif
2486 CALL ini_perturb_tile (ng, tile, &
2487 & lbi, ubi, lbj, ubj, &
2488 & imins, imaxs, jmins, jmaxs, &
2489 & nstp(ng), nnew(ng), linp, lout, &
2490# ifdef MASKING
2491 & grid(ng) % rmask, &
2492 & grid(ng) % umask, &
2493 & grid(ng) % vmask, &
2494# endif
2495# ifdef SOLVE3D
2496# if defined SEDIMENT && defined SED_MORPH
2497 & sedbed(ng) % bed_thick, &
2498# endif
2499 & grid(ng) % Hz, &
2500 & grid(ng) % h, &
2501 & grid(ng) % om_v, &
2502 & grid(ng) % on_u, &
2503# ifdef ICESHELF
2504 & grid(ng) % zice, &
2505# endif
2506 & grid(ng) % z_r, &
2507 & grid(ng) % z_w, &
2508 & coupling(ng) % Zt_avg1, &
2509 & ocean(ng) % ad_t, &
2510 & ocean(ng) % ad_u, &
2511 & ocean(ng) % ad_v, &
2512 & ocean(ng) % t, &
2513 & ocean(ng) % u, &
2514 & ocean(ng) % v, &
2515# endif
2516 & ocean(ng) % ad_ubar, &
2517 & ocean(ng) % ad_vbar, &
2518 & ocean(ng) % ad_zeta, &
2519 & ocean(ng) % ubar, &
2520 & ocean(ng) % vbar, &
2521 & ocean(ng) % zeta)
2522# ifdef PROFILE
2523 CALL wclock_off (ng, inlm, 7, __line__, myfile)
2524# endif
2525!
2526 RETURN
type(t_coupling), dimension(:), allocatable coupling
type(t_sedbed), dimension(:), allocatable sedbed
Definition sedbed_mod.h:157
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable nstp

References mod_coupling::coupling, mod_grid::grid, ini_perturb_tile(), mod_param::inlm, mod_stepping::nnew, mod_stepping::nstp, mod_ocean::ocean, mod_sedbed::sedbed, wclock_off(), and wclock_on().

Referenced by initial().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ ini_perturb_tile()

subroutine ini_adjust_mod::ini_perturb_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) nstp,
integer, intent(in) nnew,
integer, intent(in) linp,
integer, intent(in) lout,
real(r8), dimension(lbi:,lbj:), intent(in) rmask,
real(r8), dimension(lbi:,lbj:), intent(in) umask,
real(r8), dimension(lbi:,lbj:), intent(in) vmask,
real(r8), dimension(lbi:,lbj:,:), intent(in) bed_thick,
real(r8), dimension(lbi:,lbj:,:), intent(inout) hz,
real(r8), dimension(lbi:,lbj:), intent(inout) h,
real(r8), dimension(lbi:,lbj:), intent(in) om_v,
real(r8), dimension(lbi:,lbj:), intent(in) on_u,
real(r8), dimension(lbi:,lbj:), intent(in) zice,
real(r8), dimension(lbi:,lbj:,:), intent(out) z_r,
real(r8), dimension(lbi:,lbj:,0:), intent(out) z_w,
real(r8), dimension(lbi:,lbj:), intent(inout) zt_avg1,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(in) ad_t,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) ad_u,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) ad_v,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) t,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) u,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) v,
real(r8), dimension(lbi:,lbj:,:), intent(in) ad_ubar,
real(r8), dimension(lbi:,lbj:,:), intent(in) ad_vbar,
real(r8), dimension(lbi:,lbj:,:), intent(in) ad_zeta,
real(r8), dimension(lbi:,lbj:,:), intent(inout) ubar,
real(r8), dimension(lbi:,lbj:,:), intent(inout) vbar,
real(r8), dimension(lbi:,lbj:,:), intent(inout) zeta )
private

Definition at line 2530 of file ini_adjust.F.

2551!***********************************************************************
2552!
2553 USE mod_param
2554 USE mod_parallel
2555 USE mod_fourdvar
2556 USE mod_ncparam
2557 USE mod_iounits
2558 USE mod_scalars
2559!
2560 USE exchange_2d_mod
2561# ifdef SOLVE3D
2562 USE exchange_3d_mod
2563# endif
2564# ifdef DISTRIBUTE
2565 USE mp_exchange_mod, ONLY : mp_exchange2d
2566# ifdef SOLVE3D
2568# endif
2569# endif
2570# ifdef SOLVE3D
2571 USE set_depth_mod, ONLY : set_depth_tile
2572# endif
2573 USE u2dbc_mod, ONLY : u2dbc_tile
2574 USE v2dbc_mod, ONLY : v2dbc_tile
2575 USE zetabc_mod, ONLY : zetabc_tile
2576# ifdef SOLVE3D
2577 USE t3dbc_mod, ONLY : t3dbc_tile
2578 USE u3dbc_mod, ONLY : u3dbc_tile
2579 USE v3dbc_mod, ONLY : v3dbc_tile
2580# endif
2581!
2582! Imported variable declarations.
2583!
2584 integer, intent(in) :: ng, tile
2585 integer, intent(in) :: LBi, UBi, LBj, UBj
2586 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
2587 integer, intent(in) :: nstp, nnew, Linp, Lout
2588!
2589# ifdef ASSUMED_SHAPE
2590# ifdef MASKING
2591 real(r8), intent(in) :: rmask(LBi:,LBj:)
2592 real(r8), intent(in) :: umask(LBi:,LBj:)
2593 real(r8), intent(in) :: vmask(LBi:,LBj:)
2594# endif
2595# ifdef SOLVE3D
2596# if defined SEDIMENT && defined SED_MORPH
2597 real(r8), intent(in) :: bed_thick(LBi:,LBj:,:)
2598# endif
2599 real(r8), intent(in) :: om_v(LBi:,LBj:)
2600 real(r8), intent(in) :: on_u(LBi:,LBj:)
2601# ifdef ICESHELF
2602 real(r8), intent(in) :: zice(LBi:,LBj:)
2603# endif
2604 real(r8), intent(in) :: ad_t(LBi:,LBj:,:,:,:)
2605 real(r8), intent(in) :: ad_u(LBi:,LBj:,:,:)
2606 real(r8), intent(in) :: ad_v(LBi:,LBj:,:,:)
2607# endif
2608 real(r8), intent(in) :: ad_ubar(LBi:,LBj:,:)
2609 real(r8), intent(in) :: ad_vbar(LBi:,LBj:,:)
2610 real(r8), intent(in) :: ad_zeta(LBi:,LBj:,:)
2611# ifdef SOLVE3D
2612 real(r8), intent(inout) :: h(LBi:,LBj:)
2613 real(r8), intent(inout) :: Hz(LBi:,LBj:,:)
2614 real(r8), intent(inout) :: Zt_avg1(LBi:,LBj:)
2615 real(r8), intent(inout) :: t(LBi:,LBj:,:,:,:)
2616 real(r8), intent(inout) :: u(LBi:,LBj:,:,:)
2617 real(r8), intent(inout) :: v(LBi:,LBj:,:,:)
2618# endif
2619 real(r8), intent(inout) :: ubar(LBi:,LBj:,:)
2620 real(r8), intent(inout) :: vbar(LBi:,LBj:,:)
2621 real(r8), intent(inout) :: zeta(LBi:,LBj:,:)
2622# ifdef SOLVE3D
2623 real(r8), intent(out) :: z_r(LBi:,LBj:,:)
2624 real(r8), intent(out) :: z_w(LBi:,LBj:,0:)
2625# endif
2626# else
2627# ifdef MASKING
2628 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
2629 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
2630 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
2631# endif
2632# ifdef SOLVE3D
2633# if defined SEDIMENT && defined SED_MORPH
2634 real(r8), intent(in) :: bed_thick(LBi:UBi,LBj:UBj,3)
2635# endif
2636 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
2637 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
2638# ifdef ICESHELF
2639 real(r8), intent(in) :: zice(LBi:,LBj:)
2640# endif
2641 real(r8), intent(in) :: ad_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
2642 real(r8), intent(in) :: ad_u(LBi:UBi,LBj:UBj,N(ng),2)
2643 real(r8), intent(in) :: ad_v(LBi:UBi,LBj:UBj,N(ng),2)
2644# endif
2645 real(r8), intent(in) :: ad_ubar(LBi:UBi,LBj:UBj,:)
2646 real(r8), intent(in) :: ad_vbar(LBi:UBi,LBj:UBj,:)
2647 real(r8), intent(in) :: ad_zeta(LBi:UBi,LBj:UBj,:)
2648# ifdef SOLVE3D
2649 real(r8), intent(inout) :: h(LBi:UBi,LBj:UBj)
2650 real(r8), intent(inout) :: Hz(LBi:UBi,LBj:UBj,N(ng))
2651 real(r8), intent(inout) :: Zt_avg1(LBi:UBi,LBj:UBj)
2652 real(r8), intent(inout) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
2653 real(r8), intent(inout) :: u(LBi:UBi,LBj:UBj,N(ng),2)
2654 real(r8), intent(inout) :: v(LBi:UBi,LBj:UBj,N(ng),2)
2655# endif
2656 real(r8), intent(inout) :: ubar(LBi:UBi,LBj:UBj,:)
2657 real(r8), intent(inout) :: vbar(LBi:UBi,LBj:UBj,:)
2658 real(r8), intent(inout) :: zeta(LBi:UBi,LBj:UBj,:)
2659# ifdef SOLVE3D
2660 real(r8), intent(out) :: z_r(LBi:UBi,LBj:UBj,N(ng))
2661 real(r8), intent(out) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
2662# endif
2663# endif
2664!
2665! Local variable declarations.
2666!
2667 integer :: i, ic, j
2668# ifdef SOLVE3D
2669 integer :: itrc, k
2670# endif
2671 integer, dimension(NstateVar(ng)+1) :: StateVar
2672
2673 real(r8) :: p, scale
2674
2675# include "set_bounds.h"
2676!
2677!-----------------------------------------------------------------------
2678! Add a perturbation to nonlinear state variables: steepest descent
2679! direction times the perturbation amplitude "p".
2680!-----------------------------------------------------------------------
2681!
2682! Set state variable to perturb according to outer loop index.
2683!
2684# ifdef SOLVE3D
2685 statevar(1)=0
2686 statevar(2)=isfsur
2687 statevar(3)=isubar
2688 statevar(4)=isvbar
2689 statevar(5)=isuvel
2690 statevar(6)=isvvel
2691 DO i=1,nt(ng)
2692 statevar(6+i)=istvar(i)
2693 END DO
2694# else
2695 statevar(1)=0
2696 statevar(2)=isfsur
2697 statevar(3)=isubar
2698 statevar(4)=isvbar
2699# endif
2700!
2701! Set perturbation amplitude as a function of the inner loop.
2702!
2703 p=10.0_r8**real(-inner,r8)
2704 scale=1.0_r8/sqrt(addotproduct)
2705 IF (domain(ng)%SouthWest_Corner(tile)) THEN
2706 IF (master) WRITE (stdout,10)
2707 END IF
2708!
2709! Free-surface.
2710!
2711 IF ((statevar(outer).eq.0).or.(statevar(outer).eq.isfsur)) THEN
2712 DO j=jstrb,jendb
2713 DO i=istrb,iendb
2714 zeta(i,j,lout)=zeta(i,j,lout)+p*ad_zeta(i,j,linp)*scale
2715# ifdef MASKING
2716 zeta(i,j,lout)=zeta(i,j,lout)*rmask(i,j)
2717# endif
2718 END DO
2719 END DO
2720 IF (domain(ng)%SouthWest_Test(tile)) THEN
2721 IF (master) THEN
2722 WRITE (stdout,20) outer, inner, &
2723 & trim(vname(1,idsvar(isfsur)))
2724 END IF
2725 END IF
2726 END IF
2727!
2728! 2D u-momentum component.
2729!
2730 IF ((statevar(outer).eq.0).or.(statevar(outer).eq.isubar)) THEN
2731 DO j=jstrb,jendb
2732 DO i=istrm,iendb
2733 ubar(i,j,lout)=ubar(i,j,lout)+p*ad_ubar(i,j,linp)*scale
2734# ifdef MASKING
2735 ubar(i,j,lout)=ubar(i,j,lout)*umask(i,j)
2736# endif
2737 END DO
2738 END DO
2739 IF (domain(ng)%SouthWest_Test(tile)) THEN
2740 IF (master) THEN
2741 WRITE (stdout,20) outer, inner, &
2742 & trim(vname(1,idsvar(isubar)))
2743 END IF
2744 END IF
2745 END IF
2746!
2747! 2D v-momentum component.
2748!
2749 IF ((statevar(outer).eq.0).or.(statevar(outer).eq.isvbar)) THEN
2750 DO j=jstrm,jendb
2751 DO i=istrb,iendb
2752 vbar(i,j,lout)=vbar(i,j,lout)+p*ad_vbar(i,j,linp)*scale
2753# ifdef MASKING
2754 vbar(i,j,lout)=vbar(i,j,lout)*vmask(i,j)
2755# endif
2756 END DO
2757 END DO
2758 IF (domain(ng)%SouthWest_Test(tile)) THEN
2759 IF (master) THEN
2760 WRITE (stdout,20) outer, inner, &
2761 & trim(vname(1,idsvar(isvbar)))
2762 END IF
2763 END IF
2764 END IF
2765# ifdef SOLVE3D
2766!
2767! 3D u-momentum component.
2768!
2769 IF ((statevar(outer).eq.0).or.(statevar(outer).eq.isuvel)) THEN
2770 DO k=1,n(ng)
2771 DO j=jstrb,jendb
2772 DO i=istrm,iendb
2773 u(i,j,k,lout)=u(i,j,k,lout)+p*ad_u(i,j,k,linp)*scale
2774# ifdef MASKING
2775 u(i,j,k,lout)=u(i,j,k,lout)*umask(i,j)
2776# endif
2777 END DO
2778 END DO
2779 END DO
2780 IF (domain(ng)%SouthWest_Test(tile)) THEN
2781 IF (master) THEN
2782 WRITE (stdout,20) outer, inner, &
2783 & trim(vname(1,idsvar(isuvel)))
2784 END IF
2785 END IF
2786 END IF
2787!
2788! 3D v-momentum component.
2789!
2790 IF ((statevar(outer).eq.0).or.(statevar(outer).eq.isvvel)) THEN
2791 DO k=1,n(ng)
2792 DO j=jstrm,jendb
2793 DO i=istrb,iendb
2794 v(i,j,k,lout)=v(i,j,k,lout)+p*ad_v(i,j,k,linp)*scale
2795# ifdef MASKING
2796 v(i,j,k,lout)=v(i,j,k,lout)*vmask(i,j)
2797# endif
2798 END DO
2799 END DO
2800 END DO
2801 IF (domain(ng)%SouthWest_Test(tile)) THEN
2802 IF (master) THEN
2803 WRITE (stdout,20) outer, inner, &
2804 & trim(vname(1,idsvar(isvvel)))
2805 END IF
2806 END IF
2807 END IF
2808!
2809! Tracers.
2810!
2811 DO itrc=1,nt(ng)
2812 IF ((statevar(outer).eq.0).or. &
2813 & (statevar(outer).eq.istvar(itrc))) THEN
2814 DO k=1,n(ng)
2815 DO j=jstrb,jendb
2816 DO i=istrb,iendb
2817 t(i,j,k,lout,itrc)=t(i,j,k,lout,itrc)+ &
2818 & p*ad_t(i,j,k,linp,itrc)*scale
2819# ifdef MASKING
2820 t(i,j,k,lout,itrc)=t(i,j,k,lout,itrc)*rmask(i,j)
2821# endif
2822 END DO
2823 END DO
2824 END DO
2825 IF (domain(ng)%SouthWest_Test(tile)) THEN
2826 IF (master) THEN
2827 WRITE (stdout,20) outer, inner, &
2828 & trim(vname(1,idsvar(istvar(itrc))))
2829 END IF
2830 END IF
2831 END IF
2832 END DO
2833# endif
2834 IF (master) WRITE (stdout,'(/)')
2835!
2836!-----------------------------------------------------------------------
2837! Apply lateral boundary conditions to 2D fields.
2838!-----------------------------------------------------------------------
2839!
2840 CALL zetabc_tile (ng, tile, &
2841 & lbi, ubi, lbj, ubj, &
2842 & imins, imaxs, jmins, jmaxs, &
2843 & lout, lout, lout, &
2844 & zeta)
2845# ifndef SOLVE3D
2846 CALL u2dbc_tile (ng, tile, &
2847 & lbi, ubi, lbj, ubj, &
2848 & imins, imaxs, jmins, jmaxs, &
2849 & lout, lout, lout, &
2850 & ubar, vbar, zeta)
2851 CALL v2dbc_tile (ng, tile, &
2852 & lbi, ubi, lbj, ubj, &
2853 & imins, imaxs, jmins, jmaxs, &
2854 & lout, lout, lout, &
2855 & ubar, vbar, zeta)
2856# endif
2857!
2858 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
2859 CALL exchange_r2d_tile (ng, tile, &
2860 & lbi, ubi, lbj, ubj, &
2861 & zeta(:,:,lout))
2862# ifndef SOLVE3D
2863 CALL exchange_u2d_tile (ng, tile, &
2864 & lbi, ubi, lbj, ubj, &
2865 & ubar(:,:,lout))
2866 CALL exchange_v2d_tile (ng, tile, &
2867 & lbi, ubi, lbj, ubj, &
2868 & vbar(:,:,lout))
2869# endif
2870 END IF
2871
2872# ifdef DISTRIBUTE
2873!
2874 CALL mp_exchange2d (ng, tile, inlm, 1, &
2875 & lbi, ubi, lbj, ubj, &
2876 & nghostpoints, &
2877 & ewperiodic(ng), nsperiodic(ng), &
2878 & zeta(:,:,lout))
2879# ifndef SOLVE3D
2880 CALL mp_exchange2d (ng, tile, inlm, 2, &
2881 & lbi, ubi, lbj, ubj, &
2882 & nghostpoints, &
2883 & ewperiodic(ng), nsperiodic(ng), &
2884 & ubar(:,:,lout), &
2885 & vbar(:,:,lout))
2886# endif
2887# endif
2888# ifdef SOLVE3D
2889!
2890!-----------------------------------------------------------------------
2891! Compute new depths and thicknesses.
2892!-----------------------------------------------------------------------
2893!
2894 CALL set_depth_tile (ng, tile, inlm, &
2895 & lbi, ubi, lbj, ubj, &
2896 & imins, imaxs, jmins, jmaxs, &
2897 & nstp, nnew, &
2898 & h, &
2899# ifdef ICESHELF
2900 & zice, &
2901# endif
2902# if defined SEDIMENT && defined SED_MORPH
2903 & bed_thick, &
2904# endif
2905 & zt_avg1, &
2906 & hz, z_r, z_w)
2907!
2908!-----------------------------------------------------------------------
2909! Apply lateral boundary conditions to perturbed 3D fields.
2910!-----------------------------------------------------------------------
2911!
2912 CALL u3dbc_tile (ng, tile, &
2913 & lbi, ubi, lbj, ubj, n(ng), &
2914 & imins, imaxs, jmins, jmaxs, &
2915 & lout, lout, u)
2916 CALL v3dbc_tile (ng, tile, &
2917 & lbi, ubi, lbj, ubj, n(ng), &
2918 & imins, imaxs, jmins, jmaxs, &
2919 & lout, lout, v)
2920!
2921 ic=0
2922 DO itrc=1,nt(ng)
2923 IF (ltracerclm(itrc,ng).and.lnudgetclm(itrc,ng)) THEN
2924 ic=ic+1
2925 END IF
2926 CALL t3dbc_tile (ng, tile, itrc, ic, &
2927 & lbi, ubi, lbj, ubj, n(ng), nt(ng), &
2928 & imins, imaxs, jmins, jmaxs, &
2929 & lout, lout, t)
2930 END DO
2931!
2932 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
2933 CALL exchange_u3d_tile (ng, tile, &
2934 & lbi, ubi, lbj, ubj, 1, n(ng), &
2935 & u(:,:,:,lout))
2936 CALL exchange_v3d_tile (ng, tile, &
2937 & lbi, ubi, lbj, ubj, 1, n(ng), &
2938 & v(:,:,:,lout))
2939 DO itrc=1,nt(ng)
2940 CALL exchange_r3d_tile (ng, tile, &
2941 & lbi, ubi, lbj, ubj, 1, n(ng), &
2942 & t(:,:,:,lout,itrc))
2943 END DO
2944 END IF
2945
2946# ifdef DISTRIBUTE
2947!
2948 CALL mp_exchange3d (ng, tile, inlm, 2, &
2949 & lbi, ubi, lbj, ubj, 1, n(ng), &
2950 & nghostpoints, &
2951 & ewperiodic(ng), nsperiodic(ng), &
2952 & u(:,:,:,lout), &
2953 & v(:,:,:,lout))
2954 CALL mp_exchange4d (ng, tile, inlm, 1, &
2955 & lbi, ubi, lbj, ubj, 1, n(ng), 1, nt(ng), &
2956 & nghostpoints, &
2957 & ewperiodic(ng), nsperiodic(ng), &
2958 & t(:,:,:,lout,:))
2959# endif
2960# endif
2961!
2962 10 FORMAT (/,' Perturbing Nonlinear model Initial Conditions:',/)
2963 20 FORMAT (' (Outer,Inner) = ','(',i4.4,',',i4.4,')',3x, &
2964 & 'Perturbing State Variable: ',a)
2965!
2966 RETURN
subroutine exchange_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_u2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_v2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_u3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
subroutine exchange_r3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
subroutine exchange_v3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
real(r8) addotproduct
integer stdout
integer isvvel
integer isvbar
integer, dimension(:), allocatable istvar
integer isuvel
integer isfsur
character(len=maxlen), dimension(6, 0:nv) vname
integer, dimension(:), allocatable idsvar
integer isubar
logical master
integer nghostpoints
Definition mod_param.F:710
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
logical, dimension(:,:), allocatable ltracerclm
logical, dimension(:,:), allocatable lnudgetclm
integer inner
integer outer
subroutine mp_exchange4d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, lbt, ubt, nghost, ew_periodic, ns_periodic, a, b, c)
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine mp_exchange3d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine, public set_depth_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nstp, nnew, h, zice, zt_avg1, hz, z_r, z_w)
Definition set_depth.F:86
subroutine, public t3dbc_tile(ng, tile, itrc, ic, lbi, ubi, lbj, ubj, ubk, ubt, imins, imaxs, jmins, jmaxs, nstp, nout, t)
Definition t3dbc_im.F:55
subroutine, public u2dbc_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, krhs, kstp, kout, ubar, vbar, zeta)
Definition u2dbc_im.F:56
subroutine, public u3dbc_tile(ng, tile, lbi, ubi, lbj, ubj, ubk, imins, imaxs, jmins, jmaxs, nstp, nout, u)
Definition u3dbc_im.F:55
subroutine, public v2dbc_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, krhs, kstp, kout, ubar, vbar, zeta)
Definition v2dbc_im.F:57
subroutine, public v3dbc_tile(ng, tile, lbi, ubi, lbj, ubj, ubk, imins, imaxs, jmins, jmaxs, nstp, nout, v)
Definition v3dbc_im.F:55
subroutine, public zetabc_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, krhs, kstp, kout, zeta)
Definition zetabc.F:65

References mod_fourdvar::addotproduct, mod_param::domain, mod_scalars::ewperiodic, exchange_2d_mod::exchange_r2d_tile(), exchange_3d_mod::exchange_r3d_tile(), exchange_2d_mod::exchange_u2d_tile(), exchange_3d_mod::exchange_u3d_tile(), exchange_2d_mod::exchange_v2d_tile(), exchange_3d_mod::exchange_v3d_tile(), mod_ncparam::idsvar, mod_param::inlm, mod_scalars::inner, mod_ncparam::isfsur, mod_ncparam::istvar, mod_ncparam::isubar, mod_ncparam::isuvel, mod_ncparam::isvbar, mod_ncparam::isvvel, mod_scalars::lnudgetclm, mod_scalars::ltracerclm, mod_parallel::master, mp_exchange_mod::mp_exchange2d(), mp_exchange_mod::mp_exchange3d(), mp_exchange_mod::mp_exchange4d(), mod_param::nghostpoints, mod_scalars::nsperiodic, mod_scalars::outer, set_depth_mod::set_depth_tile(), mod_iounits::stdout, t3dbc_mod::t3dbc_tile(), u2dbc_mod::u2dbc_tile(), u3dbc_mod::u3dbc_tile(), v2dbc_mod::v2dbc_tile(), v3dbc_mod::v3dbc_tile(), mod_ncparam::vname, and zetabc_mod::zetabc_tile().

Referenced by ini_perturb().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ load_adtotl()

subroutine, public ini_adjust_mod::load_adtotl ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) linp,
integer, intent(in) lout,
logical, intent(in) add )

Definition at line 813 of file ini_adjust.F.

814!
815!=======================================================================
816! !
817! This routine loads or adds Linp adjoint state variables into Lout !
818! Lout tangent linear state variables. !
819! !
820! On Input: !
821! !
822! ng Nested grid number. !
823! tile Domain partition. !
824! Linp Tangent linear state time index to add. !
825! Lout Nonlinear state time index to update. !
826! add Logical switch to add to imported values. !
827! !
828!=======================================================================
829!
830 USE mod_param
831# ifdef ADJUST_BOUNDARY
832 USE mod_boundary
833# endif
834# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
835 USE mod_forces
836# endif
837 USE mod_grid
838 USE mod_ocean
839!
840! Imported variable declarations.
841!
842 logical, intent(in) :: add
843 integer, intent(in) :: ng, tile, Linp, Lout
844!
845! Local variable declarations.
846!
847 character (len=*), parameter :: MyFile = &
848 & __FILE__//", load_ADtoTL"
849
850# include "tile.h"
851!
852# ifdef PROFILE
853 CALL wclock_on (ng, itlm, 7, __line__, myfile)
854# endif
855 CALL load_adtotl_tile (ng, tile, &
856 & lbi, ubi, lbj, ubj, lbij, ubij, &
857 & imins, imaxs, jmins, jmaxs, &
858 & linp, lout, add, &
859# ifdef MASKING
860 & grid(ng) % rmask, &
861 & grid(ng) % umask, &
862 & grid(ng) % vmask, &
863# endif
864# ifdef ADJUST_BOUNDARY
865# ifdef SOLVE3D
866 & boundary(ng) % ad_t_obc, &
867 & boundary(ng) % ad_u_obc, &
868 & boundary(ng) % ad_v_obc, &
869# endif
870 & boundary(ng) % ad_ubar_obc, &
871 & boundary(ng) % ad_vbar_obc, &
872 & boundary(ng) % ad_zeta_obc, &
873# endif
874# ifdef ADJUST_WSTRESS
875 & forces(ng) % ad_ustr, &
876 & forces(ng) % ad_vstr, &
877# endif
878# ifdef SOLVE3D
879# ifdef ADJUST_STFLUX
880 & forces(ng) % ad_tflux, &
881# endif
882 & ocean(ng) % ad_t, &
883 & ocean(ng) % ad_u, &
884 & ocean(ng) % ad_v, &
885# if defined WEAK_CONSTRAINT && defined TIME_CONV
886 & ocean(ng) % ad_ubar, &
887 & ocean(ng) % ad_vbar, &
888# endif
889# else
890 & ocean(ng) % ad_ubar, &
891 & ocean(ng) % ad_vbar, &
892# endif
893 & ocean(ng) % ad_zeta, &
894# ifdef ADJUST_BOUNDARY
895# ifdef SOLVE3D
896 & boundary(ng) % tl_t_obc, &
897 & boundary(ng) % tl_u_obc, &
898 & boundary(ng) % tl_v_obc, &
899# endif
900 & boundary(ng) % tl_ubar_obc, &
901 & boundary(ng) % tl_vbar_obc, &
902 & boundary(ng) % tl_zeta_obc, &
903# endif
904# ifdef ADJUST_WSTRESS
905 & forces(ng) % tl_ustr, &
906 & forces(ng) % tl_vstr, &
907# endif
908# ifdef SOLVE3D
909# ifdef ADJUST_STFLUX
910 & forces(ng) % tl_tflux, &
911# endif
912 & ocean(ng) % tl_t, &
913 & ocean(ng) % tl_u, &
914 & ocean(ng) % tl_v, &
915# if defined WEAK_CONSTRAINT && defined TIME_CONV
916 & ocean(ng) % tl_ubar, &
917 & ocean(ng) % tl_vbar, &
918# endif
919# else
920 & ocean(ng) % tl_ubar, &
921 & ocean(ng) % tl_vbar, &
922# endif
923 & ocean(ng) % tl_zeta)
924# ifdef PROFILE
925 CALL wclock_off (ng, itlm, 7, __line__, myfile)
926# endif
927!
928 RETURN
type(t_boundary), dimension(:), allocatable boundary
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
integer, parameter itlm
Definition mod_param.F:663

References mod_boundary::boundary, mod_forces::forces, mod_grid::grid, mod_param::itlm, load_adtotl_tile(), mod_ocean::ocean, wclock_off(), and wclock_on().

Referenced by convolve_mod::convolve(), convolve_mod::error_covariance(), roms_kernel_mod::roms_run(), and convolve_mod::saddlec().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ load_adtotl_tile()

subroutine ini_adjust_mod::load_adtotl_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) lbij,
integer, intent(in) ubij,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) linp,
integer, intent(in) lout,
logical, intent(in) add,
real(r8), dimension(lbi:,lbj:), intent(in) rmask,
real(r8), dimension(lbi:,lbj:), intent(in) umask,
real(r8), dimension(lbi:,lbj:), intent(in) vmask,
real(r8), dimension(lbij:,:,:,:,:,:), intent(inout) ad_t_obc,
real(r8), dimension(lbij:,:,:,:,:), intent(inout) ad_u_obc,
real(r8), dimension(lbij:,:,:,:,:), intent(inout) ad_v_obc,
real(r8), dimension(lbij:,:,:,:), intent(inout) ad_ubar_obc,
real(r8), dimension(lbij:,:,:,:), intent(inout) ad_vbar_obc,
real(r8), dimension(lbij:,:,:,:), intent(inout) ad_zeta_obc,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) ad_ustr,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) ad_vstr,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) ad_tflux,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) ad_t,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) ad_u,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) ad_v,
real(r8), dimension(lbi:,lbj:,:), intent(inout) ad_ubar,
real(r8), dimension(lbi:,lbj:,:), intent(inout) ad_vbar,
real(r8), dimension(lbi:,lbj:,:), intent(inout) ad_zeta,
real(r8), dimension(lbij:,:,:,:,:,:), intent(inout) tl_t_obc,
real(r8), dimension(lbij:,:,:,:,:), intent(inout) tl_u_obc,
real(r8), dimension(lbij:,:,:,:,:), intent(inout) tl_v_obc,
real(r8), dimension(lbij:,:,:,:), intent(inout) tl_ubar_obc,
real(r8), dimension(lbij:,:,:,:), intent(inout) tl_vbar_obc,
real(r8), dimension(lbij:,:,:,:), intent(inout) tl_zeta_obc,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) tl_ustr,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) tl_vstr,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) tl_tflux,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) tl_t,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) tl_u,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) tl_v,
real(r8), dimension(lbi:,lbj:,:), intent(inout) tl_ubar,
real(r8), dimension(lbi:,lbj:,:), intent(inout) tl_vbar,
real(r8), dimension(lbi:,lbj:,:), intent(inout) tl_zeta )
private

Definition at line 932 of file ini_adjust.F.

983!***********************************************************************
984!
985 USE mod_param
986 USE mod_ncparam
987# if defined ADJUST_BOUNDARY || defined ADJUST_STFLUX || \
988 defined adjust_wstress
989 USE mod_scalars
990# endif
991!
993 USE state_copy_mod, ONLY : state_copy
994!
995! Imported variable declarations.
996!
997 logical, intent(in) :: add
998 integer, intent(in) :: ng, tile
999 integer, intent(in) :: LBi, UBi, LBj, UBj, LBij, UBij
1000 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
1001 integer, intent(in) :: Linp, Lout
1002!
1003# ifdef ASSUMED_SHAPE
1004# ifdef MASKING
1005 real(r8), intent(in) :: rmask(LBi:,LBj:)
1006 real(r8), intent(in) :: umask(LBi:,LBj:)
1007 real(r8), intent(in) :: vmask(LBi:,LBj:)
1008# endif
1009# ifdef ADJUST_BOUNDARY
1010# ifdef SOLVE3D
1011 real(r8), intent(inout) :: ad_t_obc(LBij:,:,:,:,:,:)
1012 real(r8), intent(inout) :: ad_u_obc(LBij:,:,:,:,:)
1013 real(r8), intent(inout) :: ad_v_obc(LBij:,:,:,:,:)
1014# endif
1015 real(r8), intent(inout) :: ad_ubar_obc(LBij:,:,:,:)
1016 real(r8), intent(inout) :: ad_vbar_obc(LBij:,:,:,:)
1017 real(r8), intent(inout) :: ad_zeta_obc(LBij:,:,:,:)
1018# endif
1019# ifdef ADJUST_WSTRESS
1020 real(r8), intent(inout) :: ad_ustr(LBi:,LBj:,:,:)
1021 real(r8), intent(inout) :: ad_vstr(LBi:,LBj:,:,:)
1022# endif
1023# ifdef SOLVE3D
1024# ifdef ADJUST_STFLUX
1025 real(r8), intent(inout) :: ad_tflux(LBi:,LBj:,:,:,:)
1026# endif
1027 real(r8), intent(inout) :: ad_t(LBi:,LBj:,:,:,:)
1028 real(r8), intent(inout) :: ad_u(LBi:,LBj:,:,:)
1029 real(r8), intent(inout) :: ad_v(LBi:,LBj:,:,:)
1030# if defined WEAK_CONSTRAINT && defined TIME_CONV
1031 real(r8), intent(inout) :: ad_ubar(LBi:,LBj:,:)
1032 real(r8), intent(inout) :: ad_vbar(LBi:,LBj:,:)
1033# endif
1034# else
1035 real(r8), intent(inout) :: ad_ubar(LBi:,LBj:,:)
1036 real(r8), intent(inout) :: ad_vbar(LBi:,LBj:,:)
1037# endif
1038 real(r8), intent(inout) :: ad_zeta(LBi:,LBj:,:)
1039# ifdef ADJUST_BOUNDARY
1040# ifdef SOLVE3D
1041 real(r8), intent(inout) :: tl_t_obc(LBij:,:,:,:,:,:)
1042 real(r8), intent(inout) :: tl_u_obc(LBij:,:,:,:,:)
1043 real(r8), intent(inout) :: tl_v_obc(LBij:,:,:,:,:)
1044# endif
1045 real(r8), intent(inout) :: tl_ubar_obc(LBij:,:,:,:)
1046 real(r8), intent(inout) :: tl_vbar_obc(LBij:,:,:,:)
1047 real(r8), intent(inout) :: tl_zeta_obc(LBij:,:,:,:)
1048# endif
1049# ifdef ADJUST_WSTRESS
1050 real(r8), intent(inout) :: tl_ustr(LBi:,LBj:,:,:)
1051 real(r8), intent(inout) :: tl_vstr(LBi:,LBj:,:,:)
1052# endif
1053# ifdef SOLVE3D
1054# ifdef ADJUST_STFLUX
1055 real(r8), intent(inout) :: tl_tflux(LBi:,LBj:,:,:,:)
1056# endif
1057 real(r8), intent(inout) :: tl_t(LBi:,LBj:,:,:,:)
1058 real(r8), intent(inout) :: tl_u(LBi:,LBj:,:,:)
1059 real(r8), intent(inout) :: tl_v(LBi:,LBj:,:,:)
1060# if defined WEAK_CONSTRAINT && defined TIME_CONV
1061 real(r8), intent(inout) :: tl_ubar(LBi:,LBj:,:)
1062 real(r8), intent(inout) :: tl_vbar(LBi:,LBj:,:)
1063# endif
1064# else
1065 real(r8), intent(inout) :: tl_ubar(LBi:,LBj:,:)
1066 real(r8), intent(inout) :: tl_vbar(LBi:,LBj:,:)
1067# endif
1068 real(r8), intent(inout) :: tl_zeta(LBi:,LBj:,:)
1069# else
1070# ifdef MASKING
1071 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
1072 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
1073 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
1074# endif
1075# ifdef ADJUST_BOUNDARY
1076# ifdef SOLVE3D
1077 real(r8), intent(inout) :: ad_t_obc(LBij:UBij,N(ng),4, &
1078 & Nbrec(ng),2,NT(ng))
1079 real(r8), intent(inout) :: ad_u_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
1080 real(r8), intent(inout) :: ad_v_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
1081# endif
1082 real(r8), intent(inout) :: ad_ubar_obc(LBij:UBij,4,Nbrec(ng),2)
1083 real(r8), intent(inout) :: ad_vbar_obc(LBij:UBij,4,Nbrec(ng),2)
1084 real(r8), intent(input) :: ad_zeta_obc(LBij:UBij,4,Nbrec(ng),2)
1085# endif
1086# ifdef ADJUST_WSTRESS
1087 real(r8), intent(inout) :: ad_ustr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
1088 real(r8), intent(inout) :: ad_vstr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
1089# endif
1090# ifdef SOLVE3D
1091# ifdef ADJUST_STFLUX
1092 real(r8), intent(inout) :: ad_tflux(LBi:UBi,LBj:UBj, &
1093 & Nfrec(ng),2,NT(ng))
1094# endif
1095 real(r8), intent(inout) :: ad_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
1096 real(r8), intent(inout) :: ad_u(LBi:UBi,LBj:UBj,N(ng),2)
1097 real(r8), intent(inout) :: ad_v(LBi:UBi,LBj:UBj,N(ng),2)
1098# if defined WEAK_CONSTRAINT && defined TIME_CONV
1099 real(r8), intent(inout) :: ad_ubar(LBi:UBi,LBj:UBj,:)
1100 real(r8), intent(inout) :: ad_vbar(LBi:UBi,LBj:UBj,:)
1101# endif
1102# else
1103 real(r8), intent(inout) :: ad_ubar(LBi:UBi,LBj:UBj,:)
1104 real(r8), intent(inout) :: ad_vbar(LBi:UBi,LBj:UBj,:)
1105# endif
1106 real(r8), intent(inout) :: ad_zeta(LBi:UBi,LBj:UBj,:)
1107# ifdef ADJUST_BOUNDARY
1108# ifdef SOLVE3D
1109 real(r8), intent(inout) :: tl_t_obc(LBij:UBij,N(ng),4, &
1110 & Nbrec(ng),2,NT(ng))
1111 real(r8), intent(inout) :: tl_u_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
1112 real(r8), intent(inout) :: tl_v_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
1113# endif
1114 real(r8), intent(inout) :: tl_ubar_obc(LBij:UBij,4,Nbrec(ng),2)
1115 real(r8), intent(inout) :: tl_vbar_obc(LBij:UBij,4,Nbrec(ng),2)
1116 real(r8), intent(inout) :: tl_zeta_obc(LBij:UBij,4,Nbrec(ng),2)
1117# endif
1118# ifdef ADJUST_WSTRESS
1119 real(r8), intent(inout) :: tl_ustr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
1120 real(r8), intent(inout) :: tl_vstr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
1121# endif
1122# ifdef SOLVE3D
1123# ifdef ADJUST_STFLUX
1124 real(r8), intent(inout) :: tl_tflux(LBi:UBi,LBj:UBj, &
1125 & Nfrec(ng),2,NT(ng))
1126# endif
1127 real(r8), intent(inout) :: tl_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
1128 real(r8), intent(inout) :: tl_u(LBi:UBi,LBj:UBj,N(ng),2)
1129 real(r8), intent(inout) :: tl_v(LBi:UBi,LBj:UBj,N(ng),2)
1130# if defined WEAK_CONSTRAINT && defined TIME_CONV
1131 real(r8), intent(inout) :: tl_ubar(LBi:UBi,LBj:UBj,:)
1132 real(r8), intent(inout) :: tl_vbar(LBi:UBi,LBj:UBj,:)
1133# endif
1134# else
1135 real(r8), intent(inout) :: tl_ubar(LBi:UBi,LBj:UBj,:)
1136 real(r8), intent(inout) :: tl_vbar(LBi:UBi,LBj:UBj,:)
1137# endif
1138 real(r8), intent(inout) :: tl_zeta(LBi:UBi,LBj:UBj,:)
1139# endif
1140!
1141! Local variable declarations.
1142!
1143# ifdef MASKING
1144 integer :: i, j, k
1145 integer :: ib, ir, it
1146# endif
1147 real(r8) :: fac1, fac2
1148
1149# ifdef MASKING
1150# include "set_bounds.h"
1151!
1152!-----------------------------------------------------------------------
1153! Multiply by land-sea mask - fail safe. Notice that if add=.FALSE. the
1154! "state_addition" routine is not called and the state arrays are not
1155! masked.
1156!-----------------------------------------------------------------------
1157!
1158! Free-surface.
1159!
1160 DO j=jstrt,jendt
1161 DO i=istrt,iendt
1162 ad_zeta(i,j,linp)=ad_zeta(i,j,linp)*rmask(i,j)
1163 END DO
1164 END DO
1165
1166# ifdef ADJUST_BOUNDARY
1167!
1168! Free-surface open boundaries.
1169!
1170 IF (any(lobc(:,isfsur,ng))) THEN
1171 DO ir=1,nbrec(ng)
1172 IF ((lobc(iwest,isfsur,ng)).and. &
1173 & domain(ng)%Western_Edge(tile)) THEN
1174 ib=iwest
1175 DO j=jstr,jend
1176 ad_zeta_obc(j,ib,ir,linp)=ad_zeta_obc(j,ib,ir,linp)* &
1177 & rmask(istr-1,j)
1178 END DO
1179 END IF
1180 IF ((lobc(ieast,isfsur,ng)).and. &
1181 & domain(ng)%Eastern_Edge(tile)) THEN
1182 ib=ieast
1183 DO j=jstr,jend
1184 ad_zeta_obc(j,ib,ir,linp)=ad_zeta_obc(j,ib,ir,linp)* &
1185 & rmask(iend+1,j)
1186 END DO
1187 END IF
1188 IF ((lobc(isouth,isfsur,ng)).and. &
1189 & domain(ng)%Southern_Edge(tile)) THEN
1190 ib=isouth
1191 DO i=istr,iend
1192 ad_zeta_obc(i,ib,ir,linp)=ad_zeta_obc(i,ib,ir,linp)* &
1193 & rmask(i,jstr-1)
1194 END DO
1195 END IF
1196 IF ((lobc(inorth,isfsur,ng)).and. &
1197 & domain(ng)%Northern_Edge(tile)) THEN
1198 ib=inorth
1199 DO i=istr,iend
1200 ad_zeta_obc(i,ib,ir,linp)=ad_zeta_obc(i,ib,ir,linp)* &
1201 & rmask(i,jend+1)
1202 END DO
1203 END IF
1204 END DO
1205 END IF
1206# endif
1207
1208# if !defined SOLVE3D || \
1209 (defined weak_constraint && defined time_conv)
1210!
1211! 2D U-momentum component.
1212!
1213 DO j=jstrt,jendt
1214 DO i=istrp,iendt
1215 ad_ubar(i,j,linp)=ad_ubar(i,j,linp)*umask(i,j)
1216 END DO
1217 END DO
1218!
1219! 2D V-momentum component.
1220!
1221 DO j=jstrp,jendt
1222 DO i=istrt,iendt
1223 ad_vbar(i,j,linp)=ad_vbar(i,j,linp)*vmask(i,j)
1224 END DO
1225 END DO
1226# endif
1227
1228# ifdef ADJUST_BOUNDARY
1229!
1230! 2D U-momentum open boundaries.
1231!
1232 IF (any(lobc(:,isubar,ng))) THEN
1233 DO ir=1,nbrec(ng)
1234 IF ((lobc(iwest,isubar,ng)).and. &
1235 & domain(ng)%Western_Edge(tile)) THEN
1236 ib=iwest
1237 DO j=jstr,jend
1238 ad_ubar_obc(j,ib,ir,linp)=ad_ubar_obc(j,ib,ir,linp)* &
1239 & umask(istr,j)
1240 END DO
1241 END IF
1242 IF ((lobc(ieast,isubar,ng)).and. &
1243 & domain(ng)%Eastern_Edge(tile)) THEN
1244 ib=ieast
1245 DO j=jstr,jend
1246 ad_ubar_obc(j,ib,ir,linp)=ad_ubar_obc(j,ib,ir,linp)* &
1247 & umask(iend+1,j)
1248 END DO
1249 END IF
1250 IF ((lobc(isouth,isubar,ng)).and. &
1251 & domain(ng)%Southern_Edge(tile)) THEN
1252 ib=isouth
1253 DO i=istru,iend
1254 ad_ubar_obc(i,ib,ir,linp)=ad_ubar_obc(i,ib,ir,linp)* &
1255 & umask(i,jstr-1)
1256 END DO
1257 END IF
1258 IF ((lobc(inorth,isubar,ng)).and. &
1259 & domain(ng)%Northern_Edge(tile)) THEN
1260 ib=inorth
1261 DO i=istru,iend
1262 ad_ubar_obc(i,ib,ir,linp)=ad_ubar_obc(i,ib,ir,linp)* &
1263 & umask(i,jend+1)
1264 END DO
1265 END IF
1266 END DO
1267 END IF
1268!
1269! 2D V-momentum open boundaries.
1270!
1271 IF (any(lobc(:,isvbar,ng))) THEN
1272 DO ir=1,nbrec(ng)
1273 IF ((lobc(iwest,isvbar,ng)).and. &
1274 & domain(ng)%Western_Edge(tile)) THEN
1275 ib=iwest
1276 DO j=jstrv,jend
1277 ad_vbar_obc(j,ib,ir,linp)=ad_vbar_obc(j,ib,ir,linp)* &
1278 & vmask(istr-1,j)
1279 END DO
1280 END IF
1281 IF ((lobc(ieast,isvbar,ng)).and. &
1282 & domain(ng)%Eastern_Edge(tile)) THEN
1283 ib=ieast
1284 DO j=jstrv,jend
1285 ad_vbar_obc(j,ib,ir,linp)=ad_vbar_obc(j,ib,ir,linp)* &
1286 & vmask(iend+1,j)
1287 END DO
1288 END IF
1289 IF ((lobc(isouth,isvbar,ng)).and. &
1290 & domain(ng)%Southern_Edge(tile)) THEN
1291 ib=isouth
1292 DO i=istr,iend
1293 ad_vbar_obc(i,ib,ir,linp)=ad_vbar_obc(i,ib,ir,linp)* &
1294 & vmask(i,jstr)
1295 END DO
1296 END IF
1297 IF ((lobc(inorth,isvbar,ng)).and. &
1298 & domain(ng)%Northern_Edge(tile)) THEN
1299 ib=inorth
1300 DO i=istr,iend
1301 ad_vbar_obc(i,ib,ir,linp)=ad_vbar_obc(i,ib,ir,linp)* &
1302 & vmask(i,jend+1)
1303 END DO
1304 END IF
1305 END DO
1306 END IF
1307# endif
1308
1309# ifdef ADJUST_WSTRESS
1310!
1311! Surface momentum stress.
1312!
1313 DO ir=1,nfrec(ng)
1314 DO j=jstrt,jendt
1315 DO i=istrp,iendt
1316 ad_ustr(i,j,ir,linp)=ad_ustr(i,j,ir,linp)*umask(i,j)
1317 END DO
1318 END DO
1319 DO j=jstrp,jendt
1320 DO i=istrt,iendt
1321 ad_vstr(i,j,ir,linp)=ad_vstr(i,j,ir,linp)*vmask(i,j)
1322 END DO
1323 END DO
1324 END DO
1325# endif
1326
1327# ifdef SOLVE3D
1328!
1329! 3D U-momentum component.
1330!
1331 DO k=1,n(ng)
1332 DO j=jstrt,jendt
1333 DO i=istrp,iendt
1334 ad_u(i,j,k,linp)=ad_u(i,j,k,linp)*umask(i,j)
1335 END DO
1336 END DO
1337 END DO
1338
1339# ifdef ADJUST_BOUNDARY
1340!
1341! 3D U-momentum open boundaries.
1342!
1343 IF (any(lobc(:,isuvel,ng))) THEN
1344 DO ir=1,nbrec(ng)
1345 IF ((lobc(iwest,isuvel,ng)).and. &
1346 & domain(ng)%Western_Edge(tile)) THEN
1347 ib=iwest
1348 DO k=1,n(ng)
1349 DO j=jstr,jend
1350 ad_u_obc(j,k,ib,ir,linp)=ad_u_obc(j,k,ib,ir,linp)* &
1351 & umask(istr,j)
1352 END DO
1353 END DO
1354 END IF
1355 IF ((lobc(ieast,isuvel,ng)).and. &
1356 & domain(ng)%Eastern_Edge(tile)) THEN
1357 ib=ieast
1358 DO k=1,n(ng)
1359 DO j=jstr,jend
1360 ad_u_obc(j,k,ib,ir,linp)=ad_u_obc(j,k,ib,ir,linp)* &
1361 & umask(iend+1,j)
1362 END DO
1363 END DO
1364 END IF
1365 IF ((lobc(isouth,isuvel,ng)).and. &
1366 & domain(ng)%Southern_Edge(tile)) THEN
1367 ib=isouth
1368 DO k=1,n(ng)
1369 DO i=istru,iend
1370 ad_u_obc(i,k,ib,ir,linp)=ad_u_obc(i,k,ib,ir,linp)* &
1371 & umask(i,jstr-1)
1372 END DO
1373 END DO
1374 END IF
1375 IF ((lobc(inorth,isuvel,ng)).and. &
1376 & domain(ng)%Northern_Edge(tile)) THEN
1377 ib=inorth
1378 DO k=1,n(ng)
1379 DO i=istru,iend
1380 ad_u_obc(i,k,ib,ir,linp)=ad_u_obc(i,k,ib,ir,linp)* &
1381 & umask(i,jend+1)
1382 END DO
1383 END DO
1384 END IF
1385 END DO
1386 END IF
1387# endif
1388!
1389! 3D V-momentum component.
1390!
1391 DO k=1,n(ng)
1392 DO j=jstrp,jendt
1393 DO i=istrt,iendt
1394 ad_v(i,j,k,linp)=ad_v(i,j,k,linp)*vmask(i,j)
1395 END DO
1396 END DO
1397 END DO
1398
1399# ifdef ADJUST_BOUNDARY
1400!
1401! 3D V-momentum open boundaries.
1402!
1403 IF (any(lobc(:,isvvel,ng))) THEN
1404 DO ir=1,nbrec(ng)
1405 IF ((lobc(iwest,isvvel,ng)).and. &
1406 & domain(ng)%Western_Edge(tile)) THEN
1407 ib=iwest
1408 DO k=1,n(ng)
1409 DO j=jstrv,jend
1410 ad_v_obc(j,k,ib,ir,linp)=ad_v_obc(j,k,ib,ir,linp)* &
1411 & vmask(istr-1,j)
1412 END DO
1413 END DO
1414 END IF
1415 IF ((lobc(ieast,isvvel,ng)).and. &
1416 & domain(ng)%Eastern_Edge(tile)) THEN
1417 ib=ieast
1418 DO k=1,n(ng)
1419 DO j=jstrv,jend
1420 ad_v_obc(j,k,ib,ir,linp)=ad_v_obc(j,k,ib,ir,linp)* &
1421 & vmask(iend+1,j)
1422 END DO
1423 END DO
1424 END IF
1425 IF ((lobc(isouth,isvvel,ng)).and. &
1426 & domain(ng)%Southern_Edge(tile)) THEN
1427 ib=isouth
1428 DO k=1,n(ng)
1429 DO i=istr,iend
1430 ad_v_obc(i,k,ib,ir,linp)=ad_v_obc(i,k,ib,ir,linp)* &
1431 & vmask(i,jstr)
1432 END DO
1433 END DO
1434 END IF
1435 IF ((lobc(inorth,isvvel,ng)).and. &
1436 & domain(ng)%Northern_Edge(tile)) THEN
1437 ib=inorth
1438 DO k=1,n(ng)
1439 DO i=istr,iend
1440 ad_v_obc(i,k,ib,ir,linp)=ad_v_obc(i,k,ib,ir,linp)* &
1441 & vmask(i,jend+1)
1442 END DO
1443 END DO
1444 END IF
1445 END DO
1446 END IF
1447# endif
1448!
1449! Tracers.
1450!
1451 DO it=1,nt(ng)
1452 DO k=1,n(ng)
1453 DO j=jstrt,jendt
1454 DO i=istrt,iendt
1455 ad_t(i,j,k,linp,it)=ad_t(i,j,k,linp,it)*rmask(i,j)
1456 END DO
1457 END DO
1458 END DO
1459 END DO
1460
1461# ifdef ADJUST_BOUNDARY
1462!
1463! Tracers open boundaries.
1464!
1465 DO it=1,nt(ng)
1466 IF (any(lobc(:,istvar(it),ng))) THEN
1467 DO ir=1,nbrec(ng)
1468 IF ((lobc(iwest,istvar(it),ng)).and. &
1469 & domain(ng)%Western_Edge(tile)) THEN
1470 ib=iwest
1471 DO k=1,n(ng)
1472 DO j=jstr,jend
1473 ad_t_obc(j,k,ib,ir,linp,it)= &
1474 & ad_t_obc(j,k,ib,ir,linp,it)*rmask(istr-1,j)
1475 END DO
1476 END DO
1477 END IF
1478 IF ((lobc(ieast,istvar(it),ng)).and. &
1479 & domain(ng)%Eastern_Edge(tile)) THEN
1480 ib=ieast
1481 DO k=1,n(ng)
1482 DO j=jstr,jend
1483 ad_t_obc(j,k,ib,ir,linp,it)= &
1484 & ad_t_obc(j,k,ib,ir,linp,it)*rmask(iend+1,j)
1485 END DO
1486 END DO
1487 END IF
1488 IF ((lobc(isouth,istvar(it),ng)).and. &
1489 & domain(ng)%Southern_Edge(tile)) THEN
1490 ib=isouth
1491 DO k=1,n(ng)
1492 DO i=istr,iend
1493 ad_t_obc(i,k,ib,ir,linp,it)= &
1494 & ad_t_obc(i,k,ib,ir,linp,it)*rmask(i,jstr-1)
1495 END DO
1496 END DO
1497 END IF
1498 IF ((lobc(inorth,istvar(it),ng)).and. &
1499 & domain(ng)%Northern_Edge(tile)) THEN
1500 ib=inorth
1501 DO k=1,n(ng)
1502 DO i=istr,iend
1503 ad_t_obc(i,k,ib,ir,linp,it)= &
1504 & ad_t_obc(i,k,ib,ir,linp,it)*rmask(i,jend+1)
1505 END DO
1506 END DO
1507 END IF
1508 END DO
1509 END IF
1510 END DO
1511# endif
1512
1513# ifdef ADJUST_STFLUX
1514!
1515! Surface tracers flux.
1516!
1517 DO it=1,nt(ng)
1518 IF (lstflux(it,ng)) THEN
1519 DO ir=1,nfrec(ng)
1520 DO j=jstrt,jendt
1521 DO i=istrt,iendt
1522 ad_tflux(i,j,ir,linp,it)=ad_tflux(i,j,ir,linp,it)* &
1523 & rmask(i,j)
1524 END DO
1525 END DO
1526 END DO
1527 END IF
1528 END DO
1529# endif
1530# endif
1531
1532# endif
1533!
1534!-----------------------------------------------------------------------
1535! Load or add adjoint state variables into tangent linear state
1536! variables.
1537!-----------------------------------------------------------------------
1538!
1539! Add adjoint state to tangent linear state.
1540!
1541! tl_var(Lout) = fac1 * tl_var(Lout) + fac2 * ad_var(Linp)
1542!
1543 IF (add) THEN
1544
1545 fac1=1.0_r8
1546 fac2=1.0_r8
1547
1548 CALL state_addition (ng, tile, &
1549 & lbi, ubi, lbj, ubj, lbij, ubij, &
1550 & lout, linp, lout, fac1, fac2, &
1551# ifdef MASKING
1552 & rmask, umask, vmask, &
1553# endif
1554# ifdef ADJUST_BOUNDARY
1555# ifdef SOLVE3D
1556 & tl_t_obc, ad_t_obc, &
1557 & tl_u_obc, ad_u_obc, &
1558 & tl_v_obc, ad_v_obc, &
1559# endif
1560 & tl_ubar_obc, ad_ubar_obc, &
1561 & tl_vbar_obc, ad_vbar_obc, &
1562 & tl_zeta_obc, ad_zeta_obc, &
1563# endif
1564# ifdef ADJUST_WSTRESS
1565 & tl_ustr, ad_ustr, &
1566 & tl_vstr, ad_vstr, &
1567# endif
1568# ifdef SOLVE3D
1569# ifdef ADJUST_STFLUX
1570 & tl_tflux, ad_tflux, &
1571# endif
1572 & tl_t, ad_t, &
1573 & tl_u, ad_u, &
1574 & tl_v, ad_v, &
1575# if defined WEAK_CONSTRAINT && defined TIME_CONV
1576 & tl_ubar, ad_ubar, &
1577 & tl_vbar, ad_vbar, &
1578# endif
1579# else
1580 & tl_ubar, ad_ubar, &
1581 & tl_vbar, ad_vbar, &
1582# endif
1583 & tl_zeta, ad_zeta)
1584!
1585! Otherwise, copy adjoint state into tangent linear state.
1586!
1587! tl_var(Lout) = ad_var(Linp)
1588!
1589 ELSE
1590
1591 CALL state_copy (ng, tile, &
1592 & lbi, ubi, lbj, ubj, lbij, ubij, &
1593 & linp, lout, &
1594# ifdef ADJUST_BOUNDARY
1595# ifdef SOLVE3D
1596 & tl_t_obc, ad_t_obc, &
1597 & tl_u_obc, ad_u_obc, &
1598 & tl_v_obc, ad_v_obc, &
1599# endif
1600 & tl_ubar_obc, ad_ubar_obc, &
1601 & tl_vbar_obc, ad_vbar_obc, &
1602 & tl_zeta_obc, ad_zeta_obc, &
1603# endif
1604# ifdef ADJUST_WSTRESS
1605 & tl_ustr, ad_ustr, &
1606 & tl_vstr, ad_vstr, &
1607# endif
1608# ifdef SOLVE3D
1609# ifdef ADJUST_STFLUX
1610 & tl_tflux, ad_tflux, &
1611# endif
1612 & tl_t, ad_t, &
1613 & tl_u, ad_u, &
1614 & tl_v, ad_v, &
1615# if defined WEAK_CONSTRAINT && defined TIME_CONV
1616 & tl_ubar, ad_ubar, &
1617 & tl_vbar, ad_vbar, &
1618# endif
1619# else
1620 & tl_ubar, ad_ubar, &
1621 & tl_vbar, ad_vbar, &
1622# endif
1623 & tl_zeta, ad_zeta)
1624 END IF
1625!
1626 RETURN
logical, dimension(:,:,:), allocatable lobc
integer, parameter iwest
logical, dimension(:,:), allocatable lstflux
integer, dimension(:), allocatable nfrec
integer, parameter isouth
integer, parameter ieast
integer, parameter inorth
integer, dimension(:), allocatable nbrec
subroutine, public state_addition(ng, tile, lbi, ubi, lbj, ubj, lbij, ubij, lin1, lin2, lout, fac1, fac2, rmask, umask, vmask, s1_t_obc, s2_t_obc, s1_u_obc, s2_u_obc, s1_v_obc, s2_v_obc, s1_ubar_obc, s2_ubar_obc, s1_vbar_obc, s2_vbar_obc, s1_zeta_obc, s2_zeta_obc, s1_sustr, s2_sustr, s1_svstr, s2_svstr, s1_tflux, s2_tflux, s1_t, s2_t, s1_u, s2_u, s1_v, s2_v, s1_ubar, s2_ubar, s1_vbar, s2_vbar, s1_zeta, s2_zeta)
subroutine, public state_copy(ng, tile, lbi, ubi, lbj, ubj, lbij, ubij, linp, lout, s1_t_obc, s2_t_obc, s1_u_obc, s2_u_obc, s1_v_obc, s2_v_obc, s1_ubar_obc, s2_ubar_obc, s1_vbar_obc, s2_vbar_obc, s1_zeta_obc, s2_zeta_obc, s1_sustr, s2_sustr, s1_svstr, s2_svstr, s1_tflux, s2_tflux, s1_t, s2_t, s1_u, s2_u, s1_v, s2_v, s1_ubar, s2_ubar, s1_vbar, s2_vbar, s1_zeta, s2_zeta)
Definition state_copy.F:57

References mod_param::domain, mod_scalars::ieast, mod_scalars::inorth, mod_ncparam::isfsur, mod_scalars::isouth, mod_ncparam::istvar, mod_ncparam::isubar, mod_ncparam::isuvel, mod_ncparam::isvbar, mod_ncparam::isvvel, mod_scalars::iwest, mod_scalars::lobc, mod_scalars::lstflux, state_addition_mod::state_addition(), and state_copy_mod::state_copy().

Referenced by load_adtotl().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ load_tltoad()

subroutine, public ini_adjust_mod::load_tltoad ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) linp,
integer, intent(in) lout,
logical, intent(in) add )

Definition at line 1629 of file ini_adjust.F.

1630!
1631!=======================================================================
1632! !
1633! This routine loads or adds Linp tangent linear state variables into !
1634! Lout adjoint state variables. !
1635! !
1636! On Input: !
1637! !
1638! ng Nested grid number. !
1639! tile Domain partition. !
1640! Linp Tangent linear state time index to add. !
1641! Lout Nonlinear state time index to update. !
1642! add Logical switch to add to imported values. !
1643! !
1644!=======================================================================
1645!
1646 USE mod_param
1647# ifdef ADJUST_BOUNDARY
1648 USE mod_boundary
1649# endif
1650# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
1651 USE mod_forces
1652# endif
1653 USE mod_grid
1654 USE mod_ocean
1655!
1656! Imported variable declarations.
1657!
1658 logical, intent(in) :: add
1659 integer, intent(in) :: ng, tile, Linp, Lout
1660!
1661! Local variable declarations.
1662!
1663 character (len=*), parameter :: MyFile = &
1664 & __FILE__//", load_TLtoAD"
1665
1666# include "tile.h"
1667!
1668# ifdef PROFILE
1669 CALL wclock_on (ng, iadm, 7, __line__, myfile)
1670# endif
1671 CALL load_tltoad_tile (ng, tile, &
1672 & lbi, ubi, lbj, ubj, lbij, ubij, &
1673 & imins, imaxs, jmins, jmaxs, &
1674 & linp, lout, add, &
1675# ifdef MASKING
1676 & grid(ng) % rmask, &
1677 & grid(ng) % umask, &
1678 & grid(ng) % vmask, &
1679# endif
1680# ifdef ADJUST_BOUNDARY
1681# ifdef SOLVE3D
1682 & boundary(ng) % tl_t_obc, &
1683 & boundary(ng) % tl_u_obc, &
1684 & boundary(ng) % tl_v_obc, &
1685# endif
1686 & boundary(ng) % tl_ubar_obc, &
1687 & boundary(ng) % tl_vbar_obc, &
1688 & boundary(ng) % tl_zeta_obc, &
1689# endif
1690# ifdef ADJUST_WSTRESS
1691 & forces(ng) % tl_ustr, &
1692 & forces(ng) % tl_vstr, &
1693# endif
1694# ifdef SOLVE3D
1695# ifdef ADJUST_STFLUX
1696 & forces(ng) % tl_tflux, &
1697# endif
1698 & ocean(ng) % tl_t, &
1699 & ocean(ng) % tl_u, &
1700 & ocean(ng) % tl_v, &
1701# if defined WEAK_CONSTRAINT && defined TIME_CONV
1702 & ocean(ng) % tl_ubar, &
1703 & ocean(ng) % tl_vbar, &
1704# endif
1705# else
1706 & ocean(ng) % tl_ubar, &
1707 & ocean(ng) % tl_vbar, &
1708# endif
1709 & ocean(ng) % tl_zeta, &
1710# ifdef ADJUST_BOUNDARY
1711# ifdef SOLVE3D
1712 & boundary(ng) % ad_t_obc, &
1713 & boundary(ng) % ad_u_obc, &
1714 & boundary(ng) % ad_v_obc, &
1715# endif
1716 & boundary(ng) % ad_ubar_obc, &
1717 & boundary(ng) % ad_vbar_obc, &
1718 & boundary(ng) % ad_zeta_obc, &
1719# endif
1720# ifdef ADJUST_WSTRESS
1721 & forces(ng) % ad_ustr, &
1722 & forces(ng) % ad_vstr, &
1723# endif
1724# ifdef SOLVE3D
1725# ifdef ADJUST_STFLUX
1726 & forces(ng) % ad_tflux, &
1727# endif
1728 & ocean(ng) % ad_t, &
1729 & ocean(ng) % ad_u, &
1730 & ocean(ng) % ad_v, &
1731# if defined WEAK_CONSTRAINT && defined TIME_CONV
1732 & ocean(ng) % ad_ubar, &
1733 & ocean(ng) % ad_vbar, &
1734# endif
1735# else
1736 & ocean(ng) % ad_ubar, &
1737 & ocean(ng) % ad_vbar, &
1738# endif
1739 & ocean(ng) % ad_zeta)
1740# ifdef PROFILE
1741 CALL wclock_off (ng, iadm, 7, __line__, myfile)
1742# endif
1743!
1744 RETURN

References mod_boundary::boundary, mod_forces::forces, mod_grid::grid, mod_param::iadm, load_tltoad_tile(), mod_ocean::ocean, wclock_off(), and wclock_on().

Referenced by convolve_mod::convolve(), convolve_mod::error_covariance(), rbl4dvar_mod::increment(), roms_kernel_mod::roms_run(), and convolve_mod::saddlec().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ load_tltoad_tile()

subroutine ini_adjust_mod::load_tltoad_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) lbij,
integer, intent(in) ubij,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) linp,
integer, intent(in) lout,
logical, intent(in) add,
real(r8), dimension(lbi:,lbj:), intent(in) rmask,
real(r8), dimension(lbi:,lbj:), intent(in) umask,
real(r8), dimension(lbi:,lbj:), intent(in) vmask,
real(r8), dimension(lbij:,:,:,:,:,:), intent(inout) tl_t_obc,
real(r8), dimension(lbij:,:,:,:,:), intent(inout) tl_u_obc,
real(r8), dimension(lbij:,:,:,:,:), intent(inout) tl_v_obc,
real(r8), dimension(lbij:,:,:,:), intent(inout) tl_ubar_obc,
real(r8), dimension(lbij:,:,:,:), intent(inout) tl_vbar_obc,
real(r8), dimension(lbij:,:,:,:), intent(inout) tl_zeta_obc,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) tl_ustr,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) tl_vstr,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) tl_tflux,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) tl_t,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) tl_u,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) tl_v,
real(r8), dimension(lbi:,lbj:,:), intent(inout) tl_ubar,
real(r8), dimension(lbi:,lbj:,:), intent(inout) tl_vbar,
real(r8), dimension(lbi:,lbj:,:), intent(inout) tl_zeta,
real(r8), dimension(lbij:,:,:,:,:,:), intent(inout) ad_t_obc,
real(r8), dimension(lbij:,:,:,:,:), intent(inout) ad_u_obc,
real(r8), dimension(lbij:,:,:,:,:), intent(inout) ad_v_obc,
real(r8), dimension(lbij:,:,:,:), intent(inout) ad_ubar_obc,
real(r8), dimension(lbij:,:,:,:), intent(inout) ad_vbar_obc,
real(r8), dimension(lbij:,:,:,:), intent(inout) ad_zeta_obc,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) ad_ustr,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) ad_vstr,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) ad_tflux,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) ad_t,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) ad_u,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) ad_v,
real(r8), dimension(lbi:,lbj:,:), intent(inout) ad_ubar,
real(r8), dimension(lbi:,lbj:,:), intent(inout) ad_vbar,
real(r8), dimension(lbi:,lbj:,:), intent(inout) ad_zeta )
private

Definition at line 1748 of file ini_adjust.F.

1799!***********************************************************************
1800!
1801 USE mod_param
1802 USE mod_ncparam
1803# if defined ADJUST_BOUNDARY || defined ADJUST_STFLUX || \
1804 defined adjust_wstress
1805 USE mod_scalars
1806# endif
1807!
1809 USE state_copy_mod, ONLY : state_copy
1810!
1811! Imported variable declarations.
1812!
1813 logical, intent(in) :: add
1814 integer, intent(in) :: ng, tile
1815 integer, intent(in) :: LBi, UBi, LBj, UBj, LBij, UBij
1816 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
1817 integer, intent(in) :: Linp, Lout
1818!
1819# ifdef ASSUMED_SHAPE
1820# ifdef MASKING
1821 real(r8), intent(in) :: rmask(LBi:,LBj:)
1822 real(r8), intent(in) :: umask(LBi:,LBj:)
1823 real(r8), intent(in) :: vmask(LBi:,LBj:)
1824# endif
1825# ifdef ADJUST_BOUNDARY
1826# ifdef SOLVE3D
1827 real(r8), intent(inout) :: tl_t_obc(LBij:,:,:,:,:,:)
1828 real(r8), intent(inout) :: tl_u_obc(LBij:,:,:,:,:)
1829 real(r8), intent(inout) :: tl_v_obc(LBij:,:,:,:,:)
1830# endif
1831 real(r8), intent(inout) :: tl_ubar_obc(LBij:,:,:,:)
1832 real(r8), intent(inout) :: tl_vbar_obc(LBij:,:,:,:)
1833 real(r8), intent(inout) :: tl_zeta_obc(LBij:,:,:,:)
1834# endif
1835# ifdef ADJUST_WSTRESS
1836 real(r8), intent(inout) :: tl_ustr(LBi:,LBj:,:,:)
1837 real(r8), intent(inout) :: tl_vstr(LBi:,LBj:,:,:)
1838# endif
1839# ifdef SOLVE3D
1840# ifdef ADJUST_STFLUX
1841 real(r8), intent(inout) :: tl_tflux(LBi:,LBj:,:,:,:)
1842# endif
1843 real(r8), intent(inout) :: tl_t(LBi:,LBj:,:,:,:)
1844 real(r8), intent(inout) :: tl_u(LBi:,LBj:,:,:)
1845 real(r8), intent(inout) :: tl_v(LBi:,LBj:,:,:)
1846# if defined WEAK_CONSTRAINT && defined TIME_CONV
1847 real(r8), intent(inout) :: tl_ubar(LBi:,LBj:,:)
1848 real(r8), intent(inout) :: tl_vbar(LBi:,LBj:,:)
1849# endif
1850# else
1851 real(r8), intent(inout) :: tl_ubar(LBi:,LBj:,:)
1852 real(r8), intent(inout) :: tl_vbar(LBi:,LBj:,:)
1853# endif
1854 real(r8), intent(inout) :: tl_zeta(LBi:,LBj:,:)
1855# ifdef ADJUST_BOUNDARY
1856# ifdef SOLVE3D
1857 real(r8), intent(inout) :: ad_t_obc(LBij:,:,:,:,:,:)
1858 real(r8), intent(inout) :: ad_u_obc(LBij:,:,:,:,:)
1859 real(r8), intent(inout) :: ad_v_obc(LBij:,:,:,:,:)
1860# endif
1861 real(r8), intent(inout) :: ad_ubar_obc(LBij:,:,:,:)
1862 real(r8), intent(inout) :: ad_vbar_obc(LBij:,:,:,:)
1863 real(r8), intent(inout) :: ad_zeta_obc(LBij:,:,:,:)
1864# endif
1865# ifdef ADJUST_WSTRESS
1866 real(r8), intent(inout) :: ad_ustr(LBi:,LBj:,:,:)
1867 real(r8), intent(inout) :: ad_vstr(LBi:,LBj:,:,:)
1868# endif
1869# ifdef SOLVE3D
1870# ifdef ADJUST_STFLUX
1871 real(r8), intent(inout) :: ad_tflux(LBi:,LBj:,:,:,:)
1872# endif
1873 real(r8), intent(inout) :: ad_t(LBi:,LBj:,:,:,:)
1874 real(r8), intent(inout) :: ad_u(LBi:,LBj:,:,:)
1875 real(r8), intent(inout) :: ad_v(LBi:,LBj:,:,:)
1876# if defined WEAK_CONSTRAINT && defined TIME_CONV
1877 real(r8), intent(inout) :: ad_ubar(LBi:,LBj:,:)
1878 real(r8), intent(inout) :: ad_vbar(LBi:,LBj:,:)
1879# endif
1880# else
1881 real(r8), intent(inout) :: ad_ubar(LBi:,LBj:,:)
1882 real(r8), intent(inout) :: ad_vbar(LBi:,LBj:,:)
1883# endif
1884 real(r8), intent(inout) :: ad_zeta(LBi:,LBj:,:)
1885# else
1886# ifdef MASKING
1887 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
1888 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
1889 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
1890# endif
1891# ifdef ADJUST_BOUNDARY
1892# ifdef SOLVE3D
1893 real(r8), intent(inout) :: tl_t_obc(LBij:UBij,N(ng),4, &
1894 & Nbrec(ng),2,NT(ng))
1895 real(r8), intent(inout) :: tl_u_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
1896 real(r8), intent(inout) :: tl_v_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
1897# endif
1898 real(r8), intent(inout) :: tl_ubar_obc(LBij:UBij,4,Nbrec(ng),2)
1899 real(r8), intent(inout) :: tl_vbar_obc(LBij:UBij,4,Nbrec(ng),2)
1900 real(r8), intent(inout) :: tl_zeta_obc(LBij:UBij,4,Nbrec(ng),2)
1901# endif
1902# ifdef ADJUST_WSTRESS
1903 real(r8), intent(inout) :: tl_ustr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
1904 real(r8), intent(inout) :: tl_vstr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
1905# endif
1906# ifdef SOLVE3D
1907# ifdef ADJUST_STFLUX
1908 real(r8), intent(inout) :: tl_tflux(LBi:UBi,LBj:UBj, &
1909 & Nfrec(ng),2,NT(ng))
1910# endif
1911 real(r8), intent(inout) :: tl_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
1912 real(r8), intent(inout) :: tl_u(LBi:UBi,LBj:UBj,N(ng),2)
1913 real(r8), intent(inout) :: tl_v(LBi:UBi,LBj:UBj,N(ng),2)
1914# if defined WEAK_CONSTRAINT && defined TIME_CONV
1915 real(r8), intent(inout) :: tl_ubar(LBi:UBi,LBj:UBj,:)
1916 real(r8), intent(inout) :: tl_vbar(LBi:UBi,LBj:UBj,:)
1917# endif
1918# else
1919 real(r8), intent(inout) :: tl_ubar(LBi:UBi,LBj:UBj,:)
1920 real(r8), intent(inout) :: tl_vbar(LBi:UBi,LBj:UBj,:)
1921# endif
1922 real(r8), intent(inout) :: tl_zeta(LBi:UBi,LBj:UBj,:)
1923# ifdef ADJUST_BOUNDARY
1924# ifdef SOLVE3D
1925 real(r8), intent(inout) :: ad_t_obc(LBij:UBij,N(ng),4, &
1926 & Nbrec(ng),2,NT(ng))
1927 real(r8), intent(inout) :: ad_u_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
1928 real(r8), intent(inout) :: ad_v_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
1929# endif
1930 real(r8), intent(inout) :: ad_ubar_obc(LBij:UBij,4,Nbrec(ng),2)
1931 real(r8), intent(inout) :: ad_vbar_obc(LBij:UBij,4,Nbrec(ng),2)
1932 real(r8), intent(inout) :: ad_zeta_obc(LBij:UBij,4,Nbrec(ng),2)
1933# endif
1934# ifdef ADJUST_WSTRESS
1935 real(r8), intent(inout) :: ad_ustr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
1936 real(r8), intent(inout) :: ad_vstr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
1937# endif
1938# ifdef SOLVE3D
1939# ifdef ADJUST_STFLUX
1940 real(r8), intent(inout) :: ad_tflux(LBi:UBi,LBj:UBj, &
1941 & Nfrec(ng),2,NT(ng))
1942# endif
1943 real(r8), intent(inout) :: ad_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
1944 real(r8), intent(inout) :: ad_u(LBi:UBi,LBj:UBj,N(ng),2)
1945 real(r8), intent(inout) :: ad_v(LBi:UBi,LBj:UBj,N(ng),2)
1946# if defined WEAK_CONSTRAINT && defined TIME_CONV
1947 real(r8), intent(inout) :: ad_ubar(LBi:UBi,LBj:UBj,:)
1948 real(r8), intent(inout) :: ad_vbar(LBi:UBi,LBj:UBj,:)
1949# endif
1950# else
1951 real(r8), intent(inout) :: ad_ubar(LBi:UBi,LBj:UBj,:)
1952 real(r8), intent(inout) :: ad_vbar(LBi:UBi,LBj:UBj,:)
1953# endif
1954 real(r8), intent(inout) :: ad_zeta(LBi:UBi,LBj:UBj,:)
1955# endif
1956!
1957! Local variable declarations.
1958!
1959# ifdef MASKING
1960 integer :: i, j, k
1961 integer :: ib, ir, it
1962# endif
1963 real(r8) :: fac1, fac2
1964# ifdef MASKING
1965# include "set_bounds.h"
1966# endif
1967
1968# ifdef MASKING
1969!
1970!-----------------------------------------------------------------------
1971! Multiply by land-sea mask - fail safe. Notice that if add=.FALSE. the
1972! "state_addition" routine is not called and the state arrays are not
1973! masked.
1974!-----------------------------------------------------------------------
1975!
1976! Free-surface.
1977!
1978 DO j=jstrt,jendt
1979 DO i=istrt,iendt
1980 tl_zeta(i,j,linp)=tl_zeta(i,j,linp)*rmask(i,j)
1981 END DO
1982 END DO
1983
1984# ifdef ADJUST_BOUNDARY
1985!
1986! Free-surface open boundaries.
1987!
1988 IF (any(lobc(:,isfsur,ng))) THEN
1989 DO ir=1,nbrec(ng)
1990 IF ((lobc(iwest,isfsur,ng)).and. &
1991 & domain(ng)%Western_Edge(tile)) THEN
1992 ib=iwest
1993 DO j=jstr,jend
1994 tl_zeta_obc(j,ib,ir,linp)=tl_zeta_obc(j,ib,ir,linp)* &
1995 & rmask(istr-1,j)
1996 END DO
1997 END IF
1998 IF ((lobc(ieast,isfsur,ng)).and. &
1999 & domain(ng)%Eastern_Edge(tile)) THEN
2000 ib=ieast
2001 DO j=jstr,jend
2002 tl_zeta_obc(j,ib,ir,linp)=tl_zeta_obc(j,ib,ir,linp)* &
2003 & rmask(iend+1,j)
2004 END DO
2005 END IF
2006 IF ((lobc(isouth,isfsur,ng)).and. &
2007 & domain(ng)%Southern_Edge(tile)) THEN
2008 ib=isouth
2009 DO i=istr,iend
2010 tl_zeta_obc(i,ib,ir,linp)=tl_zeta_obc(i,ib,ir,linp)* &
2011 & rmask(i,jstr-1)
2012 END DO
2013 END IF
2014 IF ((lobc(inorth,isfsur,ng)).and. &
2015 & domain(ng)%Northern_Edge(tile)) THEN
2016 ib=inorth
2017 DO i=istr,iend
2018 tl_zeta_obc(i,ib,ir,linp)=tl_zeta_obc(i,ib,ir,linp)* &
2019 & rmask(i,jend+1)
2020 END DO
2021 END IF
2022 END DO
2023 END IF
2024# endif
2025
2026# if !defined SOLVE3D || \
2027 (defined weak_constraint && defined time_conv)
2028!
2029! 2D U-momentum component.
2030!
2031 DO j=jstrt,jendt
2032 DO i=istrp,iendt
2033 tl_ubar(i,j,linp)=tl_ubar(i,j,linp)*umask(i,j)
2034 END DO
2035 END DO
2036!
2037! 2D V-momentum component.
2038!
2039 DO j=jstrp,jendt
2040 DO i=istrt,iendt
2041 tl_vbar(i,j,linp)=tl_vbar(i,j,linp)*vmask(i,j)
2042 END DO
2043 END DO
2044# endif
2045
2046# ifdef ADJUST_BOUNDARY
2047!
2048! 2D U-momentum open boundaries.
2049!
2050 IF (any(lobc(:,isubar,ng))) THEN
2051 DO ir=1,nbrec(ng)
2052 IF ((lobc(iwest,isubar,ng)).and. &
2053 & domain(ng)%Western_Edge(tile)) THEN
2054 ib=iwest
2055 DO j=jstr,jend
2056 tl_ubar_obc(j,ib,ir,linp)=tl_ubar_obc(j,ib,ir,linp)* &
2057 & umask(istr,j)
2058 END DO
2059 END IF
2060 IF ((lobc(ieast,isubar,ng)).and. &
2061 & domain(ng)%Eastern_Edge(tile)) THEN
2062 ib=ieast
2063 DO j=jstr,jend
2064 tl_ubar_obc(j,ib,ir,linp)=tl_ubar_obc(j,ib,ir,linp)* &
2065 & umask(iend+1,j)
2066 END DO
2067 END IF
2068 IF ((lobc(isouth,isubar,ng)).and. &
2069 & domain(ng)%Southern_Edge(tile)) THEN
2070 ib=isouth
2071 DO i=istru,iend
2072 tl_ubar_obc(i,ib,ir,linp)=tl_ubar_obc(i,ib,ir,linp)* &
2073 & umask(i,jstr-1)
2074 END DO
2075 END IF
2076 IF ((lobc(inorth,isubar,ng)).and. &
2077 & domain(ng)%Northern_Edge(tile)) THEN
2078 ib=inorth
2079 DO i=istru,iend
2080 tl_ubar_obc(i,ib,ir,linp)=tl_ubar_obc(i,ib,ir,linp)* &
2081 & umask(i,jend+1)
2082 END DO
2083 END IF
2084 END DO
2085 END IF
2086!
2087! 2D V-momentum open boundaries.
2088!
2089 IF (any(lobc(:,isvbar,ng))) THEN
2090 DO ir=1,nbrec(ng)
2091 IF ((lobc(iwest,isvbar,ng)).and. &
2092 & domain(ng)%Western_Edge(tile)) THEN
2093 ib=iwest
2094 DO j=jstrv,jend
2095 tl_vbar_obc(j,ib,ir,linp)=tl_vbar_obc(j,ib,ir,linp)* &
2096 & vmask(istr-1,j)
2097 END DO
2098 END IF
2099 IF ((lobc(ieast,isvbar,ng)).and. &
2100 & domain(ng)%Eastern_Edge(tile)) THEN
2101 ib=ieast
2102 DO j=jstrv,jend
2103 tl_vbar_obc(j,ib,ir,linp)=tl_vbar_obc(j,ib,ir,linp)* &
2104 & vmask(iend+1,j)
2105 END DO
2106 END IF
2107 IF ((lobc(isouth,isvbar,ng)).and. &
2108 & domain(ng)%Southern_Edge(tile)) THEN
2109 ib=isouth
2110 DO i=istr,iend
2111 tl_vbar_obc(i,ib,ir,linp)=tl_vbar_obc(i,ib,ir,linp)* &
2112 & vmask(i,jstr)
2113 END DO
2114 END IF
2115 IF ((lobc(inorth,isvbar,ng)).and. &
2116 & domain(ng)%Northern_Edge(tile)) THEN
2117 ib=inorth
2118 DO i=istr,iend
2119 tl_vbar_obc(i,ib,ir,linp)=tl_vbar_obc(i,ib,ir,linp)* &
2120 & vmask(i,jend+1)
2121 END DO
2122 END IF
2123 END DO
2124 END IF
2125# endif
2126
2127# ifdef ADJUST_WSTRESS
2128!
2129! Surface momentum stress.
2130!
2131 DO ir=1,nfrec(ng)
2132 DO j=jstrt,jendt
2133 DO i=istrp,iendt
2134 tl_ustr(i,j,ir,linp)=tl_ustr(i,j,ir,linp)*umask(i,j)
2135 END DO
2136 END DO
2137 DO j=jstrp,jendt
2138 DO i=istrt,iendt
2139 tl_vstr(i,j,ir,linp)=tl_vstr(i,j,ir,linp)*vmask(i,j)
2140 END DO
2141 END DO
2142 END DO
2143# endif
2144
2145# ifdef SOLVE3D
2146!
2147! 3D U-momentum component.
2148!
2149 DO k=1,n(ng)
2150 DO j=jstrt,jendt
2151 DO i=istrp,iendt
2152 tl_u(i,j,k,linp)=tl_u(i,j,k,linp)*umask(i,j)
2153 END DO
2154 END DO
2155 END DO
2156
2157# ifdef ADJUST_BOUNDARY
2158!
2159! 3D U-momentum open boundaries.
2160!
2161 IF (any(lobc(:,isuvel,ng))) THEN
2162 DO ir=1,nbrec(ng)
2163 IF ((lobc(iwest,isuvel,ng)).and. &
2164 & domain(ng)%Western_Edge(tile)) THEN
2165 ib=iwest
2166 DO k=1,n(ng)
2167 DO j=jstr,jend
2168 tl_u_obc(j,k,ib,ir,linp)=tl_u_obc(j,k,ib,ir,linp)* &
2169 & umask(istr,j)
2170 END DO
2171 END DO
2172 END IF
2173 IF ((lobc(ieast,isuvel,ng)).and. &
2174 & domain(ng)%Eastern_Edge(tile)) THEN
2175 ib=ieast
2176 DO k=1,n(ng)
2177 DO j=jstr,jend
2178 tl_u_obc(j,k,ib,ir,linp)=tl_u_obc(j,k,ib,ir,linp)* &
2179 & umask(iend+1,j)
2180 END DO
2181 END DO
2182 END IF
2183 IF ((lobc(isouth,isuvel,ng)).and. &
2184 & domain(ng)%Southern_Edge(tile)) THEN
2185 ib=isouth
2186 DO k=1,n(ng)
2187 DO i=istru,iend
2188 tl_u_obc(i,k,ib,ir,linp)=tl_u_obc(i,k,ib,ir,linp)* &
2189 & umask(i,jstr-1)
2190 END DO
2191 END DO
2192 END IF
2193 IF ((lobc(inorth,isuvel,ng)).and. &
2194 & domain(ng)%Northern_Edge(tile)) THEN
2195 ib=inorth
2196 DO k=1,n(ng)
2197 DO i=istru,iend
2198 tl_u_obc(i,k,ib,ir,linp)=tl_u_obc(i,k,ib,ir,linp)* &
2199 & umask(i,jend+1)
2200 END DO
2201 END DO
2202 END IF
2203 END DO
2204 END IF
2205# endif
2206!
2207! 3D V-momentum component.
2208!
2209 DO k=1,n(ng)
2210 DO j=jstrp,jendt
2211 DO i=istrt,iendt
2212 tl_v(i,j,k,linp)=tl_v(i,j,k,linp)*vmask(i,j)
2213 END DO
2214 END DO
2215 END DO
2216
2217# ifdef ADJUST_BOUNDARY
2218!
2219! 3D V-momentum open boundaries.
2220!
2221 IF (any(lobc(:,isvvel,ng))) THEN
2222 DO ir=1,nbrec(ng)
2223 IF ((lobc(iwest,isvvel,ng)).and. &
2224 & domain(ng)%Western_Edge(tile)) THEN
2225 ib=iwest
2226 DO k=1,n(ng)
2227 DO j=jstrv,jend
2228 tl_v_obc(j,k,ib,ir,linp)=tl_v_obc(j,k,ib,ir,linp)* &
2229 & vmask(istr-1,j)
2230 END DO
2231 END DO
2232 END IF
2233 IF ((lobc(ieast,isvvel,ng)).and. &
2234 & domain(ng)%Eastern_Edge(tile)) THEN
2235 ib=ieast
2236 DO k=1,n(ng)
2237 DO j=jstrv,jend
2238 tl_v_obc(j,k,ib,ir,linp)=tl_v_obc(j,k,ib,ir,linp)* &
2239 & vmask(iend+1,j)
2240 END DO
2241 END DO
2242 END IF
2243 IF ((lobc(isouth,isvvel,ng)).and. &
2244 & domain(ng)%Southern_Edge(tile)) THEN
2245 ib=isouth
2246 DO k=1,n(ng)
2247 DO i=istr,iend
2248 tl_v_obc(i,k,ib,ir,linp)=tl_v_obc(i,k,ib,ir,linp)* &
2249 & vmask(i,jstr)
2250 END DO
2251 END DO
2252 END IF
2253 IF ((lobc(inorth,isvvel,ng)).and. &
2254 & domain(ng)%Northern_Edge(tile)) THEN
2255 ib=inorth
2256 DO k=1,n(ng)
2257 DO i=istr,iend
2258 tl_v_obc(i,k,ib,ir,linp)=tl_v_obc(i,k,ib,ir,linp)* &
2259 & vmask(i,jend+1)
2260 END DO
2261 END DO
2262 END IF
2263 END DO
2264 END IF
2265# endif
2266!
2267! Tracers.
2268!
2269 DO it=1,nt(ng)
2270 DO k=1,n(ng)
2271 DO j=jstrt,jendt
2272 DO i=istrt,iendt
2273 tl_t(i,j,k,linp,it)=tl_t(i,j,k,linp,it)*rmask(i,j)
2274 END DO
2275 END DO
2276 END DO
2277 END DO
2278
2279# ifdef ADJUST_BOUNDARY
2280!
2281! Tracers open boundaries.
2282!
2283 DO it=1,nt(ng)
2284 IF (any(lobc(:,istvar(it),ng))) THEN
2285 DO ir=1,nbrec(ng)
2286 IF ((lobc(iwest,istvar(it),ng)).and. &
2287 & domain(ng)%Western_Edge(tile)) THEN
2288 ib=iwest
2289 DO k=1,n(ng)
2290 DO j=jstr,jend
2291 tl_t_obc(j,k,ib,ir,linp,it)= &
2292 & tl_t_obc(j,k,ib,ir,linp,it)*rmask(istr-1,j)
2293 END DO
2294 END DO
2295 END IF
2296 IF ((lobc(ieast,istvar(it),ng)).and. &
2297 & domain(ng)%Eastern_Edge(tile)) THEN
2298 ib=ieast
2299 DO k=1,n(ng)
2300 DO j=jstr,jend
2301 tl_t_obc(j,k,ib,ir,linp,it)= &
2302 & tl_t_obc(j,k,ib,ir,linp,it)*rmask(iend+1,j)
2303 END DO
2304 END DO
2305 END IF
2306 IF ((lobc(isouth,istvar(it),ng)).and. &
2307 & domain(ng)%Southern_Edge(tile)) THEN
2308 ib=isouth
2309 DO k=1,n(ng)
2310 DO i=istr,iend
2311 tl_t_obc(i,k,ib,ir,linp,it)= &
2312 & tl_t_obc(i,k,ib,ir,linp,it)*rmask(i,jstr-1)
2313 END DO
2314 END DO
2315 END IF
2316 IF ((lobc(inorth,istvar(it),ng)).and. &
2317 & domain(ng)%Northern_Edge(tile)) THEN
2318 ib=inorth
2319 DO k=1,n(ng)
2320 DO i=istr,iend
2321 tl_t_obc(i,k,ib,ir,linp,it)= &
2322 & tl_t_obc(i,k,ib,ir,linp,it)*rmask(i,jend+1)
2323 END DO
2324 END DO
2325 END IF
2326 END DO
2327 END IF
2328 END DO
2329# endif
2330
2331# ifdef ADJUST_STFLUX
2332!
2333! Surface tracers flux.
2334!
2335 DO it=1,nt(ng)
2336 IF (lstflux(it,ng)) THEN
2337 DO ir=1,nfrec(ng)
2338 DO j=jstrt,jendt
2339 DO i=istrt,iendt
2340 tl_tflux(i,j,ir,linp,it)=tl_tflux(i,j,ir,linp,it)* &
2341 & rmask(i,j)
2342 END DO
2343 END DO
2344 END DO
2345 END IF
2346 END DO
2347# endif
2348# endif
2349
2350# endif
2351!
2352!-----------------------------------------------------------------------
2353! Load or add tangent linear state variables into adjoint state
2354! variables.
2355!-----------------------------------------------------------------------
2356!
2357! Add tangent linear state to adjoint state.
2358!
2359! ad_var(Lout) = fac1 * ad_var(Lout) + fac2 * tl_var(Linp)
2360!
2361 IF (add) THEN
2362
2363 fac1=1.0_r8
2364 fac2=1.0_r8
2365
2366 CALL state_addition (ng, tile, &
2367 & lbi, ubi, lbj, ubj, lbij, ubij, &
2368 & lout, linp, lout, fac1, fac2, &
2369# ifdef MASKING
2370 & rmask, umask, vmask, &
2371# endif
2372# ifdef ADJUST_BOUNDARY
2373# ifdef SOLVE3D
2374 & ad_t_obc, tl_t_obc, &
2375 & ad_u_obc, tl_u_obc, &
2376 & ad_v_obc, tl_v_obc, &
2377# endif
2378 & ad_ubar_obc, tl_ubar_obc, &
2379 & ad_vbar_obc, tl_vbar_obc, &
2380 & ad_zeta_obc, tl_zeta_obc, &
2381# endif
2382# ifdef ADJUST_WSTRESS
2383 & ad_ustr, tl_ustr, &
2384 & ad_vstr, tl_vstr, &
2385# endif
2386# ifdef SOLVE3D
2387# ifdef ADJUST_STFLUX
2388 & ad_tflux, tl_tflux, &
2389# endif
2390 & ad_t, tl_t, &
2391 & ad_u, tl_u, &
2392 & ad_v, tl_v, &
2393# if defined WEAK_CONSTRAINT && defined TIME_CONV
2394 & ad_ubar, tl_ubar, &
2395 & ad_vbar, tl_vbar, &
2396# endif
2397# else
2398 & ad_ubar, tl_ubar, &
2399 & ad_vbar, tl_vbar, &
2400# endif
2401 & ad_zeta, tl_zeta)
2402!
2403! Otherwise, copy tangent linear state into adjoint state.
2404!
2405! ad_var(Lout) = tl_var(Linp)
2406!
2407 ELSE
2408
2409 CALL state_copy (ng, tile, &
2410 & lbi, ubi, lbj, ubj, lbij, ubij, &
2411 & linp, lout, &
2412# ifdef ADJUST_BOUNDARY
2413# ifdef SOLVE3D
2414 & ad_t_obc, tl_t_obc, &
2415 & ad_u_obc, tl_u_obc, &
2416 & ad_v_obc, tl_v_obc, &
2417# endif
2418 & ad_ubar_obc, tl_ubar_obc, &
2419 & ad_vbar_obc, tl_vbar_obc, &
2420 & ad_zeta_obc, tl_zeta_obc, &
2421# endif
2422# ifdef ADJUST_WSTRESS
2423 & ad_ustr, tl_ustr, &
2424 & ad_vstr, tl_vstr, &
2425# endif
2426# ifdef SOLVE3D
2427# ifdef ADJUST_STFLUX
2428 & ad_tflux, tl_tflux, &
2429# endif
2430 & ad_t, tl_t, &
2431 & ad_u, tl_u, &
2432 & ad_v, tl_v, &
2433# if defined WEAK_CONSTRAINT && defined TIME_CONV
2434 & ad_ubar, tl_ubar, &
2435 & ad_vbar, tl_vbar, &
2436# endif
2437# else
2438 & ad_ubar, tl_ubar, &
2439 & ad_vbar, tl_vbar, &
2440# endif
2441 & ad_zeta, tl_zeta)
2442 END IF
2443!
2444 RETURN

References mod_param::domain, mod_scalars::ieast, mod_scalars::inorth, mod_ncparam::isfsur, mod_scalars::isouth, mod_ncparam::istvar, mod_ncparam::isubar, mod_ncparam::isuvel, mod_ncparam::isvbar, mod_ncparam::isvvel, mod_scalars::iwest, mod_scalars::lobc, mod_scalars::lstflux, state_addition_mod::state_addition(), and state_copy_mod::state_copy().

Referenced by load_tltoad().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ rp_ini_adjust()

subroutine, public ini_adjust_mod::rp_ini_adjust ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) linp,
integer, intent(in) lout )

Definition at line 537 of file ini_adjust.F.

538!
539!=======================================================================
540! !
541! This routine adds weak constraint adjoint increments to nonlinear !
542! state initial conditions (background state) at the end of inner !
543! loop. The boundary condition and barotropic/baroclinic coupling, !
544! if any, are processed latter in routine "ini_fields". !
545! !
546! On Input: !
547! !
548! ng Nested grid number. !
549! Linp Tangent linear state time index to add. !
550! Lout Nonlinear state time index to update. !
551! tile Domain partition. !
552! !
553!=======================================================================
554!
555 USE mod_param
556 USE mod_grid
557 USE mod_ocean
558!
559! Imported variable declarations.
560!
561 integer, intent(in) :: ng, tile, Linp, Lout
562!
563! Local variable declarations.
564!
565 character (len=*), parameter :: MyFile = &
566 & __FILE__//", rp_ini_adjust"
567
568# include "tile.h"
569!
570# ifdef PROFILE
571 CALL wclock_on (ng, irpm, 7, __line__, myfile)
572# endif
573 CALL rp_ini_adjust_tile (ng, tile, &
574 & lbi, ubi, lbj, ubj, &
575 & imins, imaxs, jmins, jmaxs, &
576 & linp, lout, &
577# ifdef MASKING
578 & grid(ng) % rmask, &
579 & grid(ng) % umask, &
580 & grid(ng) % vmask, &
581# endif
582# ifdef SOLVE3D
583 & ocean(ng) % ad_t, &
584 & ocean(ng) % ad_u, &
585 & ocean(ng) % ad_v, &
586# else
587 & ocean(ng) % ad_ubar, &
588 & ocean(ng) % ad_vbar, &
589# endif
590 & ocean(ng) % ad_zeta, &
591# ifdef SOLVE3D
592 & ocean(ng) % t, &
593 & ocean(ng) % u, &
594 & ocean(ng) % v, &
595# else
596 & ocean(ng) % ubar, &
597 & ocean(ng) % vbar, &
598# endif
599 & ocean(ng) % zeta, &
600# ifdef SOLVE3D
601 & ocean(ng) % tl_t, &
602 & ocean(ng) % tl_u, &
603 & ocean(ng) % tl_v, &
604# else
605 & ocean(ng) % tl_ubar, &
606 & ocean(ng) % tl_vbar, &
607# endif
608 & ocean(ng) % tl_zeta)
609# ifdef PROFILE
610 CALL wclock_off (ng, irpm, 7, __line__, myfile)
611# endif
612!
613 RETURN
integer, parameter irpm
Definition mod_param.F:664

References mod_grid::grid, mod_param::irpm, mod_ocean::ocean, rp_ini_adjust_tile(), wclock_off(), and wclock_on().

Referenced by convolve_mod::error_covariance().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ rp_ini_adjust_tile()

subroutine ini_adjust_mod::rp_ini_adjust_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) linp,
integer, intent(in) lout,
real(r8), dimension(lbi:,lbj:), intent(in) rmask,
real(r8), dimension(lbi:,lbj:), intent(in) umask,
real(r8), dimension(lbi:,lbj:), intent(in) vmask,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(in) ad_t,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) ad_u,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) ad_v,
real(r8), dimension(lbi:,lbj:,:), intent(in) ad_zeta,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(in) t,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) u,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) v,
real(r8), dimension(lbi:,lbj:,:), intent(in) zeta,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(out) tl_t,
real(r8), dimension(lbi:,lbj:,:,:), intent(out) tl_u,
real(r8), dimension(lbi:,lbj:,:,:), intent(out) tl_v,
real(r8), dimension(lbi:,lbj:,:), intent(out) tl_zeta )
private

Definition at line 617 of file ini_adjust.F.

642!***********************************************************************
643!
644 USE mod_param
645!
646! Imported variable declarations.
647!
648 integer, intent(in) :: ng, tile
649 integer, intent(in) :: LBi, UBi, LBj, UBj
650 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
651 integer, intent(in) :: Linp, Lout
652!
653# ifdef ASSUMED_SHAPE
654# ifdef MASKING
655 real(r8), intent(in) :: rmask(LBi:,LBj:)
656 real(r8), intent(in) :: umask(LBi:,LBj:)
657 real(r8), intent(in) :: vmask(LBi:,LBj:)
658# endif
659# ifdef SOLVE3D
660 real(r8), intent(in) :: ad_t(LBi:,LBj:,:,:,:)
661 real(r8), intent(in) :: ad_u(LBi:,LBj:,:,:)
662 real(r8), intent(in) :: ad_v(LBi:,LBj:,:,:)
663# else
664 real(r8), intent(in) :: ad_ubar(LBi:,LBj:,:)
665 real(r8), intent(in) :: ad_vbar(LBi:,LBj:,:)
666# endif
667 real(r8), intent(in) :: ad_zeta(LBi:,LBj:,:)
668# ifdef SOLVE3D
669 real(r8), intent(in) :: t(LBi:,LBj:,:,:,:)
670 real(r8), intent(in) :: u(LBi:,LBj:,:,:)
671 real(r8), intent(in) :: v(LBi:,LBj:,:,:)
672# else
673 real(r8), intent(in) :: ubar(LBi:,LBj:,:)
674 real(r8), intent(in) :: vbar(LBi:,LBj:,:)
675# endif
676 real(r8), intent(in) :: zeta(LBi:,LBj:,:)
677# ifdef SOLVE3D
678 real(r8), intent(out) :: tl_t(LBi:,LBj:,:,:,:)
679 real(r8), intent(out) :: tl_u(LBi:,LBj:,:,:)
680 real(r8), intent(out) :: tl_v(LBi:,LBj:,:,:)
681# else
682 real(r8), intent(out) :: tl_ubar(LBi:,LBj:,:)
683 real(r8), intent(out) :: tl_vbar(LBi:,LBj:,:)
684# endif
685 real(r8), intent(out) :: tl_zeta(LBi:,LBj:,:)
686# else
687# ifdef MASKING
688 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
689 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
690 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
691# endif
692# ifdef SOLVE3D
693 real(r8), intent(in) :: ad_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
694 real(r8), intent(in) :: ad_u(LBi:UBi,LBj:UBj,N(ng),2)
695 real(r8), intent(in) :: ad_v(LBi:UBi,LBj:UBj,N(ng),2)
696# else
697 real(r8), intent(in) :: ad_ubar(LBi:UBi,LBj:UBj,:)
698 real(r8), intent(in) :: ad_vbar(LBi:UBi,LBj:UBj,:)
699# endif
700 real(r8), intent(in) :: ad_zeta(LBi:UBi,LBj:UBj,:)
701# ifdef SOLVE3D
702 real(r8), intent(in) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
703 real(r8), intent(in) :: u(LBi:UBi,LBj:UBj,N(ng),2)
704 real(r8), intent(in) :: v(LBi:UBi,LBj:UBj,N(ng),2)
705# else
706 real(r8), intent(in) :: ubar(LBi:UBi,LBj:UBj,:)
707 real(r8), intent(in) :: vbar(LBi:UBi,LBj:UBj,:)
708# endif
709 real(r8), intent(in) :: zeta(LBi:UBi,LBj:UBj,:)
710# ifdef SOLVE3D
711 real(r8), intent(out) :: tl_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
712 real(r8), intent(out) :: tl_u(LBi:UBi,LBj:UBj,N(ng),2)
713 real(r8), intent(out) :: tl_v(LBi:UBi,LBj:UBj,N(ng),2)
714# else
715 real(r8), intent(out) :: tl_ubar(LBi:UBi,LBj:UBj,:)
716 real(r8), intent(out) :: tl_vbar(LBi:UBi,LBj:UBj,:)
717# endif
718 real(r8), intent(out) :: tl_zeta(LBi:UBi,LBj:UBj,:)
719# endif
720!
721! Local variable declarations.
722!
723 integer :: i, j
724# ifdef SOLVE3D
725 integer :: itrc, k
726# endif
727
728# include "set_bounds.h"
729!
730!-----------------------------------------------------------------------
731! Adjust initial conditions for 2D state variables by adding adjoint
732! increments from weak constraint inner loop.
733!-----------------------------------------------------------------------
734!
735# ifndef SOLVE3D
736 DO j=jstrt,jendt
737 DO i=istrp,iendt
738 tl_ubar(i,j,lout)=ubar(i,j,linp)+ad_ubar(i,j,lout)
739# ifdef MASKING
740 tl_ubar(i,j,lout)=tl_ubar(i,j,lout)*umask(i,j)
741# endif
742 END DO
743 END DO
744
745 DO j=jstrp,jendt
746 DO i=istrt,iendt
747 tl_vbar(i,j,lout)=vbar(i,j,linp)+ad_vbar(i,j,lout)
748# ifdef MASKING
749 tl_vbar(i,j,lout)=tl_vbar(i,j,lout)*vmask(i,j)
750# endif
751 END DO
752 END DO
753# endif
754
755 DO j=jstrt,jendt
756 DO i=istrt,iendt
757 tl_zeta(i,j,lout)=zeta(i,j,linp)+ad_zeta(i,j,lout)
758# ifdef MASKING
759 tl_zeta(i,j,lout)=tl_zeta(i,j,lout)*rmask(i,j)
760# endif
761 END DO
762 END DO
763
764# ifdef SOLVE3D
765!
766!-----------------------------------------------------------------------
767! Adjust initial conditions for 3D state variables by adding adjoint
768! increments from weak constraint inner loop.
769!-----------------------------------------------------------------------
770!
771 DO k=1,n(ng)
772 DO j=jstrt,jendt
773 DO i=istrp,iendt
774 tl_u(i,j,k,lout)=u(i,j,k,linp)+ad_u(i,j,k,lout)
775# ifdef MASKING
776 tl_u(i,j,k,lout)=tl_u(i,j,k,lout)*umask(i,j)
777# endif
778 END DO
779 END DO
780
781 DO j=jstrp,jendt
782 DO i=istrt,iendt
783 tl_v(i,j,k,lout)=v(i,j,k,linp)+ad_v(i,j,k,lout)
784# ifdef MASKING
785 tl_v(i,j,k,lout)=tl_v(i,j,k,lout)*vmask(i,j)
786# endif
787 END DO
788 END DO
789 END DO
790
791 DO itrc=1,nt(ng)
792 DO k=1,n(ng)
793 DO j=jstrt,jendt
794 DO i=istrt,iendt
795 tl_t(i,j,k,lout,itrc)=t(i,j,k,linp,itrc)+ &
796 & ad_t(i,j,k,lout,itrc)
797# ifdef MASKING
798 tl_t(i,j,k,lout,itrc)=tl_t(i,j,k,lout,itrc)*rmask(i,j)
799# endif
800 END DO
801 END DO
802 END DO
803 END DO
804# endif
805!
806 RETURN

Referenced by rp_ini_adjust().

Here is the caller graph for this function:

◆ tl_ini_perturb()

subroutine, public ini_adjust_mod::tl_ini_perturb ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) linp,
integer, intent(in) lout )

Definition at line 3356 of file ini_adjust.F.

3357!
3358!=======================================================================
3359! !
3360! Initialize tangent linear state variable according to the outer !
3361! and inner loop iterations. The initial field is a function of !
3362! the steepest descent direction (grad(J)) times the perturbation !
3363! amplitude "p" (controlled by inner loop). !
3364! !
3365!=======================================================================
3366!
3367 USE mod_param
3368# ifdef SOLVE3D
3369 USE mod_coupling
3370# endif
3371 USE mod_grid
3372 USE mod_ocean
3373# if defined SEDIMENT && defined SED_MORPH && defined SOLVE3D
3374 USE mod_sedbed
3375# endif
3376 USE mod_stepping
3377!
3378! Imported variable declarations.
3379!
3380 integer, intent(in) :: ng, tile, Linp, Lout
3381!
3382! Local variable declarations.
3383!
3384 character (len=*), parameter :: MyFile = &
3385 & __FILE__//", tl_ini_perturb"
3386!
3387# include "tile.h"
3388!
3389# ifdef PROFILE
3390 CALL wclock_on (ng, itlm, 7, __line__, myfile)
3391# endif
3392 CALL tl_ini_perturb_tile (ng, tile, &
3393 & lbi, ubi, lbj, ubj, &
3394 & imins, imaxs, jmins, jmaxs, &
3395 & nstp(ng), nnew(ng), linp, lout, &
3396# ifdef MASKING
3397 & grid(ng) % rmask, &
3398 & grid(ng) % umask, &
3399 & grid(ng) % vmask, &
3400# endif
3401# ifdef SOLVE3D
3402# if defined SEDIMENT && defined SED_MORPH
3403 & sedbed(ng) % tl_bed_thick, &
3404# endif
3405 & grid(ng) % tl_Hz, &
3406 & grid(ng) % h, &
3407 & grid(ng) % tl_h, &
3408 & grid(ng) % om_v, &
3409 & grid(ng) % on_u, &
3410# ifdef ICESHELF
3411 & grid(ng) % zice, &
3412# endif
3413 & grid(ng) % tl_z_r, &
3414 & grid(ng) % tl_z_w, &
3415 & coupling(ng) % Zt_avg1, &
3416 & coupling(ng) % tl_Zt_avg1, &
3417# endif
3418 & ocean(ng) % ubar, &
3419 & ocean(ng) % vbar, &
3420 & ocean(ng) % zeta, &
3421# ifdef SOLVE3D
3422 & ocean(ng) % ad_t, &
3423 & ocean(ng) % ad_u, &
3424 & ocean(ng) % ad_v, &
3425 & ocean(ng) % tl_t, &
3426 & ocean(ng) % tl_u, &
3427 & ocean(ng) % tl_v, &
3428# endif
3429 & ocean(ng) % ad_ubar, &
3430 & ocean(ng) % ad_vbar, &
3431 & ocean(ng) % ad_zeta, &
3432 & ocean(ng) % tl_ubar, &
3433 & ocean(ng) % tl_vbar, &
3434 & ocean(ng) % tl_zeta)
3435# ifdef PROFILE
3436 CALL wclock_off (ng, itlm, 7, __line__, myfile)
3437# endif
3438!
3439 RETURN

References mod_coupling::coupling, mod_grid::grid, mod_param::itlm, mod_stepping::nnew, mod_stepping::nstp, mod_ocean::ocean, mod_sedbed::sedbed, tl_ini_perturb_tile(), wclock_off(), and wclock_on().

Referenced by rp_initial(), and tl_initial().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ tl_ini_perturb_tile()

subroutine ini_adjust_mod::tl_ini_perturb_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) nstp,
integer, intent(in) nnew,
integer, intent(in) linp,
integer, intent(in) lout,
real(r8), dimension(lbi:,lbj:), intent(in) rmask,
real(r8), dimension(lbi:,lbj:), intent(in) umask,
real(r8), dimension(lbi:,lbj:), intent(in) vmask,
real(r8), dimension(lbi:,lbj:,:), intent(in) tl_bed_thick,
real(r8), dimension(lbi:,lbj:,:), intent(out) tl_hz,
real(r8), dimension(lbi:,lbj:), intent(in) h,
real(r8), dimension(lbi:,lbj:), intent(inout) tl_h,
real(r8), dimension(lbi:,lbj:), intent(in) om_v,
real(r8), dimension(lbi:,lbj:), intent(in) on_u,
real(r8), dimension(lbi:,lbj:), intent(in) zice,
real(r8), dimension(lbi:,lbj:,:), intent(out) tl_z_r,
real(r8), dimension(lbi:,lbj:,0:), intent(out) tl_z_w,
real(r8), dimension(lbi:,lbj:), intent(in) zt_avg1,
real(r8), dimension(lbi:,lbj:), intent(inout) tl_zt_avg1,
real(r8), dimension(lbi:,lbj:,:), intent(in) ubar,
real(r8), dimension(lbi:,lbj:,:), intent(in) vbar,
real(r8), dimension(lbi:,lbj:,:), intent(in) zeta,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(in) ad_t,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) ad_u,
real(r8), dimension(lbi:,lbj:,:,:), intent(in) ad_v,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) tl_t,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) tl_u,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) tl_v,
real(r8), dimension(lbi:,lbj:,:), intent(in) ad_ubar,
real(r8), dimension(lbi:,lbj:,:), intent(in) ad_vbar,
real(r8), dimension(lbi:,lbj:,:), intent(in) ad_zeta,
real(r8), dimension(lbi:,lbj:,:), intent(inout) tl_ubar,
real(r8), dimension(lbi:,lbj:,:), intent(inout) tl_vbar,
real(r8), dimension(lbi:,lbj:,:), intent(inout) tl_zeta )
private

Definition at line 3443 of file ini_adjust.F.

3469!***********************************************************************
3470!
3471 USE mod_param
3472 USE mod_parallel
3473 USE mod_fourdvar
3474 USE mod_ncparam
3475 USE mod_iounits
3476 USE mod_scalars
3477!
3478 USE exchange_2d_mod
3479# ifdef SOLVE3D
3480 USE exchange_3d_mod
3481# endif
3482# ifdef DISTRIBUTE
3483 USE mp_exchange_mod, ONLY : mp_exchange2d
3484# ifdef SOLVE3D
3486# endif
3487# endif
3488# ifdef SOLVE3D
3490# endif
3491 USE tl_u2dbc_mod, ONLY : tl_u2dbc_tile
3492 USE tl_v2dbc_mod, ONLY : tl_v2dbc_tile
3493 USE tl_zetabc_mod, ONLY : tl_zetabc_tile
3494# ifdef SOLVE3D
3495 USE tl_t3dbc_mod, ONLY : tl_t3dbc_tile
3496 USE tl_u3dbc_mod, ONLY : tl_u3dbc_tile
3497 USE tl_v3dbc_mod, ONLY : tl_v3dbc_tile
3498# endif
3499!
3500! Imported variable declarations.
3501!
3502 integer, intent(in) :: ng, tile
3503 integer, intent(in) :: LBi, UBi, LBj, UBj
3504 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
3505 integer, intent(in) :: nstp, nnew, Linp, Lout
3506!
3507# ifdef ASSUMED_SHAPE
3508# ifdef MASKING
3509 real(r8), intent(in) :: rmask(LBi:,LBj:)
3510 real(r8), intent(in) :: umask(LBi:,LBj:)
3511 real(r8), intent(in) :: vmask(LBi:,LBj:)
3512# endif
3513# ifdef SOLVE3D
3514# if defined SEDIMENT && defined SED_MORPH
3515 real(r8), intent(in) :: tl_bed_thick(LBi:,LBj:,:)
3516# endif
3517 real(r8), intent(in) :: om_v(LBi:,LBj:)
3518 real(r8), intent(in) :: on_u(LBi:,LBj:)
3519# ifdef ICESHELF
3520 real(r8), intent(in) :: zice(LBi:,LBj:)
3521# endif
3522 real(r8), intent(in) :: h(LBi:,LBj:)
3523 real(r8), intent(in) :: Zt_avg1(LBi:,LBj:)
3524# endif
3525 real(r8), intent(in) :: ubar(LBi:,LBj:,:)
3526 real(r8), intent(in) :: vbar(LBi:,LBj:,:)
3527 real(r8), intent(in) :: zeta(LBi:,LBj:,:)
3528# ifdef SOLVE3D
3529 real(r8), intent(in) :: ad_t(LBi:,LBj:,:,:,:)
3530 real(r8), intent(in) :: ad_u(LBi:,LBj:,:,:)
3531 real(r8), intent(in) :: ad_v(LBi:,LBj:,:,:)
3532# endif
3533 real(r8), intent(in) :: ad_ubar(LBi:,LBj:,:)
3534 real(r8), intent(in) :: ad_vbar(LBi:,LBj:,:)
3535 real(r8), intent(in) :: ad_zeta(LBi:,LBj:,:)
3536# ifdef SOLVE3D
3537 real(r8), intent(inout) :: tl_h(LBi:,LBj:)
3538 real(r8), intent(inout) :: tl_Zt_avg1(LBi:,LBj:)
3539# endif
3540# ifdef SOLVE3D
3541 real(r8), intent(inout) :: tl_t(LBi:,LBj:,:,:,:)
3542 real(r8), intent(inout) :: tl_u(LBi:,LBj:,:,:)
3543 real(r8), intent(inout) :: tl_v(LBi:,LBj:,:,:)
3544# endif
3545 real(r8), intent(inout) :: tl_ubar(LBi:,LBj:,:)
3546 real(r8), intent(inout) :: tl_vbar(LBi:,LBj:,:)
3547 real(r8), intent(inout) :: tl_zeta(LBi:,LBj:,:)
3548# ifdef SOLVE3D
3549 real(r8), intent(out) :: tl_Hz(LBi:,LBj:,:)
3550 real(r8), intent(out) :: tl_z_r(LBi:,LBj:,:)
3551 real(r8), intent(out) :: tl_z_w(LBi:,LBj:,0:)
3552# endif
3553# else
3554# ifdef MASKING
3555 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
3556 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
3557 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
3558# endif
3559# ifdef SOLVE3D
3560# if defined SEDIMENT && defined SED_MORPH
3561 real(r8), intent(in) :: tl_bed_thick(LBi:UBi,LBj:UBj,3)
3562# endif
3563 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
3564 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
3565# ifdef ICESHELF
3566 real(r8), intent(in) :: zice(LBi:,LBj:)
3567# endif
3568 real(r8), intent(in) :: Zt_avg1(LBi:UBi,LBj:UBj)
3569# endif
3570 real(r8), intent(in) :: ubar(LBi:UBi,LBj:UBj,:)
3571 real(r8), intent(in) :: vbar(LBi:UBi,LBj:UBj,:)
3572 real(r8), intent(in) :: zeta(LBi:UBi,LBj:UBj,:)
3573# ifdef SOLVE3D
3574 real(r8), intent(in) :: ad_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
3575 real(r8), intent(in) :: ad_u(LBi:UBi,LBj:UBj,N(ng),2)
3576 real(r8), intent(in) :: ad_v(LBi:UBi,LBj:UBj,N(ng),2)
3577# endif
3578 real(r8), intent(in) :: ad_ubar(LBi:UBi,LBj:UBj,:)
3579 real(r8), intent(in) :: ad_vbar(LBi:UBi,LBj:UBj,:)
3580 real(r8), intent(in) :: ad_zeta(LBi:UBi,LBj:UBj,:)
3581# ifdef SOLVE3D
3582 real(r8), intent(inout) :: tl_Zt_avg1(LBi:UBi,LBj:UBj)
3583# endif
3584# ifdef SOLVE3D
3585 real(r8), intent(inout) :: tl_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
3586 real(r8), intent(inout) :: tl_u(LBi:UBi,LBj:UBj,N(ng),2)
3587 real(r8), intent(inout) :: tl_v(LBi:UBi,LBj:UBj,N(ng),2)
3588# endif
3589 real(r8), intent(inout) :: h(LBi:UBi,LBj:UBj)
3590 real(r8), intent(inout) :: tl_h(LBi:UBi,LBj:UBj)
3591 real(r8), intent(inout) :: tl_ubar(LBi:UBi,LBj:UBj,:)
3592 real(r8), intent(inout) :: tl_vbar(LBi:UBi,LBj:UBj,:)
3593 real(r8), intent(inout) :: tl_zeta(LBi:UBi,LBj:UBj,:)
3594# ifdef SOLVE3D
3595 real(r8), intent(out) :: tl_Hz(LBi:UBi,LBj:UBj,N(ng))
3596 real(r8), intent(out) :: tl_z_r(LBi:UBi,LBj:UBj,N(ng))
3597 real(r8), intent(out) :: tl_z_w(LBi:UBi,LBj:UBj,0:N(ng))
3598# endif
3599# endif
3600!
3601! Local variable declarations.
3602!
3603 integer :: i, ic, j
3604# ifdef SOLVE3D
3605 integer :: itrc, k
3606# endif
3607 integer, dimension(NstateVar(ng)+1) :: StateVar
3608
3609 real(r8) :: p, scale
3610!
3611# include "set_bounds.h"
3612!
3613!-----------------------------------------------------------------------
3614! Initialize tangent linear with the steepest descent direction times
3615! the perturbation amplitude "p".
3616!-----------------------------------------------------------------------
3617!
3618! Set state variable to perturb according to outer loop index.
3619!
3620# ifdef SOLVE3D
3621 statevar(1)=0
3622 statevar(2)=isfsur
3623 statevar(3)=isubar
3624 statevar(4)=isvbar
3625 statevar(5)=isuvel
3626 statevar(6)=isvvel
3627 DO i=1,nt(ng)
3628 statevar(6+i)=istvar(i)
3629 END DO
3630# else
3631 statevar(1)=0
3632 statevar(2)=isfsur
3633 statevar(3)=isubar
3634 statevar(4)=isvbar
3635# endif
3636!
3637! Set perturbation amplitude as a function of the inner loop.
3638!
3639 p=10.0_r8**real(-inner,r8)
3640 scale=1.0_r8/sqrt(addotproduct)
3641 IF (domain(ng)%SouthWest_Test(tile)) THEN
3642 IF (master) WRITE (stdout,10)
3643 END IF
3644!
3645! Free-surface.
3646!
3647 IF ((statevar(outer).eq.0).or.(statevar(outer).eq.isfsur)) THEN
3648 DO j=jstrb,jendb
3649 DO i=istrb,iendb
3650 tl_zeta(i,j,lout)=p*ad_zeta(i,j,linp)*scale
3651# ifdef MASKING
3652 tl_zeta(i,j,lout)=tl_zeta(i,j,lout)*rmask(i,j)
3653# endif
3654 END DO
3655 END DO
3656 IF (domain(ng)%SouthWest_Test(tile)) THEN
3657 IF (master) THEN
3658 WRITE (stdout,20) outer, inner, &
3659 & trim(vname(1,idsvar(isfsur)))
3660 END IF
3661 END IF
3662 END IF
3663!
3664! 2D u-momentum component.
3665!
3666 IF ((statevar(outer).eq.0).or.(statevar(outer).eq.isubar)) THEN
3667 DO j=jstrb,jendb
3668 DO i=istrm,iendb
3669 tl_ubar(i,j,lout)=p*ad_ubar(i,j,linp)*scale
3670# ifdef MASKING
3671 tl_ubar(i,j,lout)=tl_ubar(i,j,lout)*umask(i,j)
3672# endif
3673 END DO
3674 END DO
3675 IF (domain(ng)%SouthWest_Test(tile)) THEN
3676 IF (master) THEN
3677 WRITE (stdout,20) outer, inner, &
3678 & trim(vname(1,idsvar(isubar)))
3679 END IF
3680 END IF
3681 END IF
3682!
3683! 2D v-momentum component.
3684!
3685 IF ((statevar(outer).eq.0).or.(statevar(outer).eq.isvbar)) THEN
3686 DO j=jstrm,jendb
3687 DO i=istrb,iendb
3688 tl_vbar(i,j,lout)=p*ad_vbar(i,j,linp)*scale
3689# ifdef MASKING
3690 tl_vbar(i,j,lout)=tl_vbar(i,j,lout)*vmask(i,j)
3691# endif
3692 END DO
3693 END DO
3694 IF (domain(ng)%SouthWest_Test(tile)) THEN
3695 IF (master) THEN
3696 WRITE (stdout,20) outer, inner, &
3697 & trim(vname(1,idsvar(isvbar)))
3698 END IF
3699 END IF
3700 END IF
3701# ifdef SOLVE3D
3702!
3703! 3D u-momentum component.
3704!
3705 IF ((statevar(outer).eq.0).or.(statevar(outer).eq.isuvel)) THEN
3706 DO k=1,n(ng)
3707 DO j=jstrb,jendb
3708 DO i=istrm,iendb
3709 tl_u(i,j,k,lout)=p*ad_u(i,j,k,linp)*scale
3710# ifdef MASKING
3711 tl_u(i,j,k,lout)=tl_u(i,j,k,lout)*umask(i,j)
3712# endif
3713 END DO
3714 END DO
3715 END DO
3716 IF (domain(ng)%SouthWest_Test(tile)) THEN
3717 IF (master) THEN
3718 WRITE (stdout,20) outer, inner, &
3719 & trim(vname(1,idsvar(isuvel)))
3720 END IF
3721 END IF
3722 END IF
3723!
3724! 3D v-momentum component.
3725!
3726 IF ((statevar(outer).eq.0).or.(statevar(outer).eq.isvvel)) THEN
3727 DO k=1,n(ng)
3728 DO j=jstrm,jendb
3729 DO i=istrb,iendb
3730 tl_v(i,j,k,lout)=p*ad_v(i,j,k,linp)*scale
3731# ifdef MASKING
3732 tl_v(i,j,k,lout)=tl_v(i,j,k,lout)*vmask(i,j)
3733# endif
3734 END DO
3735 END DO
3736 END DO
3737 IF (domain(ng)%SouthWest_Test(tile)) THEN
3738 IF (master) THEN
3739 WRITE (stdout,20) outer, inner, &
3740 & trim(vname(1,idsvar(isvvel)))
3741 END IF
3742 END IF
3743 END IF
3744!
3745! Tracers.
3746!
3747 DO itrc=1,nt(ng)
3748 IF ((statevar(outer).eq.0).or. &
3749 & (statevar(outer).eq.istvar(itrc))) THEN
3750 DO k=1,n(ng)
3751 DO j=jstrb,jendb
3752 DO i=istrb,iendb
3753 tl_t(i,j,k,lout,itrc)=p*ad_t(i,j,k,linp,itrc)*scale
3754# ifdef MASKING
3755 tl_t(i,j,k,lout,itrc)=tl_t(i,j,k,lout,itrc)*rmask(i,j)
3756# endif
3757 END DO
3758 END DO
3759 END DO
3760 IF (domain(ng)%SouthWest_Test(tile)) THEN
3761 IF (master) THEN
3762 WRITE (stdout,20) outer, inner, &
3763 & trim(vname(1,idsvar(istvar(itrc))))
3764 END IF
3765 END IF
3766 END IF
3767 END DO
3768# endif
3769 IF (master) WRITE (stdout,'(/)')
3770!
3771!-----------------------------------------------------------------------
3772! Apply lateral boundary conditions to 2D fields.
3773!-----------------------------------------------------------------------
3774!
3775 CALL tl_zetabc_tile (ng, tile, &
3776 & lbi, ubi, lbj, ubj, &
3777 & imins, imaxs, jmins, jmaxs, &
3778 & lout, lout, lout, &
3779 & zeta, tl_zeta)
3780# ifndef SOLVE3D
3781 CALL tl_u2dbc_tile (ng, tile, &
3782 & lbi, ubi, lbj, ubj, &
3783 & imins, imaxs, jmins, jmaxs, &
3784 & lout, lout, lout, &
3785 & ubar, vbar, zeta, &
3786 & tl_ubar, tl_vbar, tl_zeta)
3787 CALL tl_v2dbc_tile (ng, tile, &
3788 & lbi, ubi, lbj, ubj, &
3789 & imins, imaxs, jmins, jmaxs, &
3790 & lout, lout, lout, &
3791 & ubar, vbar, zeta, &
3792 & tl_ubar, tl_vbar, tl_zeta)
3793# endif
3794!
3795 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
3796 CALL exchange_r2d_tile (ng, tile, &
3797 & lbi, ubi, lbj, ubj, &
3798 & tl_zeta(:,:,lout))
3799# ifndef SOLVE3D
3800 CALL exchange_u2d_tile (ng, tile, &
3801 & lbi, ubi, lbj, ubj, &
3802 & tl_ubar(:,:,lout))
3803 CALL exchange_v2d_tile (ng, tile, &
3804 & lbi, ubi, lbj, ubj, &
3805 & tl_vbar(:,:,lout))
3806# endif
3807 END IF
3808
3809# ifdef DISTRIBUTE
3810!
3811 CALL mp_exchange2d (ng, tile, itlm, 1, &
3812 & lbi, ubi, lbj, ubj, &
3813 & nghostpoints, &
3814 & ewperiodic(ng), nsperiodic(ng), &
3815 & tl_zeta(:,:,lout))
3816# ifndef SOLVE3D
3817 CALL mp_exchange2d (ng, tile, itlm, 2, &
3818 & lbi, ubi, lbj, ubj, &
3819 & nghostpoints, &
3820 & ewperiodic(ng), nsperiodic(ng), &
3821 & tl_ubar(:,:,lout), &
3822 & tl_vbar(:,:,lout))
3823# endif
3824# endif
3825# ifdef SOLVE3D
3826!
3827!-----------------------------------------------------------------------
3828! Compute new depths and thicknesses.
3829!-----------------------------------------------------------------------
3830!
3831 CALL tl_set_depth_tile (ng, tile, itlm, &
3832 & lbi, ubi, lbj, ubj, &
3833 & imins, imaxs, jmins, jmaxs, &
3834 & nstp, nnew, &
3835 & h, tl_h, &
3836# ifdef ICESHELF
3837 & zice, &
3838# endif
3839# if defined SEDIMENT && defined SED_MORPH
3840 & tl_bed_thick, &
3841# endif
3842 & zt_avg1, tl_zt_avg1, &
3843 & tl_hz, tl_z_r, tl_z_w)
3844!
3845!-----------------------------------------------------------------------
3846! Apply lateral boundary conditions to perturbed 3D fields.
3847!-----------------------------------------------------------------------
3848!
3849 CALL tl_u3dbc_tile (ng, tile, &
3850 & lbi, ubi, lbj, ubj, n(ng), &
3851 & imins, imaxs, jmins, jmaxs, &
3852 & lout, lout, tl_u)
3853 CALL tl_v3dbc_tile (ng, tile, &
3854 & lbi, ubi, lbj, ubj, n(ng), &
3855 & imins, imaxs, jmins, jmaxs, &
3856 & lout, lout, tl_v)
3857!
3858 ic=0
3859 DO itrc=1,nt(ng)
3860 IF (ltracerclm(itrc,ng).and.lnudgetclm(itrc,ng)) THEN
3861 ic=ic+1
3862 END IF
3863 CALL tl_t3dbc_tile (ng, tile, itrc, ic, &
3864 & lbi, ubi, lbj, ubj, n(ng), nt(ng), &
3865 & imins, imaxs, jmins, jmaxs, &
3866 & lout, lout, tl_t)
3867 END DO
3868!
3869 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
3870 CALL exchange_u3d_tile (ng, tile, &
3871 & lbi, ubi, lbj, ubj, 1, n(ng), &
3872 & tl_u(:,:,:,lout))
3873 CALL exchange_v3d_tile (ng, tile, &
3874 & lbi, ubi, lbj, ubj, 1, n(ng), &
3875 & tl_v(:,:,:,lout))
3876 DO itrc=1,nt(ng)
3877 CALL exchange_r3d_tile (ng, tile, &
3878 & lbi, ubi, lbj, ubj, 1, n(ng), &
3879 & tl_t(:,:,:,lout,itrc))
3880 END DO
3881 END IF
3882
3883# ifdef DISTRIBUTE
3884!
3885 CALL mp_exchange3d (ng, tile, inlm, 2, &
3886 & lbi, ubi, lbj, ubj, 1, n(ng), &
3887 & nghostpoints, &
3888 & ewperiodic(ng), nsperiodic(ng), &
3889 & tl_u(:,:,:,lout), &
3890 & tl_v(:,:,:,lout))
3891 CALL mp_exchange4d (ng, tile, inlm, 1, &
3892 & lbi, ubi, lbj, ubj, 1, n(ng), 1, nt(ng), &
3893 & nghostpoints, &
3894 & ewperiodic(ng), nsperiodic(ng), &
3895 & tl_t(:,:,:,lout,:))
3896# endif
3897# endif
3898!
3899 10 FORMAT (/,' Perturbing Tangent Linear Initial Conditions:',/)
3900 20 FORMAT (' (Outer,Inner) = ','(',i4.4,',',i4.4,')',3x, &
3901 & 'Perturbing State Variable: ',a)
3902!
3903 RETURN
subroutine, public tl_set_depth_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nstp, nnew, h, tl_h, zice, zt_avg1, tl_zt_avg1, tl_hz, tl_z_r, tl_z_w)
subroutine, public tl_t3dbc_tile(ng, tile, itrc, ic, lbi, ubi, lbj, ubj, ubk, ubt, imins, imaxs, jmins, jmaxs, nstp, nout, tl_t)
Definition tl_t3dbc_im.F:58
subroutine, public tl_u2dbc_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, krhs, kstp, kout, ubar, vbar, zeta, tl_ubar, tl_vbar, tl_zeta)
Definition tl_u2dbc_im.F:64
subroutine, public tl_u3dbc_tile(ng, tile, lbi, ubi, lbj, ubj, ubk, imins, imaxs, jmins, jmaxs, nstp, nout, tl_u)
Definition tl_u3dbc_im.F:57
subroutine, public tl_v2dbc_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, krhs, kstp, kout, ubar, vbar, zeta, tl_ubar, tl_vbar, tl_zeta)
Definition tl_v2dbc_im.F:64
subroutine, public tl_v3dbc_tile(ng, tile, lbi, ubi, lbj, ubj, ubk, imins, imaxs, jmins, jmaxs, nstp, nout, tl_v)
Definition tl_v3dbc_im.F:57
subroutine, public tl_zetabc_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, krhs, kstp, kout, zeta, tl_zeta)
Definition tl_zetabc.F:58

References mod_fourdvar::addotproduct, mod_param::domain, mod_scalars::ewperiodic, exchange_2d_mod::exchange_r2d_tile(), exchange_3d_mod::exchange_r3d_tile(), exchange_2d_mod::exchange_u2d_tile(), exchange_3d_mod::exchange_u3d_tile(), exchange_2d_mod::exchange_v2d_tile(), exchange_3d_mod::exchange_v3d_tile(), mod_ncparam::idsvar, mod_param::inlm, mod_scalars::inner, mod_ncparam::isfsur, mod_ncparam::istvar, mod_ncparam::isubar, mod_ncparam::isuvel, mod_ncparam::isvbar, mod_ncparam::isvvel, mod_param::itlm, mod_scalars::lnudgetclm, mod_scalars::ltracerclm, mod_parallel::master, mp_exchange_mod::mp_exchange2d(), mp_exchange_mod::mp_exchange3d(), mp_exchange_mod::mp_exchange4d(), mod_param::nghostpoints, mod_scalars::nsperiodic, mod_scalars::outer, mod_iounits::stdout, tl_set_depth_mod::tl_set_depth_tile(), tl_t3dbc_mod::tl_t3dbc_tile(), tl_u2dbc_mod::tl_u2dbc_tile(), tl_u3dbc_mod::tl_u3dbc_tile(), tl_v2dbc_mod::tl_v2dbc_tile(), tl_v3dbc_mod::tl_v3dbc_tile(), tl_zetabc_mod::tl_zetabc_tile(), and mod_ncparam::vname.

Referenced by tl_ini_perturb().

Here is the call graph for this function:
Here is the caller graph for this function: