juliateachingctu / scientific-programming-in-julia Goto Github PK

Abstraction was too tough, explain better, how things click it.

• The transition from theory to practice was difficult.

https://discourse.julialang.org/t/how-are-invalidated-methodinstances-found-and-how-well-does-the-process-scale/83099/4

• Simplify/split last exercise in multiple parts
• explain better why we are moving to parametric types
• should we maybe move away from the ecosystem? simpler example for parametric types would be good...?
• add full solution of lab 03 (sketched below)

EcosystemCore.mates(::Animal{S,Female}, ::Animal{S,Male}) where S<:Species = true
EcosystemCore.mates(::Animal{S,Male}, ::Animal{S,Female}) where S<:Species = true
EcosystemCore.mates(::Agent, ::Agent) = false

##########  Additional functionality  ##########################################

function simulate!(world::World, iters::Int; cb=()->())
    for i in 1:iters
        world_step!(world)
        cb()
    end
end

agent_count(p::Plant) = size(p)/EcosystemCore.max_size(p)
agent_count(::Animal) = 1
agent_count(as::Vector{<:Agent}) = sum(agent_count,as)

function agent_count(w::World)
    function op(d::Dict,a::Agent{S}) where S<:Species
        n = nameof(S)
        if n in keys(d)
            d[n] += agent_count(a)
        else
            d[n] = agent_count(a)
        end
        return d
    end
    foldl(op, w.agents |> values |> collect, init=Dict{Symbol,Real}())
end

function every_nth(f::Function, n::Int)
    i = 1
    function callback(args...)
        # display(i) # comment this out to see out the counter increases
        if i == n
            f(args...)
            i = 1
        else
            i += 1
        end
    end
end

Hi, I have a I think correct solution for the string replacement, in case you wish to add it.

sreplace_i(s) = replace(s, r"([^\w]|\b)i(?=[^\w]|\z)" => s"\1k")

It replaces "i" only if it is not preceded or followed by an alphanumeric character or there is an end/beginning of string.

example:

s = "i + i*i + y*i - sin(z) + int(i) + i"
sreplace_i(s)

"k + k*k + y*k - sin(z) + int(k) + k"

Using filter is both slow and not super straight forward in this case. We could remove hint 1 and expand on hint 2.

module A
f() = println("hello")
macro foo()
      :(f())
end
end
Main.A

julia> f() = println("goodbye")
f (generic function with 1 method)

julia> A.@foo
hello

and compare this to

module A
f() = println("hello")
macro foo()
      :(esc(f())))
end
end

https://discourse.julialang.org/t/rust-julia-comparison-post/75403/52

In 2nd lab we are provided with method agent_step!(a::Animal, w::World).

https://github.com/JuliaTeachingCTU/Scientific-Programming-in-Julia/blob/master/docs/src/lecture_02/lab.md#code-skeleton

This method assumes we have implmented kill_agent!() method which is then used here

if energy(a) <= 0 kill_agent!(a,w) return end

But this is never referenced as a task or side note, but it is present in solved solutions.

It would be worth to add task for creating this method aswell, so the simulation works in the end after adding all other methods.

We should probably use PkgTemplates.jl to simplify the lab. Orient on Vasek's workflow tutorial in lecture 7 of the julia for optimization course.

All the versions of polynomial function should share the same order of coefficients, that is consistent with how it is evaluated in Julia's evalpoly.

Currently in the lab we have coefficients ordered from a_0 to a_n, but it is probably reversed such that the order from the original polynomial a_n*x^n + a_{n-1}*x^{n-1} + ... + a_1*x + a_0 is conserved.

Affected pages

• Lab 1 (naive)
• Lab 5 (Horner)
• Lab 6 (recursive)
• Lab 7 (macro generation)

https://simeonschaub.github.io/ReverseModePluto/notebook.html

Use parametric type Sex

https://www.youtube.com/c/juliafortalentedamateurs/videos

Show them @forward macro, how to write it and what it does.

https://discourse.julialang.org/t/edited-video-of-course-gpu-programming-with-julia-at-the-swiss-national-supercomputing-centre-cscs-now-available/80290

https://github.com/JuliaDiff/ReverseDiff.jl

ODE variables/coordinates have different names in code and in latex equations...

• First, when working out the tape example, we should use the same names of variables, i.e. x,y h1, h2, h3, z...
• we should explain that the computation graph is not stored in a single place, but it is stored in a distributed fashion in nodes.
• we might break the computation to steps, for example h1 = x * y and show, how the graph is built and how edges are stored and where.

Homework: custom adjoint of tanh/map

