MODFLOW 6  version 6.7.0.dev3
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 df_export (this)
 define timeseries input variables More...
 
subroutine export_df (this, export_pkg)
 define export package More...
 
subroutine create_timeseries (this, idt, iparam, iaux, layer, export_pkg)
 create timeseries export variable 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, iaux, nc_tag, nc_fname)
 put variable internal attributes 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 314 of file MeshNCModel.f90.

315  class(MeshModelType), intent(inout) :: this
316  ! file scoped title
317  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'title', &
318  this%annotation%title), this%nc_fname)
319  ! source (MODFLOW 6)
320  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'source', &
321  this%annotation%source), this%nc_fname)
322  ! grid type (MODFLOW 6)
323  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'modflow_grid', &
324  this%annotation%grid), this%nc_fname)
325  ! mesh type (MODFLOW 6)
326  if (this%annotation%mesh /= '') then
327  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'mesh', &
328  this%annotation%mesh), this%nc_fname)
329 
330  end if
331  ! MODFLOW 6 model type
332  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'modflow_model', &
333  this%annotation%model), this%nc_fname)
334  ! generation datetime
335  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'history', &
336  this%annotation%history), this%nc_fname)
337  ! supported conventions
338  call nf_verify(nf90_put_att(this%ncid, nf90_global, 'Conventions', &
339  this%annotation%conventions), &
340  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 372 of file MeshNCModel.f90.

379  class(MeshModelType), intent(inout) :: this
380  character(LENCOMPONENTNAME) :: ptype, pname, pkgtype
381  character(len=LENMEMPATH) :: input_mempath
382  type(CharacterStringType), dimension(:), contiguous, &
383  pointer :: pkgtypes => null()
384  type(CharacterStringType), dimension(:), contiguous, &
385  pointer :: pkgnames => null()
386  type(CharacterStringType), dimension(:), contiguous, &
387  pointer :: mempaths => null()
388  type(InputParamDefinitionType), dimension(:), pointer :: param_dfns
389  character(len=LENMEMPATH) :: mempath
390  integer(I4B) :: n
391  integer(I4B), pointer :: export_arrays
392  logical(LGP) :: found
393 
394  input_mempath = create_mem_path(component=this%modelname, context=idm_context)
395 
396  ! set pointers to model path package info
397  call mem_setptr(pkgtypes, 'PKGTYPES', input_mempath)
398  call mem_setptr(pkgnames, 'PKGNAMES', input_mempath)
399  call mem_setptr(mempaths, 'MEMPATHS', input_mempath)
400 
401  allocate (export_arrays)
402 
403  do n = 1, size(mempaths)
404  ! initialize export_arrays
405  export_arrays = 0
406 
407  ! set package attributes
408  mempath = mempaths(n)
409  pname = pkgnames(n)
410  ptype = pkgtypes(n)
411 
412  ! export input arrays
413  if (mempath /= '') then
414  ! update export
415  call mem_set_value(export_arrays, 'EXPORT_NC', mempath, found)
416  if (export_arrays > 0) then
417  pkgtype = idm_subcomponent_type(this%modeltype, ptype)
418  param_dfns => param_definitions(this%modeltype, pkgtype)
419  call this%export_input_arrays(ptype, pname, mempath, param_dfns)
420  end if
421  end if
422  end do
423 
424  ! cleanup
425  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 515 of file MeshNCModel.f90.

516  class(Mesh2dModelType), intent(inout) :: this
517 
518  ! create mesh container variable
519  call nf_verify(nf90_def_var(this%ncid, this%mesh_name, nf90_int, &
520  this%var_ids%mesh), this%nc_fname)
521 
522  ! assign container variable attributes
523  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, 'cf_role', &
524  'mesh_topology'), this%nc_fname)
525  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, 'long_name', &
526  '2D mesh topology'), this%nc_fname)
527  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, &
528  'topology_dimension', 2), this%nc_fname)
529  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, 'face_dimension', &
530  'nmesh_face'), this%nc_fname)
531  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, &
532  'node_coordinates', 'mesh_node_x mesh_node_y'), &
533  this%nc_fname)
534  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, &
535  'face_coordinates', 'mesh_face_x mesh_face_y'), &
536  this%nc_fname)
537  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh, &
538  'face_node_connectivity', 'mesh_face_nodes'), &
539  this%nc_fname)
540 
541  ! create mesh x node (mesh vertex) variable
542  call nf_verify(nf90_def_var(this%ncid, 'mesh_node_x', nf90_double, &
543  (/this%dim_ids%nmesh_node/), &
544  this%var_ids%mesh_node_x), this%nc_fname)
545 
546  ! assign mesh x node variable attributes
547  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_x, &
548  'units', this%lenunits), this%nc_fname)
549  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_x, &
550  'standard_name', 'projection_x_coordinate'), &
551  this%nc_fname)
552  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_x, &
553  'long_name', 'Easting'), this%nc_fname)
554 
555  if (this%wkt /= '') then
556  ! associate with projection
557  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_x, &
558  'grid_mapping', this%gridmap_name), &
559  this%nc_fname)
560  end if
561 
562  ! create mesh y node (mesh vertex) variable
563  call nf_verify(nf90_def_var(this%ncid, 'mesh_node_y', nf90_double, &
564  (/this%dim_ids%nmesh_node/), &
565  this%var_ids%mesh_node_y), this%nc_fname)
566 
567  ! assign mesh y variable attributes
568  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_y, &
569  'units', this%lenunits), this%nc_fname)
570  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_y, &
571  'standard_name', 'projection_y_coordinate'), &
572  this%nc_fname)
573  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_y, &
574  'long_name', 'Northing'), this%nc_fname)
575 
576  if (this%wkt /= '') then
577  ! associate with projection
578  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_y, &
579  'grid_mapping', this%gridmap_name), &
580  this%nc_fname)
581  end if
582 
583  ! create mesh x face (cell vertex) variable
584  call nf_verify(nf90_def_var(this%ncid, 'mesh_face_x', nf90_double, &
585  (/this%dim_ids%nmesh_face/), &
586  this%var_ids%mesh_face_x), this%nc_fname)
587 
588  ! assign mesh x face variable attributes
589  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_x, &
590  'units', this%lenunits), this%nc_fname)
591  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_x, &
592  'standard_name', 'projection_x_coordinate'), &
593  this%nc_fname)
594  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_x, &
595  'long_name', 'Easting'), this%nc_fname)
596  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_x, 'bounds', &
597  'mesh_face_xbnds'), this%nc_fname)
598  if (this%wkt /= '') then
599  ! associate with projection
600  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_x, &
601  'grid_mapping', this%gridmap_name), &
602  this%nc_fname)
603  end if
604 
605  ! create mesh x cell bounds variable
606  call nf_verify(nf90_def_var(this%ncid, 'mesh_face_xbnds', nf90_double, &
607  (/this%dim_ids%max_nmesh_face_nodes, &
608  this%dim_ids%nmesh_face/), &
609  this%var_ids%mesh_face_xbnds), &
610  this%nc_fname)
611 
612  ! create mesh y face (cell vertex) variable
613  call nf_verify(nf90_def_var(this%ncid, 'mesh_face_y', nf90_double, &
614  (/this%dim_ids%nmesh_face/), &
615  this%var_ids%mesh_face_y), this%nc_fname)
616 
617  ! assign mesh y face variable attributes
618  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_y, &
619  'units', this%lenunits), this%nc_fname)
620  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_y, &
621  'standard_name', 'projection_y_coordinate'), &
622  this%nc_fname)
623  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_y, &
624  'long_name', 'Northing'), this%nc_fname)
625  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_y, 'bounds', &
626  'mesh_face_ybnds'), this%nc_fname)
627 
628  if (this%wkt /= '') then
629  ! associate with projection
630  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_y, &
631  'grid_mapping', this%gridmap_name), &
632  this%nc_fname)
633  end if
634 
635  ! create mesh y cell bounds variable
636  call nf_verify(nf90_def_var(this%ncid, 'mesh_face_ybnds', nf90_double, &
637  (/this%dim_ids%max_nmesh_face_nodes, &
638  this%dim_ids%nmesh_face/), &
639  this%var_ids%mesh_face_ybnds), &
640  this%nc_fname)
641 
642  ! create mesh face nodes variable
643  call nf_verify(nf90_def_var(this%ncid, 'mesh_face_nodes', nf90_int, &
644  (/this%dim_ids%max_nmesh_face_nodes, &
645  this%dim_ids%nmesh_face/), &
646  this%var_ids%mesh_face_nodes), &
647  this%nc_fname)
648 
649  ! assign variable attributes
650  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_nodes, &
651  'cf_role', 'face_node_connectivity'), &
652  this%nc_fname)
653  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_nodes, &
654  'long_name', &
655  'Vertices bounding cell (counterclockwise)'), &
656  this%nc_fname)
657  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_nodes, &
658  '_FillValue', (/nf90_fill_int/)), &
659  this%nc_fname)
660  call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_nodes, &
661  'start_index', 1), this%nc_fname)
Here is the call graph for this function:

