Fix flat face coloring and add particle interpolation test section

Still need to clean up the particle color plots. Almost there..

Author: fluidnumericsJoe

docs/user_guide/examples/tutorial_interpolation.ipynb

@@ -294,6 +294,70 @@
     "fieldset = parcels.FieldSet([U, V, UV, Tnode, Tface])"
    ]
   },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "We can now create particles to sample data using either the node registered or face registered data."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "SampleParticle = parcels.Particle.add_variable(\n",
+    "    parcels.Variable(\"Tracer\", dtype=np.float32, initial=np.nan)\n",
+    ")\n",
+    "\n",
+    "\n",
+    "def SampleTracer_Node(particles, fieldset):\n",
+    "    particles.Tracer = fieldset.T_node[particles]\n",
+    "\n",
+    "\n",
+    "def SampleTracer_Face(particles, fieldset):\n",
+    "    particles.Tracer = fieldset.T_face[particles]"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": []
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "print(fieldset.fields)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "pset = {}\n",
+    "\n",
+    "xv, yv = np.meshgrid(np.linspace(0.2, 0.8, 8), np.linspace(0.2, 0.9, 8))\n",
+    "pset[\"node\"] = parcels.ParticleSet(\n",
+    "    fieldset, pclass=SampleParticle, lon=xv.flatten(), lat=yv.flatten()\n",
+    ")\n",
+    "pset[\"node\"].execute(\n",
+    "    SampleTracer_Node, runtime=np.timedelta64(1, \"D\"), dt=np.timedelta64(1, \"D\")\n",
+    ")\n",
+    "\n",
+    "pset[\"face\"] = parcels.ParticleSet(\n",
+    "    fieldset, pclass=SampleParticle, lon=xv.flatten(), lat=yv.flatten()\n",
+    ")\n",
+    "pset[\"face\"].execute(\n",
+    "    SampleTracer_Face, runtime=np.timedelta64(1, \"D\"), dt=np.timedelta64(1, \"D\")\n",
+    ")"
+   ]
+  },
   {
    "cell_type": "code",
    "execution_count": null,
@@ -330,6 +394,17 @@
     "    x, y, c=T_node, cmap=\"viridis\", s=150, edgecolors=\"black\", vmin=-1.0, vmax=1.0\n",
     ")\n",
     "\n",
+    "ax[0].scatter(\n",
+    "    pset[\"node\"].lon,\n",
+    "    pset[\"node\"].lat,\n",
+    "    c=pset[\"node\"].Tracer,\n",
+    "    cmap=\"viridis\",\n",
+    "    edgecolors=\"k\",\n",
+    "    s=50,\n",
+    "    vmin=-1.0,\n",
+    "    vmax=1.0,\n",
+    ")\n",
+    "\n",
     "cbar = fig.colorbar(sc1, ax=ax[0])\n",
     "cbar.set_label(\"T\")\n",
     "\n",
@@ -340,9 +415,10 @@
     "\n",
     "# # Plot the face registered data\n",
     "T_face = np.squeeze(ds[\"T_face\"][0, 0, :].values)\n",
-    "tpc = ax[1].tripcolor(\n",
+    "assert triang.triangles.shape[0] == T_face.shape[0]\n",
+    "tpc2 = ax[1].tripcolor(\n",
     "    triang,\n",
-    "    T_face,\n",
+    "    facecolors=T_face,\n",
     "    shading=\"flat\",\n",
     "    cmap=\"viridis\",\n",
     "    edgecolors=\"k\",\n",
@@ -353,14 +429,22 @@
     "# Plot the element edges\n",
     "ax[1].triplot(triang, color=\"black\", linewidth=0.5)\n",
     "\n",
-    "# Plot the node registered data\n",
+    "# Plot the face registered data\n",
     "xf = ds.uxgrid.face_lon.values\n",
     "yf = ds.uxgrid.face_lat.values\n",
-    "sc2 = ax[1].scatter(\n",
-    "    xf, yf, c=T_face, cmap=\"viridis\", s=150, edgecolors=\"black\", vmin=-1.0, vmax=1.0\n",
+    "\n",
+    "ax[1].scatter(\n",
+    "    pset[\"face\"].lon,\n",
+    "    pset[\"face\"].lat,\n",
+    "    c=pset[\"face\"].Tracer,\n",
+    "    cmap=\"viridis\",\n",
+    "    edgecolors=\"k\",\n",
+    "    s=50,\n",
+    "    vmin=-1.0,\n",
+    "    vmax=1.0,\n",
     ")\n",
     "\n",
-    "cbar = fig.colorbar(sc2, ax=ax[1])\n",
+    "cbar = fig.colorbar(tpc, ax=ax[1])\n",
     "cbar.set_label(\"T\")\n",
     "\n",
     "ax[1].set_aspect(\"equal\", \"box\")\n",

src/parcels/_datasets/unstructured/generated.py

@@ -43,14 +43,14 @@ def simple_small_delaunay(nx=10, ny=10):
     Tface = np.zeros((1, nz1, uxgrid.n_face), dtype=np.float64)
 
     for i, (x, y) in enumerate(zip(uxgrid.face_lon, uxgrid.face_lat, strict=False)):
-        P[0, 0, i] = -vmax * delta * (1 - x) * (math.exp(-x / delta) - 1) * np.sin(math.pi * y)
-        U[0, 0, i] = -vmax * (1 - math.exp(-x / delta) - x) * np.cos(math.pi * y)
-        V[0, 0, i] = vmax * ((2.0 - x) * math.exp(-x / delta) - 1) * np.sin(math.pi * y)
-        Tface[0, 0, i] = np.sin(math.pi * y) * np.cos(math.pi * x)
+        P[0, :, i] = -vmax * delta * (1 - x) * (math.exp(-x / delta) - 1) * np.sin(math.pi * y)
+        U[0, :, i] = -vmax * (1 - math.exp(-x / delta) - x) * np.cos(math.pi * y)
+        V[0, :, i] = vmax * ((2.0 - x) * math.exp(-x / delta) - 1) * np.sin(math.pi * y)
+        Tface[0, :, i] = np.sin(math.pi * y) * np.cos(math.pi * x)
 
-    Tnode = np.zeros((1, nz1, uxgrid.n_node), dtype=np.float64)
+    Tnode = np.zeros((1, nz, uxgrid.n_node), dtype=np.float64)
     for i, (x, y) in enumerate(zip(uxgrid.node_lon, uxgrid.node_lat, strict=False)):
-        Tnode[0, 0, i] = np.sin(math.pi * y) * np.cos(math.pi * x)
+        Tnode[0, :, i] = np.sin(math.pi * y) * np.cos(math.pi * x)
 
     u = ux.UxDataArray(
         data=U,
@@ -124,10 +124,10 @@ def simple_small_delaunay(nx=10, ny=10):
         data=Tnode,
         name="T_node",
         uxgrid=uxgrid,
-        dims=["time", "nz1", "n_node"],
+        dims=["time", "nz", "n_node"],
         coords=dict(
             time=(["time"], [TIME[0]]),
-            nz1=(["nz1"], zc),
+            nz1=(["nz"], zf),
         ),
         attrs=dict(
             description="Tracer field sampled on face vertices",