3) Functions

3) Functions#

  1. Using Julia

  2. Introduction to floating point arithmetic
    2.1. Plotting

  3. Machine Epsilon

1. Julia#

To recap from first class, Julia is a relatively new programming language. Think of it as MATLAB done right, open source, and fast. It’s nominally general-purpose, but mostly for numerical/scientific/statistical computing. There are great learning resources. We’ll introduce concepts and language features as we go.

# The last line of a cell is output by default
x = 3
y = 4

4

println("$x + $y = $(x + y)") # string formatting/interpolation
4;  # trailing semicolon suppresses output

3 + 4 = 7

@show x + y
x * y

x + y = 7

12

1.1. Numbers#

Tip

Check this very helpful Julia documentation page.

3, 3.0, 3.0f0, big(3.0) # integers, double precision, single precision, and convert to a maximum precision representation

(3, 3.0, 3.0f0, 3.0)

typeof(3), typeof(3.0), typeof(3.0f0), typeof(big(3.0))

(Int64, Float64, Float32, BigFloat)

# automatic promotion
@show 3 + 3.0
@show 3.0 + 3.0f0
@show 3 + 3.0f0;

3 + 3.0 = 6.0
3.0 + 3.0f0 = 6.0
3 + 3.0f0 = 6.0f0

# floating and integer division
@show 4 / 2
@show -3 ÷ 2; # type `\div` and press TAB

4 / 2 = 2.0
-3 ÷ 2 = -1

1.2 Arrays#

[1, 2, 3]

3-element Vector{Int64}:
 1
 2
 3

# explicit typing
Float64[1,2,3]

3-element Vector{Float64}:
 1.0
 2.0
 3.0

# promotion rules similar to arithmetic
[1,2,3.] + [1,2,3]

3-element Vector{Float64}:
 2.0
 4.0
 6.0

x = [10., 20, 30]
x[2] # one-based indexing

20.0

x[2] = 3.5

3.5

# multi-dimensional array
A = [10 20 30; 40 50 60]

2×3 Matrix{Int64}:
 10  20  30
 40  50  60

Compare with the Python notation:

A = np.array([[10, 20, 30], [40, 50, 60]])

1.3 Functions#

# define a function called 'f' that takes three arguments, one of which has a default value
function f(x, y; z=3)
    sqrt(x*x + y*y) + z
end

# define another function also called 'f' that only takes two arguments
function f(x, y)
    sqrt(x*x + y*y)
end

f (generic function with 1 method)

What is the difference between a method and a function in Julia? Check this resource.

# if I invoke f this way, which one will the Julia compiler use?
f(3, 4, z=5)

10.0

# if instead I invoke f this way, which one will the Julia compiler use?
f(2, 3)

3.605551275463989

g(x, y) = sqrt(x^2 + y^2)
g(3, 4)

5.0

We have just seen an example of one of the most powerful features of Julia: dynamic (multiple) dispatch:

  • Functions can have multiple definitions as long each definition restricts the type of the parameters differently. It is the type of the parameters that define which “definition” (or “method” in Julia terminology) will be called.

  • This is the way Julia mimics polymorphism that some compiled languages (e.g., C++) have.

Anonymous functions are usually functions so short-lived that they do not need a name. This is done with an arrow notation. Example:

((x, y) -> sqrt(x^2 + y^2))(3, 4)

5.0

1.4 Loops#

# ranges
1:50000000

1:50000000

collect(1:5)

5-element Vector{Int64}:
 1
 2
 3
 4
 5

x = 0
for n in 1:50000000
    x += 1/n^2
end
@show x
x - pi^2/6 # Basel problem (solved by Euler; generalized by Riemann's zeta function)

x = 1.6449340467988642

-2.0049362170482254e-8

# list comprehensions
sum([1/n^2 for n in 1:1000])

1.6439345666815612

2. Floating point arithmetic#

0.1 + 0.2
(1 + 1.2e-16) - 1

2.220446049250313e-16

using Pkg
Pkg.add("Plots")