◆ create_timeseries()

subroutine meshmodelmodule::create_timeseries ( class(meshmodeltype), intent(inout)  this,
type(inputparamdefinitiontype), intent(in), pointer  idt,
integer(i4b), intent(in)  iparam,
integer(i4b), intent(in)  iaux,
integer(i4b), intent(in)  layer,
class(exportpackagetype), intent(in), pointer  export_pkg 
)

Definition at line 220 of file MeshNCModel.f90.

221  use constantsmodule, only: dnodata
223  class(MeshModelType), intent(inout) :: this
224  type(InputParamDefinitionType), pointer, intent(in) :: idt
225  integer(I4B), intent(in) :: iparam
226  integer(I4B), intent(in) :: iaux
227  integer(I4B), intent(in) :: layer
228  class(ExportPackageType), pointer, intent(in) :: export_pkg
229  character(len=LINELENGTH) :: varname, longname, nc_tag
230  integer(I4B) :: varid
231 
232  ! set variable input tag
233  nc_tag = this%input_attribute(export_pkg%mf6_input%subcomponent_name, &
234  idt)
235 
236  ! set names
237  varname = export_varname(export_pkg%mf6_input%subcomponent_name, &
238  idt%tagname, export_pkg%mf6_input%mempath, &
239  layer=layer, iaux=iaux)
240  longname = export_longname(idt%longname, &
241  export_pkg%mf6_input%subcomponent_name, &
242  idt%tagname, export_pkg%mf6_input%mempath, &
243  layer=layer, iaux=iaux)
244 
245  ! create the netcdf dependent layer variable
246  select case (idt%datatype)
247  case ('DOUBLE1D', 'DOUBLE2D')
248  call nf_verify(nf90_def_var(this%ncid, varname, nf90_double, &
249  (/this%dim_ids%nmesh_face, &
250  this%dim_ids%time/), &
251  varid), &
252  this%nc_fname)
253  call nf_verify(nf90_put_att(this%ncid, varid, &
254  '_FillValue', (/dnodata/)), &
255  this%nc_fname)
256  case ('INTEGER1D')
257  call nf_verify(nf90_def_var(this%ncid, varname, nf90_int, &
258  (/this%dim_ids%nmesh_face, &
259  this%dim_ids%time/), &
260  varid), &
261  this%nc_fname)
262  call nf_verify(nf90_put_att(this%ncid, varid, &
263  '_FillValue', (/nf90_fill_int/)), &
264  this%nc_fname)
265  end select
266 
267  ! apply chunking parameters
268  if (this%chunking_active) then
269  call nf_verify(nf90_def_var_chunking(this%ncid, &
270  varid, &
271  nf90_chunked, &
272  (/this%chunk_face, &
273  this%chunk_time/)), &
274  this%nc_fname)
275  end if
276 
277  ! deflate and shuffle
278  call ncvar_deflate(this%ncid, varid, this%deflate, &
279  this%shuffle, this%nc_fname)
280 
281  ! assign variable attributes
282  call nf_verify(nf90_put_att(this%ncid, varid, &
283  'units', this%lenunits), this%nc_fname)
284  call nf_verify(nf90_put_att(this%ncid, varid, &
285  'long_name', longname), this%nc_fname)
286  call nf_verify(nf90_put_att(this%ncid, varid, &
287  'mesh', this%mesh_name), this%nc_fname)
288  call nf_verify(nf90_put_att(this%ncid, varid, &
289  'location', 'face'), this%nc_fname)
290 
291  ! add grid mapping and mf6 attr
292  call ncvar_gridmap(this%ncid, varid, &
293  this%gridmap_name, this%nc_fname)
294  call ncvar_mf6attr(this%ncid, varid, layer, iaux, nc_tag, this%nc_fname)
295 
296  ! store variable id
297  if (idt%tagname == 'AUX') then
298  if (layer > 0) then
299  export_pkg%varids_aux(iaux, layer) = varid
300  else
301  export_pkg%varids_aux(iaux, 1) = varid
302  end if
303  else
304  if (layer > 0) then
305  export_pkg%varids_param(iparam, layer) = varid
306  else
307  export_pkg%varids_param(iparam, 1) = varid
308  end if
309  end if
This module contains simulation constants.
Definition: Constants.f90:9
real(dp), parameter dnodata
real no data constant
Definition: Constants.f90:95
This module contains the NCModelExportModule.
Definition: NCModel.f90:8
Here is the call graph for this function:

