host.F90 55.9 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 57 58
#if _FABM_BOTTOM_INDEX_==-1 && !defined(_HAS_MASK_)
!  No mask but variable bottom index. Index of depth dimension must be 1.
!  All loops over inner dimension should skip points below bottom.
#  ifdef _FABM_VERTICAL_BOTTOM_TO_SURFACE_
59 60
#    define _IMIN_ bottom_index _INDEX_HORIZONTAL_LOCATION_
#    define _IMAX_ domain_extent(1)
61
#  else
62 63
#    define _IMIN_ 1
#    define _IMAX_ bottom_index _INDEX_HORIZONTAL_LOCATION_
64 65 66
#  endif
#else
!  Loops over inner dimension should span full domain
67 68
#  define _IMIN_ 1
#  define _IMAX_ domain_extent(1)
69
#endif
70
#define _IRANGE_ _IMIN_,_IMAX_
71

72 73 74 75
#if _FABM_DIMENSION_COUNT_==0
#  define _BEGIN_GLOBAL_LOOP_
#  define _END_GLOBAL_LOOP_
#elif _FABM_DIMENSION_COUNT_==1
76 77 78
#  define _BEGIN_GLOBAL_LOOP_ do i__=_IRANGE_
#  define _END_GLOBAL_LOOP_ end do;i__=domain_extent(1)
#  ifdef _FABM_DEPTH_DIMENSION_INDEX_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
79
#    define _BEGIN_GLOBAL_HORIZONTAL_LOOP_
80
#    define _END_GLOBAL_HORIZONTAL_LOOP_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
81
#  endif
82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121
#elif _FABM_DIMENSION_COUNT_==2
#  define _BEGIN_GLOBAL_LOOP_ do j__=1,domain_extent(2);do i__=_IRANGE_
#  define _END_GLOBAL_LOOP_ end do;end do;i__=domain_extent(1);j__=domain_extent(2)
#  if _FABM_DEPTH_DIMENSION_INDEX_==1
#    define _BEGIN_GLOBAL_HORIZONTAL_LOOP_ do j__=1,domain_extent(2)
#    define _END_GLOBAL_HORIZONTAL_LOOP_ end do;j__=domain_extent(2)
#  elif _FABM_DEPTH_DIMENSION_INDEX_==2
#    define _BEGIN_GLOBAL_HORIZONTAL_LOOP_ do i__=_IRANGE_
#    define _END_GLOBAL_HORIZONTAL_LOOP_ end do;i__=domain_extent(1)
#  endif
#elif _FABM_DIMENSION_COUNT_==3
#  define _BEGIN_GLOBAL_LOOP_ do k__=1,domain_extent(3);do j__=1,domain_extent(2);do i__=_IRANGE_
#  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
#    define _BEGIN_GLOBAL_HORIZONTAL_LOOP_ do k__=1,domain_extent(3);do j__=1,domain_extent(2)
#    define _END_GLOBAL_HORIZONTAL_LOOP_ end do;end do;j__=domain_extent(2);k__=domain_extent(3)
#  elif _FABM_DEPTH_DIMENSION_INDEX_==2
#    define _BEGIN_GLOBAL_HORIZONTAL_LOOP_ do k__=1,domain_extent(3);do i__=_IRANGE_
#    define _END_GLOBAL_HORIZONTAL_LOOP_ end do;end do;i__=domain_extent(1);k__=domain_extent(3)
#  elif _FABM_DEPTH_DIMENSION_INDEX_==3
#    define _BEGIN_GLOBAL_HORIZONTAL_LOOP_ do j__=1,domain_extent(2);do i__=_IRANGE_
#    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
122
#  define _END_OUTER_INTERIOR_LOOP_
123 124
#  define _BEGIN_OUTER_HORIZONTAL_LOOP_
#  define _END_OUTER_HORIZONTAL_LOOP_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
125
#elif _FABM_DIMENSION_COUNT_==2
126 127
#  define _BEGIN_OUTER_INTERIOR_LOOP_ do j__=1,domain_extent(2)
#  define _END_OUTER_INTERIOR_LOOP_ end do;j__=domain_extent(2)
128
#  if _FABM_DEPTH_DIMENSION_INDEX_==2
129
     ! The entire horizontal is already vectorized; no outer loop necessary