Ideally we should have used deepcopy and fixed Random.seed! for the simulation profiling.
Different versions of the code might use RNG in a different way, thus the simulation is also different.
The same applies for the benchmarking. We should benchmark only what is needed: find_rand.

In homework from week 6, there semms to be a problem with the assignment itself https://github.com/JuliaTeachingCTU/Scientific-Programming-in-Julia/blob/master/docs/src/lecture_06/hw.md.
There is line:

return an array of unique alphabetically sorted symbols representing variables in an expression.

However next to it is example
[:x, :y, :z, :c]
Which is not alphabetically sorted

https://discourse.julialang.org/t/reason-behind-designing-parametric-types-as-invariant/18971

• Add more references to included images.

• Small note on NVidia/AMD framework comparison

  • SIMT/SIMD comparison (AMD calls multiple CUDA cores a SIMD unit)
  • warp is 32 threads whereas the equivalent wavefront in AMD has 64
  • on AMDGPU the requirement for warp alignment can be broken, however we cannot switch context that easily as a result

• there is a nice picture which shows the division of work over a 1D/2D array into blocks/threads etc - https://youtu.be/LG9G4aA28rU?t=368 - some analog should be available somewhere or created manually

  • x property suggest that you can partition the execution along three-dimensional cube (three nested for loops), - should be explained better add pictures
    • such as https://developer-blogs.nvidia.com/wp-content/uploads/2017/01/cuda_indexing-1024x463.png in https://developer.nvidia.com/blog/even-easier-introduction-cuda/)
    • or http://2.bp.blogspot.com/-gSMyMA7nnIU/UPbTSAgV_CI/AAAAAAAAAtE/2D942iHCg9Q/s1600/gpu2.png

• using GPU without writing kernels should mention the scalar indexing problem, which we often encounter (comparing CPU and GPU array without moving them to the same place, show methods fallback to getindex! ops)

Test and if it succeeds, put AbbreviatedStackTraces.jl to Lecture 1 on parsing the errors

the code here

https://juliateachingctu.github.io/Scientific-Programming-in-Julia/stable/lecture_02/lab/#Finding-food-for-sheep

function find_food(a::Sheep, w::World)
    as = filter(x->isa(x,Grass), w.agents |> values |> collect)
    isempty(as) ? nothing : sample(as)
end

filters objects so that only grass is present, but does not check the size of grass, so 0 sheep may try to eat 0 size grass and get 0 energy

Unit test can be setup also on GitLab, see tutorial from Tamás K. Papp.
It boils down to adding a file .gitlab-ci.yml to the root of a repo with the following few lines

image: julia:1.6                # image comes from Docker hub

default:
  script:
    - julia --project=@. -e "import Pkg; Pkg.test()"

Documentation can be hosted in a similar fashion as well.

https://arxiv.org/pdf/2109.01950.pdf

https://discourse.julialang.org/t/when-does-the-julia-compiler-step-in/81930

https://github.com/dfdx/Yota.jl/blob/main/src/grad.jl

https://discourse.julialang.org/t/c-code-much-faster-than-julia-how-can-i-optimize-it/87868/22

Lab 02

