Merge pull request #127 from jverzani/v0.19a

V0.19a

quarto/_quarto.yml

@@ -1,4 +1,4 @@
-version: "0.18"
+version: "0.19"
 engine: julia
 
 project:
@@ -84,6 +84,7 @@ book:
         - integrals/volumes_slice.qmd
         - integrals/arc_length.qmd
         - integrals/surface_area.qmd
+        - integrals/twelve-qs.qmd
 
     - part: ODEs.qmd
       chapters:

quarto/derivatives/first_second_derivatives.qmd

@@ -49,7 +49,7 @@ A parallel definition with $a < b$ implying $f(a) > f(b)$ would be used for a *s
 We can try and prove these properties for a function algebraically – we'll see both are related to the zeros of some function. However, before proceeding to that it is usually helpful to get an idea of where the answer is using exploratory graphs.
 
 
-We will use a helper function, `plotif(f, g, a, b)` that plots the function `f` over `[a,b]` coloring it red when `g` is positive (and blue otherwise). Such a function is defined for us in the accompanying `CalculusWithJulia` package, which has been previously been loaded.
+We will use a helper function, `plotif(f, g, a, b)` that plots the function `f` over `[a,b]` highlighting the regions in the domain when `g` is non-negative. Such a function is defined for us in the accompanying `CalculusWithJulia` package, which has been previously been loaded.
 
 
 To see where a function is positive, we simply pass the function object in for *both* `f` and `g` above. For example, let's look at where $f(x) = \sin(x)$ is positive:
@@ -106,7 +106,7 @@ The mean value theorem provides the reasoning behind the first statement: on $I$
 The second part, follows from the secant line equation. The derivative can be written as a limit of secant-line slopes, each of which is positive. The limit of positive things can only be non-negative, though there is no guarantee the limit will be positive.
 
 
-So, to visualize where a function is increasing, we can just pass in the derivative as the masking function in our `plotif` function, as long as we are wary about places with $0$ derivative (flat spots).
+So, to visualize where a function is increasing or flat, we can just pass in the derivative as the masking function in our `plotif` function.
 
 
 For example, here, with a more complicated function, the intervals where the function is increasing are highlighted by passing in the functions derivative to `plotif`:

quarto/integrals/center_of_mass.qmd

