Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

379 continue integrating the ground water hydrology work into cable #411

Open
wants to merge 77 commits into
base: main
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from 74 commits
Commits
Show all changes
77 commits
Select commit Hold shift + click to select a range
1b7e2e3
Done with science/soilsnow
rkutteh Sep 6, 2024
d1b5994
Done with science/canopy
rkutteh Sep 11, 2024
e9644f5
Done with science/radiation
rkutteh Sep 11, 2024
17650c5
Done with science/gw_hydro
rkutteh Sep 11, 2024
f73518c
Done with src/util
rkutteh Sep 16, 2024
10b7550
Done with src/offline
rkutteh Sep 30, 2024
ee4fb6b
Update cable_driver.F90
rkutteh Nov 5, 2024
8f354a0
Update cable_driver.F90
rkutteh Nov 5, 2024
5288555
Update src/offline/cable_initialise.F90
rkutteh Nov 5, 2024
c3b7559
Update src/offline/cable_iovars.F90
rkutteh Nov 5, 2024
17a3a16
Update cable_mpimaster.F90
rkutteh Nov 5, 2024
43af9ee
Update src/offline/cable_mpimaster.F90
rkutteh Nov 5, 2024
8d73320
Update src/offline/cable_mpimaster.F90
rkutteh Nov 5, 2024
d03d07f
Update src/offline/cable_mpiworker.F90
rkutteh Nov 5, 2024
8b98ec2
Update src/offline/cable_mpiworker.F90
rkutteh Nov 5, 2024
d9d19aa
Update src/offline/cable_mpiworker.F90
rkutteh Nov 5, 2024
b5f25af
Update cable_parameters.F90
rkutteh Nov 5, 2024
d6f92d7
Update src/offline/cable_parameters.F90
rkutteh Nov 5, 2024
68f930e
Update src/offline/cable_parameters.F90
rkutteh Nov 5, 2024
2aeeb40
Update cable_mpiworker.F90
rkutteh Nov 6, 2024
9955e81
Update cable_mpiworker.F90
rkutteh Nov 6, 2024
936acf8
Update src/offline/cable_output.F90
rkutteh Nov 6, 2024
230feec
Update cable_mpimaster.F90
rkutteh Nov 6, 2024
86f4bf4
Update cable_mpimaster.F90
rkutteh Nov 6, 2024
64dd7c3
Update cable_driver.F90
rkutteh Nov 11, 2024
9e9ca6f
Update cable_initialise.F90
rkutteh Nov 11, 2024
d9d7356
Update cable_initialise.F90
rkutteh Nov 11, 2024
a828802
Update cable_initialise.F90
rkutteh Nov 11, 2024
4b41852
Update cable_initialise.F90
rkutteh Nov 11, 2024
038d38f
Update src/offline/cable_output.F90
rkutteh Nov 11, 2024
c23ecf8
Update cable_driver.F90
rkutteh Nov 12, 2024
98b6673
Update cable_input.F90
rkutteh Nov 13, 2024
0143f33
Update cable_input.F90
rkutteh Nov 13, 2024
d3af980
Update cbl_stempv.F90
rkutteh Nov 13, 2024
1d974ae
Update cable_runtime_opts_mod.F90
rkutteh Nov 14, 2024
2b61c5e
Update src/science/soilsnow/cbl_surfbv.F90
rkutteh Nov 14, 2024
b2c4937
Update cable_checks.F90
rkutteh Nov 14, 2024
aa8997c
Update cable_output.F90
rkutteh Nov 15, 2024
e0ec557
Update cable_output.F90
rkutteh Nov 15, 2024
97d732e
Update cable_output.F90
rkutteh Nov 15, 2024
09d966a
Update cable_output.F90
rkutteh Nov 15, 2024
b0e1cd1
Update cable_output.F90
rkutteh Nov 15, 2024
f1676aa
Update cable_output.F90
rkutteh Nov 15, 2024
c5b30b1
Update cable_output.F90
rkutteh Nov 15, 2024
e4f4288
Update cable_output.F90
rkutteh Nov 15, 2024
274d8ff
Update src/offline/cable_parameters.F90
rkutteh Nov 15, 2024
ee3bb16
Update cable_output.F90
rkutteh Nov 15, 2024
30b3938
Update cable_parameters.F90
rkutteh Nov 24, 2024
fe52788
Update cable_parameters.F90
rkutteh Nov 24, 2024
49aaddc
Update cable_parameters.F90
rkutteh Nov 24, 2024
29917d2
Update src/offline/cable_parameters.F90
rkutteh Nov 24, 2024
f817a7e
Update cable_parameters.F90
rkutteh Nov 24, 2024
38ce564
Update cable_parameters.F90
rkutteh Nov 24, 2024
5679ede
Update cable_parameters.F90
rkutteh Nov 24, 2024
1c3e8a4
Update cable_parameters.F90
rkutteh Nov 24, 2024
03accf9
Update cable_parameters.F90
rkutteh Nov 24, 2024
f61c2a0
Update cable_parameters.F90
rkutteh Nov 24, 2024
6d955ed
Update cable_parameters.F90
rkutteh Nov 27, 2024
eb04b08
Update cable_parameters.F90
rkutteh Nov 27, 2024
f3a262c
Update cable_parameters.F90
rkutteh Nov 27, 2024
7bb3f32
Update cable_canopy.F90
rkutteh Dec 4, 2024
1fd7ce4
Update cable_canopy.F90
rkutteh Dec 4, 2024
878fe12
Update cable_psm.F90
rkutteh Dec 4, 2024
91a67ff
Update cable_psm.F90
rkutteh Dec 4, 2024
844fc2a
Update cbl_radiation.F90
rkutteh Dec 4, 2024
c51b50e
Update src/science/soilsnow/cbl_snowCheck.F90
rkutteh Dec 4, 2024
8014bf9
Update cable_canopy.F90
rkutteh Dec 4, 2024
30e55b6
Update cable_canopy.F90
rkutteh Dec 4, 2024
68ccad5
Update cable_canopy.F90
rkutteh Dec 4, 2024
0a79bcd
Update src/science/soilsnow/cbl_snowCheck.F90
rkutteh Dec 4, 2024
68e4658
Update src/science/soilsnow/cbl_snowCheck.F90
rkutteh Dec 4, 2024
4ca6929
Update cbl_snowCheck.F90
rkutteh Dec 4, 2024
4eb9a1a
Update cbl_snowCheck.F90
rkutteh Dec 4, 2024
34b1b77
Update cable_common.F90
rkutteh Dec 4, 2024
1fe7734
Add choice for 1D or 2D soil parameters in output
ccarouge Jan 17, 2025
d3e39e1
Fix for check_ranges() needs a rethink.
ccarouge Jan 17, 2025
a7a1c9a
Princeton forcing does not provide snowfall.
ccarouge Jan 17, 2025
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
54 changes: 50 additions & 4 deletions src/offline/cable_checks.F90
Original file line number Diff line number Diff line change
Expand Up @@ -84,7 +84,7 @@ MODULE cable_checks_module
ESoil = [-0.0015, 0.0015], &
TVeg = [-0.0003, 0.0003], &
ECanop = [-0.0003, 0.0003], &
PotEvap = [-0.005, 0.005], & !note should encompass Evap
PotEvap = [-0.005, 0.005], & !note should encompass Evap ! I have PotEvap = [-0.0006, 0.0006] if it makes a difference - rk4417
ACond = [0.0, 1.0], &
SoilWet = [-0.4, 1.2], &
Albedo = [0.0, 1.0], &
Expand Down Expand Up @@ -133,12 +133,13 @@ MODULE cable_checks_module
css = [700.0, 2200.0], &
rhosoil = [300.0, 3000.0], &
hyds = [5.0E-7, 8.5E-3], & ! vh_js ! sep14
! MMY 8.5E-3->8.5 since hyds uses m/s, but hyds_vec uses mm/s
rs20 = [0.0, 10.0], &
sand = [0.0, 1.0], &
sfc = [0.1, 0.5], &
silt = [0.0, 1.0], &
ssat = [0.35, 0.5], &
sucs = [-0.8, -0.03], &
sucs = [-0.8, -0.03], & ! MMY@23Apr2023 works for CABLE non-GW
swilt = [0.05, 0.4], &
froot = [0.0, 1.0], &
zse = [0.0, 5.0], &
Expand Down Expand Up @@ -204,6 +205,8 @@ MODULE cable_checks_module

