Add WEC-Sim power performance example

examples/data/wave/pacwave_cluster_32.csv

@@ -0,0 +1,33 @@
+,Te,Hm0,weights,Tp,J
+0,7.974491297618696,1.2539695860020375,0.05886124580653463,9.294278901653492,6031.115427233068
+1,10.79453302186441,2.6414034469857426,0.03521606218803187,12.581040818023789,37004.325097576424
+2,6.901978594273449,1.9531216543907912,0.05200095798403432,8.044264095889801,12516.214518733921
+3,12.66762764506827,7.310116435425513,0.00506999304134946,14.764134784461854,367451.94558109366
+4,12.893701222644763,2.262294477484095,0.016045936894736042,15.027623802616274,36455.136139003094
+5,10.557620939325757,4.754296605622727,0.017311356620466345,12.30491950970368,116784.36178931354
+6,8.766663943128838,2.73937982438113,0.04364638972455223,10.217557043273704,31825.667988697423
+7,6.53740270831686,1.3055783411566657,0.05074593316072677,7.619350475893777,5272.441394386226
+8,9.666290858496115,1.340694087456887,0.03707009242331086,11.266073261650485,8443.64982080737
+9,12.787307070448522,3.9203973577543567,0.016464374526706786,14.903621294229048,107680.98651108926
+10,11.605879274797648,1.8210155555767638,0.022322848153721837,13.52666582144248,19446.169168239652
+11,7.584081846459788,1.8787353416133314,0.05462400212242149,8.839256231305114,12827.143860504848
+12,10.175410691413894,6.133932424424177,0.008836092568956415,11.859453020295915,188189.68918703857
+13,9.31935707393927,4.587432323814967,0.018817182816956667,10.861721531397752,95155.06836750833
+14,7.228995507990213,1.2566906096446477,0.05769191323850427,8.425402689965283,5449.034308785218
+15,5.646966933450414,1.3391066651568038,0.03211811510063545,6.581546542483,4750.825174993035
+16,7.615979550190041,2.634620209981552,0.03649712691597074,8.876433042179535,25339.645267361964
+17,9.460406027609398,3.381146758779285,0.03382351193072047,11.026114251293006,52508.67094062201
+18,16.000440812915354,3.044223392703431,0.00429464920004312,18.64853241598526,87446.89573211693
+19,10.55034307463129,1.5636343648164452,0.03026542177097602,12.296437149919917,12622.321615858824
+20,11.817436231011273,2.9829231066038613,0.021717088244806053,13.77323570047934,53748.08982906377
+21,12.10112174926391,5.305727115709343,0.014539688517578586,14.10387150263859,177305.22043219427
+22,10.400034982693642,3.5882967612571264,0.03231245358261764,12.121252893582334,65405.27202789767
+23,9.132539594795865,2.0115676174654173,0.04813020643130128,10.643985541720124,17913.12934373824
+24,9.880295706786633,2.4629082845730914,0.05063441840994648,11.515496161755983,29157.316024586242
+25,8.321803224944599,2.0030803613144172,0.054171203145458285,9.699071357744288,16104.920042215257
+26,6.131352093023454,1.794448837319034,0.03542910867558095,7.146098010516846,9295.960215551484
+27,11.430726526098171,3.9798910657304747,0.025330550582418253,13.322525088692508,90734.82138245075
+28,14.263809369877091,2.781733273746191,0.009358726934992883,16.624486445078194,66183.17957100201
+29,13.744160648067846,5.465225456257929,0.007517584283470983,16.01883525415833,240475.8375064763
+30,8.735251181801084,1.270630004139555,0.047066226355940106,10.180945433334598,6821.344536640248
+31,8.301748078970448,3.6767674099266174,0.02206953864653263,9.675697061737118,54122.886646054554

examples/wecsim_power_performance_example.ipynb

