mantaray 0.2.0

Ray-tracing solver for ocean surface gravity waves that integrates the wave ray equations over spatially varying currents and bathymetry.
Documentation
{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c42959cc-dd0d-480f-9360-a53827641c5d",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import xarray as xr\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "import mantaray"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "53d17292-aaa2-4ee1-81f5-b00bc0deba51",
   "metadata": {},
   "outputs": [],
   "source": [
    "g = 9.8\n",
    "def period2wavenumber(T):\n",
    "    k = (2*np.pi)**2/g/T**2\n",
    "    return k\n",
    "\n",
    "def group_velocity(k):\n",
    "    cg = 0.5*(g/k)**.5\n",
    "    return cg"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ffdbcf96-597b-41ac-bcad-f2f65f187e28",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Period of incident waves in seconds\n",
    "T0 = 10\n",
    "# Direction of incident waves in radians (trig convention, going to)\n",
    "theta0 = 0\n",
    "# Convert period to wavenumber magnitude\n",
    "k0 = period2wavenumber(T0)\n",
    "# Calculate wavenumber components\n",
    "kx0 = k0*np.cos(theta0)\n",
    "ky0 = k0*np.sin(theta0)\n",
    "\n",
    "# Number of rays\n",
    "n_rays = 100\n",
    "# Initialize wavenumber for all rays\n",
    "Kx0 = kx0*np.ones(n_rays)\n",
    "Ky0 = ky0*np.ones(n_rays)\n",
    "\n",
    "# Current and bathymetry file path\n",
    "current = 'data/currents/zonal_jet.nc'\n",
    "bathymetry = 'data/bathymetry/bathy_zonal_jet.nc'\n",
    "\n",
    "# Read x and y from file to get domain size\n",
    "ds = xr.open_dataset(current)\n",
    "x = ds.x.values\n",
    "y = ds.y.values\n",
    "\n",
    "# Creates initial x position for all rays\n",
    "x0 = 10*np.ones(n_rays)\n",
    "y0 = np.linspace(0, y.max(), n_rays)\n",
    "\n",
    "# Estimates CFL\n",
    "# Computes grid smallest spacing\n",
    "dd = np.min([np.diff(x).mean(), np.diff(y).mean()])\n",
    "# Computes group velocity\n",
    "cg = group_velocity(k0)\n",
    "# Computes CFL\n",
    "cfl = dd/cg\n",
    "\n",
    "duration = round(x.max()/cg)\n",
    "step_size = cfl"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c928c5ad-edd8-46cc-91ae-733528ef22d2",
   "metadata": {},
   "outputs": [],
   "source": [
    "bando = mantaray.ray_tracing(x0, y0, Kx0, Ky0, duration, step_size, bathymetry, current)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "234b8f5f-04f4-4d00-8249-2011676d5efd",
   "metadata": {},
   "outputs": [],
   "source": [
    "bando"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "81546a85-5104-4e53-be88-6628ffbbfbd0",
   "metadata": {},
   "source": [
    "###  Plot rays and current"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "dfe5737d-ca97-4ee6-9117-c33797cfdcb6",
   "metadata": {},
   "outputs": [],
   "source": [
    "X = ds.x\n",
    "Y = ds.y\n",
    "U = (ds.u**2 + ds.v**2)**.5\n",
    "plt.figure(figsize=(12, 6))\n",
    "cs = plt.pcolormesh(X, Y, U)\n",
    "for i in range(bando.ray.size)[::2]:\n",
    "    ray = bando.isel(ray=i)\n",
    "    plt.plot(ray.x, ray.y, 'k', lw=.78)\n",
    "plt.colorbar(cs)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8e655c63-2016-4630-ae1b-50cbf8ac55bc",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.12.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}