diff --git a/docs/source/examples/multicomponent_flash.ipynb b/docs/source/examples/multicomponent_flash.ipynb new file mode 100644 index 0000000..ee5ed17 --- /dev/null +++ b/docs/source/examples/multicomponent_flash.ipynb @@ -0,0 +1,286 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Multi-Component Flash Separation (BTX)\n", + "\n", + "Simulating an isothermal flash drum for a ternary **benzene–toluene–p-xylene** (BTX) mixture.\n", + "The `MultiComponentFlash` block uses Raoult's law with Antoine correlations to compute K-values\n", + "and solves the Rachford-Rice equation via Brent's method.\n", + "\n", + "This example is inspired by [MiniSim's SimpleFlash example](https://github.com/Nukleon84/MiniSim/blob/master/doc/SimpleFlash.ipynb),\n", + "adapted to PathSim's dynamic simulation framework.\n", + "\n", + "**Feed conditions:**\n", + "- 10 mol/s total flow\n", + "- 50% benzene, 10% toluene, 40% p-xylene (molar)\n", + "- 1 atm pressure\n", + "- Temperature sweep: 340 K → 420 K" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np\n", + "import matplotlib.pyplot as plt\n", + "\n", + "from pathsim import Simulation, Connection\n", + "from pathsim.blocks import Source, Scope\n", + "\n", + "from pathsim_chem.process import MultiComponentFlash" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Flash Drum Setup\n", + "\n", + "The `MultiComponentFlash` block defaults to BTX Antoine parameters (ln form, Pa, K):\n", + "\n", + "| Component | A | B | C |\n", + "|-----------|-------|---------|--------|\n", + "| Benzene | 20.7936 | 2788.51 | -52.36 |\n", + "| Toluene | 20.9064 | 3096.52 | -53.67 |\n", + "| p-Xylene | 20.9891 | 3346.65 | -57.84 |\n", + "\n", + "We feed the drum with constant composition and pressure while ramping temperature\n", + "to observe the transition from all-liquid through two-phase to all-vapor." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Create the flash drum (3 components, BTX defaults)\n", + "flash = MultiComponentFlash(N_comp=3, holdup=100.0)\n", + "\n", + "# Feed sources\n", + "F_feed = Source(func=lambda t: 10.0) # 10 mol/s\n", + "z_benzene = Source(func=lambda t: 0.5) # 50% benzene\n", + "z_toluene = Source(func=lambda t: 0.1) # 10% toluene (p-xylene = 40% inferred)\n", + "T_sweep = Source(func=lambda t: 340.0 + t) # ramp 340 -> 420 K\n", + "P_feed = Source(func=lambda t: 101325.0) # 1 atm\n", + "\n", + "# Record all outputs: V_rate, L_rate, y_benzene, y_toluene, x_benzene, x_toluene\n", + "scp = Scope(labels=[\"V_rate\", \"L_rate\", \"y_benz\", \"y_tol\", \"x_benz\", \"x_tol\"])\n", + "\n", + "sim = Simulation(\n", + " blocks=[F_feed, z_benzene, z_toluene, T_sweep, P_feed, flash, scp],\n", + " connections=[\n", + " Connection(F_feed, flash), # F -> port 0\n", + " Connection(z_benzene, flash[1]), # z_1 (benzene) -> port 1\n", + " Connection(z_toluene, flash[2]), # z_2 (toluene) -> port 2\n", + " Connection(T_sweep, flash[3]), # T -> port 3\n", + " Connection(P_feed, flash[4]), # P -> port 4\n", + " Connection(flash, scp), # V_rate -> scope 0\n", + " Connection(flash[1], scp[1]), # L_rate -> scope 1\n", + " Connection(flash[2], scp[2]), # y_benzene -> scope 2\n", + " Connection(flash[3], scp[3]), # y_toluene -> scope 3\n", + " Connection(flash[4], scp[4]), # x_benzene -> scope 4\n", + " Connection(flash[5], scp[5]), # x_toluene -> scope 5\n", + " ],\n", + " dt=0.5,\n", + ")\n", + "\n", + "sim.run(80)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Results: Flow Rates\n", + "\n", + "As temperature increases, the vapor fraction grows. Below the bubble point the drum produces\n", + "only liquid; above the dew point it produces only vapor." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "time, signals = scp.read()\n", + "T = 340.0 + time # temperature axis\n", + "\n", + "V_rate, L_rate = signals[0], signals[1]\n", + "y_benz, y_tol = signals[2], signals[3]\n", + "x_benz, x_tol = signals[4], signals[5]\n", + "\n", + "# Infer p-xylene fractions\n", + "y_xyl = 1.0 - y_benz - y_tol\n", + "x_xyl = 1.0 - x_benz - x_tol\n", + "\n", + "fig, axes = plt.subplots(1, 3, figsize=(15, 5))\n", + "\n", + "# Flow rates\n", + "ax = axes[0]\n", + "ax.plot(T, V_rate, label=\"Vapor\")\n", + "ax.plot(T, L_rate, label=\"Liquid\")\n", + "ax.set_xlabel(\"Temperature [K]\")\n", + "ax.set_ylabel(\"Flow rate [mol/s]\")\n", + "ax.set_title(\"Flash Drum Flow Rates\")\n", + "ax.legend()\n", + "ax.grid(True, alpha=0.3)\n", + "\n", + "# Vapor compositions\n", + "ax = axes[1]\n", + "ax.plot(T, y_benz, label=\"Benzene\")\n", + "ax.plot(T, y_tol, label=\"Toluene\")\n", + "ax.plot(T, y_xyl, label=\"p-Xylene\")\n", + "ax.set_xlabel(\"Temperature [K]\")\n", + "ax.set_ylabel(\"Vapor mole fraction\")\n", + "ax.set_title(\"Vapor Composition\")\n", + "ax.legend()\n", + "ax.grid(True, alpha=0.3)\n", + "\n", + "# Liquid compositions\n", + "ax = axes[2]\n", + "ax.plot(T, x_benz, label=\"Benzene\")\n", + "ax.plot(T, x_tol, label=\"Toluene\")\n", + "ax.plot(T, x_xyl, label=\"p-Xylene\")\n", + "ax.set_xlabel(\"Temperature [K]\")\n", + "ax.set_ylabel(\"Liquid mole fraction\")\n", + "ax.set_title(\"Liquid Composition\")\n", + "ax.legend()\n", + "ax.grid(True, alpha=0.3)\n", + "\n", + "plt.tight_layout()\n", + "plt.show()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Results: VLE Diagram\n", + "\n", + "Plot the vapor vs liquid composition for each component across the temperature sweep.\n", + "The diagonal represents equal vapor and liquid composition — deviation from it shows\n", + "the separation achieved by the flash." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "fig, axes = plt.subplots(1, 3, figsize=(15, 5))\n", + "\n", + "for ax, xi, yi, name in zip(axes,\n", + " [x_benz, x_tol, x_xyl],\n", + " [y_benz, y_tol, y_xyl],\n", + " [\"Benzene\", \"Toluene\", \"p-Xylene\"]):\n", + " ax.plot(xi, yi, \".\", markersize=3)\n", + " ax.plot([0, 1], [0, 1], \"k--\", alpha=0.3)\n", + " ax.set_xlabel(f\"$x$ ({name})\")\n", + " ax.set_ylabel(f\"$y$ ({name})\")\n", + " ax.set_title(f\"{name} (x, y)-Diagram\")\n", + " ax.set_aspect(\"equal\")\n", + " ax.grid(True, alpha=0.3)\n", + "\n", + "plt.tight_layout()\n", + "plt.show()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Fixed-Temperature Flash at 380 K\n", + "\n", + "For a direct comparison with MiniSim's result (which solves at steady state),\n", + "we run a fixed-temperature flash and let the holdup reach equilibrium." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "flash2 = MultiComponentFlash(N_comp=3, holdup=100.0)\n", + "\n", + "F_src = Source(func=lambda t: 10.0)\n", + "z1_src = Source(func=lambda t: 0.5)\n", + "z2_src = Source(func=lambda t: 0.1)\n", + "T_src = Source(func=lambda t: 380.0) # fixed at 380 K (~107 °C)\n", + "P_src = Source(func=lambda t: 101325.0)\n", + "\n", + "scp2 = Scope(labels=[\"V_rate\", \"L_rate\", \"y_benz\", \"y_tol\", \"x_benz\", \"x_tol\"])\n", + "\n", + "sim2 = Simulation(\n", + " blocks=[F_src, z1_src, z2_src, T_src, P_src, flash2, scp2],\n", + " connections=[\n", + " Connection(F_src, flash2),\n", + " Connection(z1_src, flash2[1]),\n", + " Connection(z2_src, flash2[2]),\n", + " Connection(T_src, flash2[3]),\n", + " Connection(P_src, flash2[4]),\n", + " Connection(flash2, scp2),\n", + " Connection(flash2[1], scp2[1]),\n", + " Connection(flash2[2], scp2[2]),\n", + " Connection(flash2[3], scp2[3]),\n", + " Connection(flash2[4], scp2[4]),\n", + " Connection(flash2[5], scp2[5]),\n", + " ],\n", + " dt=0.5,\n", + ")\n", + "\n", + "sim2.run(100) # let it reach steady state\n", + "\n", + "time2, signals2 = scp2.read()\n", + "\n", + "# Extract final steady-state values\n", + "V_ss = signals2[0][-1]\n", + "L_ss = signals2[1][-1]\n", + "y_benz_ss = signals2[2][-1]\n", + "y_tol_ss = signals2[3][-1]\n", + "x_benz_ss = signals2[4][-1]\n", + "x_tol_ss = signals2[5][-1]\n", + "\n", + "print(\"BTX Flash at 380 K, 1 atm\")\n", + "print(\"=\" * 40)\n", + "print(f\"{'':20s} {'Vapor':>10s} {'Liquid':>10s}\")\n", + "print(f\"{'-'*40}\")\n", + "print(f\"{'Flow rate [mol/s]':20s} {V_ss:10.3f} {L_ss:10.3f}\")\n", + "print(f\"{'Benzene':20s} {y_benz_ss:10.4f} {x_benz_ss:10.4f}\")\n", + "print(f\"{'Toluene':20s} {y_tol_ss:10.4f} {x_tol_ss:10.4f}\")\n", + "print(f\"{'p-Xylene':20s} {1-y_benz_ss-y_tol_ss:10.4f} {1-x_benz_ss-x_tol_ss:10.4f}\")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The lighter component (benzene) is enriched in the vapor phase while the heavier\n", + "component (p-xylene) concentrates in the liquid — exactly the separation behaviour\n", + "expected from VLE. The dynamic formulation reaches the same steady state that\n", + "an equation-oriented solver (like MiniSim) finds directly." + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "name": "python", + "version": "3.11.0" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +} diff --git a/docs/source/examples/process_flowsheet.ipynb b/docs/source/examples/process_flowsheet.ipynb new file mode 100644 index 0000000..a5178b3 --- /dev/null +++ b/docs/source/examples/process_flowsheet.ipynb @@ -0,0 +1,246 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Process Flowsheet: Mixer → Heater → CSTR\n", + "\n", + "A simple process combining multiple unit operations into a flowsheet:\n", + "\n", + "1. **Mixer** — blends a fresh feed stream with a recycle-like second stream\n", + "2. **Heater** — brings the mixed stream to reactor inlet temperature\n", + "3. **CSTR** — exothermic first-order reaction A → products with cooling jacket\n", + "\n", + "This demonstrates how PathSim-Chem blocks compose into multi-unit simulations,\n", + "inspired by the reactor front-end of [MiniSim's Cumene Process](https://github.com/Nukleon84/MiniSim/blob/master/doc/Cumene%20Process.ipynb)." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np\n", + "import matplotlib.pyplot as plt\n", + "\n", + "from pathsim import Simulation, Connection\n", + "from pathsim.blocks import Source, Scope\n", + "\n", + "from pathsim_chem.process import Mixer, Heater, CSTR" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## System Setup\n", + "\n", + "**Fresh feed:** 0.05 m³/s at 300 K with 2.0 mol/m³ reactant A \n", + "**Second stream:** 0.05 m³/s at 310 K (simulating a warm recycle) \n", + "**Heater:** 200 kW duty to preheat the mixed stream \n", + "**CSTR:** 1 m³ reactor, first-order Arrhenius kinetics, cooled by a 290 K jacket" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Unit operations\n", + "mixer = Mixer()\n", + "heater = Heater(rho=1000.0, Cp=4184.0) # water-like fluid\n", + "cstr = CSTR(\n", + " V=1.0, # reactor volume [m³]\n", + " F=0.1, # total flow [m³/s] (will be overridden by mixer output)\n", + " k0=1e6, # pre-exponential factor [1/s]\n", + " Ea=50000.0, # activation energy [J/mol]\n", + " n=1.0, # first-order reaction\n", + " dH_rxn=-50000.0, # exothermic [J/mol]\n", + " rho=1000.0,\n", + " Cp=4184.0,\n", + " UA=500.0, # heat transfer [W/K]\n", + " C_A0=0.0, # start empty\n", + " T0=300.0, # start at 300 K\n", + ")\n", + "\n", + "# Feed sources\n", + "F1_src = Source(func=lambda t: 0.05) # fresh feed flow [m³/s]\n", + "T1_src = Source(func=lambda t: 300.0) # fresh feed temp [K]\n", + "F2_src = Source(func=lambda t: 0.05) # second stream flow [m³/s]\n", + "T2_src = Source(func=lambda t: 310.0) # second stream temp [K]\n", + "\n", + "# Heater duty [W]\n", + "Q_src = Source(func=lambda t: 200000.0)\n", + "\n", + "# Reactant concentration in combined feed [mol/m³]\n", + "C_in_src = Source(func=lambda t: 2.0)\n", + "\n", + "# Coolant temperature [K]\n", + "T_cool = Source(func=lambda t: 290.0)\n", + "\n", + "# Scopes for recording\n", + "scp_mix = Scope(labels=[\"F_mix\", \"T_mix\"])\n", + "scp_heat = Scope(labels=[\"F_heat\", \"T_heated\"])\n", + "scp_cstr = Scope(labels=[\"C_out\", \"T_reactor\"])" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Wiring the Flowsheet\n", + "\n", + "```\n", + "Feed 1 (F₁, T₁) ─┐\n", + " ├─ Mixer ─── Heater ─── CSTR ─── Products\n", + "Feed 2 (F₂, T₂) ─┘ ↑ ↑ ↑\n", + " │ Q T_coolant\n", + "```\n", + "\n", + "Note: The CSTR's internal flow rate `F` is set at construction. Here the mixer's output flow\n", + "is used for monitoring, while the CSTR uses its own `F=0.1` parameter for residence time." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "sim = Simulation(\n", + " blocks=[\n", + " F1_src, T1_src, F2_src, T2_src, Q_src, C_in_src, T_cool,\n", + " mixer, heater, cstr,\n", + " scp_mix, scp_heat, scp_cstr,\n", + " ],\n", + " connections=[\n", + " # Mixer inputs: (F_1, T_1, F_2, T_2)\n", + " Connection(F1_src, mixer),\n", + " Connection(T1_src, mixer[1]),\n", + " Connection(F2_src, mixer[2]),\n", + " Connection(T2_src, mixer[3]),\n", + "\n", + " # Mixer outputs -> Heater inputs: (F, T_in)\n", + " Connection(mixer, heater), # F_out -> F\n", + " Connection(mixer[1], heater[1]), # T_out -> T_in\n", + "\n", + " # Heater heat duty\n", + " Connection(Q_src, heater[2]), # Q\n", + "\n", + " # Heater output temperature -> CSTR inlet temperature\n", + " Connection(C_in_src, cstr), # C_in -> port 0\n", + " Connection(heater[1], cstr[1]), # T_heated -> T_in (port 1)\n", + "\n", + " # CSTR coolant\n", + " Connection(T_cool, cstr[2]), # T_c -> port 2\n", + "\n", + " # Recording\n", + " Connection(mixer, scp_mix),\n", + " Connection(mixer[1], scp_mix[1]),\n", + " Connection(heater, scp_heat),\n", + " Connection(heater[1], scp_heat[1]),\n", + " Connection(cstr, scp_cstr),\n", + " Connection(cstr[1], scp_cstr[1]),\n", + " ],\n", + " dt=0.1,\n", + ")\n", + "\n", + "sim.run(100)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Results: Temperature Profile Through the Plant\n", + "\n", + "Track the temperature at each stage: after mixing, after heating, and inside the reactor." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "t_mix, s_mix = scp_mix.read()\n", + "t_heat, s_heat = scp_heat.read()\n", + "t_cstr, s_cstr = scp_cstr.read()\n", + "\n", + "fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 5))\n", + "\n", + "# Temperature evolution\n", + "ax1.plot(t_mix, s_mix[1], label=\"After Mixer\")\n", + "ax1.plot(t_heat, s_heat[1], label=\"After Heater\")\n", + "ax1.plot(t_cstr, s_cstr[1], label=\"Reactor (CSTR)\")\n", + "ax1.axhline(290, color=\"blue\", ls=\"--\", alpha=0.3, label=\"Coolant\")\n", + "ax1.set_xlabel(\"Time [s]\")\n", + "ax1.set_ylabel(\"Temperature [K]\")\n", + "ax1.set_title(\"Temperature Profile\")\n", + "ax1.legend()\n", + "ax1.grid(True, alpha=0.3)\n", + "\n", + "# Reactor concentration\n", + "ax2.plot(t_cstr, s_cstr[0], color=\"tab:red\")\n", + "ax2.set_xlabel(\"Time [s]\")\n", + "ax2.set_ylabel(\"Concentration $C_A$ [mol/m³]\")\n", + "ax2.set_title(\"Reactor Outlet Concentration\")\n", + "ax2.grid(True, alpha=0.3)\n", + "\n", + "plt.tight_layout()\n", + "plt.show()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Steady-State Summary" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "print(\"Process Steady State\")\n", + "print(\"=\" * 45)\n", + "print(f\"{'Unit':15s} {'Flow [m³/s]':>12s} {'T [K]':>10s}\")\n", + "print(\"-\" * 45)\n", + "print(f\"{'Fresh feed':15s} {'0.05':>12s} {'300.0':>10s}\")\n", + "print(f\"{'Second stream':15s} {'0.05':>12s} {'310.0':>10s}\")\n", + "print(f\"{'Mixer outlet':15s} {s_mix[0][-1]:12.4f} {s_mix[1][-1]:10.2f}\")\n", + "print(f\"{'Heater outlet':15s} {s_heat[0][-1]:12.4f} {s_heat[1][-1]:10.2f}\")\n", + "print(f\"{'CSTR outlet':15s} {'—':>12s} {s_cstr[1][-1]:10.2f}\")\n", + "print(f\"\\nReactor C_A = {s_cstr[0][-1]:.4f} mol/m³ (feed = 2.0 mol/m³)\")\n", + "print(f\"Conversion = {(1 - s_cstr[0][-1]/2.0)*100:.1f}%\")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The mixer blends both feeds to an intermediate temperature (~305 K). The heater raises it\n", + "further. The CSTR then reaches a steady state where the exothermic reaction heat is balanced\n", + "by cooling jacket removal, and the conversion depends on the residence time and Arrhenius rate." + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "name": "python", + "version": "3.11.0" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +} diff --git a/docs/source/examples/valve_characteristics.ipynb b/docs/source/examples/valve_characteristics.ipynb new file mode 100644 index 0000000..01e8122 --- /dev/null +++ b/docs/source/examples/valve_characteristics.ipynb @@ -0,0 +1,165 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Valve Flow Characteristics\n", + "\n", + "Exploring the flow behaviour of the `Valve` block under varying pressure drop conditions.\n", + "The valve models isenthalpic flow using a standard Cv equation:\n", + "\n", + "$$F = C_v \\sqrt{|\\Delta P|} \\cdot \\text{sign}(\\Delta P)$$\n", + "\n", + "This example is inspired by [MiniSim's ValveExample](https://github.com/Nukleon84/MiniSim/blob/master/doc/ValveExample.ipynb),\n", + "adapted to PathSim's dynamic framework." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np\n", + "import matplotlib.pyplot as plt\n", + "\n", + "from pathsim import Simulation, Connection\n", + "from pathsim.blocks import Source, Scope\n", + "\n", + "from pathsim_chem.process import Valve" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Pressure Drop Sweep\n", + "\n", + "Hold the outlet pressure at 1 bar and ramp the inlet pressure from 1 bar to 10 bar.\n", + "Compare valves with different Cv coefficients." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "results = {}\n", + "\n", + "for Cv in [0.5, 1.0, 2.0, 5.0]:\n", + " valve = Valve(Cv=Cv)\n", + "\n", + " P_in = Source(func=lambda t: 100000.0 + t * 10000.0) # 1 bar -> 10 bar\n", + " P_out = Source(func=lambda t: 100000.0) # 1 bar fixed\n", + " T_in = Source(func=lambda t: 350.0) # 350 K\n", + "\n", + " scp = Scope(labels=[\"F\", \"T_out\"])\n", + "\n", + " sim = Simulation(\n", + " blocks=[P_in, P_out, T_in, valve, scp],\n", + " connections=[\n", + " Connection(P_in, valve), # P_in -> port 0\n", + " Connection(P_out, valve[1]), # P_out -> port 1\n", + " Connection(T_in, valve[2]), # T_in -> port 2\n", + " Connection(valve, scp), # F -> scope 0\n", + " Connection(valve[1], scp[1]), # T_out -> scope 1\n", + " ],\n", + " dt=0.5,\n", + " )\n", + "\n", + " sim.run(90)\n", + " time, signals = scp.read()\n", + " dP = (100000.0 + time * 10000.0) - 100000.0 # pressure drop [Pa]\n", + " results[f\"Cv={Cv}\"] = (dP / 1e5, signals[0]) # convert dP to bar" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(figsize=(10, 6))\n", + "\n", + "for label, (dP_bar, F) in results.items():\n", + " ax.plot(dP_bar, F, label=label)\n", + "\n", + "ax.set_xlabel(r\"Pressure drop $\\Delta P$ [bar]\")\n", + "ax.set_ylabel(\"Flow rate $F$\")\n", + "ax.set_title(\"Valve Characteristic Curves\")\n", + "ax.legend()\n", + "ax.grid(True, alpha=0.3)\n", + "\n", + "plt.tight_layout()\n", + "plt.show()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The square-root relationship between pressure drop and flow rate is clearly visible.\n", + "Doubling the Cv coefficient doubles the flow at the same pressure drop." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Isenthalpic Expansion\n", + "\n", + "The valve is an isenthalpic device — outlet temperature equals inlet temperature\n", + "regardless of pressure drop. Verify this across the full sweep." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "# Re-run one case and check T_out\n", + "valve = Valve(Cv=2.0)\n", + "\n", + "P_in = Source(func=lambda t: 100000.0 + t * 10000.0)\n", + "P_out = Source(func=lambda t: 100000.0)\n", + "T_in = Source(func=lambda t: 350.0)\n", + "\n", + "scp = Scope(labels=[\"F\", \"T_out\"])\n", + "\n", + "sim = Simulation(\n", + " blocks=[P_in, P_out, T_in, valve, scp],\n", + " connections=[\n", + " Connection(P_in, valve),\n", + " Connection(P_out, valve[1]),\n", + " Connection(T_in, valve[2]),\n", + " Connection(valve, scp),\n", + " Connection(valve[1], scp[1]),\n", + " ],\n", + " dt=0.5,\n", + ")\n", + "\n", + "sim.run(90)\n", + "time, signals = scp.read()\n", + "\n", + "print(f\"T_out range: {signals[1].min():.2f} K to {signals[1].max():.2f} K\")\n", + "print(f\"T_in = 350.00 K (constant) -> isenthalpic confirmed\")" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "name": "python", + "version": "3.11.0" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +}