host.F90 63.6 KB
Newer Older
Jorn Bruggeman's avatar
Jorn Bruggeman committed
1 2 3
#include "fabm_driver.h"
#include "fabm_private.h"

Jorn Bruggeman's avatar
Jorn Bruggeman committed
4 5 6 7 8
#undef _BEGIN_OUTER_INTERIOR_LOOP_
#undef _END_OUTER_INTERIOR_LOOP_
#undef _BEGIN_OUTER_HORIZONTAL_LOOP_
#undef _END_OUTER_HORIZONTAL_LOOP_

9 10
module host_hooks
   use fabm_driver
Jorn Bruggeman's avatar
Jorn Bruggeman committed
11

12
   implicit none
Jorn Bruggeman's avatar
Jorn Bruggeman committed
13

14
   type,extends(type_base_driver) :: type_test_driver
Jorn Bruggeman's avatar
Jorn Bruggeman committed
15
   contains
16 17 18
      procedure :: fatal_error => test_driver_fatal_error
      procedure :: log_message => test_driver_log_message
   end type
Jorn Bruggeman's avatar
Jorn Bruggeman committed
19

20
contains
Jorn Bruggeman's avatar
Jorn Bruggeman committed
21

22 23 24
   subroutine test_driver_fatal_error(self,location,message)
      class (type_test_driver), intent(inout) :: self
      character(len=*),         intent(in)    :: location,message
Jorn Bruggeman's avatar
Jorn Bruggeman committed
25

26 27 28
      write (*,*) trim(location)//': '//trim(message)
      stop 1
   end subroutine
Jorn Bruggeman's avatar
Jorn Bruggeman committed
29

30 31 32
   subroutine test_driver_log_message(self,message)
      class (type_test_driver), intent(inout) :: self
      character(len=*),         intent(in)    :: message
Jorn Bruggeman's avatar
Jorn Bruggeman committed
33

34 35
      write (*,*) trim(message)
   end subroutine
Jorn Bruggeman's avatar
Jorn Bruggeman committed
36

37
end module host_hooks
38

Jorn Bruggeman's avatar
Jorn Bruggeman committed
39 40
program test_host

41
   use fabm
Jorn Bruggeman's avatar
Jorn Bruggeman committed
42 43 44
   use fabm_driver
   use fabm_parameters, only: rke
   use fabm_types, only: source_do, source_do_surface, source_do_bottom, source_do_column
Jorn Bruggeman's avatar
Jorn Bruggeman committed
45

Jorn Bruggeman's avatar
Jorn Bruggeman committed
46 47
   use test_models
   use host_hooks
Jorn Bruggeman's avatar
Jorn Bruggeman committed
48

Jorn Bruggeman's avatar
Jorn Bruggeman committed
49
   implicit none
Jorn Bruggeman's avatar
Jorn Bruggeman committed
50 51

#if _FABM_DIMENSION_COUNT_>0
Jorn Bruggeman's avatar
Jorn Bruggeman committed
52
   integer :: _LOCATION_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
53 54
#endif

55 56
#ifdef _FABM_VECTORIZED_DIMENSION_INDEX_
#  define _INTERIOR_SLICE_RANGE_PLUS_1_ (_START_:_STOP_,:)
57
#else
58 59 60 61 62 63 64
#  define _INTERIOR_SLICE_RANGE_PLUS_1_ (:)
#endif

#if defined(_FABM_VECTORIZED_DIMENSION_INDEX_) &&  _FABM_VECTORIZED_DIMENSION_INDEX_!=_FABM_DEPTH_DIMENSION_INDEX_
#  define _HORIZONTAL_SLICE_RANGE_PLUS_1_ _INTERIOR_SLICE_RANGE_PLUS_1_
#else
#  define _HORIZONTAL_SLICE_RANGE_PLUS_1_ (:)
65 66
#endif

67 68 69 70
#if _FABM_DIMENSION_COUNT_==0
#  define _BEGIN_GLOBAL_LOOP_
#  define _END_GLOBAL_LOOP_
#elif _FABM_DIMENSION_COUNT_==1
71
#  define _BEGIN_GLOBAL_LOOP_ do i__=domain_start(1),domain_stop(1)
72 73
#  define _END_GLOBAL_LOOP_ end do;i__=domain_extent(1)
#  ifdef _FABM_DEPTH_DIMENSION_INDEX_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
74
#    define _BEGIN_GLOBAL_HORIZONTAL_LOOP_
75
#    define _END_GLOBAL_HORIZONTAL_LOOP_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
76
#  endif
77
#elif _FABM_DIMENSION_COUNT_==2
78
#  define _BEGIN_GLOBAL_LOOP_ do j__=domain_start(2),domain_stop(2);do i__=domain_start(1),domain_stop(1)
79 80
#  define _END_GLOBAL_LOOP_ end do;end do;i__=domain_extent(1);j__=domain_extent(2)
#  if _FABM_DEPTH_DIMENSION_INDEX_==1
81
#    define _BEGIN_GLOBAL_HORIZONTAL_LOOP_ do j__=domain_start(2),domain_stop(2)
82 83
#    define _END_GLOBAL_HORIZONTAL_LOOP_ end do;j__=domain_extent(2)
#  elif _FABM_DEPTH_DIMENSION_INDEX_==2
84
#    define _BEGIN_GLOBAL_HORIZONTAL_LOOP_ do i__=domain_start(1),domain_stop(1)
85 86 87
#    define _END_GLOBAL_HORIZONTAL_LOOP_ end do;i__=domain_extent(1)
#  endif
#elif _FABM_DIMENSION_COUNT_==3
88
#  define _BEGIN_GLOBAL_LOOP_ do k__=domain_start(3),domain_stop(3);do j__=domain_start(2),domain_stop(2);do i__=domain_start(1),domain_stop(1)
89 90
#  define _END_GLOBAL_LOOP_ end do;end do;end do;i__=domain_extent(1);j__=domain_extent(2);k__=domain_extent(3)
#  if _FABM_DEPTH_DIMENSION_INDEX_==1
91
#    define _BEGIN_GLOBAL_HORIZONTAL_LOOP_ do k__=domain_start(3),domain_stop(3);do j__=domain_start(2),domain_stop(2)
92 93
#    define _END_GLOBAL_HORIZONTAL_LOOP_ end do;end do;j__=domain_extent(2);k__=domain_extent(3)
#  elif _FABM_DEPTH_DIMENSION_INDEX_==2
94
#    define _BEGIN_GLOBAL_HORIZONTAL_LOOP_ do k__=domain_start(3),domain_stop(3);do i__=domain_start(1),domain_stop(1)
95 96
#    define _END_GLOBAL_HORIZONTAL_LOOP_ end do;end do;i__=domain_extent(1);k__=domain_extent(3)
#  elif _FABM_DEPTH_DIMENSION_INDEX_==3
97
#    define _BEGIN_GLOBAL_HORIZONTAL_LOOP_ do j__=domain_start(2),domain_stop(2);do i__=domain_start(1),domain_stop(1)
98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116
#    define _END_GLOBAL_HORIZONTAL_LOOP_ end do;end do;i__=domain_extent(1);j__=domain_extent(2)
#  endif
#endif

! If there is no depth dimension, horizontal = global
#ifndef _FABM_DEPTH_DIMENSION_INDEX_
#  define _BEGIN_GLOBAL_HORIZONTAL_LOOP_ _BEGIN_GLOBAL_LOOP_
#  define _END_GLOBAL_HORIZONTAL_LOOP_ _END_GLOBAL_LOOP_
#endif

#ifndef _FABM_VECTORIZED_DIMENSION_INDEX_
   ! No vectorization: outer loops are global loops
#  define _BEGIN_OUTER_INTERIOR_LOOP_ _BEGIN_GLOBAL_LOOP_
#  define _END_OUTER_INTERIOR_LOOP_ _END_GLOBAL_LOOP_
#  define _BEGIN_OUTER_HORIZONTAL_LOOP_ _BEGIN_GLOBAL_HORIZONTAL_LOOP_
#  define _END_OUTER_HORIZONTAL_LOOP_ _END_GLOBAL_HORIZONTAL_LOOP_
#elif _FABM_DIMENSION_COUNT_==1
   ! Entire domain is vectorized; no outer loops needed
#  define _BEGIN_OUTER_INTERIOR_LOOP_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
117
#  define _END_OUTER_INTERIOR_LOOP_
118 119
#  define _BEGIN_OUTER_HORIZONTAL_LOOP_
#  define _END_OUTER_HORIZONTAL_LOOP_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
120
#elif _FABM_DIMENSION_COUNT_==2
121
#  define _BEGIN_OUTER_INTERIOR_LOOP_ do j__=domain_start(2),domain_stop(2)
122
#  define _END_OUTER_INTERIOR_LOOP_ end do;j__=domain_extent(2)
123
#  if _FABM_DEPTH_DIMENSION_INDEX_==2
124
     ! The entire horizontal is already vectorized; no outer loop necessary
125
#    define _BEGIN_OUTER_HORIZONTAL_LOOP_
126 127 128 129 130
#    define _END_OUTER_HORIZONTAL_LOOP_
#  else
     ! No horizontal dimension vectorized; do full outer loop.
#    define _BEGIN_OUTER_HORIZONTAL_LOOP_ _BEGIN_GLOBAL_HORIZONTAL_LOOP_
#    define _END_OUTER_HORIZONTAL_LOOP_ _END_GLOBAL_HORIZONTAL_LOOP_
131 132
#  endif
#elif _FABM_DIMENSION_COUNT_==3
133
#  define _BEGIN_OUTER_INTERIOR_LOOP_ do k__=domain_start(3),domain_stop(3);do j__=domain_start(2),domain_stop(2)
134 135
#  define _END_OUTER_INTERIOR_LOOP_ end do;end do;j__=domain_extent(2);k__=domain_extent(3)
#  if _FABM_DEPTH_DIMENSION_INDEX_==2
136
#    define _BEGIN_OUTER_HORIZONTAL_LOOP_ do k__=domain_start(3),domain_stop(3)
137 138
#    define _END_OUTER_HORIZONTAL_LOOP_ end do;k__=domain_extent(3)
#  elif _FABM_DEPTH_DIMENSION_INDEX_==3
139
#    define _BEGIN_OUTER_HORIZONTAL_LOOP_ do j__=domain_start(2),domain_stop(2)
140 141 142 143 144
#    define _END_OUTER_HORIZONTAL_LOOP_ end do;j__=domain_extent(2)
#  else
     ! No horizontal dimension vectorized; do full outer loop.
#    define _BEGIN_OUTER_HORIZONTAL_LOOP_ _BEGIN_GLOBAL_HORIZONTAL_LOOP_
#    define _END_OUTER_HORIZONTAL_LOOP_ _END_GLOBAL_HORIZONTAL_LOOP_
145 146 147
#  endif
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
148
#ifdef _INTERIOR_IS_VECTORIZED_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
149
   integer :: _START_, _STOP_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
150 151
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
152 153
   real(rke), allocatable _DIMENSION_GLOBAL_            :: tmp
   real(rke), allocatable _DIMENSION_GLOBAL_HORIZONTAL_ :: tmp_hz
154

Jorn Bruggeman's avatar
Jorn Bruggeman committed
155
#ifdef _HAS_MASK_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
156
   _FABM_MASK_TYPE_, allocatable, target _DIMENSION_GLOBAL_HORIZONTAL_ :: mask_hz
Jorn Bruggeman's avatar
Jorn Bruggeman committed
157
#  ifndef _FABM_HORIZONTAL_MASK_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
158
   _FABM_MASK_TYPE_, allocatable, target _DIMENSION_GLOBAL_ :: mask
Jorn Bruggeman's avatar
Jorn Bruggeman committed
159
#  endif
160 161 162
#  ifndef _FABM_MASKED_VALUE_
#    define _FABM_MASKED_VALUE_ _FABM_UNMASKED_VALUE_+1
#  endif
163 164 165
#  ifndef _FABM_UNMASKED_VALUE_
#    define _FABM_UNMASKED_VALUE_ _FABM_MASKED_VALUE_+1
#  endif
166 167
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
168 169
   integer :: interior_count
   integer :: horizontal_count
170

171
#if _FABM_BOTTOM_INDEX_==-1
Jorn Bruggeman's avatar
Jorn Bruggeman committed
172
   integer, allocatable, target _DIMENSION_GLOBAL_HORIZONTAL_ :: bottom_index
Jorn Bruggeman's avatar
Jorn Bruggeman committed
173 174
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
175 176 177
   real(rke), allocatable, target _DIMENSION_GLOBAL_PLUS_1_            :: interior_state
   real(rke), allocatable, target _DIMENSION_GLOBAL_HORIZONTAL_PLUS_1_ :: surface_state
   real(rke), allocatable, target _DIMENSION_GLOBAL_HORIZONTAL_PLUS_1_ :: bottom_state
