From 682972b10de77c4fc034be922c511d384d2e0cad Mon Sep 17 00:00:00 2001 From: langevin-usgs Date: Fri, 28 Oct 2022 15:21:09 -0500 Subject: [PATCH] refactor(arrayreader): working on object-based array reader (#1067) --- doc/mf6io/mf6ivar/dfn/gwf-disv.dfn | 6 +- make/makefile | 50 ++- msvs/mf6core.vfproj | 7 + pymake/makefile | 299 +++++++------ src/Model/GroundWaterFlow/gwf3dis8idm.f90 | 6 +- src/Model/GroundWaterFlow/gwf3disv8.f90 | 54 +-- src/Model/GroundWaterFlow/gwf3disv8idm.f90 | 12 +- src/Utilities/ArrayRead/ArrayReaderBase.f90 | 186 ++++++++ src/Utilities/ArrayRead/Double1dReader.f90 | 116 +++++ src/Utilities/ArrayRead/Double2dReader.f90 | 123 ++++++ src/Utilities/ArrayRead/Integer1dReader.f90 | 148 +++++++ src/Utilities/ArrayRead/Integer2dReader.f90 | 123 ++++++ .../ArrayRead/LayeredArrayReader.f90 | 162 +++++++ src/Utilities/ArrayReaders.f90 | 16 +- src/Utilities/Idm/LoadMf6FileType.f90 | 403 ++++++++++-------- src/meson.build | 6 + utils/idmloader/scripts/dfn2f90.py | 1 + utils/mf5to6/make/makefile | 8 +- 18 files changed, 1366 insertions(+), 360 deletions(-) create mode 100644 src/Utilities/ArrayRead/ArrayReaderBase.f90 create mode 100644 src/Utilities/ArrayRead/Double1dReader.f90 create mode 100644 src/Utilities/ArrayRead/Double2dReader.f90 create mode 100644 src/Utilities/ArrayRead/Integer1dReader.f90 create mode 100644 src/Utilities/ArrayRead/Integer2dReader.f90 create mode 100644 src/Utilities/ArrayRead/LayeredArrayReader.f90 diff --git a/doc/mf6io/mf6ivar/dfn/gwf-disv.dfn b/doc/mf6io/mf6ivar/dfn/gwf-disv.dfn index a22fe1c9f03..c6bab33389e 100644 --- a/doc/mf6io/mf6ivar/dfn/gwf-disv.dfn +++ b/doc/mf6io/mf6ivar/dfn/gwf-disv.dfn @@ -71,7 +71,7 @@ description is the total number of (x, y) vertex pairs used to characterize the block griddata name top type double precision -shape (ncpl, 1) +shape (ncpl) reader readarray longname model top elevation description is the top elevation for each cell in the top model layer. @@ -79,7 +79,7 @@ description is the top elevation for each cell in the top model layer. block griddata name botm type double precision -shape (ncpl, 1, nlay) +shape (ncpl, nlay) reader readarray layered true longname model bottom elevation @@ -88,7 +88,7 @@ description is the bottom elevation for each cell. block griddata name idomain type integer -shape (ncpl, 1, nlay) +shape (ncpl, nlay) reader readarray layered true optional true diff --git a/make/makefile b/make/makefile index 83a6bc5c5fb..e17d8abd2bd 100644 --- a/make/makefile +++ b/make/makefile @@ -6,27 +6,28 @@ include ./makedefaults # Define the source file directories SOURCEDIR1=../src SOURCEDIR2=../src/Exchange -SOURCEDIR3=../src/Model -SOURCEDIR4=../src/Model/Connection -SOURCEDIR5=../src/Model/Geometry -SOURCEDIR6=../src/Model/GroundWaterFlow -SOURCEDIR7=../src/Model/GroundWaterTransport -SOURCEDIR8=../src/Model/ModelUtilities -SOURCEDIR9=../src/Solution -SOURCEDIR10=../src/Solution/LinearMethods -SOURCEDIR11=../src/Timing -SOURCEDIR12=../src/Utilities -SOURCEDIR13=../src/Utilities/Idm -SOURCEDIR14=../src/Utilities/Libraries -SOURCEDIR15=../src/Utilities/Libraries/blas +SOURCEDIR3=../src/Solution +SOURCEDIR4=../src/Solution/LinearMethods +SOURCEDIR5=../src/Timing +SOURCEDIR6=../src/Utilities +SOURCEDIR7=../src/Utilities/Idm +SOURCEDIR8=../src/Utilities/TimeSeries +SOURCEDIR9=../src/Utilities/Memory +SOURCEDIR10=../src/Utilities/OutputControl +SOURCEDIR11=../src/Utilities/ArrayRead +SOURCEDIR12=../src/Utilities/Libraries +SOURCEDIR13=../src/Utilities/Libraries/rcm +SOURCEDIR14=../src/Utilities/Libraries/blas +SOURCEDIR15=../src/Utilities/Libraries/sparskit2 SOURCEDIR16=../src/Utilities/Libraries/daglib -SOURCEDIR17=../src/Utilities/Libraries/rcm -SOURCEDIR18=../src/Utilities/Libraries/sparsekit -SOURCEDIR19=../src/Utilities/Libraries/sparskit2 -SOURCEDIR20=../src/Utilities/Memory -SOURCEDIR21=../src/Utilities/Observation -SOURCEDIR22=../src/Utilities/OutputControl -SOURCEDIR23=../src/Utilities/TimeSeries +SOURCEDIR17=../src/Utilities/Libraries/sparsekit +SOURCEDIR18=../src/Utilities/Observation +SOURCEDIR19=../src/Model +SOURCEDIR20=../src/Model/Connection +SOURCEDIR21=../src/Model/GroundWaterTransport +SOURCEDIR22=../src/Model/ModelUtilities +SOURCEDIR23=../src/Model/GroundWaterFlow +SOURCEDIR24=../src/Model/Geometry VPATH = \ ${SOURCEDIR1} \ @@ -51,7 +52,8 @@ ${SOURCEDIR19} \ ${SOURCEDIR20} \ ${SOURCEDIR21} \ ${SOURCEDIR22} \ -${SOURCEDIR23} +${SOURCEDIR23} \ +${SOURCEDIR24} .SUFFIXES: .f90 .F90 .o @@ -98,6 +100,7 @@ $(OBJDIR)/TimeArraySeries.o \ $(OBJDIR)/ObsOutputList.o \ $(OBJDIR)/Observe.o \ $(OBJDIR)/InputDefinition.o \ +$(OBJDIR)/ArrayReaderBase.o \ $(OBJDIR)/TimeArraySeriesLink.o \ $(OBJDIR)/ObsUtility.o \ $(OBJDIR)/ObsContainer.o \ @@ -107,6 +110,7 @@ $(OBJDIR)/gwf3npf8idm.o \ $(OBJDIR)/gwf3disv8idm.o \ $(OBJDIR)/gwf3disu8idm.o \ $(OBJDIR)/gwf3dis8idm.o \ +$(OBJDIR)/Integer2dReader.o \ $(OBJDIR)/TimeArraySeriesManager.o \ $(OBJDIR)/PackageMover.o \ $(OBJDIR)/Obs3.o \ @@ -116,11 +120,15 @@ $(OBJDIR)/BudgetFileReader.o \ $(OBJDIR)/StructVector.o \ $(OBJDIR)/IdmLogger.o \ $(OBJDIR)/InputDefinitionSelector.o \ +$(OBJDIR)/Integer1dReader.o \ +$(OBJDIR)/Double2dReader.o \ +$(OBJDIR)/Double1dReader.o \ $(OBJDIR)/BoundaryPackage.o \ $(OBJDIR)/BaseModel.o \ $(OBJDIR)/BudgetTerm.o \ $(OBJDIR)/StructArray.o \ $(OBJDIR)/ModflowInput.o \ +$(OBJDIR)/LayeredArrayReader.o \ $(OBJDIR)/NumericalModel.o \ $(OBJDIR)/mf6lists.o \ $(OBJDIR)/PackageBudget.o \ diff --git a/msvs/mf6core.vfproj b/msvs/mf6core.vfproj index 975d2eb3bd0..547d51d3ee3 100644 --- a/msvs/mf6core.vfproj +++ b/msvs/mf6core.vfproj @@ -162,6 +162,13 @@ + + + + + + + diff --git a/pymake/makefile b/pymake/makefile index d289c3ea4fd..26cb285a460 100644 --- a/pymake/makefile +++ b/pymake/makefile @@ -1,6 +1,4 @@ -# makefile created on 2021-10-12 16:13:16.130574 -# by pymake (version 1.2.1) for the 'mf6' executable -# using the 'gfortran' fortran compiler(s). +# makefile created by pymake (version 1.2.5) for the 'mf6' executable. include ./makedefaults @@ -9,18 +7,27 @@ include ./makedefaults SOURCEDIR1=../src SOURCEDIR2=../src/Exchange SOURCEDIR3=../src/Solution -SOURCEDIR4=../src/Solution/SparseMatrixSolver +SOURCEDIR4=../src/Solution/LinearMethods SOURCEDIR5=../src/Timing SOURCEDIR6=../src/Utilities -SOURCEDIR7=../src/Utilities/TimeSeries -SOURCEDIR8=../src/Utilities/Memory -SOURCEDIR9=../src/Utilities/OutputControl -SOURCEDIR10=../src/Utilities/Observation -SOURCEDIR11=../src/Model -SOURCEDIR12=../src/Model/GroundWaterTransport -SOURCEDIR13=../src/Model/ModelUtilities -SOURCEDIR14=../src/Model/GroundWaterFlow -SOURCEDIR15=../src/Model/Geometry +SOURCEDIR7=../src/Utilities/Idm +SOURCEDIR8=../src/Utilities/TimeSeries +SOURCEDIR9=../src/Utilities/Memory +SOURCEDIR10=../src/Utilities/OutputControl +SOURCEDIR11=../src/Utilities/ArrayRead +SOURCEDIR12=../src/Utilities/Libraries +SOURCEDIR13=../src/Utilities/Libraries/rcm +SOURCEDIR14=../src/Utilities/Libraries/blas +SOURCEDIR15=../src/Utilities/Libraries/sparskit2 +SOURCEDIR16=../src/Utilities/Libraries/daglib +SOURCEDIR17=../src/Utilities/Libraries/sparsekit +SOURCEDIR18=../src/Utilities/Observation +SOURCEDIR19=../src/Model +SOURCEDIR20=../src/Model/Connection +SOURCEDIR21=../src/Model/GroundWaterTransport +SOURCEDIR22=../src/Model/ModelUtilities +SOURCEDIR23=../src/Model/GroundWaterFlow +SOURCEDIR24=../src/Model/Geometry VPATH = \ ${SOURCEDIR1} \ @@ -37,153 +44,207 @@ ${SOURCEDIR11} \ ${SOURCEDIR12} \ ${SOURCEDIR13} \ ${SOURCEDIR14} \ -${SOURCEDIR15} +${SOURCEDIR15} \ +${SOURCEDIR16} \ +${SOURCEDIR17} \ +${SOURCEDIR18} \ +${SOURCEDIR19} \ +${SOURCEDIR20} \ +${SOURCEDIR21} \ +${SOURCEDIR22} \ +${SOURCEDIR23} \ +${SOURCEDIR24} .SUFFIXES: .f90 .F90 .o OBJECTS = \ -$(OBJDIR)/OpenSpec.o \ $(OBJDIR)/kind.o \ -$(OBJDIR)/HashTable.o \ -$(OBJDIR)/compilerversion.o \ $(OBJDIR)/Constants.o \ -$(OBJDIR)/BaseGeometry.o \ $(OBJDIR)/SimVariables.o \ -$(OBJDIR)/GwfNpfOptions.o \ -$(OBJDIR)/ims8reordering.o \ -$(OBJDIR)/Sparse.o \ -$(OBJDIR)/GwfNpfGridData.o \ -$(OBJDIR)/PackageBudget.o \ -$(OBJDIR)/GwfStorageUtils.o \ $(OBJDIR)/genericutils.o \ -$(OBJDIR)/defmacro.o \ -$(OBJDIR)/ims8sparsekit.o \ +$(OBJDIR)/compilerversion.o \ $(OBJDIR)/ArrayHandlers.o \ -$(OBJDIR)/Xt3dAlgorithm.o \ $(OBJDIR)/version.o \ -$(OBJDIR)/SfrCrossSectionUtils.o \ -$(OBJDIR)/SmoothingFunctions.o \ -$(OBJDIR)/List.o \ -$(OBJDIR)/Timer.o \ -$(OBJDIR)/StringList.o \ -$(OBJDIR)/TimeSeriesRecord.o \ -$(OBJDIR)/mf6lists.o \ $(OBJDIR)/Message.o \ -$(OBJDIR)/ObsOutput.o \ +$(OBJDIR)/defmacro.o \ $(OBJDIR)/Sim.o \ -$(OBJDIR)/sort.o \ -$(OBJDIR)/Budget.o \ +$(OBJDIR)/OpenSpec.o \ $(OBJDIR)/InputOutput.o \ -$(OBJDIR)/RectangularGeometry.o \ -$(OBJDIR)/BlockParser.o \ -$(OBJDIR)/Iunit.o \ -$(OBJDIR)/MemoryHelper.o \ -$(OBJDIR)/NameFile.o \ -$(OBJDIR)/ArrayReaders.o \ -$(OBJDIR)/comarg.o \ -$(OBJDIR)/SfrCrossSectionManager.o \ -$(OBJDIR)/ObsOutputList.o \ -$(OBJDIR)/HeadFileReader.o \ -$(OBJDIR)/DisvGeom.o \ -$(OBJDIR)/BudgetFileReader.o \ $(OBJDIR)/TableTerm.o \ -$(OBJDIR)/CircularGeometry.o \ -$(OBJDIR)/TimeSeries.o \ -$(OBJDIR)/PrintSaveManager.o \ -$(OBJDIR)/ims8base.o \ -$(OBJDIR)/TimeSeriesLink.o \ $(OBJDIR)/Table.o \ -$(OBJDIR)/ListReader.o \ -$(OBJDIR)/TimeSeriesFileList.o \ +$(OBJDIR)/MemoryHelper.o \ +$(OBJDIR)/CharString.o \ $(OBJDIR)/Memory.o \ +$(OBJDIR)/List.o \ $(OBJDIR)/MemoryList.o \ +$(OBJDIR)/TimeSeriesRecord.o \ +$(OBJDIR)/BlockParser.o \ $(OBJDIR)/MemoryManager.o \ +$(OBJDIR)/TimeSeries.o \ $(OBJDIR)/ats.o \ -$(OBJDIR)/BaseModel.o \ -$(OBJDIR)/PackageMover.o \ +$(OBJDIR)/TimeSeriesLink.o \ +$(OBJDIR)/TimeSeriesFileList.o \ $(OBJDIR)/tdis.o \ -$(OBJDIR)/ims8linear.o \ -$(OBJDIR)/Connections.o \ -$(OBJDIR)/MemorySetHandler.o \ -$(OBJDIR)/Mover.o \ +$(OBJDIR)/HashTable.o \ +$(OBJDIR)/Sparse.o \ +$(OBJDIR)/DisvGeom.o \ +$(OBJDIR)/ArrayReaders.o \ $(OBJDIR)/TimeSeriesManager.o \ -$(OBJDIR)/UzfCellGroup.o \ -$(OBJDIR)/BaseExchange.o \ +$(OBJDIR)/SmoothingFunctions.o \ +$(OBJDIR)/ListReader.o \ +$(OBJDIR)/Connections.o \ $(OBJDIR)/DiscretizationBase.o \ -$(OBJDIR)/gwf3dis8.o \ -$(OBJDIR)/BudgetTerm.o \ -$(OBJDIR)/Observe.o \ -$(OBJDIR)/gwf3disu8.o \ -$(OBJDIR)/BaseSolution.o \ -$(OBJDIR)/gwf3disv8.o \ -$(OBJDIR)/Xt3dInterface.o \ -$(OBJDIR)/BudgetObject.o \ $(OBJDIR)/TimeArray.o \ +$(OBJDIR)/ObsOutput.o \ +$(OBJDIR)/TimeArraySeries.o \ +$(OBJDIR)/ObsOutputList.o \ +$(OBJDIR)/Observe.o \ +$(OBJDIR)/InputDefinition.o \ +$(OBJDIR)/TimeArraySeriesLink.o \ +$(OBJDIR)/ObsUtility.o \ +$(OBJDIR)/ObsContainer.o \ +$(OBJDIR)/VectorInt.o \ +$(OBJDIR)/gwt1dspidm.o \ +$(OBJDIR)/gwf3npf8idm.o \ +$(OBJDIR)/gwf3disv8idm.o \ +$(OBJDIR)/gwf3disu8idm.o \ +$(OBJDIR)/gwf3dis8idm.o \ +$(OBJDIR)/ArrayReaderBase.o \ +$(OBJDIR)/TimeArraySeriesManager.o \ +$(OBJDIR)/PackageMover.o \ +$(OBJDIR)/Obs3.o \ $(OBJDIR)/NumericalPackage.o \ +$(OBJDIR)/Budget.o \ +$(OBJDIR)/BudgetFileReader.o \ +$(OBJDIR)/StructVector.o \ +$(OBJDIR)/IdmLogger.o \ +$(OBJDIR)/InputDefinitionSelector.o \ +$(OBJDIR)/Integer2dReader.o \ +$(OBJDIR)/Double2dReader.o \ +$(OBJDIR)/BoundaryPackage.o \ +$(OBJDIR)/BaseModel.o \ +$(OBJDIR)/BudgetTerm.o \ +$(OBJDIR)/StructArray.o \ +$(OBJDIR)/ModflowInput.o \ +$(OBJDIR)/Integer1dReader.o \ +$(OBJDIR)/Double1dReader.o \ +$(OBJDIR)/NumericalModel.o \ +$(OBJDIR)/mf6lists.o \ +$(OBJDIR)/PackageBudget.o \ +$(OBJDIR)/HeadFileReader.o \ +$(OBJDIR)/BudgetObject.o \ +$(OBJDIR)/sort.o \ +$(OBJDIR)/SfrCrossSectionUtils.o \ +$(OBJDIR)/PrintSaveManager.o \ +$(OBJDIR)/Xt3dAlgorithm.o \ +$(OBJDIR)/gwf3tvbase8.o \ +$(OBJDIR)/LoadMf6FileType.o \ +$(OBJDIR)/DistributedModel.o \ +$(OBJDIR)/BaseExchange.o \ +$(OBJDIR)/UzfCellGroup.o \ +$(OBJDIR)/gwt1fmi1.o \ +$(OBJDIR)/SfrCrossSectionManager.o \ +$(OBJDIR)/dag_module.o \ $(OBJDIR)/OutputControlData.o \ -$(OBJDIR)/ObsContainer.o \ -$(OBJDIR)/ObsUtility.o \ -$(OBJDIR)/gwf3mvr8.o \ -$(OBJDIR)/gwf3hfb8.o \ -$(OBJDIR)/SolutionGroup.o \ -$(OBJDIR)/TimeArraySeries.o \ $(OBJDIR)/gwf3ic8.o \ -$(OBJDIR)/Obs3.o \ +$(OBJDIR)/Xt3dInterface.o \ +$(OBJDIR)/gwf3tvk8.o \ +$(OBJDIR)/MemoryManagerExt.o \ +$(OBJDIR)/IdmMf6FileLoader.o \ +$(OBJDIR)/GwfNpfOptions.o \ +$(OBJDIR)/CellWithNbrs.o \ +$(OBJDIR)/NumericalExchange.o \ +$(OBJDIR)/Iunit.o \ +$(OBJDIR)/gwf3uzf8.o \ +$(OBJDIR)/gwt1apt1.o \ +$(OBJDIR)/GwtSpc.o \ +$(OBJDIR)/gwf3sfr8.o \ $(OBJDIR)/OutputControl.o \ -$(OBJDIR)/gwf3tvbase8.o \ -$(OBJDIR)/TimeArraySeriesLink.o \ -$(OBJDIR)/gwf3csub8.o \ -$(OBJDIR)/gwt1oc1.o \ $(OBJDIR)/gwt1ic1.o \ -$(OBJDIR)/gwt1obs1.o \ -$(OBJDIR)/TimeArraySeriesManager.o \ -$(OBJDIR)/gwf3tvk8.o \ -$(OBJDIR)/gwf3obs8.o \ +$(OBJDIR)/gwf3maw8.o \ +$(OBJDIR)/gwf3lak8.o \ +$(OBJDIR)/gwt1mst1.o \ +$(OBJDIR)/GwtDspOptions.o \ $(OBJDIR)/gwf3npf8.o \ -$(OBJDIR)/GwtSpc.o \ +$(OBJDIR)/GwtAdvOptions.o \ $(OBJDIR)/gwf3tvs8.o \ -$(OBJDIR)/BoundaryPackage.o \ -$(OBJDIR)/gwf3oc8.o \ -$(OBJDIR)/gwf3lak8.o \ -$(OBJDIR)/gwf3riv8.o \ -$(OBJDIR)/gwf3drn8.o \ -$(OBJDIR)/gwf3ghb8.o \ +$(OBJDIR)/GwfStorageUtils.o \ +$(OBJDIR)/Mover.o \ +$(OBJDIR)/GwfMvrPeriodData.o \ +$(OBJDIR)/ims8misc.o \ +$(OBJDIR)/GwfBuyInputData.o \ +$(OBJDIR)/InterfaceMap.o \ +$(OBJDIR)/gwf3disu8.o \ +$(OBJDIR)/GridSorting.o \ +$(OBJDIR)/DisConnExchange.o \ +$(OBJDIR)/CsrUtils.o \ +$(OBJDIR)/MappedVariable.o \ +$(OBJDIR)/TransportModel.o \ +$(OBJDIR)/NameFile.o \ +$(OBJDIR)/gwt1uzt1.o \ +$(OBJDIR)/gwt1ssm1.o \ +$(OBJDIR)/gwt1src1.o \ +$(OBJDIR)/gwt1sft1.o \ +$(OBJDIR)/gwt1oc1.o \ +$(OBJDIR)/gwt1obs1.o \ +$(OBJDIR)/gwt1mwt1.o \ +$(OBJDIR)/gwt1mvt1.o \ +$(OBJDIR)/gwt1lkt1.o \ +$(OBJDIR)/gwt1ist1.o \ +$(OBJDIR)/gwt1dsp.o \ $(OBJDIR)/gwt1cnc1.o \ -$(OBJDIR)/gwf3uzf8.o \ +$(OBJDIR)/gwt1adv1.o \ +$(OBJDIR)/gwf3disv8.o \ +$(OBJDIR)/gwf3dis8.o \ $(OBJDIR)/gwf3api8.o \ $(OBJDIR)/gwf3wel8.o \ -$(OBJDIR)/gwf3sto8.o \ +$(OBJDIR)/gwf3riv8.o \ $(OBJDIR)/gwf3rch8.o \ -$(OBJDIR)/gwf3maw8.o \ -$(OBJDIR)/gwt1fmi1.o \ -$(OBJDIR)/gwf3evt8.o \ -$(OBJDIR)/gwt1src1.o \ -$(OBJDIR)/gwf3chd8.o \ -$(OBJDIR)/NumericalModel.o \ -$(OBJDIR)/gwf3sfr8.o \ -$(OBJDIR)/gwt1mvt1.o \ -$(OBJDIR)/gwt1mst1.o \ -$(OBJDIR)/gwt1adv1.o \ -$(OBJDIR)/gwt1ist1.o \ +$(OBJDIR)/gwf3sto8.o \ +$(OBJDIR)/gwf3oc8.o \ +$(OBJDIR)/gwf3obs8.o \ +$(OBJDIR)/gwf3mvr8.o \ +$(OBJDIR)/gwf3hfb8.o \ +$(OBJDIR)/gwf3csub8.o \ $(OBJDIR)/gwf3buy8.o \ -$(OBJDIR)/NumericalExchange.o \ -$(OBJDIR)/gwt1apt1.o \ -$(OBJDIR)/gwt1dsp.o \ -$(OBJDIR)/gwt1uzt1.o \ -$(OBJDIR)/gwt1ssm1.o \ $(OBJDIR)/GhostNode.o \ -$(OBJDIR)/gwt1lkt1.o \ -$(OBJDIR)/NumericalSolution.o \ -$(OBJDIR)/gwt1mwt1.o \ -$(OBJDIR)/DisConnExchange.o \ -$(OBJDIR)/gwt1sft1.o \ +$(OBJDIR)/gwf3ghb8.o \ +$(OBJDIR)/gwf3evt8.o \ +$(OBJDIR)/gwf3drn8.o \ +$(OBJDIR)/gwf3chd8.o \ +$(OBJDIR)/ims8reordering.o \ +$(OBJDIR)/GridConnection.o \ +$(OBJDIR)/DistributedData.o \ +$(OBJDIR)/gwt1.o \ $(OBJDIR)/gwf3.o \ +$(OBJDIR)/ims8base.o \ +$(OBJDIR)/SpatialModelConnection.o \ +$(OBJDIR)/GwtInterfaceModel.o \ +$(OBJDIR)/GwtGwtExchange.o \ +$(OBJDIR)/GwfInterfaceModel.o \ $(OBJDIR)/GwfGwfExchange.o \ -$(OBJDIR)/gwt1.o \ +$(OBJDIR)/BaseSolution.o \ +$(OBJDIR)/Timer.o \ +$(OBJDIR)/ims8linear.o \ +$(OBJDIR)/GwtGwtConnection.o \ +$(OBJDIR)/GwfGwfConnection.o \ +$(OBJDIR)/SolutionGroup.o \ +$(OBJDIR)/NumericalSolution.o \ $(OBJDIR)/GwfGwtExchange.o \ $(OBJDIR)/SimulationCreate.o \ +$(OBJDIR)/ConnectionBuilder.o \ +$(OBJDIR)/comarg.o \ $(OBJDIR)/mf6core.o \ -$(OBJDIR)/mf6.o +$(OBJDIR)/BaseGeometry.o \ +$(OBJDIR)/mf6.o \ +$(OBJDIR)/StringList.o \ +$(OBJDIR)/MemorySetHandler.o \ +$(OBJDIR)/ilut.o \ +$(OBJDIR)/sparsekit.o \ +$(OBJDIR)/rcm.o \ +$(OBJDIR)/blas1_d.o \ +$(OBJDIR)/RectangularGeometry.o \ +$(OBJDIR)/CircularGeometry.o # Define the objects that make up the program $(PROGRAM) : $(OBJECTS) diff --git a/src/Model/GroundWaterFlow/gwf3dis8idm.f90 b/src/Model/GroundWaterFlow/gwf3dis8idm.f90 index 69df46d4ee3..fd4bbd51b3a 100644 --- a/src/Model/GroundWaterFlow/gwf3dis8idm.f90 +++ b/src/Model/GroundWaterFlow/gwf3dis8idm.f90 @@ -175,7 +175,7 @@ module GwfDisInputModule 'TOP', & ! tag name 'TOP', & ! fortran variable 'DOUBLE2D', & ! type - 'NCOL, NROW', & ! shape + 'NCOL NROW', & ! shape .true., & ! required .false., & ! multi-record .false., & ! preserve case @@ -191,7 +191,7 @@ module GwfDisInputModule 'BOTM', & ! tag name 'BOTM', & ! fortran variable 'DOUBLE3D', & ! type - 'NCOL, NROW, NLAY', & ! shape + 'NCOL NROW NLAY', & ! shape .true., & ! required .false., & ! multi-record .false., & ! preserve case @@ -207,7 +207,7 @@ module GwfDisInputModule 'IDOMAIN', & ! tag name 'IDOMAIN', & ! fortran variable 'INTEGER3D', & ! type - 'NCOL, NROW, NLAY', & ! shape + 'NCOL NROW NLAY', & ! shape .false., & ! required .false., & ! multi-record .false., & ! preserve case diff --git a/src/Model/GroundWaterFlow/gwf3disv8.f90 b/src/Model/GroundWaterFlow/gwf3disv8.f90 index 72d3c164ceb..e1581abff11 100644 --- a/src/Model/GroundWaterFlow/gwf3disv8.f90 +++ b/src/Model/GroundWaterFlow/gwf3disv8.f90 @@ -26,9 +26,9 @@ module GwfDisvModule real(DP), dimension(:, :), pointer, contiguous :: cellxy => null() ! cell center stored as 2d array of x and y integer(I4B), dimension(:), pointer, contiguous :: iavert => null() ! cell vertex pointer ia array integer(I4B), dimension(:), pointer, contiguous :: javert => null() ! cell vertex pointer ja array - real(DP), dimension(:, :), pointer, contiguous :: top2d => null() ! top elevations for each cell at top of model (ncpl, 1) - real(DP), dimension(:, :, :), pointer, contiguous :: bot3d => null() ! bottom elevations for each cell (ncpl, 1, nlay) - integer(I4B), dimension(:, :, :), pointer, contiguous :: idomain => null() ! idomain (ncpl, 1, nlay) + real(DP), dimension(:), pointer, contiguous :: top1d => null() ! top elevations for each cell at top of model (ncpl) + real(DP), dimension(:, :), pointer, contiguous :: bot2d => null() ! bottom elevations for each cell (ncpl, nlay) + integer(I4B), dimension(:, :), pointer, contiguous :: idomain => null() ! idomain (ncpl, nlay) contains procedure :: dis_df => disv_df procedure :: dis_da => disv_da @@ -200,8 +200,8 @@ subroutine disv_da(this) call mem_deallocate(this%cellxy) call mem_deallocate(this%iavert) call mem_deallocate(this%javert) - call mem_deallocate(this%top2d) - call mem_deallocate(this%bot3d) + call mem_deallocate(this%top1d) + call mem_deallocate(this%bot2d) call mem_deallocate(this%idomain) ! ! -- Return @@ -338,10 +338,10 @@ subroutine source_dimensions(this) this%nodesuser = this%nlay * this%ncpl ! ! -- Allocate non-reduced vectors for disv - call mem_allocate(this%idomain, this%ncpl, 1, this%nlay, 'IDOMAIN', & + call mem_allocate(this%idomain, this%ncpl, this%nlay, 'IDOMAIN', & this%memoryPath) - call mem_allocate(this%top2d, this%ncpl, 1, 'TOP2D', this%memoryPath) - call mem_allocate(this%bot3d, this%ncpl, 1, this%nlay, 'BOT3D', & + call mem_allocate(this%top1d, this%ncpl, 'TOP1D', this%memoryPath) + call mem_allocate(this%bot2d, this%ncpl, this%nlay, 'BOT2D', & this%memoryPath) ! ! -- Allocate vertices array @@ -351,7 +351,7 @@ subroutine source_dimensions(this) ! -- initialize all cells to be active (idomain = 1) do k = 1, this%nlay do j = 1, this%ncpl - this%idomain(j, 1, k) = 1 + this%idomain(j, k) = 1 end do end do ! @@ -405,8 +405,8 @@ subroutine source_griddata(this) idmMemoryPath = create_mem_path(this%name_model, 'DISV', idm_context) ! ! -- update defaults with idm sourced values - call mem_set_value(this%top2d, 'TOP', idmMemoryPath, afound(1)) - call mem_set_value(this%bot3d, 'BOTM', idmMemoryPath, afound(2)) + call mem_set_value(this%top1d, 'TOP', idmMemoryPath, afound(1)) + call mem_set_value(this%bot2d, 'BOTM', idmMemoryPath, afound(2)) call mem_set_value(this%idomain, 'IDOMAIN', idmMemoryPath, afound(3)) ! ! -- log simulation values @@ -472,7 +472,7 @@ subroutine grid_finalize(this) this%nodes = 0 do k = 1, this%nlay do j = 1, this%ncpl - if (this%idomain(j, 1, k) > 0) this%nodes = this%nodes + 1 + if (this%idomain(j, k) > 0) this%nodes = this%nodes + 1 end do end do ! @@ -487,16 +487,16 @@ subroutine grid_finalize(this) ! -- Check cell thicknesses do k = 1, this%nlay do j = 1, this%ncpl - if (this%idomain(j, 1, k) == 0) cycle - if (this%idomain(j, 1, k) > 0) then + if (this%idomain(j, k) == 0) cycle + if (this%idomain(j, k) > 0) then if (k > 1) then - top = this%bot3d(j, 1, k - 1) + top = this%bot2d(j, k - 1) else - top = this%top2d(j, 1) + top = this%top1d(j) end if - dz = top - this%bot3d(j, 1, k) + dz = top - this%bot2d(j, k) if (dz <= DZERO) then - write (ermsg, fmt=fmtdz) k, j, top, this%bot3d(j, 1, k) + write (ermsg, fmt=fmtdz) k, j, top, this%bot2d(j, k) call store_error(ermsg) end if end if @@ -523,10 +523,10 @@ subroutine grid_finalize(this) noder = 1 do k = 1, this%nlay do j = 1, this%ncpl - if (this%idomain(j, 1, k) > 0) then + if (this%idomain(j, k) > 0) then this%nodereduced(node) = noder noder = noder + 1 - elseif (this%idomain(j, 1, k) < 0) then + elseif (this%idomain(j, k) < 0) then this%nodereduced(node) = -1 else this%nodereduced(node) = 0 @@ -542,7 +542,7 @@ subroutine grid_finalize(this) noder = 1 do k = 1, this%nlay do j = 1, this%ncpl - if (this%idomain(j, 1, k) > 0) then + if (this%idomain(j, k) > 0) then this%nodeuser(noder) = node noder = noder + 1 end if @@ -551,7 +551,7 @@ subroutine grid_finalize(this) end do end if ! - ! -- Move top2d and bot3d into top and bot + ! -- Move top1d and bot2d into top and bot ! and set x and y center coordinates node = 0 do k = 1, this%nlay @@ -561,12 +561,12 @@ subroutine grid_finalize(this) if (this%nodes < this%nodesuser) noder = this%nodereduced(node) if (noder <= 0) cycle if (k > 1) then - top = this%bot3d(j, 1, k - 1) + top = this%bot2d(j, k - 1) else - top = this%top2d(j, 1) + top = this%top1d(j) end if this%top(noder) = top - this%bot(noder) = this%bot3d(j, 1, k) + this%bot(noder) = this%bot2d(j, k) this%xc(noder) = this%cellxy(1, j) this%yc(noder) = this%cellxy(2, j) end do @@ -937,8 +937,8 @@ subroutine write_grb(this, icelltype) write (iunit) this%xorigin ! xorigin write (iunit) this%yorigin ! yorigin write (iunit) this%angrot ! angrot - write (iunit) this%top2d ! top - write (iunit) this%bot3d ! botm + write (iunit) this%top1d ! top + write (iunit) this%bot2d ! botm write (iunit) this%vertices ! vertices write (iunit) (this%cellxy(1, i), i=1, this%ncpl) ! cellx write (iunit) (this%cellxy(2, i), i=1, this%ncpl) ! celly diff --git a/src/Model/GroundWaterFlow/gwf3disv8idm.f90 b/src/Model/GroundWaterFlow/gwf3disv8idm.f90 index 31259ddddf6..81bb4db9635 100644 --- a/src/Model/GroundWaterFlow/gwf3disv8idm.f90 +++ b/src/Model/GroundWaterFlow/gwf3disv8idm.f90 @@ -142,8 +142,8 @@ module GwfDisvInputModule 'GRIDDATA', & ! block 'TOP', & ! tag name 'TOP', & ! fortran variable - 'DOUBLE2D', & ! type - 'NCPL, 1', & ! shape + 'DOUBLE1D', & ! type + 'NCPL', & ! shape .true., & ! required .false., & ! multi-record .false., & ! preserve case @@ -158,8 +158,8 @@ module GwfDisvInputModule 'GRIDDATA', & ! block 'BOTM', & ! tag name 'BOTM', & ! fortran variable - 'DOUBLE3D', & ! type - 'NCPL, 1, NLAY', & ! shape + 'DOUBLE2D', & ! type + 'NCPL NLAY', & ! shape .true., & ! required .false., & ! multi-record .false., & ! preserve case @@ -174,8 +174,8 @@ module GwfDisvInputModule 'GRIDDATA', & ! block 'IDOMAIN', & ! tag name 'IDOMAIN', & ! fortran variable - 'INTEGER3D', & ! type - 'NCPL, 1, NLAY', & ! shape + 'INTEGER2D', & ! type + 'NCPL NLAY', & ! shape .false., & ! required .false., & ! multi-record .false., & ! preserve case diff --git a/src/Utilities/ArrayRead/ArrayReaderBase.f90 b/src/Utilities/ArrayRead/ArrayReaderBase.f90 new file mode 100644 index 00000000000..e387178f087 --- /dev/null +++ b/src/Utilities/ArrayRead/ArrayReaderBase.f90 @@ -0,0 +1,186 @@ +module ArrayReaderBaseModule + + use KindModule, only: DP, I4B, LGP + use ConstantsModule, only: MAXCHARLEN + use BlockParserModule, only: BlockParserType + use SimVariablesModule, only: errmsg + use SimModule, only: store_error + use InputOutputModule, only: openfile + + implicit none + private + public :: ArrayReaderBaseType + + type ArrayReaderBaseType + + type(BlockParserType), pointer :: parser => null() + integer(I4B) :: iout = 0 + integer(I4B) :: input_unit = 0 + character(len=:), allocatable :: array_name + character(len=:), allocatable :: filename + integer(I4B) :: iprn = 0 + logical(LGP) :: isConstant = .false. + logical(LGP) :: isInternal = .false. + logical(LGP) :: isOpenClose = .false. + logical(LGP) :: isBinary = .false. + + contains + + procedure :: read_array + procedure :: reset_reader + procedure :: read_control_record + procedure :: set_constant ! must be overriden + procedure :: fill_constant ! must be overriden + procedure :: fill_internal + procedure :: fill_open_close + procedure :: read_ascii ! must be overriden + procedure :: read_binary ! must be overriden + procedure :: set_factor ! must be overriden + procedure :: apply_factor ! must be overriden + procedure :: open_file + + end type ArrayReaderBaseType + +contains + + subroutine read_array(this) + class(ArrayReaderBaseType) :: this + + ! read control record + call this%read_control_record() + + ! fill array + if (this%isConstant) then + call this%fill_constant() + else if (this%isInternal) then + call this%fill_internal() + else if (this%isOpenClose) then + call this%fill_open_close() + end if + + end subroutine read_array + + subroutine reset_reader(this) + class(ArrayReaderBaseType) :: this + this%iprn = 0 + this%isConstant = .false. + this%isInternal = .false. + this%isOpenClose = .false. + this%isBinary = .false. + end subroutine reset_reader + + subroutine read_control_record(this) + class(ArrayReaderBaseType) :: this + logical(LGP) :: endOfBlock + character(len=100) :: keyword + character(len=MAXCHARLEN) :: string + + ! read the array input style + call this%parser%GetNextLine(endOfBlock) + call this%parser%GetStringCaps(keyword) + + ! load array based on the different styles + select case (keyword) + case ('CONSTANT') + this%isConstant = .true. + call this%set_constant() + case ('INTERNAL') + this%isInternal = .true. + case ('OPEN/CLOSE') + this%isOpenClose = .true. + call this%parser%GetString(string) + this%filename = trim(string) + case default + write (errmsg, *) 'Error reading control record for '// & + trim(adjustl(this%array_name))//'. & + & Use CONSTANT, INTERNAL, or OPEN/CLOSE.' + call store_error(errmsg) + call this%parser%StoreErrorUnit() + end select + + ! if INTERNAL or OPEN/CLOSE then look for FACTOR and IPRN + if (this%isInternal .or. this%isOpenClose) then + do + call this%parser%GetStringCaps(keyword) + if (keyword == '') exit + select case (keyword) + case ('FACTOR') + call this%set_factor() + case ('IPRN') + this%iprn = this%parser%GetInteger() + case ('(BINARY)') + this%isBinary = .true. + end select + end do + end if + + end subroutine read_control_record + + subroutine set_constant(this) + class(ArrayReaderBaseType) :: this + errmsg = 'Programming error in ArrayReader' + call store_error(errmsg, terminate=.true.) + end subroutine set_constant + + subroutine fill_constant(this) + class(ArrayReaderBaseType) :: this + errmsg = 'Programming error in ArrayReader' + call store_error(errmsg, terminate=.true.) + end subroutine fill_constant + + subroutine fill_internal(this) + class(ArrayReaderBaseType) :: this + this%input_unit = this%parser%iuactive + call this%read_ascii() + call this%apply_factor() + end subroutine fill_internal + + subroutine fill_open_close(this) + class(ArrayReaderBaseType) :: this + this%input_unit = 0 + call this%open_file() + if (this%isBinary) then + call this%read_binary() + else + call this%read_ascii() + end if + close (this%input_unit) + call this%apply_factor() + end subroutine fill_open_close + + subroutine read_ascii(this) + class(ArrayReaderBaseType) :: this + errmsg = 'Programming error in ArrayReader' + call store_error(errmsg, terminate=.true.) + end subroutine read_ascii + + subroutine read_binary(this) + class(ArrayReaderBaseType) :: this + errmsg = 'Programming error in ArrayReader' + call store_error(errmsg, terminate=.true.) + end subroutine read_binary + + subroutine set_factor(this) + class(ArrayReaderBaseType) :: this + errmsg = 'Programming error in ArrayReader' + call store_error(errmsg, terminate=.true.) + end subroutine set_factor + + subroutine apply_factor(this) + class(ArrayReaderBaseType) :: this + errmsg = 'Programming error in ArrayReader' + call store_error(errmsg, terminate=.true.) + end subroutine apply_factor + + subroutine open_file(this) + use OpenSpecModule, only: FORM, ACCESS + class(ArrayReaderBaseType) :: this + if (this%isBinary) then + call openfile(this%input_unit, this%iout, this%filename, & + 'OPEN/CLOSE', fmtarg_opt=FORM, accarg_opt=ACCESS) + else + call openfile(this%input_unit, this%iout, this%filename, 'OPEN/CLOSE') + end if + end subroutine open_file + +end module ArrayReaderBaseModule diff --git a/src/Utilities/ArrayRead/Double1dReader.f90 b/src/Utilities/ArrayRead/Double1dReader.f90 new file mode 100644 index 00000000000..a34202a5d5a --- /dev/null +++ b/src/Utilities/ArrayRead/Double1dReader.f90 @@ -0,0 +1,116 @@ +module Double1dReaderModule + + use KindModule, only: DP, I4B, LGP + use ConstantsModule, only: DZERO, DONE + use BlockParserModule, only: BlockParserType + use SimVariablesModule, only: errmsg + use SimModule, only: store_error, store_error_unit + use ArrayReadersModule, only: read_binary_header + use ArrayReaderBaseModule, only: ArrayReaderBaseType + + implicit none + private + public :: read_dbl1d + + type, extends(ArrayReaderBaseType) :: Double1dReaderType + + real(DP) :: constant_array_value = DZERO + real(DP) :: factor = DONE + real(DP), dimension(:), contiguous, pointer :: dbl1d => null() + + contains + + procedure :: reset_reader + procedure :: set_constant ! must be overriden + procedure :: fill_constant ! must be overriden + procedure :: read_ascii ! must be overriden + procedure :: read_binary ! must be overriden + procedure :: set_factor ! must be overriden + procedure :: apply_factor ! must be overriden + + end type Double1dReaderType + +contains + + subroutine read_dbl1d(parser, dbl1d, aname) + ! -- dummy + type(BlockParserType), intent(in), target :: parser + real(DP), dimension(:), contiguous, target :: dbl1d + character(len=*), intent(in) :: aname + ! -- local + type(Double1dReaderType) :: this + + this%parser => parser + this%dbl1d => dbl1d + this%array_name = aname + + call this%read_array() + + end subroutine read_dbl1d + + subroutine reset_reader(this) + class(Double1dReaderType) :: this + call this%ArrayReaderBaseType%reset_reader() + this%constant_array_value = DZERO + this%factor = DONE + end subroutine reset_reader + + subroutine set_constant(this) + class(Double1dReaderType) :: this + this%constant_array_value = this%parser%GetDouble() + end subroutine set_constant + + subroutine fill_constant(this) + class(Double1dReaderType) :: this + integer(I4B) :: i + do i = 1, size(this%dbl1d) + this%dbl1d(i) = this%constant_array_value + end do + end subroutine fill_constant + + subroutine read_ascii(this) + class(Double1dReaderType) :: this + integer(I4B) :: i + integer(I4B) :: istat + read (this%input_unit, *, iostat=istat, iomsg=errmsg) & + (this%dbl1d(i), i=1, size(this%dbl1d)) + if (istat /= 0) then + errmsg = 'Error reading data for array '//trim(this%array_name)// & + '. '//trim(errmsg) + call store_error(errmsg) + call store_error_unit(this%input_unit) + end if + end subroutine read_ascii + + subroutine read_binary(this) + class(Double1dReaderType) :: this + integer(I4B) :: i + integer(I4B) :: nvals + integer(I4B) :: istat + call read_binary_header(this%input_unit, this%iout, this%array_name, nvals) + read (this%input_unit, iostat=istat, iomsg=errmsg) & + (this%dbl1d(i), i=1, size(this%dbl1d)) + if (istat /= 0) then + errmsg = 'Error reading data for array '//trim(this%array_name)// & + '. '//trim(errmsg) + call store_error(errmsg) + call store_error_unit(this%input_unit) + end if + end subroutine read_binary + + subroutine set_factor(this) + class(Double1dReaderType) :: this + this%factor = this%parser%GetDouble() + end subroutine set_factor + + subroutine apply_factor(this) + class(Double1dReaderType) :: this + integer(I4B) :: i + if (this%factor /= DZERO) then + do i = 1, size(this%dbl1d) + this%dbl1d(i) = this%dbl1d(i) * this%factor + end do + end if + end subroutine apply_factor + +end module Double1dReaderModule diff --git a/src/Utilities/ArrayRead/Double2dReader.f90 b/src/Utilities/ArrayRead/Double2dReader.f90 new file mode 100644 index 00000000000..0df3d7bce15 --- /dev/null +++ b/src/Utilities/ArrayRead/Double2dReader.f90 @@ -0,0 +1,123 @@ +module Double2dReaderModule + + use KindModule, only: DP, I4B, LGP + use ConstantsModule, only: DZERO, DONE + use BlockParserModule, only: BlockParserType + use SimVariablesModule, only: errmsg + use SimModule, only: store_error, store_error_unit + use ArrayReadersModule, only: read_binary_header + use ArrayReaderBaseModule, only: ArrayReaderBaseType + + implicit none + private + public :: read_dbl2d + + type, extends(ArrayReaderBaseType) :: Double2dReaderType + + real(DP) :: constant_array_value = DZERO + real(DP) :: factor = DONE + real(DP), dimension(:, :), contiguous, pointer :: dbl2d => null() + + contains + + procedure :: reset_reader + procedure :: set_constant ! must be overriden + procedure :: fill_constant ! must be overriden + procedure :: read_ascii ! must be overriden + procedure :: read_binary ! must be overriden + procedure :: set_factor ! must be overriden + procedure :: apply_factor ! must be overriden + + end type Double2dReaderType + +contains + + subroutine read_dbl2d(parser, dbl2d, aname) + ! -- dummy + type(BlockParserType), intent(in), target :: parser + real(DP), dimension(:, :), contiguous, target :: dbl2d + character(len=*), intent(in) :: aname + ! -- local + type(Double2dReaderType) :: this + + this%parser => parser + this%dbl2d => dbl2d + this%array_name = aname + + call this%read_array() + + end subroutine read_dbl2d + + subroutine reset_reader(this) + class(Double2dReaderType) :: this + call this%ArrayReaderBaseType%reset_reader() + this%constant_array_value = DZERO + this%factor = DONE + end subroutine reset_reader + + subroutine set_constant(this) + class(Double2dReaderType) :: this + this%constant_array_value = this%parser%GetDouble() + end subroutine set_constant + + subroutine fill_constant(this) + class(Double2dReaderType) :: this + integer(I4B) :: i, j + do i = 1, size(this%dbl2d, dim=2) + do j = 1, size(this%dbl2d, dim=1) + this%dbl2d(j, i) = this%constant_array_value + end do + end do + end subroutine fill_constant + + subroutine read_ascii(this) + class(Double2dReaderType) :: this + integer(I4B) :: i, j + integer(I4B) :: istat + do i = 1, size(this%dbl2d, dim=2) + read (this%input_unit, *, iostat=istat, iomsg=errmsg) & + (this%dbl2d(j, i), j=1, size(this%dbl2d, dim=1)) + end do + if (istat /= 0) then + errmsg = 'Error reading data for array '//trim(this%array_name)// & + '. '//trim(errmsg) + call store_error(errmsg) + call store_error_unit(this%input_unit) + end if + end subroutine read_ascii + + subroutine read_binary(this) + class(Double2dReaderType) :: this + integer(I4B) :: i, j + integer(I4B) :: nvals + integer(I4B) :: istat + call read_binary_header(this%input_unit, this%iout, this%array_name, nvals) + read (this%input_unit, iostat=istat, iomsg=errmsg) & + ((this%dbl2d(j, i), j=1, size(this%dbl2d, dim=1)), & + i=1, size(this%dbl2d, dim=2)) + if (istat /= 0) then + errmsg = 'Error reading data for array '//trim(this%array_name)// & + '. '//trim(errmsg) + call store_error(errmsg) + call store_error_unit(this%input_unit) + end if + end subroutine read_binary + + subroutine set_factor(this) + class(Double2dReaderType) :: this + this%factor = this%parser%GetDouble() + end subroutine set_factor + + subroutine apply_factor(this) + class(Double2dReaderType) :: this + integer(I4B) :: i, j + if (this%factor /= DZERO) then + do i = 1, size(this%dbl2d, dim=2) + do j = 1, size(this%dbl2d, dim=1) + this%dbl2d(j, i) = this%dbl2d(j, i) * this%factor + end do + end do + end if + end subroutine apply_factor + +end module Double2dReaderModule diff --git a/src/Utilities/ArrayRead/Integer1dReader.f90 b/src/Utilities/ArrayRead/Integer1dReader.f90 new file mode 100644 index 00000000000..60b279e9fee --- /dev/null +++ b/src/Utilities/ArrayRead/Integer1dReader.f90 @@ -0,0 +1,148 @@ +module Integer1dReaderModule + + use KindModule, only: DP, I4B, LGP + use BlockParserModule, only: BlockParserType + use SimVariablesModule, only: errmsg + use SimModule, only: store_error, store_error_unit + use ArrayReadersModule, only: read_binary_header + use ArrayReaderBaseModule, only: ArrayReaderBaseType + + implicit none + private + public :: read_int1d, read_int1d_layered + + type, extends(ArrayReaderBaseType) :: Integer1dReaderType + + integer(I4B) :: constant_array_value = 0 + integer(I4B) :: factor = 1 + integer(I4B), dimension(:), contiguous, pointer :: int1d => null() + + contains + + procedure :: reset_reader + procedure :: set_constant ! must be overriden + procedure :: fill_constant ! must be overriden + procedure :: read_ascii ! must be overriden + procedure :: read_binary ! must be overriden + procedure :: set_factor ! must be overriden + procedure :: apply_factor ! must be overriden + + end type Integer1dReaderType + +contains + + subroutine read_int1d(parser, int1d, aname) + ! -- dummy + type(BlockParserType), intent(in), target :: parser + integer(I4B), dimension(:), contiguous, target :: int1d + character(len=*), intent(in) :: aname + ! -- local + type(Integer1dReaderType) :: this + + this%parser => parser + this%int1d => int1d + this%array_name = aname + + call this%read_array() + + end subroutine read_int1d + + subroutine read_int1d_layered(parser, int1d, aname, nlay, layer_shape) + use Integer2dReaderModule, only: read_int2d + ! -- dummy + type(BlockParserType), intent(in), target :: parser + integer(I4B), dimension(:), contiguous, target :: int1d + character(len=*), intent(in) :: aname + integer(I4B), intent(in) :: nlay + integer(I4B), dimension(:), intent(in) :: layer_shape + ! -- local + integer(I4B) :: k + integer(I4B) :: ncpl, nrow, ncol + integer(I4B) :: index_start, index_stop + integer(I4B), dimension(:, :), contiguous, pointer :: int2d_ptr + + ncpl = product(layer_shape) + index_start = 1 + do k = 1, nlay + index_stop = index_start + ncpl - 1 + if (size(layer_shape) == 2) then + ncol = layer_shape(1) + nrow = layer_shape(2) + int2d_ptr(1:ncol, 1:nrow) => int1d(index_start:index_stop) + call read_int2d(parser, int2d_ptr, aname) + else + call read_int1d(parser, int1d(index_start:index_stop), aname) + end if + index_start = index_stop + 1 + end do + + end subroutine read_int1d_layered + + subroutine reset_reader(this) + class(Integer1dReaderType) :: this + call this%ArrayReaderBaseType%reset_reader() + this%constant_array_value = 0 + this%factor = 1 + end subroutine reset_reader + + subroutine set_constant(this) + class(Integer1dReaderType) :: this + this%constant_array_value = this%parser%GetInteger() + end subroutine set_constant + + subroutine fill_constant(this) + class(Integer1dReaderType) :: this + integer(I4B) :: i + do i = 1, size(this%int1d) + this%int1d(i) = this%constant_array_value + end do + end subroutine fill_constant + + subroutine read_ascii(this) + class(Integer1dReaderType) :: this + integer(I4B) :: i + integer(I4B) :: nvals + integer(I4B) :: istat + nvals = size(this%int1d) + read (this%input_unit, *, iostat=istat, iomsg=errmsg) & + (this%int1d(i), i=1, size(this%int1d)) + if (istat /= 0) then + errmsg = 'Error reading data for array '//trim(this%array_name)// & + '. '//trim(errmsg) + call store_error(errmsg) + call store_error_unit(this%input_unit) + end if + end subroutine read_ascii + + subroutine read_binary(this) + class(Integer1dReaderType) :: this + integer(I4B) :: i + integer(I4B) :: nvals + integer(I4B) :: istat + call read_binary_header(this%input_unit, this%iout, this%array_name, nvals) + read (this%input_unit, iostat=istat, iomsg=errmsg) & + (this%int1d(i), i=1, size(this%int1d)) + if (istat /= 0) then + errmsg = 'Error reading data for array '//trim(this%array_name)// & + '. '//trim(errmsg) + call store_error(errmsg) + call store_error_unit(this%input_unit) + end if + end subroutine read_binary + + subroutine set_factor(this) + class(Integer1dReaderType) :: this + this%factor = this%parser%GetInteger() + end subroutine set_factor + + subroutine apply_factor(this) + class(Integer1dReaderType) :: this + integer(I4B) :: i + if (this%factor /= 0) then + do i = 1, size(this%int1d) + this%int1d(i) = this%int1d(i) * this%factor + end do + end if + end subroutine apply_factor + +end module Integer1dReaderModule diff --git a/src/Utilities/ArrayRead/Integer2dReader.f90 b/src/Utilities/ArrayRead/Integer2dReader.f90 new file mode 100644 index 00000000000..2fe4364d359 --- /dev/null +++ b/src/Utilities/ArrayRead/Integer2dReader.f90 @@ -0,0 +1,123 @@ +module Integer2dReaderModule + + use KindModule, only: DP, I4B, LGP + use ConstantsModule, only: DZERO, DONE + use BlockParserModule, only: BlockParserType + use SimVariablesModule, only: errmsg + use SimModule, only: store_error, store_error_unit + use ArrayReadersModule, only: read_binary_header + use ArrayReaderBaseModule, only: ArrayReaderBaseType + + implicit none + private + public :: read_int2d + + type, extends(ArrayReaderBaseType) :: Integer2dReaderType + + integer(I4B) :: constant_array_value = DZERO + integer(I4B) :: factor = DONE + integer(I4B), dimension(:, :), contiguous, pointer :: int2d => null() + + contains + + procedure :: reset_reader + procedure :: set_constant ! must be overriden + procedure :: fill_constant ! must be overriden + procedure :: read_ascii ! must be overriden + procedure :: read_binary ! must be overriden + procedure :: set_factor ! must be overriden + procedure :: apply_factor ! must be overriden + + end type Integer2dReaderType + +contains + + subroutine read_int2d(parser, int2d, aname) + ! -- dummy + type(BlockParserType), intent(in), target :: parser + integer(I4B), dimension(:, :), contiguous, target :: int2d + character(len=*), intent(in) :: aname + ! -- local + type(Integer2dReaderType) :: this + + this%parser => parser + this%int2d => int2d + this%array_name = aname + + call this%read_array() + + end subroutine read_int2d + + subroutine reset_reader(this) + class(Integer2dReaderType) :: this + call this%ArrayReaderBaseType%reset_reader() + this%constant_array_value = 0 + this%factor = 1 + end subroutine reset_reader + + subroutine set_constant(this) + class(Integer2dReaderType) :: this + this%constant_array_value = this%parser%GetInteger() + end subroutine set_constant + + subroutine fill_constant(this) + class(Integer2dReaderType) :: this + integer(I4B) :: i, j + do i = 1, size(this%int2d, dim=2) + do j = 1, size(this%int2d, dim=1) + this%int2d(j, i) = this%constant_array_value + end do + end do + end subroutine fill_constant + + subroutine read_ascii(this) + class(Integer2dReaderType) :: this + integer(I4B) :: i, j + integer(I4B) :: istat + do i = 1, size(this%int2d, dim=2) + read (this%input_unit, *, iostat=istat, iomsg=errmsg) & + (this%int2d(j, i), j=1, size(this%int2d, dim=1)) + end do + if (istat /= 0) then + errmsg = 'Error reading data for array '//trim(this%array_name)// & + '. '//trim(errmsg) + call store_error(errmsg) + call store_error_unit(this%input_unit) + end if + end subroutine read_ascii + + subroutine read_binary(this) + class(Integer2dReaderType) :: this + integer(I4B) :: i, j + integer(I4B) :: nvals + integer(I4B) :: istat + call read_binary_header(this%input_unit, this%iout, this%array_name, nvals) + read (this%input_unit, iostat=istat, iomsg=errmsg) & + ((this%int2d(j, i), j=1, size(this%int2d, dim=1)), & + i=1, size(this%int2d, dim=2)) + if (istat /= 0) then + errmsg = 'Error reading data for array '//trim(this%array_name)// & + '. '//trim(errmsg) + call store_error(errmsg) + call store_error_unit(this%input_unit) + end if + end subroutine read_binary + + subroutine set_factor(this) + class(Integer2dReaderType) :: this + this%factor = this%parser%GetInteger() + end subroutine set_factor + + subroutine apply_factor(this) + class(Integer2dReaderType) :: this + integer(I4B) :: i, j + if (this%factor /= DZERO) then + do i = 1, size(this%int2d, dim=2) + do j = 1, size(this%int2d, dim=1) + this%int2d(j, i) = this%int2d(j, i) * this%factor + end do + end do + end if + end subroutine apply_factor + +end module Integer2dReaderModule diff --git a/src/Utilities/ArrayRead/LayeredArrayReader.f90 b/src/Utilities/ArrayRead/LayeredArrayReader.f90 new file mode 100644 index 00000000000..26170fd7939 --- /dev/null +++ b/src/Utilities/ArrayRead/LayeredArrayReader.f90 @@ -0,0 +1,162 @@ +module LayeredArrayReaderModule + + use KindModule, only: DP, I4B, LGP + use BlockParserModule, only: BlockParserType + use Double1dReaderModule, only: read_dbl1d + use Double2dReaderModule, only: read_dbl2d + use Integer1dReaderModule, only: read_int1d + use Integer2dReaderModule, only: read_int2d + + implicit none + public :: read_dbl1d_layered + public :: read_dbl2d_layered + public :: read_dbl3d_layered + public :: read_int1d_layered + public :: read_int2d_layered + public :: read_int3d_layered + +contains + + subroutine read_dbl1d_layered(parser, dbl1d, aname, nlay, layer_shape) + ! -- dummy + type(BlockParserType), intent(in), target :: parser + real(DP), dimension(:), contiguous, target :: dbl1d + character(len=*), intent(in) :: aname + integer(I4B), intent(in) :: nlay + integer(I4B), dimension(:), intent(in) :: layer_shape + ! -- local + integer(I4B) :: k + integer(I4B) :: ncpl, nrow, ncol + integer(I4B) :: index_start, index_stop + real(DP), dimension(:, :), contiguous, pointer :: dbl2d_ptr + + ncpl = product(layer_shape) + index_start = 1 + do k = 1, nlay + index_stop = index_start + ncpl - 1 + if (size(layer_shape) == 2) then + ncol = layer_shape(1) + nrow = layer_shape(2) + dbl2d_ptr(1:ncol, 1:nrow) => dbl1d(index_start:index_stop) + call read_dbl2d(parser, dbl2d_ptr, aname) + else + call read_dbl1d(parser, dbl1d(index_start:index_stop), aname) + end if + index_start = index_stop + 1 + end do + + end subroutine read_dbl1d_layered + + subroutine read_dbl2d_layered(parser, dbl2d, aname, nlay, layer_shape) + ! -- dummy + type(BlockParserType), intent(in), target :: parser + real(DP), dimension(:, :), contiguous, target :: dbl2d + character(len=*), intent(in) :: aname + integer(I4B), intent(in) :: nlay + integer(I4B), dimension(:), intent(in) :: layer_shape + ! -- local + integer(I4B) :: k + integer(I4B) :: ncpl + real(DP), dimension(:), contiguous, pointer :: dbl1d_ptr + + ncpl = layer_shape(1) + do k = 1, nlay + dbl1d_ptr(1:ncpl) => dbl2d(1:ncpl, k) + call read_dbl1d(parser, dbl1d_ptr, aname) + end do + + end subroutine read_dbl2d_layered + + subroutine read_dbl3d_layered(parser, dbl3d, aname, nlay, layer_shape) + ! -- dummy + type(BlockParserType), intent(in), target :: parser + real(DP), dimension(:, :, :), contiguous, target :: dbl3d + character(len=*), intent(in) :: aname + integer(I4B), intent(in) :: nlay + integer(I4B), dimension(:), intent(in) :: layer_shape + ! -- local + integer(I4B) :: k + integer(I4B) :: ncol, nrow + real(DP), dimension(:, :), contiguous, pointer :: dbl2d_ptr + + ncol = layer_shape(1) + nrow = layer_shape(2) + do k = 1, nlay + dbl2d_ptr(1:ncol, 1:nrow) => dbl3d(:, :, k:k) + call read_dbl2d(parser, dbl2d_ptr, aname) + end do + + end subroutine read_dbl3d_layered + + subroutine read_int1d_layered(parser, int1d, aname, nlay, layer_shape) + ! -- dummy + type(BlockParserType), intent(in), target :: parser + integer(I4B), dimension(:), contiguous, target :: int1d + character(len=*), intent(in) :: aname + integer(I4B), intent(in) :: nlay + integer(I4B), dimension(:), intent(in) :: layer_shape + ! -- local + integer(I4B) :: k + integer(I4B) :: ncpl, nrow, ncol + integer(I4B) :: index_start, index_stop + integer(I4B), dimension(:, :), contiguous, pointer :: int2d_ptr + + ncpl = product(layer_shape) + index_start = 1 + do k = 1, nlay + index_stop = index_start + ncpl - 1 + if (size(layer_shape) == 2) then + ncol = layer_shape(1) + nrow = layer_shape(2) + int2d_ptr(1:ncol, 1:nrow) => int1d(index_start:index_stop) + call read_int2d(parser, int2d_ptr, aname) + else + call read_int1d(parser, int1d(index_start:index_stop), aname) + end if + index_start = index_stop + 1 + end do + + end subroutine read_int1d_layered + + subroutine read_int2d_layered(parser, int2d, aname, nlay, layer_shape) + ! -- dummy + type(BlockParserType), intent(in), target :: parser + integer(I4B), dimension(:, :), contiguous, target :: int2d + character(len=*), intent(in) :: aname + integer(I4B), intent(in) :: nlay + integer(I4B), dimension(:), intent(in) :: layer_shape + ! -- local + integer(I4B) :: k + integer(I4B) :: ncpl + integer(I4B), dimension(:), contiguous, pointer :: int1d_ptr + + ncpl = layer_shape(1) + do k = 1, nlay + int1d_ptr(1:ncpl) => int2d(1:ncpl, k) + call read_int1d(parser, int1d_ptr, aname) + end do + + end subroutine read_int2d_layered + + subroutine read_int3d_layered(parser, int3d, aname, nlay, layer_shape) + ! -- dummy + type(BlockParserType), intent(in), target :: parser + integer(I4B), dimension(:, :, :), contiguous, target :: int3d + character(len=*), intent(in) :: aname + integer(I4B), intent(in) :: nlay + integer(I4B), dimension(:), intent(in) :: layer_shape + ! -- local + integer(I4B) :: k + integer(I4B) :: ncol, nrow + integer(I4B), dimension(:, :), contiguous, pointer :: int2d_ptr + + ncol = layer_shape(1) + nrow = layer_shape(2) + do k = 1, nlay + int2d_ptr(1:ncol, 1:nrow) => int3d(:, :, k:k) + call read_int2d(parser, int2d_ptr, aname) + end do + + end subroutine read_int3d_layered + +end module LayeredArrayReaderModule diff --git a/src/Utilities/ArrayReaders.f90 b/src/Utilities/ArrayReaders.f90 index 973e4d5087e..c3a98961285 100644 --- a/src/Utilities/ArrayReaders.f90 +++ b/src/Utilities/ArrayReaders.f90 @@ -14,12 +14,20 @@ module ArrayReadersModule private public :: ReadArray + public :: read_binary_header interface ReadArray - module procedure read_array_int1d, read_array_int2d, read_array_int3d, & - read_array_dbl1d, read_array_dbl2d, read_array_dbl3d, & - read_array_dbl1d_layered, read_array_int1d_layered, & - read_array_dbl3d_all, read_array_int3d_all + module procedure & + read_array_int1d, & + read_array_int2d, & + read_array_int3d, & + read_array_dbl1d, & + read_array_dbl2d, & + read_array_dbl3d, & + read_array_dbl1d_layered, & + read_array_int1d_layered, & + read_array_dbl3d_all, & + read_array_int3d_all end interface ReadArray ! Integer readers diff --git a/src/Utilities/Idm/LoadMf6FileType.f90 b/src/Utilities/Idm/LoadMf6FileType.f90 index dd7080482b0..ce6c6943427 100644 --- a/src/Utilities/Idm/LoadMf6FileType.f90 +++ b/src/Utilities/Idm/LoadMf6FileType.f90 @@ -12,7 +12,16 @@ module LoadMf6FileTypeModule use SimVariablesModule, only: errmsg use SimModule, only: store_error use BlockParserModule, only: BlockParserType - use ArrayReadersModule, only: ReadArray + use LayeredArrayReaderModule, only: read_dbl1d_layered, & + read_dbl2d_layered, & + read_dbl3d_layered, & + read_int1d_layered, & + read_int2d_layered, & + read_int3d_layered + use Double1dReaderModule, only: read_dbl1d + use Double2dReaderModule, only: read_dbl2d + use Integer1dReaderModule, only: read_int1d + use Integer2dReaderModule, only: read_int2d use InputOutputModule, only: parseline use InputDefinitionModule, only: InputParamDefinitionType use InputDefinitionSelectorModule, only: get_param_definition_type, & @@ -73,7 +82,7 @@ subroutine idm_load_from_blockparser(parser, filetype, & do iblock = 1, size(mf6_input%p_block_dfns) call parse_block(parser, mf6_input, iblock, mshape, iout) ! - ! -- set model shape if discretion dimensions have been read + ! -- set model shape if discretization dimensions have been read if (mf6_input%p_block_dfns(iblock)%blockname == 'DIMENSIONS' .and. & filetype(1:3) == 'DIS') then call set_model_shape(mf6_input%file_type, componentMemPath, & @@ -234,6 +243,8 @@ recursive subroutine parse_tag(parser, mf6_input, iblock, mshape, iout, & call load_integer_type(parser, idt, mf6_input%memoryPath, iout) case ('INTEGER1D') call load_integer1d_type(parser, idt, mf6_input%memoryPath, mshape, iout) + case ('INTEGER2D') + call load_integer2d_type(parser, idt, mf6_input%memoryPath, mshape, iout) case ('INTEGER3D') call load_integer3d_type(parser, idt, mf6_input%memoryPath, mshape, iout) case ('DOUBLE') @@ -378,72 +389,132 @@ subroutine load_integer1d_type(parser, idt, memoryPath, mshape, iout) integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape !< model shape integer(I4B), intent(in) :: iout !< unit number for output integer(I4B), dimension(:), pointer, contiguous :: int1d - integer(I4B), pointer :: nsize1 + !integer(I4B), pointer :: nsize1 + integer(I4B) :: nlay integer(I4B) :: nvals + integer(I4B), dimension(:), allocatable :: array_shape + integer(I4B), dimension(:), allocatable :: layer_shape + character(len=LINELENGTH) :: keyword + ! Check if it is a full grid sized array (NODES) if (idt%shape == 'NODES') then nvals = product(mshape) - call mem_allocate(int1d, nvals, idt%mf6varname, memoryPath) else - call mem_setptr(nsize1, idt%shape, memoryPath) - call mem_allocate(int1d, nsize1, idt%mf6varname, memoryPath) + call get_shape_from_string(idt%shape, array_shape, memoryPath) + nvals = array_shape(1) end if - call read_grid_array(parser, mshape, idt%tagname, idt%layered, intarray=int1d) + ! allocate memory for the array + call mem_allocate(int1d, nvals, idt%mf6varname, memoryPath) + + ! check to see if the user specified "LAYERED" input + keyword = '' + if (idt%layered) then + call parser%GetStringCaps(keyword) + end if + ! read the array from the input file + if (keyword == 'LAYERED' .and. idt%layered) then + call get_layered_shape(mshape, nlay, layer_shape) + call read_int1d_layered(parser, int1d, idt%mf6varname, nlay, layer_shape) + else + call read_int1d(parser, int1d, idt%mf6varname) + end if + + ! log information on the loaded array to the list file call idm_log_var(int1d, idt%mf6varname, memoryPath, iout) return end subroutine load_integer1d_type - !> @brief load type 3d integer + !> @brief load type 2d integer !< - subroutine load_integer3d_type(parser, idt, memoryPath, mshape, iout) + subroutine load_integer2d_type(parser, idt, memoryPath, mshape, iout) type(BlockParserType), intent(inout) :: parser !< block parser type(InputParamDefinitionType), intent(in) :: idt !< input data type object describing this record character(len=*), intent(in) :: memoryPath !< memorypath to put loaded information integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape !< model shape integer(I4B), intent(in) :: iout !< unit number for output - integer(I4B), dimension(:, :, :), pointer, contiguous :: int3d - integer(I4B) :: ndim - integer(I4B) :: nsize1, nsize2, nsize3 + integer(I4B), dimension(:, :), pointer, contiguous :: int2d + integer(I4B) :: nlay + integer(I4B) :: nsize1, nsize2 + integer(I4B), dimension(:), allocatable :: array_shape + integer(I4B), dimension(:), allocatable :: layer_shape character(len=LINELENGTH) :: keyword - ndim = size(mshape) + ! determine the array shape from the input data defintion (idt%shape), + ! which looks like "NCOL, NROW, NLAY" + call get_shape_from_string(idt%shape, array_shape, memoryPath) + nsize1 = array_shape(1) + nsize2 = array_shape(2) - ! set sizes - if (ndim == 2) then - nsize1 = mshape(2) ! NCPL - nsize2 = 1 - nsize3 = mshape(1) - elseif (ndim == 3) then - nsize1 = mshape(3) ! NCOL - nsize2 = mshape(2) ! NROW - nsize3 = mshape(1) ! NLAY + ! create a new 3d memory managed variable + call mem_allocate(int2d, nsize1, nsize2, idt%mf6varname, memoryPath) + + ! check to see if the user specified "LAYERED" input + keyword = '' + if (idt%layered) then + call parser%GetStringCaps(keyword) end if - ! allocate the array using the memory manager - call mem_allocate(int3d, nsize1, nsize2, nsize3, idt%mf6varname, memoryPath) + ! read the array from the input file + if (keyword == 'LAYERED' .and. idt%layered) then + call get_layered_shape(mshape, nlay, layer_shape) + call read_int2d_layered(parser, int2d, idt%mf6varname, nlay, layer_shape) + else + call read_int2d(parser, int2d, idt%mf6varname) + end if + + ! log information on the loaded array to the list file + call idm_log_var(int2d, idt%mf6varname, memoryPath, iout) + return + end subroutine load_integer2d_type - ! fill the array from the file - if (idt%blockname == 'GRIDDATA') then + !> @brief load type 3d integer + !< + subroutine load_integer3d_type(parser, idt, memoryPath, mshape, iout) + type(BlockParserType), intent(inout) :: parser !< block parser + type(InputParamDefinitionType), intent(in) :: idt !< input data type object describing this record + character(len=*), intent(in) :: memoryPath !< memorypath to put loaded information + integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape !< model shape + integer(I4B), intent(in) :: iout !< unit number for output + integer(I4B), dimension(:, :, :), pointer, contiguous :: int3d + integer(I4B) :: nlay + integer(I4B) :: nsize1, nsize2, nsize3 + integer(I4B), dimension(:), allocatable :: array_shape + integer(I4B), dimension(:), allocatable :: layer_shape + character(len=LINELENGTH) :: keyword + integer(I4B), dimension(:), pointer, contiguous :: int1d_ptr + + ! determine the array shape from the input data defintion (idt%shape), + ! which looks like "NCOL, NROW, NLAY" + call get_shape_from_string(idt%shape, array_shape, memoryPath) + nsize1 = array_shape(1) + nsize2 = array_shape(2) + nsize3 = array_shape(3) + + ! create a new 3d memory managed variable + call mem_allocate(int3d, nsize1, nsize2, nsize3, idt%mf6varname, & + memoryPath) + + ! check to see if the user specified "LAYERED" input + keyword = '' + if (idt%layered) then call parser%GetStringCaps(keyword) - if (keyword == 'LAYERED') then - ! read by layer - call ReadArray(parser%iuactive, int3d(:, :, :), & - idt%mf6varname, ndim, nsize1, nsize2, & - nsize3, iout, 1, nsize3) - else - ! read full 3d array - call ReadArray(parser%iuactive, int3d(:, :, :), idt%mf6varname, & - ndim, nsize1 * nsize2 * nsize3, iout) - end if + end if + + ! read the array from the input file + if (keyword == 'LAYERED' .and. idt%layered) then + call get_layered_shape(mshape, nlay, layer_shape) + call read_int3d_layered(parser, int3d, idt%mf6varname, nlay, & + layer_shape) else - ! read full 3d array - call ReadArray(parser%iuactive, int3d(:, :, :), idt%mf6varname, & - ndim, nsize1 * nsize2 * nsize3, iout) + int1d_ptr(1:nsize1 * nsize2 * nsize3) => int3d(:, :, :) + call read_int1d(parser, int1d_ptr, idt%mf6varname) end if + ! log information on the loaded array to the list file call idm_log_var(int3d, idt%mf6varname, memoryPath, iout) + return end subroutine load_integer3d_type @@ -470,18 +541,39 @@ subroutine load_double1d_type(parser, idt, memoryPath, mshape, iout) integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape !< model shape integer(I4B), intent(in) :: iout !< unit number for output real(DP), dimension(:), pointer, contiguous :: dbl1d - integer(I4B), pointer :: nsize1 + !integer(I4B), pointer :: nsize1 + integer(I4B) :: nlay integer(I4B) :: nvals + integer(I4B), dimension(:), allocatable :: array_shape + integer(I4B), dimension(:), allocatable :: layer_shape + character(len=LINELENGTH) :: keyword + ! Check if it is a full grid sized array (NODES) if (idt%shape == 'NODES') then nvals = product(mshape) - call mem_allocate(dbl1d, nvals, idt%mf6varname, memoryPath) else - call mem_setptr(nsize1, idt%shape, memoryPath) - call mem_allocate(dbl1d, nsize1, idt%mf6varname, memoryPath) + call get_shape_from_string(idt%shape, array_shape, memoryPath) + nvals = array_shape(1) + end if + + ! allocate memory for the array + call mem_allocate(dbl1d, nvals, idt%mf6varname, memoryPath) + + ! check to see if the user specified "LAYERED" input + keyword = '' + if (idt%layered) then + call parser%GetStringCaps(keyword) end if - call read_grid_array(parser, mshape, idt%tagname, idt%layered, dbl1d) + ! read the array from the input file + if (keyword == 'LAYERED' .and. idt%layered) then + call get_layered_shape(mshape, nlay, layer_shape) + call read_dbl1d_layered(parser, dbl1d, idt%mf6varname, nlay, layer_shape) + else + call read_dbl1d(parser, dbl1d, idt%mf6varname) + end if + + ! log information on the loaded array to the list file call idm_log_var(dbl1d, idt%mf6varname, memoryPath, iout) return end subroutine load_double1d_type @@ -495,27 +587,36 @@ subroutine load_double2d_type(parser, idt, memoryPath, mshape, iout) integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape !< model shape integer(I4B), intent(in) :: iout !< unit number for output real(DP), dimension(:, :), pointer, contiguous :: dbl2d - integer(I4B) :: ndim + integer(I4B) :: nlay integer(I4B) :: nsize1, nsize2 + integer(I4B), dimension(:), allocatable :: array_shape + integer(I4B), dimension(:), allocatable :: layer_shape + character(len=LINELENGTH) :: keyword - ndim = size(mshape) - - ! set sizes - if (ndim == 2) then - nsize1 = mshape(2) ! NCPL - nsize2 = 1 - elseif (ndim == 3) then - nsize1 = mshape(3) ! NCOL - nsize2 = mshape(2) ! NROW - end if + ! determine the array shape from the input data defintion (idt%shape), + ! which looks like "NCOL, NROW, NLAY" + call get_shape_from_string(idt%shape, array_shape, memoryPath) + nsize1 = array_shape(1) + nsize2 = array_shape(2) - ! allocate the array using the memory manager + ! create a new 3d memory managed variable call mem_allocate(dbl2d, nsize1, nsize2, idt%mf6varname, memoryPath) - ! fill the array from the file - call ReadArray(parser%iuactive, dbl2d, idt%mf6varname, & - ndim, nsize1, nsize2, iout, 0) + ! check to see if the user specified "LAYERED" input + keyword = '' + if (idt%layered) then + call parser%GetStringCaps(keyword) + end if + + ! read the array from the input file + if (keyword == 'LAYERED' .and. idt%layered) then + call get_layered_shape(mshape, nlay, layer_shape) + call read_dbl2d_layered(parser, dbl2d, idt%mf6varname, nlay, layer_shape) + else + call read_dbl2d(parser, dbl2d, idt%mf6varname) + end if + ! log information on the loaded array to the list file call idm_log_var(dbl2d, idt%mf6varname, memoryPath, iout) return end subroutine load_double2d_type @@ -529,46 +630,43 @@ subroutine load_double3d_type(parser, idt, memoryPath, mshape, iout) integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape !< model shape integer(I4B), intent(in) :: iout !< unit number for output real(DP), dimension(:, :, :), pointer, contiguous :: dbl3d - integer(I4B) :: ndim + integer(I4B) :: nlay integer(I4B) :: nsize1, nsize2, nsize3 + integer(I4B), dimension(:), allocatable :: array_shape + integer(I4B), dimension(:), allocatable :: layer_shape character(len=LINELENGTH) :: keyword - - ndim = size(mshape) - - ! set sizes - if (ndim == 2) then - nsize1 = mshape(2) ! NCPL - nsize2 = 1 - nsize3 = mshape(1) - elseif (ndim == 3) then - nsize1 = mshape(3) ! NCOL - nsize2 = mshape(2) ! NROW - nsize3 = mshape(1) ! NLAY + real(DP), dimension(:), pointer, contiguous :: dbl1d_ptr + + ! determine the array shape from the input data defintion (idt%shape), + ! which looks like "NCOL, NROW, NLAY" + call get_shape_from_string(idt%shape, array_shape, memoryPath) + nsize1 = array_shape(1) + nsize2 = array_shape(2) + nsize3 = array_shape(3) + + ! create a new 3d memory managed variable + call mem_allocate(dbl3d, nsize1, nsize2, nsize3, idt%mf6varname, & + memoryPath) + + ! check to see if the user specified "LAYERED" input + keyword = '' + if (idt%layered) then + call parser%GetStringCaps(keyword) end if - ! allocate the array using the memory manager - call mem_allocate(dbl3d, nsize1, nsize2, nsize3, idt%mf6varname, memoryPath) - - ! fill the array from the file - if (idt%blockname == 'GRIDDATA') then - call parser%GetStringCaps(keyword) - if (keyword == 'LAYERED') then - ! read by layer - call ReadArray(parser%iuactive, dbl3d(:, :, :), & - idt%mf6varname, ndim, nsize1, nsize2, & - nsize3, iout, 1, nsize3) - else - ! read full 3d array - call ReadArray(parser%iuactive, dbl3d(:, :, :), idt%mf6varname, & - ndim, nsize1 * nsize2 * nsize3, iout) - end if + ! read the array from the input file + if (keyword == 'LAYERED' .and. idt%layered) then + call get_layered_shape(mshape, nlay, layer_shape) + call read_dbl3d_layered(parser, dbl3d, idt%mf6varname, nlay, & + layer_shape) else - ! read full 3d array - call ReadArray(parser%iuactive, dbl3d(:, :, :), idt%mf6varname, & - ndim, nsize1 * nsize2 * nsize3, iout) + dbl1d_ptr(1:nsize1 * nsize2 * nsize3) => dbl3d(:, :, :) + call read_dbl1d(parser, dbl1d_ptr, idt%mf6varname) end if + ! log information on the loaded array to the list file call idm_log_var(dbl3d, idt%mf6varname, memoryPath, iout) + return end subroutine load_double3d_type @@ -611,94 +709,53 @@ subroutine set_model_shape(ftype, model_mempath, dis_mempath, model_shape) return end subroutine set_model_shape - !> @brief read an array that is the size of the model grid - !< - subroutine read_grid_array(parser, mshape, array_name, layered, dblarray, & - intarray) - type(BlockParserType), intent(inout) :: parser !< block parser - integer(I4B), dimension(:), intent(in) :: mshape !< model shape - character(len=*), intent(in) :: array_name - logical(LGP), intent(in) :: layered - real(DP), dimension(:), optional, intent(inout) :: dblarray - integer(I4B), dimension(:), optional, intent(inout) :: intarray - integer(I4B) :: nvals + subroutine get_layered_shape(mshape, nlay, layer_shape) + integer(I4B), dimension(:), intent(in) :: mshape + integer(I4B), intent(out) :: nlay + integer(I4B), dimension(:), allocatable, intent(out) :: layer_shape integer(I4B) :: ndim - integer(I4B) :: ndim1 - integer(I4B) :: ndim2 - integer(I4B) :: ndim3 - integer(I4B) :: k1 - integer(I4B) :: k2 - integer(I4B) :: iout !< unit number for output - character(len=LINELENGTH) :: keyword ndim = size(mshape) - if (present(dblarray)) then - nvals = size(dblarray) - end if - if (present(intarray)) then - nvals = size(intarray) - end if - iout = 0 - - ! disu - if (ndim == 1) then - ndim1 = mshape(1) ! nodesuser - ndim2 = 1 ! none - ndim3 = 1 ! none - k1 = 0 - k2 = 0 - - ! disv - else if (ndim == 2) then - ndim1 = mshape(1) ! nlay - ndim2 = 1 ! none - ndim3 = mshape(2) ! ncpl - k1 = 1 - k2 = ndim1 - - ! dis - else if (ndim == 3) then - ndim1 = mshape(1) ! nlay - ndim2 = mshape(2) ! nrow - ndim3 = mshape(3) ! ncol - k1 = 1 - k2 = ndim1 + nlay = 0 + + if (ndim == 1) then ! disu + nlay = 1 + allocate (layer_shape(1)) + layer_shape(1) = mshape(1) + else if (ndim == 2) then ! disv + nlay = mshape(1) + allocate (layer_shape(1)) + layer_shape(1) = mshape(2) + else if (ndim == 3) then ! disu + nlay = mshape(1) + allocate (layer_shape(2)) + layer_shape(1) = mshape(3) ! ncol + layer_shape(2) = mshape(2) ! nrow end if - call parser%GetStringCaps(keyword) - if (keyword == 'LAYERED' .and. layered) then + end subroutine get_layered_shape - ! float array - if (present(dblarray)) then - call ReadArray(parser%iuactive, dblarray, & - array_name, ndim, ndim3, ndim2, & - ndim1, nvals, iout, k1, k2) - end if - - ! integer array - if (present(intarray)) then - call ReadArray(parser%iuactive, intarray, & - array_name, ndim, ndim3, ndim2, & - ndim1, nvals, iout, k1, k2) - end if - - else - - ! float array - if (present(dblarray)) then - call ReadArray(parser%iuactive, dblarray, array_name, & - ndim, nvals, iout, 0) - end if - - ! integer array - if (present(intarray)) then - call ReadArray(parser%iuactive, intarray, array_name, & - ndim, nvals, iout, 0) - end if - - end if + subroutine get_shape_from_string(shape_string, array_shape, memoryPath) + character(len=*), intent(in) :: shape_string + integer(I4B), dimension(:), allocatable, intent(inout) :: array_shape + character(len=*), intent(in) :: memoryPath !< memorypath to put loaded information + integer(I4B) :: ndim + integer(I4B) :: i + integer(I4B), pointer :: int_ptr + character(len=16), dimension(:), allocatable :: array_shape_string + character(len=:), allocatable :: shape_string_copy + + ! parse the string into multiple words + shape_string_copy = trim(shape_string)//' ' + call ParseLine(shape_string_copy, ndim, array_shape_string) + allocate (array_shape(ndim)) + + ! find shape in memory manager and put into array_shape + do i = 1, ndim + call mem_setptr(int_ptr, array_shape_string(i), memoryPath) + array_shape(i) = int_ptr + end do - return - end subroutine read_grid_array + end subroutine get_shape_from_string end module LoadMf6FileTypeModule diff --git a/src/meson.build b/src/meson.build index f5df942650a..0532308f8cd 100644 --- a/src/meson.build +++ b/src/meson.build @@ -113,6 +113,12 @@ modflow_sources = files( 'Solution' / 'SolutionGroup.f90', 'Timing' / 'ats.f90', 'Timing' / 'tdis.f90', + 'Utilities' / 'ArrayRead' / 'ArrayReaderBase.f90', + 'Utilities' / 'ArrayRead' / 'Double1dReader.f90', + 'Utilities' / 'ArrayRead' / 'Double2dReader.f90', + 'Utilities' / 'ArrayRead' / 'Integer1dReader.f90', + 'Utilities' / 'ArrayRead' / 'Integer2dReader.f90', + 'Utilities' / 'ArrayRead' / 'LayeredArrayReader.f90', 'Utilities' / 'Idm' / 'IdmLogger.f90', 'Utilities' / 'Idm' / 'IdmMf6FileLoader.f90', 'Utilities' / 'Idm' / 'ModflowInput.f90', diff --git a/utils/idmloader/scripts/dfn2f90.py b/utils/idmloader/scripts/dfn2f90.py index eb063e89fff..1863439397b 100644 --- a/utils/idmloader/scripts/dfn2f90.py +++ b/utils/idmloader/scripts/dfn2f90.py @@ -243,6 +243,7 @@ def _set_blk_param_strs(self, blockname, component, subcomponent): shape = v["shape"] shape = shape.replace("(", "") shape = shape.replace(")", "") + shape = shape.replace(",", "") shape = shape.upper() shapelist = shape.strip().split() ndim = len(shapelist) diff --git a/utils/mf5to6/make/makefile b/utils/mf5to6/make/makefile index 9708dc02d99..161a83f072d 100644 --- a/utils/mf5to6/make/makefile +++ b/utils/mf5to6/make/makefile @@ -5,10 +5,10 @@ include ./makedefaults # Define the source file directories SOURCEDIR1=../src -SOURCEDIR2=../src/LGR -SOURCEDIR3=../src/MF2005 -SOURCEDIR4=../src/NWT -SOURCEDIR5=../src/Preproc +SOURCEDIR2=../src/NWT +SOURCEDIR3=../src/LGR +SOURCEDIR4=../src/Preproc +SOURCEDIR5=../src/MF2005 SOURCEDIR6=../../../src/Utilities/Memory SOURCEDIR7=../../../src/Utilities/TimeSeries SOURCEDIR8=../../../src/Utilities