MODFLOW 6  version 6.7.0.dev1
USGS Modular Hydrologic Model
meshmodelmodule Module Reference

This module contains the MeshModelModule. More...

Data Types

type  meshncdimidtype
 type for storing model export dimension ids More...
 
type  meshncvaridtype
 type for storing model export variable ids More...
 
type  meshmodeltype
 base ugrid netcdf export type More...
 
interface  nc_array_export_if
 abstract interfaces for derived ugrid netcd export types More...
 
type  mesh2dmodeltype
 

Functions/Subroutines

subroutine mesh_init (this, modelname, modeltype, modelfname, nc_fname, disenum, nctype, lenuni, iout)
 initialize More...
 
subroutine mesh_destroy (this)
 initialize More...
 
subroutine add_global_att (this)
 create file (group) attributes More...
 
subroutine export_input_arrays (this, pkgtype, pkgname, mempath, param_dfns)
 write package gridded input data More...
 
subroutine add_pkg_data (this)
 determine packages to write gridded input More...
 
subroutine define_dependent (this)
 create the model layer dependent variables More...
 
subroutine define_gridmap (this)
 create the file grid mapping container variable More...
 
subroutine create_mesh (this)
 create the file mesh container variable More...
 
subroutine, public ncvar_chunk (ncid, varid, chunk_face, nc_fname)
 define variable chunking More...
 
subroutine, public ncvar_deflate (ncid, varid, deflate, shuffle, nc_fname)
 define variable compression More...
 
subroutine, public ncvar_gridmap (ncid, varid, gridmap_name, nc_fname)
 put variable gridmap attributes More...
 
subroutine, public ncvar_mf6attr (ncid, varid, layer, iper, iaux, nc_tag, nc_fname)
 put variable internal attributes More...
 
character(len=linelength) function, public export_varname (varname, layer, iper, iaux)
 build netcdf variable name More...
 

Detailed Description

This module defines a base class for UGRID based (mesh) model netcdf exports. It is dependent on external netcdf libraries.

Function/Subroutine Documentation

◆ add_global_att()

subroutine meshmodelmodule::add_global_att ( class(meshmodeltype), intent(inout)  this)

Definition at line 160 of file MeshNCModel.f90.

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)
Here is the call graph for this function:

◆ add_pkg_data()

subroutine meshmodelmodule::add_pkg_data ( class(meshmodeltype), intent(inout)  this)

Definition at line 212 of file MeshNCModel.f90.

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)
type(inputparamdefinitiontype) function, dimension(:), pointer, public param_definitions(component, subcomponent)
logical function, public idm_multi_package(component, subcomponent)
character(len=lenmempath) function create_mem_path(component, subcomponent, context)
returns the path to the memory object
This module contains simulation variables.
Definition: SimVariables.f90:9
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
Here is the call graph for this function:

◆ create_mesh()

subroutine meshmodelmodule::create_mesh ( class(mesh2dmodeltype), intent(inout)  this)
private

Definition at line 355 of file MeshNCModel.f90.

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)
Here is the call graph for this function:

◆ define_dependent()

subroutine meshmodelmodule::define_dependent ( class(meshmodeltype), intent(inout)  this)

Definition at line 270 of file MeshNCModel.f90.

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
Here is the call graph for this function:

◆ define_gridmap()

subroutine meshmodelmodule::define_gridmap ( class(meshmodeltype), intent(inout)  this)
private

Definition at line 331 of file MeshNCModel.f90.

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
Here is the call graph for this function:

◆ export_input_arrays()

subroutine meshmodelmodule::export_input_arrays ( class(meshmodeltype), intent(inout)  this,
character(len=*), intent(in)  pkgtype,
character(len=*), intent(in)  pkgname,
character(len=*), intent(in)  mempath,
type(inputparamdefinitiontype), dimension(:), intent(in), pointer  param_dfns 
)
private

Definition at line 185 of file MeshNCModel.f90.

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
subroutine, public get_isize(name, mem_path, isize)
@ brief Get the number of elements for this variable
Here is the call graph for this function:

◆ export_varname()

character(len=linelength) function, public meshmodelmodule::export_varname ( character(len=*), intent(in)  varname,
integer(i4b), intent(in), optional  layer,
integer(i4b), intent(in), optional  iper,
integer(i4b), intent(in), optional  iaux 
)

Definition at line 577 of file MeshNCModel.f90.

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
subroutine, public lowcase(word)
Convert to lower case.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ mesh_destroy()

subroutine meshmodelmodule::mesh_destroy ( class(meshmodeltype), intent(inout)  this)

Definition at line 150 of file MeshNCModel.f90.

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)
Here is the call graph for this function:

◆ mesh_init()

subroutine meshmodelmodule::mesh_init ( class(meshmodeltype), intent(inout)  this,
character(len=*), intent(in)  modelname,
character(len=*), intent(in)  modeltype,
character(len=*), intent(in)  modelfname,
character(len=*), intent(in)  nc_fname,
integer(i4b), intent(in)  disenum,
integer(i4b), intent(in)  nctype,
integer(i4b), intent(in)  lenuni,
integer(i4b), intent(in)  iout 
)
private

Definition at line 98 of file MeshNCModel.f90.

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)
Here is the call graph for this function:

◆ ncvar_chunk()

subroutine, public meshmodelmodule::ncvar_chunk ( integer(i4b), intent(in)  ncid,
integer(i4b), intent(in)  varid,
integer(i4b), intent(in)  chunk_face,
character(len=*), intent(in)  nc_fname 
)

Definition at line 506 of file MeshNCModel.f90.

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
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ncvar_deflate()

subroutine, public meshmodelmodule::ncvar_deflate ( integer(i4b), intent(in)  ncid,
integer(i4b), intent(in)  varid,
integer(i4b), intent(in)  deflate,
integer(i4b), intent(in)  shuffle,
character(len=*), intent(in)  nc_fname 
)

Definition at line 519 of file MeshNCModel.f90.

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
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ncvar_gridmap()

subroutine, public meshmodelmodule::ncvar_gridmap ( integer(i4b), intent(in)  ncid,
integer(i4b), intent(in)  varid,
character(len=*), intent(in)  gridmap_name,
character(len=*), intent(in)  nc_fname 
)

Definition at line 534 of file MeshNCModel.f90.

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
Here is the call graph for this function:
Here is the caller graph for this function:

◆ ncvar_mf6attr()

subroutine, public meshmodelmodule::ncvar_mf6attr ( integer(i4b), intent(in)  ncid,
integer(i4b), intent(in)  varid,
integer(i4b), intent(in)  layer,
integer(i4b), intent(in)  iper,
integer(i4b), intent(in)  iaux,
character(len=*), intent(in)  nc_tag,
character(len=*), intent(in)  nc_fname 
)

Definition at line 549 of file MeshNCModel.f90.

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
Here is the call graph for this function:
Here is the caller graph for this function: