From 79cfb7d596e80c7822ef0cfa2f964127c2f4f087 Mon Sep 17 00:00:00 2001 From: uturuncoglu Date: Sat, 9 Mar 2024 23:07:21 -0700 Subject: [PATCH 1/2] initial pass to fix issue in debug mode --- drivers/nuopc/lnd_comp_domain.F90 | 35 +++-- drivers/nuopc/lnd_comp_driver.F90 | 2 +- drivers/nuopc/lnd_comp_io.F90 | 150 ++++++++++++++----- drivers/nuopc/lnd_comp_types.F90 | 235 +++++++++++++++--------------- 4 files changed, 255 insertions(+), 167 deletions(-) diff --git a/drivers/nuopc/lnd_comp_domain.F90 b/drivers/nuopc/lnd_comp_domain.F90 index 92c9be8..db6bd8b 100644 --- a/drivers/nuopc/lnd_comp_domain.F90 +++ b/drivers/nuopc/lnd_comp_domain.F90 @@ -60,7 +60,8 @@ subroutine lnd_set_decomp_and_domain_from_mosaic(gcomp, noahmp, rc) integer, intent(out) :: rc ! local variables - real(r4), target, allocatable :: tmpr4(:) + real(r4), target, allocatable :: tmpr4(:) + real(r8), target, allocatable :: tmpr8(:) integer :: n integer :: decomptile(2,6) integer :: maxIndex(2) @@ -183,6 +184,10 @@ subroutine lnd_set_decomp_and_domain_from_mosaic(gcomp, noahmp, rc) allocate(tmpr4(noahmp%domain%begl:noahmp%domain%endl)) tmpr4(:) = 0.0 end if + if (.not. allocated(tmpr8)) then + allocate(tmpr8(noahmp%domain%begl:noahmp%domain%endl)) + tmpr8(:) = 0.0 + end if ! --------------------- ! Get fraction from orography file @@ -207,18 +212,27 @@ subroutine lnd_set_decomp_and_domain_from_mosaic(gcomp, noahmp, rc) ! following link: https://github.com/ufs-community/ufs-weather-model/issues/1423 ! --------------------- - ! read field - write(filename, fmt="(A,I0,A)") trim(noahmp%nmlist%input_dir)//'C',noahmp%domain%ni, '.vegetation_type.tile*.nc' - flds(1)%short_name = 'vegetation_type' - flds(1)%ptr1r4 => tmpr4 - call read_tiled_file(noahmp, filename, flds, rc=rc) - if (ChkErr(rc,__LINE__,u_FILE_u)) return - ! allocate data if (.not. allocated(vegtype)) then allocate(vegtype(noahmp%domain%begl:noahmp%domain%endl)) end if - vegtype(:) = int(tmpr4) + + ! read field + if (trim(noahmp%nmlist%ic_type) == 'sfc') then + filename = trim(noahmp%nmlist%input_dir)//'sfc_data.tile*.nc' + flds(1)%short_name = 'vtype' + flds(1)%ptr1r8 => tmpr8 + call read_tiled_file(noahmp, filename, flds, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + vegtype(:) = int(tmpr8) + else + write(filename, fmt="(A,I0,A)") trim(noahmp%nmlist%input_dir)//'C',noahmp%domain%ni, '.vegetation_type.tile*.nc' + flds(1)%short_name = 'vegetation_type' + flds(1)%ptr1r4 => tmpr4 + call read_tiled_file(noahmp, filename, flds, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + vegtype(:) = int(tmpr4) + end if ! --------------------- ! Calculate mask from land-sea fraction @@ -228,7 +242,7 @@ subroutine lnd_set_decomp_and_domain_from_mosaic(gcomp, noahmp, rc) allocate(noahmp%domain%mask(noahmp%domain%begl:noahmp%domain%endl)) end if - where (noahmp%domain%frac(:) > 0.0_r8 .and. vegtype(:) >= 0) + where (noahmp%domain%frac(:) > 0.0_r8 .and. vegtype(:) > 0) noahmp%domain%mask(:) = 1 elsewhere noahmp%domain%mask(:) = 0 @@ -319,6 +333,7 @@ subroutine lnd_set_decomp_and_domain_from_mosaic(gcomp, noahmp, rc) ! --------------------- if (allocated(tmpr4)) deallocate(tmpr4) + if (allocated(tmpr8)) deallocate(tmpr8) call ESMF_LogWrite(subname//' done', ESMF_LOGMSG_INFO) diff --git a/drivers/nuopc/lnd_comp_driver.F90 b/drivers/nuopc/lnd_comp_driver.F90 index 7ba04db..eec9fe3 100644 --- a/drivers/nuopc/lnd_comp_driver.F90 +++ b/drivers/nuopc/lnd_comp_driver.F90 @@ -283,7 +283,7 @@ subroutine drv_run(gcomp, noahmp, rc) end if ! initialize model variables - call noahmp%InitializeStates(noahmp%nmlist, noahmp%static, month) + call noahmp%InitializeStates(noahmp%nmlist, noahmp%static, noahmp%domain, month) end if end if diff --git a/drivers/nuopc/lnd_comp_io.F90 b/drivers/nuopc/lnd_comp_io.F90 index dfa08d9..60a53bc 100644 --- a/drivers/nuopc/lnd_comp_io.F90 +++ b/drivers/nuopc/lnd_comp_io.F90 @@ -467,7 +467,9 @@ subroutine read_static(noahmp, rc) ! local variables type(field_type), allocatable :: flds(:) real(r4), target, allocatable :: tmpr4(:) + real(r8), target, allocatable :: tmpr8(:) real(r4), target, allocatable :: tmp2r4(:,:) + real(r8), target, allocatable :: tmp2r8(:,:) character(len=CL) :: filename real(ESMF_KIND_R8), parameter :: pi_8 = 3.14159265358979323846_r8 character(len=*), parameter :: subname=trim(modName)//':(read_static) ' @@ -485,11 +487,21 @@ subroutine read_static(noahmp, rc) tmpr4(:) = 0.0 end if + if (.not. allocated(tmpr8)) then + allocate(tmpr8(noahmp%domain%begl:noahmp%domain%endl)) + tmpr8(:) = 0.0_r8 + end if + if (.not. allocated(tmp2r4)) then allocate(tmp2r4(noahmp%domain%begl:noahmp%domain%endl,12)) tmp2r4(:,:) = 0.0 end if + if (.not. allocated(tmp2r8)) then + allocate(tmp2r8(noahmp%domain%begl:noahmp%domain%endl,1)) + tmp2r8(:,:) = 0.0_r8 + end if + !---------------------- ! Read latitude, we could also retrive from ESMF mesh object !---------------------- @@ -509,66 +521,122 @@ subroutine read_static(noahmp, rc) ! Read soil type !---------------------- - allocate(flds(1)) - write(filename, fmt="(A,I0,A)") trim(noahmp%nmlist%input_dir)//'C', noahmp%domain%ni, '.soil_type.tile*.nc' - flds(1)%short_name = 'soil_type' - flds(1)%ptr1r4 => tmpr4 - call read_tiled_file(noahmp, filename, flds, rc=rc) - if (ChkErr(rc,__LINE__,u_FILE_u)) return - noahmp%model%soiltyp = int(tmpr4) - deallocate(flds) + if (trim(noahmp%nmlist%ic_type) == 'sfc') then + allocate(flds(1)) + write(filename, fmt="(A)") trim(noahmp%nmlist%input_dir)//'sfc_data.tile*.nc' + flds(1)%short_name = 'stype' + flds(1)%ptr1r8 => tmpr8 + call read_tiled_file(noahmp, filename, flds, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + noahmp%model%soiltyp = int(tmpr8) + deallocate(flds) + else + allocate(flds(1)) + write(filename, fmt="(A,I0,A)") trim(noahmp%nmlist%input_dir)//'C', noahmp%domain%ni, '.soil_type.tile*.nc' + flds(1)%short_name = 'soil_type' + flds(1)%ptr1r4 => tmpr4 + call read_tiled_file(noahmp, filename, flds, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + noahmp%model%soiltyp = int(tmpr4) + deallocate(flds) + end if !---------------------- ! Read vegetation type !---------------------- - allocate(flds(1)) - write(filename, fmt="(A,I0,A)") trim(noahmp%nmlist%input_dir)//'C', noahmp%domain%ni, '.vegetation_type.tile*.nc' - flds(1)%short_name = 'vegetation_type' - flds(1)%ptr1r4 => tmpr4 - call read_tiled_file(noahmp, filename, flds, rc=rc) - if (ChkErr(rc,__LINE__,u_FILE_u)) return - noahmp%model%vegtype = int(tmpr4) - deallocate(flds) + if (trim(noahmp%nmlist%ic_type) == 'sfc') then + allocate(flds(1)) + write(filename, fmt="(A)") trim(noahmp%nmlist%input_dir)//'sfc_data.tile*.nc' + flds(1)%short_name = 'vtype' + flds(1)%ptr1r8 => tmpr8 + call read_tiled_file(noahmp, filename, flds, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + noahmp%model%vegtype = int(tmpr8) + deallocate(flds) + else + allocate(flds(1)) + write(filename, fmt="(A,I0,A)") trim(noahmp%nmlist%input_dir)//'C', noahmp%domain%ni, '.vegetation_type.tile*.nc' + flds(1)%short_name = 'vegetation_type' + flds(1)%ptr1r4 => tmpr4 + call read_tiled_file(noahmp, filename, flds, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + noahmp%model%vegtype = int(tmpr4) + deallocate(flds) + end if !---------------------- ! Read slope type !---------------------- - allocate(flds(1)) - write(filename, fmt="(A,I0,A)") trim(noahmp%nmlist%input_dir)//'C', noahmp%domain%ni, '.slope_type.tile*.nc' - flds(1)%short_name = 'slope_type' - flds(1)%ptr1r4 => tmpr4 - call read_tiled_file(noahmp, filename, flds, rc=rc) - if (ChkErr(rc,__LINE__,u_FILE_u)) return - noahmp%model%slopetyp = int(tmpr4) - deallocate(flds) + if (trim(noahmp%nmlist%ic_type) == 'sfc') then + allocate(flds(1)) + write(filename, fmt="(A)") trim(noahmp%nmlist%input_dir)//'sfc_data.tile*.nc' + flds(1)%short_name = 'slope' + flds(1)%ptr1r8 => tmpr8 + call read_tiled_file(noahmp, filename, flds, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + noahmp%model%slopetyp = int(tmpr8) + deallocate(flds) + else + allocate(flds(1)) + write(filename, fmt="(A,I0,A)") trim(noahmp%nmlist%input_dir)//'C', noahmp%domain%ni, '.slope_type.tile*.nc' + flds(1)%short_name = 'slope_type' + flds(1)%ptr1r4 => tmpr4 + call read_tiled_file(noahmp, filename, flds, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + noahmp%model%slopetyp = int(tmpr4) + deallocate(flds) + end if !---------------------- ! Read deep soil temperature !---------------------- - allocate(flds(1)) - write(filename, fmt="(A,I0,A)") trim(noahmp%nmlist%input_dir)//'C', noahmp%domain%ni, '.substrate_temperature.tile*.nc' - flds(1)%short_name = 'substrate_temperature' - flds(1)%ptr1r4 => tmpr4 - call read_tiled_file(noahmp, filename, flds, rc=rc) - if (ChkErr(rc,__LINE__,u_FILE_u)) return - noahmp%model%tg3 = dble(tmpr4) - deallocate(flds) + if (trim(noahmp%nmlist%ic_type) == 'sfc') then + allocate(flds(1)) + write(filename, fmt="(A)") trim(noahmp%nmlist%input_dir)//'sfc_data.tile*.nc' + flds(1)%short_name = 'tg3' + flds(1)%nrec = 1; flds(1)%ptr2r8 => tmp2r8 + call read_tiled_file(noahmp, filename, flds, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + print*, 'tg3, lb, ub = ', lbound(tmp2r8, dim=1), ubound(tmp2r8, dim=1) + noahmp%model%tg3 = tmp2r8(:,1) + deallocate(flds) + else + allocate(flds(1)) + write(filename, fmt="(A,I0,A)") trim(noahmp%nmlist%input_dir)//'C', noahmp%domain%ni, '.substrate_temperature.tile*.nc' + flds(1)%short_name = 'substrate_temperature' + flds(1)%ptr1r4 => tmpr4 + call read_tiled_file(noahmp, filename, flds, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + noahmp%model%tg3 = dble(tmpr4) + deallocate(flds) + end if !---------------------- ! Read maximum snow albedo !---------------------- - allocate(flds(1)) - write(filename, fmt="(A,I0,A)") trim(noahmp%nmlist%input_dir)//'C', noahmp%domain%ni, '.maximum_snow_albedo.tile*.nc' - flds(1)%short_name = 'maximum_snow_albedo' - flds(1)%ptr1r4 => tmpr4 - call read_tiled_file(noahmp, filename, flds, rc=rc) - if (ChkErr(rc,__LINE__,u_FILE_u)) return - noahmp%model%snoalb = dble(tmpr4) - deallocate(flds) + if (trim(noahmp%nmlist%ic_type) == 'sfc') then + allocate(flds(1)) + write(filename, fmt="(A)") trim(noahmp%nmlist%input_dir)//'sfc_data.tile*.nc' + flds(1)%short_name = 'snoalb' + flds(1)%nrec = 1; flds(1)%ptr2r8 => tmp2r8 + call read_tiled_file(noahmp, filename, flds, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + noahmp%model%snoalb = tmp2r8(:,1) + deallocate(flds) + else + allocate(flds(1)) + write(filename, fmt="(A,I0,A)") trim(noahmp%nmlist%input_dir)//'C', noahmp%domain%ni, '.maximum_snow_albedo.tile*.nc' + flds(1)%short_name = 'maximum_snow_albedo' + flds(1)%ptr1r4 => tmpr4 + call read_tiled_file(noahmp, filename, flds, rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + noahmp%model%snoalb = dble(tmpr4) + deallocate(flds) + end if !---------------------- ! Read vegetation greenness, monthly average @@ -610,7 +678,9 @@ subroutine read_static(noahmp, rc) !---------------------- if (allocated(tmpr4)) deallocate(tmpr4) + if (allocated(tmpr8)) deallocate(tmpr8) if (allocated(tmp2r4)) deallocate(tmp2r4) + if (allocated(tmp2r8)) deallocate(tmp2r8) call ESMF_LogWrite(subname//' done for '//trim(filename), ESMF_LOGMSG_INFO) diff --git a/drivers/nuopc/lnd_comp_types.F90 b/drivers/nuopc/lnd_comp_types.F90 index d263322..bec9fbd 100644 --- a/drivers/nuopc/lnd_comp_types.F90 +++ b/drivers/nuopc/lnd_comp_types.F90 @@ -675,7 +675,7 @@ subroutine InitializeDefault(this) end subroutine InitializeDefault - subroutine InitializeStates(this, namelist, static, month) + subroutine InitializeStates(this, namelist, static, domain, month) ! out of necessity, this is extracted from FV3GFS_io.F90 @@ -685,140 +685,143 @@ subroutine InitializeStates(this, namelist, static, month) class(noahmp_type) :: this type(namelist_type) :: namelist type(static_type) :: static + type(domain_type) :: domain integer :: iloc, ilevel, isnow real :: masslai, masssai, depthm real :: dzsno(static%lsnowl:static%km) do iloc = 1, static%im - this%model%tvxy (iloc) = this%model%tskin(iloc) - this%model%tgxy (iloc) = this%model%tskin(iloc) - this%model%tahxy(iloc) = this%model%tskin(iloc) + if (domain%mask(iloc) > 0) then + this%model%tvxy (iloc) = this%model%tskin(iloc) + this%model%tgxy (iloc) = this%model%tskin(iloc) + this%model%tahxy(iloc) = this%model%tskin(iloc) - if (this%model%snwdph(iloc) > 0.01 .and. this%model%tskin(iloc) > 273.15 ) then - this%model%tvxy (iloc) = 273.15 - this%model%tgxy (iloc) = 273.15 - this%model%tahxy(iloc) = 273.15 - end if + if (this%model%snwdph(iloc) > 0.01 .and. this%model%tskin(iloc) > 273.15 ) then + this%model%tvxy (iloc) = 273.15 + this%model%tgxy (iloc) = 273.15 + this%model%tahxy(iloc) = 273.15 + end if - this%model%canicexy(iloc) = 0.0 - this%model%canliqxy(iloc) = this%model%canopy(iloc) - this%model%eahxy(iloc) = 2000.0 - this%model%cmxy (iloc) = 0.0 - this%model%chxy (iloc) = 0.0 - this%model%fwetxy (iloc) = 0.0 - this%model%sneqvoxy(iloc) = this%model%weasd(iloc) - this%model%alboldxy(iloc) = 0.65 - this%model%qsnowxy (iloc) = 0.0 - this%model%wslakexy(iloc) = 0.0 - this%model%taussxy (iloc) = 0.0 - this%model%waxy (iloc) = 4900.0 ! this assumes water table is at 2.5m - this%model%wtxy (iloc) = this%model%waxy(iloc) - this%model%zwtxy (iloc) = (25.0 + 2.0)-this%model%waxy(iloc)/1000.0/0.2 + this%model%canicexy(iloc) = 0.0 + this%model%canliqxy(iloc) = this%model%canopy(iloc) + this%model%eahxy(iloc) = 2000.0 + this%model%cmxy (iloc) = 0.0 + this%model%chxy (iloc) = 0.0 + this%model%fwetxy (iloc) = 0.0 + this%model%sneqvoxy(iloc) = this%model%weasd(iloc) + this%model%alboldxy(iloc) = 0.65 + this%model%qsnowxy (iloc) = 0.0 + this%model%wslakexy(iloc) = 0.0 + this%model%taussxy (iloc) = 0.0 + this%model%waxy (iloc) = 4900.0 ! this assumes water table is at 2.5m + this%model%wtxy (iloc) = this%model%waxy(iloc) + this%model%zwtxy (iloc) = (25.0 + 2.0)-this%model%waxy(iloc)/1000.0/0.2 - if ((this%model%vegtype(iloc) == isbarren_table) .or. (this%model%vegtype(iloc) == isice_table) .or. & - (this%model%vegtype(iloc) == isurban_table) .or. (this%model%vegtype(iloc) == iswater_table) .or. & - (this%model%vegtype(iloc) < 0)) then - this%model%xlaixy (iloc) = 0.0 - this%model%xsaixy (iloc) = 0.0 - this%model%lfmassxy(iloc) = 0.0 - this%model%stmassxy(iloc) = 0.0 - this%model%rtmassxy(iloc) = 0.0 - this%model%woodxy (iloc) = 0.0 - this%model%stblcpxy(iloc) = 0.0 - this%model%fastcpxy(iloc) = 0.0 - else - this%model%xlaixy (iloc) = max(laim_table(this%model%vegtype(iloc), month),0.05) - this%model%xsaixy (iloc) = max(this%model%xlaixy(iloc)*0.1,0.05) - masslai = 1000.0/max(sla_table(this%model%vegtype(iloc)),1.0) - this%model%lfmassxy(iloc) = this%model%xlaixy(iloc)*masslai - masssai = 1000.0/3.0 - this%model%stmassxy(iloc) = this%model%xsaixy(iloc)*masssai - this%model%rtmassxy(iloc) = 500.0 - this%model%woodxy (iloc) = 500.0 - this%model%stblcpxy(iloc) = 1000.0 - this%model%fastcpxy(iloc) = 1000.0 - end if + if ((this%model%vegtype(iloc) == isbarren_table) .or. (this%model%vegtype(iloc) == isice_table) .or. & + (this%model%vegtype(iloc) == isurban_table) .or. (this%model%vegtype(iloc) == iswater_table) .or. & + (this%model%vegtype(iloc) < 0)) then + this%model%xlaixy (iloc) = 0.0 + this%model%xsaixy (iloc) = 0.0 + this%model%lfmassxy(iloc) = 0.0 + this%model%stmassxy(iloc) = 0.0 + this%model%rtmassxy(iloc) = 0.0 + this%model%woodxy (iloc) = 0.0 + this%model%stblcpxy(iloc) = 0.0 + this%model%fastcpxy(iloc) = 0.0 + else + this%model%xlaixy (iloc) = max(laim_table(this%model%vegtype(iloc), month),0.05) + this%model%xsaixy (iloc) = max(this%model%xlaixy(iloc)*0.1,0.05) + masslai = 1000.0/max(sla_table(this%model%vegtype(iloc)),1.0) + this%model%lfmassxy(iloc) = this%model%xlaixy(iloc)*masslai + masssai = 1000.0/3.0 + this%model%stmassxy(iloc) = this%model%xsaixy(iloc)*masssai + this%model%rtmassxy(iloc) = 500.0 + this%model%woodxy (iloc) = 500.0 + this%model%stblcpxy(iloc) = 1000.0 + this%model%fastcpxy(iloc) = 1000.0 + end if - if (this%model%vegtype(iloc) == isice_table ) then - do ilevel = 1, namelist%num_soil_levels - this%model%stc(iloc,ilevel) = min(this%model%stc(iloc,ilevel),min(this%model%tg3(iloc),263.15)) - end do - this%model%smc(iloc,:) = 1.0 - this%model%slc(iloc,:) = 0.0 - end if + if (this%model%vegtype(iloc) == isice_table ) then + do ilevel = 1, namelist%num_soil_levels + this%model%stc(iloc,ilevel) = min(this%model%stc(iloc,ilevel),min(this%model%tg3(iloc),263.15)) + end do + this%model%smc(iloc,:) = 1.0 + this%model%slc(iloc,:) = 0.0 + end if - depthm = this%model%snwdph(iloc)/1000.0 ! snow depth in meters + depthm = this%model%snwdph(iloc)/1000.0 ! snow depth in meters - if (this%model%weasd(iloc) /= 0.0 .and. depthm == 0.0 ) then - depthm = this%model%weasd(iloc)/1000.0*10.0 ! assume 10/1 ratio - endif + if (this%model%weasd(iloc) /= 0.0 .and. depthm == 0.0 ) then + depthm = this%model%weasd(iloc)/1000.0*10.0 ! assume 10/1 ratio + endif - if (this%model%vegtype(iloc) == isice_table) then ! land ice in MODIS/IGBP - if (this%model%weasd(iloc) < 0.1) then - this%model%weasd(iloc) = 0.1 - depthm = this%model%weasd(iloc)/1000.0*10.0 - end if - end if + if (this%model%vegtype(iloc) == isice_table) then ! land ice in MODIS/IGBP + if (this%model%weasd(iloc) < 0.1) then + this%model%weasd(iloc) = 0.1 + depthm = this%model%weasd(iloc)/1000.0*10.0 + end if + end if - dzsno = 0.0 - if (depthm < 0.025 ) then - this%model%snowxy(iloc) = 0.0 - dzsno(-2:0) = 0.0 - elseif (depthm >= 0.025 .and. depthm <= 0.05 ) then - this%model%snowxy(iloc) = -1.0 - dzsno(0) = depthm - elseif (depthm > 0.05 .and. depthm <= 0.10 ) then - this%model%snowxy(iloc) = -2.0 - dzsno(-1) = 0.5*depthm - dzsno(0) = 0.5*depthm - elseif (depthm > 0.10 .and. depthm <= 0.25 ) then - this%model%snowxy(iloc) = -2.0 - dzsno(-1) = 0.05 - dzsno(0) = depthm - 0.05 - elseif (depthm > 0.25 .and. depthm <= 0.45 ) then - this%model%snowxy(iloc) = -3.0 - dzsno(-2) = 0.05 - dzsno(-1) = 0.5*(depthm-0.05) - dzsno(0) = 0.5*(depthm-0.05) - elseif (depthm > 0.45) then - this%model%snowxy(iloc) = -3.0 - dzsno(-2) = 0.05 - dzsno(-1) = 0.20 - dzsno(0) = depthm - 0.05 - 0.20 - endif + dzsno = 0.0 + if (depthm < 0.025 ) then + this%model%snowxy(iloc) = 0.0 + dzsno(-2:0) = 0.0 + elseif (depthm >= 0.025 .and. depthm <= 0.05 ) then + this%model%snowxy(iloc) = -1.0 + dzsno(0) = depthm + elseif (depthm > 0.05 .and. depthm <= 0.10 ) then + this%model%snowxy(iloc) = -2.0 + dzsno(-1) = 0.5*depthm + dzsno(0) = 0.5*depthm + elseif (depthm > 0.10 .and. depthm <= 0.25 ) then + this%model%snowxy(iloc) = -2.0 + dzsno(-1) = 0.05 + dzsno(0) = depthm - 0.05 + elseif (depthm > 0.25 .and. depthm <= 0.45 ) then + this%model%snowxy(iloc) = -3.0 + dzsno(-2) = 0.05 + dzsno(-1) = 0.5*(depthm-0.05) + dzsno(0) = 0.5*(depthm-0.05) + elseif (depthm > 0.45) then + this%model%snowxy(iloc) = -3.0 + dzsno(-2) = 0.05 + dzsno(-1) = 0.20 + dzsno(0) = depthm - 0.05 - 0.20 + endif - ! Now we have the snowxy field - ! snice + snliq + tsno allocation and compute them from what we have - this%model%tsnoxy (iloc,-2:0) = 0.0 - this%model%snicexy(iloc,-2:0) = 0.0 - this%model%snliqxy(iloc,-2:0) = 0.0 - this%model%zsnsoxy(iloc,-2:namelist%num_soil_levels) = 0.0 + ! Now we have the snowxy field + ! snice + snliq + tsno allocation and compute them from what we have + this%model%tsnoxy (iloc,-2:0) = 0.0 + this%model%snicexy(iloc,-2:0) = 0.0 + this%model%snliqxy(iloc,-2:0) = 0.0 + this%model%zsnsoxy(iloc,-2:namelist%num_soil_levels) = 0.0 - isnow = nint(this%model%snowxy(iloc))+1 ! snowxy <=0.0, dzsno >= 0.0 + isnow = nint(this%model%snowxy(iloc))+1 ! snowxy <=0.0, dzsno >= 0.0 - do ilevel = isnow, 0 - this%model%tsnoxy (iloc,ilevel) = this%model%tgxy(iloc) - this%model%snliqxy(iloc,ilevel) = 0.0 - this%model%snicexy(iloc,ilevel) = dzsno(ilevel)/depthm*this%model%weasd(iloc) - end do + do ilevel = isnow, 0 + this%model%tsnoxy (iloc,ilevel) = this%model%tgxy(iloc) + this%model%snliqxy(iloc,ilevel) = 0.0 + this%model%snicexy(iloc,ilevel) = dzsno(ilevel)/depthm*this%model%weasd(iloc) + end do - this%model%zsnsoxy(iloc,isnow) = -1.0*dzsno(isnow) - do ilevel = isnow+1, 0 - this%model%zsnsoxy(iloc,ilevel) = this%model%zsnsoxy(iloc,ilevel-1)-dzsno(ilevel) - end do - do ilevel = 1, namelist%num_soil_levels - this%model%zsnsoxy(iloc,ilevel) = this%model%zsnsoxy(iloc,ilevel-1)-namelist%soil_level_thickness(ilevel) - end do + this%model%zsnsoxy(iloc,isnow) = -1.0*dzsno(isnow) + do ilevel = isnow+1, 0 + this%model%zsnsoxy(iloc,ilevel) = this%model%zsnsoxy(iloc,ilevel-1)-dzsno(ilevel) + end do + do ilevel = 1, namelist%num_soil_levels + this%model%zsnsoxy(iloc,ilevel) = this%model%zsnsoxy(iloc,ilevel-1)-namelist%soil_level_thickness(ilevel) + end do - ! TODO: Should not need to initialize these - this%model%smoiseq(iloc,:) = 0.0 - this%model%smcwtdxy(iloc) = 0.0 - this%model%deeprechxy(iloc) = 0.0 - this%model%rechxy(iloc) = 0.0 - - ! TODO: Fixed number given. This could be coupling field - this%model%rhonewsn1(iloc)= 200.0 + ! TODO: Should not need to initialize these + this%model%smoiseq(iloc,:) = 0.0 + this%model%smcwtdxy(iloc) = 0.0 + this%model%deeprechxy(iloc) = 0.0 + this%model%rechxy(iloc) = 0.0 + + ! TODO: Fixed number given. This could be coupling field + this%model%rhonewsn1(iloc)= 200.0 + end if end do ! iloc end subroutine InitializeStates From c0b9ac3a545fc956e2a4064440b007a689e8a6cf Mon Sep 17 00:00:00 2001 From: Ufuk Turuncoglu Date: Tue, 12 Mar 2024 22:55:05 -0500 Subject: [PATCH 2/2] clean debug print --- drivers/nuopc/lnd_comp_io.F90 | 1 - 1 file changed, 1 deletion(-) diff --git a/drivers/nuopc/lnd_comp_io.F90 b/drivers/nuopc/lnd_comp_io.F90 index 60a53bc..2bd3299 100644 --- a/drivers/nuopc/lnd_comp_io.F90 +++ b/drivers/nuopc/lnd_comp_io.F90 @@ -600,7 +600,6 @@ subroutine read_static(noahmp, rc) flds(1)%nrec = 1; flds(1)%ptr2r8 => tmp2r8 call read_tiled_file(noahmp, filename, flds, rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return - print*, 'tg3, lb, ub = ', lbound(tmp2r8, dim=1), ubound(tmp2r8, dim=1) noahmp%model%tg3 = tmp2r8(:,1) deallocate(flds) else