JanssonDimensions2023.lagda.tex

% -*-Latex-*-
\RequirePackage{amsmath}
\documentclass[aspectratio=169]{beamer}
\usetheme{Madrid}
% Hide navigation symbols
\setbeamertemplate{navigation symbols}{}
\usepackage[utf8]{inputenc}
\usepackage[T1]{fontenc}
\usepackage{graphicx}
\usepackage{hyperref}
\hypersetup{pdfpagemode={FullScreen}}
\usepackage{newunicodechar}
\input{newunicodedefs}
\author[Jansson]{Patrik Jansson}
\date{\today}
\institute[FP unit, Chalmers]{Functional Programming unit, Chalmers University of Technology}
\date{2023-02-28}
%TODO perhaps \titlegraphic{\vspace{-0.5cm}\includegraphics[width=0.6\textwidth]{}}
\title{Dimensional Analysis and Graded Algebras}
\hypersetup{
 pdfauthor={Patrik Jansson},
 pdftitle={Dimensional Analysis and Graded Algebras},
 pdfkeywords={},
 pdfsubject={},
 pdflang={English}}
\setlength{\parindent}{0pt}
\setlength{\parskip}{6pt plus 2pt minus 1pt}
%include lhs2TeX.fmt
%include lhs2TeX.sty
%include polycode.fmt
% %format / = "\mathbin{/}"
%format : = "\mathop{:}"
%format -> = "\,\to\,"
%format Type1 = Type "_1"
\newcommand{\un}{\mbox{\texttt{\_}}}   % short underscore
\newcommand{\unopun}[1]{\ensuremath{\,\un{#1}\un\,}}
%format _≡_ = "\unopun{" ≡ "}"
%format _==_ = "\unopun{" == "}"
%format ^ = "\mathbin{\!\hat{}\!}"
%format _*d_ = "\unopun{" *d "}"
%format *d = "\mathbin{{*}_d}"
%format _/d_ = "\unopun{" /d "}"
%format /d = "\mathbin{{/\!}_d}"
%format _^d_ = "\unopun{" ^d "}"
%format ^d = "\mathbin{\hat{\;}^d}"
%format _+q_ = "\unopun{" +q "}"
%format +q = "\mathbin{{+\!}_q}"
%format _*q_ = "\unopun{" *q "}"
%format *q = "\mathbin{{*}_q}"
%format _/q_ = "\unopun{" /q "}"
%format /q = "\mathbin{{/\!}_q}"
%format _^q_ = "\unopun{" ^q "}"
%format ^q = "\mathbin{\hat{\;}^q}"
%format scaleq = scale "_q"
%format _+v_ = "\unopun{" +v "}"
%format +v = "\mathbin{{+\!}_v}"
%format _-v_ = "\unopun{" -v "}"
%format -v = "\mathbin{{-\!}_v}"
%format _+s_ = "\unopun{" +s "}"
%format +s = "\mathbin{{+\!}_s}"
%format _*s_ = "\unopun{" *s "}"
%format *s = "\mathbin{{*\!}_s}"
%format _/s_ = "\unopun{" /s "}"
%format /s = "\mathbin{{/\!}_s}"
%format d1 = "d_1"
%format d2 = "d_2"
%format q0 = 0 "_{" q "}"
%format q1 = 1 "_{" q "}"
%format s0 = 0 "_{" s "}"
%format s1 = 1 "_{" s "}"
%format ℤ0 = 0 "_{" ℤ "}"
%format ℤ1 = 1 "_{" ℤ "}"
%format ℤ2 = 2 "_{" ℤ "}"
% the `doubleequals' macro is due to Jeremy Gibbons
\def\doubleequals{\mathrel{\unitlength 0.01em
  \begin{picture}(78,40)
    \put(7,34){\line(1,0){25}} \put(45,34){\line(1,0){25}}
    \put(7,14){\line(1,0){25}} \put(45,14){\line(1,0){25}}
  \end{picture}}}
%% If you remove the %format == command the lhs2TeX default yields ≡, which can be a problem
%format ==   = "\doubleequals"
\begin{document}
\begin{frame}
  \maketitle
  Joint work with Nicola Botta and Guilherme Silva.
\end{frame}
\begin{frame}
\frametitle{Dimension analysis and graded algebras: Intro}
\begin{itemize}
\item Talk at the \href{https://ifipwg21wiki.cs.kuleuven.be/IFIP21/OnlineFeb23}{2023-02 online meeting} of WG 2.1 on Algorithmic Languages and Calculi.
\item Abstract:
This talk describes work in progress aimed at understanding dimension analysis
though coding it up using dependent types (in Agda). I will talk you through
definitions of physical quantities, systems of units, dimensions, dimensional
analysis, and an example of applying it to modelling a simple pendulum
experiment.
\item Some history: I did an Agda tutorial at the WG2.1 meeting \#63 in Kyoto
(in 2007) and I found mentions of Agda in 8 other physical meetings after
that.
\item Source code: \url{https://github.com/DSLsofMath/dimensionanalysis}
\end{itemize}
%if allCode
\subsection{Agda version etc.}
\begin{itemize}
%Old: This file loads in Agda 2.6.2.2 with Agda stdlib-1.7
\item This file loads in Agda 2.6.4.1 with Agda stdlib-2.0
\item It is a literate Agda file with embedded LaTeX code.
\end{itemize}
%endif
%if allCode
\subsection{Skip some imports and setup}
\subsubsection{Basic imports: equality, integers}
\begin{code}
module JanssonDimensions2023 where
open Agda.Primitive using (lzero)
open import Relation.Binary.PropositionalEquality renaming (_≡_ to _==_)
open import Relation.Nullary renaming (¬_ to Not)
Type = Set
Type1 = Set1

import Data.Integer using (ℤ; +0; +_; -[1+_])
                    renaming (0ℤ to ℤ0; 1ℤ to ℤ1)
open Data.Integer
Integer : Type
Integer = ℤ
ℤ2 : ℤ
ℤ2 = + 2
\end{code}
\subsubsection{Lift |Ring| operations to vectors}
\begin{code}
import Algebra
Ring = Algebra.Ring lzero lzero
open import Data.Nat hiding (_^_; _/_) renaming (ℕ to Nat; _+_ to _+n_; _*_ to _*n_)
open import Data.Vec renaming (foldr to depFoldr)
foldr : ∀ {A : Type} {B : Type} {m} →
        (A → B → B) → B → Vec A m → B
foldr {A} {B} op = depFoldr (\ _ -> B) (\{m} -> op)
module VectScope (r : Ring) where

  A = Algebra.Ring.Carrier r
  open Algebra.Ring r

  0v : {n : Nat} -> Vec A n
  0v {n} = replicate n 0#

  _+v_  :  {n : Nat} ->
           Vec A n -> Vec A n -> Vec A n
  _+v_  =  zipWith _+_

  _-v_  :  {n : Nat} ->
           Vec A n -> Vec A n -> Vec A n
  _-v_  =  zipWith _-_

  _*v_  :  {n : Nat} ->
           A -> Vec A n -> Vec A n
  x *v v =  map (x *_) v
\end{code}
\subsubsection{Postulate a Ring of reals}
\begin{itemize}
\item implemented using |module Data.Float| and some postulates
\item \textbf{NOTE} - this is not true but used for pragmatic reasons
\end{itemize}
\begin{code}
import Data.Float using (Float)
Real : Type
Real = Data.Float.Float

import Data.Float.Base

module RealRingPostulates where
  Carrier = Real
  _≈_ = _==_
  _+_ = Data.Float.Base._+_
  _*_ = Data.Float.Base._*_
  -_  = Data.Float.Base.-_
  0#  = 0.0
  1#  = 1.0
  postulate isRing : Algebra.IsRing _≈_ _+_ _*_ -_ 0# 1#

RealRing : Ring
RealRing = record { RealRingPostulates }
open Data.Float.Base renaming (_*_ to _*r_; _+_ to _+r_; _-_ to _-r_; _÷_ to _/r_)
\end{code}
\subsubsection{Raw Field record type + Real field record}
\begin{code}
record RawField (S : Type) : Type where
  -- TODO redo with reuse of Ring structure?
  field
    s0 s1 : S
    _+s_ : S -> S -> S
    _*s_ : S -> S -> S
    _/s_ : S -> S -> S

  _^_ : S → Nat → S
  s ^ zero   = s1
  s ^ suc n  = (s ^ n) *s s

  pow : S -> ℤ -> S
  pow s  ( + n )  =        s ^ n
  pow s -[1+ n ]  = s1 /s (s ^ (suc n))

rfReal : RawField Real
rfReal = record {s0 = 0.0; s1 = 1.0; _+s_ = _+r_; _*s_ = _*r_; _/s_ = _/r_}

RealPlus3 : Type
RealPlus3 = Vec Real 3
\end{code}
%endif
\end{frame}
\begin{frame}
  \frametitle{Physical quantities and systems of units}
\begin{itemize}
\item We can assign a ``dimension'' to each physical quantity:
informally |dim : Q -> D|
\item Physical quantities like length, speed, force, etc. are usually
measured with respect to a system of units of measurement, one unit
per ``base dimension''.

\pause
\item These systems can be grouped into classes
\begin{itemize}
\item For geometry (the class L) we have just one base dimension of length.
\item For kinematics (the class LT) we have a length and a time.
\item For mechancics (the class LTM) we have length, time, and mass.
\item Etc. (current, temperature, etc.)
\end{itemize}
% The seven SI base units:
% | Quantity            | base unit |
% |---------------------+-----------|
% | time                | s         |
% | length              | m         |
% | mass                | kg        |
% | electric current    | A         |
% | temperature         | K         |
% | amount of substance | mol       |
% | luminous intensity  | cd        |
\pause
\item To make the most of dimension analysis, stay in as small a class as you can.
\end{itemize}
\end{frame}
\begin{frame}
\frametitle{Dimensions and dimension functions}
\begin{itemize}
\item Informal example: acceleration ``has dimension \(L /
  T^2\)''. Formal reading: the \textbf{dimension function}
  is |\ L T M -> L / T ^ 2|. The function describes how the measured
  value of the acceleration changes when we change the units of
  measurements.
\pause
\item My height (of dimension L) is measured to 1.78 in meters but 178 in
  cm.
\item In general, if we make the unit of measurement L times smaller,
  we make the measured height L times bigger.
\item The height is actually invariant, but the measured value changes
  in the opposite direction of the measuring rod.
\pause
\item This simplest (linear scaling) relation holds for quantities of
  one of the base dimensions, but changes if we move to derived
  dimensions like that for area or acceleration.
\item In general, if we make the unit of length |L| times smaller, the
  unit of time |T| times smaller and the unit of mass |M| times
  smaller, the measuread acceleration becomes |L / T ^ 2| times bigger
  and the measured force |M * L / T ^ 2| times bigger.
\end{itemize}
\end{frame}
\begin{frame}
  \frametitle{Physics and dimensions}
\begin{itemize}
\item In an equation like |F = m * a| %(force equals mass times acceleration)
  in physics, the dimensions must match:
\begin{spec}
dim F = dim (m * a)
\end{spec}
\item Similarly for addition: the dimensions of the operands must match.
\pause
\item For multiplication we don't need to require matching dimension:
  |dim| is a homomorphism: |dim (m * a) = dim m *d dim a|.
\item Now, how do we type this? We clearly need to be able to multiply
  and divide dimensions, and we also need a value to denote
  ``dimensionless''.
\end{itemize}
\end{frame}
\begin{frame}
\frametitle{DimStuff}
We can capture the signature of what we need in a record type |DimStuff|.
\begin{code}
record DimStuff : Type1 where
  field
    Dim : Type
    Dimless : Dim
    _*d_  : Dim -> Dim -> Dim
    _/d_  : Dim -> Dim -> Dim

  _^d_ : Dim -> ℤ -> Dim
\end{code}
(We will get back to laws later, and to concrete implementations.)
%if allCode
\subsubsection{Skip details (exponentiation)}
\begin{code}
  _^dn_ : Dim -> Nat -> Dim
  d ^dn zero   = Dimless
  d ^dn suc n  = (d ^dn n) *d d

  d ^d  ( + n )  = d ^dn n
  d ^d -[1+ n ]  = Dimless /d (d ^dn suc n)
\end{code}
%endif
\end{frame}
\begin{frame}
\frametitle{Quantities indexed by dimensions}
Then we need to make sure the type for ``quantities'' is indexed by
dimensions. We assume a type for dimensions (|Dim|), and a type of
scalars |S|.
\begin{code}
module Quantities (dimstuff : DimStuff) (S : Type) where
  open DimStuff dimstuff

  variable d d1 d2 : Dim
  record QStuff : Type1 where
    field
      Q     : Dim -> Type
      dim   : Q d -> Dim

      q0      : Q d
      _+q_    : Q d  -> Q d  -> Q d
      scaleq  : S    -> Q d  -> Q d
\end{code}
\end{frame}
\begin{frame}
\frametitle{|QStuff| continued: a few more operations}
Quantities form almost a field:
\begin{code}
      q1    : Q Dimless
      _*q_  : Q d1 -> Q d2 -> Q (d1 *d d2)
      _/q_  : Q d1 -> Q d2 -> Q (d1 /d d2)
\end{code}

Here we see that the value of the dimension index tracks the
operations performed on the quantities.

%if allCode
\subsubsection{Skip (fixities, exponentiation)}
\begin{code}
    infixl 7 _*q_
    infixl 7 _/q_
    infixl 6 _+q_

    _^qn_ : Q d -> (n : Nat) -> Q (d ^dn n)
    s ^qn zero   = q1
    s ^qn suc n  = (s ^qn n) *q s

    _^q_ : Q d -> (i : ℤ) -> Q (d ^d i)
    s ^q  ( + n )  = s ^qn n
    s ^q -[1+ n ]  = q1 /q (s ^qn (suc n))
  open QStuff
\end{code}
%endif
\end{frame}
\begin{frame}
\frametitle{Dimension structure}
The algebraic structures at play for the dimensions part is a group.
And we can use a vector of integers to keep track of the exponents of
the base dimensions:

\frametitle{Concrete example: the LTM class}
%if allCode
\subsubsection{Skip module header / imports}
\begin{code}
module LTM (S : Type) (rf : RawField S) where
  open RawField rf public
  module LTMDim where
    open import Data.Integer.Properties using (+-*-ring) public
    open module V = VectScope +-*-ring public
\end{code}
%endif

%\frametitle{LTM Dim implementation example}
\begin{code}
    Dim = Vec ℤ 3
    _*d_ = _+v_
    _/d_ = _-v_
    Dimless  : Dim
    Dimless  = ℤ0 ∷ ℤ0 ∷ ℤ0 ∷ []  -- zero vector
    L T M    : Dim                -- base vectors / dimensions
    L        = ℤ1 ∷ ℤ0 ∷ ℤ0 ∷ []
    T        = ℤ0 ∷ ℤ1 ∷ ℤ0 ∷ []
    M        = ℤ0 ∷ ℤ0 ∷ ℤ1 ∷ []
\end{code}
\end{frame}
\begin{frame}
\frametitle{Implementing physical quantities}
%if allCode
\subsubsection{Skip details (module end + opens)}
\begin{code}
  dimstuff : DimStuff
  dimstuff = record { LTMDim }

  open Quantities dimstuff S public
  open LTMDim
\end{code}
%endif
The |Q| record type: just a wrapper around a scalar value:
\begin{code}
  module LTMQ where
    record Q (d : Dim) : Type where
      constructor Val
      field
        val : S
    open Q public
    dim : Q d -> Dim
    dim {d} _ = d
\end{code}
\end{frame}
\begin{frame}
\frametitle{|Q d| is a 1D vector space for each |d : Dim|}
\begin{code}
    q0 : Q d
    q0 = Val s0

    _+q_ : Q d  -> Q d  -> Q d
    Val x +q Val y = Val (x +s y)

    scaleq : S -> Q d -> Q d
    scaleq s (Val x) = Val (s *s x)
\end{code}

\pause
\ldots and |Q| is a graded field (informally - no proofs here)
\begin{code}
    q1  : Q Dimless
    q1 = Val s1
    _*q_ : Q d1 -> Q d2 -> Q (d1 *d d2)
    Val x *q Val y = Val (x *s y)
    _/q_ : Q d1 -> Q d2 -> Q (d1 /d d2)
    Val x /q Val y = Val (x /s y)
\end{code}
%if allCode
\subsubsection{Skip (module end, open)}
\begin{code}
  qstuff : QStuff
  qstuff = record { LTMQ }
  open QStuff qstuff
\end{code}
%endif

\end{frame}
\begin{frame}
\frametitle{Discussion: field or not a field?}

For the scalars we need an underlying field |S| and we get a field
back for |Q Dimless| but not for |Q d| for other |d|. But we get an
example of an indexed structure: an indexed, or graded, field.

\pause
Example of why |Q L| is not a field: the multiplication operation
does not stay within the type:

\begin{code}
  mul : Q L -> Q L -> Q (L *d L) -- |Q Area|
  mul x y = x *q y
\end{code}

\end{frame}
\begin{frame}
\frametitle{Measuring quantities}

There is one operation missing from |Q d|: measuring the quantity in a
particular system of units. A system of units consists of one scale
factor for each base dimension.
\subsubsection{System of units}
\begin{code}
  SysUnit : Type
  SysUnit = Vec S 3   -- for the |LTM| class with |n = 3|
\end{code}
We will later implement these two examples:
\begin{spec}
  si   : SysUnit  -- m, kg, second
  cgs  : SysUnit  -- cm, gram, second
\end{spec}

Now we can interpret a ``syntactic dimension'' (a vector of exponents of
the base dimensions) as its ``dimension function'' - the corresponding
monomial which computes the change in measured value for a change in
units.
\end{frame}

\begin{frame}
\frametitle{Dimension function and |measure| implemented}
\begin{code}
  dimfun : Dim -> (SysUnit -> S)
  dimfun d su = foldr _*s_ s1 (zipWith pow su d)
\end{code}

\pause

Finally we can define the |measure| of a quantity in a system of
units:

%{
%format LTMQ.val = val
\begin{code}
  measure : Q d -> SysUnit -> S
  measure q su = dimfun (dim q) su *s (LTMQ.val q)
\end{code}
(Remember that |val : Q d -> S| extracts the scalar inside a quantity.)
%}
\end{frame}
\begin{frame}
\frametitle{Examples: systems of units and quantities}
%if allCode
\subsubsection{Skip module / opens}
\begin{code}
module Examples where
  open LTM Real rfReal
  open LTMDim using (L; T; M)
  open DimStuff dimstuff
  open QStuff qstuff using (_^q_)
  open LTMQ
\end{code}
%endif
%\subsubsection{Length and Acceleration}
\begin{columns}
\begin{column}{0.49\textwidth}
\begin{code}
  si   : SysUnit   -- m,  second, kg
  si   =   1.0 ∷ 1.0 ∷    1.0 ∷ []
  cgs  : SysUnit   -- cm, second, gram
  cgs  = 100.0 ∷ 1.0 ∷ 1000.0 ∷ []
\end{code}
\pause
\vspace*{-1cm}
\begin{code}
  mass   : Q M
  mass   = Val 76.0
  hei    : Q L
  hei    = Val 1.78

  test1   : Real
  test1   = measure hei cgs
  check1  : test1 == 178.0
  check1  = refl
\end{code}
\end{column}
\begin{column}{0.49\textwidth}
\pause
% \frametitle{Mass and Force}
\begin{code}
  Acceleration : Dim
  Acceleration = L /d (T ^d ℤ2)
  g : Q Acceleration
  g = Val 9.82
  Force  : Dim
  Force  = M *d Acceleration
  f  : Q Force
  f  = mass *q g
\end{code}
\pause
Unit conversion (measure)
\begin{code}
  test2   : Real
  test2   = measure f cgs
  check2  : test2 == 76000.0 *s 982.0
  check2  = refl
\end{code}
\end{column}
\end{columns}
\end{frame}
\begin{frame}
\frametitle{Dimension analysis (Pi theorem example)}
OK, now we have dimensions, quantities, the dimension function and
some examples. Time to get to the core of dimension analysis: the Pi
theorem.

\subsubsection{Pendulum example intro}
\begin{itemize}
\item Assume we are experimenting with an ideal pendulum:
  \begin{itemize}
  \item a point mass |m| hanging from
  \item a piece of string of length |x|,
  \item the time |t| for one period
  \item when starting from an angle |theta|.
  \end{itemize}
\item We want to find how the gravitational constant |g| can be computed
from the other four parameters.
\end{itemize}

\begin{spec}
  gravity : Q M -> Q L -> Q T -> Q Dimless -> Q Acceleration
  gravity m x t theta = ?
\end{spec}

\end{frame}
\begin{frame}
\frametitle{Pi theorem for this example}
\begin{itemize}
\item The core theorem of dimension analysis says that such a function
  can be simplified significantly.
\item First, pick three quantities of independent dimension: here |m|,
  |x|, |t|.
\pause
\item Second, make the remaining quantities dimensionless by combining
  them with monomials in the first three.
\item Here |theta| is already dimensionless, but |g = gravity m x
  theta t| is replaced by |a = g *q (t ^q ℤ2) /q x : Q Dimless|.
\pause
\item Finally, the function |gravity| is factored into a function
  |acc| from just |theta| to |a|, and a monomial factor used to
  translate back to a quantity of the correct dimension.
\end{itemize}
\begin{code}
  acc : Q Dimless -> Q Dimless
  gravity : Q M -> Q L -> Q T -> Q Dimless -> Q Acceleration
  gravity m x t theta = monomial  *q  acc theta
    where  monomial : Q Acceleration
           monomial = x /q (t ^q ℤ2) -- The type dictates the value.
\end{code}
\end{frame}
\begin{frame}
\frametitle{Pi theorem implications}
\begin{itemize}
\item If we want to figure out the \textbf{4-argument} function
  |gravity| from experiments, we actually only need to figure out the
  \textbf{1-argument} function |acc|.
\item Experiments (or numerical simulations with the proper
  differential equations) further show that the remaining function is
  constant for small angles |theta|.
\end{itemize}

\begin{code}
  twoPiSq : Q Dimless
  acc theta = twoPiSq        -- for small |theta|

  twoPiSq = scaleq 39.5 q1
\end{code}

\end{frame}
\begin{frame}
\frametitle{Summary \& future work}
\begin{itemize}
\item Dimension analysis can be usefully expressed done with dep. types
\item Physical quantities form a field graded by an abelian group of
  dimensions.
\item Other example: tensors are graded by their order
\item Library code and paper is work in progress (mostly in Idris).
\item Multiple grading: |T ord (Q dim s)| or |Q dim (T ord s)|?
  \begin{itemize}
  \item tensors with rank containing quantities with dimensions?
  \item or quantities with dimensions containing tensors with rank?
  \item or \ldots{}
  \end{itemize}
\end{itemize}
\end{frame}
\end{document}

Some experiments with numeric literals.
\begin{code}
module NumLitExperiment where
  -- Natural literals work "out of the box"
  n : Nat
  n = 3

  -- And they can be used with the (smart) constructor +_ to build
  -- integers.
  z5 : ℤ
  z5 = + 5

  -- With the "FromNat" builtin, the number literals can also be used
  -- for integers directly.
  open import Agda.Builtin.FromNat
  import Data.Integer.Literals
  import Data.Unit -- needed for the (trivial) constraint
  instance
    _ = Data.Integer.Literals.number

  z6 : ℤ
  z6 = 6

  -- And we can also do the same for (embed naturals as) rationals:
  open import Data.Rational

  import Data.Rational.Literals
  instance
    _ = Data.Rational.Literals.number

  q3 : ℚ
  q3 = 3

  -- To express proper rationals we also need the FromNat machinery
  -- for the naturals (I'm not sure why).
  import Data.Nat.Literals
  instance
    _ = Data.Nat.Literals.number

  pointnine : ℚ
  pointnine = 9 / 10
\end{code}