@@ -0,0 +1,268 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# MHKiT WEC-Sim Example\n",
+    "\n",
+    "This notebook demonstrates using the MHKiT wave module and [WEC-Sim](http://wec-sim.github.io/WEC-Sim/index.html) together to perform a resource characterization study in MHKiT, simulate representative cases with WEC-Sim, and visualize the results in MHKiT to estimate MAEP (Mean Annual Energy Production).\n",
+    "\n",
+    "  1. Characterize the available resource at a location\n",
+    "     - See the PacWave example notebook\n",
+    "  2. Write a WEC-Sim batch file for the given clusters\n",
+    "  3. Simulate the device *in WEC-Sim*.\n",
+    "     - Ensure that the spectra used in WEC-Sim is identical to the one used in MHKiT.\n",
+    "  4. Load WEC-Sim batch results\n",
+    "  5. Assess results and visualize quantities of interest\n",
+    "\n",
+    "This example uses WEC-Sim to simulate the [Oscillating Surge Wave Energy Converter (OSWEC)](https://wec-sim.github.io/WEC-Sim/main/user/tutorials.html#oscillating-surge-wec-oswec), a flap-type device.\n",
+    "\n",
+    "Start by importing MHKiT and the necessary python packages (e.g.`scipy.io`, `matplotlib.pyplot`, `pandas`, `numpy`)."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 1,
+   "metadata": {},
+   "outputs": [
+    {
+     "ename": "ModuleNotFoundError",
+     "evalue": "No module named 'mhkit'",
+     "output_type": "error",
+     "traceback": [
+      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
+      "\u001b[1;31mModuleNotFoundError\u001b[0m                       Traceback (most recent call last)",
+      "Cell \u001b[1;32mIn[1], line 1\u001b[0m\n\u001b[1;32m----> 1\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mmhkit\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m wave\n\u001b[0;32m      2\u001b[0m \u001b[38;5;28;01mimport\u001b[39;00m \u001b[38;5;21;01mscipy\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mio\u001b[39;00m \u001b[38;5;28;01mas\u001b[39;00m \u001b[38;5;21;01msio\u001b[39;00m\n\u001b[0;32m      3\u001b[0m \u001b[38;5;28;01mimport\u001b[39;00m \u001b[38;5;21;01mmatplotlib\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mpyplot\u001b[39;00m \u001b[38;5;28;01mas\u001b[39;00m \u001b[38;5;21;01mplt\u001b[39;00m\n",
+      "\u001b[1;31mModuleNotFoundError\u001b[0m: No module named 'mhkit'"
+     ]
+    }
+   ],
+   "source": [
+    "from mhkit import wave\n",
+    "import scipy.io as sio\n",
+    "import matplotlib.pyplot as plt\n",
+    "import pandas as pd\n",
+    "import numpy as np"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## 1. Characterize the available resource at a location\n",
+    "\n",
+    "This example will use an abbreviated version of `PacWave_resource_characterization_example.ipynb`. \n",
+    "For full details on downloading, calculating, and visualizing the k-means clusters representation of the site's wave resouce, see that example.\n",
+    "\n",
+    "We select the N=32 cluster as it's total energy flux is closet to the total energy flux of the site considering all wave conditions. We will load the PacWave example output, which can be easily saved after running the example with the command `results[32].to_csv(\"pacwave_cluster_32.csv\")`. We will start this example by reading in that csv output and formatting it for WEC-Sim."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 2,
+   "metadata": {},
+   "outputs": [
+    {
+     "ename": "NameError",
+     "evalue": "name 'pd' is not defined",
+     "output_type": "error",
+     "traceback": [
+      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
+      "\u001b[1;31mNameError\u001b[0m                                 Traceback (most recent call last)",
+      "Cell \u001b[1;32mIn[2], line 1\u001b[0m\n\u001b[1;32m----> 1\u001b[0m results \u001b[38;5;241m=\u001b[39m \u001b[43mpd\u001b[49m\u001b[38;5;241m.\u001b[39mread_csv(\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mdata/wave/pacwave_cluster_32.csv\u001b[39m\u001b[38;5;124m\"\u001b[39m, index_col\u001b[38;5;241m=\u001b[39m\u001b[38;5;241m0\u001b[39m)\n\u001b[0;32m      2\u001b[0m results\n",
+      "\u001b[1;31mNameError\u001b[0m: name 'pd' is not defined"
+     ]
+    }
+   ],
+   "source": [
+    "results = pd.read_csv(\"data/wave/pacwave_cluster_32.csv\", index_col=0)\n",
+    "results"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## 2. Write a WEC-Sim batch file for the given clusters\n",
+    "\n",
+    "[WEC-Sim MCR](https://wec-sim.github.io/WEC-Sim/main/user/advanced_features.html#multiple-condition-runs-mcr) (multiple condition run) files should contain a structure `mcr` that contains two variables: `header` and `cases`. Each column of `header` and `cases` denotes a variable and it's value respectively. Each row is another simulation. WEC-Sim defines waves using the significant wave height and peak period. We will isolate these values from the results of the cluster analysis and create a dictionary that is written to the `.mat` file."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 3,
+   "metadata": {
+    "scrolled": true
+   },
+   "outputs": [
+    {
+     "ename": "NameError",
+     "evalue": "name 'results' is not defined",
+     "output_type": "error",
+     "traceback": [
+      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
+      "\u001b[1;31mNameError\u001b[0m                                 Traceback (most recent call last)",
+      "Cell \u001b[1;32mIn[3], line 1\u001b[0m\n\u001b[1;32m----> 1\u001b[0m ws_mcr_cases \u001b[38;5;241m=\u001b[39m \u001b[43mresults\u001b[49m[[\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mHm0\u001b[39m\u001b[38;5;124m\"\u001b[39m,\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mTp\u001b[39m\u001b[38;5;124m\"\u001b[39m]]\n\u001b[0;32m      2\u001b[0m ws_mcr_header \u001b[38;5;241m=\u001b[39m np\u001b[38;5;241m.\u001b[39marray([\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mwaves.height\u001b[39m\u001b[38;5;124m\"\u001b[39m,\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mwaves.period\u001b[39m\u001b[38;5;124m\"\u001b[39m], dtype\u001b[38;5;241m=\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mobject\u001b[39m\u001b[38;5;124m'\u001b[39m)\n\u001b[0;32m      3\u001b[0m ws_mcr_out \u001b[38;5;241m=\u001b[39m {\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mmcr\u001b[39m\u001b[38;5;124m'\u001b[39m: {\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mheader\u001b[39m\u001b[38;5;124m'\u001b[39m: ws_mcr_header, \u001b[38;5;124m'\u001b[39m\u001b[38;5;124mcases\u001b[39m\u001b[38;5;124m'\u001b[39m: ws_mcr_cases}}\n",
+      "\u001b[1;31mNameError\u001b[0m: name 'results' is not defined"
+     ]
+    }
+   ],
+   "source": [
+    "ws_mcr_cases = results[[\"Hm0\",\"Tp\"]]\n",
+    "ws_mcr_header = np.array([\"waves.height\",\"waves.period\"], dtype='object')\n",
+    "ws_mcr_out = {'mcr': {'header': ws_mcr_header, 'cases': ws_mcr_cases}}\n",
+    "sio.savemat('mcr_mhkit.mat', ws_mcr_out)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## 3. Simulate the device *in WEC-Sim*\n",
+    "\n",
+    "Now that the MCR file is created, we need to go simulate WEC performance in these wave conditions using WEC-Sim. To recreate the data used in the next step, use the created MCR file with WEC-Sim's [OSWEC example](https://github.com/WEC-Sim/WEC-Sim/tree/main/examples/OSWEC). For an accurate comparison to the power calculated in the resource characterization, we should ensure that the WEC-Sim cases use irregular JONSWAP wave spectra as in the PacWave example.\n",
+    "\n",
+    "For convenience in this demonstration, we enforce OSWEC model stability in the extreme wave conditions by arbitrarily applying a large PTO stiffness and damping:\n",
+    "```\n",
+    "pto(1).stiffness = 1e5;\n",
+    "pto(1).damping = 5e7;\n",
+    "```\n",
+    "\n",
+    "To reduce the amount of extranenous data saved for this example, we limit the WEC-Sim output to the PTO's power output in the `userDefinedFunctions.m` script:\n",
+    "```\n",
+    "if exist('imcr','var')\n",
+    "    if imcr == 1\n",
+    "        nmcr = size(mcr.cases,1);\n",
+    "        power = nan(1, nmcr);\n",
+    "    end\n",
+    "\n",
+    "    iRampEnd = simu.rampTime./simu.dtOut + 1;\n",
+    "    power(imcr) = -mean(output.ptos(1).powerInternalMechanics(iRampEnd:end,5));\n",
+    "\n",
+    "    if imcr == nmcr\n",
+    "        % Save output\n",
+    "        save('mcr_mhkit_power.mat', 'power');\n",
+    "    end\n",
+    "end\n",
+    "```"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## 4. Load WEC-Sim batch results\n",
+    "\n",
+    "Note that in this example we do not save the entire WEC-Sim `output` structure for each case. See the `wecsim_example.ipynb` for information on loading that WEC-Sim data. Here the output is one array of average power output that we will load and compare to the resource characterization.\n",
+    "\n",
+    "Note that the power output \\[W\\] is significantly larger than the energy flux \\[W/m\\] due to the width of the OSWEC. "
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 4,
+   "metadata": {
+    "scrolled": true
+   },
+   "outputs": [
+    {
+     "ename": "NameError",
+     "evalue": "name 'sio' is not defined",
+     "output_type": "error",
+     "traceback": [
+      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
+      "\u001b[1;31mNameError\u001b[0m                                 Traceback (most recent call last)",
+      "Cell \u001b[1;32mIn[4], line 5\u001b[0m\n\u001b[0;32m      2\u001b[0m filename \u001b[38;5;241m=\u001b[39m \u001b[38;5;124m\"\u001b[39m\u001b[38;5;124m./data/wave/mcr_mhkit_power.mat\u001b[39m\u001b[38;5;124m\"\u001b[39m\n\u001b[0;32m      4\u001b[0m \u001b[38;5;66;03m# Load data using the `wecsim.read_output` function which returns a dictionary of Datasets.\u001b[39;00m\n\u001b[1;32m----> 5\u001b[0m wecsim_data \u001b[38;5;241m=\u001b[39m \u001b[43msio\u001b[49m\u001b[38;5;241m.\u001b[39mloadmat(filename)\n\u001b[0;32m      6\u001b[0m results[\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mP\u001b[39m\u001b[38;5;124m'\u001b[39m] \u001b[38;5;241m=\u001b[39m wecsim_data[\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mpower\u001b[39m\u001b[38;5;124m'\u001b[39m][\u001b[38;5;241m0\u001b[39m]\n\u001b[0;32m      7\u001b[0m results\n",
+      "\u001b[1;31mNameError\u001b[0m: name 'sio' is not defined"
+     ]
+    }
+   ],
+   "source": [
+    "# Relative location and filename of simulated WEC-Sim data (run with mooring)\n",
+    "filename = \"./data/wave/mcr_mhkit_power.mat\"\n",
+    "\n",
+    "# Load data using the `wecsim.read_output` function which returns a dictionary of Datasets.\n",
+    "wecsim_data = sio.loadmat(filename)\n",
+    "results['P'] = wecsim_data['power'][0]\n",
+    "results"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "## 5. Assess results and visualize quantities of interest\n",
+    "Now that we have loaded the OSWEC's modeled power, we can assess it's performance relative to the incoming wave and calculate the mean annual energy production (MAEP) using MHKiT."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 5,
+   "metadata": {},
+   "outputs": [
+    {
+     "ename": "NameError",
+     "evalue": "name 'wave' is not defined",
+     "output_type": "error",
+     "traceback": [
+      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
+      "\u001b[1;31mNameError\u001b[0m                                 Traceback (most recent call last)",
+      "Cell \u001b[1;32mIn[5], line 1\u001b[0m\n\u001b[1;32m----> 1\u001b[0m results[\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mCW\u001b[39m\u001b[38;5;124m'\u001b[39m] \u001b[38;5;241m=\u001b[39m \u001b[43mwave\u001b[49m\u001b[38;5;241m.\u001b[39mperformance\u001b[38;5;241m.\u001b[39mcapture_width(results[\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mP\u001b[39m\u001b[38;5;124m'\u001b[39m], results[\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mJ\u001b[39m\u001b[38;5;124m'\u001b[39m])\n\u001b[0;32m      2\u001b[0m oswec_width \u001b[38;5;241m=\u001b[39m \u001b[38;5;241m18\u001b[39m\n\u001b[0;32m      3\u001b[0m results[\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mCWR\u001b[39m\u001b[38;5;124m'\u001b[39m] \u001b[38;5;241m=\u001b[39m results[\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mCW\u001b[39m\u001b[38;5;124m'\u001b[39m] \u001b[38;5;241m/\u001b[39m oswec_width\n",
+      "\u001b[1;31mNameError\u001b[0m: name 'wave' is not defined"
+     ]
+    }
+   ],
+   "source": [
+    "results['CW'] = wave.performance.capture_width(results['P'], results['J'])\n",
+    "oswec_width = 18\n",
+    "results['CWR'] = results['CW'] / oswec_width\n",
+    "results"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 6,
+   "metadata": {},
+   "outputs": [
+    {
+     "ename": "NameError",
+     "evalue": "name 'wave' is not defined",
+     "output_type": "error",
+     "traceback": [
+      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
+      "\u001b[1;31mNameError\u001b[0m                                 Traceback (most recent call last)",
+      "Cell \u001b[1;32mIn[6], line 1\u001b[0m\n\u001b[1;32m----> 1\u001b[0m MAEP \u001b[38;5;241m=\u001b[39m \u001b[43mwave\u001b[49m\u001b[38;5;241m.\u001b[39mperformance\u001b[38;5;241m.\u001b[39mmean_annual_energy_production_matrix(results[\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mCW\u001b[39m\u001b[38;5;124m'\u001b[39m], results[\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mJ\u001b[39m\u001b[38;5;124m'\u001b[39m], results[\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mweights\u001b[39m\u001b[38;5;124m'\u001b[39m]) \u001b[38;5;241m/\u001b[39m \u001b[38;5;241m1000\u001b[39m \u001b[38;5;66;03m# kWh\u001b[39;00m\n\u001b[0;32m      2\u001b[0m MAEP \u001b[38;5;241m=\u001b[39m np\u001b[38;5;241m.\u001b[39mround(MAEP, \u001b[38;5;241m0\u001b[39m)\u001b[38;5;241m.\u001b[39mitem()\n\u001b[0;32m      3\u001b[0m MAEP\n",
+      "\u001b[1;31mNameError\u001b[0m: name 'wave' is not defined"
+     ]
+    }
+   ],
+   "source": [
+    "MAEP = wave.performance.mean_annual_energy_production_matrix(results['CW'], results['J'], results['weights']) / 1000 # kWh\n",
+    "MAEP = np.round(MAEP, 0).item()\n",
+    "MAEP"
+   ]
+  }
+ ],
+ "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.10.9"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}

mhkit/tests/wave/io/test_wecsim.py

@@ -56,7 +56,7 @@ def test_read_wecSim_cable(self):
         )
         self.assertEqual(ws_output["wave"]["elevation"].name, "elevation")
         self.assertEqual(
-            ws_output["bodies"]["body1"]["position_dof1"].name, "position_dof1"
+            ws_output["bodies"].sel(body=1, dof=1)["position"].name, "position"
         )
         self.assertEqual(len(ws_output["mooring"]), 0)
         self.assertEqual(len(ws_output["moorDyn"]), 0)