TYPE(ranges_type), SAVE :: ranges

IF (cable_user%GW_MODEL) ranges%sucs = [30., 800.] ! MMY the range [-0.8, -0.03] doesn't suit Mark Decker's version
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Can't be done in the module definition part. This will require some more in-depth changes that will be done later.


INTERFACE check_range
MODULE PROCEDURE :: check_range_d1
MODULE PROCEDURE :: check_range_d2
Expand Down Expand Up @@ -486,6 +489,8 @@ SUBROUTINE mass_balance(dels,ktau, ssnow,soil,canopy,met,

! Local variables
REAL(r_2), DIMENSION(:,:,:),POINTER, SAVE :: bwb ! volumetric soil moisture
! line below inserted by rk4417 - phase2
REAL(r_2), DIMENSION(:,:),POINTER, SAVE :: bwb_gw ! volumetric gw soil moisture ! MMY
REAL(r_2), DIMENSION(mp) :: delwb ! change in soilmoisture
! b/w tsteps
REAL, DIMENSION(mp) :: canopy_wbal !canopy water balance
Expand All @@ -494,12 +499,39 @@ SUBROUTINE mass_balance(dels,ktau, ssnow,soil,canopy,met,

IF(ktau==1) THEN
ALLOCATE( bwb(mp,ms,2) )
! initial vlaue of soil moisture
ALLOCATE( bwb_gw(mp,2) ) ! MMY ! inserted by rk4417 - phase2
! initial value of soil moisture
bwb(:,:,1)=ssnow%wb
bwb(:,:,2)=ssnow%wb
bwb_gw(:,1)=ssnow%GWwb ! MMY ! 2 lines inserted by rk4417 - phase2
bwb_gw(:,2)=ssnow%GWwb ! MMY
delwb(:) = 0.
ELSE
! Calculate change in soil moisture b/w timesteps:

! ________________ MMY, Water Balance Equation for GW_ON _______________
! MMY to fix water balance bug when gw-on
IF ( cable_user%GW_MODEL) THEN ! MMY

IF(MOD(REAL(ktau),2.0)==1.0) THEN ! if odd timestep
bwb(:,:,1)=ssnow%wb
bwb_gw(:,1)=ssnow%GWwb
DO k=1,mp ! current smoist - prev tstep smoist
delwb(k) = ( SUM( (bwb(k,:,1) - bwb(k,:,2))*soil%zse ) + &
(bwb_gw(k,1) - bwb_gw(k,2))*soil%GWdz(k) ) *1000.0
END DO
ELSE IF(MOD(REAL(ktau),2.0)==0.0) THEN ! if even timestep
bwb(:,:,2)=ssnow%wb
bwb_gw(:,2)=ssnow%GWwb
DO k=1,mp ! current smoist - prev tstep smoist
delwb(k) = ( SUM( (bwb(k,:,2) - bwb(k,:,1))*soil%zse ) + &
(bwb_gw(k,2) -bwb_gw(k,1))*soil%GWdz(k) ) *1000.0
END DO
END IF
! ______________________________________________________________________

ELSE ! MMY ! IF part above inserted by rk4417 - phase2

IF(MOD(REAL(ktau),2.0)==1.0) THEN ! if odd timestep
bwb(:,:,1)=ssnow%wb
DO k=1,mp ! current smoist - prev tstep smoist
Expand All @@ -513,7 +545,10 @@ SUBROUTINE mass_balance(dels,ktau, ssnow,soil,canopy,met,
- (bwb(k,:,1)))*soil%zse)*1000.0
END DO
END IF
END IF

END IF ! MMY ! inserted by rk4417 - phase2

END IF



Expand All @@ -523,10 +558,21 @@ SUBROUTINE mass_balance(dels,ktau, ssnow,soil,canopy,met,
! it's included in change in canopy storage calculation))
! rml 28/2/11 ! BP changed rnof1+rnof2 to ssnow%runoff which also included rnof5
! which is used when nglacier=2 in soilsnow routines (BP feb2011)

IF ( cable_user%GW_MODEL) THEN ! MMY
! ________________ MMY, Water Balance Equation for GW_ON _______________
bal%wbal = REAL(met%precip - canopy%delwc - ssnow%snowd+ssnow%osnowd &
- ssnow%runoff-(canopy%fevw+canopy%fevc &
+ canopy%fes/ssnow%cls)*dels/air%rlam - delwb) ! remove qrecharge
! ______________________________________________________________________
ELSE ! MMY ! IF part above inserted by rk4417 - phase2

bal%wbal = REAL(met%precip - canopy%delwc - ssnow%snowd+ssnow%osnowd &
- ssnow%runoff-(canopy%fevw+canopy%fevc &
+ canopy%fes/ssnow%cls)*dels/air%rlam - delwb - ssnow%qrecharge)

END IF ! MMY ! inserted by rk4417 - phase2

! Canopy water balance: precip-change.can.storage-throughfall-evap+dew
canopy_wbal = REAL(met%precip-canopy%delwc-canopy%through &
- (canopy%fevw+MIN(canopy%fevc,0.0_r_2))*dels/air%rlam)
Expand Down
5 changes: 3 additions & 2 deletions src/offline/cable_define_types.F90
Original file line number Diff line number Diff line change
Expand Up @@ -210,8 +210,9 @@ MODULE cable_def_types_mod

! Soil and snow variables:
TYPE soil_snow_type

INTEGER, DIMENSION(:), POINTER :: isflag ! 0 => no snow 1 => snow
!> isflag 0 => one snow layer 1 => three snow layer
INTEGER, DIMENSION(:), POINTER :: isflag
! 0 => one snow layer 1 => three snow layer ! MMY ! above and here inserted by rk4417 - phase2

REAL, DIMENSION(:), POINTER :: &
iantrct, & ! pointer to Antarctic land points
Expand Down
97 changes: 91 additions & 6 deletions src/offline/cable_driver.F90
Original file line number Diff line number Diff line change
Expand Up @@ -62,7 +62,7 @@ PROGRAM cable_offline_driver
USE cable_def_types_mod
USE cable_IO_vars_module, ONLY: logn,gswpfile,ncciy,leaps, &
verbose, fixedCO2,output,check,patchout, &
patch_type,landpt,soilparmnew,&
patch_type,landpt,soilparmnew, &
defaultLAI, sdoy, smoy, syear, timeunits, exists, calendar, set_group_output_values, &
NO_CHECK
USE casa_ncdf_module, ONLY: is_casa_time
Expand Down Expand Up @@ -245,7 +245,7 @@ PROGRAM cable_offline_driver
NAMELIST/CABLE/ &
filename, & ! TYPE, containing input filenames
vegparmnew, & ! use new soil param. method
soilparmnew, & ! use new soil param. method
soilparmnew, & ! use new soil param. method
calcsoilalbedo, & ! albedo considers soil color Ticket #27
spinup, & ! spinup model (soil) to steady state
delsoilM,delsoilT,& !
Expand Down Expand Up @@ -409,7 +409,8 @@ PROGRAM cable_offline_driver
! STOP 'casaCNP required to get prognostic LAI or Vcmax'
IF( l_vcmaxFeedbk .AND. icycle < 1 ) &
STOP 'icycle must be 2 to 3 to get prognostic Vcmax'
IF( icycle > 0 .AND. ( .NOT. soilparmnew ) ) &

IF( icycle > 0 .AND. ( .NOT. soilparmnew ) ) &
STOP 'casaCNP must use new soil parameters'

NRRRR = MERGE(MAX(CABLE_USER%CASA_NREP,1), 1, CASAONLY)
Expand Down Expand Up @@ -512,6 +513,38 @@ PROGRAM cable_offline_driver
kend = ktauday * LOY
ENDIF

! __________________________ MMY using Princeton _______________________________
ELSE IF ( TRIM(cable_user%MetType) .EQ. 'prin' ) THEN
ncciy = CurYear

CALL prepareFiles_princeton(ncciy) ! MMY
IF ( RRRR .EQ. 1 ) THEN
CALL open_met_file( dels, koffset, kend, spinup, CTFRZ )
IF (leaps.and.is_leapyear(YYYY).and.kend.eq.2920) THEN
STOP 'LEAP YEAR INCOMPATIBILITY WITH INPUT MET !'
ENDIF
IF ( NRRRR .GT. 1 ) THEN
GSWP_MID(1,YYYY) = ncid_rain
! GSWP_MID(2,YYYY) = ncid_snow MMY
GSWP_MID(3,YYYY) = ncid_lw
GSWP_MID(4,YYYY) = ncid_sw
GSWP_MID(5,YYYY) = ncid_ps
GSWP_MID(6,YYYY) = ncid_qa
GSWP_MID(7,YYYY) = ncid_ta
GSWP_MID(8,YYYY) = ncid_wd
ENDIF
ELSE
ncid_rain = GSWP_MID(1,YYYY)
! ncid_snow = GSWP_MID(2,YYYY) MMY
ncid_lw = GSWP_MID(3,YYYY)
ncid_sw = GSWP_MID(4,YYYY)
ncid_ps = GSWP_MID(5,YYYY)
ncid_qa = GSWP_MID(6,YYYY)
ncid_ta = GSWP_MID(7,YYYY)
ncid_wd = GSWP_MID(8,YYYY)
kend = ktauday * LOY ! MMY
ENDIF

ELSE IF ( TRIM(cable_user%MetType) .EQ. 'plum' ) THEN
! PLUME experiment setup using WATCH
IF ( CALL1 ) THEN
Expand Down Expand Up @@ -1049,9 +1082,17 @@ PROGRAM cable_offline_driver
YYYY.EQ. CABLE_USER%YearEnd ) THEN

! evaluate spinup
! =================== MMY_phase2 commented out this region in favor of the one below - rk4417 =====================
IF( ANY( ABS(ssnow%wb-soilMtemp)>delsoilM).OR. &
ANY( ABS(ssnow%tgg-soilTtemp)>delsoilT) .OR. &
MAXVAL(ABS(ssnow%GWwb-GWtemp),dim=1) > delgwM) THEN
ANY( ABS(ssnow%tgg-soilTtemp)>delsoilT) .OR. &
MAXVAL(ABS(ssnow%GWwb-GWtemp),dim=1) > delgwM) THEN