using Plots
default(linewidth=4)
plot(x -> (1 + x) - 1, xlims=(-1e-15, 1e-15),
xlabel="x", ylabel="y")
plot!(x -> x)

   Resolving package versions...

   Installed mtdev_jll ──────────────────── v1.1.7+0
   Installed libpng_jll ─────────────────── v1.6.49+0

   Installed Xorg_xcb_util_wm_jll ───────── v0.4.2+0
   Installed Xorg_xcb_util_image_jll ────── v0.4.1+0
   Installed Unitful ────────────────────── v1.23.1
   Installed Xorg_xcb_util_jll ──────────── v0.4.1+0
   Installed Xorg_xcb_util_keysyms_jll ──── v0.4.1+0

   Installed HTTP ───────────────────────── v1.10.17

   Installed libevdev_jll ───────────────── v1.13.4+0
   Installed ColorVectorSpace ───────────── v0.11.0
   Installed DocStringExtensions ────────── v0.9.5
   Installed libinput_jll ───────────────── v1.28.1+0

   Installed eudev_jll ──────────────────── v3.2.14+0
   Installed Scratch ────────────────────── v1.3.0
   Installed Wayland_protocols_jll ──────── v1.44.0+0
   Installed Plots ──────────────────────── v1.40.14

   Installed Xorg_xcb_util_renderutil_jll ─ v0.3.10+0
   Installed URIs ───────────────────────── v1.6.0

    Updating `~/.julia/environments/v1.10/Project.toml`
  [91a5bcdd] + Plots v1.40.14
    Updating `~/.julia/environments/v1.10/Manifest.toml`
 

 [66dad0bd] + AliasTables v1.1.3
  [d1d4a3ce] + BitFlags v0.1.9
  [944b1d66] + CodecZlib v0.7.8
  [35d6a980] + ColorSchemes v3.29.0
  [3da002f7] + ColorTypes v0.12.1
  [c3611d14] + ColorVectorSpace v0.11.0
  [5ae59095] + Colors v0.13.1
  [34da2185] + Compat v4.16.0
  [f0e56b4a] + ConcurrentUtilities v2.5.0
  [d38c429a] + Contour v0.6.3
  [9a962f9c] + DataAPI v1.16.0
  [864edb3b] + DataStructures v0.18.22
  [8bb1440f] + DelimitedFiles v1.9.1
  [ffbed154] + DocStringExtensions v0.9.5
  [460bff9d] + ExceptionUnwrapping v0.1.11
  [c87230d0] + FFMPEG v0.4.2
  [53c48c17] + FixedPointNumbers v0.8.5
  [1fa38f19] + Format v1.3.7
  [28b8d3ca] + GR v0.73.16
  [42e2da0e] + Grisu v1.0.2
  [cd3eb016] + HTTP v1.10.17
  [92d709cd] + IrrationalConstants v0.2.4
  [1019f520] + JLFzf v0.1.11
  [b964fa9f] + LaTeXStrings v1.4.0
  [23fbe1c1] + Latexify v0.16.8
  [2ab3a3ac] + LogExpFunctions v0.3.29
  [e6f89c97] + LoggingExtras v1.1.0
  [1914dd2f] + MacroTools v0.5.16
  [442fdcdd] + Measures v0.3.2
  [e1d29d7a] + Missings v1.2.0
  [77ba4419] + NaNMath v1.1.3
  [4d8831e6] + OpenSSL v1.5.0
  [bac558e1] + OrderedCollections v1.8.1
  [ccf2f8ad] + PlotThemes v3.3.0
  [995b91a9] + PlotUtils v1.4.3
  [91a5bcdd] + Plots v1.40.14
  [43287f4e] + PtrArrays v1.3.0
  [3cdcf5f2] + RecipesBase v1.3.4
  [01d81517] + RecipesPipeline v0.6.12
  [189a3867] + Reexport v1.2.2
  [05181044] + RelocatableFolders v1.0.1
  [ae029012] + Requires v1.3.1
  [6c6a2e73] + Scratch v1.3.0
  [992d4aef] + Showoff v1.0.3
  [777ac1f9] + SimpleBufferStream v1.2.0
  [a2af1166] + SortingAlgorithms v1.2.1
  [860ef19b] + StableRNGs v1.0.3
  [82ae8749] + StatsAPI v1.7.1
  [2913bbd2] + StatsBase v0.34.5
  [62fd8b95] + TensorCore v0.1.1
  [3bb67fe8] + TranscodingStreams v0.11.3
  [5c2747f8] + URIs v1.6.0
  [1cfade01] + UnicodeFun v0.4.1
  [1986cc42] + Unitful v1.23.1
  [45397f5d] + UnitfulLatexify v1.7.0
  [41fe7b60] + Unzip v0.2.0
  [6e34b625] + Bzip2_jll v1.0.9+0
  [83423d85] + Cairo_jll v1.18.5+0
  [ee1fde0b] + Dbus_jll v1.16.2+0
  [2702e6a9] + EpollShim_jll v0.0.20230411+1
  [2e619515] + Expat_jll v2.6.5+0
