(0.95.6) Correctly account for field names in ContinuousBoundaryFunctions

[simone-silvestri]

continuous boundary functions with field dependencies are broken at the moment because of the recent changes with split explicit that changed the order of the fields outputted by fields(model).

This PR corrects this bug and introduces
(1) some tests
(2) a field_names function to be able to change the field names in a HydrostaticFreeSurfaceModel depending on different inputs

[simone-silvestri]

I am a bit confused about this function. I think that the biogeochemical auxiliary fields have been added to the auxiliary_fields during model construction, so we shouldn't need an extra biogeochemical_auxiliary_fields(model.biogeochemistry) here right? @jagoosw is it true?

[glwagner]

no, the biogeochemistry struct is allowed to have extra auxiliary fields that are not attached to the model. Basically model.auxiliary_fields are purely for the user, whereas biogeochemistry can have its own auxiliary fields, much like the closure has auxiliary fields that are stored in diffusivity_fields (which could be called closure_auxiliary_fields).

[simone-silvestri]

Then I guess the implementation is bugged because It actually looks like the fields are added to the auxiliary_fields, and, thus they are accounted for twice in the fields function

Oceananigans.jl/src/Models/HydrostaticFreeSurfaceModels/hydrostatic_free_surface_model.jl

Lines 131 to 133 in d8273e7

	biogeochemical_fields = merge(auxiliary_fields, biogeochemical_auxiliary_fields(biogeochemistry))
	tracers, auxiliary_fields = validate_biogeochemistry(tracers, biogeochemical_fields, biogeochemistry, grid, clock)

[glwagner]

I think this is a bad idea. I think that fields and prognostic_fields should be the source of truth. Those functions return NamedTuple --- both names, and fields. Thus the names can be derived from fields or prognostic_fields.

[simone-silvestri]

we need something to regularize the boundary conditions, so we need some names that correspond to keys(fields(model)) before the model is built. I am open to a different implementation though.
At the moment there is a (wrong) assumption that fields(model) will always have the same structure. We could change fields(model) to return nothing in the positions where we would have a field, in that way the regularization would work

[glwagner]

I see

[simone-silvestri]

I have changed the design for the fields(model) to always output a predictable namedtuple with nothing in place of fields that don't exist. In this way we don't need a new function

[glwagner]

(2) a field_names function to be able to change the field names in a HydrostaticFreeSurfaceModel depending on different inputs

I think it is a bad design to have multiplicitous sources of truth. This PR does not solve the underlying issue, which is that there are two conflicting sources of truth. Instead, it creates a third source of truth which is definitely not learning the lesson that the bug should have given.

[glwagner]

Get rid of field_names, and instead create a design where there is only one source of truth about the field names --- either from fields or prognostic_fields. We may not need fields independently from prognostic_fields.

[glwagner]

Can you add a test that exposes the issue with ContinuousBoundaryFunction directly?

[simone-silvestri]

done in a different way

[glwagner]

I see the issue now that the source of truth was hardcoded int he constructor... hmm..

[glwagner]

We want this for the test?

[simone-silvestri]

It seems that I^{PAR} is repeated in fields. So, indeed, it looks like the bialogical_auxiliary_fields are added to the model's auxiliary fields upon model construction here:

    biogeochemical_fields = merge(auxiliary_fields, biogeochemical_auxiliary_fields(biogeochemistry))
    tracers, auxiliary_fields = validate_biogeochemistry(tracers, biogeochemical_fields, biogeochemistry, grid, clock)

which leads to having repeated I^{PAR} here

@inline fields(model::HydrostaticFreeSurfaceModel) = 
    merge(hydrostatic_fields(model.velocities, model.free_surface, model.tracers), 
          model.auxiliary_fields,
          biogeochemical_auxiliary_fields(model.biogeochemistry))

So I think we have a decision to make
(1) leave the biogeochemical fields only in the model.biogeochemistry
(2) copy them in model.auxiliary_fields but modify the fields function

I guess from previous comments, the intended behavior is option (1), so I will remove the merging of the biogeochemical auxiliaries to the model auxiliaries at model construction.

[glwagner]

required_biogeochemical_auxiliary_fields are the model.auxiliary_fields. Then there are also biogeochemical_auxiliary_fields(model.biogeochemistry), which are distinct? Notice that model.auxiliary_fields are not accessible from biogeochemical_auxiliary_fields(model.biogeochemistry). So it seems like there are two different kinds of auxiliary fields with bgc.

We can change the design.