130
#    define _BEGIN_OUTER_HORIZONTAL_LOOP_
131 132 133 134 135
#    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_
136 137
#  endif
#elif _FABM_DIMENSION_COUNT_==3
138 139 140 141 142 143
#  define _BEGIN_OUTER_INTERIOR_LOOP_ do k__=1,domain_extent(3);do j__=1,domain_extent(2)
#  define _END_OUTER_INTERIOR_LOOP_ end do;end do;j__=domain_extent(2);k__=domain_extent(3)
#  if _FABM_DEPTH_DIMENSION_INDEX_==2
#    define _BEGIN_OUTER_HORIZONTAL_LOOP_ do k__=1,domain_extent(3)
#    define _END_OUTER_HORIZONTAL_LOOP_ end do;k__=domain_extent(3)
#  elif _FABM_DEPTH_DIMENSION_INDEX_==3
144
#    define _BEGIN_OUTER_HORIZONTAL_LOOP_ do j__=1,domain_extent(2)
145 146 147 148 149
#    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_
150 151 152
#  endif
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
153
#ifdef _INTERIOR_IS_VECTORIZED_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
154
   integer :: _START_, _STOP_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
155 156
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
157 158
   real(rke), allocatable _DIMENSION_GLOBAL_            :: tmp
   real(rke), allocatable _DIMENSION_GLOBAL_HORIZONTAL_ :: tmp_hz
159

Jorn Bruggeman's avatar
Jorn Bruggeman committed
160
#ifdef _HAS_MASK_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
161
   _FABM_MASK_TYPE_, allocatable, target _DIMENSION_GLOBAL_HORIZONTAL_ :: mask_hz
Jorn Bruggeman's avatar
Jorn Bruggeman committed
162
#  ifndef _FABM_HORIZONTAL_MASK_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
163
   _FABM_MASK_TYPE_, allocatable, target _DIMENSION_GLOBAL_ :: mask
Jorn Bruggeman's avatar
Jorn Bruggeman committed
164
#  endif
165 166 167
#  ifndef _FABM_MASKED_VALUE_
#    define _FABM_MASKED_VALUE_ _FABM_UNMASKED_VALUE_+1
#  endif
168 169 170
#  ifndef _FABM_UNMASKED_VALUE_
#    define _FABM_UNMASKED_VALUE_ _FABM_MASKED_VALUE_+1
#  endif
171 172
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
173 174
   integer :: interior_count
   integer :: horizontal_count
175

176
#if _FABM_BOTTOM_INDEX_==-1
Jorn Bruggeman's avatar
Jorn Bruggeman committed
177
   integer, allocatable, target _DIMENSION_GLOBAL_HORIZONTAL_ :: bottom_index
Jorn Bruggeman's avatar
Jorn Bruggeman committed
178 179
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
180 181 182
   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
183

Jorn Bruggeman's avatar
Jorn Bruggeman committed
184 185 186 187 188
   real(rke), allocatable _DIMENSION_SLICE_PLUS_1_            :: dy
   real(rke), allocatable _DIMENSION_SLICE_PLUS_1_            :: w
   real(rke), allocatable _DIMENSION_HORIZONTAL_SLICE_PLUS_1_ :: flux
   real(rke), allocatable _DIMENSION_HORIZONTAL_SLICE_PLUS_1_ :: sms_sf
   real(rke), allocatable _DIMENSION_HORIZONTAL_SLICE_PLUS_1_ :: sms_bt
Jorn Bruggeman's avatar
Jorn Bruggeman committed
189

Jorn Bruggeman's avatar
Jorn Bruggeman committed
190 191 192
   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
193

Jorn Bruggeman's avatar
Jorn Bruggeman committed
194
   class (type_fabm_model), pointer :: model
Jorn Bruggeman's avatar
Jorn Bruggeman committed
195

Jorn Bruggeman's avatar
Jorn Bruggeman committed
196
   class (type_test_model), pointer :: test_model
Jorn Bruggeman's avatar
Jorn Bruggeman committed
197

Jorn Bruggeman's avatar
Jorn Bruggeman committed
198
   integer :: domain_extent(_FABM_DIMENSION_COUNT_)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
199

Jorn Bruggeman's avatar
Jorn Bruggeman committed
200 201 202 203 204 205
   character(len=20) :: arg
   integer :: ivar
   integer :: i
   integer :: mode = 1
   integer :: ntest = -1
   logical :: no_mask = .false.
206

207
#if _FABM_DIMENSION_COUNT_>0
Jorn Bruggeman's avatar
Jorn Bruggeman committed
208
   i__ = 50
209 210
#endif
#if _FABM_DIMENSION_COUNT_>1
Jorn Bruggeman's avatar
Jorn Bruggeman committed
211
   j__ = 40
212 213
#endif
#if _FABM_DIMENSION_COUNT_>2
Jorn Bruggeman's avatar
Jorn Bruggeman committed
214
   k__ = 45
215 216
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
217 218 219 220
   ! Parse command line arguments
   call start_test('parsing command line arguments')
   i = 1
   do
221
      call get_command_argument(i, arg)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
222 223 224 225 226 227
      if (arg == '') exit
      select case (arg)
      case ('-s', '--simulate')
         mode = 2
      case ('--nomask')
         no_mask = .true.