⌅ [b22a6f82] + FFMPEG_jll v4.4.4+1
  [a3f928ae] + Fontconfig_jll v2.16.0+0
  [d7e528f0] + FreeType2_jll v2.13.4+0
  [559328eb] + FriBidi_jll v1.0.17+0
  [0656b61e] + GLFW_jll v3.4.0+2
  [d2c73de3] + GR_jll v0.73.16+0
  [78b55507] + Gettext_jll v0.21.0+0
  [7746bdde] + Glib_jll v2.84.0+0
  [3b182d85] + Graphite2_jll v1.3.15+0
  [2e76f6c2] + HarfBuzz_jll v8.5.1+0
  [aacddb02] + JpegTurbo_jll v3.1.1+0
  [c1c5ebd0] + LAME_jll v3.100.2+0
  [88015f11] + LERC_jll v4.0.1+0
  [1d63c593] + LLVMOpenMP_jll v18.1.8+0
  [dd4b983a] + LZO_jll v2.10.3+0
  [e9f186c6] + Libffi_jll v3.4.7+0
  [7e76a0d4] + Libglvnd_jll v1.7.1+1
  [94ce4f54] + Libiconv_jll v1.18.0+0
  [4b2f31a3] + Libmount_jll v2.41.0+0
  [89763e89] + Libtiff_jll v4.7.1+0
  [38a345b3] + Libuuid_jll v2.41.0+0
  [e7412a2a] + Ogg_jll v1.3.5+1
  [458c3c95] + OpenSSL_jll v3.5.0+0
  [91d4177d] + Opus_jll v1.3.3+0
  [36c8627f] + Pango_jll v1.56.3+0
⌅ [30392449] + Pixman_jll v0.44.2+0
  [c0090381] + Qt6Base_jll v6.8.2+1
  [629bc702] + Qt6Declarative_jll v6.8.2+1
  [ce943373] + Qt6ShaderTools_jll v6.8.2+1
  [e99dba38] + Qt6Wayland_jll v6.8.2+0
  [a44049a8] + Vulkan_Loader_jll v1.3.243+0
  [a2964d1f] + Wayland_jll v1.23.1+0
  [2381bf8a] + Wayland_protocols_jll v1.44.0+0
⌅ [02c8fc9c] + XML2_jll v2.13.6+1
  [ffd25f8a] + XZ_jll v5.8.1+0
  [f67eecfb] + Xorg_libICE_jll v1.1.2+0
  [c834827a] + Xorg_libSM_jll v1.2.6+0
  [4f6342f7] + Xorg_libX11_jll v1.8.12+0
  [0c0b7dd1] + Xorg_libXau_jll v1.0.13+0
  [935fb764] + Xorg_libXcursor_jll v1.2.4+0
  [a3789734] + Xorg_libXdmcp_jll v1.1.6+0
  [1082639a] + Xorg_libXext_jll v1.3.7+0
  [d091e8ba] + Xorg_libXfixes_jll v6.0.1+0
  [a51aa0fd] + Xorg_libXi_jll v1.8.3+0
  [d1454406] + Xorg_libXinerama_jll v1.1.6+0
  [ec84b674] + Xorg_libXrandr_jll v1.5.5+0
  [ea2f1a96] + Xorg_libXrender_jll v0.9.12+0
  [c7cfdc94] + Xorg_libxcb_jll v1.17.1+0
  [cc61e674] + Xorg_libxkbfile_jll v1.1.3+0
  [e920d4aa] + Xorg_xcb_util_cursor_jll v0.1.4+0
  [12413925] + Xorg_xcb_util_image_jll v0.4.1+0
  [2def613f] + Xorg_xcb_util_jll v0.4.1+0
  [975044d2] + Xorg_xcb_util_keysyms_jll v0.4.1+0
  [0d47668e] + Xorg_xcb_util_renderutil_jll v0.3.10+0
  [c22f9ab0] + Xorg_xcb_util_wm_jll v0.4.2+0
  [35661453] + Xorg_xkbcomp_jll v1.4.7+0
  [33bec58e] + Xorg_xkeyboard_config_jll v2.44.0+0
  [c5fb5394] + Xorg_xtrans_jll v1.6.0+0
  [3161d3a3] + Zstd_jll v1.5.7+1
  [35ca27e7] + eudev_jll v3.2.14+0
  [214eeab7] + fzf_jll v0.61.1+0
  [a4ae2306] + libaom_jll v3.11.0+0
  [0ac62f75] + libass_jll v0.15.2+0
  [1183f4f0] + libdecor_jll v0.2.2+0
  [2db6ffa8] + libevdev_jll v1.13.4+0
  [f638f0a6] + libfdk_aac_jll v2.0.3+0
  [36db933b] + libinput_jll v1.28.1+0
  [b53b4c65] + libpng_jll v1.6.49+0
  [f27f6e37] + libvorbis_jll v1.3.7+2
  [009596ad] + mtdev_jll v1.1.7+0
