WIP: Prototype water state breakout

[billsacks]

This PR is not meant to be merged: it is for discussion only.

This records some (evolving) thoughts on how to handle the separation of WaterStateType for water isotopes and water tracers in general (#358).

[billsacks]

@mathewvrothstein @martynpclark @mvertens I'm recording some evolving thoughts here. There are a couple more iterations (at least) that I still plan to commit.

[billsacks]

Here are some diagrams of the three designs I have prototyped (in the three commits currently on this branch) as well as two other designs that I have not yet prototyped:

WaterType possible designs.pdf

Below, I'll refer to the different designs as (D1), (D2), etc., and the enumerated list of benefits as (1), (2), etc.

The first three designs (D1-D3) correspond to the 3 commits on this branch. They start fairly simple and then get a bit more complex via extra layers – though the extra complexity could result in some extra simplicity / clarity in the science code:

Benefits of (D2) over (D1):

1. This groups together instances relating to the same tracer - e.g.,
   the state and diagnostic variables for a given tracer (and later the
   fluxes for that tracer, too).

2. This potentially will simplify subroutine argument lists, by allowing
   us to pass a single container instance, e.g., for the water bulk
   variables.

Benefits of (D3) over (D2) - and also over (D1):

3. It reduces the number of arguments to some calls, because
   water_bulk_inst and water_bulk_only_inst are packaged together now into
   a single instance. Also, this packaging is now more apparent.

4. These communicate intention better in subroutine argument lists:
   subroutines can expect an argument of bulk type, rather than the
   generic type, ensuring that you're not passing a tracer instance into
   a routine that expects bulk.

5. This should make it easier to migrate variables between bulk-only and
   bulk-and-tracer.

There's an important addition to (5): A related point here is that science modules don't need to know whether a given variable comes from (for example) waterstate_type or water_bulk_only_state_type.

To see the code for each of these designs, you can check out the corresponding commits on this branch:

git clone https://github.com/billsacks/ctsm.git ctsm_bill
cd ctsm_bill
# For design 1:
git checkout d3c26106
# For design 2:
git checkout 9db4f12c
# For design 3:
git checkout f5d2149b

I ran designs (D1) and (D3) by @mvertens and @mathewvrothstein yesterday. @mathewvrothstein raised the very valid critique that the extra layers in (D3) - and in particular the use of type extension - might make it hard to understand. On the other hand, they both liked that (D3) allows science routines to explicitly declare that they expect the bulk water state, not just any water state (since the latter could, in principle, be some tracer state, which would not be valid for the science of that routine) - i.e., benefit (4) above. Personally, I actually feel that benefit (5) - which is connected with benefit (3) - will be even more important in terms of making the code easier to write and maintain.

In an attempt to keep some of the benefits of (D3) while simplifying the design, I came up with (D5). I feel like this is a reasonable design as long as we don't break up the state/diag/flux types much more: If we go with design (D5), then we won't have containers like waterbulk_inst that we can pass to science routines. I don't really like the idea of passing the whole water_inst (which includes all the tracers) to all science routines, because it makes it hard to track where water tracers are dealt with in the code. So this would mean that we'd need to pass the separate instances of waterstate_bulk_inst, waterdiag_bulk_inst, etc. Again, this is okay if we just have a few such instances. But if we break these down further, argument lists to science routines will become unmanageable - though we could always re-introduce one or more container types at that point.

I also tried out a different inheritance hierarchy in (D4), but I don't like this as much as (D3) or (D5). This gives benefits (3) and (4), but not benefit (5) – and we further lose some of these benefits if we want to explicitly pass (e.g.) just the state variables into a routine, as opposed to passing the whole waterbulk_inst.

My current preferences are for (3) or (5), but I can live with any of these designs – or possibly others that I haven't thought of which other people might see.

[billsacks]

Here are the different options for what the subroutine calls from the driver to science routines could look like under design options (D3) and (D5):

! ========================================================================
! For design (D3) we have two options for how to pass arguments from the driver to the
! science routines.
! ========================================================================

! ------------------------------------------------------------------------
! Option 1: Pass in the whole waterbulk_inst container. This keeps argument lists shorter
! and less volatile, and probably simplifies the initial refactoring. On the other hand,
! it makes it harder to track the flow of individual types (e.g., waterstate_bulk_type)
! through the system, and you can't (for example) declare that state variables are
! intent(in) but diagnostic variables are intent(inout) for a given subroutine.
! ------------------------------------------------------------------------

! Driver call
call canopy_hydrology(..., water_inst%waterbulk_inst, ...)

! Subroutine header
subroutine canopy_hydrology(..., waterbulk_inst, ...)
  type(waterbulk_type), intent(inout) :: waterbulk_inst

  associate( &
       foo => waterbulk_inst%waterstate_bulk_inst%foo, &  ! Input
       bar => waterbulk_inst%waterdiag_bulk_inst%bar   &  ! Output
       )
  end associate
end subroutine canopy_hydrology


! ------------------------------------------------------------------------
! Option 2: Pass in the individual instances from waterbulk_inst. This has the opposite
! pros and cons of option 1.
! ------------------------------------------------------------------------

! Driver call
call canopy_hydrology(..., water_inst%waterbulk_inst%waterstate_bulk_inst, &
     water_inst%waterbulk_inst%waterdiag_bulk_inst, ...)

! Subroutine header
subroutine canopy_hydrology(..., waterstate_bulk_inst, waterdiag_bulk_inst, ...)
  type(waterstate_bulk_type), intent(in) :: waterstate_bulk_inst
  type(waterdiag_bulk_type), intent(inout) :: waterdiag_bulk_inst

  associate( &
       foo => waterstate_bulk_inst%foo, &  ! Input
       bar => waterdiag_bulk_inst%bar   &  ! Output
       )
  end associate
end subroutine canopy_hydrology

! ========================================================================
! For designs (D1) and (D5), option 1 is not possible. Here is what option 2 would look
! like if we go with design (D5).
! ========================================================================

! Driver call
call canopy_hydrology(..., water_inst%waterstate_bulk_inst, &
     water_inst%waterdiag_bulk_inst, ...)

! Subroutine header
! Same as above

If people prefer option (2) then there is little benefit to (D3) over (D5). If people prefer option (1) then designs (D1) and (D5) are ruled out.

[billsacks]

Okay, I'm mostly done commenting on this for now. @martynpclark @mvertens @mathewvrothstein and others: I'd welcome any thoughts you have.

[billsacks]

Oh, and if you review this: I don't really recommend looking at diffs in the different commits: instead, I recommend looking at the individual commits, as I mention in #395 (comment)

[billsacks]

As I discussed with @mvertens yesterday: I've been struggling to see how to apply the water isotope / tracer design to more modularized fluxes, where the flux variable is contained by the science module that computes it. An example of this is the irrigation flux in IrrigationMod; there are other examples in the unified_land_model branch. The problem is: it's likely that some of the variables in a science-oriented class (such as the fluxes themselves) will need a separate instance for each isotope/tracer, but some of the variables (e.g., various auxiliary variables that are used to compute the bulk water fluxes) should only have a single instance for the bulk water.

This could be easier to handle if we ditched the idea of handling isotopes and tracers via multiple instances, and instead went back to the idea of having an extra array dimension for the tracers. But this carries its own set of issues – and even then, this would mean that all of the water-related science classes need to know something about the different water tracers, making it harder to write (and read) the code for each science class.

One possible solution is to have two types for each science module: e.g., an irrigation_type that contains things that exist for tracers, too; and an irrigation_bulk_type that contains things that are just relevant for the bulk. But this adds complexity to the code.

Another solution is to not worry about the fact that some variables are going to be duplicated for tracers when they don't need to be... but I think that's messy (for restart and other reasons).

@mvertens and I feel that the best solution is to move fluxes back into the centralized flux types, giving up on the idea of declaring flux variables in the science module responsible for computing the flux. I feel like there was some value in the idea of having the science modules own these fluxes, but not enough value to warrant the extra complexity it will cause in the face of water isotopes / tracers. I'd still like to keep some variables local to the science modules – especially variables that can be made private to those modules. But the fluxes were never handled in a completely (textbook) object-oriented fashion anyway, since they were accessed publicly by other modules, so I'm not too sad to move these public flux variables back out to the generic flux types. And while this has the disadvantage of losing the attachment of a variable to the module that computes that variable, it has some advantages, too, in terms of understanding the suite of water flux variables, and reducing the number of instances that need to get passed all around (e.g., we won't need to pass irrigation_inst to as many places).

[billsacks]

In a discussion with @mvertens and @mathewvrothstein we agreed that (D5) is probably best. I'll make another commit to this branch that prototypes that design, then @mathewvrothstein will move ahead with refactoring the code to implement it.

[billsacks]

@mathewvrothstein okay I have committed the prototype for (D5). This is ready for you to take over now, I think.

[billsacks]

@martynpclark @swensosc @dlawrenncar - We're planning to go with the design as illustrated by the latest commit in this PR (illustrated by design (D5) in the pdf linked from #395 (comment)). If you're interested in seeing what this looks like, easiest thing to do is probably to check out this branch and look at these files:

src/biogeophys/WaterBalanceType.F90
src/biogeophys/WaterDiagnosticBulkType.F90
src/biogeophys/WaterDiagnosticType.F90
src/biogeophys/WaterStateBulkType.F90
src/biogeophys/WaterStateType.F90
src/biogeophys/WaterType.F90

Let us know very soon if you have any problems with this, because @mathewvrothstein is planning to start implementing this soon.

[billsacks]

Closing this, since it has served its purpose.