MODFLOW 6  version 6.6.0.dev0
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
12  use simvariablesmodule, only: errmsg
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, disenum, &
99  nctype, 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  integer(I4B), intent(in) :: disenum
106  integer(I4B), intent(in) :: nctype
107  integer(I4B), intent(in) :: iout
108  logical(LGP) :: found
109  !
110  ! -- initialize base class
111  call this%NCModelExportType%init(modelname, modeltype, modelfname, disenum, &
112  nctype, iout)
113  !
114  ! -- allocate and initialize
115  allocate (this%chunk_face)
116  this%chunk_face = -1
117  !
118  ! -- update values from input context
119  if (this%ncf_mempath /= '') then
120  call mem_set_value(this%chunk_face, 'CHUNK_FACE', this%ncf_mempath, found)
121  end if
122  !
123  if (this%chunk_time > 0 .and. this%chunk_face > 0) then
124  this%chunking_active = .true.
125  end if
126  !
127  ! -- create the netcdf file
128  call nf_verify(nf90_create(this%nc_fname, &
129  iand(nf90_clobber, nf90_netcdf4), this%ncid), &
130  this%nc_fname)
131  end subroutine mesh_init
132 
133  !> @brief initialize
134  !<
135  subroutine mesh_destroy(this)
137  class(meshmodeltype), intent(inout) :: this
138  !
139  call nf_verify(nf90_close(this%ncid), this%nc_fname)
140  !
141  deallocate (this%chunk_face)
142  nullify (this%chunk_face)
143  end subroutine mesh_destroy
144 
145  !> @brief create file (group) attributes
146  !<
147  subroutine add_global_att(this)
148  class(meshmodeltype), intent(inout) :: this
149  ! -- file scoped title
150  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'title', &
151  this%annotation%title), this%nc_fname)
152  ! -- source (MODFLOW 6)
153  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'source', &
154  this%annotation%source), this%nc_fname)
155  ! -- export type (MODFLOW 6)
156  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'modflow6_grid', &
157  this%annotation%grid), this%nc_fname)
158  ! -- MODFLOW 6 model type
159  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'modflow6_model', &
160  this%annotation%model), this%nc_fname)
161  ! -- generation datetime
162  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'history', &
163  this%annotation%history), this%nc_fname)
164  ! -- supported conventions
165  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'Conventions', &
166  this%annotation%conventions), &
167  this%nc_fname)
168  end subroutine add_global_att
169 
170  !> @brief write package gridded input data
171  !<
172  subroutine export_input_arrays(this, pkgtype, pkgname, mempath, param_dfns)
173  use memorymanagermodule, only: get_isize
174  class(meshmodeltype), intent(inout) :: this
175  character(len=*), intent(in) :: pkgtype
176  character(len=*), intent(in) :: pkgname
177  character(len=*), intent(in) :: mempath
178  type(inputparamdefinitiontype), dimension(:), pointer, &
179  intent(in) :: param_dfns
180  type(inputparamdefinitiontype), pointer :: idt
181  integer(I4B) :: iparam, isize
182  !
183  ! -- export griddata block parameters
184  do iparam = 1, size(param_dfns)
185  ! -- assign param definition pointer
186  idt => param_dfns(iparam)
187  ! -- for now
188  if (idt%blockname == 'GRIDDATA') then
189  ! -- veriy variable is allocated
190  call get_isize(idt%mf6varname, mempath, isize)
191  !
192  if (isize > 0) then
193  call this%export_input_array(pkgtype, pkgname, mempath, idt)
194  end if
195  end if
196  end do
197  end subroutine export_input_arrays
198 
199  !> @brief determine packages to write gridded input
200  !<
201  subroutine add_pkg_data(this)
208  class(meshmodeltype), intent(inout) :: this
209  character(LENCOMPONENTNAME) :: ptype, pname, pkgtype
210  character(len=LENMEMPATH) :: input_mempath
211  type(characterstringtype), dimension(:), contiguous, &
212  pointer :: pkgtypes => null()
213  type(characterstringtype), dimension(:), contiguous, &
214  pointer :: pkgnames => null()
215  type(characterstringtype), dimension(:), contiguous, &
216  pointer :: mempaths => null()
217  type(inputparamdefinitiontype), dimension(:), pointer :: param_dfns
218  character(len=LENMEMPATH) :: mempath
219  integer(I4B) :: n
220  integer(I4B), pointer :: export_arrays
221  logical(LGP) :: found
222  !
223  input_mempath = create_mem_path(component=this%modelname, context=idm_context)
224  !
225  ! -- set pointers to model path package info
226  call mem_setptr(pkgtypes, 'PKGTYPES', input_mempath)
227  call mem_setptr(pkgnames, 'PKGNAMES', input_mempath)
228  call mem_setptr(mempaths, 'MEMPATHS', input_mempath)
229  !
230  allocate (export_arrays)
231  !
232  do n = 1, size(mempaths)
233  !
234  ! -- initialize export_arrays
235  export_arrays = 0
236  !
237  ! -- set package attributes
238  mempath = mempaths(n)
239  pname = pkgnames(n)
240  ptype = pkgtypes(n)
241  !
242  ! -- export input arrays
243  if (mempath /= '') then
244  ! -- update export
245  call mem_set_value(export_arrays, 'EXPORT_NC', mempath, found)
246  !
247  if (export_arrays > 0) then
248  pkgtype = idm_subcomponent_type(this%modeltype, ptype)
249  param_dfns => param_definitions(this%modeltype, pkgtype)
250  call this%export_input_arrays(ptype, pname, mempath, param_dfns)
251  end if
252  end if
253  end do
254  !
255  ! -- cleanup
256  deallocate (export_arrays)
257  end subroutine add_pkg_data
258 
259  !> @brief create the model layer dependent variables
260  !<
261  subroutine define_dependent(this)
262  class(meshmodeltype), intent(inout) :: this
263  character(len=LINELENGTH) :: varname, longname
264  integer(I4B) :: k
265  !
266  ! -- create a dependent variable for each layer
267  do k = 1, this%nlay
268  !
269  ! -- initialize names
270  varname = ''
271  longname = ''
272  !
273  ! -- set layer variable and longnames
274  write (varname, '(a,i0)') trim(this%xname)//'_l', k
275  write (longname, '(a,i0,a)') trim(this%annotation%longname)// &
276  ' (layer ', k, ')'
277  !
278  ! -- create the netcdf dependent layer variable
279  call nf_verify(nf90_def_var(this%ncid, varname, nf90_double, &
280  (/this%dim_ids%nmesh_face, &
281  this%dim_ids%time/), &
282  this%var_ids%dependent(k)), &
283  this%nc_fname)
284  !
285  ! -- apply chunking parameters
286  if (this%chunking_active) then
287  call nf_verify(nf90_def_var_chunking(this%ncid, &
288  this%var_ids%dependent(k), &
289  nf90_chunked, &
290  (/this%chunk_face, &
291  this%chunk_time/)), &
292  this%nc_fname)
293  end if
294  ! -- deflate and shuffle
295  call ncvar_deflate(this%ncid, this%var_ids%dependent(k), this%deflate, &
296  this%shuffle, this%nc_fname)
297  !
298  ! -- assign variable attributes
299  call nf_verify(nf90_put_att(this%ncid, this%var_ids%dependent(k), &
300  'units', 'm'), this%nc_fname)
301  call nf_verify(nf90_put_att(this%ncid, this%var_ids%dependent(k), &
302  'standard_name', this%annotation%stdname), &
303  this%nc_fname)
304  call nf_verify(nf90_put_att(this%ncid, this%var_ids%dependent(k), &
305  'long_name', longname), this%nc_fname)
306  call nf_verify(nf90_put_att(this%ncid, this%var_ids%dependent(k), &
307  '_FillValue', (/dhnoflo/)), &
308  this%nc_fname)
309  call nf_verify(nf90_put_att(this%ncid, this%var_ids%dependent(k), &
310  'mesh', this%mesh_name), this%nc_fname)
311  call nf_verify(nf90_put_att(this%ncid, this%var_ids%dependent(k), &
312  'location', 'face'), this%nc_fname)
313  !
314  ! -- add grid mapping
315  call ncvar_gridmap(this%ncid, this%var_ids%dependent(k), &
316  this%gridmap_name, this%nc_fname)
317  end do
318  end subroutine define_dependent
319 
320  !> @brief create the file grid mapping container variable
321  !<
322  subroutine define_gridmap(this)
323  class(meshmodeltype), intent(inout) :: this
324  integer(I4B) :: var_id
325  !
326  ! -- was projection info provided
327  if (this%ogc_wkt /= '') then
328  !
329  ! -- create projection variable
330  call nf_verify(nf90_redef(this%ncid), this%nc_fname)
331  call nf_verify(nf90_def_var(this%ncid, this%gridmap_name, nf90_int, &
332  var_id), this%nc_fname)
333  ! -- cf-conventions prefers 'crs_wkt'
334  !call nf_verify(nf90_put_att(this%ncid, var_id, 'crs_wkt', this%ogc_wkt), &
335  ! this%nc_fname)
336  ! -- QGIS recognizes 'wkt'
337  call nf_verify(nf90_put_att(this%ncid, var_id, 'wkt', this%ogc_wkt), &
338  this%nc_fname)
339  call nf_verify(nf90_enddef(this%ncid), this%nc_fname)
340  call nf_verify(nf90_put_var(this%ncid, var_id, 1), &
341  this%nc_fname)
342  end if
343  end subroutine define_gridmap
344 
345  !> @brief create the file mesh container variable
346  !<
347  subroutine create_mesh(this)
348  class(mesh2dmodeltype), intent(inout) :: this
349  !
350  ! -- create mesh container variable
351  call nf_verify(nf90_def_var(this%ncid, this%mesh_name, nf90_int, &
352  this%var_ids%mesh), this%nc_fname)
353  !
354  ! -- assign container variable attributes
355  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, 'cf_role', &
356  'mesh_topology'), this%nc_fname)
357  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, 'long_name', &
358  '2D mesh topology'), this%nc_fname)
359  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, &
360  'topology_dimension', 2), this%nc_fname)
361  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, 'face_dimension', &
362  'nmesh_face'), this%nc_fname)
363  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, &
364  'node_coordinates', 'mesh_node_x mesh_node_y'), &
365  this%nc_fname)
366  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, &
367  'face_coordinates', 'mesh_face_x mesh_face_y'), &
368  this%nc_fname)
369  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, &
370  'face_node_connectivity', 'mesh_face_nodes'), &
371  this%nc_fname)
372 
373  ! -- create mesh x node (mesh vertex) variable
374  call nf_verify(nf90_def_var(this%ncid, 'mesh_node_x', nf90_double, &
375  (/this%dim_ids%nmesh_node/), &
376  this%var_ids%mesh_node_x), this%nc_fname)
377  !
378  ! -- assign mesh x node variable attributes
379  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_x, &
380  'units', 'm'), this%nc_fname)
381  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_x, &
382  'standard_name', 'projection_x_coordinate'), &
383  this%nc_fname)
384  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_x, &
385  'long_name', 'Easting'), this%nc_fname)
386  !
387  if (this%ogc_wkt /= '') then
388  ! -- associate with projection
389  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_x, &
390  'grid_mapping', this%gridmap_name), &
391  this%nc_fname)
392  end if
393 
394  ! -- create mesh y node (mesh vertex) variable
395  call nf_verify(nf90_def_var(this%ncid, 'mesh_node_y', nf90_double, &
396  (/this%dim_ids%nmesh_node/), &
397  this%var_ids%mesh_node_y), this%nc_fname)
398  !
399  ! -- assign mesh y variable attributes
400  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_y, &
401  'units', 'm'), this%nc_fname)
402  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_y, &
403  'standard_name', 'projection_y_coordinate'), &
404  this%nc_fname)
405  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_y, &
406  'long_name', 'Northing'), this%nc_fname)
407  !
408  if (this%ogc_wkt /= '') then
409  ! -- associate with projection
410  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_y, &
411  'grid_mapping', this%gridmap_name), &
412  this%nc_fname)
413  end if
414 
415  ! -- create mesh x face (cell vertex) variable
416  call nf_verify(nf90_def_var(this%ncid, 'mesh_face_x', nf90_double, &
417  (/this%dim_ids%nmesh_face/), &
418  this%var_ids%mesh_face_x), this%nc_fname)
419  !
420  ! -- assign mesh x face variable attributes
421  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_x, &
422  'units', 'm'), this%nc_fname)
423  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_x, &
424  'standard_name', 'projection_x_coordinate'), &
425  this%nc_fname)
426  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_x, &
427  'long_name', 'Easting'), this%nc_fname)
428  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_x, 'bounds', &
429  'mesh_face_xbnds'), this%nc_fname)
430  if (this%ogc_wkt /= '') then
431  ! -- associate with projection
432  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_x, &
433  'grid_mapping', this%gridmap_name), &
434  this%nc_fname)
435  end if
436 
437  ! -- create mesh x cell bounds variable
438  call nf_verify(nf90_def_var(this%ncid, 'mesh_face_xbnds', nf90_double, &
439  (/this%dim_ids%max_nmesh_face_nodes, &
440  this%dim_ids%nmesh_face/), &
441  this%var_ids%mesh_face_xbnds), &
442  this%nc_fname)
443  !
444  ! -- create mesh y face (cell vertex) variable
445  call nf_verify(nf90_def_var(this%ncid, 'mesh_face_y', nf90_double, &
446  (/this%dim_ids%nmesh_face/), &
447  this%var_ids%mesh_face_y), this%nc_fname)
448  !
449  ! -- assign mesh y face variable attributes
450  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_y, &
451  'units', 'm'), this%nc_fname)
452  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_y, &
453  'standard_name', 'projection_y_coordinate'), &
454  this%nc_fname)
455  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_y, &
456  'long_name', 'Northing'), this%nc_fname)
457  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_y, 'bounds', &
458  'mesh_face_ybnds'), this%nc_fname)
459  !
460  if (this%ogc_wkt /= '') then
461  ! -- associate with projection
462  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_y, &
463  'grid_mapping', this%gridmap_name), &
464  this%nc_fname)
465  end if
466 
467  ! -- create mesh y cell bounds variable
468  call nf_verify(nf90_def_var(this%ncid, 'mesh_face_ybnds', nf90_double, &
469  (/this%dim_ids%max_nmesh_face_nodes, &
470  this%dim_ids%nmesh_face/), &
471  this%var_ids%mesh_face_ybnds), &
472  this%nc_fname)
473  !
474  ! -- create mesh face nodes variable
475  call nf_verify(nf90_def_var(this%ncid, 'mesh_face_nodes', nf90_int, &
476  (/this%dim_ids%max_nmesh_face_nodes, &
477  this%dim_ids%nmesh_face/), &
478  this%var_ids%mesh_face_nodes), &
479  this%nc_fname)
480  !
481  ! -- assign variable attributes
482  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_nodes, &
483  'cf_role', 'face_node_connectivity'), &
484  this%nc_fname)
485  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_nodes, &
486  'long_name', &
487  'Vertices bounding cell (counterclockwise)'), &
488  this%nc_fname)
489  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_nodes, &
490  '_FillValue', (/nf90_fill_int/)), &
491  this%nc_fname)
492  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_nodes, &
493  'start_index', 1), this%nc_fname)
494  end subroutine create_mesh
495 
496  !> @brief define variable chunking
497  !<
498  subroutine ncvar_chunk(ncid, varid, chunk_face, nc_fname)
499  integer(I4B), intent(in) :: ncid
500  integer(I4B), intent(in) :: varid
501  integer(I4B), intent(in) :: chunk_face
502  character(len=*), intent(in) :: nc_fname
503  !
504  ! -- apply chunking parameters
505  if (chunk_face > 0) then
506  call nf_verify(nf90_def_var_chunking(ncid, varid, nf90_chunked, &
507  (/chunk_face/)), nc_fname)
508  end if
509  end subroutine ncvar_chunk
510 
511  !> @brief define variable compression
512  !<
513  subroutine ncvar_deflate(ncid, varid, deflate, shuffle, nc_fname)
514  integer(I4B), intent(in) :: ncid
515  integer(I4B), intent(in) :: varid
516  integer(I4B), intent(in) :: deflate
517  integer(I4B), intent(in) :: shuffle
518  character(len=*), intent(in) :: nc_fname
519  ! -- deflate and shuffle
520  if (deflate >= 0) then
521  call nf_verify(nf90_def_var_deflate(ncid, varid, shuffle=shuffle, &
522  deflate=1, deflate_level=deflate), &
523  nc_fname)
524  end if
525  end subroutine ncvar_deflate
526 
527  !> @brief put variable gridmap attributes
528  !<
529  subroutine ncvar_gridmap(ncid, varid, gridmap_name, nc_fname)
530  integer(I4B), intent(in) :: ncid
531  integer(I4B), intent(in) :: varid
532  character(len=*), intent(in) :: gridmap_name
533  character(len=*), intent(in) :: nc_fname
534  !
535  if (gridmap_name /= '') then
536  call nf_verify(nf90_put_att(ncid, varid, 'coordinates', &
537  'mesh_face_x mesh_face_y'), nc_fname)
538  call nf_verify(nf90_put_att(ncid, varid, 'grid_mapping', &
539  gridmap_name), nc_fname)
540  end if
541  end subroutine ncvar_gridmap
542 
543  !> @brief put variable internal attributes
544  !<
545  subroutine ncvar_mf6attr(ncid, varid, layer, iper, iaux, nc_tag, nc_fname)
546  integer(I4B), intent(in) :: ncid
547  integer(I4B), intent(in) :: varid
548  integer(I4B), intent(in) :: layer
549  integer(I4B), intent(in) :: iper
550  integer(I4B), intent(in) :: iaux
551  character(len=*), intent(in) :: nc_tag
552  character(len=*), intent(in) :: nc_fname
553  !
554  if (nc_tag /= '') then
555  call nf_verify(nf90_put_att(ncid, varid, 'modflow6_input', &
556  nc_tag), nc_fname)
557  if (layer > 0) then
558  call nf_verify(nf90_put_att(ncid, varid, 'modflow6_layer', &
559  layer), nc_fname)
560  end if
561  !
562  if (iper > 0) then
563  call nf_verify(nf90_put_att(ncid, varid, 'modflow6_iper', &
564  iper), nc_fname)
565  end if
566  !
567  if (iaux > 0) then
568  call nf_verify(nf90_put_att(ncid, varid, 'modflow6_iaux', &
569  iaux), nc_fname)
570  end if
571  end if
572  end subroutine ncvar_mf6attr
573 
574  !> @brief build netcdf variable name
575  !<
576  function export_varname(varname, layer, iper, iaux) result(vname)
577  use inputoutputmodule, only: lowcase
578  character(len=*), intent(in) :: varname
579  integer(I4B), optional, intent(in) :: layer
580  integer(I4B), optional, intent(in) :: iper
581  integer(I4B), optional, intent(in) :: iaux
582  character(len=LINELENGTH) :: vname
583  !
584  vname = ''
585  !
586  if (varname /= '') then
587  vname = varname
588  call lowcase(vname)
589  if (present(layer)) then
590  if (layer > 0) then
591  write (vname, '(a,i0)') trim(vname)//'_l', layer
592  end if
593  end if
594  if (present(iper)) then
595  if (iper > 0) then
596  write (vname, '(a,i0)') trim(vname)//'_p', iper
597  end if
598  end if
599  if (present(iaux)) then
600  if (iaux > 0) then
601  write (vname, '(a,i0)') trim(vname)//'a', iaux
602  end if
603  end if
604  end if
605  end function export_varname
606 
607 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, disenum, nctype, 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_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
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:100