Kd work and buoyancy flux diagnostics (#853)

* Adds new buoyancy flux diagnostics to diabatic_ALE

- Adds the ability to diagnose the implied work and buoyancy flux due to given diffusivities inside of diabatic_ALE
- Adds outputs of N2 after applying diffusion, the implied bouyancy flux, the integrated buoyancy flux, the salinity components of buoyancy flux and N2, the temperature components of buoyancy flux and N2, the ePBL contribution to the buoyancy flux, and double diffusion contribution to the buoyancy flux,
- Only considers the contribution of KPP to the total diffusive flux, the diagnostics should not be used w/ KPP until the non-local term and its contribution are properly accounted for.
- Only considers the contribution of convection to the total diffusive flux, since diagnose its implied work would require accounting for the convective contributions to diffusivity (it would be easily doable following this template in the future, but it is not done here).
- More work is needed to apply these methods to diffusivities that are diagnosed and stored in set_diffusivity.
- This code is not repeated for diabatic_ALE_legacy since the double diffusion component to diffusivity is applied separately, thus complicating the interpretation of these diagnostics in the legacy driver.

* Renaming diagnostics for buoyancy flux dz in diabatic_driver

* Fix bad logic (>0 checks) in MOM_diabatic_driver id_ tests

* Style updates in MOM_diagnose_KdWork

* Correct bug in previous style commit

* Atempt to fix unclear issue related to MOM_error in MOM_diagnose_KdWork

* Fix some potential diagnostic logic issue in MOM_diabatic_driver.

* Fix scaling and update units for Buoyancy Flux diagnostics

- Change buoyancy flux units to W/m3 (W/m2 for vertical integrated terms).
- Correct issues in unit scaling for Boussinesq and non-Boussinesq code.
- Compute N2 from specific volume for non-Boussinesq mode.

* Update units/scaling for non-Boussinesq in buoyancy flux diagnostic

* Reorder spatial do-loops to k,j,i in MOM_diagnose_KdWork for efficiency.

* Move Kd Work/buoyancy flux diagnostics

- Kd work diagnostic calculations migrated from MOM_diabatic_driver to MOM_diagnose_KdWork
- Kd work diagnostics stored in new type which contains copies of all requested diffusivities to facilitate storing and moving data around so it is accessible at the conclusion of diabatic_driver.
- All Kd arrays are allocatable that won't be allocated unless needed.
- 2d Vertical and volume integral quantities added.
- Many more Kd work diagnostics added including many from MOM_set_diffusivity, still not exhaustive.

* Undo reordering control structures in MOM_set_diffusivity

* Add KdWork diagnostics for legacy driver and add warning about not using in layer mode

* Fix KdWork issue in legacy driver and whitespace

* Reduce line lengths in MOM_diagnose_KdWork

* Missing factor of 0.5 in buoyancy flux

* Addresses failure when KdWork from ePBL diagnostic is on, but ePBL is off

* Remove white space around ' = ' for optional arguments

* Remove all array syntax addition and replace with explicit do-loops

* Fix typo mismatch of bdif_dz diagnostics

---------

Co-authored-by: brandon.reichl <Brandon.Reichl@noaa.gov>

Author: breichl