0
# Geostrophic Functions
1

2
Ocean circulation and dynamic oceanography calculations for analyzing geostrophic currents, ocean stability, and potential vorticity. These functions are essential for understanding ocean dynamics and circulation patterns.
3

4
## Capabilities
5

6
### Brünt-Väisälä Frequency
7

8
Ocean stability analysis through buoyancy frequency calculations.
9

10
```python { .api }
11
def bfrq(s, t, p, lat=None):
12
    """
13
    Brünt-Väisälä frequency squared and potential vorticity.
14
    
15
    Calculates the buoyancy frequency (also known as the Brünt-Väisälä frequency)
16
    which indicates ocean stability and stratification strength.
17
    
18
    Parameters:
19
    - s: array_like, salinity profile [psu (PSS-78)]
20
    - t: array_like, temperature profile [°C (ITS-90)]
21
    - p: array_like, pressure profile [db]
22
    - lat: array_like, latitude [decimal degrees], optional
23
    
24
    Returns:
25
    - tuple: (n2, q, p_ave)
26
        - n2: array_like, Brünt-Väisälä frequency squared [rad²/s²]
27
        - q: array_like, potential vorticity [(m·s)⁻¹]
28
        - p_ave: array_like, mid-point pressures [db]
29
    """
30
```
31

32
### Specific Volume Anomaly
33

34
Specific volume calculations for geostrophic velocity analysis.
35

36
```python { .api }
37
def svan(s, t, p=0):
38
    """
39
    Specific volume anomaly of seawater.
40
    
41
    Calculates the specific volume anomaly which is used in geostrophic
42
    velocity calculations and dynamic height computations.
43
    
44
    Parameters:
45
    - s: array_like, salinity [psu (PSS-78)]
46
    - t: array_like, temperature [°C (ITS-90)]
47
    - p: array_like, pressure [db], default 0
48
    
49
    Returns:
50
    - array_like: specific volume anomaly [m³/kg]
51
    """
52
```
53

54
### Geopotential Anomaly
55

56
Geopotential calculations for dynamic oceanography.
57

58
```python { .api }
59
def gpan(s, t, p):
60
    """
61
    Geopotential anomaly of seawater.
62
    
63
    Calculates geopotential anomaly which represents the work done against
64
    gravity in lifting a unit mass from a reference level.
65
    
66
    Parameters:
67
    - s: array_like, salinity [psu (PSS-78)]
68
    - t: array_like, temperature [°C (ITS-90)]
69
    - p: array_like, pressure [db]
70
    
71
    Returns:
72
    - array_like: geopotential anomaly [m²/s²]
73
    """
74
```
75

76
### Geostrophic Velocity
77

78
Geostrophic current calculations from geopotential gradients.
79

80
```python { .api }
81
def gvel(ga, lat, lon):
82
    """
83
    Geostrophic velocity from geopotential anomaly gradients.
84
    
85
    Calculates geostrophic velocities from horizontal gradients of
86
    geopotential anomaly, accounting for Earth's rotation (Coriolis effect).
87
    
88
    Parameters:
89
    - ga: array_like, geopotential anomaly [m²/s²]
90
    - lat: array_like, latitude [decimal degrees]
91
    - lon: array_like, longitude [decimal degrees]
92
    
93
    Returns:
94
    - array_like: geostrophic velocity [m/s]
95
    """
96
```
97

98
## Usage Examples
99

100
### Ocean Stability Analysis
101

