MODFLOW 6  version 6.7.0.dev1
USGS Modular Hydrologic Model
MeshNCModel.f90
Go to the documentation of this file.
1 !> @brief This module contains the MeshModelModule
2 !!
3 !! This module defines a base class for UGRID based
4 !! (mesh) model netcdf exports. It is dependent on
5 !! external netcdf libraries.
6 !<
8 
9  use kindmodule, only: dp, i4b, lgp
18  use netcdfcommonmodule, only: nf_verify
19  use netcdf
20 
21  implicit none
22  private
24  public :: mesh2dmodeltype
25  public :: ncvar_chunk
26  public :: ncvar_deflate
27  public :: ncvar_gridmap
28  public :: ncvar_mf6attr
29  public :: export_varname
30 
31  !> @brief type for storing model export dimension ids
32  !<
34  integer(I4B) :: nmesh_node !< number of nodes in mesh
35  integer(I4B) :: nmesh_face !< number of faces in mesh
36  integer(I4B) :: max_nmesh_face_nodes !< max number of nodes in a single face
37  integer(I4B) :: nlay !< number of layers
38  integer(I4B) :: time !< number of steps
39  contains
40  end type meshncdimidtype
41 
42  !> @brief type for storing model export variable ids
43  !<
45  integer(I4B) :: mesh !< mesh container variable
46  integer(I4B) :: mesh_node_x !< mesh nodes x array
47  integer(I4B) :: mesh_node_y !< mesh nodes y array
48  integer(I4B) :: mesh_face_x !< mesh faces x location array
49  integer(I4B) :: mesh_face_y !< mesh faces y location array
50  integer(I4B) :: mesh_face_xbnds !< mesh faces 2D x bounds array
51  integer(I4B) :: mesh_face_ybnds !< mesh faces 2D y bounds array
52  integer(I4B) :: mesh_face_nodes !< mesh faces 2D nodes array
53  integer(I4B) :: time !< time coordinate variable
54  integer(I4B), dimension(:), allocatable :: dependent !< layered dependent variables array
55  contains
56  end type meshncvaridtype
57 
58  !> @brief base ugrid netcdf export type
59  !<
60  type, abstract, extends(ncbasemodelexporttype) :: meshmodeltype
61  type(meshncdimidtype) :: dim_ids !< dimension ids
62  type(meshncvaridtype) :: var_ids !< variable ids
63  integer(I4B) :: nlay !< number of layers
64  integer(I4B), pointer :: chunk_face !< chunking parameter for face dimension
65  contains
66  procedure :: mesh_init
67  procedure :: mesh_destroy
68  procedure :: add_global_att
69  procedure(nc_array_export_if), deferred :: export_input_array
70  procedure :: export_input_arrays
71  procedure :: add_pkg_data
72  procedure :: define_dependent
73  procedure :: define_gridmap
74  end type meshmodeltype
75 
76  !> @brief abstract interfaces for derived ugrid netcd export types
77  !<
78  abstract interface
79  subroutine nc_array_export_if(this, pkgtype, pkgname, mempath, idt)
81  class(meshmodeltype), intent(inout) :: this
82  character(len=*), intent(in) :: pkgtype
83  character(len=*), intent(in) :: pkgname
84  character(len=*), intent(in) :: mempath
85  type(inputparamdefinitiontype), pointer, intent(in) :: idt
86  end subroutine
87  end interface
88 
89  type, abstract, extends(meshmodeltype) :: mesh2dmodeltype
90  contains
91  procedure :: create_mesh
92  end type mesh2dmodeltype
93 
94 contains
95 
96  !> @brief initialize
97  !<
98  subroutine mesh_init(this, modelname, modeltype, modelfname, nc_fname, &
99  disenum, nctype, lenuni, iout)
101  class(meshmodeltype), intent(inout) :: this
102  character(len=*), intent(in) :: modelname
103  character(len=*), intent(in) :: modeltype
104  character(len=*), intent(in) :: modelfname
105  character(len=*), intent(in) :: nc_fname
106  integer(I4B), intent(in) :: disenum
107  integer(I4B), intent(in) :: nctype
108  integer(I4B), intent(in) :: lenuni
109  integer(I4B), intent(in) :: iout
110  logical(LGP) :: found
111 
112  ! initialize base class
113  call this%NCModelExportType%init(modelname, modeltype, modelfname, nc_fname, &
114  disenum, nctype, iout)
115 
116  ! allocate and initialize
117  allocate (this%chunk_face)
118  this%chunk_face = -1
119 
120  ! update values from input context
121  if (this%ncf_mempath /= '') then
122  call mem_set_value(this%chunk_face, 'CHUNK_FACE', this%ncf_mempath, found)
123  end if
124 
125  if (this%chunk_time > 0 .and. this%chunk_face > 0) then
126  this%chunking_active = .true.
127  else if (this%chunk_time > 0 .or. this%chunk_face > 0) then
128  this%chunk_face = -1
129  this%chunk_time = -1
130  write (warnmsg, '(a)') 'Ignoring user provided NetCDF chunking parameter. &
131  &Define chunk_time and chunk_face input parameters to see an effect in &
132  &file "'//trim(nc_fname)//'".'
133  call store_warning(warnmsg)
134  end if
135 
136  if (lenuni == 1) then
137  this%lenunits = 'ft'
138  else
139  this%lenunits = 'm'
140  end if
141 
142  ! create the netcdf file
143  call nf_verify(nf90_create(this%nc_fname, &
144  ior(nf90_clobber, nf90_netcdf4), this%ncid), &
145  this%nc_fname)
146  end subroutine mesh_init
147 
148  !> @brief initialize
149  !<
150  subroutine mesh_destroy(this)
152  class(meshmodeltype), intent(inout) :: this
153  call nf_verify(nf90_close(this%ncid), this%nc_fname)
154  deallocate (this%chunk_face)
155  nullify (this%chunk_face)
156  end subroutine mesh_destroy
157 
158  !> @brief create file (group) attributes
159  !<
160  subroutine add_global_att(this)
161  class(meshmodeltype), intent(inout) :: this
162  ! file scoped title
163  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'title', &
164  this%annotation%title), this%nc_fname)
165  ! source (MODFLOW 6)
166  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'source', &
167  this%annotation%source), this%nc_fname)
168  ! export type (MODFLOW 6)
169  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'modflow_grid', &
170  this%annotation%grid), this%nc_fname)
171  ! MODFLOW 6 model type
172  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'modflow_model', &
173  this%annotation%model), this%nc_fname)
174  ! generation datetime
175  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'history', &
176  this%annotation%history), this%nc_fname)
177  ! supported conventions
178  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'Conventions', &
179  this%annotation%conventions), &
180  this%nc_fname)
181  end subroutine add_global_att
182 
183  !> @brief write package gridded input data
184  !<
185  subroutine export_input_arrays(this, pkgtype, pkgname, mempath, param_dfns)
186  use memorymanagermodule, only: get_isize
187  class(meshmodeltype), intent(inout) :: this
188  character(len=*), intent(in) :: pkgtype
189  character(len=*), intent(in) :: pkgname
190  character(len=*), intent(in) :: mempath
191  type(inputparamdefinitiontype), dimension(:), pointer, &
192  intent(in) :: param_dfns
193  type(inputparamdefinitiontype), pointer :: idt
194  integer(I4B) :: iparam, isize
195  ! export griddata block parameters
196  do iparam = 1, size(param_dfns)
197  ! assign param definition pointer
198  idt => param_dfns(iparam)
199  ! for now only griddata is exported
200  if (idt%blockname == 'GRIDDATA') then
201  ! veriy variable is allocated
202  call get_isize(idt%mf6varname, mempath, isize)
203  if (isize > 0) then
204  call this%export_input_array(pkgtype, pkgname, mempath, idt)
205  end if
206  end if
207  end do
208  end subroutine export_input_arrays
209 
210  !> @brief determine packages to write gridded input
211  !<
212  subroutine add_pkg_data(this)
219  class(meshmodeltype), intent(inout) :: this
220  character(LENCOMPONENTNAME) :: ptype, pname, pkgtype
221  character(len=LENMEMPATH) :: input_mempath
222  type(characterstringtype), dimension(:), contiguous, &
223  pointer :: pkgtypes => null()
224  type(characterstringtype), dimension(:), contiguous, &
225  pointer :: pkgnames => null()
226  type(characterstringtype), dimension(:), contiguous, &
227  pointer :: mempaths => null()
228  type(inputparamdefinitiontype), dimension(:), pointer :: param_dfns
229  character(len=LENMEMPATH) :: mempath
230  integer(I4B) :: n
231  integer(I4B), pointer :: export_arrays
232  logical(LGP) :: found
233 
234  input_mempath = create_mem_path(component=this%modelname, context=idm_context)
235 
236  ! set pointers to model path package info
237  call mem_setptr(pkgtypes, 'PKGTYPES', input_mempath)
238  call mem_setptr(pkgnames, 'PKGNAMES', input_mempath)
239  call mem_setptr(mempaths, 'MEMPATHS', input_mempath)
240 
241  allocate (export_arrays)
242 
243  do n = 1, size(mempaths)
244  ! initialize export_arrays
245  export_arrays = 0
246 
247  ! set package attributes
248  mempath = mempaths(n)
249  pname = pkgnames(n)
250  ptype = pkgtypes(n)
251 
252  ! export input arrays
253  if (mempath /= '') then
254  ! update export
255  call mem_set_value(export_arrays, 'EXPORT_NC', mempath, found)
256  if (export_arrays > 0) then
257  pkgtype = idm_subcomponent_type(this%modeltype, ptype)
258  param_dfns => param_definitions(this%modeltype, pkgtype)
259  call this%export_input_arrays(ptype, pname, mempath, param_dfns)
260  end if
261  end if
262  end do
263 
264  ! cleanup
265  deallocate (export_arrays)
266  end subroutine add_pkg_data
267 
268  !> @brief create the model layer dependent variables
269  !<
270  subroutine define_dependent(this)
271  class(meshmodeltype), intent(inout) :: this
272  character(len=LINELENGTH) :: varname, longname
273  integer(I4B) :: k
274 
275  ! create a dependent variable for each layer
276  do k = 1, this%nlay
277  ! initialize names
278  varname = ''
279  longname = ''
280 
281  ! set layer variable and longnames
282  write (varname, '(a,i0)') trim(this%xname)//'_l', k
283  write (longname, '(a,i0,a)') trim(this%annotation%longname)// &
284  ' (layer ', k, ')'
285 
286  ! create the netcdf dependent layer variable
287  call nf_verify(nf90_def_var(this%ncid, varname, nf90_double, &
288  (/this%dim_ids%nmesh_face, &
289  this%dim_ids%time/), &
290  this%var_ids%dependent(k)), &
291  this%nc_fname)
292 
293  ! apply chunking parameters
294  if (this%chunking_active) then
295  call nf_verify(nf90_def_var_chunking(this%ncid, &
296  this%var_ids%dependent(k), &
297  nf90_chunked, &
298  (/this%chunk_face, &
299  this%chunk_time/)), &
300  this%nc_fname)
301  end if
302 
303  ! deflate and shuffle
304  call ncvar_deflate(this%ncid, this%var_ids%dependent(k), this%deflate, &
305  this%shuffle, this%nc_fname)
306 
307  ! assign variable attributes
308  call nf_verify(nf90_put_att(this%ncid, this%var_ids%dependent(k), &
309  'units', this%lenunits), this%nc_fname)
310  call nf_verify(nf90_put_att(this%ncid, this%var_ids%dependent(k), &
311  'standard_name', this%annotation%stdname), &
312  this%nc_fname)
313  call nf_verify(nf90_put_att(this%ncid, this%var_ids%dependent(k), &
314  'long_name', longname), this%nc_fname)
315  call nf_verify(nf90_put_att(this%ncid, this%var_ids%dependent(k), &
316  '_FillValue', (/dhnoflo/)), &
317  this%nc_fname)
318  call nf_verify(nf90_put_att(this%ncid, this%var_ids%dependent(k), &
319  'mesh', this%mesh_name), this%nc_fname)
320  call nf_verify(nf90_put_att(this%ncid, this%var_ids%dependent(k), &
321  'location', 'face'), this%nc_fname)
322 
323  ! add grid mapping
324  call ncvar_gridmap(this%ncid, this%var_ids%dependent(k), &
325  this%gridmap_name, this%nc_fname)
326  end do
327  end subroutine define_dependent
328 
329  !> @brief create the file grid mapping container variable
330  !<
331  subroutine define_gridmap(this)
332  class(meshmodeltype), intent(inout) :: this
333  integer(I4B) :: var_id
334 
335  ! was projection info provided
336  if (this%wkt /= '') then
337  ! create projection variable
338  call nf_verify(nf90_redef(this%ncid), this%nc_fname)
339  call nf_verify(nf90_def_var(this%ncid, this%gridmap_name, nf90_int, &
340  var_id), this%nc_fname)
341  ! cf-conventions prefers 'crs_wkt'
342  !call nf_verify(nf90_put_att(this%ncid, var_id, 'crs_wkt', this%wkt), &
343  ! this%nc_fname)
344  ! QGIS recognizes 'wkt'
345  call nf_verify(nf90_put_att(this%ncid, var_id, 'wkt', this%wkt), &
346  this%nc_fname)
347  call nf_verify(nf90_enddef(this%ncid), this%nc_fname)
348  call nf_verify(nf90_put_var(this%ncid, var_id, 1), &
349  this%nc_fname)
350  end if
351  end subroutine define_gridmap
352 
353  !> @brief create the file mesh container variable
354  !<
355  subroutine create_mesh(this)
356  class(mesh2dmodeltype), intent(inout) :: this
357 
358  ! create mesh container variable
359  call nf_verify(nf90_def_var(this%ncid, this%mesh_name, nf90_int, &
360  this%var_ids%mesh), this%nc_fname)
361 
362  ! assign container variable attributes
363  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, 'cf_role', &
364  'mesh_topology'), this%nc_fname)
365  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, 'long_name', &
366  '2D mesh topology'), this%nc_fname)
367  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, &
368  'topology_dimension', 2), this%nc_fname)
369  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, 'face_dimension', &
370  'nmesh_face'), this%nc_fname)
371  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, &
372  'node_coordinates', 'mesh_node_x mesh_node_y'), &
373  this%nc_fname)
374  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, &
375  'face_coordinates', 'mesh_face_x mesh_face_y'), &
376  this%nc_fname)
377  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, &
378  'face_node_connectivity', 'mesh_face_nodes'), &
379  this%nc_fname)
380 
381  ! create mesh x node (mesh vertex) variable
382  call nf_verify(nf90_def_var(this%ncid, 'mesh_node_x', nf90_double, &
383  (/this%dim_ids%nmesh_node/), &
384  this%var_ids%mesh_node_x), this%nc_fname)
385 
386  ! assign mesh x node variable attributes
387  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_x, &
388  'units', this%lenunits), this%nc_fname)
389  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_x, &
390  'standard_name', 'projection_x_coordinate'), &
391  this%nc_fname)
392  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_x, &
393  'long_name', 'Easting'), this%nc_fname)
394 
395  if (this%wkt /= '') then
396  ! associate with projection
397  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_x, &
398  'grid_mapping', this%gridmap_name), &
399  this%nc_fname)
400  end if
401 
402  ! create mesh y node (mesh vertex) variable
403  call nf_verify(nf90_def_var(this%ncid, 'mesh_node_y', nf90_double, &
404  (/this%dim_ids%nmesh_node/), &
405  this%var_ids%mesh_node_y), this%nc_fname)
406 
407  ! assign mesh y variable attributes
408  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_y, &
409  'units', this%lenunits), this%nc_fname)
410  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_y, &
411  'standard_name', 'projection_y_coordinate'), &
412  this%nc_fname)
413  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_y, &
414  'long_name', 'Northing'), this%nc_fname)
415 
416  if (this%wkt /= '') then
417  ! associate with projection
418  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_y, &
419  'grid_mapping', this%gridmap_name), &
420  this%nc_fname)
421  end if
422 
423  ! create mesh x face (cell vertex) variable
424  call nf_verify(nf90_def_var(this%ncid, 'mesh_face_x', nf90_double, &
425  (/this%dim_ids%nmesh_face/), &
426  this%var_ids%mesh_face_x), this%nc_fname)
427 
428  ! assign mesh x face variable attributes
429  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_x, &
430  'units', this%lenunits), this%nc_fname)
431  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_x, &
432  'standard_name', 'projection_x_coordinate'), &
433  this%nc_fname)
434  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_x, &
435  'long_name', 'Easting'), this%nc_fname)
436  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_x, 'bounds', &
437  'mesh_face_xbnds'), this%nc_fname)
438  if (this%wkt /= '') then
439  ! associate with projection
440  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_x, &
441  'grid_mapping', this%gridmap_name), &
442  this%nc_fname)
443  end if
444 
445  ! create mesh x cell bounds variable
446  call nf_verify(nf90_def_var(this%ncid, 'mesh_face_xbnds', nf90_double, &
447  (/this%dim_ids%max_nmesh_face_nodes, &
448  this%dim_ids%nmesh_face/), &
449  this%var_ids%mesh_face_xbnds), &
450  this%nc_fname)
451 
452  ! create mesh y face (cell vertex) variable
453  call nf_verify(nf90_def_var(this%ncid, 'mesh_face_y', nf90_double, &
454  (/this%dim_ids%nmesh_face/), &
455  this%var_ids%mesh_face_y), this%nc_fname)
456 
457  ! assign mesh y face variable attributes
458  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_y, &
459  'units', this%lenunits), this%nc_fname)
460  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_y, &
461  'standard_name', 'projection_y_coordinate'), &
462  this%nc_fname)
463  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_y, &
464  'long_name', 'Northing'), this%nc_fname)
465  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_y, 'bounds', &
466  'mesh_face_ybnds'), this%nc_fname)
467 
468  if (this%wkt /= '') then
469  ! associate with projection
470  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_y, &
471  'grid_mapping', this%gridmap_name), &
472  this%nc_fname)
473  end if
474 
475  ! create mesh y cell bounds variable
476  call nf_verify(nf90_def_var(this%ncid, 'mesh_face_ybnds', nf90_double, &
477  (/this%dim_ids%max_nmesh_face_nodes, &
478  this%dim_ids%nmesh_face/), &
479  this%var_ids%mesh_face_ybnds), &
480  this%nc_fname)
481 
482  ! create mesh face nodes variable
483  call nf_verify(nf90_def_var(this%ncid, 'mesh_face_nodes', nf90_int, &
484  (/this%dim_ids%max_nmesh_face_nodes, &
485  this%dim_ids%nmesh_face/), &
486  this%var_ids%mesh_face_nodes), &
487  this%nc_fname)
488 
489  ! assign variable attributes
490  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_nodes, &
491  'cf_role', 'face_node_connectivity'), &
492  this%nc_fname)
493  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_nodes, &
494  'long_name', &
495  'Vertices bounding cell (counterclockwise)'), &
496  this%nc_fname)
497  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_nodes, &
498  '_FillValue', (/nf90_fill_int/)), &
499  this%nc_fname)
500  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_nodes, &
501  'start_index', 1), this%nc_fname)
502  end subroutine create_mesh
503 
504  !> @brief define variable chunking
505  !<
506  subroutine ncvar_chunk(ncid, varid, chunk_face, nc_fname)
507  integer(I4B), intent(in) :: ncid
508  integer(I4B), intent(in) :: varid
509  integer(I4B), intent(in) :: chunk_face
510  character(len=*), intent(in) :: nc_fname
511  if (chunk_face > 0) then
512  call nf_verify(nf90_def_var_chunking(ncid, varid, nf90_chunked, &
513  (/chunk_face/)), nc_fname)
514  end if
515  end subroutine ncvar_chunk
516 
517  !> @brief define variable compression
518  !<
519  subroutine ncvar_deflate(ncid, varid, deflate, shuffle, nc_fname)
520  integer(I4B), intent(in) :: ncid
521  integer(I4B), intent(in) :: varid
522  integer(I4B), intent(in) :: deflate
523  integer(I4B), intent(in) :: shuffle
524  character(len=*), intent(in) :: nc_fname
525  if (deflate >= 0) then
526  call nf_verify(nf90_def_var_deflate(ncid, varid, shuffle=shuffle, &
527  deflate=1, deflate_level=deflate), &
528  nc_fname)
529  end if
530  end subroutine ncvar_deflate
531 
532  !> @brief put variable gridmap attributes
533  !<
534  subroutine ncvar_gridmap(ncid, varid, gridmap_name, nc_fname)
535  integer(I4B), intent(in) :: ncid
536  integer(I4B), intent(in) :: varid
537  character(len=*), intent(in) :: gridmap_name
538  character(len=*), intent(in) :: nc_fname
539  if (gridmap_name /= '') then
540  call nf_verify(nf90_put_att(ncid, varid, 'coordinates', &
541  'mesh_face_x mesh_face_y'), nc_fname)
542  call nf_verify(nf90_put_att(ncid, varid, 'grid_mapping', &
543  gridmap_name), nc_fname)
544  end if
545  end subroutine ncvar_gridmap
546 
547  !> @brief put variable internal attributes
548  !<
549  subroutine ncvar_mf6attr(ncid, varid, layer, iper, iaux, nc_tag, nc_fname)
550  integer(I4B), intent(in) :: ncid
551  integer(I4B), intent(in) :: varid
552  integer(I4B), intent(in) :: layer
553  integer(I4B), intent(in) :: iper
554  integer(I4B), intent(in) :: iaux
555  character(len=*), intent(in) :: nc_tag
556  character(len=*), intent(in) :: nc_fname
557  if (nc_tag /= '') then
558  call nf_verify(nf90_put_att(ncid, varid, 'modflow6_input', &
559  nc_tag), nc_fname)
560  if (layer > 0) then
561  call nf_verify(nf90_put_att(ncid, varid, 'modflow6_layer', &
562  layer), nc_fname)
563  end if
564  if (iper > 0) then
565  call nf_verify(nf90_put_att(ncid, varid, 'modflow6_iper', &
566  iper), nc_fname)
567  end if
568  if (iaux > 0) then
569  call nf_verify(nf90_put_att(ncid, varid, 'modflow6_iaux', &
570  iaux), nc_fname)
571  end if
572  end if
573  end subroutine ncvar_mf6attr
574 
575  !> @brief build netcdf variable name
576  !<
577  function export_varname(varname, layer, iper, iaux) result(vname)
578  use inputoutputmodule, only: lowcase
579  character(len=*), intent(in) :: varname
580  integer(I4B), optional, intent(in) :: layer
581  integer(I4B), optional, intent(in) :: iper
582  integer(I4B), optional, intent(in) :: iaux
583  character(len=LINELENGTH) :: vname
584  vname = ''
585  if (varname /= '') then
586  vname = varname
587  call lowcase(vname)
588  if (present(layer)) then
589  if (layer > 0) then
590  write (vname, '(a,i0)') trim(vname)//'_l', layer
591  end if
592  end if
593  if (present(iper)) then
594  if (iper > 0) then
595  write (vname, '(a,i0)') trim(vname)//'_p', iper
596  end if
597  end if
598  if (present(iaux)) then
599  if (iaux > 0) then
600  write (vname, '(a,i0)') trim(vname)//'a', iaux
601  end if
602  end if
603  end if
604  end function export_varname
605 
606 end module meshmodelmodule
abstract interfaces for derived ugrid netcd export types
Definition: MeshNCModel.f90:79
This module contains simulation constants.
Definition: Constants.f90:9
integer(i4b), parameter linelength
maximum length of a standard line
Definition: Constants.f90:45
integer(i4b), parameter lencomponentname
maximum length of a component name
Definition: Constants.f90:18
real(dp), parameter dnodata
real no data constant
Definition: Constants.f90:95
real(dp), parameter dhnoflo
real no flow constant
Definition: Constants.f90:93
integer(i4b), parameter lenmempath
maximum length of the memory path
Definition: Constants.f90:27
type(inputparamdefinitiontype) function, dimension(:), pointer, public param_definitions(component, subcomponent)
logical function, public idm_multi_package(component, subcomponent)
This module contains the InputDefinitionModule.
subroutine, public lowcase(word)
Convert to lower case.
This module defines variable data types.
Definition: kind.f90:8
character(len=lenmempath) function create_mem_path(component, subcomponent, context)
returns the path to the memory object
subroutine, public get_isize(name, mem_path, isize)
@ brief Get the number of elements for this variable
This module contains the MeshModelModule.
Definition: MeshNCModel.f90:7
subroutine define_gridmap(this)
create the file grid mapping container variable
character(len=linelength) function, public export_varname(varname, layer, iper, iaux)
build netcdf variable name
subroutine, public ncvar_gridmap(ncid, varid, gridmap_name, nc_fname)
put variable gridmap attributes
subroutine, public ncvar_chunk(ncid, varid, chunk_face, nc_fname)
define variable chunking
subroutine mesh_init(this, modelname, modeltype, modelfname, nc_fname, disenum, nctype, lenuni, iout)
initialize
subroutine mesh_destroy(this)
initialize
subroutine add_global_att(this)
create file (group) attributes
subroutine add_pkg_data(this)
determine packages to write gridded input
subroutine export_input_arrays(this, pkgtype, pkgname, mempath, param_dfns)
write package gridded input data
subroutine define_dependent(this)
create the model layer dependent variables
subroutine, public ncvar_deflate(ncid, varid, deflate, shuffle, nc_fname)
define variable compression
subroutine, public ncvar_mf6attr(ncid, varid, layer, iper, iaux, nc_tag, nc_fname)
put variable internal attributes
subroutine create_mesh(this)
create the file mesh container variable
This module contains the NCModelExportModule.
Definition: NCModel.f90:8
This module contains the NetCDFCommonModule.
Definition: NetCDFCommon.f90:6
subroutine, public nf_verify(res, nc_fname)
error check a netcdf-fortran interface call
This module contains simulation methods.
Definition: Sim.f90:10
subroutine, public store_warning(msg, substring)
Store warning message.
Definition: Sim.f90:236
subroutine, public store_error(msg, terminate)
Store an error message.
Definition: Sim.f90:92
subroutine, public store_error_filename(filename, terminate)
Store the erroring file name.
Definition: Sim.f90:203
This module contains simulation variables.
Definition: SimVariables.f90:9
character(len=maxcharlen) errmsg
error message string
character(len=linelength) idm_context
character(len=maxcharlen) warnmsg
warning message string
This module contains the SourceCommonModule.
Definition: SourceCommon.f90:7
character(len=lencomponentname) function, public idm_subcomponent_type(component, subcomponent)
component from package or model type
This class is used to store a single deferred-length character string. It was designed to work in an ...
Definition: CharString.f90:23
base ugrid netcdf export type
Definition: MeshNCModel.f90:60
type for storing model export dimension ids
Definition: MeshNCModel.f90:33
type for storing model export variable ids
Definition: MeshNCModel.f90:44
abstract type for model netcdf export type
Definition: NCModel.f90:102