Add bib entry for Pauluis (2008) + few other fixes/code alignment (#51)

* add Pauluis2008 bib entry

* alignment

* lowercase

* use bib citation

* code alignment; order imports

* [ after (

* fix rendering

* code -> math

* fix typo in p(T)

* fix latex rendering and add definitions

Author: navidcy

docs/Project.toml

@@ -1,4 +1,5 @@
 [deps]
+Breeze = "660aa2fb-d4c8-4359-a52c-9c057bc511da"
 CairoMakie = "13f3f980-e62b-5c42-98c6-ff1f3baf88f0"
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
 DocumenterCitations = "daee34ce-89f3-4625-b898-19384cb65244"

docs/src/breeze.bib

@@ -9,11 +9,22 @@ @article{Chammas2023
     year = {2023}
 }
 
+@article{Pauluis2008,
+    author = {Olivier Pauluis},
+    title = {Thermodynamic consistency of the anelastic approximation for a moist atmosphere},
+    journal = {Journal of the Atmospheric Sciences},
+    year = {2008},
+    volume = {65},
+    number = {8},
+    doi = {10.1175/2007JAS2475.1},
+    pages = {2719-2729},
+}
+
 @book{Straka2009,
-    place={Cambridge},
-    title={Cloud and precipitation microphysics: Principles and parameterizations},
-    publisher={Cambridge University Press},
-    author={Straka, Jerry M.},
-    year={2009},
-    doi={10.1017/CBO9780511581168}
+    place = {Cambridge},
+    title = {Cloud and precipitation microphysics: Principles and parameterizations},
+    publisher = {Cambridge University Press},
+    author = {Straka, Jerry M.},
+    year = {2009},
+    doi = {10.1017/CBO9780511581168}
 }

docs/src/index.md

@@ -11,7 +11,7 @@ Breeze provides two ways to simulate moist atmospheres:
 
 * A [`MoistAirBuoyancy`](@ref) that can be used with [Oceananigans](https://clima.github.io/OceananigansDocumentation/stable/)' [`NonhydrostaticModel`](https://clima.github.io/OceananigansDocumentation/stable/appendix/library/#Oceananigans.Models.NonhydrostaticModels.NonhydrostaticModel-Tuple{}) to simulate atmospheric flows with the Boussinesq approximation.
 
-* A prototype [`AtmosphereModel`](@ref Breeze.AtmosphereModels.AtmosphereModel) that uses the anelastic approximation following [Pauluis (2008)](https://journals.ametsoc.org/view/journals/atsc/65/8/200s7jas2475.1.xml).
+* A prototype [`AtmosphereModel`](@ref Breeze.AtmosphereModels.AtmosphereModel) that uses the anelastic approximation following [Pauluis2008](@citet).
 
 ## Installation
 
@@ -26,7 +26,7 @@ Pkg.add("https://github.com/NumericalEarth/Breeze.jl.git")
 
 A basic free convection simulation:
 
-```@example intro
+```julia intro
 using Oceananigans
 using Oceananigans.Units
 using CairoMakie

src/AtmosphereModels/anelastic_formulation.jl

@@ -1,9 +1,3 @@
-using Oceananigans.Architectures: architecture
-using Oceananigans.Grids: inactive_cell
-using Oceananigans.Utils: prettysummary
-using Oceananigans.Operators: Δzᵃᵃᶜ, Δzᵃᵃᶠ, divᶜᶜᶜ
-using Oceananigans.Solvers: solve!
-
 using ..Thermodynamics:
     AtmosphereThermodynamics,
     ReferenceStateConstants,
@@ -13,6 +7,12 @@ using ..Thermodynamics:
     mixture_heat_capacity,
     dry_air_gas_constant
 
+using Oceananigans.Architectures: architecture
+using Oceananigans.Grids: inactive_cell
+using Oceananigans.Operators: Δzᵃᵃᶜ, Δzᵃᵃᶠ, divᶜᶜᶜ
+using Oceananigans.Solvers: solve!
+using Oceananigans.Utils: prettysummary
+
 using KernelAbstractions: @kernel, @index
 
 import Oceananigans.Solvers: tridiagonal_direction, compute_main_diagonal!, compute_lower_diagonal!

src/AtmosphereModels/atmosphere_model.jl

@@ -9,12 +9,12 @@ using ..Thermodynamics:
 
 using Oceananigans
 using Oceananigans.Advection: Centered, adapt_advection_order
-using Oceananigans.Grids: ZDirection
-using Oceananigans.Models: AbstractModel
 using Oceananigans.Architectures: AbstractArchitecture
-using Oceananigans.TimeSteppers: TimeStepper
 using Oceananigans.BoundaryConditions: FieldBoundaryConditions, regularize_field_boundary_conditions
+using Oceananigans.Grids: ZDirection
+using Oceananigans.Models: AbstractModel
 using Oceananigans.Solvers: FourierTridiagonalPoissonSolver, HomogeneousNeumannFormulation
+using Oceananigans.TimeSteppers: TimeStepper
 using Oceananigans.Utils: launch!
 
 using KernelAbstractions: @kernel, @index
@@ -81,10 +81,11 @@ end
                     timestepper = :RungeKutta3)
 
 Return an AtmosphereModel that uses the anelastic approximation following
-[Pauluis (2008)](https://journals.ametsoc.org/view/journals/atsc/65/8/200s7jas2475.1.xml).
+[Pauluis2008](@citet).
 
 Example
 =======
+
 ```jldoctest
 julia> using Breeze, Breeze.AtmosphereModels, Oceananigans
 
@@ -100,6 +101,11 @@ AtmosphereModel{CPU, RectilinearGrid}(time = 0 seconds, iteration = 0)
 ├── coriolis: Nothing
 └── microphysics: WarmPhaseSaturationAdjustment
 ```
+
+References
+==========
+Pauluis, O. (2008). Thermodynamic consistency of the anelastic approximation for a moist atmosphere.
+  Journal of the Atmospheric Sciences 65, 2719–2729.
 """
 function AtmosphereModel(grid;
                          clock = Clock(grid),

src/AtmosphereModels/set_atmosphere_model.jl

@@ -1,7 +1,8 @@
-import Oceananigans.Fields: set!
-using Oceananigans.TimeSteppers: update_state!
 using Oceananigans.BoundaryConditions: fill_halo_regions!
 using Oceananigans.Models.NonhydrostaticModels: compute_pressure_correction!, make_pressure_correction!
+using Oceananigans.TimeSteppers: update_state!
+
+import Oceananigans.Fields: set!
 
 function set!(model::AtmosphereModel; enforce_mass_conservation=true, kw...)
     for (name, value) in kw
@@ -43,14 +44,14 @@ function set!(model::AtmosphereModel; enforce_mass_conservation=true, kw...)
             value = ρʳ * u
         end
 
-        set!(ϕ, value)                
+        set!(ϕ, value)
         fill_halo_regions!(ϕ, model.clock, fields(model))
     end
 
     # Apply a mask
     # foreach(mask_immersed_field!, prognostic_fields(model))
     update_state!(model, compute_tendencies=false)
-    
+
     if enforce_mass_conservation
         FT = eltype(model.grid)
         Δt = one(FT)

src/AtmosphereModels/update_atmosphere_model_state.jl

@@ -3,12 +3,12 @@ using ..Thermodynamics:
     mixture_heat_capacity,
     mixture_gas_constant
 
+using Oceananigans.Architectures: architecture
 using Oceananigans.BoundaryConditions: fill_halo_regions!, compute_x_bcs!, compute_y_bcs!, compute_z_bcs!
 using Oceananigans.ImmersedBoundaries: mask_immersed_field!
-using Oceananigans.Architectures: architecture
 
-import Oceananigans.TimeSteppers: update_state!, compute_flux_bc_tendencies!
 import Oceananigans: fields, prognostic_fields
+import Oceananigans.TimeSteppers: update_state!, compute_flux_bc_tendencies!
 
 const AnelasticModel = AtmosphereModel{<:AnelasticFormulation}
 
@@ -113,7 +113,7 @@ function compute_tendencies!(model::AnelasticModel)
                    model.momentum,
                    model.coriolis,
                    model.clock,
-                   fields(model))    
+                   fields(model))
 
     pₕ′ = model.hydrostatic_pressure_anomaly
     ρᵣ = model.formulation.reference_density
@@ -178,7 +178,7 @@ end
                                      hydrostatic_pressure_anomaly)
 
     # Note: independent of x
-    ρᵣ = @inbounds reference_density[i, j, k]    
+    ρᵣ = @inbounds reference_density[i, j, k]
 
     return ( - div_𝐯u(i, j, k, grid, advection, velocities, momentum.ρu)
              - x_f_cross_U(i, j, k, grid, coriolis, momentum)
@@ -198,7 +198,7 @@ end
                                      hydrostatic_pressure_anomaly)
 
     # Note: independent of y
-    ρᵣ = @inbounds reference_density[i, j, k]    
+    ρᵣ = @inbounds reference_density[i, j, k]
 
     return ( - div_𝐯v(i, j, k, grid, advection, velocities, momentum.ρv)
              - y_f_cross_U(i, j, k, grid, coriolis, momentum)
@@ -222,7 +222,7 @@ end
 
     ρᵣᶜᶜᶠ = ℑzᵃᵃᶠ(i, j, k, grid, reference_density)
     bᶜᶜᶠ = ℑzᵃᵃᶠ(i, j, k, grid, buoyancy,
-                 formulation, temperature, specific_humidity, thermo)    
+                 formulation, temperature, specific_humidity, thermo)
 
     return ( - div_𝐯w(i, j, k, grid, advection, velocities, momentum.ρw)
              - z_f_cross_U(i, j, k, grid, coriolis, momentum)
@@ -259,10 +259,10 @@ end
              + forcing(i, j, k, grid, clock, model_fields))
 end
 =#
-                                        
+
 """ Apply boundary conditions by adding flux divergences to the right-hand-side. """
 function compute_flux_bc_tendencies!(model::AtmosphereModel)
-    
+
     Gⁿ    = model.timestepper.Gⁿ
     arch  = model.architecture
     clock = model.clock

src/AtmosphereModels/update_hydrostatic_pressure.jl

@@ -2,10 +2,9 @@
 ##### update pressure
 #####
 
-using Oceananigans.Operators: Δzᶜᶜᶜ, Δzᶜᶜᶠ, ℑzᵃᵃᶠ
+using Oceananigans.Grids: topology, XFlatGrid, YFlatGrid
 using Oceananigans.ImmersedBoundaries: PartialCellBottom, ImmersedBoundaryGrid
-using Oceananigans.Grids: topology
-using Oceananigans.Grids: XFlatGrid, YFlatGrid
+using Oceananigans.Operators: Δzᶜᶜᶜ, Δzᶜᶜᶠ, ℑzᵃᵃᶠ
 using Oceananigans.Utils: KernelParameters
 
 const c = Center()

src/Thermodynamics/Thermodynamics.jl

@@ -1,6 +1,8 @@
 module Thermodynamics
 
-export AtmosphereThermodynamics, ReferenceStateConstants, IdealGas, PhaseTransitionConstants, CondensedPhase, mixture_gas_constant, mixture_heat_capacity
+export AtmosphereThermodynamics, ReferenceStateConstants, IdealGas,
+       PhaseTransitionConstants, CondensedPhase,
+       mixture_gas_constant, mixture_heat_capacity
 
 include("atmosphere_thermodynamics.jl")
 include("vapor_saturation.jl")

src/Thermodynamics/atmosphere_thermodynamics.jl

@@ -160,27 +160,30 @@ end
 
 Create `AtmosphereThermodynamics` with parameters that represent gaseous mixture of dry "air"
 and vapor, as well as condensed liquid and solid phases.
-The `triple_point_temperature` and `triple_point_pressure` may be combined with 
+The `triple_point_temperature` and `triple_point_pressure` may be combined with
 internal energy parameters for condensed phases to compute the vapor pressure
 at the boundary between vapor and a homogeneous sample of the condensed phase.
 The `gravitational_acceleration` parameter is included to compute reference_state
 quantities associated with hydrostatic balance.
 
-The Clausius-Clapeyron relation describes the pressure-temperature relationship during phase transitions:
+The Clausius-Clapeyron relation describes the pressure-temperature relationship during phase
+transitions:
 
-    d/dT log(pⁱ⁺) = ℒⁱ / (Rⁱ * T²)
+```math
+d[\\log(pⁱ⁺)] / dT = ℒⁱ / (Rⁱ T²)
+```
 
 where:
 
-- `pⁱ⁺` is the saturation vapor pressure for a transition between vapor and the `ⁱ`th phase
-- `T` is temperature
-- `ℒⁱ` is the latent heat of vaporization
-- `Rⁱ` is the specific gas constant for the `ⁱ`th phase
+- ``pⁱ⁺`` is the saturation vapor pressure for a transition between vapor and the ``i``-th phase
+- ``T`` is temperature
+- ``ℒⁱ`` is the latent heat of vaporization
+- ``Rⁱ`` is the specific gas constant for the ``i``-th phase
 
 For water vapor, this integrates to:
 
 ```math
-    pⁱ⁺ = pᵗʳ * exp( ℒⁱ (1/Tᵗʳ - 1/T) / Rⁱ )
+pⁱ⁺ = pᵗʳ \\exp[ ℒⁱ (1/Tᵗʳ - 1/T) / Rⁱ ]
 ```
 
 where
@@ -191,7 +194,6 @@ where
 Note: any reference values for pressure and temperature can be used in principle.
 The advantage of using reference values at the triple point is that the same values
 can then be used for both condensation (vapor → liquid) and deposition (vapor → ice).
-
 """
 function AtmosphereThermodynamics(FT = Oceananigans.defaults.FloatType;
                                   molar_gas_constant = 8.314462618,
@@ -234,7 +236,7 @@ const NonCondensingAtmosphereThermodynamics{FT} = AtmosphereThermodynamics{FT, N
 """
     mixture_gas_constant(q, thermo)
 
-Compute the gas constant of moist air given the specific humidity `q` and 
+Compute the gas constant of moist air given the specific humidity `q` and
 thermodynamic parameters `thermo`.
 
 The mixture gas constant is calculated as a weighted average of the dry air
@@ -246,7 +248,7 @@ R_m = R_d (1 - q) + R_v q
 
 where:
 - `R_d` is the dry air gas constant
-- `R_v` is the water vapor gas constant  
+- `R_v` is the water vapor gas constant
 - `q` is the specific humidity (mass fraction of water vapor)
 
 # Arguments