102
```python
103
import seawater as sw
104
import numpy as np
105
import matplotlib.pyplot as plt
106

107
# Create a stratified ocean profile
108
pressure = np.arange(0, 1000, 10)  # 0 to 1000 db
109
# Typical thermocline structure
110
temperature = 20 * np.exp(-pressure/200) + 2
111
salinity = 34.5 + 0.5 * (pressure/1000)
112
latitude = 30.0  # degrees N
113

114
# Calculate Brünt-Väisälä frequency
115
n2, q, p_mid = sw.bfrq(salinity, temperature, pressure, latitude)
116

117
# Convert to buoyancy frequency in cycles per hour
118
n_cph = np.sqrt(np.maximum(n2, 0)) * 3600 / (2 * np.pi)
119

120
# Plot stratification
121
plt.figure(figsize=(10, 6))
122
plt.subplot(1, 2, 1)
123
plt.plot(temperature, pressure, 'b-', label='Temperature')
124
plt.plot(salinity, pressure, 'r-', label='Salinity')
125
plt.gca().invert_yaxis()
126
plt.xlabel('Temperature (°C) / Salinity (psu)')
127
plt.ylabel('Pressure (db)')
128
plt.legend()
129
plt.title('Ocean Profile')
130

131
plt.subplot(1, 2, 2)
132
plt.plot(n_cph, p_mid, 'g-', linewidth=2)
133
plt.gca().invert_yaxis()
134
plt.xlabel('Buoyancy Frequency (cph)')
135
plt.ylabel('Pressure (db)')
136
plt.title('Ocean Stability')
137
plt.grid(True)
138

139
plt.tight_layout()
140
plt.show()
141
```
142

143
### Geostrophic Current Calculation
144

145
```python
146
import seawater as sw
147
import numpy as np
148

149
# Create two ocean stations for geostrophic calculation
150
pressure = np.arange(0, 500, 25)  # pressure levels
151

152
# Station 1 (warmer, less saline)
153
s1 = np.full_like(pressure, 34.8) - 0.001 * pressure
154
t1 = 15.0 * np.exp(-pressure/150) + 4
155

156
# Station 2 (cooler, more saline)  
157
s2 = np.full_like(pressure, 35.0) - 0.0005 * pressure
158
t2 = 12.0 * np.exp(-pressure/180) + 3
159

160
# Calculate geopotential anomaly at each station
161
gpa1 = sw.gpan(s1, t1, pressure)
162
gpa2 = sw.gpan(s2, t2, pressure)
163

164
# Calculate geopotential difference (proxy for geostrophic velocity)
165
dgpa = gpa2 - gpa1
166
print("Geopotential anomaly difference:")
167
for i in range(0, len(pressure), 4):
168
    print(f"P: {pressure[i]:3.0f} db, Δφ: {dgpa[i]:8.3f} m²/s²")
169

170
# Calculate specific volume anomaly
171
sva1 = sw.svan(s1, t1, pressure)
172
sva2 = sw.svan(s2, t2, pressure)
173

174
print(f"\\nSpecific volume anomaly range:")
175
print(f"Station 1: {sva1.min():.2e} to {sva1.max():.2e} m³/kg")
176
print(f"Station 2: {sva2.min():.2e} to {sva2.max():.2e} m³/kg")
177
```
178

179
### Dynamic Height Analysis
180

181
```python
182
import seawater as sw
183
import numpy as np
184

185
# Typical subtropical ocean profile
186
depth = np.arange(0, 2000, 50)
187
lat = 25.0  # degrees N
188

189
# Convert depth to pressure
190
pressure = sw.pres(depth, lat)
191

192
# Create realistic T-S profile
193
temperature = 22 * np.exp(-depth/300) + 2.5
194
salinity = 36.0 - 1.0 * np.exp(-depth/500)
195

196
# Calculate dynamic properties
197
n2, q, p_mid = sw.bfrq(salinity, temperature, pressure, lat)
198
sva = sw.svan(salinity, temperature, pressure)
199
gpa = sw.gpan(salinity, temperature, pressure)
200

201
# Calculate dynamic height (integrated specific volume anomaly)
202
# This is a simplified calculation - proper dynamic height requires
203
# integration between pressure levels
204
dyn_height = np.cumsum(sva * np.gradient(pressure)) / 10  # convert to dynamic meters
205

206
print("Dynamic oceanography summary:")
207
print(f"Maximum N²: {np.max(n2):.2e} rad²/s²")
208
print(f"Surface specific volume anomaly: {sva[0]:.2e} m³/kg")
209
print(f"Deep specific volume anomaly: {sva[-1]:.2e} m³/kg")
210
print(f"Total dynamic height: {dyn_height[-1]:.3f} dynamic meters")
211

212
# Find thermocline depth (maximum stratification)
213
max_strat_idx = np.argmax(n2)
214
thermocline_depth = depth[max_strat_idx//2]  # approximate due to mid-point pressures
215
print(f"Thermocline depth: {thermocline_depth:.0f} m")
216
```