*+Add ePBL bottom boundary layer mixing option

[Hallberg-NOAA]

This series of two commits adds the option to do energetically consistent bottom boundary layer mixing with
the new routine ePBL_BBL_column(). ePBL_BBL_column() is closely based on the surface-focused ePBL mixing in ePBL_column(), but without adding convective instability driven mixing or mean-TKE driven mixing to avoid possible
double-counting. This new option is enabled by setting the new runtime parameter EPBL_BBL_EFFIC to be positive. The changes that actually adds ePBL_BBL_column is in the second of the two commits.

The first commit corrects a single-line bug that causes ePBL_column() to set the wrong variable and then use an uninitialized variable when EPBL_ORIGINAL_PE_CALC is set to false. This bug was present when the EPBL_ORIGINAL_PE_CALC was first added on Sept. 30, 2016, but it was not detected because only the default case with EPBL_ORIGINAL_PE_CALC = True appears to being used or tested. Any runs that used this code with debugging compile options would have trapped it immediately. This will change answers when EPBL_ORIGINAL_PE_CALC is false.

If both EPBL_BBL_EFFIC and BBL_EFFIC are set to positive values, there is a risk of double-counting, but this case is not being trapped for now.

Within ePBL_column(), the hp_a variable was converted from a scalar into a 1-d array to emphasize where these variables are located and to bring ePBL_column() closer to what would be required if the new diffusivities were being applied on top of diffusivities from other processes, as is the case in ePBL_BBL_column(). Several comments in find_PE_change() and ePBL_column() were corrected. Two allocatable three-dimensional diagnostic arrays that had been in the energetic_PBL control structure were converted into local arrays to reduce the dedicated memory requirements when these diagnostics are enabled.

The changes include the addition of a new mandatory vertvisc_type argument to the publicly visible routine energetic_PBL().

When this new ePBL bottom boundary layer mixing option is enabled, there are several new diagnostics available that are related to bottom boundary layer mixing. Several new checksum calls were also added with this new option when DEBUG = True. The MOM_parameter_doc files are altered by the addition of two new runtime parameters, and by the correction of several spelling errors in the descriptions of other ePBL parameters. By default, all answers are bitwise identical.

This commits in this PR include:

• 3034aed0b +(*)Add ePBL bottom boundary mixing option
• 8a2a7c99f *Fix a bug when EPBL_ORIGINAL_PE_CALC is false

[Hallberg-NOAA]

I am converting this to a draft PR , because while I believe that it is correct (and worth assessing), I think that it could be made much more efficient. We are not using released mean kinetic energy to drive added bottom boundary layer mixing, which means that we can solve directly for the diffusivity that can be supported by a given energy increment without having to do so iteratively with Newton's method.