SABER
bump_interpolation_mod.F90
Go to the documentation of this file.
1 ! (C) Copyright 2019-2020 UCAR
2 !
3 ! This software is licensed under the terms of the Apache Licence Version 2.0
4 ! which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
5 
6 !> \brief Bump Interpolation module
7 !!
8 !!
9 !! \date Jan, 2020: Created by M. Miesch (JCSDA/UCAR) based on previously exising type_obsop.F90
10 !! written by Benjamin Menetrier (JCSDA/UCAR, CERFACS, METEO-FRANCE, IRIT)
11 !!
13 
14 use atlas_module
15 use atlas_metadata_module, only: atlas_metadata
16 use fckit_configuration_module, only : fckit_configuration
17 use fckit_log_module, only: fckit_log
18 use fckit_mpi_module, only: fckit_mpi_comm, fckit_mpi_sum,fckit_mpi_min,&
19  fckit_mpi_max,fckit_mpi_status
20 use netcdf
21 use tools_atlas
22 use tools_const, only: pi,deg2rad,rad2deg
24 use tools_kinds, only: kind_real
25 use tools_qsort, only: qsort
26 use tools_repro, only: rth
27 use type_bump, only : bump_type
28 use type_com, only: com_type
30 use type_geom, only: geom_type
31 use type_linop, only: linop_type
32 use type_mesh, only: mesh_type
33 use type_mpl, only: mpl_type
34 use type_nam, only: nam_type
35 use type_rng, only: rng_type
36 
37 implicit none
38 private
39 
40 public :: bump_interpolator
42 
43 integer, parameter :: max_string = 1024
44 
45 ! ------------------------------------------------------------------------------
46 !!
47 !! bint = shorthand for bump interpolator
48 !!
50 private
51  type(bump_type), public :: bump
52 
53  !> Grid information
54  type(geom_type) :: outgeom !<< output grid - unstructured
55  type(atlas_functionspace) :: in_funcspace !<< atlas functionspace for input grid
56  type(atlas_functionspace) :: out_funcspace !<< atlas functionspace for output grid
57 
58  !> Number of points
59  integer :: nc0b !<< Halo B size
60  integer, public :: nlev !<< number of levels
61 
62  integer, public :: nout !<< global number of output grid points
63  integer, public :: nout_local !<< local number of output grid points
64 
65  !> Interpolation data (operator)
66  type(linop_type) :: h
67 
68  !> Communication data
69  type(com_type) :: com
70 
71 contains
72  private
73 
74  procedure, private :: driver => bint_driver
75 
76  procedure, public :: init => bint_init
77 
78  procedure, public :: apply => bint_apply
79  procedure, public :: apply_ad => bint_apply_ad
80 
81  procedure :: deallocate_outgrid => bint_deallocate_outgrid
82  procedure, public :: delete => bint_delete
83 
84  ! low-level apply methods
85  procedure, private :: apply_interp, apply_interp_ad
86  final :: dummy
87 
88 end type bump_interpolator
89 
90 ! ------------------------------------------------------------------------------
91 !> Registry for bump_interpolator objects
92 
93 #define LISTED_TYPE bump_interpolator
94 
95 !> Linked list interface - defines registry_t type
96 #include "saber/util/linkedList_i.f"
97 
98 !> Global registry
99 type(registry_t) :: bump_interpolator_registry
100 
101 ! ------------------------------------------------------------------------------
102 
103 contains
104 !-------------------------------------------------------------------------------
105 !> Linked list implementation
106 #include "saber/util/linkedList_c.f"
107 
108 ! ------------------------------------------------------------------------------
109 !> Initialize interpolation object
110 !!
111 !! The input and output fields are atlas_FieldSet objects that are assumed
112 !! to be created from atlas functionspaces. So, they have the grid and
113 !! mesh information built in.
114 !!
115 !! \param[in] config = configuration
116 !!
117 !! \param[in] in_afs = This is the input grid, rendered as an atlas
118 !! functionspace, so it includes information about the mesh (parallel
119 !! connectivity) as well.
120 !!
121 !! \param[in] out_afs = This is the output grid, also represented as
122 !! an altas_functionspace.
123 !!
124 !! \param[in] masks = This contains metadata needed for the interpolation.
125 !! It is rendered as an atlas FieldSet with the following named fields.
126 !! Each of these named fields is optional; if omitted default values will
127 !! be provided
128 !! * area : cell area
129 !! * vunit : vertical unit
130 !! * gmask : geometry mask
131 !! * smask : sampling mask
132 !!
133 subroutine bint_init(self, config, comm, in_funcspace, out_funcspace, masks)
134  use fckit_configuration_module, only : fckit_configuration
135  use iso_c_binding, only : c_char
136  class(bump_interpolator), intent(inout) :: self
137  class(fckit_configuration), intent(in) :: config
138  type(fckit_mpi_comm) , intent(in) :: comm
139  class(atlas_functionspace), intent(in) :: in_funcspace
140  class(atlas_functionspace), intent(in) :: out_funcspace
141  type(fieldset_type) , intent(in), optional :: masks
142 
143  ! local variables
144  type(fckit_configuration) :: bump_config
145  integer :: msvali, j
146  real(kind_real) :: msvalr
147  integer, allocatable :: levels(:)
148  character(len=max_string) :: msg
149  character(len=max_string) :: myname = "saber::interpolation::bump_interpolation_mod::bint_init "
150 
151  !--------------------------------------------------------------------------------
152  ! set bump namelist parameters.
153 
154  call config%get_or_die("bump",bump_config)
155 
156  ! bump defaults
157  call self%bump%nam%init(comm%size()) ! set up defaults
158 
159  ! modified defaults for bump interpolator
160  self%bump%nam%default_seed = .true.
161  self%bump%nam%prefix = 'bump_interpolator'
162  self%bump%nam%datadir = 'testdata'
163 
164  ! optionally override defaults with config
165  call self%bump%nam%from_conf(bump_config)
166 
167  ! optional missing values for integers and reals
168  msvali = -999
169  msvalr = -999.0
170  If (config%has("missingvalue_int")) &
171  call config%get_or_die("missingvalue_int",msvali)
172  If (config%has("missingvalue_real")) &
173  call config%get_or_die("missingvalue_real",msvalr)
174 
175  ! save these for future use
176  self%in_funcspace = in_funcspace
177  self%out_funcspace = out_funcspace
178 
179  !--------------------------------------------------------------------------------
180  ! determine the number of vertical levels
181 
182  self%nlev = 1
183  If (config%has("nlevels")) call config%get_or_die("nlevels",self%nlev)
184 
185  self%bump%nam%nl = self%nlev
186 
187  allocate(levels(1:self%nlev))
188  if (config%has("levels")) then
189  call config%get_or_die("levels",levels)
190  else
191  do j = 1, self%nlev
192  levels(j) = j
193  enddo
194  endif
195  self%bump%nam%levs(1:self%nlev) = levels
196  deallocate(levels)
197 
198  !--------------------------------------------------------------------------------
199  ! Initialize BUMP
200  ! ---------------
201 
202  ! clear any pre-existing data
203  call self%delete
204 
205  ! basic bump setup
206  if (present(masks)) then
207  call self%bump%setup(self%bump%mpl%f_comm, in_funcspace, masks, &
208  msvali=msvali, msvalr=msvalr)
209  else
210  call self%bump%setup(self%bump%mpl%f_comm, in_funcspace, &
211  msvali=msvali, msvalr=msvalr)
212  endif
213 
214  !--------------------------------------------------------------------------------
215  ! initialize output grid
216 
217  call fckit_log%info('-------------------------------------------------------------------')
218  call fckit_log%info('--- Initialize output grid')
219 
220  self%outgeom%nl0 = self%nlev
221 
222  ! unpack output grid from output functionspace
223  call self%outgeom%from_atlas(self%bump%mpl, out_funcspace)
224  self%nout_local = size(self%outgeom%lon_mga)
225 
226  !--------------------------------------------------------------------------------
227  ! Run basic BUMP drivers
228  ! The only thing here that may be needed is to possibly run the nicas driver
229  ! after we get things to work, try commenting this out and see if it still works
230  call self%bump%run_drivers
231 
232  ! Run interpolation driver
233  call fckit_log%info('-------------------------------------------------------------------')
234  call fckit_log%info('--- Run bump_interpolation driver')
235  call self%driver(self%bump%mpl,self%bump%rng,self%bump%nam,self%bump%geom)
236  if (self%bump%nam%default_seed) call self%bump%rng%reseed(self%bump%mpl)
237 
238  !--------------------------------------------------------------------------------
239  ! more initializations and checks that the setup is correct
240 
241  if (self%nout < 1) then
242  write(msg,'(a)') trim(myname) // ": " // "Output grid not properly defined"
243  call abor1_ftn(msg)
244  endif
245 
246  if (self%nlev /= self%bump%geom%nl0) then
247  write(msg,'(a)') trim(myname) // ": " // &
248  "number of levels does not match bump geom"
249  call abor1_ftn(msg)
250  endif
251 
252  if (self%nlev /= self%outgeom%nl0) then
253  write(msg,'(a)') trim(myname) // ": " // &
254  "ERROR: Output grid has different number of levels than input grid!"
255  call abor1_ftn(msg)
256  endif
257 
258  !--------------------------------------------------------------------------------
259 
260 end subroutine bint_init
261 
262 ! ------------------------------------------------------------------------------
263 ! ------------------------------------------------------------------------------
264 !> Initialize bump to perform interpolation.
265 !!
266 !----------------------------------------------------------------------
267 subroutine bint_driver(self,mpl,rng,nam,geom)
268  class(bump_interpolator), intent(inout) :: self
269  type(mpl_type),intent(inout) :: mpl !<< MPI data
270  type(rng_type),intent(inout) :: rng !<< Random number generator
271  type(nam_type),intent(in) :: nam !<< Namelist
272  type(geom_type),intent(in) :: geom !<< Geometry
273 
274  ! Local variables
275  integer :: iouta,iproc,i_s,ic0,ic0u,jc0u,ic0b,ic0a,nouta_eff
276  integer :: nout_eff,nn_index(1),proc_to_nouta(mpl%nproc),proc_to_nouta_eff(mpl%nproc)
277  integer :: c0u_to_c0b(geom%nc0u)
278  integer,allocatable :: c0b_to_c0(:)
279  real(kind_real) :: nn_dist(1),N_max,C_max
280  logical :: maskouta(self%nout_local),lcheck_nc0b(geom%nc0)
281  character(len=1024),parameter :: subr = 'bint_driver'
282 
283  ! Check that universe is global
284  if (any(.not.geom%myuniverse)) call mpl%abort(subr,'universe should be global for interpolation')
285 
286  ! Check whether output grid points are inside the mesh
287  if (self%nout_local > 0) then
288  do iouta=1,self%nout_local
289  call geom%mesh_c0u%inside(mpl,self%outgeom%lon_mga(iouta),self%outgeom%lat_mga(iouta),maskouta(iouta))
290  if (.not.maskouta(iouta)) then
291  ! Check for very close points
292  call geom%tree_c0u%find_nearest_neighbors(self%outgeom%lon_mga(iouta),self%outgeom%lat_mga(iouta), &
293  & 1,nn_index,nn_dist)
294  if (nn_dist(1)<rth) maskouta(iouta) = .true.
295  end if
296  end do
297  nouta_eff = count(maskouta)
298  else
299  nouta_eff = 0
300  endif
301 
302  ! Get global number of output grid points
303  call mpl%f_comm%allgather(self%nout_local,proc_to_nouta)
304  call mpl%f_comm%allgather(nouta_eff,proc_to_nouta_eff)
305  self%nout = sum(proc_to_nouta)
306  nout_eff = sum(proc_to_nouta_eff)
307 
308  ! Print input
309  write(mpl%info,'(a7,a)') '','Number of points in output grid / valid points per MPI task:'
310  call mpl%flush
311  do iproc=1,mpl%nproc
312  write(mpl%info,'(a10,a,i3,a,i8,a,i8)') '','Task ',iproc,': ', &
313  proc_to_nouta(iproc),' / ',proc_to_nouta_eff(iproc)
314  call mpl%flush
315  end do
316  write(mpl%info,'(a10,a,i8,a,i8)') '','Total : ',self%nout,' / ',nout_eff
317  call mpl%flush
318 
319  ! Compute interpolation
320  self%h%prefix = 'o'
321  write(mpl%info,'(a7,a)') '','Single level:'
322  call mpl%flush
323  call self%h%interp(mpl,rng,nam,geom,0,geom%nc0u,geom%lon_c0u,geom%lat_c0u,&
324  geom%gmask_hor_c0u,self%nout_local,self%outgeom%lon_mga,&
325  self%outgeom%lat_mga,maskouta,10)
326 
327  ! Define halo B
328  lcheck_nc0b = .false.
329  do ic0a=1,geom%nc0a
330  ic0u = geom%c0a_to_c0u(ic0a)
331  lcheck_nc0b(ic0u) = .true.
332  end do
333  do iouta=1,self%nout_local
334  do i_s=1,self%h%n_s
335  jc0u = self%h%col(i_s)
336  lcheck_nc0b(jc0u) = .true.
337  end do
338  end do
339  self%nc0b = count(lcheck_nc0b)
340 
341  ! Allocation
342  allocate(c0b_to_c0(self%nc0b))
343 
344  ! Global-local conversion for halo B
345  c0u_to_c0b = mpl%msv%vali
346  ic0b = 0
347  do ic0u=1,geom%nc0u
348  if (lcheck_nc0b(ic0u)) then
349  ic0b = ic0b+1
350  ic0 = geom%c0u_to_c0(ic0u)
351  c0b_to_c0(ic0b) = ic0
352  c0u_to_c0b(ic0u) = ic0b
353  end if
354  end do
355 
356  ! Local interpolation source
357  self%h%n_src = self%nc0b
358  do i_s=1,self%h%n_s
359  self%h%col(i_s) = c0u_to_c0b(self%h%col(i_s))
360  end do
361 
362  ! Setup communications
363  call self%com%setup(mpl,'com',geom%nc0a,self%nc0b,geom%nc0,geom%c0a_to_c0,c0b_to_c0)
364 
365  ! Compute scores
366  if (nout_eff > 0) then
367  call mpl%f_comm%allreduce(real(self%com%nhalo,kind_real),c_max,fckit_mpi_max())
368  c_max = c_max/(3.0*real(nout_eff,kind_real)/real(mpl%nproc,kind_real))
369  n_max = real(maxval(proc_to_nouta_eff),kind_real)/(real(nout_eff,kind_real)/real(mpl%nproc,kind_real))
370 
371  ! Print results
372  write(mpl%info,'(a7,a,f5.1,a)') '','Output grid repartition imbalance: ',100.0*real(maxval(proc_to_nouta_eff) &
373  & -minval(proc_to_nouta_eff),kind_real)/(real(sum(proc_to_nouta_eff),kind_real)/real(mpl%nproc,kind_real)),' %'
374  call mpl%flush
375  write(mpl%info,'(a7,a,i8,a,i8,a,i8)') '','Number of grid points / halo size / number of received values: ', &
376  & self%com%nred,' / ',self%com%next,' / ',self%com%nhalo
377  call mpl%flush
378  write(mpl%info,'(a7,a,f10.2,a,f10.2)') '','Scores (N_max / C_max):',n_max,' / ',c_max
379  call mpl%flush
380  endif
381 
382  ! Release memory
383  deallocate(c0b_to_c0)
384 
385 end subroutine bint_driver
386 
387 ! ------------------------------------------------------------------------------
388 ! ------------------------------------------------------------------------------
389 !> Apply interpolation
390 !!
391 !! \param[in] instate = input fields represented as an atlas fieldset,
392 !! created from a functionspace
393 !!
394 !! \param[out] outstate = output fields represented as an atlas fieldset.
395 !! If the fields that constitudte the fieldset are not already allocated by
396 !! the caller, then they will be created and allocated by this method.
397 !! So, the user can optionally pass this routine an empty output fieldset.
398 !!
399 subroutine bint_apply(self, infields, outfields)
400  class(bump_interpolator), intent(inout) :: self
401  type(fieldset_type) , intent(in) :: infields
402  type(fieldset_type) , intent(inout) :: outfields
403 
404  ! local variables
405  type(atlas_field) :: infield, outfield
406  real(kind_real), allocatable :: infld_mga(:,:), infld_c0a(:,:)
407  real(kind_real), allocatable :: outfld(:,:)
408  integer :: ifield
409  character(len=max_string) :: fieldname
410 
411  ! allocate bump arrays
412  allocate(infld_mga(self%bump%geom%nmga,self%nlev))
413  allocate(outfld(self%nout_local,self%nlev))
414 
415  if (.not.self%bump%geom%same_grid) then
416  allocate(infld_c0a(self%bump%geom%nc0a, self%nlev))
417  endif
418 
419  do ifield = 1, infields%size()
420 
421  infield = infields%field(ifield)
422  fieldname = infield%name()
423 
424  !--------------------------------------------
425  ! allocate output field if necessary
426 
427  if (.not. outfields%has_field(fieldname)) then
428  outfield = self%out_funcspace%create_field(name=fieldname, &
429  kind=atlas_real(kind_real),levels=self%nlev)
430  else
431  outfield = outfields%field(name=fieldname)
432  endif
433 
434  !--------------------------------------------
435  ! compute interpolation
436 
437  ! atlas field to fortran array
438  call field_to_array(self%bump%mpl, infield, infld_mga)
439 
440  if (self%bump%geom%same_grid) then
441  call self%apply_interp(infld_mga,outfld)
442  else
443  ! Model grid to subset Sc0
444  infld_c0a = 0.0_kind_real
445  call self%bump%geom%copy_mga_to_c0a(self%bump%mpl,infld_mga,infld_c0a)
446 
447  ! Apply interpolation
448  call self%apply_interp(infld_c0a,outfld)
449  end if
450 
451  ! fortran array to atlas field
452  call field_from_array(self%bump%mpl, outfld, outfield)
453 
454  !--------------------------------------------
455  ! Add output field to output fields
456  if (.not. outfields%has_field(fieldname)) then
457  call outfields%add(outfield)
458  endif
459 
460  ! release pointers
461  call infield%final()
462  call outfield%final()
463 
464  enddo
465 
466  ! clean up
467  if (allocated(infld_mga)) deallocate(infld_mga)
468  if (allocated(infld_c0a)) deallocate(infld_c0a)
469  if (allocated(outfld)) deallocate(outfld)
470 
471 end subroutine bint_apply
472 
473 !----------------------------------------------------------------------
474 !> Subroutine: apply_interp
475 !! Purpose: low-level routine to apply the interpolation to a single
476 !! field on a single level
477 !! \param[in] infield: input field
478 !! \param[out] outfield: output field
479 !!
480 subroutine apply_interp(self,infield,outfield)
481  class(bump_interpolator), intent(inout) :: self
482  real(kind_real) , intent(in) :: infield(self%bump%geom%nc0a,self%nlev)
483  real(kind_real) , intent(out) :: outfield(self%nout_local,self%nlev)
484 
485  integer :: ilev
486  real(kind_real), allocatable :: infield_ext(:,:)
487 
488  allocate(infield_ext(self%nc0b,self%nlev))
489 
490  ! Halo extension
491  call self%com%ext(self%bump%mpl, self%nlev, infield ,infield_ext)
492 
493  if (self%nout_local > 0) then
494  ! Horizontal interpolation
495  !$omp parallel do schedule(static) private(ilev)
496  do ilev = 1, self%nlev
497  call self%h%apply(self%bump%mpl,infield_ext(:,ilev),outfield(:,ilev))
498  end do
499  !$omp end parallel do
500  end if
501 
502  deallocate(infield_ext)
503 
504 end subroutine apply_interp
505 
506 !----------------------------------------------------------------------
507 !----------------------------------------------------------------------
508 !> Apply interpolator operator adjoint
509 !!
510 !! \param[in] fields_outgrid = These are the fields on the second grid,
511 !! i.e. the target grid of the original interpolation. For the adjoint,
512 !! these fields are treated as an input.
513 !!
514 !! \param[inout] fields_ingrid = These are the fields defined on the
515 !! first grid, i.e. the source grid of the original interpolation. For
516 !! the adjoint, these are treated as an output. The caller can optionally
517 !! pass this arguement as an empty fieldset and the routine will create
518 !! and allocate each component of the fieldset. Or, if the field components
519 !! of the fieldset are already allocated by the caller, then this routine
520 !! will merely replace the field values with the result of the computation.
521 !!
522 !----------------------------------------------------------------------
523 subroutine bint_apply_ad(self, fields_outgrid, fields_ingrid)
524  class(bump_interpolator), intent(inout) :: self
525  type(fieldset_type) , intent(in) :: fields_outgrid
526  type(fieldset_type) , intent(inout) :: fields_ingrid
527 
528  ! Local variables
529  type(atlas_field) :: field_ingrid, field_outgrid
530  real(kind_real), allocatable :: fld_ingrid_mga(:,:), fld_ingrid_c0a(:,:)
531  real(kind_real), allocatable :: fld_outgrid(:,:)
532  integer :: ifield
533  character(len=max_string) :: fieldname
534 
535  ! allocate bump arrays
536  allocate(fld_ingrid_mga(self%bump%geom%nmga,self%nlev))
537  allocate(fld_outgrid(self%nout_local,self%nlev))
538 
539  if (.not.self%bump%geom%same_grid) then
540  allocate(fld_ingrid_c0a(self%bump%geom%nc0a, self%nlev))
541  endif
542 
543  do ifield = 1, fields_outgrid%size()
544 
545  field_outgrid = fields_outgrid%field(ifield)
546  fieldname = field_outgrid%name()
547 
548  !--------------------------------------------
549  ! allocate output field if necessary
550 
551  if (.not. fields_ingrid%has_field(fieldname)) then
552  field_ingrid = self%in_funcspace%create_field(name=fieldname, &
553  kind=atlas_real(kind_real),levels=self%nlev)
554  else
555  field_ingrid = fields_ingrid%field(name=fieldname)
556  endif
557 
558  !--------------------------------------------
559  ! compute interpolation adjoint
560 
561  ! atlas field to bump fld
562  call field_to_array(self%bump%mpl, field_outgrid, fld_outgrid)
563 
564  if (self%bump%geom%same_grid) then
565  ! Apply observation operator adjoint
566  call self%apply_interp_ad(fld_outgrid,fld_ingrid_mga)
567  else
568  ! Apply observation operator adjoint
569  call self%apply_interp_ad(fld_outgrid,fld_ingrid_c0a)
570 
571  ! Subset Sc0 to model grid
572  call self%bump%geom%copy_c0a_to_mga(self%bump%mpl,fld_ingrid_c0a,fld_ingrid_mga)
573  end if
574 
575  ! bump fld to atlas field
576  call field_from_array(self%bump%mpl, fld_ingrid_mga, field_ingrid)
577 
578  !--------------------------------------------
579 
580  ! add field to result
581  if (.not. fields_ingrid%has_field(fieldname)) then
582  call fields_ingrid%add(field_ingrid)
583  endif
584 
585  ! release pointers
586  call field_ingrid%final()
587  call field_outgrid%final()
588 
589  enddo
590 
591  ! clean up
592  if (allocated(fld_ingrid_mga)) deallocate(fld_ingrid_mga)
593  if (allocated(fld_ingrid_c0a)) deallocate(fld_ingrid_c0a)
594  if (allocated(fld_outgrid)) deallocate(fld_outgrid)
595 
596 end subroutine bint_apply_ad
597 
598 !----------------------------------------------------------------------
599 !> Subroutine: apply_interp_ad
600 !! Purpose: low-level routine to apply the adjoint of the interpolation operator
601 !! to a single field on a single level
602 !!
603 !! \param[in] mpl bump MPI data
604 !! \param[in] geom bump geometry data
605 !! \param[in] fld_outgrid field on output grid
606 !! \param[out] fld_ingrid field on input grid
607 !!
608 
609 subroutine apply_interp_ad(self,fld_outgrid,fld_ingrid)
610  class(bump_interpolator), intent(inout) :: self
611  real(kind_real),intent(in) :: fld_outgrid(self%nout_local,self%nlev)
612  real(kind_real),intent(out) :: fld_ingrid(self%bump%geom%nc0a,self%nlev)
613 
614  integer :: ilev
615  real(kind_real) :: fld_ingrid_ext(self%nc0b,self%nlev)
616 
617  if (self%nout_local > 0) then
618  ! Horizontal interpolation
619  !$omp parallel do schedule(static) private(ilev)
620  do ilev = 1, self%nlev
621  call self%h%apply_ad(self%bump%mpl,fld_outgrid(:,ilev),fld_ingrid_ext(:,ilev))
622  end do
623  !$omp end parallel do
624  else
625  ! No observation on this task
626  fld_ingrid_ext = 0.0
627  end if
628 
629  ! Halo reduction
630  call self%com%red(self%bump%mpl,self%nlev,fld_ingrid_ext,fld_ingrid)
631 
632 end subroutine apply_interp_ad
633 
634 !----------------------------------------------------------------------
635 !> Release memory (partial) by deallocating output grid
636 !----------------------------------------------------------------------
637 subroutine bint_deallocate_outgrid(self)
638  class(bump_interpolator), intent(inout) :: self
639 
640  ! Release memory
641  call self%outgeom%dealloc()
642 
643 end subroutine bint_deallocate_outgrid
644 
645 ! ------------------------------------------------------------------------------
646 !> Release all memory
647 ! ------------------------------------------------------------------------------
648 
649 subroutine bint_delete(self)
650  class(bump_interpolator), intent(inout) :: self
651 
652  call self%bump%dealloc()
653 
654  call self%deallocate_outgrid()
655  call self%h%dealloc()
656  call self%com%dealloc()
657 
658 end subroutine bint_delete
659 
660 !----------------------------------------------------------------------
661 ! Subroutine: dummy
662 !> Dummy finalization
663 !----------------------------------------------------------------------
664 subroutine dummy(bump)
665 
666 implicit none
667 
668 ! Passed variables
669 type(bump_interpolator),intent(inout) :: bump !< BUMP
670 
671 end subroutine dummy
672 
673 end module bump_interpolation_mod
tools_func::sphere_dist
subroutine, public sphere_dist(lon_i, lat_i, lon_f, lat_f, dist)
Compute the great-circle distance between two points.
Definition: tools_func.F90:126
bump_interpolation_mod::bint_delete
subroutine bint_delete(self)
Release all memory.
Definition: bump_interpolation_mod.F90:650
type_com
Communications derived type.
Definition: type_com.F90:8
type_rng::rng_type
Definition: type_rng.F90:22
bump_interpolation_mod
Bump Interpolation module.
Definition: bump_interpolation_mod.F90:12
bump_interpolation_mod::apply_interp_ad
subroutine apply_interp_ad(self, fld_outgrid, fld_ingrid)
Subroutine: apply_interp_ad Purpose: low-level routine to apply the adjoint of the interpolation oper...
Definition: bump_interpolation_mod.F90:610
tools_const
Define usual constants and missing values.
Definition: tools_const.F90:8
bump_interpolation_mod::bint_deallocate_outgrid
subroutine bint_deallocate_outgrid(self)
Release memory (partial) by deallocating output grid.
Definition: bump_interpolation_mod.F90:638
type_bump
BUMP derived type.
Definition: type_bump.F90:8
bump_interpolation_mod::bint_driver
subroutine bint_driver(self, mpl, rng, nam, geom)
Initialize bump to perform interpolation.
Definition: bump_interpolation_mod.F90:268
type_mesh
Mesh derived type.
Definition: type_mesh.F90:8
tools_func
Usual functions.
Definition: tools_func.F90:8
type_nam
Namelist derived type.
Definition: type_nam.F90:8
bump_interpolation_mod::bint_apply_ad
subroutine bint_apply_ad(self, fields_outgrid, fields_ingrid)
Apply interpolator operator adjoint.
Definition: bump_interpolation_mod.F90:524
tools_qsort
Qsort routines.
Definition: tools_qsort.F90:12
type_geom::geom_type
Definition: type_geom.F90:31
type_nam::nam_type
Definition: type_nam.F90:29
type_rng
Random numbers generator derived type.
Definition: type_rng.F90:8
bump_interpolation_mod::bump_interpolator
Definition: bump_interpolation_mod.F90:49
tools_repro
Reproducibility functions.
Definition: tools_repro.F90:8
type_linop::linop_type
Definition: type_linop.F90:39
tools_atlas::field_from_array
Definition: tools_atlas.F90:24
tools_atlas
Random numbers generator derived type.
Definition: tools_atlas.F90:9
bump_interpolation_mod::apply_interp
subroutine apply_interp(self, infield, outfield)
Subroutine: apply_interp Purpose: low-level routine to apply the interpolation to a single field on a...
Definition: bump_interpolation_mod.F90:481
type_fieldset
Random numbers generator derived type.
Definition: type_fieldset.F90:9
type_mesh::mesh_type
Definition: type_mesh.F90:24
type_bump::bump_type
Definition: type_bump.F90:38
bump_interpolation_mod::bint_init
subroutine bint_init(self, config, comm, in_funcspace, out_funcspace, masks)
Linked list implementation.
Definition: bump_interpolation_mod.F90:134
tools_repro::rth
real(kind_real), parameter, public rth
Definition: tools_repro.F90:17
tools_qsort::qsort
Definition: tools_qsort.F90:19
type_fieldset::fieldset_type
Definition: type_fieldset.F90:18
tools_atlas::field_to_array
Definition: tools_atlas.F90:19
tools_kinds
Kinds definition.
Definition: tools_kinds.F90:8
bump_interpolation_mod::dummy
subroutine dummy(bump)
Dummy finalization.
Definition: bump_interpolation_mod.F90:665
tools_const::rad2deg
real(kind_real), parameter, public rad2deg
Definition: tools_const.F90:16
tools_const::deg2rad
real(kind_real), parameter, public deg2rad
Definition: tools_const.F90:15
type_geom
Geometry derived type.
Definition: type_geom.F90:8
type_linop
Linear operator derived type.
Definition: type_linop.F90:8
type_mpl
MPI parameters derived type.
Definition: type_mpl.F90:8
bump_interpolation_mod::bint_apply
subroutine bint_apply(self, infields, outfields)
Apply interpolation.
Definition: bump_interpolation_mod.F90:400
type_com::com_type
Definition: type_com.F90:21
bump_interpolation_mod::bump_interpolator_registry
type(registry_t), public bump_interpolator_registry
Registry for bump_interpolator objects.
Definition: bump_interpolation_mod.F90:99
tools_func::lonlatmod
subroutine, public lonlatmod(lon, lat)
Set latitude between -pi/2 and pi/2 and longitude between -pi and pi.
Definition: tools_func.F90:66
tools_const::pi
real(kind_real), parameter, public pi
Definition: tools_const.F90:14
bump_interpolation_mod::max_string
integer, parameter max_string
Definition: bump_interpolation_mod.F90:43
type_mpl::mpl_type
Definition: type_mpl.F90:24
tools_kinds::kind_real
integer, parameter, public kind_real
Definition: tools_kinds.F90:18