◆ define_dependent()

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

Definition at line 430 of file MeshNCModel.f90.

431  class(MeshModelType), intent(inout) :: this
432  character(len=LINELENGTH) :: varname, longname
433  integer(I4B) :: k
434 
435  ! create a dependent variable for each layer
436  do k = 1, this%nlay
437  ! initialize names
438  varname = ''
439  longname = ''
440 
441  ! set layer variable and longnames
442  write (varname, '(a,i0)') trim(this%xname)//'_l', k
443  write (longname, '(a,i0,a)') trim(this%annotation%longname)// &
444  ' (layer ', k, ')'
445 
446  ! create the netcdf dependent layer variable
447  call nf_verify(nf90_def_var(this%ncid, varname, nf90_double, &
448  (/this%dim_ids%nmesh_face, &
449  this%dim_ids%time/), &
450  this%var_ids%dependent(k)), &
451  this%nc_fname)
452 
453  ! apply chunking parameters
454  if (this%chunking_active) then
455  call nf_verify(nf90_def_var_chunking(this%ncid, &
456  this%var_ids%dependent(k), &
457  nf90_chunked, &
458  (/this%chunk_face, &
459  this%chunk_time/)), &
460  this%nc_fname)
461  end if
462 
463  ! deflate and shuffle
464  call ncvar_deflate(this%ncid, this%var_ids%dependent(k), this%deflate, &
465  this%shuffle, this%nc_fname)
466 
467  ! assign variable attributes
468  call nf_verify(nf90_put_att(this%ncid, this%var_ids%dependent(k), &
469  'units', this%lenunits), this%nc_fname)
470  call nf_verify(nf90_put_att(this%ncid, this%var_ids%dependent(k), &
471  'standard_name', this%annotation%stdname), &
472  this%nc_fname)
473  call nf_verify(nf90_put_att(this%ncid, this%var_ids%dependent(k), &
474  'long_name', longname), this%nc_fname)
475  call nf_verify(nf90_put_att(this%ncid, this%var_ids%dependent(k), &
476  '_FillValue', (/dhnoflo/)), &
477  this%nc_fname)
478  call nf_verify(nf90_put_att(this%ncid, this%var_ids%dependent(k), &
479  'mesh', this%mesh_name), this%nc_fname)
480  call nf_verify(nf90_put_att(this%ncid, this%var_ids%dependent(k), &
481  'location', 'face'), this%nc_fname)
482 
483  ! add grid mapping
484  call ncvar_gridmap(this%ncid, this%var_ids%dependent(k), &
485  this%gridmap_name, this%nc_fname)
486  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 491 of file MeshNCModel.f90.