• implement rational numbers (do we have a real world use case?)
• implement ecosystem without world
  • do it bottom up with simple examples instead of top down
  • implement: types (non-parametric Wolf/Sheep/Grass, show, eats

Improvement suggestions from last year

• get rid of getters/setters (#16)

Lab 03

• quick, theoretical walk through of ⚥Sheep forward method including overloading of getproperty (and setproperty!) for ⚥Sheep:

function Base.getproperty(s::⚥Sheep, x::Symbol)
    if x == :sex || x == :sheep
        getfield(s,x)
    else
        getfield(s.sheep,x)
    end
end

• introduce parametric type for species/sex
• implement rest of ecosystem (immutable version?!!)

Improvement suggestions from last year

• Simplify/split last exercise in multiple parts
• explain better why we are moving to parametric types
• Previous #13: Using filter is both slow and not super straight forward in this case. We could remove hint 1 and expand on hint 2.

  Scientific-Programming-in-Julia/docs/src/lecture_02/lab.md

  Lines 372 to 377 in c77e2a5

  	1. Hint: For the functional programming way of coding this can use `filter` and
  	`isa` to filter for a certain type and `StatsBase.sample` to choose a random
  	element from a vector. You can get an `Iterator` of values of your dictionary
  	via `values` which you might want to `collect` before passing it to `sample`.
  	2. Hint: You could also program this with a for-loop that iterates over your
  	agents in a random order.

• add full solution of lab 03 (sketched below)

EcosystemCore.mates(::Animal{S,Female}, ::Animal{S,Male}) where S<:Species = true
EcosystemCore.mates(::Animal{S,Male}, ::Animal{S,Female}) where S<:Species = true
EcosystemCore.mates(::Agent, ::Agent) = false

##########  Additional functionality  ##########################################

function simulate!(world::World, iters::Int; cb=()->())
    for i in 1:iters
        world_step!(world)
        cb()
    end
end

agent_count(p::Plant) = size(p)/EcosystemCore.max_size(p)
agent_count(::Animal) = 1
agent_count(as::Vector{<:Agent}) = sum(agent_count,as)

function agent_count(w::World)
    function op(d::Dict,a::Agent{S}) where S<:Species
        n = nameof(S)
        if n in keys(d)
            d[n] += agent_count(a)
        else
            d[n] = agent_count(a)
        end
        return d
    end
    foldl(op, w.agents |> values |> collect, init=Dict{Symbol,Real}())
end

function every_nth(f::Function, n::Int)
    i = 1
    function callback(args...)
        # display(i) # comment this out to see out the counter increases
        if i == n
            f(args...)
            i = 1
        else
            i += 1
        end
    end
end

Checking type stability on function with optional arguments does not show the whole picture

@code_warntype find_root(p, Newton(), -5.0, 5.0) # shows only the wrapper function which is stable
@code_warntype find_root(p, Newton(), -5.0, 5.0, 100, 0.95, 1e-12)  # fully specified (shows type instabilities - also used in the eval)

The problem here is that the function is type stable in the way that the compiler knows the output type before execution, however the eval pipeline checks for all ::Any entries in the fully specified function, i.e. for runtime dispatches of which there are few. This leads to confusion based on what has been said in the lab.

We should have a naive explanation for a macro when we first use it (probably already in the second lecture).

https://juliadiff.org/ChainRulesOverloadGeneration.jl/dev/examples/reverse_mode.html

type anotated global

https://discourse.julialang.org/t/rust-julia-comparison-post/75403

https://www.youtube.com/watch?v=bXHe7Kj3Xxg

• https://www.youtube.com/watch?v=an6_N6RM8Zo

https://discourse.julialang.org/t/ann-interactivecodesearch-jl-interactively-search-julia-code/17657/12

Implement / test / benchmark a parallelized version of the Ecosystem which would be based on multiple phases of an agent step:

1. create a list of all potential interactions (e.g.: "wolf eats sheep", "sheep eats grass", ...)
2. find collisions like "sheep1 eats grass1" / "sheep2 eats grass1" and resolve them
3. execute interactions
4. remove / spawn agents

• Mention Peter Norvig's claim that in functional programming most design patterns are not needed?
• Make a difference between closure and partial evaluation
• Do we mention Inner Constructor to enforce type invariants?
• Note that the structure ##1 is not directly accessible. Try f.x and g.x. is I think incorrect as you can access it. Try in repl

function call_counter(f)
	called = 0
	(args...;kwargs...) -> begin
		called += 1
		f(args...;kwargs...)
	end
end

counted_println = call_counter(println);
counted_println("hello world"

But somehow the above does not trigger that
https://discourse.julialang.org/t/accessing-data-in-a-closure/54930/4

• Should we mention the issue of Discoveribility of Interfaces? A recently recommended way (adopted from ?Haskell?) is to isolate interface to a package defining them, example is ArrayInterfaces
• Regarding making comparison to OOP, there is now this package https://discourse.julialang.org/t/ann-objectoriented-jl-renamed-from-tyoop-jl-complete-oop-support/83911

The getters and setters we currently use in the Ecosystem are only there because of the ⚥Sheep type which should demonstrate the forwarding of fields. We could get rid of all of them by overloading getproperty (and setproperty!) for ⚥Sheep instead:

function Base.getproperty(s::⚥Sheep, x::Symbol)
    if x == :sex || x == :sheep
        getfield(s,x)
    else
        getfield(s.sheep,x)
    end
end

In labs, we talk a lot about polynomials as an example where meta programming buys you speed due to loop unrolling.
But students might not appreciate why polynomials are such a big deal.

• they are important for computing complicated functions, which are mostly computed from taylor expansion.
• We have found out polynomials to be super important for private machine learning, as secret sharing schemes rely on that https://mortendahl.github.io/

https://juliagraphics.github.io/ColorSchemes.jl/stable/images/

• Precompilation has improved by a great mean
• As mentioned in notes to 3, Interface -> Separate(sub)Module for discoverability
• Might be good place to talk about PackageCompiler, which is just an awesome way to avoid TTFT.