{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "aa423240-3daf-4189-8d12-7c251b4beb19",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "\"\"\"\n",
    "    ans = mytrapezoids(fun, a, b, n)\n",
    "\n",
    "Numerically evaluate the integral int_a^b fun(x) dx using the\n",
    "trapezoidal rule: I = h/2*(f_1 + 2f_2 + ... + 2f_{n-1} + f_n),\n",
    "where h = (b - a)/(n - 1), x_i = a + (i - 1)*h, f_i = fun(x_i).\n",
    "\"\"\"\n",
    "function mytrapezoids(fun, a, b, n)\n",
    "    h = (b - a)/(n - 1)\n",
    "    s1 = fun(a) + fun(b)\n",
    "    s2 = 0.0\n",
    "    for i = 2:(n-1)\n",
    "        s2 += fun(a + h*(i-1))\n",
    "    end\n",
    "    return h/2*(s1 + 2*s2)\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "209f238c-e1fa-4d12-86f6-ac9a49630a53",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "?mytrapezoids"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a9c89f13-53ef-4c21-951d-75006d16f43a",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "fun1(x) = 1/x\n",
    "a1 = 1.0\n",
    "b1 = 2.0\n",
    "exact1 = log(2)\n",
    "np1 = 20"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4ccdd2af-ce84-4c3e-8bb9-38a57b13aaa6",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "res1 = mytrapezoids(fun1, a1, b1, np1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "26aea22b-0436-4c18-9302-2564897f4d6a",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "round(abs(res1 - exact1), sigdigits=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fb8c7bd7-60be-4965-aa1e-d65c7b3d6273",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "fun2(x) = exp(x)\n",
    "a2 = 0.0\n",
    "b2 = 1.0\n",
    "exact2 = exp(1) - 1.0\n",
    "np2 = 100"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a248f450-480b-46c3-b787-4378e81b2730",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "res2 = mytrapezoids(fun2, a2, b2, np2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8bf0821d-22e5-4ec2-af6d-04434ef4cbfe",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "round(abs(res2 - exact2), sigdigits=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5d8e6517-3338-42c5-992f-73142ce6321a",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "fun3(x) = sin(x)\n",
    "a3 = 0.0\n",
    "b3 = pi\n",
    "exact3 = 2.0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "02a3d11f-8f8c-413a-814c-6952172271e9",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "ndp = 10\n",
    "hh = zeros(ndp)       # create an array of ndp elements\n",
    "abserrs = zeros(ndp);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "90d532ba-1353-4ee8-abec-dc1f53c83249",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "for i = 1:ndp\n",
    "    np = 2^i + 1\n",
    "    hh[i] = (b3 - a3)/(np - 1)\n",
    "    abserrs[i] = abs(mytrapezoids(fun3, a3, b3, np) - exact3)\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "19786643-2489-4e70-9baa-eb7b7e5181df",
   "metadata": {},
   "outputs": [],
   "source": [
    "abserrs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4a1727ad-bd43-4ab8-aa37-558ac982c2b1",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "using PyPlot"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a7404e10-9273-4bca-bd12-2a7291202920",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "loglog(hh, abserrs, marker=\".\", label=\"Numerical experiment\")\n",
    "loglog(hh, hh.^2, linestyle=\"dashed\", label=L\"h^2\")\n",
    "grid(true)\n",
    "xlabel(\"h\")\n",
    "ylabel(\"absolute error\")\n",
    "title(\"Error of trapezoidal formula\")\n",
    "legend();"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "08f4e6d0-0cfe-48a2-9b5d-d35511b88333",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "function mytrapezoids2(fv, h)\n",
    "    return h/2*(fv[1] + fv[end] + 2*sum(fv[2:(end-1)]))\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "247c78e9-2c7b-4432-963b-96135cb48634",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "x1 = range(a1, b1, np1)\n",
    "fv1 = fun1.(x1)\n",
    "h1 = x1[2] - x1[1]\n",
    "round(abs(mytrapezoids2(fv1, h1) - exact1), sigdigits=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c240c996-df28-42a9-9e25-7fcb3338ef11",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "x2 = range(a2, b2, np2)\n",
    "fv2 = fun2.(x2)\n",
    "h2 = x2[2] - x2[1]\n",
    "round(abs(mytrapezoids2(fv2, h2) - exact2), sigdigits=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "37e948ad-60bb-4161-b54f-0758ef087cce",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "function simpsons2(fun, a, b, n)\n",
    "    x = range(a, b, n)\n",
    "    fv = fun.(x)\n",
    "    h = x[2] - x[1]\n",
    "    th = mytrapezoids2(fv, h)\n",
    "    t2h = mytrapezoids2(fv[1:2:end], 2*h)\n",
    "    simp = 1/3*(4*th - t2h)\n",
    "    return simp\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "bbe4b2df-4047-436c-964c-ef404b306eb4",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "nps1 = 19\n",
    "round(abs(simpsons2(fun1, a1, b1, nps1) - exact1), sigdigits=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7ada2b52-15cc-43cd-8384-0a224f0a083a",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "nps2 = 99\n",
    "round(abs(simpsons2(fun2, a2, b2, nps2) - exact2), sigdigits=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3c6d5416-e6be-4a62-9523-08a89db8a906",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "nn = 500\n",
    "round(abs(mytrapezoids(fun1, a1, b1, nn) - exact1), sigdigits=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2b7051c4-799b-447b-a4ee-155257efb2ed",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "function test(n)\n",
    "    if iseven(n) || n <= 3\n",
    "        return Inf\n",
    "    end\n",
    "    return 0\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fe6762e9-36b3-4c61-9b79-068322530c9e",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "test(7)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b331817c-813e-42bf-97d2-ed55c6f03a88",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Julia 1.11.6",
   "language": "julia",
   "name": "julia-1.11"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "1.11.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}