492  class(MeshModelType), intent(inout) :: this
493  integer(I4B) :: var_id
494 
495  ! was projection info provided
496  if (this%wkt /= '') then
497  ! create projection variable
498  call nf_verify(nf90_redef(this%ncid), this%nc_fname)
499  call nf_verify(nf90_def_var(this%ncid, this%gridmap_name, nf90_int, &
500  var_id), this%nc_fname)
501  ! cf-conventions prefers 'crs_wkt'
502  !call nf_verify(nf90_put_att(this%ncid, var_id, 'crs_wkt', this%wkt), &
503  ! this%nc_fname)
504  ! QGIS recognizes 'wkt'
505  call nf_verify(nf90_put_att(this%ncid, var_id, 'wkt', this%wkt), &
506  this%nc_fname)
507  call nf_verify(nf90_enddef(this%ncid), this%nc_fname)
508  call nf_verify(nf90_put_var(this%ncid, var_id, 1), &
509  this%nc_fname)
510  end if
Here is the call graph for this function:

◆ df_export()

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

Definition at line 163 of file MeshNCModel.f90.

165  class(MeshModelType), intent(inout) :: this
166  class(ExportPackageType), pointer :: export_pkg
167  integer(I4B) :: idx
168  do idx = 1, this%pkglist%Count()
169  export_pkg => this%get(idx)
170  call this%export_df(export_pkg)
171  end do