@@ -344,6 +344,18 @@ Here the center of mass is below $1/2$ as the bulk of the area is. (The bottom a
 As seen, the computation of the center of mass in the $y$ direction has an identical formula, though may be more involved if an inverse function must be computed.
 
 
+::: {.callout-note}
+#### An alternative formula
+
+An alternative formula, which is easily derived once double integrals are introduced, to find the center of mass in the $y$ direction is
+
+$$
+\text{cm}_y = \frac{\int_a^b \frac{1}{2}(f(x)^2 - g(x)^2) dx}{\int_a^b (f(x) -g(x)) dx}.
+$$
+
+:::
+
+
 ##### Example
 
 

quarto/integrals/twelve-qs.qmd

@@ -0,0 +1,431 @@
+#  A dozen minima for a parabola
+
+This section uses these packages:
+
+```{julia}
+using SymPy
+using Plots
+```
+
+----
+
+In the March 2003 issue of the College Mathematics Journal, Leon M Hall posed 12 questions related to the following figure:
+
+```{julia}
+#| echo: false
+f(x) = x^2
+fp(x) = 2x
+a₀ = 7/8
+q₀ = -a₀ - 1/(2a₀)
+tangent(x) = f(a₀) + fp(a₀) * (x - a₀)
+normal(x) = f(a₀) - (1 / fp(a₀)) * (x - a₀)
+function make_plot(a₀=7/8, q₀=-a₀ - 1/2a₀)
+    plt = plot(; xlim=(-2,2), ylim=(-1, (1.5)^2),
+               xticks=nothing, yticks=nothing,
+               aspect_ratio=:equal, border=:none, legend=false)
+    f(x) = x^2
+    fp(x) = 2x
+    plot!(f, -1.5, 1.5)
+    plot!(zero)
+
+    tl = x -> f(a₀) + fp(a₀) * (x-a₀)
+    nl = x -> f(a₀) - 1/(fp(a₀)) * (x-a₀)
+    plot!(tl, -0.02, 1.6; linecolor=:black)
+    plot!(nl, -1.6, 1; linecolor=:black)
+    # add in right triangle
+    scatter!([a₀, q₀], f.([a₀, q₀]), markersize=5)
+    Δ = 0.01
+    annotate!([(a₀ + Δ, nl(a₀+Δ), "P", :bottom),
+               (q₀ - Δ, nl(q₀-Δ), "Q", :top)])
+    plt
+end
+make_plot()
+```
+
+The figure shows $f(x) = x^2$, the tangent line at $(a, f(a))$ (for $a > 0$), and the *normal* line at $(a, f(a))$. The  questions all involve finding the value $a$ which minimizes a related quantity.
+
+
+We set up some variables to work symbolically:
+
+```{julia}
+@syms a::positive x::real
+f(x) = x^2
+fp(x) = 2x
+m = fp(a)
+mᴵ = - 1/m
+tl = f(a) + m * (x - a)
+nl = f(a) + mᴵ * (x - a)
+zs = solve(f(x) ~ nl, x)
+q = only(filter(!=(a), zs))
+```
+
+----
+
+The first question is simply:
+
+> 1a. The $y$ coordinate of $Q$
+
+```{julia}
+#| echo: false
+let
+    p = make_plot()
+    plot!(p, [q₀,q₀], [0,normal(q₀)], linewidth=5)
+end
+```
+
+The value is $f(q)$
+
+```{julia}
+yvalue = f(q)
+```
+
+To minimize we solve for critical points:
+
+```{julia}
+cps = solve(diff(yvalue, a), a)
+```
+
+The lone critical point must be at a minimum. (Given the geometry of the problem, as $a$ goes to $\infty$ the height does too, and as $a$ goes to $0$ the height will also go to $\infty$. This can also be seen analytically, as $q = -a - 1/(2a)$ which goes to $-\infty$ when $a$ heads to $0$ or $\infty$.)
+
+::: {.callout-note}
+## We hide the code
+
+In the remaining examples we don't show the code by default.
+
+:::
+
+----
+
+> 1b. The length of the line segment $PQ$
+
+```{julia}
+#| echo: false
+p = make_plot()
+plot!([q₀, a₀], [f(q₀), f(a₀)], linewidth=5)
+```
+
+
+```{julia}
+#| code-fold: true
+#| code-summary: "Show the code"
+lseg = sqrt((f(a) - f(q))^2 + (a - q)^2);
+```
+
+
+----
+
+> 2a. The horizontal distance between $P$ and $Q$
+
+
+```{julia}
+#| echo: false
+p = make_plot()
+plot!([q₀, a₀], [f(a₀), f(a₀)], linewidth=5)
+```
+
+
+```{julia}
+#| code-fold: true
+#| code-summary: "Show the code"
+
+hd = a - q;
+```
+----
+
+> 2b. The area of the parabolic segment
+
+
+```{julia}
+#| echo: false
+p = make_plot()
+xs′ = range(q₀, a₀, 50)
+xs = vcat(xs′, reverse(xs′), first(xs′))
+ys = vcat(f.(xs′), normal.(reverse(xs′)), f(first(xs′)))
+plot!(xs, ys, fill=(:green, 0.25, 0))
+```
+
+
+```{julia}
+#| code-fold: true
+#| code-summary: "Show the code"
+A = simplify(integrate(nl - f(x), (x, q, a)));
+```
+
+----
+
+> 2c. The volume of the rotated solid formed by revolving the parabolic segment around the vertical line $k$ units to the right of $P$ or to the left of $Q$ where $k > 0$.
+
+
+
+```{julia}
+#| code-fold: true
+#| code-summary: "Show the code"
+@syms k::nonnegative
+V = simplify(integrate(PI * (nl - f(x) - k)^2, (x, q, a)));
+```
+
+----
+
+> 3. The $y$ coordinate of the centroid of the parabolic segment
+
+
+```{julia}
+#| echo: false
+p = make_plot()
+scatter!(p, [-1/(4a₀)], [1], marker=(10, :diamond))
+p
+```
+
+
+We warm up with the $x$ coordinate, given by:
+
+```{julia}
+xₘ = integrate(x * (nl - f(x)), (x, q, a)) / A
+simplify(xₘ)
+```
+
+a fact noted by the author.
+
+```{julia}
+#| code-fold: true
+#| code-summary: "Show the code"
+yₘ = integrate( (1//2) * (nl^2 - f(x)^2), (x, q, a)) / A
+yₘ = simplify(yₘ);
+```
+
+----
+
+> 4. The length of the arc of the parabola between $P$ and $Q$
+
+
+```{julia}
+#| echo: false
+p = make_plot()
+xs = range(q₀, a₀, 100)
+plot!(xs, f.(xs), linewidth=5)
+p
+```
+
+```{julia}
+#| code-fold: true
+#| code-summary: "Show the code"
+L = integrate(sqrt(1 + fp(x)^2), (x, q, a));
+```
+
+----
+
+> 5. The $y$ coordinate of the midpoint ofthe line segment $PQ$
+
+
+```{julia}
+#| echo: false
+p = make_plot()
+mid = a₀/2 + q₀/2
+vline!([mid])
+scatter!([mid], [normal(mid)], markersize=5)
+p
+```
+
+
+```{julia}
+#| code-fold: true
+#| code-summary: "Show the code"
+mp = nl(x => (a + q)/2);
+```
+
+----
+
+> 6. The area of the trapezoid bound by the normal line, the $x$-axis, and the vertical lines through $P$ and $Q$.
+
+
+```{julia}
+#| echo: false
+p = make_plot()
+plot!([q₀, a₀, a₀, q₀, q₀],
+      [0,  0, f(a₀), f(q₀), 0]; fill=(:green, 0.25, 0))
+p
+```
+
+```{julia}
+#| code-fold: true
+#| code-summary: "Show the code"
+trap = 1//2 * (f(q) + f(a)) * (a - q);
+```
+
+----
+
+> 7. The area bounded by the parabola and the $x$ axis and the vertical lines through $P$ and $Q$
+
+```{julia}
+#| echo: false
+p = make_plot()
+ts = range(q₀, a₀, 100)
+xs = vcat(ts, a₀, q₀, q₀)
+ys = vcat(f.(ts), 0, 0, f(q₀))
+plot!(xs, ys, fill=(:green, 0.25, 0))
+p
+```
+
+
+
+```{julia}
+#| code-fold: true
+#| code-summary: "Show the code"
+pa = integrate(x^2, (x, q, a));
+```
+
+----
+
+> 8. The area of the surface formed by revolving the arc of the parabola between $P$ and $Q$ around the vertical line through $P$
+
+```{julia}
+#| echo: false
+let
+    using SplitApplyCombine
+    S(x, θ) = (a₀ - (a₀-x)*cos(θ), (a₀-x)*sin(θ), f(x))
+    xs = range(q₀, a₀, 50)
+    θs = range(0, 2pi, 50)
+    surface(SplitApplyCombine.invert(S.(xs,θs'))...)
+end
+```
+
+
+
+```{julia}
+#| code-fold: true
+#| code-summary: "Show the code"
+# use parametric and  2π ∫ u(t) √(u'(t)^2 + v'(t)^2) dt
+uu(x) = a - x
+vv(x) = f(uu(x))
+SA = 2PI * integrate(uu(x) * sqrt(diff(uu(x),x)^2 + diff(vv(x),x)^2), (x, q, a));
+```
+
+----
+
+> 9. The height of the parabolic segment (i.e. the distance between the normal line and the tangent line to the parabola that is parallel to the normal line)
+
+```{julia}
+#| echo: false
+# distance point to a line
+# > solve(diff(x^2,x) ~ -1/fp(a), x) # mean value theorem
+#  -1/(4*a)
+
+b₀ = -1/(4a₀) # mean value approach
+# f(b0) + fp(a0)*(x-b0) = f(a0) - 1/fp(a0)*(x-a0)
+x₀ =  a₀/2 - 1/(8*a₀)
+make_plot()
+plot!(x -> f(b₀) + (-1/fp(a₀))*(x - b₀), -1, 1/2)
+plot!([b₀,x₀], [f(b₀), normal(x₀)]; linewidth=5)
+```
+
+```{julia}
+#| code-fold: true
+#| code-summary: "Show the code"
+# find b through mean value theorem,
+# then solve for point of intersection
+b = only(solve(diff(f(x),x) ~ -1/fp(a), x))
+b′ = only(solve(f(b) + fp(a)*(x-b) ~ nl, x))
+segment_height = sqrt((b-b′)^2 + (f(b) - nl(x=>b′))^2);
+```
+
+
+----
+
+> 10. The volume of the solid formed by revolving the parabolic segment around the $x$-axis
+
+```{julia}
+#| echo: false
+let
+    using SplitApplyCombine
+    S(x,θ) = (x, ((f(a₀) - 1/fp(a₀)*(x - a₀)) .* (sin(θ), cos(θ)))...)
+    xs = range(q₀, a₀, 50)
+    θs = range(0, 2pi, 50)
+    surface(SplitApplyCombine.invert(S.(xs,θs'))...)
+end
+```
+
+```{julia}
+#| code-fold: true
+#| code-summary: "Show the code"
+Vₓ = integrate(pi * (nl^2 - f(x)^2), (x, q, a));
+```
+
+----
+
+> 11. The area of the triangle bound by the normal line, the vertical line through $Q$ and the $x$-axis
+
+```{julia}
+#| echo: false
+make_plot()
+# solve 0 = y0 + m * (x-x0) --> x0 - y0/m
+p₀ = a₀ - f(a₀) / (-1/fp(a₀))
+xlims!((-2,  p₀ + 0.2))
+plot!([p₀,q₀,q₀,p₀], [0,f(q₀),0,0];
+      fill=(:green, 0.25,0))
+```
+
+```{julia}
+#| code-fold: true
+#| code-summary: "Show the code"
+
+triangle = 1/2 * f(q) * (a - f(a)/(-1/fp(a)) - q);
+```
+
+----
+
+> 12. The area of the quadrilateral bound by the normal line, the tangent line, the vertical line through $Q$ and the $x$-axis
+
+```{julia}
+#| echo: false
+make_plot()
+plot!([a₀,q₀,q₀,a₀-f(a₀)/fp(a₀),a₀],
+      [f(a₀), f(q₀), 0, 0,f(a₀)], fill=(:green, 0.25, 0))
+```
+
+
+```{julia}
+#| code-fold: true
+#| code-summary: "Show the code"
+# @syms x[1:4], y[1:4]
+# v1, v2, v3 = [[x[i]-x[1],y[i]-y[1], 0] for i in 2:4]
+# area = 1//2 * last(cross(v3,v2) + cross(v2, v1)) # 1/2 area of parallelogram
+# print(simplify(area))
+# -(x₁ - x₂)*(y₁ - y₃)/2 + (x₁ - x₃)*(y₁ - y₂)/2 - (x₁ - x₃)*(y₁ - y₄)/2 + (x₁ - x₄)*(y₁ - y₃)/2
+tl₀ = a - f(a) / fp(a)
+x₁,x₂,x₃,x₄ = (a,q,q,tl₀)
+y₁, y₂, y₃, y₄ = (f(a), f(q), 0, 0)
+quadrilateral = -(x₁ - x₂)*(y₁ - y₃)/2 + (x₁ - x₃)*(y₁ - y₂)/2 - (x₁ - x₃)*(y₁ - y₄)/2 + (x₁ - x₄)*(y₁ - y₃)/2;
+```
+----
+
+The answers appear here in sorted order, some given as approximate floating point values:
+
+
+```{julia}
+article_answers = (1/(2sqrt(2)), 1/2, sqrt(3/10), 0.558480, 0.564641,
+                   0.569723, 0.574646,
+                   1/sqrt(3), 1/8^(1/4), 1/6^(1/4), .644004, 1/sqrt(2))
+```
+
+```{julia}
+#| echo: false
+# check
+problems = ("1a"=>yvalue, "1b"=>lseg, "1c"=>hd,
+            "2a" => A, "2b" => V,
+            "3" => yₘ,
+            "4" => L,
+            "5" => mp,
+            "6" => trap,
+            "7" => pa,
+            "8" => SA,
+            "9" => segment_height,
+            "10" => Vₓ,
+            "11" => triangle,
+            "12" => quadrilateral
+            )
+≈ₒ(a,b) = isapprox(a, b; atol=1e-5, rtol=sqrt(eps()))
+∂ = Differential(a)
+solutions = [k => only(solve(∂(p) ~ 0, a)) for (k,p) in problems]
+[(sol=k, correct=(any(isapprox.(s, article_answers; atol=1e-5)))) for (k,s) ∈ solutions]
+nothing
+```

quarto/limits/limits.qmd

@@ -528,7 +528,7 @@ For example:
 lim(f, c)
 ```
 
-The numbers are displayed in decreasing order so the values on the left are read from top to bottom:
+The numbers are displayed in decreasing order so the values on the left side of $c$ are read from bottom to top:
 
 ```{julia}
 lim(f, c; dir="-")  # or even lim(f, c, -)