The TELEMAC-MASCARET system  trunk
fv_zones.F90
Go to the documentation of this file.
1 !###############################################################################
2 !# #
3 !# fv_zones.F90 #
4 !# #
5 !# Interface for FV (Finite Volume) Model to AED2 modules. #
6 !# Designed for TUFLOW-FV, released by BMT-WBM: #
7 !# http://www.tuflow.com/Tuflow%20FV.aspx #
8 !# #
9 !# This is a support module to allow ability for benthic/sediment zones in #
10 !# AED2, including zone-averaging #
11 !# #
12 !# ----------------------------------------------------------------------- #
13 !# #
14 !# Developed by : #
15 !# AquaticEcoDynamics (AED) Group #
16 !# School of Earth & Environment #
17 !# (C) The University of Western Australia #
18 !# #
19 !# Copyright by the AED-team @ UWA under the GNU Public License - www.gnu.org #
20 !# #
21 !# ----------------------------------------------------------------------- #
22 !# #
23 !# Created Apr 2015 #
24 !# #
25 !###############################################################################
26 
27 #if defined HAVE_AED2
28 #include "aed2.h"
29 !#include "aed2_common.h"
30 #endif
31 
32 #ifndef DEBUG
33 #define DEBUG 0
34 #endif
35 
36 !###############################################################################
37 MODULE fv_zones
38 !-------------------------------------------------------------------------------
39 #if defined HAVE_AED2
40  USE aed2_common
42 
43  IMPLICIT NONE
44 
46  PUBLIC compute_zone_benthic_fluxes !, distribute_flux_from_zone
47  PUBLIC zone
48 
49  !#--------------------------------------------------------------------------#
50  !# Module Data
51 
52  !# Arrays for environmental variables not supplied externally.
53  aed_real,DIMENSION(:), ALLOCATABLE,TARGET :: zone
54  aed_real,DIMENSION(:,:),ALLOCATABLE,TARGET :: zone_cc
55  aed_real,DIMENSION(:,:),ALLOCATABLE,TARGET :: zone_cc_diag
56  aed_real,DIMENSION(:), ALLOCATABLE,TARGET :: zone_area
57  aed_real,DIMENSION(:), ALLOCATABLE,TARGET :: zone_temp
58  aed_real,DIMENSION(:), ALLOCATABLE,TARGET :: zone_salt
59  aed_real,DIMENSION(:), ALLOCATABLE,TARGET :: zone_rho
60  aed_real,DIMENSION(:), ALLOCATABLE,TARGET :: zone_height
61  aed_real,DIMENSION(:), ALLOCATABLE,TARGET :: zone_extc
62  aed_real,DIMENSION(:), ALLOCATABLE,TARGET :: zone_tss
63  aed_real,DIMENSION(:), ALLOCATABLE,TARGET :: zone_par
64  aed_real,DIMENSION(:), ALLOCATABLE,TARGET :: zone_wind
65  aed_real,DIMENSION(:), ALLOCATABLE,TARGET :: zone_rain ! JC not sure it makes sense
66  aed_real,DIMENSION(:), ALLOCATABLE,TARGET :: zone_i_0
67  aed_real,DIMENSION(:), ALLOCATABLE,TARGET :: zone_taub
68  INTEGER, DIMENSION(:), ALLOCATABLE :: zone_count, zm
69 
70  aed_real,DIMENSION(:,:),ALLOCATABLE,TARGET :: flux_pelz
71  aed_real,DIMENSION(:,:),ALLOCATABLE,TARGET :: flux_benz
72 
73  INTEGER :: nzones, nwq_var, nben_var
74 
75 !#####################################################
76 
77 CONTAINS
78 !===============================================================================
79 
80 
81 !###############################################################################
82 SUBROUTINE init_zones(nCols, mat_id, n_aed2_vars, n_vars, n_vars_ben)
83 !-------------------------------------------------------------------------------
84 !ARGUMENTS
85  INTEGER,INTENT(in) :: nCols
86  INTEGER,DIMENSION(:,:),INTENT(in) :: mat_id
87  INTEGER,INTENT(in) :: n_aed2_vars, n_vars, n_vars_ben
88 !
89 !LOCALS
90  INTEGER :: cType, nTypes
91  INTEGER :: col, zon
92  INTEGER,DIMENSION(:),ALLOCATABLE :: mat_t
93 !
94 !-------------------------------------------------------------------------------
95 !BEGIN
96  ALLOCATE(mat_t(ncols)) ; ALLOCATE(zm(ncols))
97  !# The new form of zones
98  ctype = mat_id(1,1) ; ntypes = 1 ; mat_t(ntypes) = mat_id(1,1)
99  DO col=1, ubound(mat_id,2)
100  !# use the bottom index to fill the array
101  IF ( ctype /= mat_id(1,col) ) THEN
102  DO zon=1,ntypes
103  IF ( mat_t(zon) .eq. mat_id(1,col) ) THEN
104  ctype = mat_id(1,col)
105  EXIT
106  ENDIF
107  ENDDO
108  ENDIF
109  IF ( ctype /= mat_id(1,col) ) THEN
110  ntypes = ntypes + 1
111  mat_t(ntypes) = mat_id(1,col)
112  ctype = mat_id(1,col)
113  zon = ntypes
114  ENDIF
115  zm(col) = zon
116  ENDDO
117  WRITE(lu, *) "Number of mats = ", ntypes, " = ", mat_t(1:ntypes)
118  DEALLOCATE(mat_t)
119  nzones = ntypes
120 
121  ALLOCATE(zone(nzones))
122  ALLOCATE(zone_area(nzones))
123  ALLOCATE(zone_temp(nzones))
124  ALLOCATE(zone_salt(nzones))
125  ALLOCATE(zone_rho(nzones))
126  ALLOCATE(zone_height(nzones))
127 
128  ALLOCATE(zone_extc(nzones))
129  ALLOCATE(zone_tss(nzones))
130  ALLOCATE(zone_par(nzones))
131  ALLOCATE(zone_wind(nzones))
132  ALLOCATE(zone_rain(nzones)) ! JC
133  ALLOCATE(zone_i_0(nzones))
134  ALLOCATE(zone_taub(nzones))
135 
136  ALLOCATE(zone_count(nzones))
137 
138  ALLOCATE(zone_cc(n_vars, nzones))
139  ALLOCATE(zone_cc_diag(n_aed2_vars-n_vars, nzones))
140 
141  ALLOCATE(flux_pelz(n_vars, nzones))
142  ALLOCATE(flux_benz(n_vars, nzones))
143 
144  nwq_var = n_vars
145  nben_var = n_vars_ben
146 END SUBROUTINE init_zones
147 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
148 
149 
150 !###############################################################################
151 SUBROUTINE calc_zone_areas(nCols, temp, salt, h, area, wnd, rho, extcoeff, I_0, par, tss, active, rain)!, bathy ,rain JC
152 !-------------------------------------------------------------------------------
153 !ARGUMENTS
154  INTEGER,INTENT(in) :: nCols
155  aed_real,DIMENSION(:),INTENT(in) :: temp, salt, h, area, wnd, rain
156  aed_real,DIMENSION(:),INTENT(in) :: rho, extcoeff, i_0, par, tss
157  LOGICAL,DIMENSION(:), INTENT(in) :: active
158 !
159 !LOCALS
160  INTEGER :: col, zon
161 !
162 !-------------------------------------------------------------------------------
163 !BEGIN
164  zone_area = zero_
165  zone_temp = zero_
166  zone_salt = zero_
167  zone_height = zero_
168  zone_wind = zero_
169  zone_rain = zero_ ! JC
170  zone_count = 0
171 
172  DO col=1, ncols
173  IF (.NOT. active(col)) cycle
174 
175  zon = zm(col)
176  zone_area(zon) = zone_area(zon) + area(col)
177  zone_temp(zon) = zone_temp(zon) + temp(col)
178  zone_salt(zon) = zone_salt(zon) + salt(col)
179  zone_rho(zon) = zone_rho(zon) + rho(col)
180  zone_height(zon) = zone_height(zon) + h(col)
181  zone_extc(zon) = zone_extc(zon) + extcoeff(col)
182  zone_tss(zon) = zone_tss(zon) + tss(col)
183  zone_par(zon) = zone_par(zon) + par(col)
184  zone_wind(zon) = zone_wind(zon) + wnd(col)
185  zone_rain(zon) = zone_rain(zon) + rain(col) !JC
186  zone_i_0(zon) = zone_i_0(zon) + i_0(col)
187  ! zone_taub(zon) = zone_taub(zon) + col_taub
188 
189  zone_count(zon) = zone_count(zon) + 1
190  ENDDO
191 
192  zone_temp = zone_temp / zone_count
193  zone_salt = zone_salt / zone_count
194 
195  zone_rho = zone_rho / zone_count
196  zone_extc = zone_extc / zone_count
197  zone_tss = zone_tss / zone_count
198  zone_par = zone_par / zone_count
199  zone_wind = zone_wind / zone_count
200  zone_rain = zone_rain / zone_count ! JC
201  zone_i_0 = zone_i_0 / zone_count
202  zone_taub = zone_taub / zone_count
203 
204 END SUBROUTINE calc_zone_areas
205 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
206 
207 
208 !###############################################################################
209 SUBROUTINE copy_to_zone(nCols, cc, area, active, benth_map)
210 !-------------------------------------------------------------------------------
211 !ARGUMENTS
212  INTEGER,INTENT(in) :: nCols
213  aed_real,INTENT(in) :: cc(:,:) !# (n_vars, n_layers)
214  aed_real,DIMENSION(:),INTENT(in) :: area
215  LOGICAL,DIMENSION(:), INTENT(in) :: active
216  INTEGER,DIMENSION(:), INTENT(in) :: benth_map
217 !
218 !LOCALS
219  INTEGER :: col, zon, bot
220 !
221 !-------------------------------------------------------------------------------
222 !BEGIN
223  DO col=1, ncols
224  IF (.NOT. active(col)) cycle
225 
226  bot = benth_map(col)
227  zon = zm(col)
228 
229  zone_cc(1:nwq_var,zon) = zone_cc(1:nwq_var,zon) + (cc(1:nwq_var,bot) * area(col) / zone_area(zon))
230  ENDDO
231 END SUBROUTINE copy_to_zone
232 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
233 
234 
235 !###############################################################################
236 SUBROUTINE copy_from_zone(nCols, cc, area, active, benth_map)
237 !-------------------------------------------------------------------------------
238 !ARGUMENTS
239  INTEGER,INTENT(in) :: nCols
240  aed_real,INTENT(inout) :: cc(:,:) !# (n_vars, n_layers)
241  aed_real,DIMENSION(:),INTENT(in) :: area
242  LOGICAL,DIMENSION(:), INTENT(in) :: active
243  INTEGER,DIMENSION(:), INTENT(in) :: benth_map
244 !
245 !LOCALS
246  INTEGER :: col, zon, bot
247 !
248 !-------------------------------------------------------------------------------
249 !BEGIN
250  DO col=1, ncols
251  IF (.NOT. active(col)) cycle
252 
253  bot = benth_map(col)
254  zon = zm(col)
255 
256  cc(nwq_var+1:nwq_var+nben_var,bot) = zone_cc(nwq_var+1:nwq_var+nben_var,zon)
257  ENDDO
258 END SUBROUTINE copy_from_zone
259 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
260 
261 
262 !###############################################################################
263 SUBROUTINE stopit(message)
264 !-------------------------------------------------------------------------------
265 !ARGUMENTS
266  CHARACTER(*) :: message
267 !-------------------------------------------------------------------------------
268  print *,message
269  stop
270 END SUBROUTINE stopit
271 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
272 
273 
274 !###############################################################################
275 SUBROUTINE define_column_zone(column, zon, n_aed2_vars)
276 !-------------------------------------------------------------------------------
277 !ARGUMENTS
278  type(aed2_column_t), INTENT(inout) :: column(:)
279  INTEGER, INTENT(in) :: zon, n_aed2_vars
280 !
281 !LOCALS
282  INTEGER :: av
283  INTEGER :: v, d, sv, sd, ev
284  TYPE(aed2_variable_t),POINTER :: tvar
285 !
286 !-------------------------------------------------------------------------------
287 !BEGIN
288  v = 0 ; d = 0; sv = 0; sd = 0 ; ev = 0
289  DO av=1,n_aed2_vars
290 
291  IF ( .NOT. aed2_get_var(av, tvar) ) stop "Error getting variable info"
292 
293  IF ( tvar%extern ) THEN !# global variable
294  ev = ev + 1
295  SELECT CASE (tvar%name)
296  CASE ( 'temperature' ) ; column(av)%cell => zone_temp
297  CASE ( 'salinity' ) ; column(av)%cell => zone_salt
298  CASE ( 'density' ) ; column(av)%cell => zone_rho
299  CASE ( 'layer_ht' ) ; column(av)%cell => zone_height
300  CASE ( 'layer_area' ) ; column(av)%cell_sheet => zone_area(zon)
301  CASE ( 'rain' ) ; column(av)%cell_sheet => zone_rain(zon) ! JC
302  ! CASE ( 'rainloss' ) ; column(av)%cell_sheet => zone_rainloss(zon) ! JC
303  CASE ( 'material' ) ; column(av)%cell_sheet => zone(zon)
304  ! CASE ( 'bathy' ) ; column(av)%cell_sheet => zone_bathy(zon)! Does it mean this var is already pointed to the correct variable of TFFV? JC
305  CASE ( 'extc_coef' ) ; column(av)%cell => zone_extc
306  CASE ( 'tss' ) ; column(av)%cell => zone_tss
307  ! CASE ( 'par' ) ; column(av)%cell => zone_par
308  ! CASE ( 'par' ) ; IF (link_ext_par) THEN
309  ! column(av)%cell => lpar(top:bot)
310  ! ELSE
311  ! column(av)%cell => par(top:bot)
312  ! ENDIF
313  ! CASE ( 'nir' ) ; column(av)%cell => zone_nir
314  ! CASE ( 'uva' ) ; column(av)%cell => zone_uva
315  ! CASE ( 'uvb' ) ; column(av)%cell => zone_uvb
316  CASE ( 'sed_zone' ) ; column(av)%cell_sheet => zone(zon)
317  CASE ( 'wind_speed' ) ; column(av)%cell_sheet => zone_wind(zon)
318  CASE ( 'par_sf' ) ; column(av)%cell_sheet => zone_i_0(zon)
319  ! CASE ( 'taub' ) ; column(av)%cell_sheet => zone_taub
320  ! CASE ( 'air_temp' ) ; column(av)%cell_sheet => zone_air_temp(zon)
321  CASE DEFAULT ; CALL stopit("ERROR: external variable "//trim(tvar%name)//" not found.")
322  END SELECT
323  ELSEIF ( tvar%diag ) THEN !# Diagnostic variable
324  d = d + 1
325  IF ( tvar%sheet ) THEN
326  column(av)%cell_sheet => zone_cc_diag(d,zon)
327  ELSE
328  column(av)%cell => zone_cc_diag(d,zon:zon)
329  ENDIF
330  ELSE !# state variable
331  IF ( tvar%sheet ) THEN
332  sv = sv + 1
333  IF ( tvar%bot ) THEN
334  column(av)%cell_sheet => zone_cc(nwq_var+sv, zon)
335  ELSEIF ( tvar%top ) THEN
336  ! column(av)%cell_sheet => zone_cc(nwq_var+sv, top)
337  ENDIF
338  column(av)%flux_ben => flux_benz(nwq_var+sv, zon)
339  ! column(av)%flux_atm => flux_atm(nwq_var+sv)
340  ELSE
341  v = v + 1
342  column(av)%cell => zone_cc(v, 1:1)
343  column(av)%flux_pel => flux_pelz(v, 1:1)
344  column(av)%flux_ben => flux_benz(v, zon)
345  ! column(av)%flux_atm => flux_atm(v)
346  ENDIF
347  ENDIF
348  ENDDO
349 END SUBROUTINE define_column_zone
350 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
351 
352 
353 !###############################################################################
354 SUBROUTINE compute_zone_benthic_fluxes(n_aed2_vars, dt)
355 !-------------------------------------------------------------------------------
356 !ARGUMENTS
357  INTEGER,INTENT(in) :: n_aed2_vars
358  aed_real,INTENT(in):: dt
359 !
360 !LOCALS
361  INTEGER :: zon, v
362  type(aed2_column_t) :: column(n_aed2_vars)
363 !
364 !-------------------------------------------------------------------------------
365 !BEGIN
366  DO zon=1, nzones
367  CALL define_column_zone(column, zon, n_aed2_vars)!, nwq_var)
368 
369  CALL aed2_calculate_benthic(column, 1)
370  DO v=nwq_var+1,nwq_var+nben_var
371  zone_cc(v, 1) = zone_cc(v, 1) + dt*flux_benz(v, zon);
372  ENDDO
373  ENDDO
374 END SUBROUTINE compute_zone_benthic_fluxes
375 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
376 
377 
378 #endif
379 !===============================================================================
380 END MODULE fv_zones
integer nzones
Definition: fv_zones.F90:73
integer, dimension(:), allocatable zm
Definition: fv_zones.F90:68
integer nwq_var
Definition: fv_zones.F90:73
integer nben_var
Definition: fv_zones.F90:73
integer, dimension(:), allocatable zone_count
Definition: fv_zones.F90:68