◆ export_df()

subroutine meshmodelmodule::export_df ( class(meshmodeltype), intent(inout)  this,
class(exportpackagetype), intent(in), pointer  export_pkg 
)

Definition at line 176 of file MeshNCModel.f90.

179  class(MeshModelType), intent(inout) :: this
180  class(ExportPackageType), pointer, intent(in) :: export_pkg
181  type(InputParamDefinitionType), pointer :: idt
182  integer(I4B) :: iparam, iaux, layer
183 
184  ! export defined period input
185  do iparam = 1, export_pkg%nparam
186  ! initialize
187  iaux = 0
188  layer = 0
189  ! set input definition
190  idt => &
191  get_param_definition_type(export_pkg%mf6_input%param_dfns, &
192  export_pkg%mf6_input%component_type, &
193  export_pkg%mf6_input%subcomponent_type, &
194  'PERIOD', export_pkg%param_names(iparam), '')
195 
196  select case (idt%shape)
197  case ('NCPL')
198  call this%create_timeseries(idt, iparam, iaux, layer, export_pkg)
199  case ('NODES')
200  do layer = 1, this%nlay
201  call this%create_timeseries(idt, iparam, iaux, layer, export_pkg)
202  end do
203  case ('NAUX NCPL')
204  do iaux = 1, export_pkg%naux
205  call this%create_timeseries(idt, iparam, iaux, layer, export_pkg)
206  end do
207  case ('NAUX NODES')
208  do iaux = 1, export_pkg%naux
209  do layer = 1, this%nlay
210  call this%create_timeseries(idt, iparam, iaux, layer, export_pkg)
211  end do
212  end do
213  case default
214  end select
215  end do
This module contains the DefinitionSelectModule.
type(inputparamdefinitiontype) function, pointer, public get_param_definition_type(input_definition_types, component_type, subcomponent_type, blockname, tagname, filename)
Return parameter definition.
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 345 of file MeshNCModel.f90.

346  use memorymanagermodule, only: get_isize
347  class(MeshModelType), intent(inout) :: this
348  character(len=*), intent(in) :: pkgtype
349  character(len=*), intent(in) :: pkgname
350  character(len=*), intent(in) :: mempath
351  type(InputParamDefinitionType), dimension(:), pointer, &
352  intent(in) :: param_dfns
353  type(InputParamDefinitionType), pointer :: idt
354  integer(I4B) :: iparam, isize
355  ! export griddata block parameters
356  do iparam = 1, size(param_dfns)
357  ! assign param definition pointer
358  idt => param_dfns(iparam)
359  ! for now only griddata is exported
360  if (idt%blockname == 'GRIDDATA') then
361  ! veriy variable is allocated
362  call get_isize(idt%mf6varname, mempath, isize)
363  if (isize > 0) then
364  call this%export_input_array(pkgtype, pkgname, mempath, idt)
365  end if
366  end if
367  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:

◆ mesh_destroy()

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

Definition at line 153 of file MeshNCModel.f90.

155  class(MeshModelType), intent(inout) :: this
156  call nf_verify(nf90_close(this%ncid), this%nc_fname)
157  deallocate (this%chunk_face)
158  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 101 of file MeshNCModel.f90.

