diff --git a/Makefile b/Makefile
index 47a7c0b14..735f5acdd 100644
--- a/Makefile
+++ b/Makefile
@@ -16,11 +16,11 @@ gnu: # BUILDTARGET GNU Fortran, C, and C++ compilers
"CC_SERIAL = gcc" \
"CXX_SERIAL = g++" \
"FFLAGS_PROMOTION = -fdefault-real-8 -fdefault-double-8" \
- "FFLAGS_OPT = -std=f2008 -O3 -ffree-line-length-none -fconvert=big-endian -ffree-form" \
+ "FFLAGS_OPT = -std=f2008 -O3 -fallow-argument-mismatch -ffree-line-length-none -fconvert=big-endian -ffree-form" \
"CFLAGS_OPT = -O3" \
"CXXFLAGS_OPT = -O3" \
"LDFLAGS_OPT = -O3" \
- "FFLAGS_DEBUG = -std=f2008 -g -ffree-line-length-none -fconvert=big-endian -ffree-form -fcheck=all -fbacktrace -ffpe-trap=invalid,zero,overflow" \
+ "FFLAGS_DEBUG = -g -ffree-line-length-none -fallow-argument-mismatch -fconvert=big-endian -ffree-form -fcheck=all -fbacktrace -ffpe-trap=invalid,zero,overflow" \
"CFLAGS_DEBUG = -g" \
"CXXFLAGS_DEBUG = -g" \
"LDFLAGS_DEBUG = -g" \
@@ -154,7 +154,7 @@ nvhpc: # BUILDTARGET NVIDIA HPC SDK
"FFLAGS_DEBUG = -O0 -g -Mbounds -Mchkptr -byteswapio -Mfree -Ktrap=divz,fp,inv,ovf -traceback" \
"CFLAGS_DEBUG = -O0 -g -traceback" \
"CXXFLAGS_DEBUG = -O0 -g -traceback" \
- "LDFLAGS_DEBUG = -O0 -g -Mbounds -Ktrap=divz,fp,inv,ovf -traceback" \
+ "LDFLAGS_DEBUG = -O0 -g -Mbounds -Mchkptr -Ktrap=divz,fp,inv,ovf -traceback" \
"FFLAGS_OMP = -mp" \
"CFLAGS_OMP = -mp" \
"FFLAGS_ACC = -Mnofma -acc -gpu=cc70,cc80 -Minfo=accel" \
@@ -184,7 +184,7 @@ pgi: # BUILDTARGET PGI compiler suite
"FFLAGS_DEBUG = -O0 -g -Mbounds -Mchkptr -byteswapio -Mfree -Ktrap=divz,fp,inv,ovf -traceback" \
"CFLAGS_DEBUG = -O0 -g -traceback" \
"CXXFLAGS_DEBUG = -O0 -g -traceback" \
- "LDFLAGS_DEBUG = -O0 -g -Mbounds -Ktrap=divz,fp,inv,ovf -traceback" \
+ "LDFLAGS_DEBUG = -O0 -g -Mbounds -Mchkptr -Ktrap=divz,fp,inv,ovf -traceback" \
"FFLAGS_OMP = -mp" \
"CFLAGS_OMP = -mp" \
"FFLAGS_ACC = -Mnofma -acc -Minfo=accel" \
@@ -216,7 +216,7 @@ pgi-summit: # BUILDTARGET PGI compiler suite w/OpenACC options for ORNL Summit
"FFLAGS_DEBUG = -O0 -g -Mbounds -Mchkptr -byteswapio -Mfree -Ktrap=divz,fp,inv,ovf -traceback" \
"CFLAGS_DEBUG = -O0 -g -traceback" \
"CXXFLAGS_DEBUG = -O0 -g -traceback" \
- "LDFLAGS_DEBUG = -O0 -g -Mbounds -Ktrap=divz,fp,inv,ovf -traceback" \
+ "LDFLAGS_DEBUG = -O0 -g -Mbounds -Mchkptr -Ktrap=divz,fp,inv,ovf -traceback" \
"FFLAGS_OMP = -mp" \
"CFLAGS_OMP = -mp" \
"PICFLAG = -fpic" \
@@ -667,10 +667,10 @@ intel: # BUILDTARGET Intel oneAPI Fortran, C, and C++ compiler suite
"CFLAGS_OPT = -O3" \
"CXXFLAGS_OPT = -O3" \
"LDFLAGS_OPT = -O3" \
- "FFLAGS_DEBUG = -g -convert big_endian -free -check bounds,pointers,arg_temp_created,format,shape,contiguous -fpe0 -traceback" \
+ "FFLAGS_DEBUG = -g -convert big_endian -free -check all -fpe0 -traceback" \
"CFLAGS_DEBUG = -g -traceback" \
"CXXFLAGS_DEBUG = -g -traceback" \
- "LDFLAGS_DEBUG = -g -traceback" \
+ "LDFLAGS_DEBUG = -g -fpe0 -traceback" \
"FFLAGS_OMP = -qopenmp" \
"CFLAGS_OMP = -qopenmp" \
"PICFLAG = -fpic" \
@@ -681,6 +681,86 @@ intel: # BUILDTARGET Intel oneAPI Fortran, C, and C++ compiler suite
"OPENMP = $(OPENMP)" \
"CPPFLAGS = $(MODEL_FORMULATION) -D_MPI" )
+intel-xd2000:
+ ( $(MAKE) all \
+ "FC_PARALLEL = ftn" \
+ "CC_PARALLEL = cc" \
+ "CXX_PARALLEL = CC" \
+ "FC_SERIAL = ftn" \
+ "CC_SERIAL = cc" \
+ "CXX_SERIAL = CC" \
+ "FFLAGS_PROMOTION = -real-size 64" \
+ "FFLAGS_OPT = -O3 -traceback -convert big_endian -FR -march=core-avx2 -mtune=core-avx2" \
+ "CFLAGS_OPT = -O3 -traceback -std=gnu90" \
+ "CXXFLAGS_OPT = -O3 -traceback" \
+ "LDFLAGS_OPT = -O3 -traceback" \
+ "FFLAGS_DEBUG = -g -convert big_endian -FR -CU -CB -check all -fpe0 -traceback" \
+ "CFLAGS_DEBUG = -g -traceback" \
+ "CXXFLAGS_DEBUG = -g -traceback" \
+ "LDFLAGS_DEBUG = -g -fpe0 -traceback" \
+ "FFLAGS_OMP = -qopenmp" \
+ "CFLAGS_OMP = -qopenmp" \
+ "CORE = $(CORE)" \
+ "DEBUG = $(DEBUG)" \
+ "USE_PAPI = $(USE_PAPI)" \
+ "OPENMP = $(OPENMP)" \
+ "CPPFLAGS = $(MODEL_FORMULATION) -D_MPI -DUNDERSCORE" )
+
+intel2-xd2000:
+ ( $(MAKE) all \
+ "FC_PARALLEL = ftn" \
+ "CC_PARALLEL = cc" \
+ "CXX_PARALLEL = CC" \
+ "FC_SERIAL = ftn" \
+ "CC_SERIAL = cc" \
+ "CXX_SERIAL = CC" \
+ "FFLAGS_PROMOTION = -real-size 64" \
+ "FFLAGS_OPT = -O2 -traceback -convert big_endian -FR -march=core-avx2 -mtune=core-avx2" \
+ "CFLAGS_OPT = -O2 -traceback -std=gnu90" \
+ "CXXFLAGS_OPT = -O2 -traceback" \
+ "LDFLAGS_OPT = -O2 -traceback" \
+ "FFLAGS_DEBUG = -g -convert big_endian -FR -CU -CB -check all -fpe0 -traceback" \
+ "CFLAGS_DEBUG = -g -traceback" \
+ "CXXFLAGS_DEBUG = -g -traceback" \
+ "LDFLAGS_DEBUG = -g -fpe0 -traceback" \
+ "FFLAGS_OMP = -qopenmp" \
+ "CFLAGS_OMP = -qopenmp" \
+ "CORE = $(CORE)" \
+ "DEBUG = $(DEBUG)" \
+ "USE_PAPI = $(USE_PAPI)" \
+ "OPENMP = $(OPENMP)" \
+ "CPPFLAGS = $(MODEL_FORMULATION) -D_MPI -DUNDERSCORE" )
+
+gfortran-xd2000: # BUILDTARGET GNU Fortran, C, and C++ compilers
+ ( $(MAKE) all \
+ "FC_PARALLEL = ftn" \
+ "CC_PARALLEL = cc" \
+ "CXX_PARALLEL = CC" \
+ "FC_SERIAL = ftn" \
+ "CC_SERIAL = cc" \
+ "CXX_SERIAL = CC" \
+ "FFLAGS_PROMOTION = -fdefault-real-8 -fdefault-double-8" \
+ "FFLAGS_OPT = -O3 -ffree-line-length-none -fconvert=big-endian -ffree-form -fallow-argument-mismatch" \
+ "CFLAGS_OPT = -O3" \
+ "CXXFLAGS_OPT = -O3" \
+ "LDFLAGS_OPT = -O3" \
+ "FFLAGS_DEBUG = -g -ffree-line-length-none -fconvert=big-endian -ffree-form -fallow-argument-mismatch -fcheck=all -fbacktrace -ffpe-trap=invalid,zero,overflow" \
+ "CFLAGS_DEBUG = -g" \
+ "CXXFLAGS_DEBUG = -g" \
+ "LDFLAGS_DEBUG = -g" \
+ "FFLAGS_OMP = -fopenmp" \
+ "CFLAGS_OMP = -fopenmp" \
+ "FFLAGS_ACC =" \
+ "CFLAGS_ACC =" \
+ "PICFLAG = -fPIC" \
+ "BUILD_TARGET = $(@)" \
+ "CORE = $(CORE)" \
+ "DEBUG = $(DEBUG)" \
+ "USE_PAPI = $(USE_PAPI)" \
+ "OPENMP = $(OPENMP)" \
+ "OPENACC = $(OPENACC)" \
+ "CPPFLAGS = $(MODEL_FORMULATION) -D_MPI" )
+
CPPINCLUDES =
FCINCLUDES =
LIBS =
@@ -690,6 +770,9 @@ ifneq "$(PIO)" ""
# Regardless of PIO library version, look for a lib subdirectory of PIO path
# NB: PIO_LIB is used later, so we don't just set LIBS directly
#
+#ifneq ($(wildcard $(PIO)/lib64), )
+# PIO_LIB = $(PIO)/lib64
+#else
ifneq ($(wildcard $(PIO)/lib), )
PIO_LIB = $(PIO)/lib
else
@@ -737,10 +820,10 @@ else # Not using PIO, using SMIOL
endif
ifneq "$(NETCDF)" ""
-ifneq ($(wildcard $(NETCDF)/lib/libnetcdf.*), )
+ifneq ($(wildcard $(NETCDF)/lib), )
NETCDFLIBLOC = lib
endif
-ifneq ($(wildcard $(NETCDF)/lib64/libnetcdf.*), )
+ifneq ($(wildcard $(NETCDF)/lib64), )
NETCDFLIBLOC = lib64
endif
CPPINCLUDES += -I$(NETCDF)/include
@@ -761,10 +844,10 @@ endif
ifneq "$(PNETCDF)" ""
-ifneq ($(wildcard $(PNETCDF)/lib/libpnetcdf.*), )
+ifneq ($(wildcard $(PNETCDF)/lib), )
PNETCDFLIBLOC = lib
endif
-ifneq ($(wildcard $(PNETCDF)/lib64/libpnetcdf.*), )
+ifneq ($(wildcard $(PNETCDF)/lib64), )
PNETCDFLIBLOC = lib64
endif
CPPINCLUDES += -I$(PNETCDF)/include
@@ -958,6 +1041,28 @@ else
OPENACC_MESSAGE="MPAS was built without OpenACC accelerator support."
endif
+ifneq ($(wildcard .mpas_core_*), ) # CHECK FOR BUILT CORE
+
+ifneq ($(wildcard .mpas_core_$(CORE)), ) # CHECK FOR SAME CORE AS ATTEMPTED BUILD.
+ override AUTOCLEAN=false
+ CONTINUE=true
+else
+ LAST_CORE=`cat .mpas_core_*`
+
+ifeq "$(AUTOCLEAN)" "true" # CHECK FOR CLEAN PRIOR TO BUILD OF A NEW CORE.
+ CONTINUE=true
+ AUTOCLEAN_MESSAGE="Infrastructure was cleaned prior to building ."
+else
+ CONTINUE=false
+endif # END OF AUTOCLEAN CHECK
+
+endif # END OF CORE=LAST_CORE CHECK
+
+else
+
+ override AUTOCLEAN=false
+ CONTINUE=true
+endif # END IF BUILT CORE CHECK
ifneq ($(wildcard namelist.$(NAMELIST_SUFFIX)), ) # Check for generated namelist file.
NAMELIST_MESSAGE="A default namelist file (namelist.$(NAMELIST_SUFFIX).defaults) has been generated, but namelist.$(NAMELIST_SUFFIX) has not been modified."
@@ -1014,119 +1119,12 @@ report_builds:
@echo "CORE=$(CORE)"
endif
+ifeq "$(CONTINUE)" "true"
all: mpas_main
-
+else
+all: clean_core
endif
-#
-# The rebuild_check target determines whether the shared framework or $(CORE) were
-# previously compiled with incompatible options, and stops the build with an error
-# message if so.
-#
-rebuild_check:
- @#
- @# Write current build options to a file .build_opts.tmp, to later be
- @# compared with build options use for the shared framework or core.
- @# Only build options that affect compatibility are written, while options
- @# like $(RM), $(BUILD_TARGET), and $(CORE) are not.
- @#
- $(shell printf "FC=$(FC)\n$\
- CC=$(CC)\n$\
- CXX=$(CXX)\n$\
- SFC=$(SFC)\n$\
- SCC=$(SCC)\n$\
- CFLAGS=$(CFLAGS)\n$\
- CXXFLAGS=$(CXXFLAGS)\n$\
- FFLAGS=$(FFLAGS)\n$\
- LDFLAGS=$(LDFLAGS)\n$\
- CPPFLAGS=$(CPPFLAGS)\n$\
- LIBS=$(LIBS)\n$\
- CPPINCLUDES=$(CPPINCLUDES)\n$\
- OPENMP=$(OPENMP)\n$\
- OPENMP_OFFLOAD=$(OPENMP_OFFLOAD)\n$\
- OPENACC=$(OPENACC)\n$\
- TAU=$(TAU)\n$\
- PICFLAG=$(PICFLAG)\n$\
- TIMER_LIB=$(TIMER_LIB)\n$\
- GEN_F90=$(GEN_F90)\n" | sed 's/-DMPAS_EXE_NAME=[^[:space:]]*//' | sed 's/-DMPAS_NAMELIST_SUFFIX=[^[:space:]]*//' | sed 's/-DCORE_[^[:space:]]*//' | sed 's/-DMPAS_GIT_VERSION=[^[:space:]]*//' > .build_opts.tmp )
-
- @#
- @# PREV_BUILD is set to "OK" if the shared framework and core are either
- @# clean or were previously compiled with compatible options. Otherwise,
- @# PREV_BUILD is set to "shared framework" if the shared framework was
- @# built with incompatible options, or "$(CORE) core" if the core was
- @# built with incompatible options.
- @#
- $(eval PREV_BUILD := $(shell $\
- if [ -f ".build_opts.framework" ]; then $\
- cmp -s .build_opts.tmp .build_opts.framework; $\
- if [ $$? -eq 0 ]; then $\
- stat=0; $\
- else $\
- stat=1; $\
- x="shared framework"; $\
- if [ "$(AUTOCLEAN)" = "true" ]; then $\
- cp .build_opts.tmp .build_opts.framework; $\
- fi; $\
- fi $\
- else $\
- stat=0; $\
- cp .build_opts.tmp .build_opts.framework; $\
- fi; $\
- : ; $\
- : At this this point, stat is already set, and we should only ; $\
- : set it to 1 but never to 0, as that might mask an incompatibility ; $\
- : in the framework build. ; $\
- : ; $\
- if [ -f ".build_opts.$(CORE)" ]; then $\
- cmp -s .build_opts.tmp .build_opts.$(CORE); $\
- if [ $$? -ne 0 ]; then $\
- stat=1; $\
- if [ "$$x" = "" ]; then $\
- x="$(CORE) core"; $\
- else $\
- x="$$x and $(CORE) core"; $\
- fi; $\
- if [ "$(AUTOCLEAN)" = "true" ]; then $\
- cp .build_opts.tmp .build_opts.$(CORE); $\
- fi; $\
- fi; $\
- else $\
- if [ $$stat -eq 0 ]; then $\
- cp .build_opts.tmp .build_opts.$(CORE); $\
- fi; $\
- fi; $\
- rm -f .build_opts.tmp; $\
- if [ $$stat -eq 1 ]; then $\
- printf "$$x"; $\
- else $\
- printf "OK"; $\
- fi; $\
- ))
-
- $(if $(findstring and,$(PREV_BUILD)),$(eval VERB=were),$(eval VERB=was))
-ifeq "$(AUTOCLEAN)" "true"
- $(if $(findstring framework,$(PREV_BUILD)),$(eval AUTOCLEAN_DEPS+=clean_shared))
- $(if $(findstring core,$(PREV_BUILD)),$(eval AUTOCLEAN_DEPS+=clean_core))
- $(if $(findstring OK,$(PREV_BUILD)), $(eval override AUTOCLEAN=false), )
- $(eval AUTOCLEAN_MESSAGE=The $(PREV_BUILD) $(VERB) cleaned and re-compiled.)
-else
- $(if $(findstring OK,$(PREV_BUILD)), \
- , \
- $(info ************************************************************************) \
- $(info The $(PREV_BUILD) $(VERB) previously compiled with ) \
- $(info incompatible options. Please do one of the following:) \
- $(info ) \
- $(info - Clean the $(CORE) core, which will also cause the shared) \
- $(info framework to be cleaned; then compile the $(CORE) core.) \
- $(info ) \
- $(info or)\
- $(info ) \
- $(info - Add AUTOCLEAN=true to the build command to automatically clean) \
- $(info and re-compile the $(PREV_BUILD).) \
- $(info ) \
- $(info ************************************************************************) \
- $(error ))
endif
@@ -1266,7 +1264,7 @@ ifeq "$(OPENACC)" "true"
endif # OPENACC eq true
-pio_test: openmp_test openacc_test pnetcdf_test
+pio_test: openmp_test openacc_test
@#
@# PIO_VERS will be set to:
@# 0 if no working PIO library was detected (and .piotest.log will contain error messages)
@@ -1352,10 +1350,9 @@ mpi_f08_test:
$(info Checking for mpi_f08 support...)
$(eval MPAS_MPI_F08 := $(shell $\
printf "program main\n$\
- & use mpi_f08, only : MPI_Init, MPI_Comm, MPI_INTEGER, MPI_Datatype\n$\
+ & use mpi_f08, only : MPI_Init, MPI_Comm\n$\
& integer :: ierr\n$\
& type (MPI_Comm) :: comm\n$\
- & type (MPI_Datatype), parameter :: MPI_INTEGERKIND = MPI_INTEGER\n$\
& call MPI_Init(ierr)\n$\
end program main\n" | sed 's/&/ /' > mpi_f08.f90; $\
$\
@@ -1374,57 +1371,20 @@ mpi_f08_test:
$(if $(findstring 1,$(MPAS_MPI_F08)), $(eval MPI_F08_MESSAGE = "Using the mpi_f08 module."), )
$(if $(findstring 1,$(MPAS_MPI_F08)), $(info mpi_f08 module detected.))
-
-pnetcdf_test:
- @#
- @# Create test C programs that look for PNetCDF header file and some symbols in it
- @#
-ifneq "$(PNETCDF)" ""
- @echo "Checking for a working PnetCDF library..."
- @printf "#include \"pnetcdf.h\"\n\
- &#include \"mpi.h\"\n\
- &int main(){\n\
- & int err, ncid;\n\
- & err = ncmpi_create(MPI_COMM_WORLD, \"foo.nc\", NC_NOCLOBBER, MPI_INFO_NULL, &ncid);\n\
- & return 0;\n\
- &}\n" | sed 's/&/ /' > pnetcdf.c
- @( $(CC) pnetcdf.c $(CPPINCLUDES) $(CFLAGS) $(LDFLAGS) -L$(PNETCDF)/$(PNETCDFLIBLOC) -lpnetcdf -o pnetcdf.out > pnetcdf.log 2>&1; \
- if [ $$? -eq 0 ] ; then \
- echo "$(CC) can compile test PnetCDF C program."; \
- else \
- echo "*********************************************************"; \
- echo "ERROR: Test PnetCDF C program could not be compiled by $(CC)."; \
- echo "Please ensure you have a working PnetCDF library installed."; \
- echo ""; \
- echo "The following compilation command failed with errors:" ; \
- echo "$(CC) pnetcdf.c $(CPPINCLUDES) $(CFLAGS) $(LDFLAGS) -L$(PNETCDF)/$(PNETCDFLIBLOC) -lpnetcdf -o pnetcdf.out"; \
- echo ""; \
- echo "Test program pnetcdf.c and output pnetcdf.log have been left"; \
- echo "in the top-level MPAS directory for further debugging"; \
- echo "*********************************************************"; \
- rm -f pnetcdf.out; exit 1; \
- fi )
-
- @rm -f pnetcdf.c pnetcdf.out pnetcdf.log
-else
- @echo "*********************************************************"; \
- echo "ERROR: The PNETCDF environment variable isn't set."; \
- echo "Please set this variable to where PnetCDF is installed."; \
- echo "*********************************************************"; \
- exit 1
-endif
-
-
ifneq "$(PIO)" ""
-MAIN_DEPS = rebuild_check openmp_test openacc_test pnetcdf_test pio_test mpi_f08_test
+MAIN_DEPS = openmp_test openacc_test pio_test mpi_f08_test
override CPPFLAGS += "-DMPAS_PIO_SUPPORT"
else
-MAIN_DEPS = rebuild_check openmp_test openacc_test pnetcdf_test mpi_f08_test
+MAIN_DEPS = openmp_test openacc_test mpi_f08_test
IO_MESSAGE = "Using the SMIOL library."
override CPPFLAGS += "-DMPAS_SMIOL_SUPPORT"
endif
+
mpas_main: $(MAIN_DEPS)
+ifeq "$(AUTOCLEAN)" "true"
+ $(RM) .mpas_core_*
+endif
cd src; $(MAKE) FC="$(FC)" \
CC="$(CC)" \
CXX="$(CXX)" \
@@ -1443,11 +1403,11 @@ mpas_main: $(MAIN_DEPS)
FCINCLUDES="$(FCINCLUDES)" \
CORE="$(CORE)"\
AUTOCLEAN="$(AUTOCLEAN)" \
- AUTOCLEAN_DEPS="$(AUTOCLEAN_DEPS)" \
GEN_F90="$(GEN_F90)" \
NAMELIST_SUFFIX="$(NAMELIST_SUFFIX)" \
EXE_NAME="$(EXE_NAME)"
+ @echo "$(EXE_NAME)" > .mpas_core_$(CORE)
if [ -e src/$(EXE_NAME) ]; then mv src/$(EXE_NAME) .; fi
( cd src/core_$(CORE); $(MAKE) ROOT_DIR="$(PWD)" post_build )
@echo "*******************************************************************************"
@@ -1469,13 +1429,11 @@ endif
@echo $(IO_MESSAGE)
@echo "*******************************************************************************"
clean:
- cd src; $(MAKE) clean RM="$(RM)" CORE="$(CORE)" AUTOCLEAN="$(AUTOCLEAN)"
+ cd src; $(MAKE) clean RM="$(RM)" CORE="$(CORE)"
+ $(RM) .mpas_core_*
$(RM) $(EXE_NAME)
$(RM) namelist.$(NAMELIST_SUFFIX).defaults
$(RM) streams.$(NAMELIST_SUFFIX).defaults
- if [ -f .build_opts.framework ]; then $(RM) .build_opts.framework; fi
- if [ -f .build_opts.$(CORE) ]; then $(RM) .build_opts.$(CORE); fi
-
core_error:
@echo ""
@echo "*******************************************************************************"
@@ -1486,6 +1444,26 @@ core_error:
exit 1
error: errmsg
+clean_core:
+ @echo ""
+ @echo "*******************************************************************************"
+ @echo " The MPAS infrastructure is currently built for the $(LAST_CORE) core."
+ @echo " Before building the $(CORE) core, please do one of the following."
+ @echo ""
+ @echo ""
+ @echo " To remove the $(LAST_CORE)_model executable and clean the MPAS infrastructure, run:"
+ @echo " make clean CORE=$(LAST_CORE)"
+ @echo ""
+ @echo " To preserve all executables except $(CORE)_model and clean the MPAS infrastructure, run:"
+ @echo " make clean CORE=$(CORE)"
+ @echo ""
+ @echo " Alternatively, AUTOCLEAN=true can be appended to the make command to force a clean,"
+ @echo " build a new $(CORE)_model executable, and preserve all other executables."
+ @echo ""
+ @echo "*******************************************************************************"
+ @echo ""
+ exit 1
+
else # CORE IF
all: error
@@ -1513,7 +1491,7 @@ errmsg:
@echo " DEBUG=true - builds debug version. Default is optimized version."
@echo " USE_PAPI=true - builds version using PAPI for timers. Default is off."
@echo " TAU=true - builds version using TAU hooks for profiling. Default is off."
- @echo " AUTOCLEAN=true - Enables automatic cleaning and re-compilation of code as needed."
+ @echo " AUTOCLEAN=true - forces a clean of infrastructure prior to build new core."
@echo " GEN_F90=true - Generates intermediate .f90 files through CPP, and builds with them."
@echo " TIMER_LIB=opt - Selects the timer library interface to be used for profiling the model. Options are:"
@echo " TIMER_LIB=native - Uses native built-in timers in MPAS"
@@ -1530,4 +1508,3 @@ errmsg:
ifdef CORE
exit 1
endif
-
diff --git a/src/core_atmosphere/Registry.xml b/src/core_atmosphere/Registry.xml
index 1dae83257..4f0c274d0 100644
--- a/src/core_atmosphere/Registry.xml
+++ b/src/core_atmosphere/Registry.xml
@@ -508,8 +508,10 @@
#endif
#ifdef DO_PHYSICS
+
+
@@ -984,6 +986,7 @@
+
@@ -4462,7 +4465,10 @@
+ description="dominant soil texture category"/>
+
+
@@ -4470,6 +4476,9 @@
+
+
diff --git a/src/core_atmosphere/mpas_atm_core.F b/src/core_atmosphere/mpas_atm_core.F
index 7f9be1196..7acf760ed 100644
--- a/src/core_atmosphere/mpas_atm_core.F
+++ b/src/core_atmosphere/mpas_atm_core.F
@@ -409,6 +409,7 @@ subroutine atm_mpas_init_block(dminfo, stream_manager, block, mesh, dt)
integer :: nCells, nEdges, nVertices, nVertLevels
integer :: thread
character(len=StrKIND), pointer :: mminlu
+ character(len=StrKIND), pointer :: mminsc
integer, pointer :: nThreads
integer, dimension(:), pointer :: cellThreadStart, cellThreadEnd
@@ -562,6 +563,16 @@ subroutine atm_mpas_init_block(dminfo, stream_manager, block, mesh, dt)
write(mminlu,'(a)') 'USGS'
end if
+ ! Before calling physics_init, ensure that mminsc contains the name of the soil colour dataset
+ call mpas_pool_get_array(sfc_input, 'mminsc', mminsc)
+ if (len_trim(mminsc) == 0) then
+ call mpas_log_write('****************************************************************')
+ call mpas_log_write('No information on soil colour dataset is available.')
+ call mpas_log_write('Assume that we are using ''DEFAULT_RAD_NOAH''.')
+ call mpas_log_write('****************************************************************')
+ write(mminsc,'(a)') 'DEFAULT_RAD_NOAH'
+ end if
+
if (moist_physics) then
!initialization of some input variables in registry:
diff --git a/src/core_atmosphere/physics/mpas_atmphys_driver_lsm.F b/src/core_atmosphere/physics/mpas_atmphys_driver_lsm.F
index 0116dcf56..ba5faf286 100644
--- a/src/core_atmosphere/physics/mpas_atmphys_driver_lsm.F
+++ b/src/core_atmosphere/physics/mpas_atmphys_driver_lsm.F
@@ -139,6 +139,7 @@ subroutine allocate_lsm
if(.not.allocated(cpm_p) ) allocate(cpm_p(ims:ime,jms:jme) )
if(.not.allocated(cqs2_p) ) allocate(cqs2_p(ims:ime,jms:jme) )
if(.not.allocated(isltyp_p) ) allocate(isltyp_p(ims:ime,jms:jme) )
+ if(.not.allocated(isctyp_p) ) allocate(isctyp_p(ims:ime,jms:jme) )
if(.not.allocated(ivgtyp_p) ) allocate(ivgtyp_p(ims:ime,jms:jme) )
if(.not.allocated(glw_p) ) allocate(glw_p(ims:ime,jms:jme) )
if(.not.allocated(grdflx_p) ) allocate(grdflx_p(ims:ime,jms:jme) )
@@ -216,6 +217,7 @@ subroutine deallocate_lsm
if(allocated(gsw_p) ) deallocate(gsw_p )
if(allocated(hfx_p) ) deallocate(hfx_p )
if(allocated(isltyp_p) ) deallocate(isltyp_p )
+ if(allocated(isctyp_p) ) deallocate(isctyp_p )
if(allocated(ivgtyp_p) ) deallocate(ivgtyp_p )
if(allocated(lai_p) ) deallocate(lai_p )
if(allocated(lh_p) ) deallocate(lh_p )
@@ -279,7 +281,7 @@ subroutine lsm_from_MPAS(configs,mesh,diag_physics,sfc_input,its,ite)
character(len=StrKIND),pointer:: config_microp_scheme, &
config_convection_scheme
- integer,dimension(:),pointer:: isltyp,ivgtyp
+ integer,dimension(:),pointer:: isltyp,isctyp,ivgtyp
real(kind=RKIND),dimension(:),pointer :: acsnom,acsnow,canwat,chs,chs2,chklowq,cpm,cqs2,glw, &
grdflx,gsw,hfx,lai,lh,noahres,potevp,qfx,qgh,qsfc, &
@@ -336,6 +338,7 @@ subroutine lsm_from_MPAS(configs,mesh,diag_physics,sfc_input,its,ite)
call mpas_pool_get_array(diag_physics,'znt' ,znt )
call mpas_pool_get_array(sfc_input,'isltyp' ,isltyp )
+ call mpas_pool_get_array(sfc_input,'isctyp' ,isctyp )
call mpas_pool_get_array(sfc_input,'ivgtyp' ,ivgtyp )
call mpas_pool_get_array(sfc_input,'shdmin' ,shdmin )
call mpas_pool_get_array(sfc_input,'shdmax' ,shdmax )
@@ -409,6 +412,7 @@ subroutine lsm_from_MPAS(configs,mesh,diag_physics,sfc_input,its,ite)
znt_p(i,j) = znt(i)
isltyp_p(i,j) = isltyp(i)
+ isctyp_p(i,j) = isctyp(i)
ivgtyp_p(i,j) = ivgtyp(i)
shdmin_p(i,j) = shdmin(i)
shdmax_p(i,j) = shdmax(i)
@@ -490,7 +494,7 @@ subroutine lsm_to_MPAS(configs,mesh,diag_physics,sfc_input,its,ite)
!local pointers:
character(len=StrKIND),pointer:: config_microp_scheme
- integer,dimension(:),pointer:: isltyp,ivgtyp
+ integer,dimension(:),pointer:: isltyp,isctyp,ivgtyp
real(kind=RKIND),dimension(:),pointer :: acsnom,acsnow,canwat,chs,chs2,chklowq,cpm,cqs2,glw, &
grdflx,gsw,hfx,lai,lh,noahres,potevp,qfx,qgh,qsfc, &
@@ -544,6 +548,7 @@ subroutine lsm_to_MPAS(configs,mesh,diag_physics,sfc_input,its,ite)
call mpas_pool_get_array(diag_physics,'znt' ,znt )
call mpas_pool_get_array(sfc_input,'isltyp' ,isltyp )
+ call mpas_pool_get_array(sfc_input,'isctyp' ,isctyp )
call mpas_pool_get_array(sfc_input,'ivgtyp' ,ivgtyp )
call mpas_pool_get_array(sfc_input,'shdmin' ,shdmin )
call mpas_pool_get_array(sfc_input,'shdmax' ,shdmax )
diff --git a/src/core_atmosphere/physics/mpas_atmphys_driver_lsm_noahmp.F b/src/core_atmosphere/physics/mpas_atmphys_driver_lsm_noahmp.F
index 7b93e7cf6..b2ba23c85 100644
--- a/src/core_atmosphere/physics/mpas_atmphys_driver_lsm_noahmp.F
+++ b/src/core_atmosphere/physics/mpas_atmphys_driver_lsm_noahmp.F
@@ -58,7 +58,7 @@ subroutine lsm_noahmp_fromMPAS(configs,mesh,diag,diag_physics,diag_physics_noahm
integer:: i,its,ite
integer:: n,ns,nsoil,nsnow,nzsnow
- integer,dimension(:),pointer:: isltyp,ivgtyp
+ integer,dimension(:),pointer:: isltyp,isctyp,ivgtyp
real(kind=RKIND),dimension(:),pointer:: latCell,lonCell
real(kind=RKIND),dimension(:),pointer:: shdmax,shdmin,vegfra,tmn,xice,xland
diff --git a/src/core_atmosphere/physics/mpas_atmphys_lsm_noahinit.F b/src/core_atmosphere/physics/mpas_atmphys_lsm_noahinit.F
index 07481fd6f..bce73fe67 100644
--- a/src/core_atmosphere/physics/mpas_atmphys_lsm_noahinit.F
+++ b/src/core_atmosphere/physics/mpas_atmphys_lsm_noahinit.F
@@ -91,7 +91,7 @@ subroutine lsminit(dminfo,mesh,configs,diag_physics,sfc_input)
character(len=StrKIND),pointer:: mminlu,mminsl
integer,pointer:: nCells,nSoilLevels
- integer,dimension(:),pointer:: ivgtyp,isltyp
+ integer,dimension(:),pointer:: ivgtyp,isltyp,isctyp
real(kind=RKIND),dimension(:),pointer:: snoalb,snow,snowh
real(kind=RKIND),dimension(:,:),pointer:: tslb,smois,sh2o
@@ -118,6 +118,7 @@ subroutine lsminit(dminfo,mesh,configs,diag_physics,sfc_input)
call mpas_pool_get_dimension(mesh,'nSoilLevels',nSoilLevels)
call mpas_pool_get_array(sfc_input,'isltyp', isltyp)
+ call mpas_pool_get_array(sfc_input,'isctyp', isctyp)
call mpas_pool_get_array(sfc_input,'ivgtyp', ivgtyp)
call mpas_pool_get_array(sfc_input,'sh2o' , sh2o )
call mpas_pool_get_array(sfc_input,'smois' , smois )
@@ -149,6 +150,20 @@ subroutine lsminit(dminfo,mesh,configs,diag_physics,sfc_input)
call physics_error_fatal("module_sf_noahlsm.F: lsminit: out of range value "// &
"of ISLTYP. Is this field in the input?" )
+ ! Make sure all cells have reasonable soil colour indices.
+ errflag = 0
+ do iCell = 1, nCells
+ if(isctyp(iCell) < 1) then
+ errflag = 1
+ write(err_message,*) "module_sf_noahlsm.F: lsminit: out of range ISCTYP ", &
+ iCell,isctyp(iCell),isltyp(iCell),ivgtyp(iCell)
+ call physics_message(err_message)
+ endif
+ end do
+ if (errflag == 1) &
+ call physics_error_fatal("module_sf_noahlsm.F: lsminit: out of range value "// &
+ "of ISCTYP. Is this field in the input?" )
+
!initializes soil liquid water content SH2O:
do iCell = 1, nCells
bx = bb(isltyp(iCell))
diff --git a/src/core_atmosphere/physics/mpas_atmphys_lsm_noahmpinit.F b/src/core_atmosphere/physics/mpas_atmphys_lsm_noahmpinit.F
index 19ad287d0..05c6c8b65 100644
--- a/src/core_atmosphere/physics/mpas_atmphys_lsm_noahmpinit.F
+++ b/src/core_atmosphere/physics/mpas_atmphys_lsm_noahmpinit.F
@@ -46,6 +46,7 @@ subroutine init_lsm_noahmp(configs,mesh,clock,diag_physics,diag_physics_noahmp,o
!--- local variables and arrays:
character(len=StrKIND),pointer:: mminlu
+ character(len=StrKIND),pointer:: mminsc
integer:: ns
@@ -73,6 +74,10 @@ subroutine init_lsm_noahmp(configs,mesh,clock,diag_physics,diag_physics_noahmp,o
call mpas_pool_get_array(sfc_input,'mminlu',mminlu)
mpas_noahmp%llanduse = mminlu
+!--- read NoahmpTable.TBL:
+ call mpas_pool_get_array(sfc_input,'mminsc',mminsc)
+ mpas_noahmp%lsoilcol = mminsc
+
!call mpas_log_write(' ')
!call mpas_log_write('--- enter subroutine NoahmpReadTable:')
call NoahmpReadTable(mpas_noahmp)
@@ -234,7 +239,7 @@ subroutine noahmp_init(configs,mesh,clock,diag_physics,diag_physics_noahmp,outpu
logical,pointer:: urban_physics
integer,pointer:: nsoilcomps
- integer,dimension(:),pointer:: isltyp,ivgtyp
+ integer,dimension(:),pointer:: isltyp,isctyp,ivgtyp
integer,dimension(:),pointer:: isnowxy
integer,dimension(:),pointer:: irnumsi,irnummi,irnumfi
@@ -322,11 +327,14 @@ subroutine noahmp_init(configs,mesh,clock,diag_physics,diag_physics_noahmp,outpu
!--- initialization of time-invariant surface variables needed in subroutine NoahmpInitMain:
call mpas_pool_get_array(sfc_input,'dzs' ,dzs )
call mpas_pool_get_array(sfc_input,'isltyp',isltyp)
+ call mpas_pool_get_array(sfc_input,'isctyp',isctyp)
call mpas_pool_get_array(sfc_input,'ivgtyp',ivgtyp)
do i = its, ite
mpas_noahmp%isltyp(i) = isltyp(i)
+ mpas_noahmp%isctyp(i) = isctyp(i)
mpas_noahmp%ivgtyp(i) = ivgtyp(i)
+
enddo
do ns = 1, nsoil
mpas_noahmp%dzs(ns) = dzs(ns,its)
diff --git a/src/core_atmosphere/physics/mpas_atmphys_update_surface.F b/src/core_atmosphere/physics/mpas_atmphys_update_surface.F
index 6e2057d5c..4d6a59e2d 100644
--- a/src/core_atmosphere/physics/mpas_atmphys_update_surface.F
+++ b/src/core_atmosphere/physics/mpas_atmphys_update_surface.F
@@ -151,7 +151,7 @@ subroutine physics_update_sst(dminfo,config_frac_seaice,mesh,sfc_input,diag_phys
call mpas_pool_get_array(sfc_input,'isice' ,isice )
call mpas_pool_get_array(sfc_input,'iswater' ,iswater )
- call mpas_pool_get_array(sfc_input,'isltyp' ,isltyp )
+ call mpas_pool_get_array(sfc_input,'isltyp' ,isltyp )
call mpas_pool_get_array(sfc_input,'ivgtyp' ,ivgtyp )
call mpas_pool_get_array(sfc_input,'landmask' ,landmask )
call mpas_pool_get_array(sfc_input,'vegfra' ,vegfra )
diff --git a/src/core_atmosphere/physics/mpas_atmphys_vars.F b/src/core_atmosphere/physics/mpas_atmphys_vars.F
index c6ebb80e8..2e8a26aef 100644
--- a/src/core_atmosphere/physics/mpas_atmphys_vars.F
+++ b/src/core_atmosphere/physics/mpas_atmphys_vars.F
@@ -882,6 +882,7 @@ module mpas_atmphys_vars
integer,dimension(:,:),allocatable:: &
isltyp_p, &!dominant soil type category [-]
+ isctyp_p, &!dominant soil colour category [-]
ivgtyp_p !dominant vegetation category [-]
real(kind=RKIND),dimension(:),allocatable:: &
diff --git a/src/core_atmosphere/physics/physics_noahmp/RELEASE_NOTES.md b/src/core_atmosphere/physics/physics_noahmp/RELEASE_NOTES.md
index ae45c39f5..6427e6b8f 100644
--- a/src/core_atmosphere/physics/physics_noahmp/RELEASE_NOTES.md
+++ b/src/core_atmosphere/physics/physics_noahmp/RELEASE_NOTES.md
@@ -400,7 +400,7 @@
- but in the file listed in the namelist as: HRLDAS_SETUP_FILE = "
- The initialization fields are: SNOW,CANWAT,TSK,TSLB,SMOIS
- - This file also contains the static grid/domain information: XLAT,XLONG,TMN,HGT,SEAICE,MAPFAC_MX,MAPFAC_MY,SHDMAX,SHDMIN,XLAND,IVGTYP,ISLTYP,DZS,ZS
+ - This file also contains the static grid/domain information: XLAT,XLONG,TMN,HGT,SEAICE,MAPFAC_MX,MAPFAC_MY,SHDMAX,SHDMIN,XLAND,IVGTYP,ISLTYP,ISCTYP,DZS,ZS
- This file can also contains some optional fields: LAI
- NOTE: a WRF input file can be used as a HRLDAS_SETUP_FILE
diff --git a/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/ConfigVarInTransferMod.F90 b/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/ConfigVarInTransferMod.F90
index 2de35ed9c..c78fee27a 100644
--- a/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/ConfigVarInTransferMod.F90
+++ b/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/ConfigVarInTransferMod.F90
@@ -62,7 +62,7 @@ subroutine ConfigVarInTransfer(noahmp, NoahmpIO)
! config domain variable
noahmp%config%domain%SurfaceType = 1
noahmp%config%domain%NumSwRadBand = 2
- noahmp%config%domain%SoilColor = 4
+ noahmp%config%domain%SoilColor = NoahmpIO%ISCTYP(I)
noahmp%config%domain%NumCropGrowStage = 8
noahmp%config%domain%FlagSoilProcess = NoahmpIO%calculate_soil
noahmp%config%domain%NumSoilTimeStep = NoahmpIO%soil_update_steps
@@ -75,6 +75,7 @@ subroutine ConfigVarInTransfer(noahmp, NoahmpIO)
noahmp%config%domain%SoilTimeStep = NoahmpIO%DTBL * NoahmpIO%soil_update_steps
noahmp%config%domain%GridSize = NoahmpIO%DX
noahmp%config%domain%LandUseDataName = NoahmpIO%LLANDUSE
+ noahmp%config%domain%SoilColDataName = NoahmpIO%LSOILCOL
noahmp%config%domain%VegType = NoahmpIO%IVGTYP(I)
noahmp%config%domain%CropType = NoahmpIO%CROPCAT(I)
noahmp%config%domain%IndicatorIceSfc = NoahmpIO%ICE
diff --git a/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/NoahmpIOVarFinalizeMod.F90 b/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/NoahmpIOVarFinalizeMod.F90
index 66458d638..49a709f12 100644
--- a/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/NoahmpIOVarFinalizeMod.F90
+++ b/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/NoahmpIOVarFinalizeMod.F90
@@ -37,6 +37,7 @@ subroutine NoahmpIOVarFinalizeDefault(NoahmpIO)
if ( allocated (NoahmpIO%zsoil) ) deallocate ( NoahmpIO%zsoil ) ! depth to soil interfaces [m]
if ( allocated (NoahmpIO%ivgtyp) ) deallocate ( NoahmpIO%ivgtyp ) ! vegetation type
if ( allocated (NoahmpIO%isltyp) ) deallocate ( NoahmpIO%isltyp ) ! soil type
+ if ( allocated (NoahmpIO%isctyp) ) deallocate ( NoahmpIO%isctyp ) ! soil colour class
if ( allocated (NoahmpIO%vegfra) ) deallocate ( NoahmpIO%vegfra ) ! vegetation fraction []
if ( allocated (NoahmpIO%tmn) ) deallocate ( NoahmpIO%tmn ) ! deep soil temperature [K]
if ( allocated (NoahmpIO%xland) ) deallocate ( NoahmpIO%xland ) ! =2 ocean; =1 land/seaice
diff --git a/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/NoahmpIOVarInitMod.F90 b/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/NoahmpIOVarInitMod.F90
index 4f3c3f4f2..cfd91bb91 100644
--- a/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/NoahmpIOVarInitMod.F90
+++ b/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/NoahmpIOVarInitMod.F90
@@ -41,6 +41,7 @@ subroutine NoahmpIOVarInitDefault(NoahmpIO)
if ( .not. allocated (NoahmpIO%zsoil) ) allocate ( NoahmpIO%zsoil (1:nsoil ) ) ! depth to soil interfaces [m]
if ( .not. allocated (NoahmpIO%ivgtyp) ) allocate ( NoahmpIO%ivgtyp (its:ite ) ) ! vegetation type
if ( .not. allocated (NoahmpIO%isltyp) ) allocate ( NoahmpIO%isltyp (its:ite ) ) ! soil type
+ if ( .not. allocated (NoahmpIO%isctyp) ) allocate ( NoahmpIO%isctyp (its:ite ) ) ! soil colour class
if ( .not. allocated (NoahmpIO%vegfra) ) allocate ( NoahmpIO%vegfra (its:ite ) ) ! vegetation fraction []
if ( .not. allocated (NoahmpIO%tmn) ) allocate ( NoahmpIO%tmn (its:ite ) ) ! deep soil temperature [K]
if ( .not. allocated (NoahmpIO%xland) ) allocate ( NoahmpIO%xland (its:ite ) ) ! =2 ocean; =1 land/seaice
@@ -469,6 +470,7 @@ subroutine NoahmpIOVarInitDefault(NoahmpIO)
NoahmpIO%ice = undefined_int
NoahmpIO%ivgtyp = undefined_int
NoahmpIO%isltyp = undefined_int
+ NoahmpIO%isctyp = undefined_int
NoahmpIO%isnowxy = undefined_int
NoahmpIO%coszen = undefined_real
NoahmpIO%xlat = undefined_real
diff --git a/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/NoahmpIOVarType.F90 b/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/NoahmpIOVarType.F90
index 05a29d703..2abf99670 100644
--- a/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/NoahmpIOVarType.F90
+++ b/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/NoahmpIOVarType.F90
@@ -68,6 +68,7 @@ module NoahmpIOVarType
integer :: soil_update_steps ! number of model time steps to update soil process
integer, allocatable, dimension(:) :: ivgtyp ! vegetation type
integer, allocatable, dimension(:) :: isltyp ! soil type
+ integer, allocatable, dimension(:) :: isctyp ! soil colour class
real(kind=kind_noahmp), allocatable, dimension(:) :: coszen ! cosine zenith angle
real(kind=kind_noahmp), allocatable, dimension(:) :: xlat ! latitude [rad]
real(kind=kind_noahmp), allocatable, dimension(:,:) :: dz8w ! thickness of atmo layers [m]
@@ -516,6 +517,7 @@ module NoahmpIOVarType
CHARACTER(LEN=256) :: mminsl = 'STAS' ! soil classification
CHARACTER(LEN=256) :: llanduse ! (=USGS, using USGS landuse classification)
+ CHARACTER(LEN=256) :: lsoilcol ! (=Noah, using Noah's native soil colour classification)
!------------------------------------------------------------------------
! Timing:
diff --git a/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/NoahmpReadTableMod.F90 b/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/NoahmpReadTableMod.F90
index eb01ceb2f..9ae895780 100644
--- a/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/NoahmpReadTableMod.F90
+++ b/src/core_atmosphere/physics/physics_noahmp/drivers/mpas/NoahmpReadTableMod.F90
@@ -28,7 +28,7 @@ subroutine NoahmpReadTable(NoahmpIO)
!-------------------------------------------------------
integer, parameter :: MVT = 27 ! number of vegetation types
integer, parameter :: MBAND = 2 ! number of radiation bands
- integer, parameter :: MSC = 8 ! number of soil texture
+ integer, parameter :: MSC = 21 ! number of soil colour types
integer, parameter :: MAX_SOILTYP = 30 ! max number of soil types
integer, parameter :: NCROP = 5 ! number of crop types
integer, parameter :: NSTAGE = 8 ! number of crop growth stages
@@ -98,12 +98,18 @@ subroutine NoahmpReadTable(NoahmpIO)
CZIL_DATA
! radiation parameters
+ character(len=256) :: RAD_DATASET_DESCRIPTION
+ integer :: NSC
real(kind=kind_noahmp) :: BETADS, BETAIS, EICE
real(kind=kind_noahmp), dimension(MBAND) :: ALBICE, ALBLAK, OMEGAS
real(kind=kind_noahmp), dimension(2) :: EG
real(kind=kind_noahmp), dimension(MSC) :: ALBSAT_VIS, ALBSAT_NIR, ALBDRY_VIS, ALBDRY_NIR
+ namelist / noahmp_rad_categories / RAD_DATASET_DESCRIPTION, NSC
namelist / noahmp_rad_parameters / ALBSAT_VIS, ALBSAT_NIR, ALBDRY_VIS, ALBDRY_NIR, ALBICE, ALBLAK, OMEGAS, &
BETADS, BETAIS, EG, EICE
+ namelist / noahmp_clm_rad_categories / RAD_DATASET_DESCRIPTION, NSC
+ namelist / noahmp_clm_rad_parameters / ALBSAT_VIS, ALBSAT_NIR, ALBDRY_VIS, ALBDRY_NIR, ALBICE, ALBLAK, OMEGAS, &
+ BETADS, BETAIS, EG, EICE
! global parameters
real(kind=kind_noahmp) :: CO2, O2, TIMEAN, FSATMX, Z0SNO, SSI, SNOW_RET_FAC ,SNOW_EMIS, SWEMX, TAU0, &
@@ -847,21 +853,54 @@ subroutine NoahmpReadTable(NoahmpIO)
if (ierr /= 0) then
write(*,'("WARNING: Cannot find file NoahmpTable.TBL")')
endif
- read(15,noahmp_rad_parameters)
+
+
+
+ DATASET_IDENTIFIER = NoahmpIO%LSOILCOL
+
+ inquire( file='NoahmpTable.TBL', exist=file_named )
+ if ( file_named ) then
+ open(15, file="NoahmpTable.TBL", status='old', form='formatted', action='read', iostat=ierr)
+ else
+ open(15, status='old', form='formatted', action='read', iostat=ierr)
+ end if
+ if ( ierr /= 0 ) then
+ write(*,'("WARNING: Cannot find file NoahmpTable.TBL")')
+ endif
+
+ select case (trim(DATASET_IDENTIFIER))
+ case ("DEFAULT_RAD_NOAH")
+ read(15,noahmp_rad_categories,iostat=ierr)
+ if ( ierr /= 0 ) then
+ write(*,'("WARNING: Namelist ''noahmp_rad_categories'' not found; using the default values.")')
+ write(*,'("WARNING: Update your NoahmpTable.TBL. This will eventually become an error.")')
+
+ RAD_DATASET_DESCRIPTION = "DEFAULT_RAD_NOAH" ! radiation (soil colour) type dataset
+ NSC = 8 ! total number of soil colour categories in Noah
+
+ end if
+ read(15,noahmp_rad_parameters)
+ case ("MODIFIED_RAD_CLM_NOAH")
+ read(15,noahmp_clm_rad_categories)
+ read(15,noahmp_clm_rad_parameters)
+ case default
+ write(*,'("WARNING: Unrecognised DATASET_IDENTIFIER in subroutine ReadNoahmpTable")')
+ write(*,'("WARNING: DATASET_IDENTIFIER = ''", A, "''")') trim(DATASET_IDENTIFIER)
+ end select
close(15)
! assign values
- NoahmpIO%ALBSAT_TABLE(:,1) = ALBSAT_VIS ! saturated soil albedos: 1=vis, 2=nir
- NoahmpIO%ALBSAT_TABLE(:,2) = ALBSAT_NIR ! saturated soil albedos: 1=vis, 2=nir
- NoahmpIO%ALBDRY_TABLE(:,1) = ALBDRY_VIS ! dry soil albedos: 1=vis, 2=nir
- NoahmpIO%ALBDRY_TABLE(:,2) = ALBDRY_NIR ! dry soil albedos: 1=vis, 2=nir
- NoahmpIO%ALBICE_TABLE = ALBICE
- NoahmpIO%ALBLAK_TABLE = ALBLAK
- NoahmpIO%OMEGAS_TABLE = OMEGAS
- NoahmpIO%BETADS_TABLE = BETADS
- NoahmpIO%BETAIS_TABLE = BETAIS
- NoahmpIO%EG_TABLE = EG
- NoahmpIO%EICE_TABLE = EICE
+ NoahmpIO%ALBSAT_TABLE(1:NSC,1) = ALBSAT_VIS ! saturated soil albedos: 1=vis, 2=nir
+ NoahmpIO%ALBSAT_TABLE(1:NSC,2) = ALBSAT_NIR ! saturated soil albedos: 1=vis, 2=nir
+ NoahmpIO%ALBDRY_TABLE(1:NSC,1) = ALBDRY_VIS ! dry soil albedos: 1=vis, 2=nir
+ NoahmpIO%ALBDRY_TABLE(1:NSC,2) = ALBDRY_NIR ! dry soil albedos: 1=vis, 2=nir
+ NoahmpIO%ALBICE_TABLE = ALBICE
+ NoahmpIO%ALBLAK_TABLE = ALBLAK
+ NoahmpIO%OMEGAS_TABLE = OMEGAS
+ NoahmpIO%BETADS_TABLE = BETADS
+ NoahmpIO%BETAIS_TABLE = BETAIS
+ NoahmpIO%EG_TABLE = EG
+ NoahmpIO%EICE_TABLE = EICE
!---------------- NoahmpTable.TBL global parameters
inquire( file='NoahmpTable.TBL', exist=file_named )
diff --git a/src/core_atmosphere/physics/physics_noahmp/parameters/NoahmpTable.TBL b/src/core_atmosphere/physics/physics_noahmp/parameters/NoahmpTable.TBL
index c9d37c5b4..e25c233ff 100644
--- a/src/core_atmosphere/physics/physics_noahmp/parameters/NoahmpTable.TBL
+++ b/src/core_atmosphere/physics/physics_noahmp/parameters/NoahmpTable.TBL
@@ -404,6 +404,11 @@
LAI_DEC = 4.0, 4.5, 0.0, 0.0, 2.0, 0.0, 0.0, 0.2, 0.3, 0.4, 0.2, 0.0, 0.0, 0.2, 0.0, 0.0, 0.0, 1.0, 0.6, 0.0,
/
+&noahmp_rad_categories
+ RAD_DATASET_DESCRIPTION = "DEFAULT_RAD_NOAH" ! radiation (soil colour) type dataset
+ NSC = 8 ! total number of soil colour categories in Noah
+/
+
&noahmp_rad_parameters
!------------------------------------------------------------------------------
! soil color: 1 2 3 4 5 6 7 8 soil color index for soil albedo
@@ -421,6 +426,34 @@
EICE = 0.98 ! emissivity ice surface
/
+&noahmp_clm_rad_categories
+ RAD_DATASET_DESCRIPTION = "MODIFIED_RAD_CLM_NOAH" ! radiation (soil colour) type dataset
+ NSC = 21 ! total number of soil colour categories in CLM
+/
+
+&noahmp_clm_rad_parameters
+ !------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
+ ! soil color index for soil albedo
+ ! soil color: 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21
+ !------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
+ ! saturated soil albedo at visible band
+ ALBSAT_VIS = 0.25, 0.23, 0.21, 0.20, 0.19, 0.18, 0.17, 0.16, 0.15, 0.14, 0.13, 0.12, 0.11, 0.10, 0.09, 0.08, 0.07, 0.06, 0.05, 0.04, 0.00
+ ! saturated soil albedo at NIR band
+ ALBSAT_NIR = 0.50, 0.46, 0.42, 0.40, 0.38, 0.36, 0.34, 0.32, 0.30, 0.28, 0.26, 0.24, 0.22, 0.20, 0.18, 0.16, 0.14, 0.12, 0.10, 0.08, 0.00
+ ! dry soil albedo at visible band
+ ALBDRY_VIS = 0.36, 0.34, 0.32, 0.31, 0.30, 0.29, 0.28, 0.27, 0.26, 0.25, 0.24, 0.23, 0.22, 0.20, 0.18, 0.16, 0.14, 0.12, 0.10, 0.08, 0.00
+ ! dry soil albedo at NIR band
+ ALBDRY_NIR = 0.61, 0.57, 0.53, 0.51, 0.49, 0.48, 0.45, 0.43, 0.41, 0.39, 0.37, 0.35, 0.33, 0.31, 0.29, 0.27, 0.25, 0.23, 0.21, 0.16, 0.00
+ !------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
+ ALBICE = 0.80, 0.55 ! albedo land ice: 1=vis, 2=nir
+ ALBLAK = 0.60, 0.40 ! albedo frozen lakes: 1=vis, 2=nir
+ OMEGAS = 0.8 , 0.4 ! two-stream parameter omega for snow
+ BETADS = 0.5 ! two-stream parameter betad for snow
+ BETAIS = 0.5 ! two-stream parameter betaI for snow
+ EG = 0.97, 0.98 ! emissivity soil surface 1-soil;2-lake
+ EICE = 0.98 ! emissivity ice surface
+/
+
&noahmp_global_parameters
! atmospheric constituants
CO2 = 395.0e-06 ! CO2 partial pressure
diff --git a/src/core_atmosphere/physics/physics_noahmp/src/ConfigVarInitMod.F90 b/src/core_atmosphere/physics/physics_noahmp/src/ConfigVarInitMod.F90
index 5c8af537b..2f9c5e9e2 100644
--- a/src/core_atmosphere/physics/physics_noahmp/src/ConfigVarInitMod.F90
+++ b/src/core_atmosphere/physics/physics_noahmp/src/ConfigVarInitMod.F90
@@ -49,6 +49,7 @@ subroutine ConfigVarInitDefault(noahmp)
! config domain variable
noahmp%config%domain%LandUseDataName = "MODIFIED_IGBP_MODIS_NOAH"
+ noahmp%config%domain%SoilColDataName = "DEFAULT_RAD_NOAH"
noahmp%config%domain%FlagUrban = .false.
noahmp%config%domain%FlagCropland = .false.
noahmp%config%domain%FlagDynamicCrop = .false.
diff --git a/src/core_atmosphere/physics/physics_noahmp/src/ConfigVarType.F90 b/src/core_atmosphere/physics/physics_noahmp/src/ConfigVarType.F90
index dc7979f3c..1ea5013d6 100644
--- a/src/core_atmosphere/physics/physics_noahmp/src/ConfigVarType.F90
+++ b/src/core_atmosphere/physics/physics_noahmp/src/ConfigVarType.F90
@@ -128,6 +128,7 @@ module ConfigVarType
type :: domain_type
character(len=256) :: LandUseDataName ! landuse dataset name (USGS or MODIFIED_IGBP_MODIS_NOAH)
+ character(len=256) :: SoilColDataName ! soil colour dataset name (DEFAULT_RAD_NOAH or MODIFIED_RAD_CLM_NOAH)
logical :: FlagUrban ! flag for urban grid
logical :: FlagCropland ! flag to identify croplands
logical :: FlagDynamicCrop ! flag to activate dynamic crop model
diff --git a/src/core_init_atmosphere/Registry.xml b/src/core_init_atmosphere/Registry.xml
index ceec72efd..f23d9fd7e 100644
--- a/src/core_init_atmosphere/Registry.xml
+++ b/src/core_init_atmosphere/Registry.xml
@@ -171,6 +171,11 @@
description="The soil category classification to use"
possible_values="`STATSGO' or `BNU'"/>
+
+
+
+
@@ -543,10 +550,12 @@
+
+
@@ -844,9 +853,15 @@
+
+
+
+
diff --git a/src/core_init_atmosphere/mpas_geotile_manager.F b/src/core_init_atmosphere/mpas_geotile_manager.F
index 04f9c60d0..89748f2c2 100644
--- a/src/core_init_atmosphere/mpas_geotile_manager.F
+++ b/src/core_init_atmosphere/mpas_geotile_manager.F
@@ -98,6 +98,7 @@ module mpas_geotile_manager
!> * isice = 24
!> * isurban = 1
!> * isoilwater = 14
+ !> * islcolwater = 9
!
!-----------------------------------------------------------------------
function mpas_geotile_mgr_init(mgr, path) result(ierr)
@@ -119,7 +120,7 @@ function mpas_geotile_mgr_init(mgr, path) result(ierr)
integer, pointer :: tile_z_start, tile_z_end
integer, pointer :: signed
integer, pointer :: tile_bdr
- integer, pointer :: iswater, islake, isice, isurban, isoilwater
+ integer, pointer :: iswater, islake, isice, isurban, isoilwater, islcolwater
integer, pointer :: category_min, category_max
integer :: err_level
real (kind=RKIND), pointer :: dx ! Grid spacing in the x-direction
@@ -188,7 +189,7 @@ function mpas_geotile_mgr_init(mgr, path) result(ierr)
!
! If this is a categorical field, then check to see if it has category_max and category_min,
- ! and then set the defaults of iswater, islake, isice, isurban and isoilwater
+ ! and then set the defaults of iswater, islake, isice, isurban, isoilwater and islcolwater
!
call mpas_pool_get_config(mgr % pool, 'type', fieldType)
if (fieldType == 'categorical') then
@@ -219,12 +220,14 @@ function mpas_geotile_mgr_init(mgr, path) result(ierr)
isice => null()
isurban => null()
isoilwater => null()
+ islcolwater => null()
call mpas_pool_get_config(mgr % pool, 'iswater', iswater)
call mpas_pool_get_config(mgr % pool, 'islake', islake)
call mpas_pool_get_config(mgr % pool, 'isice', isice)
call mpas_pool_get_config(mgr % pool, 'isurban', isurban)
call mpas_pool_get_config(mgr % pool, 'isoilwater', isoilwater)
+ call mpas_pool_get_config(mgr % pool, 'islcolwater', islcolwater)
if (.not. associated(iswater)) then
call mpas_pool_add_config(mgr % pool, 'iswater', 16)
@@ -245,6 +248,10 @@ function mpas_geotile_mgr_init(mgr, path) result(ierr)
if (.not. associated(isoilwater)) then
call mpas_pool_add_config(mgr % pool, 'isoilwater', 14)
endif
+
+ if (.not. associated(islcolwater)) then
+ call mpas_pool_add_config(mgr % pool, 'islcolwater', 1)
+ end if
endif
!
diff --git a/src/core_init_atmosphere/mpas_init_atm_cases.F b/src/core_init_atmosphere/mpas_init_atm_cases.F
index 673ebfc52..2f0331738 100644
--- a/src/core_init_atmosphere/mpas_init_atm_cases.F
+++ b/src/core_init_atmosphere/mpas_init_atm_cases.F
@@ -74,6 +74,7 @@ subroutine init_atm_setup_case(domain, stream_manager)
character (len=StrKIND), pointer :: config_specified_zeta_levels
character(len=StrKIND), pointer :: mminlu
+ character(len=StrKIND), pointer :: mminsc
character(len=StrKIND), pointer :: xtime
real (kind=RKIND) :: dt
real (kind=RKIND), pointer :: Time
@@ -257,6 +258,20 @@ subroutine init_atm_setup_case(domain, stream_manager)
write(mminlu,'(a)') 'USGS'
end if
+ !
+ ! If at this point the mminsc variable is blank, we assume that the static interp step was
+ ! not run, and that we are working with a static file created before there was a choice
+ ! of soil colour datasets; in this case, the dataset was almost necessarily DEFAULT_RAD_NOAH
+ !
+ call mpas_pool_get_array(mesh, 'mminsc', mminsc)
+ if (len_trim(mminsc) == 0) then
+ call mpas_log_write('****************************************************************')
+ call mpas_log_write('No information on soil colour dataset is available.')
+ call mpas_log_write('Assume that we are using ''DEFAULT_RAD_NOAH''.')
+ call mpas_log_write('****************************************************************')
+ write(mminsc,'(a)') 'DEFAULT_RAD_NOAH'
+ end if
+
call init_atm_case_gfs(block_ptr, mesh, nCells, nEdges, nVertLevels, fg, state, &
diag, diag_physics, block_ptr % dimensions, block_ptr % configs)
diff --git a/src/core_init_atmosphere/mpas_init_atm_static.F b/src/core_init_atmosphere/mpas_init_atm_static.F
index 0c540094d..993691e12 100644
--- a/src/core_init_atmosphere/mpas_init_atm_static.F
+++ b/src/core_init_atmosphere/mpas_init_atm_static.F
@@ -93,6 +93,7 @@ end subroutine interp_accumulation_function
real (kind=RKIND) :: soilcomp_msgval = 255.0_RKIND ! Modified later based on index file for soilcomp
integer, dimension(:), pointer :: lu_index
integer, dimension(:), pointer :: soilcat_top
+ integer, dimension(:), pointer :: soilcol_idx
integer, dimension(:), pointer :: nhs
integer, dimension(:,:), allocatable:: ncat
! Landmask is used by the accumulation function for maxsnoalb and soilcomp,
@@ -122,6 +123,7 @@ subroutine init_atm_static(mesh, dims, configs)
character(kind=c_char), dimension(StrKIND+1) :: c_fname
character(len=StrKIND), pointer :: config_geog_data_path
character(len=StrKIND), pointer :: config_landuse_data
+ character(len=StrKIND), pointer :: config_soilcol_data
character(len=StrKIND), pointer :: config_topo_data
character(len=StrKIND), pointer :: config_vegfrac_data
character(len=StrKIND), pointer :: config_albedo_data
@@ -155,6 +157,8 @@ subroutine init_atm_static(mesh, dims, configs)
integer, pointer :: isice_lu, iswater_lu
integer :: iswater_soil
+ integer :: iswater_slcol
+ integer :: isdefault_slcol
integer, pointer :: nCells, nCellsSolve, nEdges, nVertices, maxEdges
logical, pointer :: on_a_sphere
real (kind=RKIND), pointer :: sphere_radius
@@ -185,9 +189,11 @@ subroutine init_atm_static(mesh, dims, configs)
real (kind=RKIND), pointer :: missing_value
integer, dimension(:), pointer :: lu_index
integer, dimension(:), pointer :: soilcat_top
+ integer, dimension(:), pointer :: soilcol_idx
integer, dimension(:), pointer :: landmask
integer, dimension(:), pointer :: bdyMaskCell
character(len=StrKIND), pointer :: mminlu
+ character(len=StrKIND), pointer :: mminsc
real (kind=RKIND) :: xPixel, yPixel, zPixel
@@ -217,6 +223,7 @@ subroutine init_atm_static(mesh, dims, configs)
call mpas_pool_get_config(configs, 'config_geog_data_path', config_geog_data_path)
call mpas_pool_get_config(configs, 'config_landuse_data', config_landuse_data)
call mpas_pool_get_config(configs, 'config_soilcat_data', config_soilcat_data)
+ call mpas_pool_get_config(configs, 'config_soilcol_data', config_soilcol_data)
call mpas_pool_get_config(configs, 'config_topo_data', config_topo_data)
call mpas_pool_get_config(configs, 'config_vegfrac_data', config_vegfrac_data)
call mpas_pool_get_config(configs, 'config_albedo_data', config_albedo_data)
@@ -267,9 +274,11 @@ subroutine init_atm_static(mesh, dims, configs)
call mpas_pool_get_array(mesh, 'lu_index', lu_index)
call mpas_pool_get_array(mesh, 'mminlu', mminlu)
+ call mpas_pool_get_array(mesh, 'mminsc', mminsc)
call mpas_pool_get_array(mesh, 'isice_lu', isice_lu)
call mpas_pool_get_array(mesh, 'iswater_lu', iswater_lu)
call mpas_pool_get_array(mesh, 'soilcat_top', soilcat_top)
+ call mpas_pool_get_array(mesh, 'soilcol_idx', soilcol_idx)
call mpas_pool_get_array(mesh, 'landmask', landmask)
call mpas_pool_get_array(mesh, 'snoalb', snoalb)
call mpas_pool_get_array(mesh, 'greenfrac', greenfrac)
@@ -366,6 +375,23 @@ subroutine init_atm_static(mesh, dims, configs)
call mpas_log_write('Please correct the namelist.', messageType=MPAS_LOG_CRIT)
end select surface_input_select0
+!
+! Set soil colour parameters for soil albedo
+!
+ slcol_input_select0: select case(trim(config_soilcol_data))
+ case('DEFAULT_RAD_NOAH')
+ write(mminsc,'(a)') 'DEFAULT_RAD_NOAH'
+ case('MODIFIED_CLM_NOAH')
+ write(mminsc,'(a)') 'MODIFIED_CLM_NOAH'
+ case default
+ call mpas_log_write('*****************************************************************', messageType=MPAS_LOG_ERR)
+ call mpas_log_write('Invalid soil colour dataset '''//trim(config_soilcol_data) &
+ //''' selected for config_soilcol_data', messageType=MPAS_LOG_ERR)
+ call mpas_log_write(' Possible options are: ''DEFAULT_RAD_NOAH'', ''MODIFIED_CLM_NOAHM''', messageType=MPAS_LOG_ERR)
+ call mpas_log_write('*****************************************************************', messageType=MPAS_LOG_ERR)
+ call mpas_log_write('Please correct the namelist.', messageType=MPAS_LOG_CRIT)
+ end select slcol_input_select0
+
!
! Interpolate HGT
!
@@ -443,17 +469,55 @@ subroutine init_atm_static(mesh, dims, configs)
supersample_fac=supersample_fac_30s)
call mpas_log_write('--- end interpolate SOILCAT_TOP')
+!
+! Interpolate SOILCOL_IDX
+!
+ slcol_input_select1: select case(trim(config_soilcol_data))
+ case('DEFAULT_RAD_NOAH')
+ call mpas_log_write('Using 9-class default NOAH soil colour classes')
+ geog_sub_path = 'soilcolour_30s/'
+ ! Set default category and default water category for NOAH soil colour.
+ ! These will be overwritten if soil colour data sets exist.
+ isdefault_slcol = 4
+ iswater_slcol = 1
+
+ case('MODIFIED_CLM_NOAHM')
+ call mpas_log_write('Using 21-class MODIS 30-arc-second land cover dataset')
+ geog_sub_path = 'clm_soilcolour_21class_30s/'
+
+ ! Set default category and default water category for NOAH soil colour
+ ! These will be overwritten if soil colour data sets exist.
+ isdefault_slcol = 14
+ iswater_slcol = 21
+ case default
+ call mpas_log_write('*****************************************************************', messageType=MPAS_LOG_ERR)
+ call mpas_log_write('Invalid soil colour dataset '''//trim(config_soilcol_data) &
+ //''' selected for config_landuse_data', messageType=MPAS_LOG_ERR)
+ call mpas_log_write(' Possible options are: ''DEFAULT_RAD_NOAH'', ''MODIFIED_CLM_NOAHM''', messageType=MPAS_LOG_ERR)
+ call mpas_log_write('*****************************************************************', messageType=MPAS_LOG_ERR)
+ call mpas_log_write('Please correct the namelist.', messageType=MPAS_LOG_CRIT)
+ end select slcol_input_select1
+
+ ! Initialise soil colour with default index. Gridded data sets may or may not be
+ ! available.
+ soilcol_idx(1:iCell) = isdefault_slcol
+ call mpas_log_write('--- start interpolate SOILCOL_IDX')
+ call interp_soilcol(mesh, tree, trim(geog_data_path)//trim(geog_sub_path), iswater_slcol, &
+ supersample_fac=supersample_fac_30s)
+ call mpas_log_write('--- end interpolate SOILCOL_IDX')
+
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! KLUDGE TO FIX SOIL TYPE OVER ANTARCTICA
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
where (lu_index == isice_lu) soilcat_top = 16
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
-! CORRECT INCONSISTENT SOIL AND LAND USE DATA
+! CORRECT INCONSISTENT SOIL TYPE, SOIL COLOUR AND LAND USE DATA
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
do iCell = 1,nCells
- if (lu_index(iCell) == iswater_lu .or. &
- soilcat_top(iCell) == iswater_soil) then
+ if (lu_index(iCell) == iswater_lu .or. &
+ soilcat_top(iCell) == iswater_soil .or. &
+ soilcol_idx(iCell) == iswater_slcol ) then
if (lu_index(iCell) /= iswater_lu) then
call mpas_log_write('Turning lu_index into water at $i', intArgs=(/iCell/))
lu_index(iCell) = iswater_lu
@@ -462,6 +526,10 @@ subroutine init_atm_static(mesh, dims, configs)
call mpas_log_write('Turning soilcat_top into water at $i', intArgs=(/iCell/))
soilcat_top(iCell) = iswater_soil
end if
+ if (soilcol_idx(iCell) /= iswater_slcol) then
+ call mpas_log_write('Turning soilcol_idx into water at $i', intArgs=(/iCell/))
+ soilcol_idx(iCell) = iswater_slcol
+ end if
end if
end do
@@ -1936,6 +2004,89 @@ subroutine interp_soilcat(mesh, kdtree, geog_data_path, iswater_soil, supersampl
end subroutine interp_soilcat
+ !***********************************************************************
+ !
+ ! routine interp_soilcol
+ !
+ !> \brief Interpolate soil colour indices
+ !> \author Marcos Longo
+ !> \date 15 October 2025
+ !> \details
+ !> This sub-routine is based on interp_soilcat. The procedure interpolates soil
+ !> colour indices by using the init_atm_map_static_data routine and then by
+ !> accumulating the pixel values into each cell using
+ !> category_interp_accumulation.
+ !>
+ !> Variable mesh should be an mpas_pool that contains nCells and soilcol_idx;
+ !> kdtree should be an initialised mpas_kd_type tree with dimensions
+ !> (xCell, yCell, zCell); and geog_data_path should be the path to the
+ !> soil colour index data set. The values used by this data set to flag
+ !> pixels with no data (e.g., pixels fully covered with water or ice). Some
+ !> data sets may not have the no data flag, but we harmonise this data set
+ !> with the land use class mapping.
+ !-----------------------------------------------------------------------
+ subroutine interp_soilcol(mesh, kdtree, geog_data_path, iswater_slcol, supersample_fac)
+
+ implicit none
+
+ ! Input variables
+ type (mpas_pool_type), intent(inout) :: mesh
+ type (mpas_kd_type), pointer, intent(in) :: kdtree
+ character (len=*), intent(in) :: geog_data_path
+ integer, intent(inout) :: iswater_slcol
+ integer, intent(in), optional :: supersample_fac
+
+ ! Local variables
+ type (mpas_geotile_mgr_type) :: mgr
+ integer, pointer :: nCells
+ integer, pointer :: iswater_slcol_ptr
+
+ real (kind=RKIND), pointer :: scalefactor
+
+ integer :: iCell
+ integer :: ierr
+
+ ierr = mgr % init(trim(geog_data_path))
+ if (ierr /= 0) then
+ call mpas_log_write("Error occured initalising interpolation for "//trim(geog_data_path), messageType=MPAS_LOG_CRIT)
+ return
+ end if
+
+ call mpas_pool_get_dimension(mesh, 'nCells', nCells)
+ call mpas_pool_get_array(mesh, 'soilcol_idx', soilcol_idx)
+ call mpas_pool_get_config(mgr % pool, 'scale_factor', scalefactor )
+ call mpas_pool_get_config(mgr % pool, 'category_min', category_min )
+ call mpas_pool_get_config(mgr % pool, 'category_max', category_max )
+ call mpas_pool_get_config(mgr % pool, 'islcolwater' , iswater_slcol_ptr)
+
+ iswater_slcol = iswater_slcol_ptr
+
+ allocate(ncat(category_min:category_max, nCells))
+ ncat(:,:) = 0
+
+ call init_atm_map_static_data(mesh, mgr, kdtree, categorical_interp_criteria, categorical_interp_accumulation, &
+ supersample_fac=supersample_fac)
+
+ do iCell = 1, nCells
+ ! Variable ncat will tally the classes occurring in the grid cell. We select the class that has the
+ ! highest count as the value representing the grid cell. We use maxloc for this, however, maxloc returns
+ ! the index corresponding to the maximum value starting from 1. We shift the result to make sure the
+ ! category value correctly maps between category_min and category_max.
+ soilcol_idx(iCell) = maxloc(ncat(:,iCell), dim=1) - 1 + category_min
+ end do
+ deallocate(ncat)
+
+ ierr = mgr % finalize()
+ if (ierr /= 0) then
+ call mpas_log_write("Error occured finalising interpolation for "//trim(geog_data_path), messageType=MPAS_LOG_CRIT)
+ return
+ end if
+
+ nullify(category_min)
+ nullify(category_max)
+
+ end subroutine interp_soilcol
+
!***********************************************************************
!
! routine derive_landmask
diff --git a/src/core_init_atmosphere/mpas_parse_geoindex.F b/src/core_init_atmosphere/mpas_parse_geoindex.F
index 753ed4ee8..0be2c8098 100644
--- a/src/core_init_atmosphere/mpas_parse_geoindex.F
+++ b/src/core_init_atmosphere/mpas_parse_geoindex.F
@@ -131,7 +131,8 @@ function mpas_parse_index(path, geo_pool) result(ierr)
.or. trim(lhs) == 'description' &
.or. trim(lhs) == 'row_order' &
.or. trim(lhs) == 'endian' &
- .or. trim(lhs) == 'mminlu' ) then
+ .or. trim(lhs) == 'mminlu' &
+ .or. trim(lhs) == 'mminsc' ) then
char_t = rhs
call mpas_pool_add_config(geo_pool, trim(lhs), char_t)
@@ -171,6 +172,7 @@ function mpas_parse_index(path, geo_pool) result(ierr)
.or. trim(lhs) == 'isice' &
.or. trim(lhs) == 'isurban' &
.or. trim(lhs) == 'isoilwater' &
+ .or. trim(lhs) == 'islcolwater' &
.or. trim(lhs) == 'filename_digits' ) then
! Because each compiler handles reporting type errors when transferring
diff --git a/src/framework/mpas_io.F b/src/framework/mpas_io.F
index 09514a366..bf5a89804 100644
--- a/src/framework/mpas_io.F
+++ b/src/framework/mpas_io.F
@@ -868,7 +868,7 @@ subroutine MPAS_io_inq_var(handle, fieldname, fieldtype, ndims, dimnames, dimsiz
if (present(ierr)) ierr = MPAS_IO_ERR_BACKEND
deallocate(new_fieldlist_node % fieldhandle)
deallocate(new_fieldlist_node)
- call mpas_log_write('Variable ' // trim(fieldname) // ' not in input file.', MPAS_LOG_WARN)
+ call mpas_log_write('PIO_SUPPORT: Variable ' // trim(fieldname) // ' not in input file.', MPAS_LOG_WARN)
return
end if
!call mpas_log_write('Inquired about variable ID $i.', intArgs=(/new_fieldlist_node % fieldhandle % fieldid/) )
@@ -907,7 +907,7 @@ subroutine MPAS_io_inq_var(handle, fieldname, fieldtype, ndims, dimnames, dimsiz
if (present(ierr)) ierr = MPAS_IO_ERR_BACKEND
deallocate(new_fieldlist_node % fieldhandle)
deallocate(new_fieldlist_node)
- call mpas_log_write('Variable ' // trim(fieldname) // ' not in input file.', MPAS_LOG_WARN)
+ call mpas_log_write('SMIOL_SUPPORT_A: Variable ' // trim(fieldname) // ' not in input file.', MPAS_LOG_WARN)
return
end if
new_fieldlist_node % fieldhandle % field_type = smiol_type
@@ -955,7 +955,7 @@ subroutine MPAS_io_inq_var(handle, fieldname, fieldtype, ndims, dimnames, dimsiz
deallocate(new_fieldlist_node % fieldhandle)
deallocate(new_fieldlist_node)
- call mpas_log_write('Variable ' // trim(fieldname) // ' not in input file.', MPAS_LOG_WARN)
+ call mpas_log_write('SMIOL_SUPPORT_B: Variable ' // trim(fieldname) // ' not in input file.', MPAS_LOG_WARN)
return
end if
new_fieldlist_node% fieldhandle % ndims = smiol_ndims