Jorn Bruggeman's avatar
Jorn Bruggeman committed
178

179 180
   real(rke), allocatable _DIMENSION_SLICE_PLUS_1_            :: dy, w, total_int
   real(rke), allocatable _DIMENSION_HORIZONTAL_SLICE_PLUS_1_ :: flux, sms_sf, sms_bt, total_hz
Jorn Bruggeman's avatar
Jorn Bruggeman committed
181

Jorn Bruggeman's avatar
Jorn Bruggeman committed
182 183 184
   real(rke), allocatable, target _DIMENSION_GLOBAL_            :: temperature
   real(rke), allocatable, target _DIMENSION_GLOBAL_            :: depth
   real(rke), allocatable, target _DIMENSION_GLOBAL_HORIZONTAL_ :: wind_speed
Jorn Bruggeman's avatar
Jorn Bruggeman committed
185

Jorn Bruggeman's avatar
Jorn Bruggeman committed
186
   class (type_fabm_model), pointer :: model
Jorn Bruggeman's avatar
Jorn Bruggeman committed
187

Jorn Bruggeman's avatar
Jorn Bruggeman committed
188
   class (type_test_model), pointer :: test_model
Jorn Bruggeman's avatar
Jorn Bruggeman committed
189

Jorn Bruggeman's avatar
Jorn Bruggeman committed
190
   integer :: domain_extent(_FABM_DIMENSION_COUNT_)
191 192
   integer :: domain_start(_FABM_DIMENSION_COUNT_)
   integer :: domain_stop(_FABM_DIMENSION_COUNT_)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
193

Jorn Bruggeman's avatar
Jorn Bruggeman committed
194 195 196 197 198 199
   character(len=20) :: arg
   integer :: ivar
   integer :: i
   integer :: mode = 1
   integer :: ntest = -1
   logical :: no_mask = .false.
200

201
#if _FABM_DIMENSION_COUNT_>0
Jorn Bruggeman's avatar
Jorn Bruggeman committed
202
   i__ = 50
203 204
#endif
#if _FABM_DIMENSION_COUNT_>1
Jorn Bruggeman's avatar
Jorn Bruggeman committed
205
   j__ = 40
206 207
#endif
#if _FABM_DIMENSION_COUNT_>2
Jorn Bruggeman's avatar
Jorn Bruggeman committed
208
   k__ = 45
209 210
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
211 212 213 214
   ! Parse command line arguments
   call start_test('parsing command line arguments')
   i = 1
   do
215
      call get_command_argument(i, arg)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
216 217 218 219 220 221
      if (arg == '') exit
      select case (arg)
      case ('-s', '--simulate')
         mode = 2
      case ('--nomask')
         no_mask = .true.
222
#if _FABM_DIMENSION_COUNT_>0
Jorn Bruggeman's avatar
Jorn Bruggeman committed
223 224 225 226
      case ('--nx')
         i = i + 1
         call get_command_argument(i, arg)
         read (arg,*) i__
227 228
#endif
#if _FABM_DIMENSION_COUNT_>1
Jorn Bruggeman's avatar
Jorn Bruggeman committed
229 230 231 232
      case ('--ny')
         i = i + 1
         call get_command_argument(i, arg)
         read (arg,*) j__
233 234
#endif
#if _FABM_DIMENSION_COUNT_>2
Jorn Bruggeman's avatar
Jorn Bruggeman committed
235 236 237 238
      case ('--nz')
         i = i + 1
         call get_command_argument(i, arg)
         read (arg,*) k__
239
#endif
Jorn Bruggeman's avatar
Jorn Bruggeman committed
240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256
      case ('-n')
         i = i + 1
         call get_command_argument(i, arg)
         read (arg,*) ntest
      case ('-h')
         write (*,'(a)') ''
         write (*,'(a)') ''
         write (*,'(a)') 'FABM host emulator'
         write (*,'(a)') ''
         write (*,'(a)') 'Accepted arguments:'
         write (*,'(a)') '-s/--simulate: simulate using provided fabm.yaml/environment.yaml'
         write (*,'(a)') '-n:            number of replicates when simulating'
         stop 0
      case default
         write (*,'(a)') 'Unknown command line argument: ' // trim(arg)
         stop 2
      end select
257
      i = i + 1
Jorn Bruggeman's avatar
Jorn Bruggeman committed
258 259
   end do
   call report_test_result()
260

Jorn Bruggeman's avatar
Jorn Bruggeman committed
261
#if _FABM_DIMENSION_COUNT_>0
Jorn Bruggeman's avatar
Jorn Bruggeman committed
262
   domain_extent = (/ _LOCATION_ /)
263 264 265
   domain_start(:) = 1
   domain_stop(:) = domain_extent
   interior_count = product(domain_stop - domain_start + 1)
266
#else
Jorn Bruggeman's avatar
Jorn Bruggeman committed
267
   interior_count = 1
Jorn Bruggeman's avatar
Jorn Bruggeman committed
268 269
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
270 271 272 273 274 275 276
   ! Set defaults
   if (ntest == -1) then
      if (mode == 1) then
         ntest = 1
      else
         ntest = 50000000 / interior_count
      end if
277 278
   end if

Jorn Bruggeman's avatar
Jorn Bruggeman committed
279
#ifdef _INTERIOR_IS_VECTORIZED_
280 281
   _START_ = domain_start(_FABM_VECTORIZED_DIMENSION_INDEX_)
   _STOP_ = domain_stop(_FABM_VECTORIZED_DIMENSION_INDEX_)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
282 283
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304
   allocate(tmp _INDEX_LOCATION_)
   allocate(tmp_hz _INDEX_HORIZONTAL_LOCATION_)

   allocate(type_test_driver::driver)

   call start_test('fabm_initialize_library')
   call fabm_initialize_library()
   call report_test_result()

   call start_test('building model tree')
   select case (mode)
   case (1)
       ! Unit testing with built-in model
       allocate(model)
       allocate(test_model)
       call model%root%add_child(test_model, 'test_model', 'test model', configunit=-1)
   case (2)
       ! Test with user-provided fabm.yaml
       model => fabm_create_model(initialize=.false.)
   end select
   call report_test_result()
305

Jorn Bruggeman's avatar
Jorn Bruggeman committed
306 307
   call start_test('initialize')
   call model%initialize()