! =================== MMY_phase2 uncomment =====================
! IF( (ANY( ABS(ssnow%wb-soilMtemp)>delsoilM).OR. &
! ANY( ABS(ssnow%tgg-soilTtemp)>delsoilT) .or. &
! maxval(ABS(ssnow%GWwb-GWtemp),dim=1) > delgwM) .and. &
! ( (int(ktau_tot/kend) .lt. cable_user%max_spins) .and.&
! (cable_user%max_spins .gt. 0) ) ) THEN

! No complete convergence yet
PRINT *, 'ssnow%wb : ', ssnow%wb
Expand Down Expand Up @@ -1185,7 +1226,7 @@ PROGRAM cable_offline_driver

IF (icycle > 0.and. .not.l_landuse) THEN

!CALL casa_poolout( ktau, veg, soil, casabiome, &
!CALL casa_poolout( ktau, veg, soil, casabiome, &
! casapool, casaflux, casamet, casabal, phen )
CALL write_casa_restart_nc ( casamet, casapool,casaflux,phen, CASAONLY )

Expand Down Expand Up @@ -1307,6 +1348,50 @@ SUBROUTINE renameFiles(logn,inFile,ncciy,inName)

END SUBROUTINE renameFiles


! 2 subroutines below inserted by rk4417 - phase2
! _______________________ MMY for Princeton input ______________________________
SUBROUTINE prepareFiles_princeton(ncciy)
USE cable_IO_vars_module, ONLY: logn,gswpfile
IMPLICIT NONE
INTEGER, INTENT(IN) :: ncciy

WRITE(logn,*) 'CABLE offline global run using princeton forcing for ', ncciy
PRINT *, 'CABLE offline global run using princeton forcing for ', ncciy

CALL renameFiles_princeton(logn,gswpfile%rainf,ncciy,'rainf')
CALL renameFiles_princeton(logn,gswpfile%LWdown,ncciy,'LWdown')
CALL renameFiles_princeton(logn,gswpfile%SWdown,ncciy,'SWdown')
CALL renameFiles_princeton(logn,gswpfile%PSurf,ncciy,'PSurf')
CALL renameFiles_princeton(logn,gswpfile%Qair,ncciy,'Qair')
CALL renameFiles_princeton(logn,gswpfile%Tair,ncciy,'Tair')
CALL renameFiles_princeton(logn,gswpfile%wind,ncciy,'wind')

END SUBROUTINE prepareFiles_princeton

SUBROUTINE renameFiles_princeton(logn,inFile,ncciy,inName)
IMPLICIT NONE
INTEGER, INTENT(IN) :: logn,ncciy
INTEGER:: nn
CHARACTER(LEN=200), INTENT(INOUT) :: inFile
CHARACTER(LEN=*), INTENT(IN) :: inName
INTEGER :: idummy

nn = INDEX(inFile,'19')
READ(inFile(nn:nn+3),'(i4)') idummy
WRITE(inFile(nn:nn+3),'(i4.4)') ncciy
nn = INDEX(inFile,'19')
READ(inFile(nn:nn+3),'(i4)') idummy
WRITE(inFile(nn:nn+3),'(i4.4)') ncciy
Comment on lines +1380 to +1385
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

why is it doing twice the same thing?

Is there any explanation anywhere of what this does?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I traced this back to the code from MMY (phase 1). The only explanation is what you see "MMY for Princeton input". It could be intentional or a simple duplication error. We can ping MMY about it of you like

READ(inFile(nn-5:nn-2),'(i4)') idummy
WRITE(inFile(nn-5:nn-2),'(i4.4)') ncciy
WRITE(logn,*) TRIM(inName), ' global data from ', TRIM(inFile)

END SUBROUTINE renameFiles_princeton

! ______________________________________________________________________________


!==============================================================================
! subroutine for reading LU input data, zeroing biomass in empty secondary forest tiles
! and tranferring LUC-based age weights for secondary forest to POP structure
Expand Down
Loading
Loading