Skip to content

Commit 1c47c67

Browse files
rajeejarajeeja
authored
Merge pull request #1387 from UXARRAY/rajeeja/vector_calc_div
Divergence
2 parents f7b3ad4 + 3013f29 commit 1c47c67

File tree

5 files changed

+643
-145
lines changed

5 files changed

+643
-145
lines changed

‎docs/user-guide/divergence.ipynb‎

Lines changed: 255 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,255 @@
1+
{
2+
"cells": [
3+
{
4+
"cell_type": "markdown",
5+
"id": "c51bdab8875e6859",
6+
"metadata": {},
7+
"source": [
8+
"# Divergence"
9+
]
10+
},
11+
{
12+
"cell_type": "markdown",
13+
"id": "table_of_contents",
14+
"metadata": {},
15+
"source": [
16+
"Computing the divergence of vector fields is a fundamental operation in vector calculus with applications in fluid dynamics, electromagnetism, and many other fields. This user-guide notebook showcases how to compute the divergence of a vector field."
17+
]
18+
},
19+
{
20+
"cell_type": "code",
21+
"execution_count": null,
22+
"id": "initial_id",
23+
"metadata": {},
24+
"outputs": [],
25+
"source": [
26+
"import numpy as np\n",
27+
"\n",
28+
"import uxarray as ux"
29+
]
30+
},
31+
{
32+
"cell_type": "markdown",
33+
"id": "divergence_computation",
34+
"metadata": {},
35+
"source": [
36+
"## Divergence Computation\n",
37+
"\n",
38+
"### Background\n",
39+
"\n",
40+
"The **divergence** of a vector field **V** = (u, v) is a scalar field that measures the \"outflow\" of the vector field from each point:\n",
41+
"\n",
42+
"$$\\nabla \\cdot \\mathbf{V} = \\frac{\\partial u}{\\partial x} + \\frac{\\partial v}{\\partial y}$$\n",
43+
"\n",
44+
"**Physical Interpretation:**\n",
45+
"- **Positive divergence**: Indicates a \"source\" - the vector field is flowing outward from that point\n",
46+
"- **Negative divergence**: Indicates a \"sink\" - the vector field is flowing inward to that point\n",
47+
"- **Zero divergence**: Indicates incompressible flow - no net outflow or inflow\n",
48+
"\n",
49+
"### Implementation\n",
50+
"\n",
51+
"In UXarray, the divergence is computed using the existing gradient infrastructure. The method leverages the finite-volume discretization to compute gradients of each vector component and then sums the relevant partial derivatives.\n",
52+
"\n",
53+
"| **Input** | **Usage** | **Output** |\n",
54+
"| ----------------------------- | :-----------------------------: | --------------------------- |\n",
55+
"| Vector field (u, v) | `u.divergence(v)` | Scalar field $\\nabla \\cdot \\mathbf{V}$ |"
56+
]
57+
},
58+
{
59+
"cell_type": "markdown",
60+
"id": "data_section",
61+
"metadata": {},
62+
"source": [
63+
"## Data\n",
64+
"\n",
65+
"This notebook uses a subset of a 30km MPAS atmosphere grid, taken centered at 45 degrees longitude and 0 degrees latitude with a radius of 2 degrees.\n",
66+
"- `face_lon`: Longitude at cell-centers\n",
67+
"- `face_lat`: Latitude at cell-centers\n",
68+
"- `gaussian`: Gaussian initialized at the center of the grid\n",
69+
"- `inverse_gaussian`: Inverse of the gaussian above."
70+
]
71+
},
72+
{
73+
"cell_type": "code",
74+
"execution_count": null,
75+
"id": "load_data",
76+
"metadata": {},
77+
"outputs": [],
78+
"source": [
79+
"base_path = \"../../test/meshfiles/mpas/dyamond-30km/\"\n",
80+
"grid_path = base_path + \"gradient_grid_subset.nc\"\n",
81+
"data_path = base_path + \"gradient_data_subset.nc\"\n",
82+
"\n",
83+
"uxds = ux.open_dataset(grid_path, data_path)\n",
84+
"print(f\"Grid has {uxds.uxgrid.n_face} faces\")\n",
85+
"print(f\"Available variables: {list(uxds.data_vars.keys())}\")\n",
86+
"uxds"
87+
]
88+
},
89+
{
90+
"cell_type": "markdown",
91+
"id": "usage_section",
92+
"metadata": {},
93+
"source": [
94+
"## Usage\n",
95+
"\n",
96+
"The divergence method is available on `UxDataArray` objects and follows this signature:\n",
97+
"\n",
98+
"```python\n",
99+
"div_result = u_component.divergence(v_component)\n",
100+
"```\n",
101+
"\n",
102+
"The method returns a `UxDataArray` containing the divergence values with the same shape and grid as the input components."
103+
]
104+
},
105+
{
106+
"cell_type": "markdown",
107+
"id": "gaussian_subsection",
108+
"metadata": {},
109+
"source": [
110+
"### Gaussian Fields"
111+
]
112+
},
113+
{
114+
"cell_type": "code",
115+
"execution_count": null,
116+
"id": "gaussian_example",
117+
"metadata": {},
118+
"outputs": [],
119+
"source": [
120+
"# Use Gaussian fields as vector components\n",
121+
"u_gauss = uxds[\"gaussian\"]\n",
122+
"v_gauss = uxds[\"inverse_gaussian\"]\n",
123+
"\n",
124+
"# Compute divergence\n",
125+
"div_gauss = u_gauss.divergence(v_gauss)\n",
126+
"\n",
127+
"# Handle NaN values from boundary faces\n",
128+
"finite_mask = np.isfinite(div_gauss.values)\n",
129+
"finite_values = div_gauss.values[finite_mask]\n",
130+
"\n",
131+
"print(f\"Total faces: {len(div_gauss.values)}\")\n",
132+
"print(f\"Interior faces: {len(finite_values)}\")\n",
133+
"print(f\"Boundary faces: {np.isnan(div_gauss.values).sum()}\")\n",
134+
"\n",
135+
"if len(finite_values) > 0:\n",
136+
" print(\n",
137+
" f\"Finite divergence range: [{finite_values.min():.6f}, {finite_values.max():.6f}]\"\n",
138+
" )\n",
139+
" print(f\"Mean divergence (finite): {finite_values.mean():.6f}\")"
140+
]
141+
},
142+
{
143+
"cell_type": "markdown",
144+
"id": "constant_subsection",
145+
"metadata": {},
146+
"source": [
147+
"### Constant Fields (Mathematical Validation)\n",
148+
"\n",
149+
"The divergence of a constant vector field should be zero everywhere (within numerical precision). This provides an important validation of our implementation."
150+
]
151+
},
152+
{
153+
"cell_type": "code",
154+
"execution_count": null,
155+
"id": "constant_example",
156+
"metadata": {},
157+
"outputs": [],
158+
"source": [
159+
"# The divergence of a constant vector field should be zero everywhere\n",
160+
"u_constant = uxds[\"face_lat\"] * 0 + 1.0 # Constant = 1\n",
161+
"v_constant = uxds[\"face_lat\"] * 0 + 2.0 # Constant = 2\n",
162+
"\n",
163+
"# Compute divergence\n",
164+
"div_constant = u_constant.divergence(v_constant)\n",
165+
"\n",
166+
"# Handle NaN values from boundary faces\n",
167+
"finite_mask = np.isfinite(div_constant.values)\n",
168+
"finite_values = div_constant.values[finite_mask]\n",
169+
"\n",
170+
"print(f\"Total faces: {len(div_constant.values)}\")\n",
171+
"print(f\"Finite values: {len(finite_values)}\")\n",
172+
"print(f\"NaN values (boundary faces): {np.isnan(div_constant.values).sum()}\")\n",
173+
"\n",
174+
"if len(finite_values) > 0:\n",
175+
" print(\n",
176+
" f\"Finite divergence range: [{finite_values.min():.10f}, {finite_values.max():.10f}]\"\n",
177+
" )\n",
178+
" print(f\"Maximum absolute divergence (finite): {np.abs(finite_values).max():.2e}\")\n",
179+
" print(f\"Mean absolute divergence (finite): {np.abs(finite_values).mean():.2e}\")\n",
180+
"\n",
181+
" # Should be close to zero for constant field\n",
182+
" max_abs_div = np.abs(finite_values).max()\n",
183+
" if max_abs_div < 1e-10:\n",
184+
" print(\"✓ Divergence is approximately zero as expected\")\n",
185+
" else:\n",
186+
" print(f\"⚠ Divergence is {max_abs_div:.2e} (may be due to discretization)\")"
187+
]
188+
},
189+
{
190+
"cell_type": "markdown",
191+
"id": "laplacian_computation",
192+
"metadata": {},
193+
"source": [
194+
"## Computing the Laplacian\n",
195+
"\n",
196+
"The Laplacian of a scalar field can be computed as the divergence of its gradient: $\\nabla^2 \\phi = \\nabla \\cdot (\\nabla \\phi)$\n",
197+
"\n",
198+
"This demonstrates the power of combining UXarray's vector calculus operations."
199+
]
200+
},
201+
{
202+
"cell_type": "code",
203+
"execution_count": null,
204+
"id": "laplacian_example",
205+
"metadata": {},
206+
"outputs": [],
207+
"source": [
208+
"# Compute gradient of the Gaussian field\n",
209+
"grad_gaussian = uxds[\"gaussian\"].gradient()\n",
210+
"\n",
211+
"# Extract gradient components\n",
212+
"grad_x = grad_gaussian[\"zonal_gradient\"]\n",
213+
"grad_y = grad_gaussian[\"meridional_gradient\"]\n",
214+
"\n",
215+
"# Compute Laplacian as divergence of gradient\n",
216+
"laplacian = grad_x.divergence(grad_y)\n",
217+
"\n",
218+
"# Analyze results\n",
219+
"finite_mask = np.isfinite(laplacian.values)\n",
220+
"finite_laplacian = laplacian.values[finite_mask]\n",
221+
"\n",
222+
"print(\"Laplacian computed successfully!\")\n",
223+
"print(f\"Interior faces: {len(finite_laplacian)}\")\n",
224+
"print(f\"Boundary faces: {np.isnan(laplacian.values).sum()}\")\n",
225+
"\n",
226+
"if len(finite_laplacian) > 0:\n",
227+
" print(\n",
228+
" f\"Laplacian range: [{finite_laplacian.min():.6f}, {finite_laplacian.max():.6f}]\"\n",
229+
" )\n",
230+
" print(f\"Mean Laplacian: {finite_laplacian.mean():.6f}\")"
231+
]
232+
}
233+
],
234+
"metadata": {
235+
"kernelspec": {
236+
"display_name": "Python 3 (ipykernel)",
237+
"language": "python",
238+
"name": "python3"
239+
},
240+
"language_info": {
241+
"codemirror_mode": {
242+
"name": "ipython",
243+
"version": 3
244+
},
245+
"file_extension": ".py",
246+
"mimetype": "text/x-python",
247+
"name": "python",
248+
"nbconvert_exporter": "python",
249+
"pygments_lexer": "ipython3",
250+
"version": "3.12.9"
251+
}
252+
},
253+
"nbformat": 4,
254+
"nbformat_minor": 5
255+
}

‎docs/userguide.rst‎

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -70,6 +70,9 @@ These user guides provide detailed explanations of the core functionality in UXa
7070
`Gradients <user-guide/gradients.ipynb>`_
7171
Compute the gradient of a data variable
7272

73+
`Divergence <user-guide/divergence.ipynb>`_
74+
Compute the divergence of a vector field
75+
7376
`Tree Structures <user-guide/tree_structures.ipynb>`_
7477
Data structures for nearest neighbor queries
7578

@@ -118,6 +121,7 @@ These user guides provide additional details about specific features in UXarray.
118121
user-guide/topological-aggregations.ipynb
119122
user-guide/weighted_mean.ipynb
120123
user-guide/gradients.ipynb
124+
user-guide/divergence.ipynb
121125
user-guide/tree_structures.ipynb
122126
user-guide/area_calc.ipynb
123127
user-guide/dual-mesh.ipynb

0 commit comments

Comments
 (0)