ROMS
Loading...
Searching...
No Matches
set_3dfldr.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#if defined ADJOINT && defined SOLVE3D
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This routine time-interpolates backwards in time requested 3D field !
13! from snapshots of input data. !
14! !
15!=======================================================================
16!
17 implicit none
18
19 CONTAINS
20!
21!***********************************************************************
22 SUBROUTINE set_3dfldr_tile (ng, tile, model, ifield, &
23 & LBi, UBi, LBj, UBj, LBk, UBk, &
24 & Finp, Fout, update, &
25 & SetBC)
26!***********************************************************************
27!
28 USE mod_param
29 USE mod_parallel
30 USE mod_iounits
31 USE mod_ncparam
32 USE mod_scalars
33!
35# ifdef DISTRIBUTE
37# endif
38 USE dateclock_mod, ONLY : time_string
39!
40! Imported variable declarations.
41!
42 logical, intent(in), optional :: SetBC
43
44 logical, intent(out) :: update
45
46 integer, intent(in) :: ng, tile, model, ifield
47 integer, intent(in) :: LBi, UBi, LBj, UBj, LBk, UBk
48!
49# ifdef ASSUMED_SHAPE
50 real(r8), intent(in) :: Finp(LBi:,LBj:,LBk:,:)
51 real(r8), intent(out) :: Fout(LBi:,LBj:,LBk:)
52# else
53 real(r8), intent(in) :: Finp(LBi:UBi,LBj:UBj,LBk:UBk,2)
54 real(r8), intent(out) :: Fout(LBi:UBi,LBj:UBj,LBk:UBk)
55# endif
56!
57! Local variable declarations.
58!
59 logical :: LapplyBC, Lgrided, Lonerec
60!
61 integer :: Tindex, gtype, i, it1, it2, j, k
62!
63 real(dp) :: SecScale, fac, fac1, fac2
64 real(r8) :: Fval
65!
66 character (len=22) :: DateString1, DateString2
67
68# include "set_bounds.h"
69!
70!--------------------------------------------------------------------
71! Set-up requested field for current tile.
72!--------------------------------------------------------------------
73!
74! Set switch to apply boundary conditions.
75!
76 IF (PRESENT(setbc)) THEN
77 lapplybc=setbc
78 ELSE
79 lapplybc=.true.
80 END IF
81!
82! Get requested field information from global storage.
83!
84 lgrided=linfo(1,ifield,ng)
85 lonerec=linfo(3,ifield,ng)
86 gtype =iinfo(1,ifield,ng)
87 tindex =iinfo(8,ifield,ng)
88 update=.true.
89!
90! Set linear, time interpolation factors. Fractional seconds are
91! rounded to the nearest milliseconds integer towards zero in the
92! time interpolation weights.
93!
94 secscale=1000.0_dp ! seconds to milliseconds
95 it1=3-tindex
96 it2=tindex
97 fac1=anint((time(ng)-tintrp(it2,ifield,ng))*secscale,r8)
98 fac2=anint((tintrp(it1,ifield,ng)-time(ng))*secscale,r8)
99!! CALL time_string (Tintrp(it1,ifield,ng), DateString1)
100!! CALL time_string (Tintrp(it2,ifield,ng), DateString2)
101!
102! Load time-invariant data. Time interpolation is not necessary.
103!
104 IF (lonerec) THEN
105 IF (lgrided) THEN
106 DO k=lbk,ubk
107 DO j=jstrr,jendr
108 DO i=istrr,iendr
109 fout(i,j,k)=finp(i,j,k,tindex)
110 END DO
111 END DO
112 END DO
113 ELSE
114 fval=fpoint(tindex,ifield,ng)
115 DO k=lbk,ubk
116 DO j=jstrr,jendr
117 DO i=istrr,iendr
118 fout(i,j,k)=fval
119 END DO
120 END DO
121 END DO
122 END IF
123!
124! Time-interpolate from gridded or point data.
125!
126 ELSE IF (((fac1*fac2).ge.0.0_dp).and. &
127 & ((fac1+fac2).gt.0.0_dp)) THEN
128 fac=1.0_dp/(fac1+fac2)
129 fac1=fac*fac1 ! nondimensional
130 fac2=fac*fac2 ! nondimensional
131 IF (lgrided) THEN
132 DO k=lbk,ubk
133 DO j=jstrr,jendr
134 DO i=istrr,iendr
135 fout(i,j,k)=fac1*finp(i,j,k,it1)+fac2*finp(i,j,k,it2)
136 END DO
137 END DO
138 END DO
139 ELSE
140 fval=fac1*fpoint(it1,ifield,ng)+fac2*fpoint(it2,ifield,ng)
141 DO k=lbk,ubk
142 DO j=jstrr,jendr
143 DO i=istrr,iendr
144 fout(i,j,k)=fval
145 END DO
146 END DO
147 END DO
148 END IF
149!
150! Activate synchronization flag if a new time record needs to be
151! read in at the next time step.
152!
153 IF ((time(ng)-dt(ng)).lt.tintrp(it2,ifield,ng)) THEN
154 IF (domain(ng)%SouthWest_Test(tile)) synchro_flag(ng)=.true.
155 END IF
156!
157! Unable to set-up requested field. Activate error flag to quit.
158!
159 ELSE
160 IF (domain(ng)%SouthWest_Test(tile)) THEN
161 IF (master) THEN
162 WRITE (stdout,10) trim(vname(1,ifield)), tdays(ng), &
163 & finfo(1,ifield,ng), finfo(2,ifield,ng), &
164 & finfo(3,ifield,ng), finfo(4,ifield,ng), &
165 & tintrp(it1,ifield,ng)*sec2day, &
166 & tintrp(it2,ifield,ng)*sec2day, &
167 & fac1*sec2day/secscale, &
168 & fac2*sec2day/secscale
169 END IF
170 10 FORMAT (/,' SET_3DFLDR - current model time', &
171 & ' exceeds ending value for variable: ',a, &
172 & /,14x,'TDAYS = ',f15.4, &
173 & /,14x,'Data Tmin = ',f15.4,2x,'Data Tmax = ',f15.4, &
174 & /,14x,'Data Tstr = ',f15.4,2x,'Data Tend = ',f15.4, &
175 & /,14x,'TINTRP1 = ',f15.4,2x,'TINTRP2 = ',f15.4, &
176 & /,14x,'FAC1 = ',f15.4,2x,'FAC2 = ',f15.4)
177 exit_flag=2
178 update=.false.
179 END IF
180 END IF
181!
182! Exchange boundary data.
183!
184 IF (update) THEN
185 IF (lapplybc.and.(ewperiodic(ng).or.nsperiodic(ng))) THEN
186 IF (gtype.eq.r3dvar) THEN
187 CALL exchange_r3d_tile (ng, tile, &
188 & lbi, ubi, lbj, ubj, lbk, ubk, &
189 & fout)
190 ELSE IF (gtype.eq.u3dvar) THEN
191 CALL exchange_u3d_tile (ng, tile, &
192 & lbi, ubi, lbj, ubj, lbk, ubk, &
193 & fout)
194 ELSE IF (gtype.eq.v3dvar) THEN
195 CALL exchange_v3d_tile (ng, tile, &
196 & lbi, ubi, lbj, ubj, lbk, ubk, &
197 & fout)
198 ELSE IF (gtype.eq.w3dvar) THEN
199 CALL exchange_w3d_tile (ng, tile, &
200 & lbi, ubi, lbj, ubj, lbk, ubk, &
201 & fout)
202 END IF
203 END IF
204
205# ifdef DISTRIBUTE
206 IF (.not.lapplybc) THEN
207 CALL mp_exchange3d (ng, tile, model, 1, &
208 & lbi, ubi, lbj, ubj, lbk, ubk, &
209 & nghostpoints, &
210 & .false., .false., &
211 & fout)
212 ELSE
213 CALL mp_exchange3d (ng, tile, model, 1, &
214 & lbi, ubi, lbj, ubj, lbk, ubk, &
215 & nghostpoints, &
216 & ewperiodic(ng), nsperiodic(ng), &
217 & fout)
218 END IF
219# endif
220 END IF
221
222 RETURN
223 END SUBROUTINE set_3dfldr_tile
224#endif
225 END MODULE set_3dfldr_mod
subroutine, public time_string(mytime, date_string)
Definition dateclock.F:1272
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)
subroutine exchange_w3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
integer stdout
logical, dimension(:,:,:), allocatable linfo
real(dp), dimension(:,:,:), allocatable tintrp
real(dp), dimension(:,:,:), allocatable fpoint
real(dp), dimension(:,:,:), allocatable finfo
character(len=maxlen), dimension(6, 0:nv) vname
integer, dimension(:,:,:), allocatable iinfo
logical master
integer, parameter r3dvar
Definition mod_param.F:721
integer nghostpoints
Definition mod_param.F:710
integer, parameter u3dvar
Definition mod_param.F:722
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, parameter w3dvar
Definition mod_param.F:724
integer, parameter v3dvar
Definition mod_param.F:723
real(dp), dimension(:), allocatable dt
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
logical, dimension(:), allocatable synchro_flag
real(dp), dimension(:), allocatable tdays
real(dp), parameter sec2day
integer exit_flag
real(dp), dimension(:), allocatable time
subroutine mp_exchange3d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine set_3dfldr_tile(ng, tile, model, ifield, lbi, ubi, lbj, ubj, lbk, ubk, finp, fout, update, setbc)
Definition set_3dfldr.F:26