228
#if _FABM_DIMENSION_COUNT_>0
Jorn Bruggeman's avatar
Jorn Bruggeman committed
229 230 231 232
      case ('--nx')
         i = i + 1
         call get_command_argument(i, arg)
         read (arg,*) i__
233 234
#endif
#if _FABM_DIMENSION_COUNT_>1
Jorn Bruggeman's avatar
Jorn Bruggeman committed
235 236 237 238
      case ('--ny')
         i = i + 1
         call get_command_argument(i, arg)
         read (arg,*) j__
239 240
#endif
#if _FABM_DIMENSION_COUNT_>2
Jorn Bruggeman's avatar
Jorn Bruggeman committed
241 242 243 244
      case ('--nz')
         i = i + 1
         call get_command_argument(i, arg)
         read (arg,*) k__
245
#endif
Jorn Bruggeman's avatar
Jorn Bruggeman committed
246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262
      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
263
      i = i + 1
Jorn Bruggeman's avatar
Jorn Bruggeman committed
264 265
   end do
   call report_test_result()
266

Jorn Bruggeman's avatar
Jorn Bruggeman committed
267
#if _FABM_DIMENSION_COUNT_>0
Jorn Bruggeman's avatar
Jorn Bruggeman committed
268 269
   domain_extent = (/ _LOCATION_ /)
   interior_count = product(domain_extent)
270
#else
Jorn Bruggeman's avatar
Jorn Bruggeman committed
271
   interior_count = 1
Jorn Bruggeman's avatar
Jorn Bruggeman committed
272 273
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
274 275 276 277 278 279 280
   ! Set defaults
   if (ntest == -1) then
      if (mode == 1) then
         ntest = 1
      else
         ntest = 50000000 / interior_count
      end if
281 282
   end if

283
#ifdef _FABM_DEPTH_DIMENSION_INDEX_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
284
   horizontal_count = interior_count / domain_extent(_FABM_DEPTH_DIMENSION_INDEX_)
285
#else
Jorn Bruggeman's avatar
Jorn Bruggeman committed
286
   horizontal_count = interior_count
287 288
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
289
#ifdef _INTERIOR_IS_VECTORIZED_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
290 291
   _START_ = 1
   _STOP_ = domain_extent(_FABM_VECTORIZED_DIMENSION_INDEX_)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
292 293
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314
   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()
315

Jorn Bruggeman's avatar
Jorn Bruggeman committed
316 317 318
   call start_test('initialize')
   call model%initialize()
   call report_test_result()
Jorn Bruggeman's avatar
Jorn Bruggeman committed
319

Jorn Bruggeman's avatar
Jorn Bruggeman committed
320 321 322
   ! ======================================================================
   ! Provide extents of the spatial domain.
   ! ======================================================================
Jorn Bruggeman's avatar
Jorn Bruggeman committed
323

Jorn Bruggeman's avatar
Jorn Bruggeman committed
324 325 326
   call start_test('set_domain')
   call model%set_domain(_LOCATION_)
   call report_test_result()
Jorn Bruggeman's avatar
Jorn Bruggeman committed
327

Jorn Bruggeman's avatar
Jorn Bruggeman committed
328 329 330
   ! ======================================================================
   ! Set up spatial mask.
   ! ======================================================================
Jorn Bruggeman's avatar
Jorn Bruggeman committed
331 332

#ifdef _HAS_MASK_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
333 334
   allocate(mask_hz _INDEX_HORIZONTAL_LOCATION_)
   call start_test('set_mask')
Jorn Bruggeman's avatar
Jorn Bruggeman committed
335
#  ifdef _FABM_HORIZONTAL_MASK_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
336
   call model%set_mask(mask_hz)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
337
#  else
Jorn Bruggeman's avatar
Jorn Bruggeman committed
338 339
   allocate(mask _INDEX_LOCATION_)
   call model%set_mask(mask, mask_hz)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
340
#  endif
Jorn Bruggeman's avatar
Jorn Bruggeman committed
341
   call report_test_result()
Jorn Bruggeman's avatar
Jorn Bruggeman committed
342 343
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
344 345 346
   ! ======================================================================
   ! Specify vertical indices of surface and bottom.
   ! ======================================================================
Jorn Bruggeman's avatar
Jorn Bruggeman committed
347

348
#if defined(_FABM_DEPTH_DIMENSION_INDEX_)&&_FABM_BOTTOM_INDEX_==-1
Jorn Bruggeman's avatar
Jorn Bruggeman committed
349 350 351 352
   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
353 354
#endif

355
   allocate(interior_state(_PREARG_LOCATION_ size(model%interior_state_variables)))
Jorn Bruggeman's avatar
Jorn Bruggeman committed
356 357 358 359 360 361 362 363
   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')
364
   do ivar = 1, size(model%interior_state_variables)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
365 366 367 368 369 370 371 372 373 374 375 376 377 378 379 380 381 382 383 384 385 386 387 388 389 390 391
      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')
392 393
       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
394 395 396
       call report_test_result()

       call start_test('link_horizontal_data')
397
       call model%link_horizontal_data(fabm_standard_variables%wind_speed, wind_speed)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
398 399 400 401 402 403 404 405 406 407 408 409 410
       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
411 412

#ifdef _FABM_VECTORIZED_DIMENSION_INDEX_
413 414
   allocate(dy(_START_:_STOP_, size(model%interior_state_variables)))
   allocate(w(_START_:_STOP_, size(model%interior_state_variables)))
Jorn Bruggeman's avatar
Jorn Bruggeman committed
415
#else
416 417
   allocate(dy(size(model%interior_state_variables)))
   allocate(w(size(model%interior_state_variables)))
Jorn Bruggeman's avatar
Jorn Bruggeman committed
418 419 420
#endif

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

Jorn Bruggeman's avatar
Jorn Bruggeman committed
430 431 432 433 434 435 436 437
   select case (mode)
   case (1)
      do i=1,ntest
         call test_update
      end do
   case(2)
      call simulate(ntest)
   end select
Jorn Bruggeman's avatar
Jorn Bruggeman committed
438

439 440 441 442
contains

   subroutine read_environment
      use yaml, only: yaml_parse => parse, yaml_error_length => error_length
443 444
      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
445 446 447

      integer, parameter :: yaml_unit = 100
      character(yaml_error_length) :: yaml_error
448
      class (type_node),pointer :: yaml_root
449
      type (type_yaml_key_value_pair), pointer :: yaml_pair
450
      real(rke) :: value
451 452
      logical :: success
      type type_input
453 454 455
         type (type_fabm_interior_variable_id)                :: interior_id
         type (type_fabm_horizontal_variable_id)              :: horizontal_id
         type (type_fabm_scalar_variable_id)                  :: scalar_id
456 457 458
         real(rke), allocatable _DIMENSION_GLOBAL_            :: interior_data
         real(rke), allocatable _DIMENSION_GLOBAL_HORIZONTAL_ :: horizontal_data
         real(rke)                                            :: scalar_data
459 460 461 462 463 464 465 466 467 468 469 470 471 472 473
      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))
474
                  value = node%to_real(0._yaml_real_kind, success)
475 476 477 478 479
                  if (.not. success) then
                     call driver%log_message('Cannot parse '//trim(node%string)//' as real.')
                     stop 2
                  end if
                  allocate(input)
480
                  input%interior_id = model%get_interior_variable_id(trim(yaml_pair%key))
481
                  if (model%is_variable_used(input%interior_id)) then
482 483 484 485 486
                      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))
487
                      if (model%is_variable_used(input%horizontal_id)) then
488 489 490 491 492
                         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))
493
                         if (model%is_variable_used(input%scalar_id)) then
494 495 496
                            input%scalar_data = value
                            call model%link_scalar(input%scalar_id, input%scalar_data)
                         else