308 309 310
   call assert(size(model%interior_state_variables) == test_model%nstate, 'model%initialize', 'Incorrect number of interior state variables.')
   call assert(size(model%bottom_state_variables) == test_model%nbottom_state, 'model%initialize', 'Incorrect number of bottom state variables.')
   call assert(size(model%surface_state_variables) == test_model%nsurface_state, 'model%initialize', 'Incorrect number of surface state variables.')
311
   call assert(size(model%conserved_quantities) == 1, 'model%initialize', 'Incorrect number of conserved quantities.')
Jorn Bruggeman's avatar
Jorn Bruggeman committed
312
   call report_test_result()
Jorn Bruggeman's avatar
Jorn Bruggeman committed
313

Jorn Bruggeman's avatar
Jorn Bruggeman committed
314 315 316
   ! ======================================================================
   ! Provide extents of the spatial domain.
   ! ======================================================================
Jorn Bruggeman's avatar
Jorn Bruggeman committed
317

Jorn Bruggeman's avatar
Jorn Bruggeman committed
318 319 320
   call start_test('set_domain')
   call model%set_domain(_LOCATION_)
   call report_test_result()
Jorn Bruggeman's avatar
Jorn Bruggeman committed
321

Jorn Bruggeman's avatar
Jorn Bruggeman committed
322 323 324
   ! ======================================================================
   ! Set up spatial mask.
   ! ======================================================================
Jorn Bruggeman's avatar
Jorn Bruggeman committed
325 326

#ifdef _HAS_MASK_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
327 328
   allocate(mask_hz _INDEX_HORIZONTAL_LOCATION_)
   call start_test('set_mask')
Jorn Bruggeman's avatar
Jorn Bruggeman committed
329
#  ifdef _FABM_HORIZONTAL_MASK_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
330
   call model%set_mask(mask_hz)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
331
#  else
Jorn Bruggeman's avatar
Jorn Bruggeman committed
332 333
   allocate(mask _INDEX_LOCATION_)
   call model%set_mask(mask, mask_hz)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
334
#  endif
Jorn Bruggeman's avatar
Jorn Bruggeman committed
335
   call report_test_result()
Jorn Bruggeman's avatar
Jorn Bruggeman committed
336 337
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
338 339 340
   ! ======================================================================
   ! Specify vertical indices of surface and bottom.
   ! ======================================================================
Jorn Bruggeman's avatar
Jorn Bruggeman committed
341

342
#if defined(_FABM_DEPTH_DIMENSION_INDEX_)&&_FABM_BOTTOM_INDEX_==-1
Jorn Bruggeman's avatar
Jorn Bruggeman committed
343 344 345 346
   call start_test('set_bottom_index')
   allocate(bottom_index _INDEX_HORIZONTAL_LOCATION_)
   call model%set_bottom_index(bottom_index)
   call report_test_result()
Jorn Bruggeman's avatar
Jorn Bruggeman committed
347 348
#endif

349
   allocate(interior_state(_PREARG_LOCATION_ size(model%interior_state_variables)))
Jorn Bruggeman's avatar
Jorn Bruggeman committed
350 351 352 353 354 355 356 357
   allocate(surface_state(_PREARG_HORIZONTAL_LOCATION_ size(model%surface_state_variables)))
   allocate(bottom_state(_PREARG_HORIZONTAL_LOCATION_ size(model%bottom_state_variables)))

   ! ======================================================================
   ! Send pointers to state variable data to FABM.
   ! ======================================================================

   call start_test('link_interior_state_data')
358
   do ivar = 1, size(model%interior_state_variables)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
359 360 361 362 363 364 365 366 367 368 369 370 371 372 373 374 375 376 377 378 379 380 381 382 383 384 385
      call model%link_interior_state_data(ivar, interior_state(_PREARG_LOCATION_DIMENSIONS_ ivar))
   end do
   call report_test_result()

   call start_test('link_surface_state_data')
   do ivar = 1, size(model%surface_state_variables)
      call model%link_surface_state_data(ivar, surface_state(_PREARG_HORIZONTAL_LOCATION_DIMENSIONS_ ivar))
   end do
   call report_test_result()

   call start_test('link_bottom_state_data')
   do ivar = 1, size(model%bottom_state_variables)
      call model%link_bottom_state_data(ivar, bottom_state(_PREARG_HORIZONTAL_LOCATION_DIMENSIONS_ ivar))
   end do
   call report_test_result()

   ! ======================================================================
   ! Transfer pointers to environmental data
   ! ======================================================================

   select case (mode)
   case (1)
       allocate(depth _INDEX_LOCATION_)
       allocate(temperature _INDEX_LOCATION_)
       allocate(wind_speed _INDEX_HORIZONTAL_LOCATION_)

       call start_test('link_interior_data')
386 387
       call model%link_interior_data(fabm_standard_variables%temperature, temperature)
       call model%link_interior_data(fabm_standard_variables%depth, depth)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
388 389 390
       call report_test_result()

       call start_test('link_horizontal_data')
391
       call model%link_horizontal_data(fabm_standard_variables%wind_speed, wind_speed)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
392 393 394 395 396 397 398 399 400 401 402 403 404
       call report_test_result()
   case (2)
       call read_environment
   end select

   ! ======================================================================
   ! Check whether FABM has all dependencies fulfilled
   ! (i.e., whether all required calls for link_*_data have been made)
   ! ======================================================================

   call start_test('start')
   call model%start()
   call report_test_result()
Jorn Bruggeman's avatar
Jorn Bruggeman committed
405 406

#ifdef _FABM_VECTORIZED_DIMENSION_INDEX_
407 408
   allocate(dy(_START_:_STOP_, size(model%interior_state_variables)))
   allocate(w(_START_:_STOP_, size(model%interior_state_variables)))
409
   allocate(total_int(_START_:_STOP_, size(model%conserved_quantities)))
Jorn Bruggeman's avatar
Jorn Bruggeman committed
410
#else
411 412
   allocate(dy(size(model%interior_state_variables)))
   allocate(w(size(model%interior_state_variables)))
413
   allocate(total_int(size(model%conserved_quantities)))
Jorn Bruggeman's avatar
Jorn Bruggeman committed
414 415 416
#endif

#if defined(_FABM_VECTORIZED_DIMENSION_INDEX_)&&(_FABM_DEPTH_DIMENSION_INDEX_!=_FABM_VECTORIZED_DIMENSION_INDEX_)
417
   allocate(flux(_START_:_STOP_, size(model%interior_state_variables)))