104  class(MeshModelType), intent(inout) :: this
105  character(len=*), intent(in) :: modelname
106  character(len=*), intent(in) :: modeltype
107  character(len=*), intent(in) :: modelfname
108  character(len=*), intent(in) :: nc_fname
109  integer(I4B), intent(in) :: disenum
110  integer(I4B), intent(in) :: nctype
111  integer(I4B), intent(in) :: lenuni
112  integer(I4B), intent(in) :: iout
113  logical(LGP) :: found
114 
115  ! initialize base class
116  call this%NCModelExportType%init(modelname, modeltype, modelfname, nc_fname, &
117  disenum, nctype, iout)
118 
119  ! allocate and initialize
120  allocate (this%chunk_face)
121  this%chunk_face = -1
122 
123  ! update values from input context
124  if (this%ncf_mempath /= '') then
125  call mem_set_value(this%chunk_face, 'CHUNK_FACE', this%ncf_mempath, found)
126  end if
127 
128  if (this%chunk_time > 0 .and. this%chunk_face > 0) then
129  this%chunking_active = .true.
130  else if (this%chunk_time > 0 .or. this%chunk_face > 0) then
131  this%chunk_face = -1
132  this%chunk_time = -1
133  write (warnmsg, '(a)') 'Ignoring user provided NetCDF chunking parameter. &
134  &Define chunk_time and chunk_face input parameters to see an effect in &
135  &file "'//trim(nc_fname)//'".'
136  call store_warning(warnmsg)
137  end if
138 
139  if (lenuni == 1) then
140  this%lenunits = 'ft'
141  else
142  this%lenunits = 'm'
143  end if
144 
145  ! create the netcdf file
146  call nf_verify(nf90_create(this%nc_fname, &
147  ior(nf90_clobber, nf90_netcdf4), this%ncid), &
148  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 666 of file MeshNCModel.f90.

667  integer(I4B), intent(in) :: ncid
668  integer(I4B), intent(in) :: varid
669  integer(I4B), intent(in) :: chunk_face
670  character(len=*), intent(in) :: nc_fname
671  if (chunk_face > 0) then
672  call nf_verify(nf90_def_var_chunking(ncid, varid, nf90_chunked, &
673  (/chunk_face/)), nc_fname)
674  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 679 of file MeshNCModel.f90.

680  integer(I4B), intent(in) :: ncid
681  integer(I4B), intent(in) :: varid
682  integer(I4B), intent(in) :: deflate
683  integer(I4B), intent(in) :: shuffle
684  character(len=*), intent(in) :: nc_fname
685  if (deflate >= 0) then
686  call nf_verify(nf90_def_var_deflate(ncid, varid, shuffle=shuffle, &
687  deflate=1, deflate_level=deflate), &
688  nc_fname)
689  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 694 of file MeshNCModel.f90.

695  integer(I4B), intent(in) :: ncid
696  integer(I4B), intent(in) :: varid
697  character(len=*), intent(in) :: gridmap_name
698  character(len=*), intent(in) :: nc_fname
699  if (gridmap_name /= '') then
700  call nf_verify(nf90_put_att(ncid, varid, 'coordinates', &
701  'mesh_face_x mesh_face_y'), nc_fname)
702  call nf_verify(nf90_put_att(ncid, varid, 'grid_mapping', &
703  gridmap_name), nc_fname)
704  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)  iaux,
character(len=*), intent(in)  nc_tag,
character(len=*), intent(in)  nc_fname 
)

Definition at line 709 of file MeshNCModel.f90.

710  integer(I4B), intent(in) :: ncid
711  integer(I4B), intent(in) :: varid
712  integer(I4B), intent(in) :: layer
713  integer(I4B), intent(in) :: iaux
714  character(len=*), intent(in) :: nc_tag
715  character(len=*), intent(in) :: nc_fname
716  if (nc_tag /= '') then
717  call nf_verify(nf90_put_att(ncid, varid, 'modflow_input', &
718  nc_tag), nc_fname)
719  if (layer > 0) then
720  call nf_verify(nf90_put_att(ncid, varid, 'layer', &
721  layer), nc_fname)
722  end if
723  if (iaux > 0) then
724  call nf_verify(nf90_put_att(ncid, varid, 'modflow_iaux', &
725  iaux), nc_fname)
726  end if
727  end if
Here is the call graph for this function:
Here is the caller graph for this function: