OOPS
qg_gom_mod.F90
Go to the documentation of this file.
1 ! (C) Copyright 2009-2016 ECMWF.
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 ! In applying this licence, ECMWF does not waive the privileges and immunities
6 ! granted to it by virtue of its status as an intergovernmental organisation nor
7 ! does it submit to any jurisdiction.
8 
9 module qg_gom_mod
10 
11 use atlas_module, only: atlas_field
12 use fckit_configuration_module, only: fckit_configuration
13 use fckit_log_module, only: fckit_log
14 use iso_c_binding
15 use kinds
16 use netcdf
19 use qg_geom_mod
20 use qg_locs_mod
22 use qg_tools_mod
23 use random_mod
24 
25 implicit none
26 private
27 public :: qg_gom
28 public :: qg_gom_registry
32 ! ------------------------------------------------------------------------------
33 type :: qg_gom
34  integer :: nobs !< Number of observations
35  integer :: used !< Index of used observation
36  integer :: ix !< Streamfunction index
37  integer :: iq !< Potential vorticity index
38  integer :: iu !< Zonal wind index
39  integer :: iv !< Meridian wind index
40  integer :: nv !< Number of variables
41  real(kind_real), allocatable :: values(:,:) !< Observations values
42  logical :: lalloc !< Allocation flag
43 end type qg_gom
44 
45 #define LISTED_TYPE qg_gom
46 
47 !> Linked list interface - defines registry_t type
48 #include "oops/util/linkedList_i.f"
49 
50 !> Global registry
51 type(registry_t) :: qg_gom_registry
52 ! ------------------------------------------------------------------------------
53 contains
54 ! ------------------------------------------------------------------------------
55 ! Public
56 ! ------------------------------------------------------------------------------
57 !> Linked list implementation
58 #include "oops/util/linkedList_c.f"
59 ! ------------------------------------------------------------------------------
60 !> Setup GOM
61 subroutine qg_gom_setup(self,nobs,vars)
62 
63 implicit none
64 
65 ! Passed variables
66 type(qg_gom),intent(inout) :: self !< GOM
67 integer, intent(in) :: nobs !< Number of observations
68 type(oops_variables),intent(in) :: vars !< Variables
69 
70 ! Local variables
71 integer :: ivar
72 
73 ! Set attributes
74 self%nobs = nobs
75 self%used = 0
76 self%nv = 0
77 self%ix = 0; self%iq = 0; self%iu = 0; self%iv = 0
78 do ivar = 1, vars%nvars()
79  if (vars%variable(ivar) == 'x') then
80  self%nv = self%nv+1
81  self%ix = self%nv
82  endif
83  if (vars%variable(ivar) == 'q') then
84  self%nv = self%nv+1
85  self%iq = self%nv
86  endif
87  if (vars%variable(ivar) == 'u') then
88  self%nv = self%nv+1
89  self%iu = self%nv
90  endif
91  if (vars%variable(ivar) == 'v') then
92  self%nv = self%nv+1
93  self%iv = self%nv
94  endif
95 enddo
96 
97 ! Allocation
98 allocate(self%values(self%nv,self%nobs))
99 self%lalloc = .true.
100 
101 end subroutine qg_gom_setup
102 ! ------------------------------------------------------------------------------
103 !> Create GOM
104 subroutine qg_gom_create(self)
105 
106 implicit none
107 
108 ! Passed variables
109 type(qg_gom),intent(inout) :: self !< GOM
110 
111 ! Set allocation flag
112 self%lalloc = .false.
113 
114 end subroutine qg_gom_create
115 ! ------------------------------------------------------------------------------
116 !> Delete GOM
117 subroutine qg_gom_delete(self)
118 
119 implicit none
120 
121 ! Passed variables
122 type(qg_gom),intent(inout) :: self !< GOM
123 
124 ! Release memory
125 if (self%lalloc) then
126  deallocate(self%values)
127 endif
128 
129 end subroutine qg_gom_delete
130 ! ------------------------------------------------------------------------------
131 !> Copy GOM
132 subroutine qg_gom_copy(self,other)
133 
134 implicit none
135 
136 ! Passed variables
137 type(qg_gom),intent(inout) :: self !< GOM
138 type(qg_gom),intent(in) :: other !< Other GOM
139 
140 ! Copy attribues
141 self%nobs = other%nobs
142 self%ix = other%ix
143 self%iq = other%iq
144 self%iu = other%iu
145 self%iv = other%iv
146 self%nv = other%nv
147 self%used = other%used
148 ! Allocation
149 if (.not.self%lalloc) then
150  allocate(self%values(self%nv,self%nobs))
151  self%lalloc = .true.
152 endif
153 
154 ! Copy
155 self%values = other%values
156 
157 end subroutine qg_gom_copy
158 ! ------------------------------------------------------------------------------
159 !> Set GOM to zero
160 subroutine qg_gom_zero(self)
161 
162 implicit none
163 
164 ! Passed variables
165 type(qg_gom),intent(inout) :: self !< GOM
166 
167 ! Set to zero
168 self%values = 0.0
169 
170 end subroutine qg_gom_zero
171 ! ------------------------------------------------------------------------------
172 !> Get GOM absolute value
173 subroutine qg_gom_abs(self)
174 
175 implicit none
176 
177 ! Passed variables
178 type(qg_gom),intent(inout) :: self !< GOM
179 
180 ! Get absolute value
181 self%values = abs(self%values)
182 
183 end subroutine qg_gom_abs
184 ! ------------------------------------------------------------------------------
185 !> Generate random GOM values
186 subroutine qg_gom_random(self)
187 
188 implicit none
189 
190 ! Passed variables
191 type(qg_gom),intent(inout) :: self !< GOM
192 
193 ! Generate random GOM values
194 call normal_distribution(self%values,0.0_kind_real,1.0_kind_real)
195 
196 end subroutine qg_gom_random
197 ! ------------------------------------------------------------------------------
198 !> Multiply GOM with a scalar
199 subroutine qg_gom_mult(self,zz)
200 
201 implicit none
202 
203 ! Passed variables
204 type(qg_gom),intent(inout) :: self !< GOM
205 real(kind_real),intent(in) :: zz !< Multiplier
206 
207 ! Local variables
208 integer :: jo,jv
209 
210 ! Multiply GOM with a scalar
211 do jo=1,self%nobs
212  do jv=1,self%nv
213  self%values(jv,jo) = zz*self%values(jv,jo)
214  enddo
215 enddo
216 
217 end subroutine qg_gom_mult
218 ! ------------------------------------------------------------------------------
219 !> Add GOM
220 subroutine qg_gom_add(self,other)
221 
222 implicit none
223 
224 ! Passed variables
225 type(qg_gom),intent(inout) :: self !< GOM
226 type(qg_gom),intent(in) :: other !< Other GOM
227 
228 ! Local variables
229 integer :: jo,jv
230 
231 ! Add GOM
232 do jo=1,self%nobs
233  do jv=1,self%nv
234  self%values(jv,jo) = self%values(jv,jo)+other%values(jv,jo)
235  enddo
236 enddo
237 
238 end subroutine qg_gom_add
239 ! ------------------------------------------------------------------------------
240 !> Subtract GOM
241 subroutine qg_gom_diff(self,other)
242 
243 implicit none
244 
245 ! Passed variables
246 type(qg_gom),intent(inout) :: self !< GOM
247 type(qg_gom),intent(in) :: other !< Other GOM
248 
249 ! Local variables
250 integer :: jo,jv
251 
252 ! Subtract GOM
253 do jo=1,self%nobs
254  do jv=1,self%nv
255  self%values(jv,jo) = self%values(jv,jo)-other%values(jv,jo)
256  enddo
257 enddo
258 
259 end subroutine qg_gom_diff
260 ! ------------------------------------------------------------------------------
261 !> Schur product for GOM
262 subroutine qg_gom_schurmult(self,other)
263 
264 implicit none
265 
266 ! Passed variables
267 type(qg_gom),intent(inout) :: self !< GOM
268 type(qg_gom),intent(in) :: other !< Other GOM
269 
270 ! Local variables
271 integer :: jo,jv
272 
273 ! Add GOM
274 do jo=1,self%nobs
275  do jv=1,self%nv
276  self%values(jv,jo) = self%values(jv,jo)*other%values(jv,jo)
277  enddo
278 enddo
279 
280 end subroutine qg_gom_schurmult
281 ! ------------------------------------------------------------------------------
282 !> Schur division for GOM
283 subroutine qg_gom_divide(self,other)
284 
285 implicit none
286 
287 ! Passed variables
288 type(qg_gom),intent(inout) :: self !< GOM
289 type(qg_gom),intent(in) :: other !< Other GOM
290 
291 ! Local variables
292 real(kind_real) :: tol
293 integer :: jloc,jvar
294 
295 ! Set tolerance
296 tol = epsilon(tol)
297 
298 ! Conditional division
299 do jvar=1,self%nv
300  do jloc=1,self%nobs
301  if (abs(other%values(jvar,jloc))>tol) then
302  self%values(jvar,jloc) = self%values(jvar,jloc)/other%values(jvar,jloc)
303  else
304  self%values(jvar,jloc) = 0.0
305  endif
306  enddo
307 enddo
308 
309 end subroutine qg_gom_divide
310 ! ------------------------------------------------------------------------------
311 !> Compute GOM RMS
312 subroutine qg_gom_rms(self,rms)
313 
314 implicit none
315 
316 ! Passed variables
317 type(qg_gom),intent(inout) :: self !< GOM
318 real(kind_real),intent(inout) :: rms !< RMS
319 
320 ! Local variables
321 integer :: jo,jv
322 
323 ! Initialization
324 rms = 0.0
325 
326 ! Loop over values
327 do jo=1,self%nobs
328  do jv=1,self%nv
329  rms = rms+self%values(jv,jo)**2
330  enddo
331 enddo
332 
333 ! Normalize and take square-root
334 rms = sqrt(rms/(self%nobs*self%nv))
335 
336 end subroutine qg_gom_rms
337 ! ------------------------------------------------------------------------------
338 !> GOM dot product
339 subroutine qg_gom_dotprod(gom1,gom2,prod)
340 
341 implicit none
342 
343 ! Passed variables
344 type(qg_gom),intent(in) :: gom1 !< GOM 1
345 type(qg_gom),intent(in) :: gom2 !< GOM 2
346 real(kind_real),intent(inout) :: prod !< Dot product
347 
348 ! Local variables
349 integer :: jo,jv
350 
351 ! Check
352 if ((gom1%nv/=gom2%nv).or.(gom1%nobs/=gom2%nobs)) call abor1_ftn('qg_gom_dotprod: inconsistent GOM sizes')
353 
354 ! Initialization
355 prod = 0.0
356 
357 ! Loop over values
358 do jo=1,gom1%nobs
359  do jv=1,gom1%nv
360  prod = prod+gom1%values(jv,jo)*gom2%values(jv,jo)
361  enddo
362 enddo
363 
364 end subroutine qg_gom_dotprod
365 ! ------------------------------------------------------------------------------
366 !> Compute GOM stats
367 subroutine qg_gom_stats(self,kobs,scaling,pmin,pmax,prms)
368 
369 implicit none
370 
371 ! Passed variables
372 type(qg_gom),intent(inout) :: self !< GOM
373 integer,intent(inout) :: kobs !< Number of observations
374 real(kind_real),intent(inout) :: scaling !< Scaling value
375 real(kind_real),intent(inout) :: pmin !< Minimum value
376 real(kind_real),intent(inout) :: pmax !< Maximum value
377 real(kind_real),intent(inout) :: prms !< RMS
378 
379 ! Local variables
380 real(kind_real) :: expo,scalinginv,svalues(self%nv,self%nobs)
381 character(len=1024) :: msg
382 
383 ! Scaling
384 if (maxval(abs(self%values))>0.0) then
385  expo = aint(log(maxval(abs(self%values)))/log(10.0_kind_real))-1
386  scaling = 10.0**expo
387 else
388  scaling = 1.0
389 endif
390 scalinginv = 1.0/scaling
391 
392 ! Scale values
393 svalues = self%values*scalinginv
394 
395 ! Compute GOM stats
396 kobs = self%nobs
397 pmin = minval(svalues)
398 pmax = maxval(svalues)
399 prms = sqrt(sum(svalues**2)/real(self%nobs*self%nv,kind_real))
400 
401 end subroutine qg_gom_stats
402 ! ------------------------------------------------------------------------------
403 !> Find and locate GOM max. value
404 subroutine qg_gom_maxloc(self,mxval,iloc,ivar)
405 
406 implicit none
407 
408 ! Passed variables
409 type(qg_gom),intent(inout) :: self !< GOM
410 real(kind_real),intent(inout) :: mxval !< Maximum value
411 integer,intent(inout) :: iloc !< Location of maximum value
412 integer,intent(inout) :: ivar !< Variable with maximum value
413 
414 ! Local variables
415 integer :: mxloc(2)
416 
417 ! Find GOM max. value
418 mxval = maxval(self%values)
419 mxloc = maxloc(self%values)
420 
421 ! Locate GOM max. value
422 ivar = mxloc(1)
423 iloc = mxloc(2)
424 
425 end subroutine qg_gom_maxloc
426 ! ------------------------------------------------------------------------------
427 !> Read GOM from file
428 subroutine qg_gom_read_file(self,f_conf)
429 
430 implicit none
431 
432 ! Passed variables
433 type(qg_gom),intent(inout) :: self !< GOM
434 type(fckit_configuration),intent(in) :: f_conf !< FCKIT configuration
435 
436 ! Local variables
437 integer :: ncid,nobs_id,nv_id,values_id
438 character(len=1024) :: filename
439 character(len=:),allocatable :: str
440 
441 ! Check allocation
442 if (self%lalloc) call abor1_ftn('qg_gom_read_file: gom alredy allocated')
443 
444 ! Get filename
445 call f_conf%get_or_die("filename",str)
446 filename = str
447 call fckit_log%info('qg_gom_read_file: reading '//trim(filename))
448 
449 ! Open NetCDF file
450 call ncerr(nf90_open(trim(filename)//'.nc',nf90_nowrite,ncid))
451 
452 ! Get dimensions ids
453 call ncerr(nf90_inq_dimid(ncid,'nobs',nobs_id))
454 call ncerr(nf90_inq_dimid(ncid,'nv',nv_id))
455 
456 ! Get dimensions
457 call ncerr(nf90_inquire_dimension(ncid,nobs_id,len=self%nobs))
458 call ncerr(nf90_inquire_dimension(ncid,nv_id,len=self%nv))
459 
460 ! Get attributes
461 call ncerr(nf90_get_att(ncid,nf90_global,'used',self%used))
462 call ncerr(nf90_get_att(ncid,nf90_global,'ix',self%ix))
463 call ncerr(nf90_get_att(ncid,nf90_global,'iq',self%iq))
464 call ncerr(nf90_get_att(ncid,nf90_global,'iu',self%iu))
465 call ncerr(nf90_get_att(ncid,nf90_global,'iv',self%iv))
466 call ncerr(nf90_get_att(ncid,nf90_global,'nv',self%nv))
467 
468 ! Allocation
469 allocate(self%values(self%nv,self%nobs))
470 self%lalloc = .true.
471 
472 ! Initialization
473 call qg_gom_zero(self)
474 
475 ! Get variables ids
476 call ncerr(nf90_inq_varid(ncid,'values',values_id))
477 
478 ! Get variables
479 call ncerr(nf90_get_var(ncid,values_id,self%values))
480 
481 ! Close NetCDF file
482 call ncerr(nf90_close(ncid))
483 
484 end subroutine qg_gom_read_file
485 ! ------------------------------------------------------------------------------
486 !> Write GOM to file
487 subroutine qg_gom_write_file(self,f_conf)
488 
489 implicit none
490 
491 ! Passed variables
492 type(qg_gom),intent(inout) :: self !< GOM
493 type(fckit_configuration),intent(in) :: f_conf !< FCKIT configuration
494 
495 ! Local variables
496 integer :: ncid,nobs_id,nv_id,values_id
497 character(len=1024) :: filename
498 character(len=:),allocatable :: str
499 
500 ! Check allocation
501 if (.not.self%lalloc) call abor1_ftn('qg_gom_write_file: gom not allocated')
502 
503 ! Set filename
504 call f_conf%get_or_die("filename",str)
505 filename = str
506 call fckit_log%info('qg_gom_write_file: writing '//trim(filename))
507 
508 ! Create NetCDF file
509 call ncerr(nf90_create(trim(filename)//'.nc',or(nf90_clobber,nf90_64bit_offset),ncid))
510 
511 ! Define dimensions
512 call ncerr(nf90_def_dim(ncid,'nobs',self%nobs,nobs_id))
513 call ncerr(nf90_def_dim(ncid,'nv',self%nv,nv_id))
514 
515 ! Define attributes
516 call ncerr(nf90_put_att(ncid,nf90_global,'used',self%used))
517 call ncerr(nf90_put_att(ncid,nf90_global,'ix',self%ix))
518 call ncerr(nf90_put_att(ncid,nf90_global,'iq',self%iq))
519 call ncerr(nf90_put_att(ncid,nf90_global,'iu',self%iu))
520 call ncerr(nf90_put_att(ncid,nf90_global,'iv',self%iv))
521 call ncerr(nf90_put_att(ncid,nf90_global,'nv',self%nv))
522 
523 ! Define variables
524 call ncerr(nf90_def_var(ncid,'values',nf90_double,(/nv_id,nobs_id/),values_id))
525 
526 ! End definitions
527 call ncerr(nf90_enddef(ncid))
528 
529 ! Put variables
530 call ncerr(nf90_put_var(ncid,values_id,self%values))
531 
532 ! Close NetCDF file
533 call ncerr(nf90_close(ncid))
534 
535 end subroutine qg_gom_write_file
536 ! ------------------------------------------------------------------------------
537 !> GOM analytic initialization
538 subroutine qg_gom_analytic_init(self,locs,f_conf)
539 
540 implicit none
541 
542 ! Passed variables
543 type(qg_gom),intent(inout) :: self !< GOM
544 type(qg_locs),intent(inout) :: locs !< Locations
545 type(fckit_configuration),intent(in) :: f_conf !< FCKIT configuration
546 
547 ! Local variables
548 integer :: iloc
549 real(kind_real) :: x,y
550 character(len=30) :: ic
551 character(len=:),allocatable :: str
552 real(kind_real), pointer :: lonlat(:,:), z(:)
553 type(atlas_field) :: lonlat_field, z_field
554 
555 ! get locations
556 lonlat_field = locs%lonlat()
557 call lonlat_field%data(lonlat)
558 
559 z_field = locs%altitude()
560 call z_field%data(z)
561 
562 ! Check allocation
563 if (.not. self%lalloc) call abor1_ftn('qg_gom_analytic init: gom not allocated')
564 
565 ! Get analytic configuration
566 call f_conf%get_or_die("analytic_init",str)
567 ic = str
568 call fckit_log%info('qg_gom_analytic_init: ic = '//trim(ic))
569 do iloc=1,locs%nlocs()
570  select case (trim(ic))
571  case ('baroclinic-instability')
572  ! Go to cartesian coordinates
573  call lonlat_to_xy(lonlat(1,iloc),lonlat(2,iloc),x,y)
574 
575  ! Compute values for baroclinic instability
576  if (self%ix>0) call baroclinic_instability(x,y,z(iloc),'x',self%values(self%ix,iloc))
577  if (self%iq>0) call baroclinic_instability(x,y,z(iloc),'q',self%values(self%iq,iloc))
578  if (self%iu>0) call baroclinic_instability(x,y,z(iloc),'u',self%values(self%iu,iloc))
579  if (self%iv>0) call baroclinic_instability(x,y,z(iloc),'v',self%values(self%iv,iloc))
580  case ('large-vortices')
581  ! Go to cartesian coordinates
582  call lonlat_to_xy(lonlat(1,iloc),lonlat(2,iloc),x,y)
583 
584  ! Compute values for large vortices
585  if (self%ix>0) call large_vortices(x,y,z(iloc),'x',self%values(self%ix,iloc))
586  if (self%iq>0) call large_vortices(x,y,z(iloc),'q',self%values(self%iq,iloc))
587  if (self%iu>0) call large_vortices(x,y,z(iloc),'u',self%values(self%iu,iloc))
588  if (self%iv>0) call large_vortices(x,y,z(iloc),'v',self%values(self%iv,iloc))
589  case default
590  call abor1_ftn('qg_gom_analytic_init: unknown initialization')
591  endselect
592 enddo
593 
594 call lonlat_field%final()
595 call z_field%final()
596 
597 end subroutine qg_gom_analytic_init
598 ! ------------------------------------------------------------------------------
599 end module qg_gom_mod
qg_geom_mod
Definition: qg_geom_mod.F90:9
qg_projection_mod
Definition: qg_projection_mod.F90:9
qg_tools_mod::large_vortices
subroutine, public large_vortices(x, y, z, var, res)
Generate values for large vortices.
Definition: qg_tools_mod.F90:148
qg_gom_mod::qg_gom_zero
subroutine, public qg_gom_zero(self)
Set GOM to zero.
Definition: qg_gom_mod.F90:161
qg_gom_mod::qg_gom_copy
subroutine, public qg_gom_copy(self, other)
Copy GOM.
Definition: qg_gom_mod.F90:133
qg_gom_mod::qg_gom_setup
subroutine, public qg_gom_setup(self, nobs, vars)
Linked list implementation.
Definition: qg_gom_mod.F90:62
oops_variables_mod
Fortran interface to Variables.
Definition: variables_mod.F90:9
qg_locs_mod
Definition: qg_locs_mod.F90:9
qg_projection_mod::lonlat_to_xy
subroutine, public lonlat_to_xy(lon, lat, x, y)
Convert lon/lat to x/y.
Definition: qg_projection_mod.F90:57
qg_gom_mod::qg_gom_dotprod
subroutine, public qg_gom_dotprod(gom1, gom2, prod)
GOM dot product.
Definition: qg_gom_mod.F90:340
qg_gom_mod::qg_gom_stats
subroutine, public qg_gom_stats(self, kobs, scaling, pmin, pmax, prms)
Compute GOM stats.
Definition: qg_gom_mod.F90:368
qg_gom_mod::qg_gom_registry
type(registry_t), public qg_gom_registry
Linked list interface - defines registry_t type.
Definition: qg_gom_mod.F90:51
qg_locs_mod::qg_locs
Definition: qg_locs_mod.F90:22
qg_gom_mod::qg_gom
Definition: qg_gom_mod.F90:33
qg_gom_mod::qg_gom_delete
subroutine, public qg_gom_delete(self)
Delete GOM.
Definition: qg_gom_mod.F90:118
qg_gom_mod::qg_gom_create
subroutine, public qg_gom_create(self)
Create GOM.
Definition: qg_gom_mod.F90:105
qg_constants_mod
Definition: qg_constants_mod.F90:9
qg_gom_mod::qg_gom_diff
subroutine, public qg_gom_diff(self, other)
Subtract GOM.
Definition: qg_gom_mod.F90:242
qg_gom_mod::qg_gom_read_file
subroutine, public qg_gom_read_file(self, f_conf)
Read GOM from file.
Definition: qg_gom_mod.F90:429
qg_gom_mod::qg_gom_abs
subroutine, public qg_gom_abs(self)
Get GOM absolute value.
Definition: qg_gom_mod.F90:174
qg_gom_mod::qg_gom_rms
subroutine, public qg_gom_rms(self, rms)
Compute GOM RMS.
Definition: qg_gom_mod.F90:313
qg_tools_mod::baroclinic_instability
subroutine, public baroclinic_instability(x, y, z, var, res)
Generate values for baroclinic instability.
Definition: qg_tools_mod.F90:112
qg_gom_mod::qg_gom_maxloc
subroutine, public qg_gom_maxloc(self, mxval, iloc, ivar)
Find and locate GOM max. value.
Definition: qg_gom_mod.F90:405
qg_gom_mod::qg_gom_write_file
subroutine, public qg_gom_write_file(self, f_conf)
Write GOM to file.
Definition: qg_gom_mod.F90:488
oops_variables_mod::oops_variables
Definition: variables_mod.F90:16
qg_gom_mod::qg_gom_random
subroutine, public qg_gom_random(self)
Generate random GOM values.
Definition: qg_gom_mod.F90:187
qg_tools_mod
Definition: qg_tools_mod.F90:9
qg_gom_mod
Definition: qg_gom_mod.F90:9
qg_gom_mod::qg_gom_add
subroutine, public qg_gom_add(self, other)
Add GOM.
Definition: qg_gom_mod.F90:221
qg_gom_mod::qg_gom_schurmult
subroutine, public qg_gom_schurmult(self, other)
Schur product for GOM.
Definition: qg_gom_mod.F90:263
qg_tools_mod::ncerr
subroutine, public ncerr(info)
Check NetCDF status.
Definition: qg_tools_mod.F90:99
qg_gom_mod::qg_gom_analytic_init
subroutine, public qg_gom_analytic_init(self, locs, f_conf)
GOM analytic initialization.
Definition: qg_gom_mod.F90:539
qg_gom_mod::qg_gom_mult
subroutine, public qg_gom_mult(self, zz)
Multiply GOM with a scalar.
Definition: qg_gom_mod.F90:200
qg_gom_mod::qg_gom_divide
subroutine, public qg_gom_divide(self, other)
Schur division for GOM.
Definition: qg_gom_mod.F90:284