Jorn Bruggeman's avatar
Jorn Bruggeman committed
418 419
   allocate(sms_sf(_START_:_STOP_, size(model%surface_state_variables)))
   allocate(sms_bt(_START_:_STOP_, size(model%bottom_state_variables)))
420
   allocate(total_hz(_START_:_STOP_, size(model%conserved_quantities)))
Jorn Bruggeman's avatar
Jorn Bruggeman committed
421
#else
422
   allocate(flux(size(model%interior_state_variables)))
Jorn Bruggeman's avatar
Jorn Bruggeman committed
423 424
   allocate(sms_sf(size(model%surface_state_variables)))
   allocate(sms_bt(size(model%bottom_state_variables)))
425
   allocate(total_hz(size(model%conserved_quantities)))
Jorn Bruggeman's avatar
Jorn Bruggeman committed
426 427
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
428 429
   select case (mode)
   case (1)
430 431
      write (*,'(a)') 'Testing without mask and unrestricted domain.'
      call test_update(apply_mask=.false., restrict_range=.false.)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
432
#if defined(_HAS_MASK_) || _FABM_DIMENSION_COUNT_>0
433
      write (*,'(a,i0,a)') 'Running ', ntest, ' tests with randomized domain settings.'
Jorn Bruggeman's avatar
Jorn Bruggeman committed
434
      do i=1,ntest
435 436 437 438 439 440 441 442
#  ifdef _HAS_MASK_
         write (*,'(a)') 'Random mask (unrestricted domain).'
         call test_update(apply_mask=.true., restrict_range=.false.)
#  endif
#  if _FABM_DIMENSION_COUNT_>0
         write (*,'(a)') 'Randomly restricted domain (no mask).'
         call test_update(apply_mask=.false., restrict_range=.true.)
#  endif
Jorn Bruggeman's avatar
Jorn Bruggeman committed
443
#  if defined(_HAS_MASK_) && _FABM_DIMENSION_COUNT_>0
444 445 446
         write (*,'(a)') 'Random mask and randomly restricted domain.'
         call test_update(apply_mask=.true., restrict_range=.true.)
#  endif
Jorn Bruggeman's avatar
Jorn Bruggeman committed
447
      end do
448
#endif
Jorn Bruggeman's avatar
Jorn Bruggeman committed
449 450 451
   case(2)
      call simulate(ntest)
   end select
Jorn Bruggeman's avatar
Jorn Bruggeman committed
452

453 454 455 456
contains

   subroutine read_environment
      use yaml, only: yaml_parse => parse, yaml_error_length => error_length
457 458
      use yaml_types, only: type_node, type_yaml_dictionary => type_dictionary, type_yaml_scalar => type_scalar, &
         type_yaml_key_value_pair => type_key_value_pair, yaml_real_kind => real_kind
459 460 461

      integer, parameter :: yaml_unit = 100
      character(yaml_error_length) :: yaml_error
462
      class (type_node),pointer :: yaml_root
463
      type (type_yaml_key_value_pair), pointer :: yaml_pair
464
      real(rke) :: value
465 466
      logical :: success
      type type_input
467 468 469
         type (type_fabm_interior_variable_id)                :: interior_id
         type (type_fabm_horizontal_variable_id)              :: horizontal_id
         type (type_fabm_scalar_variable_id)                  :: scalar_id
470 471 472
         real(rke), allocatable _DIMENSION_GLOBAL_            :: interior_data
         real(rke), allocatable _DIMENSION_GLOBAL_HORIZONTAL_ :: horizontal_data
         real(rke)                                            :: scalar_data