497 498
                            call driver%log_message('WARNING: environment variable '//trim(yaml_pair%key) &
                               //' is not used by FABM model and will be ignored.')
499 500 501 502 503 504 505 506 507 508 509
                         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
510

Jorn Bruggeman's avatar
Jorn Bruggeman committed
511
   subroutine randomize_mask
Jorn Bruggeman's avatar
Jorn Bruggeman committed
512
      if (.not.no_mask) then
513
#if _FABM_BOTTOM_INDEX_==-1
514
      ! Depth index of bottom varies in the horizontal
515
      call random_number(tmp_hz)
516 517
#  ifdef _HAS_MASK_
      !  Pick random numbers between 0 (land) and maximum index
518
      bottom_index = floor(tmp_hz * (1 + domain_extent(_FABM_DEPTH_DIMENSION_INDEX_)))
519
#    ifdef _FABM_VERTICAL_BOTTOM_TO_SURFACE_
520 521
     ! Ensure invalid bottom indices [land points] are set such that vertical loops have 0 iterations.
     where (bottom_index == 0) bottom_index = domain_extent(_FABM_DEPTH_DIMENSION_INDEX_) + 1
522 523 524
#    endif
#  else
      ! Pick random numbers between 1 and maximum index
525
      bottom_index = 1 + floor(tmp_hz * domain_extent(_FABM_DEPTH_DIMENSION_INDEX_))
526 527 528
#  endif
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
529
#ifdef _HAS_MASK_
530 531 532 533
#  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_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
534
      where (tmp_hz>0.5_rke) mask_hz = _FABM_MASKED_VALUE_
535 536 537
#  else
      ! Apply random mask across interior domain (half of grid cells masked)
      call random_number(tmp)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
538
      mask = _FABM_UNMASKED_VALUE_
539
      where (tmp>0.5_rke) mask = _FABM_MASKED_VALUE_
540 541

#    if _FABM_BOTTOM_INDEX_==-1
542
      ! Bottom index varies in the horizontal. Ensure the bottom cell itself is unmasked, and anything deeper is masked.
543
      _BEGIN_GLOBAL_HORIZONTAL_LOOP_
544
         ! Valid bottom index - unmask associated cell, then mask all deeper ones
545
#      ifdef _FABM_VERTICAL_BOTTOM_TO_SURFACE_
546
         if (bottom_index _INDEX_HORIZONTAL_LOCATION_ <= domain_extent(_FABM_DEPTH_DIMENSION_INDEX_)) mask _INDEX_GLOBAL_VERTICAL_(bottom_index _INDEX_HORIZONTAL_LOCATION_) = _FABM_UNMASKED_VALUE_
547
         mask _INDEX_GLOBAL_VERTICAL_(:bottom_index _INDEX_HORIZONTAL_LOCATION_ - 1) = _FABM_MASKED_VALUE_
548
#      else
549
         if (bottom_index _INDEX_HORIZONTAL_LOCATION_ >= 1) mask _INDEX_GLOBAL_VERTICAL_(bottom_index _INDEX_HORIZONTAL_LOCATION_) = _FABM_UNMASKED_VALUE_
550
         mask _INDEX_GLOBAL_VERTICAL_(bottom_index _INDEX_HORIZONTAL_LOCATION_ + 1:) = _FABM_MASKED_VALUE_
551
#      endif
552 553 554 555 556 557 558 559 560 561 562 563 564 565 566 567
      _END_GLOBAL_HORIZONTAL_LOOP_
#    endif

      ! Infer horizontal mask (mask points that have all column layers masked)
      mask_hz = _FABM_UNMASKED_VALUE_
      where (.not.any(_IS_UNMASKED_(mask),_FABM_DEPTH_DIMENSION_INDEX_)) mask_hz = _FABM_MASKED_VALUE_

      ! 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
            mask _INDEX_GLOBAL_VERTICAL_(1) = _FABM_UNMASKED_VALUE_
#    if _FABM_BOTTOM_INDEX_!=-1
            mask _INDEX_GLOBAL_VERTICAL_(domain_extent(_FABM_DEPTH_DIMENSION_INDEX_)) = _FABM_UNMASKED_VALUE_
#    endif
         end if
      _END_GLOBAL_HORIZONTAL_LOOP_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
568 569 570 571 572 573 574 575 576 577 578 579 580 581 582 583 584 585 586 587 588 589 590 591 592 593 594 595 596 597 598
#  endif
#endif
      else
#if _FABM_BOTTOM_INDEX_==-1
#  ifdef _FABM_VERTICAL_BOTTOM_TO_SURFACE_
      bottom_index = 1
#  else
      bottom_index = domain_extent(_FABM_DEPTH_DIMENSION_INDEX_)
#  endif
#endif
#ifdef _HAS_MASK_
#  ifndef _FABM_HORIZONTAL_MASK_
      mask = _FABM_UNMASKED_VALUE_
#  endif
      mask_hz = _FABM_UNMASKED_VALUE_
#endif
      end if

      call count_active_points()
   end subroutine randomize_mask

   subroutine count_active_points()
#ifdef _HAS_MASK_
#  ifdef _FABM_HORIZONTAL_MASK_
      horizontal_count = count(_IS_UNMASKED_(mask_hz))
#    ifdef _FABM_DEPTH_DIMENSION_INDEX_
      interior_count = horizontal_count * domain_extent(_FABM_DEPTH_DIMENSION_INDEX_)
#    else
      interior_count = horizontal_count
#    endif
#  else
599 600
      horizontal_count = count(_IS_UNMASKED_(mask_hz))
      interior_count = count(_IS_UNMASKED_(mask))
601
#  endif
602
#elif _FABM_BOTTOM_INDEX_==-1
603
#  ifdef _FABM_VERTICAL_BOTTOM_TO_SURFACE_
604
      horizontal_count = count(bottom_index <= domain_extent(_FABM_DEPTH_DIMENSION_INDEX_))
605 606
      interior_count = sum(domain_extent(_FABM_DEPTH_DIMENSION_INDEX_) - bottom_index + 1)
#  else
607
      horizontal_count = count(bottom_index >= 1)
608 609
      interior_count = sum(bottom_index)
#  endif
Jorn Bruggeman's avatar
Jorn Bruggeman committed
610
#endif
Jorn Bruggeman's avatar
Jorn Bruggeman committed
611
   end subroutine
612

613 614
   subroutine simulate(n)
      integer, intent(in) :: n
615
      real(rke) :: time_begin, time_end
616 617 618 619 620 621 622 623 624
      integer :: nseed
      integer, allocatable :: seed(:)

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

      call randomize_mask
625

626
      call start_test('initialize_interior_state')
627
      _BEGIN_OUTER_INTERIOR_LOOP_
628
         call model%initialize_interior_state(_ARG_INTERIOR_IN_)
629 630 631
      _END_OUTER_INTERIOR_LOOP_
      call report_test_result()

632
      call start_test('initialize_bottom_state')
633
      _BEGIN_OUTER_HORIZONTAL_LOOP_
634
         call model%initialize_bottom_state(_ARG_HORIZONTAL_IN_)
635 636 637
      _END_OUTER_HORIZONTAL_LOOP_
      call report_test_result()

638
      call start_test('initialize_surface_state')
639
      _BEGIN_OUTER_HORIZONTAL_LOOP_
640
         call model%initialize_surface_state(_ARG_HORIZONTAL_IN_)
641 642 643
      _END_OUTER_HORIZONTAL_LOOP_
      call report_test_result()

644
      write (*,'(a,i0,a)') 'Simulating with ', interior_count, ' wet cells...'
645 646 647

      call cpu_time(time_begin)

648 649
      do i = 1, n
         call model%prepare_inputs()
650 651 652 653

         _BEGIN_OUTER_HORIZONTAL_LOOP_
            flux = 0
            sms_bt = 0
Jorn Bruggeman's avatar
Jorn Bruggeman committed
654
            call model%get_bottom_sources(_PREARG_HORIZONTAL_IN_ flux, sms_bt)
655 656 657 658 659
         _END_OUTER_HORIZONTAL_LOOP_

         _BEGIN_OUTER_HORIZONTAL_LOOP_
            flux = 0
            sms_sf = 0
Jorn Bruggeman's avatar
Jorn Bruggeman committed
660
            call model%get_surface_sources(_PREARG_HORIZONTAL_IN_ flux, sms_sf)
661 662 663 664
         _END_OUTER_HORIZONTAL_LOOP_

         _BEGIN_OUTER_INTERIOR_LOOP_
            dy = 0
Jorn Bruggeman's avatar
Jorn Bruggeman committed
665
            call model%get_interior_sources(_PREARG_INTERIOR_IN_ dy)
666
         _END_OUTER_INTERIOR_LOOP_
667

668
         call model%finalize_outputs()
669

670
         if (mod(i, 100) == 0) write (*,'(i0,a)') int(100*i/real(n, rke)), ' % complete'
671 672 673 674
      end do

      call cpu_time(time_end)

675 676
      write (*,'(a)') 'Simulation complete.'
      write (*,'(a,f8.3,a)') 'Total time: ', time_end - time_begin, ' s'
677 678
   end subroutine

Jorn Bruggeman's avatar
Jorn Bruggeman committed
679
   subroutine test_update
680 681 682
      real(rke), pointer _DIMENSION_GLOBAL_            :: pdata
      real(rke), pointer _DIMENSION_GLOBAL_HORIZONTAL_ :: pdata_hz
      logical                                          :: valid
Jorn Bruggeman's avatar
Jorn Bruggeman committed
683

Jorn Bruggeman's avatar
Jorn Bruggeman committed
684 685 686 687 688
      call randomize_mask

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

690
      call start_test('initialize_interior_state')
691
      _BEGIN_OUTER_INTERIOR_LOOP_
692
         call model%initialize_interior_state(_ARG_INTERIOR_IN_)
693
      _END_OUTER_INTERIOR_LOOP_
694
      do ivar = 1, size(model%interior_state_variables)
695
         call check_interior(interior_state(_PREARG_LOCATION_DIMENSIONS_ ivar), &
696
            model%interior_state_variables(ivar)%missing_value, ivar+interior_state_offset+1._rke)
697
      end do
Jorn Bruggeman's avatar
Jorn Bruggeman committed
698
      call report_test_result()
Jorn Bruggeman's avatar
Jorn Bruggeman committed
699

700
      call start_test('initialize_bottom_state')
701
      _BEGIN_OUTER_HORIZONTAL_LOOP_
702
         call model%initialize_bottom_state(_ARG_HORIZONTAL_IN_)
703
      _END_OUTER_HORIZONTAL_LOOP_
704
      do ivar = 1, size(model%bottom_state_variables)
705 706
         call check_horizontal(bottom_state(_PREARG_HORIZONTAL_LOCATION_DIMENSIONS_ ivar), &
            model%bottom_state_variables(ivar)%missing_value, ivar+bottom_state_offset+1._rke)
707
      end do
Jorn Bruggeman's avatar
Jorn Bruggeman committed
708
      call report_test_result()
Jorn Bruggeman's avatar
Jorn Bruggeman committed
709

710
      call start_test('initialize_surface_state')
711
      _BEGIN_OUTER_HORIZONTAL_LOOP_
712
         call model%initialize_surface_state(_ARG_HORIZONTAL_IN_)
713
      _END_OUTER_HORIZONTAL_LOOP_
714
      do ivar = 1, size(model%surface_state_variables)
715 716
         call check_horizontal(surface_state(_PREARG_HORIZONTAL_LOCATION_DIMENSIONS_ ivar), &
            model%surface_state_variables(ivar)%missing_value, ivar+surface_state_offset+1._rke)
717
      end do
Jorn Bruggeman's avatar
Jorn Bruggeman committed
718
      call report_test_result()
Jorn Bruggeman's avatar
Jorn Bruggeman committed
719

Jorn Bruggeman's avatar
Jorn Bruggeman committed
720 721 722 723
      ! ======================================================================
      ! Initialize environmental dependencies
      ! ======================================================================

724 725
      temperature = 1 + interior_dependency_offset
      call apply_mask_3d(temperature, -999._rke - interior_dependency_offset)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
726

727 728
      wind_speed = 1 + horizontal_dependency_offset
      call apply_mask_2d(wind_speed, -999._rke - horizontal_dependency_offset)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
729 730 731

#ifdef _FABM_DEPTH_DIMENSION_INDEX_
      ! Model has depth dimension: make sure depth varies from 0 at the surface till 1 at the bottom
732
      _BEGIN_GLOBAL_LOOP_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
733
#  ifdef _FABM_VERTICAL_BOTTOM_TO_SURFACE_
734 735 736 737
#    if _FABM_BOTTOM_INDEX_==-1
         if (bottom_index _INDEX_HORIZONTAL_LOCATION_ == domain_extent(_FABM_DEPTH_DIMENSION_INDEX_)) then
            depth _INDEX_LOCATION_ = 2
         else
738
            depth _INDEX_LOCATION_ = real(domain_extent(_FABM_DEPTH_DIMENSION_INDEX_)-_VERTICAL_ITERATOR_,rke)/(domain_extent(_FABM_DEPTH_DIMENSION_INDEX_)-bottom_index _INDEX_HORIZONTAL_LOCATION_)
739 740
         end if
#    else
741
         depth _INDEX_LOCATION_ = real(domain_extent(_FABM_DEPTH_DIMENSION_INDEX_)-_VERTICAL_ITERATOR_,rke)/(domain_extent(_FABM_DEPTH_DIMENSION_INDEX_)-1)
742
#    endif
Jorn Bruggeman's avatar
Jorn Bruggeman committed
743
#  else
744
#    if _FABM_BOTTOM_INDEX_==-1
745 746 747
         if (bottom_index _INDEX_HORIZONTAL_LOCATION_==1) then
            depth _INDEX_LOCATION_ = 2
         else
748
            depth _INDEX_LOCATION_ = real(_VERTICAL_ITERATOR_-1,rke)/(bottom_index _INDEX_HORIZONTAL_LOCATION_-1)
749
         end if
750
#    else
751
         depth _INDEX_LOCATION_ = real(_VERTICAL_ITERATOR_-1,rke)/(domain_extent(_FABM_DEPTH_DIMENSION_INDEX_)-1)
752
#    endif
Jorn Bruggeman's avatar
Jorn Bruggeman committed
753
#  endif
754
      _END_GLOBAL_LOOP_
Jorn Bruggeman's avatar
Jorn Bruggeman committed
755 756 757 758
#else
      ! No depth dimension
      depth = 2
#endif
759
      call apply_mask_3d(depth, -999._rke - interior_dependency_offset)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
760

761
      do ivar = 1, size(model%interior_state_variables)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
762
         interior_state(_PREARG_LOCATION_DIMENSIONS_ ivar) = ivar+interior_state_offset
763
         call apply_mask_3d(interior_state(_PREARG_LOCATION_DIMENSIONS_ ivar),model%interior_state_variables(ivar)%missing_value)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
764
      end do
765
      do ivar = 1, size(model%surface_state_variables)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
766
         surface_state(_PREARG_HORIZONTAL_LOCATION_DIMENSIONS_ ivar) = ivar+surface_state_offset
767
         call apply_mask_2d(surface_state(_PREARG_HORIZONTAL_LOCATION_DIMENSIONS_ ivar),model%surface_state_variables(ivar)%missing_value)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
768
      end do
769
      do ivar = 1, size(model%bottom_state_variables)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
770
         bottom_state(_PREARG_HORIZONTAL_LOCATION_DIMENSIONS_ ivar) = ivar+bottom_state_offset
771
         call apply_mask_2d(bottom_state(_PREARG_HORIZONTAL_LOCATION_DIMENSIONS_ ivar),model%bottom_state_variables(ivar)%missing_value)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
772 773
      end do

Jorn Bruggeman's avatar
Jorn Bruggeman committed
774 775 776 777 778 779
      ! ======================================================================
      ! Preprocessing
      ! ======================================================================

      column_loop_count = 0

780 781
      call start_test('prepare_inputs')
      call model%prepare_inputs()
Jorn Bruggeman's avatar
Jorn Bruggeman committed
782 783
      call report_test_result()

Jorn Bruggeman's avatar
Jorn Bruggeman committed
784 785 786 787
      ! ======================================================================
      ! Retrieve source terms of interior state variables.
      ! ======================================================================

Jorn Bruggeman's avatar
Jorn Bruggeman committed
788
      call start_test('get_interior_sources')
789
      interior_loop_count = 0
790 791
      _BEGIN_OUTER_INTERIOR_LOOP_
         dy = 0
792 793 794
#if _FABM_BOTTOM_INDEX_==-1 && !defined(_HAS_MASK_) && _FABM_VECTORIZED_DIMENSION_INDEX_==_FABM_DEPTH_DIMENSION_INDEX_ && defined(_FABM_DEPTH_DIMENSION_INDEX_)
         ! We are looping over depth, but as we have a non-constant bottom index (yet no mask), we need to skip everything below bottom
         if (bottom_index _INDEX_HORIZONTAL_LOCATION_ >= 1 .and. bottom_index _INDEX_HORIZONTAL_LOCATION_ <= domain_extent(_FABM_DEPTH_DIMENSION_INDEX_)) &
795
            call model%get_interior_sources(_IMIN_, _IMAX_ _ARG_INTERIOR_FIXED_LOCATION_, dy(_IMIN_:_IMAX_, :))
796
#else
Jorn Bruggeman's avatar
Jorn Bruggeman committed
797
         call model%get_interior_sources(_PREARG_INTERIOR_IN_ dy)
798
#endif
799
         do ivar = 1, size(model%interior_state_variables)
800
            call check_interior_slice_plus_1(dy, ivar, 0.0_rke, -real(ivar + interior_state_offset, rke) _POSTARG_INTERIOR_IN_)
801 802
         end do
      _END_OUTER_INTERIOR_LOOP_
803 804
      call assert(interior_loop_count == interior_count, 'get_interior_sources', &
         'call count does not match number of (unmasked) interior points')
Jorn Bruggeman's avatar
Jorn Bruggeman committed
805

806 807
      do ivar = 1, size(model%interior_diagnostic_variables)
         if (model%interior_diagnostic_variables(ivar)%save .and. model%interior_diagnostic_variables(ivar)%target%source == source_do) then
808
            pdata => model%get_interior_diagnostic_data(ivar)
809 810
            call check_interior(pdata, model%interior_diagnostic_variables(ivar)%missing_value, &
               -model%interior_diagnostic_variables(ivar)%missing_value)
Jorn Bruggeman's avatar
Jorn Bruggeman committed
811 812 813
         end if
      end do

Jorn Bruggeman's avatar
Jorn Bruggeman committed
814 815 816 817 818 819
      call report_test_result()

      ! ======================================================================
      ! Retrieve surface fluxes of interior state variables, source terms of surface-associated state variables.
      ! ======================================================================

820 821
#if _FABM_BOTTOM_INDEX_==-1
      _BEGIN_GLOBAL_HORIZONTAL_LOOP_
822 823 824 825 826
#  ifdef _FABM_VERTICAL_BOTTOM_TO_SURFACE_
         if (bottom_index _INDEX_HORIZONTAL_LOCATION_ == domain_extent(_FABM_DEPTH_DIMENSION_INDEX_)) depth _INDEX_GLOBAL_VERTICAL_(domain_extent(_FABM_DEPTH_DIMENSION_INDEX_)) = 0
#  else
         if (bottom_index _INDEX_HORIZONTAL_LOCATION_ == 1) depth _INDEX_GLOBAL_VERTICAL_(1) = 0
#  endif
827 828 829
      _END_GLOBAL_HORIZONTAL_LOOP_
#endif

Jorn Bruggeman's avatar
Jorn Bruggeman committed
830
      call start_test('get_surface_sources')
831
      surface_loop_count = 0
832
      _BEGIN_OUTER_HORIZONTAL_LOOP_
833 834 835
#if _FABM_BOTTOM_INDEX_==-1 && !defined(_HAS_MASK_)
         if (bottom_index _INDEX_HORIZONTAL_LOCATION_ >= 1 .and. bottom_index _INDEX_HORIZONTAL_LOCATION_ <= domain_extent(_FABM_DEPTH_DIMENSION_INDEX_)) then
#endif
836 837
         flux = 0
         sms_sf = 0
Jorn Bruggeman's avatar
Jorn Bruggeman committed
838
         call model%get_surface_sources(_PREARG_HORIZONTAL_IN_ flux, sms_sf)
839
         do ivar = 1, size(model%interior_state_variables)