{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "# BEE 4750 Lab 1: Convergence and Discretization\n",
        "\n",
        "**Name**:\n",
        "\n",
        "**ID**:\n",
        "\n",
        "> **Due Date**\n",
        ">\n",
        "> TBD\n",
        "\n",
        "## Setup"
      ],
      "id": "16a471a6-de10-4c49-af2c-2293e1a903b7"
    },
    {
      "cell_type": "code",
      "execution_count": 1,
      "metadata": {},
      "outputs": [],
      "source": [
        "import Pkg\n",
        "Pkg.activate(\".\")\n",
        "Pkg.instantiate()"
      ],
      "id": "2"
    },
    {
      "cell_type": "code",
      "execution_count": 1,
      "metadata": {},
      "outputs": [],
      "source": [
        "using Plots\n",
        "using LaTeXStrings"
      ],
      "id": "4"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Introduction\n",
        "\n",
        "In lecture, we discretized a first-order ODE with **forward Euler** and\n",
        "used a **convergence study** (halve $\\Delta t$, watch a quantity of\n",
        "interest, plot error against $\\Delta t$ on log-log axes) to decide\n",
        "whether a step size was small enough. We did this for a box model of an\n",
        "airshed, where the pollutant concentration $C(t)$ satisfied\n",
        "\n",
        "$$\\frac{dC}{dt} = P(t) - lC(t),$$\n",
        "\n",
        "with $P(t)$ a (possibly time-varying) source term and $l$ a first-order\n",
        "loss rate. Once $P(t)$ varied in time, there was no closed-form\n",
        "solution, so we built a trustworthy “truth” by running forward Euler at\n",
        "a very fine $\\Delta t$, and checked coarser $\\Delta t$ against it.\n",
        "\n",
        "This lab applies exactly that workflow to a new setting: a **landfill\n",
        "gas generation model**. If you get stuck on the numerics, it’s worth\n",
        "pulling up that lecture’s code side-by-side with this notebook.\n",
        "\n",
        "### The Model\n",
        "\n",
        "Landfilled waste decomposes and generates gas (mostly CH<sub>4</sub> and\n",
        "CO<sub>2</sub>) at a rate that responds, with some lag, to how much\n",
        "waste is currently in place. We’ll model the landfill gas generation\n",
        "rate $G(t)$ with\n",
        "\n",
        "$$\\frac{dG}{dt} = k\\left(L_0 \\, W(t) - G(t)\\right),$$\n",
        "\n",
        "| Symbol | Meaning | Units |\n",
        "|:--:|:---|:---|\n",
        "| $G(t)$ | landfill gas generation rate (the state) | m<sup>3</sup> CH<sub>4</sub>/yr |\n",
        "| $W(t)$ | waste acceptance rate | Mg/yr |\n",
        "| $L_0$ | methane generation potential per unit waste | m<sup>3</sup> CH<sub>4</sub>/Mg |\n",
        "| $k$ | first-order decay/generation-response rate | 1/yr |\n",
        "\n",
        "This has the same structure as the airshed model from lecture: a source\n",
        "term driving the state up, and a first-order loss term (here, $-kG$)\n",
        "pulling it back down, so that $G(t)$ tracks $L_0 W(t)$ with a lag set by\n",
        "$k$.\n",
        "\n",
        "We’ll use the following scenario: a landfill accepts waste at a constant\n",
        "rate $W_\\text{max}$ for a 20-year active life, then closes (waste\n",
        "acceptance drops to zero), and we track gas generation for another 20\n",
        "years after that ($T = 40$ yr total)."
      ],
      "id": "0f091de1-7b03-4992-843c-69a6e6167089"
    },
    {
      "cell_type": "code",
      "execution_count": 1,
      "metadata": {},
      "outputs": [],
      "source": [
        "landfill_params = let\n",
        "    k = 0.15            # first-order decay/generation-response rate  [1/yr]\n",
        "    L0 = 100.0          # methane generation potential                [m³ CH₄ / Mg]\n",
        "    Wmax = 1.0e5        # waste acceptance rate while open             [Mg/yr]\n",
        "    active_life = 20.0  # years the landfill accepts waste             [yr]\n",
        "    (; k, L0, Wmax, active_life)\n",
        "end\n",
        "\n",
        "T_total = 40.0   # total simulation horizon [yr]\n",
        "G0 = 0.0         # initial gas generation rate [m³ CH₄/yr]\n",
        "\n",
        "# Waste acceptance rate: constant while the landfill is open, zero after closure.\n",
        "waste_rate(t, p) = t < p.active_life ? p.Wmax : 0.0"
      ],
      "id": "6"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "> **Two Ways to Work Through This Lab**\n",
        ">\n",
        "> You can do the coding either in this notebook (`lab01.ipynb`) or in\n",
        "> the plain script `lab01.jl`, whichever you prefer — both contain the\n",
        "> same problems. Either way, the written/derivation parts (marked below)\n",
        "> go on the separate **Activity Sheet** (handed out in class), which\n",
        "> you’ll hand in or submit by the end of the day to Gradescope.\n",
        "\n",
        "## Problem 1: Discretize and Implement\n",
        "\n",
        "Before writing any code, complete Problem 1 on the Activity Sheet:\n",
        "derive the forward Euler update for this model by hand, starting from\n",
        "$\\frac{dG}{dt} = k(L_0 W(t) - G(t))$.\n",
        "\n",
        "Then complete the skeleton below, using your derivation.\n",
        "`landfill_gas_rate` should return $dG/dt$ given the current state, time,\n",
        "and parameters; you shouldn’t need to change it. The one line you need\n",
        "to fill in is the forward Euler update inside `landfill_gas_simulate`,\n",
        "using your derivation from Problem 1."
      ],
      "id": "a0502cb2-4375-4cf0-bca1-c7c6621b9ac3"
    },
    {
      "cell_type": "code",
      "execution_count": 0,
      "metadata": {},
      "outputs": [],
      "source": [
        "landfill_gas_rate(G, t, p) = p.k * (p.L0 * waste_rate(t, p) - G)\n",
        "\n",
        "function landfill_gas_simulate(G_ic, T, Δt, p)\n",
        "    steps = Int(round(T / Δt))\n",
        "    G = zeros(steps + 1)   # index 1 holds the initial condition\n",
        "    G[1] = G_ic\n",
        "    for i in 1:steps\n",
        "        t = (i - 1) * Δt\n",
        "        # TODO: fill in the forward Euler update from Problem 1\n",
        "        G[i+1] = missing\n",
        "    end\n",
        "    return G\n",
        "end"
      ],
      "id": "8"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Sanity-check your implementation by running it at a moderate step size\n",
        "and plotting $G(t)$ over the full 40-year horizon. Does the shape make\n",
        "sense (rising while the landfill is open, declining after closure)?"
      ],
      "id": "485862bc-3a57-4ee8-a8f6-7d94c9c0e69a"
    },
    {
      "cell_type": "code",
      "execution_count": 0,
      "metadata": {},
      "outputs": [],
      "source": [
        "# insert your sanity-check simulation and plot here"
      ],
      "id": "10"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Problem 2: How Small Does $\\Delta t$ Need to Be?\n",
        "\n",
        "Before running anything, complete Problem 2 on the Activity Sheet:\n",
        "compute the model’s characteristic timescale, $1/|\\text{rate}|$, and use\n",
        "it to propose a starting $\\Delta t$ that should be “well inside” it.\n",
        "\n",
        "Now check that prediction empirically. Build a reference solution at a\n",
        "very fine step size (finer than anything you’ll test), and use the peak\n",
        "gas generation rate, $\\max_t G(t)$, as your quantity of interest — this\n",
        "is the number you’d use to size flare or collection-system capacity."
      ],
      "id": "89871b4c-55ce-47a8-979a-95a0878af522"
    },
    {
      "cell_type": "code",
      "execution_count": 0,
      "metadata": {},
      "outputs": [],
      "source": [
        "Δt_ref = 0 # replace this with a very fine reference value\n",
        "G_ref = landfill_gas_simulate(G0, T_total, Δt_ref, landfill_params)\n",
        "peak_ref = maximum(G_ref)"
      ],
      "id": "12"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Choose your **own** sequence of step sizes to test, starting from\n",
        "wherever your Problem 2a estimate suggested and successively halving.\n",
        "For each, compute the peak gas generation rate and its error relative to\n",
        "the reference. Then make a log-log plot of error vs. $\\Delta t$ (include\n",
        "a reference line of slope 1, as in lecture, to compare against)."
      ],
      "id": "92c49d35-4681-4d22-8255-34252e16afe6"
    },
    {
      "cell_type": "code",
      "execution_count": 0,
      "metadata": {},
      "outputs": [],
      "source": [
        "# TODO: choose your own step sizes here (successively halved)\n",
        "Δts = []\n",
        "\n",
        "# TODO: compute the peak G at each Δt, and the error against peak_ref\n",
        "\n",
        "# TODO: make the log-log convergence plot"
      ],
      "id": "14"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Estimate the relative error improvement between successive step-size\n",
        "refinements (: $\\log_2(\\text{err}_i / \\text{err}_{i+1})$)."
      ],
      "id": "cc6c36db-1356-4e45-ac39-a59aec431d03"
    },
    {
      "cell_type": "code",
      "execution_count": 0,
      "metadata": {},
      "outputs": [],
      "source": [
        "# TODO: compute empirical orders between successive Δts"
      ],
      "id": "16"
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Transcribe your Δt sequence, peak values, errors, and empirical orders\n",
        "into Table 1 on the Activity Sheet, then answer the interpretation\n",
        "question there (is this consistent with forward Euler being first order,\n",
        "and why).\n",
        "\n",
        "## Problem 3: Make and Justify a Recommendation\n",
        "\n",
        "Complete Problem 3 on the Activity Sheet: a short written recommendation\n",
        "for what $\\Delta t$ you’d actually use to size a flare system from this\n",
        "model, justified using both the characteristic timescale and your\n",
        "empirical convergence check."
      ],
      "id": "50df1786-5c3b-48b8-bd15-9f3880f8dd2b"
    }
  ],
  "nbformat": 4,
  "nbformat_minor": 5,
  "metadata": {
    "kernel_info": {
      "name": "julia"
    },
    "kernelspec": {
      "name": "julia",
      "display_name": "Julia",
      "language": "julia"
    },
    "language_info": {
      "name": "julia",
      "codemirror_mode": "julia",
      "version": "1.11.5"
    }
  }
}