473 474 475 476 477 478 479 480 481 482 483 484 485 486 487
      end type
      type (type_input), pointer :: input

      yaml_root => yaml_parse('environment.yaml', yaml_unit, yaml_error)
      if (yaml_error /= '') then
         call driver%log_message(yaml_error)
         stop 2
      end if
      select type (yaml_root)
      class is (type_yaml_dictionary)
          yaml_pair => yaml_root%first
          do while (associated(yaml_pair))
              select type (node => yaml_pair%value)
              class is (type_yaml_scalar)
                  call driver%log_message('Setting '//trim(yaml_pair%key)//' to '//trim(node%string))
488
                  value = node%to_real(0._yaml_real_kind, success)
489 490 491 492 493
                  if (.not. success) then
                     call driver%log_message('Cannot parse '//trim(node%string)//' as real.')
                     stop 2
                  end if
                  allocate(input)
494
                  input%interior_id = model%get_interior_variable_id(trim(yaml_pair%key))
495
                  if (model%is_variable_used(input%interior_id)) then
496 497 498 499 500
                      allocate(input%interior_data _INDEX_LOCATION_)
                      input%interior_data = value
                      call model%link_interior_data(input%interior_id, input%interior_data)
                  else
                      input%horizontal_id = model%get_horizontal_variable_id(trim(yaml_pair%key))
501
                      if (model%is_variable_used(input%horizontal_id)) then
502 503 504 505 506
                         allocate(input%horizontal_data _INDEX_HORIZONTAL_LOCATION_)
                         input%horizontal_data = value
                         call model%link_horizontal_data(input%horizontal_id, input%horizontal_data)
                      else
                         input%scalar_id = model%get_scalar_variable_id(trim(yaml_pair%key))
507
                         if (model%is_variable_used(input%scalar_id)) then
508 509 510
                            input%scalar_data = value
                            call model%link_scalar(input%scalar_id, input%scalar_data)
                         else
511 512
                            call driver%log_message('WARNING: environment variable '//trim(yaml_pair%key) &
                               //' is not used by FABM model and will be ignored.')
513 514 515 516 517 518 519 520 521 522 523
                         end if
                      end if
                  end if
              end select
              yaml_pair => yaml_pair%next
          end do
      class default
         call driver%log_message('environment.yaml should contain a dictionary at root level')
         stop 2
      end select
   end subroutine read_environment
Jorn Bruggeman's avatar
Jorn Bruggeman committed
524

525 526 527 528 529 530 531 532 533 534 535 536 537 538 539 540 541 542 543 544 545 546 547 548 549 550 551 552 553 554 555 556 557 558 559 560 561 562 563 564 565 566 567 568 569
   subroutine configure_range(randomize)
      logical, intent(in) :: randomize

      real(rke), parameter :: excluded_fraction = 0.5_rke

#if _FABM_DIMENSION_COUNT_ > 0
      integer :: i
      real(rke) :: rnd(2)
      real(rke), parameter :: cut_fraction = 0.5_rke * (1._rke - excluded_fraction ** (1._rke / _FABM_DIMENSION_COUNT_))

      if (.not. randomize) then
         domain_start(:) = 1
         domain_stop(:) = domain_extent
      else
         do i = 1, _FABM_DIMENSION_COUNT_
            call random_number(rnd)
            domain_start(i) = 1 + int(domain_extent(i) * cut_fraction * rnd(1))
            domain_stop(i) = domain_start(i) + int((domain_extent(i) - domain_start(i) + 1) * (1._rke - cut_fraction) * rnd(2))
            write (*,'(A,I0,A,I0,A,I0)') 'Dimension ', i, ': ', domain_start(i), ' - ', domain_stop(i)
         end do
      end if

#  if _FABM_DIMENSION_COUNT_ == 1
      call model%set_domain_start(domain_start(1))
      call model%set_domain_stop(domain_stop(1))
#  elif _FABM_DIMENSION_COUNT_ == 2
      call model%set_domain_start(domain_start(1), domain_start(2))
      call model%set_domain_stop(domain_stop(1), domain_stop(2))
#  elif _FABM_DIMENSION_COUNT_ == 3
      call model%set_domain_start(domain_start(1), domain_start(2), domain_start(3))
      call model%set_domain_stop(domain_stop(1), domain_stop(2), domain_stop(3))
#  endif
#  ifdef _FABM_VECTORIZED_DIMENSION_INDEX_
      _START_ = domain_start(_FABM_VECTORIZED_DIMENSION_INDEX_)
      _STOP_ = domain_stop(_FABM_VECTORIZED_DIMENSION_INDEX_)
#  endif
#endif
   end subroutine configure_range

   subroutine configure_mask(randomize)
      logical, intent(in) :: randomize

      real(rke), parameter :: masked_fraction = 0.5_rke

      if (randomize) then
570
#if _FABM_BOTTOM_INDEX_==-1
571
      ! Depth index of bottom varies in the horizontal
572
      call random_number(tmp_hz)
573
#  ifdef _HAS_MASK_
574 575
      !  Pick random numbers between start-1 and stop index [inclusive]. start-1 means invalid (land)
      bottom_index = domain_start(_FABM_DEPTH_DIMENSION_INDEX_) - 1 + floor(tmp_hz * (2 + domain_stop(_FABM_DEPTH_DIMENSION_INDEX_) - domain_start(_FABM_DEPTH_DIMENSION_INDEX_)))
576
#    ifdef _FABM_VERTICAL_BOTTOM_TO_SURFACE_
577 578
      ! Use stop+1 as invalid bottom index to ensure vertical loops will have 0 iterations
      bottom_index = bottom_index + 1
579 580
#    endif
#  else
581 582
      ! Pick random numbers between start and stop index [inclusive]
      bottom_index = domain_start(_FABM_DEPTH_DIMENSION_INDEX_) + floor(tmp_hz * (1 + domain_stop(_FABM_DEPTH_DIMENSION_INDEX_) - domain_start(_FABM_DEPTH_DIMENSION_INDEX_)))
583 584 585
#  endif
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
586
#ifdef _HAS_MASK_
587 588 589 590
#  ifdef _FABM_HORIZONTAL_MASK_
      ! Apply random mask across horizontal domain (half of grid cells masked)
      call random_number(tmp_hz)
      mask_hz = _FABM_UNMASKED_VALUE_
591
      where (tmp_hz > masked_fraction) mask_hz = _FABM_MASKED_VALUE_
592 593 594
#  else
      ! Apply random mask across interior domain (half of grid cells masked)
      call random_number(tmp)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
595
      mask = _FABM_UNMASKED_VALUE_
596
      where (tmp > masked_fraction) mask = _FABM_MASKED_VALUE_
597 598

#    if _FABM_BOTTOM_INDEX_==-1
599
      ! Bottom index varies in the horizontal. Ensure the bottom cell itself is unmasked, and anything deeper is masked.
600
      _BEGIN_GLOBAL_HORIZONTAL_LOOP_
601
         ! Valid bottom index - unmask associated cell, then mask all deeper ones
602
#      ifdef _FABM_VERTICAL_BOTTOM_TO_SURFACE_
603
         if (bottom_index _INDEX_HORIZONTAL_LOCATION_ <= domain_stop(_FABM_DEPTH_DIMENSION_INDEX_)) mask _INDEX_GLOBAL_VERTICAL_(bottom_index _INDEX_HORIZONTAL_LOCATION_) = _FABM_UNMASKED_VALUE_
604
         mask _INDEX_GLOBAL_VERTICAL_(:bottom_index _INDEX_HORIZONTAL_LOCATION_ - 1) = _FABM_MASKED_VALUE_
605
#      else
606
         if (bottom_index _INDEX_HORIZONTAL_LOCATION_ >= domain_start(_FABM_DEPTH_DIMENSION_INDEX_)) mask _INDEX_GLOBAL_VERTICAL_(bottom_index _INDEX_HORIZONTAL_LOCATION_) = _FABM_UNMASKED_VALUE_
607
         mask _INDEX_GLOBAL_VERTICAL_(bottom_index _INDEX_HORIZONTAL_LOCATION_ + 1:) = _FABM_MASKED_VALUE_
608
#      endif
609 610 611 612
      _END_GLOBAL_HORIZONTAL_LOOP_
#    endif

      ! Infer horizontal mask (mask points that have all column layers masked)
613 614 615 616
      mask_hz = _FABM_MASKED_VALUE_
      _BEGIN_GLOBAL_LOOP_
         if (_IS_UNMASKED_(mask _INDEX_LOCATION_)) mask_hz _INDEX_HORIZONTAL_LOCATION_ = _FABM_UNMASKED_VALUE_
      _END_GLOBAL_LOOP_
617 618 619 620

      ! For valid points in the horizontal, make sure that index 1 (bottom or surface) is unmasked
      _BEGIN_GLOBAL_HORIZONTAL_LOOP_
         if (_IS_UNMASKED_(mask_hz _INDEX_HORIZONTAL_LOCATION_)) then
621
            mask _INDEX_GLOBAL_VERTICAL_(domain_start(_FABM_DEPTH_DIMENSION_INDEX_)) = _FABM_UNMASKED_VALUE_
622
#    if _FABM_BOTTOM_INDEX_!=-1
623
            mask _INDEX_GLOBAL_VERTICAL_(domain_stop(_FABM_DEPTH_DIMENSION_INDEX_)) = _FABM_UNMASKED_VALUE_
624 625 626
#    endif
         end if
      _END_GLOBAL_HORIZONTAL_LOOP_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
627 628 629 630 631
#  endif
#endif
      else
#if _FABM_BOTTOM_INDEX_==-1
#  ifdef _FABM_VERTICAL_BOTTOM_TO_SURFACE_
632
      bottom_index = domain_start(_FABM_DEPTH_DIMENSION_INDEX_)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
633
#  else
634
      bottom_index = domain_stop(_FABM_DEPTH_DIMENSION_INDEX_)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
635 636 637 638 639 640 641 642 643 644
#  endif
#endif
#ifdef _HAS_MASK_
#  ifndef _FABM_HORIZONTAL_MASK_
      mask = _FABM_UNMASKED_VALUE_
#  endif
      mask_hz = _FABM_UNMASKED_VALUE_
#endif
      end if

645 646 647 648 649 650 651 652 653 654 655
#if _FABM_BOTTOM_INDEX_==-1 && defined(_HAS_MASK_)
      ! Bottom index varies in the horizontal. Ensure the bottom cell itself is unmasked, and anything deeper is masked.
      _BEGIN_GLOBAL_HORIZONTAL_LOOP_
         if (_IS_UNMASKED_(mask_hz _INDEX_HORIZONTAL_LOCATION_)) then
            call assert(bottom_index _INDEX_HORIZONTAL_LOCATION_ >= domain_start(_FABM_DEPTH_DIMENSION_INDEX_) &
              .and. bottom_index _INDEX_HORIZONTAL_LOCATION_ <= domain_stop(_FABM_DEPTH_DIMENSION_INDEX_), 'randomize_mask', &
               'BUG: unmaked horizontal location has invalid bottom index')
         end if
      _END_GLOBAL_HORIZONTAL_LOOP_
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
656
      call count_active_points()
657
   end subroutine configure_mask
Jorn Bruggeman's avatar
Jorn Bruggeman committed
658 659

   subroutine count_active_points()
660 661 662 663 664 665 666
      logical :: active
      integer :: i, nhz, nhz_active

      active = .true.

      interior_count = 0
      _BEGIN_GLOBAL_LOOP_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
667 668
#ifdef _HAS_MASK_
#  ifdef _FABM_HORIZONTAL_MASK_
669
         active = _IS_UNMASKED_(mask_hz _INDEX_HORIZONTAL_LOCATION_)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
670
#  else
671
         active = _IS_UNMASKED_(mask _INDEX_LOCATION_)
672
#  endif
673
#elif _FABM_BOTTOM_INDEX_==-1 && _FABM_VECTORIZED_DIMENSION_INDEX_==_FABM_DEPTH_DIMENSION_INDEX_ && defined(_FABM_DEPTH_DIMENSION_INDEX_)
674
#  ifdef _FABM_VERTICAL_BOTTOM_TO_SURFACE_
675
         active = _ITERATOR_ >= bottom_index _INDEX_HORIZONTAL_LOCATION_
676
#  else
677
         active = _ITERATOR_ <= bottom_index _INDEX_HORIZONTAL_LOCATION_
678
#  endif
Jorn Bruggeman's avatar
Jorn Bruggeman committed
679
#endif
680 681 682 683 684 685 686 687 688 689 690 691 692 693 694 695 696 697 698 699 700 701 702 703 704
         if (active) interior_count = interior_count + 1
      _END_GLOBAL_LOOP_

      horizontal_count = 0
      _BEGIN_GLOBAL_HORIZONTAL_LOOP_
#ifdef _HAS_MASK_
         active = _IS_UNMASKED_(mask_hz _INDEX_HORIZONTAL_LOCATION_)
#endif
         if (active) horizontal_count = horizontal_count + 1
      _END_GLOBAL_HORIZONTAL_LOOP_

      nhz = 1
      nhz_active = 1
      do i = 1, _FABM_DIMENSION_COUNT_
#ifdef _FABM_DEPTH_DIMENSION_INDEX_
         if (i /= _FABM_DEPTH_DIMENSION_INDEX_) then
#endif
         nhz = nhz * domain_extent(i)
         nhz_active = nhz_active * (domain_stop(i) - domain_start(i) + 1)
#ifdef _FABM_DEPTH_DIMENSION_INDEX_
         end if
#endif
      end do
      write (*,'(a,i0,a,i0,a,i0,a)') 'Interior: ', product(domain_extent), ' points, ', product(domain_stop - domain_start + 1), ' in active range, ', interior_count, ' unmasked'
      write (*,'(a,i0,a,i0,a,i0,a)') 'Horizontal: ', nhz, ' points, ', nhz_active, ' in active range, ', horizontal_count, ' unmasked'
Jorn Bruggeman's avatar
Jorn Bruggeman committed
705
   end subroutine
706

707 708
   subroutine simulate(n)
      integer, intent(in) :: n
709
      real(rke) :: time_begin, time_end
710 711 712 713 714 715 716 717
      integer :: nseed
      integer, allocatable :: seed(:)

      call random_seed(size=nseed)
      allocate(seed(nseed))
      seed(:) = 1
      call random_seed(put=seed)

718 719
      call configure_range(.false.)
      call configure_mask(.not. no_mask)
720

721
      call start_test('initialize_interior_state')
722
      _BEGIN_OUTER_INTERIOR_LOOP_
723
         call model%initialize_interior_state(_ARG_INTERIOR_IN_)
724 725 726
      _END_OUTER_INTERIOR_LOOP_
      call report_test_result()

727
      call start_test('initialize_bottom_state')
728
      _BEGIN_OUTER_HORIZONTAL_LOOP_
729
         call model%initialize_bottom_state(_ARG_HORIZONTAL_IN_)
730 731 732
      _END_OUTER_HORIZONTAL_LOOP_
      call report_test_result()

733
      call start_test('initialize_surface_state')
734
      _BEGIN_OUTER_HORIZONTAL_LOOP_
735
         call model%initialize_surface_state(_ARG_HORIZONTAL_IN_)
736 737 738
      _END_OUTER_HORIZONTAL_LOOP_
      call report_test_result()

739
      write (*,'(a,i0,a)') 'Simulating with ', interior_count, ' wet cells...'
740 741 742

      call cpu_time(time_begin)

743 744
      do i = 1, n
         call model%prepare_inputs()
745 746 747 748

         _BEGIN_OUTER_HORIZONTAL_LOOP_
            flux = 0
            sms_bt = 0
Jorn Bruggeman's avatar
Jorn Bruggeman committed
749
            call model%get_bottom_sources(_PREARG_HORIZONTAL_IN_ flux, sms_bt)
750 751 752 753 754
         _END_OUTER_HORIZONTAL_LOOP_

         _BEGIN_OUTER_HORIZONTAL_LOOP_
            flux = 0
            sms_sf = 0
Jorn Bruggeman's avatar
Jorn Bruggeman committed
755
            call model%get_surface_sources(_PREARG_HORIZONTAL_IN_ flux, sms_sf)
756 757 758 759
         _END_OUTER_HORIZONTAL_LOOP_

         _BEGIN_OUTER_INTERIOR_LOOP_
            dy = 0
Jorn Bruggeman's avatar
Jorn Bruggeman committed
760
            call model%get_interior_sources(_PREARG_INTERIOR_IN_ dy)
761
         _END_OUTER_INTERIOR_LOOP_
762

763
         call model%finalize_outputs()
764

765
         if (mod(i, 100) == 0) write (*,'(i0,a)') int(100*i/real(n, rke)), ' % complete'
766 767 768 769
      end do

      call cpu_time(time_end)

770 771
      write (*,'(a)') 'Simulation complete.'
      write (*,'(a,f8.3,a)') 'Total time: ', time_end - time_begin, ' s'
772 773
   end subroutine

774 775 776
   subroutine test_update(apply_mask, restrict_range)
      logical, intent(in) :: apply_mask, restrict_range

777 778 779
      real(rke), pointer _DIMENSION_GLOBAL_            :: pdata
      real(rke), pointer _DIMENSION_GLOBAL_HORIZONTAL_ :: pdata_hz
      logical                                          :: valid
Jorn Bruggeman's avatar
Jorn Bruggeman committed
780

781 782 783 784
      call configure_range(randomize=restrict_range)
      call configure_mask(randomize=apply_mask)

      call model%start()
Jorn Bruggeman's avatar
Jorn Bruggeman committed
785 786 787 788

      ! ======================================================================
      ! Initialize all state variables
      ! ======================================================================
Jorn Bruggeman's avatar
Jorn Bruggeman committed
789

790
      call start_test('initialize_interior_state')
791
      _BEGIN_OUTER_INTERIOR_LOOP_
792
         call model%initialize_interior_state(_ARG_INTERIOR_IN_)
793
      _END_OUTER_INTERIOR_LOOP_
794
      do ivar = 1, size(model%interior_state_variables)
795
         call check_interior(interior_state(_PREARG_LOCATION_DIMENSIONS_ ivar), &
796
            model%interior_state_variables(ivar)%missing_value, ivar+interior_state_offset+1._rke)
797
      end do
Jorn Bruggeman's avatar
Jorn Bruggeman committed
798
      call report_test_result()
Jorn Bruggeman's avatar
Jorn Bruggeman committed
799

800
      call start_test('initialize_bottom_state')
801
      _BEGIN_OUTER_HORIZONTAL_LOOP_
802
         call model%initialize_bottom_state(_ARG_HORIZONTAL_IN_)
803
      _END_OUTER_HORIZONTAL_LOOP_
804
      do ivar = 1, size(model%bottom_state_variables)
805 806
         call check_horizontal(bottom_state(_PREARG_HORIZONTAL_LOCATION_DIMENSIONS_ ivar), &
            model%bottom_state_variables(ivar)%missing_value, ivar+bottom_state_offset+1._rke)
807
      end do
Jorn Bruggeman's avatar
Jorn Bruggeman committed
808
      call report_test_result()
Jorn Bruggeman's avatar
Jorn Bruggeman committed
809

810
      call start_test('initialize_surface_state')
811
      _BEGIN_OUTER_HORIZONTAL_LOOP_
812
         call model%initialize_surface_state(_ARG_HORIZONTAL_IN_)
813
      _END_OUTER_HORIZONTAL_LOOP_
814
      do ivar = 1, size(model%surface_state_variables)
815 816
         call check_horizontal(surface_state(_PREARG_HORIZONTAL_LOCATION_DIMENSIONS_ ivar), &
            model%surface_state_variables(ivar)%missing_value, ivar+surface_state_offset+1._rke)
817
      end do
Jorn Bruggeman's avatar
Jorn Bruggeman committed
818
      call report_test_result()
Jorn Bruggeman's avatar
Jorn Bruggeman committed
819

Jorn Bruggeman's avatar
Jorn Bruggeman committed
820 821 822 823
      ! ======================================================================
      ! Initialize environmental dependencies
      ! ======================================================================

824 825
      temperature = 1 + interior_dependency_offset
      call apply_mask_3d(temperature, -999._rke - interior_dependency_offset)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
826

827 828
      wind_speed = 1 + horizontal_dependency_offset
      call apply_mask_2d(wind_speed, -999._rke - horizontal_dependency_offset)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
829 830 831

#ifdef _FABM_DEPTH_DIMENSION_INDEX_
      ! Model has depth dimension: make sure depth varies from 0 at the surface till 1 at the bottom
832
      _BEGIN_GLOBAL_LOOP_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
833
#  ifdef _FABM_VERTICAL_BOTTOM_TO_SURFACE_
834
#    if _FABM_BOTTOM_INDEX_==-1
835
         if (bottom_index _INDEX_HORIZONTAL_LOCATION_ == domain_stop(_FABM_DEPTH_DIMENSION_INDEX_)) then
836 837
            depth _INDEX_LOCATION_ = 2
         else
838
            depth _INDEX_LOCATION_ = real(domain_stop(_FABM_DEPTH_DIMENSION_INDEX_) - _VERTICAL_ITERATOR_, rke) / (domain_stop(_FABM_DEPTH_DIMENSION_INDEX_) - bottom_index _INDEX_HORIZONTAL_LOCATION_)
839 840
         end if
#    else
841
         depth _INDEX_LOCATION_ = real(domain_stop(_FABM_DEPTH_DIMENSION_INDEX_) - _VERTICAL_ITERATOR_, rke) / (domain_stop(_FABM_DEPTH_DIMENSION_INDEX_) - domain_start(_FABM_DEPTH_DIMENSION_INDEX_))
842
#    endif
Jorn Bruggeman's avatar
Jorn Bruggeman committed