Bug fix of mo_restart on reading the correct dimension of Landcover (years), LAI (time) and Soil (horizons)

Closes: #225 (closed) #226 (closed)

Bug fix: mo_restart now reads the correct dimension for the number of dimension slices of landcover, soil and LAI

Minor: fixed incorrect restart out path for test domain 2 in namelist

assigned to @muellese

mentioned in issue #225 (closed)

mentioned in issue #226 (closed)

Thanks @shresthp for putting this together. I think it would be better to use the solution from @ottor, since my "transpose" idea was really just a hotfix. It could/should look like this:

    ! get the dimensions
    var = nc%getVariable(trim(soilHorizonsVarName)//'_bnds')
    call var%getData(dummyD2)
    if ( read_old_style_restart_bounds ) dummyD2 = transpose(dummyD2)
    nSoilHorizons_temp = size(dummyD2, 2)
    allocate(soilHorizonBoundaries_temp(nSoilHorizons_temp+1))
    soilHorizonBoundaries_temp(1:nSoilHorizons_temp) = dummyD2(1, :)
    soilHorizonBoundaries_temp(nSoilHorizons_temp+1) = dummyD2(2, nSoilHorizons_temp)

    ! get the landcover dimension
    var = nc%getVariable(trim(landCoverPeriodsVarName)//'_bnds')
    call var%getData(dummyD2)
    if ( read_old_style_restart_bounds ) dummyD2 = transpose(dummyD2)
    nLandCoverPeriods_temp = size(dummyD2, 2)
    allocate(landCoverPeriodBoundaries_temp(nLandCoverPeriods_temp+1))
    landCoverPeriodBoundaries_temp(1:nLandCoverPeriods_temp) = dummyD2(1, :)
    landCoverPeriodBoundaries_temp(nLandCoverPeriods_temp+1) = dummyD2(2, nLandCoverPeriods_temp)

    ! get the LAI dimension
    if (nc%hasVariable(trim(LAIVarName)//'_bnds')) then
      var = nc%getVariable(trim(LAIVarName)//'_bnds')
      call var%getData(dummyD2)
      if ( read_old_style_restart_bounds ) dummyD2 = transpose(dummyD2)
      nLAIs_temp = size(dummyD2, 2)
      allocate(LAIBoundaries_temp(nLAIs_temp+1))
      LAIBoundaries_temp(1:nLAIs_temp) = dummyD2(1, :)
      LAIBoundaries_temp(nLAIs_temp+1) = dummyD2(2, nLAIs_temp)
    else if (nc%hasDimension('L1_LAITimesteps')) then
      nc_dim = nc%getDimension('L1_LAITimesteps')
      nLAIs_temp = nc_dim%getLength()
      allocate(LAIBoundaries_temp(nLAIs_temp+1))
      LAIBoundaries_temp = [(ii, ii=1, nLAIs_temp+1)]
    end if

Still some problems:

1. check case 1 is failing now, since the bounds are stored transposed in the "old-style" restart files. We would need to add a switch to the mainconfig_mhm_mrm namelist (read_old_style_restart_bounds maybe) to state, that you want to use an old-style restart file (mhm<=v5.11)
2. written restart files now use "new-style" bounds, that are transposed. I don't think that we would want to add a switch here, since that would be only needed if you write a restart file with mhm v5.12 and want to read it with an older version.

So my roadmap proposal here would be:

1. add read_old_style_restart_bounds as bool to mainconfig_mhm_mrm namelist, with default value of .false..
2. use read_old_style_restart_bounds in mo_restart::read_restart_states to transpose read bounds if wanted
3. set read_old_style_restart_bounds = .true. in check/case_01/mhm.nml

changed the description

Thanks @shresthp for putting this together. I agree with @muellese that we should switch to new bounds (this needs to be properly communicate to all users to prevent confusion/frustration). That said, I don't fully understand how you would use your switch here. I would expect that there would be an if statement in your above code to make use of the switch and do the transform as necessary. Is my understanding correct?

The example code only demonstrated the correct reading. I updated the code-snippet with the following line added at three places:

if ( read_old_style_restart_bounds ) dummyD2 = transpose(dummyD2)

added 1 commit

• 60228781 - Bug fix contd... : mo_restart and the correct dimension for the number of...

Compare with previous version

OK I made the prescribed changes from @muellese and pushed. Locally I made a simple restart test with test domains and it worked. Haven't checked check case 01, I leave it to our beloved CI, its running as I type this.

@shresthp thank you! Didn't expect this that fast. Looks good to me. I hope all compilers can deal with the auto-reallocation of dummyD2, so we can prevent introducing a temporary variable.

Sigh. In debug-mode the compilers don't like implicit reallocation (which is a good thing of course). OK, than we need a temporary variable (since I wouldn't change compiler flags because of that):

    call var%getData(dummyD2_tmp)
    if ( read_old_style_restart_bounds ) then
      dummyD2 = transpose(dummyD2_tmp)
    else
      dummyD2 = dummyD2_tmp
    end if
    deallocate(dummyD2_tmp)

I wouldn't change this in the standard nml.

changed this line in version 4 of the diff

added 1 commit

• 35761dcc - Bug fix contd... : mo_restart and avoiding implicit reallocation with dummyD2_tmp

Compare with previous version

added 1 commit

• 790a6601 - Reverting the read_restart flag in the standard mhm.nml

Compare with previous version

Made couple of commits for avoiding implicit re-allocation and reinstated the default flag for read_restart.

For gfortran 7.3 we need additional allocate statements for dummyD2 (and maybe also allocated checks with deallocate before that)

added 1 commit

• 53bd6dfc - Allocation and checks for dummyD2

Compare with previous version

added 1 commit

• c88cc796 - Allocation and checks for dummyD2 ... redone

Compare with previous version

Nice! LGTM!

marked this merge request as ready

approved this merge request

changed milestone to %v5.12

requested review from @muellese