⌅ [1270edf5] + x264_jll v2021.5.5+0
⌅ [dfaa095f] + x265_jll v3.5.0+0
  [d8fb68d0] + xkbcommon_jll v1.8.1+0
  [37e2e46d] + LinearAlgebra
  [2f01184e] + SparseArrays v1.10.0
  [10745b16] + Statistics v1.10.0
  [8dfed614] + Test
  [e66e0078] + CompilerSupportLibraries_jll v1.1.1+0
  [4536629a] + OpenBLAS_jll v0.3.23+4
  [05823500] + OpenLibm_jll v0.8.1+4
  [efcefdf7] + PCRE2_jll v10.42.0+1
  [bea87d4a] + SuiteSparse_jll v7.2.1+1
  [8e850b90] + libblastrampoline_jll v5.11.0+0
        Info Packages marked with ⌅ have new versions available but compatibility constraints restrict them from upgrading. To see why use `status --outdated -m`

Precompiling 

packages...

    422.1 ms  ✓ Scratch

    524.6 ms  ✓ eudev_jll
    529.5 ms  ✓ libpng_jll

    609.7 ms  ✓ libevdev_jll
    898.1 ms  ✓ DocStringExtensions
    652.7 ms  ✓ Wayland_protocols_jll
    661.8 ms  ✓ Xorg_xcb_util_jll

    734.0 ms  ✓ mtdev_jll

    884.5 ms  ✓ URIs
    461.6 ms  ✓ RelocatableFolders

    733.1 ms  ✓ Cairo_jll
    616.8 ms  ✓ Xorg_xcb_util_keysyms_jll

    559.0 ms  ✓ libinput_jll
    632.9 ms  ✓ Xorg_xcb_util_image_jll
    649.7 ms  ✓ xkbcommon_jll
    712.9 ms  ✓ LogExpFunctions
    666.9 ms  ✓ Xorg_xcb_util_wm_jll
    669.1 ms  ✓ Xorg_xcb_util_renderutil_jll

    516.8 ms  ✓ Vulkan_Loader_jll
    600.0 ms  ✓ HarfBuzz_jll

    556.4 ms  ✓ Xorg_xcb_util_cursor_jll

    584.6 ms  ✓ libass_jll
   2464.0 ms  ✓ ColorVectorSpace

    695.2 ms  ✓ Pango_jll

    834.9 ms  ✓ Qt6Base_jll

    587.6 ms  ✓ libdecor_jll
    726.1 ms  ✓ FFMPEG_jll

   1978.2 ms  ✓ StatsBase
    647.8 ms  ✓ Qt6ShaderTools_jll

    457.8 ms  ✓ FFMPEG

    572.6 ms  ✓ GLFW_jll

   4176.4 ms  ✓ Colors

    743.1 ms  ✓ GR_jll

   1445.4 ms  ✓ Qt6Declarative_jll

   2794.8 ms  ✓ ColorSchemes

   9562.6 ms  ✓ HTTP

3. Machine epsilon#

We approximate real numbers with floating point arithmetic, which can only represent discrete values. In particular, there exists a largest number, which we call \(\epsilon_{\text{machine}}\), such that

\[ 1 \oplus x = 1 \quad \text{for all}\ \lvert x \rvert < \epsilon_{\text{machine}}.\]

The notation \(\oplus, \ominus, \odot, \oslash\) represent the elementary operation carried out in floating point arithmetic.

eps = 1
while 1 + eps != 1
    eps = eps / 2
end
eps

1.1102230246251565e-16

eps = 1.f0
while 1 + eps != 1
    eps = eps / 2
end
eps

5.9604645f-8

3.1 Approximating exp#

\[e^x = 1 + x + \frac{x^2}{2!} + \frac{x^3}{3!} + \dotsb\]

Suppose we want to compute \(f(x) = e^x - 1\) for small values of \(x\).

f1(x) = exp(x) - 1
y1 = f1(1e-8)

9.99999993922529e-9

f2(x) = x + x^2/2 + x^3/6
y2 = f2(1e-8)

1.000000005e-8

Which answer is more accurate?

@show (y1 - y2)        # Absolute difference
@show (y1 - y2) / y2;  # Relative difference

y1 - y2 = -1.1077470910720506e-16
(y1 - y2) / y2 = -1.1077470855333152e-8

0 / 0

NaN
Contents