# Syntax: Kritik und Vorschläge

Wir schreiben hier unsere Gedanken zum Vorschlag [syntax_wishlist_OJ.jl](https://gitlab.lrz.de/software1/gaio.jl/-/blob/818be6e95bf7573593a36913ff9edcf130fec7dd/protocols/syntax_wishlist_OJ.jl) auf.

## relative_attractor!

function relative_attractor!(B::BoxSet, f::DynamicalSystem, depth::Int)
    for k = 1:depth
        B = subdivide(B)
        B = intersect(B,f(B))   # or f.(B) ?
    end
end

1. Warum `BoxSet`?
Der eigentliche Algorithmus startet mit einer einzelnen Box und wird anschließend mithilfe eines Suchbaums realisiert.
Ist hier gemeint, dass dieser Algorithmus für jede Box in `BoxSet` aufgerufen werden soll?

In Julia ist es so, dass Methoden standardmäßig *nicht* vektorisiert implementiert werden, sondern durch Broadcasting automatisch vektorisiert werden.

**Vorschlag:**

Die Methode sollte nur eine einzelne Box nehmen.

2. Was ist `DynamicalSystem`?

In Julia gilt es zu unterscheiden zwischen einer `struct` die `Callable` ist: `(::Bla)(x) = 2*x` und Methoden die von Haus aus `Callable` sind. Was ist hier gemeint?

3. Warum wird `subdivide` nicht *inplace* aufgerufen: `subdivide!(B)`?

4. Wie wird die Diskretisierung einer Box definiert?

So wie der aktuelle Vorschlag aussieht, müsste das entweder in `B` oder `f` definiert sein.
Das sind beides keine guten Optionen, da das anhand der Namen `BoxSet` oder `DynamicalSystem` nicht klar werden würde, dass da noch eine Diskretisierungsstrategie drinsteckt.

**Vorschlag 1:**

Wir bezeichnen das `f`-Argument nicht mehr als `DynamicalSystem`, sondern verändern die Voraussetzungen an `f`:

We assume that `f` takes a box and returns an iterator over points.
Example: `g` maps points to points and `points(box)` is an iterator over points.
Then `f = box -> (g(point) for point in points(box))` statisfies the assumptions for `f`.
That is, the strategy is included in `f`.

**Vorschlag 2:**

`Box` durchgehend so definieren, dass klar wird, dass da gleichzeitig auch eine Diskretisierungsstrategie mit drinsteckt. Großer Nachteil: Wenn wir fürs plotten Millionen von Boxen collecten, dann muss die Box `struct` so schlank wie möglich sein. Zusätzlich die Strategie in jeder Box zu speichern wäre ein großer Overhead.

Würde man dagegen die Strategie nur einmal pro `BoxSet` speichern gibt es ein anderes Problem: Dann kann der User nicht mehr einfach `BoxSet([box1, box2, ...])` schreiben, sondern muss bereits hier eine Strategie mitgeben.

5. Was ist `f(B)`?

Vorweg: Falls `f` eine Methode ist, die nur einen Punkt nimmt, wird dieser Aufruf gar nicht funktionieren.
Wir müssen hier die `f.(B)`-Syntax nehmen und dann `broadcast(::Any, ::BoxSet)` definieren.
Jetzt ist aber immer noch nicht klar, was `f.(B)` zurückgeben soll.

Allgemeineres Problem: `f(B)` hängt von mehr ab als `f` und `B`.
Das direkte Objekt, die Bilder aller Boxen `B` unter der Abbildung `f`, wird per se nicht direkt von dem Algorithmus benötigt. Darüber hinaus ist dieses Objekt abhängig von der Partitionierung, was aus der Notation nicht direkt klar wird. (z.B. ein `BoxSet`, das aus verschiedenen Verfeinerungen/Tiefen besteht)
Wenn man `f(B)` explizit mit liefern möchte, muss man die ganze Partition mit geben/berechnen, das ist ein höherer, unnötiger Arbeitsaufwand.
Wollen wir diesen künstlichen Aufsatz für die Userfreundlichkeit in Kauf nehmen?

**Vorschlag 1:**

Sei $S$ eine Menge von Boxen und `f` bilde Punkte auf Punkte ab.
Dann ist `f.(S)` die Teilmenge der Boxen von $S$, die `f` zumindest einmal trifft.
Der `intersect` Aufruf ist dadurch überflüssig, da dann `f.(S) == intersect(B,f(B))`.

Problem mit diesem Vorschlag: Das würde für jeden Aufruf von `f.` einen neuen Baum allokieren.
Es wäre besser, wenn mit nur einen Baum arbeiten würden, den wir verändern.
Da aber `f.(S)` in Julia niemals einen `view` zurückgibt, lässt sich das mit dieser Syntax nicht umsetzen.

**Vorschlag 2:**

`f.(tree)` erzeugt einen Iterator über `f(box)` wobei `box` ein Blatt des Baumes ist.

    # this does not seem very idiomatic...
    # Base.Broadcast(f, tree::Tree) = (f(box) for box in leaves(tree))
    Base.Broadcast(f, tree::Tree) = (f(box) for box in collect(leaves(tree)))

    function relative_attractor(f, domain::Box, depth)
        tree = Tree(domain)

        for k = 1:depth
            subdivide!(tree)
            tree = tree[f.(tree)]
        end

        return leaves(tree)
    end

Anmerkung zur Implementation von `f.`:
Das kann nicht komplett *Lazy* sein (wie im auskommentieren Code zu sehen wäre), da der Nutzer später den Baum ja wieder verändern kann. Daher muss diese Funktion alle Blätter im Baum bestimmen und extra abspeichern. Wenn man statt `tree1[f.(tree2)]` zum Beispiel sowas hätte wie `select(f, tree1, tree2)`, dann würde diese Allokation wegfallen.

Anmerkung zu `tree[f.(tree)]`
Die tatsächliche *Arbeit* (also die Suche) steckt hier im Aufruf `tree[(...)]`.

## unstable_set!

function unstable_set!(B::BoxSet, f::DynamicalSystem)
    C = B
    fB = f(B)
    while fB ≠ ∅             # or !isempty(fB)
        B = setdiff(fB,C)
        C = union(C,fB)
        fB = f(B)
    end
end

1. Selbes Problem: Was ist `f(B)`?

2. Wie kann der Algorithmus mit einem allgemeinen `BoxSet` funktionieren?
Die Realisierung scheint schwierig, da wir eigentlich einen Suchbaum verwenden wollen. (siehe Punkt 3.)

3. In diesem Algorithmus starten wir mit einer großen `Domain` (also dem Bereich innerhalb dem wir arbeiten) und fügen dann in dieser Domain eine kleine Box ein. Die Syntax in Matlab-GAIO dafür ist `tree.insert(x, depth)`. Wir schlagen dagegen vor, das ein bisschen expliziter zu machen: Wir nehmen ein Blatt von einem perfekten Binärbaum einer gewissen Tiefe und fügen dieses Blatt in einen leeren Baum ein. (`C = perfect_tree[(startpoint,)]`)
Gleichzeitig ist dieser perfekte Binärbaum in anderer Hinsicht nützlich: Er macht klar, bezüglich welcher Partition die Menge der Boxen bestimmt wird.
In Matlab-GAIO gibt es den eigentlich auch (siehe die Befehle mit einem `depth`-Argument); wir machen das aber expliziter.
Vorteil: Man muss die Tiefe nur einmal am Anfang festlegen. Bei Matlab muss man aufpassen, das die Tiefe immer richtig übergeben wird.

4. Wie funktioniert `union` und `setdiff`?

**Vorschlag:**

    struct Tree
        # good old tree
    end

    struct PerfectTree
        # efficiently represents the perfect binary tree with given depth
        # nodes are not saved; just the domain and the depth
    end

    function unstable_set(f, domain::Box, startpoint, depth) 
        perfect_tree = PerfectTree(domain, depth)

        C = perfect_tree[(startpoint,)]
        B = perfect_tree[f.(C)]

        while !isempty(B)
            setdiff!(B, C)
            union!(C, B)
            B = perfect_tree[f.(B)]
        end

        return leaves(C)
    end 

Note: We have to wrap `startpoint` into an iterator over points, since `perfect_tree[startpoint]` would only return a Box, not a Tree.

## transition_matrix!

function transition_matrix!(B::BoxSet, f::DynamicalSystem)
    # topological transition matrix between the boxes in B
    L = Lipschitz(f)   # estimates the Lipschitz constant of f
    n = size(B)
    E = zeros(n*L, 3)
    
    for b in B
        fb = f(b), m = size(fb)
        append!(E, [ones(m)*no(b), no.(fb), ones(m)]) # no(b) is in 1:n
    end
    
    I, J, V = E[1], E[2], E[3]
    T = sparse(I, J, V, n, n, combine)
end

1. Hier macht es tatsächlich Sinn, eine Menge von Boxen als Argument zu übergeben. Wenn wir aber eine effiziente Suche wollen, sollten diese Boxen bereits als Baum repräsentiert sein. Im Weiteren betrachten wir die Version des Algorithmus, die von einem Baum als Argument ausgeht.

2. Die Diskretisierungsstrategie ist sehr relevant und wird vom Algorithmus direkt benötigt. Zum Beispiel muss der Algorithmus für die Bestimmung der Übergangswahrscheinlichkeiten ein gewisses Maß für die Häufigkeit der Treffer haben, wie es zum Beispiel durch Zählen von Testpunkten gegeben ist.
Deswegen: Wir fordern, dass `f(box)` ein Iterator ist mit
`Base.IteratorSize(...) = Base.HasLength()`
Jetzt können wir `length(f(box))` aufrufen.

function transition_matrix!(tree, f)
    n = length(leaves(tree))
    E = Vector{Float64}[] # make append! work... might still be the wrong choice
    
    for box in leaves(tree)
        fb = f(box)
        hits = count_hits(fb, tree)
        m = length(hits)
        no_box = length(fb) # number of points in box
        append!(E, [ones(m)*no_box, hits, ones(m)])
    end
    
    I, J, V = E[1], E[2], E[3]
    T = sparse(I, J, V, n, n, combine)
end

This last section is a work in progress.