# Taking advantage of Apple M1?

**URL:** <https://discourse.julialang.org/t/taking-advantage-of-apple-m1/93310>\
**Category:** Performance\
**Tags:** mac-m1, hardware\
**Created:** [January 21, 2023, 1:17pm UTC](https://discourse.julialang.org/t/taking-advantage-of-apple-m1/93310 "2023-01-21T13:17:25Z")\
**Posts on this page:** 1\
**Showing post:** 9

<div class="post-metadata">

**Author:** ![LaurentPlagne](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/laurentplagne/32/10103_2.png) [@LaurentPlagne](https://discourse.julialang.org/u/LaurentPlagne)\
**Post date:** [January 22, 2023, 4:05pm UTC](https://discourse.julialang.org/t/taking-advantage-of-apple-m1/93310/9 "2023-01-22T16:05:03Z")

</div>

Here is a memory bound example (wave propag with low order F.D. scheme) with an Makie (GLMakie) animation running during the computations.

> **SingleSpring.jl**
>
> ```julia
> using GLMakie
> using PlotUtils #for colorant and color gradient (cgrad)
> 
> Base.@kwdef struct SpringParam
> ls::Float64 # Spring length
> ms::Float64 # Mass
> ks::Float64 # Spring Strength
> ns::Int # Mass number
> end
> 
> Base.@kwdef struct InitParam
> λ::Float64 # Wavelength of the initial deformation
> shift::Float64 # Maximal initial node displacement
> pos::Float64 # Central position of the node displacement
> end
> 
> Base.@kwdef struct AnimParam
> δt::Float64 # Time step for the dynamic simulation
> nδt::Int # Number of timesteps
> nδtperframe::Int # Number of timesteps per animation frame
> end
> 
> steady_position(i, ls) = (i - 1) * ls
> function initial_position(i,j,ls,λ,shift,posx,posy)
> xs=steady_position(i,ls)
> ys=steady_position(j,ls)
> dx=xs-posx 
> dy=ys-posy
> dr=sqrt(dx^2+dy^2)
> ds=(-1/λ^2)*exp(-(dr/λ)^2/2)*shift
> xs+dx*ds,ys+dy*ds
> end 
> 
> function display_perf(td,xc,nδt)
> nflops=20 # Manual count floating point ops
> nRW=12 # Manual count Array Reads and Writes
> T=eltype(xc) # i.e. Float64
> floatsize_inbytes=sizeof(T) #Float64 -> 8 Bytes
> nbytes=nRW*floatsize_inbytes # Float64 -> 96 Bytes
>     
> mupns=length(xc)*nδt/(td*1.e9) # mass pos updates per ns
> println("$(round((mupns*nflops),sigdigits=3))\t GFLOPS")
> println("$(round((mupns*nbytes),sigdigits=3))\t GB/s")
> end
> 
> #V is the type of the returned arrays (type are first class citizens !)
> function initialize_arrays(sp,ip,V) 
> (;ls,ns) = sp
> (;λ,shift,pos) = ip
> 
> #2D Array comprehensions 
> # -> cxs,xys,cxc,cyc are 2D arrays of size (ns,ns)
> cxs=[steady_position(i,ls) for i ∈ 1:ns, _ ∈ 1:ns] 
> cys=[steady_position(j,ls) for _ ∈ 1:ns, j ∈ 1:ns]
> cxc=[initial_position(i,j,ls,λ,shift,pos,pos)[1] for i ∈ 1:ns, j ∈ 1:ns]
> cyc=[initial_position(i,j,ls,λ,shift,pos,pos)[2] for i ∈ 1:ns, j ∈ 1:ns]
> 
> #Broadcast of V Ctor equivalent to xs=V(cxs);ys=V(cys);xc=V(cxc);yc=V(cyc)
> xs,ys,xc,yc=V.((cxs,cys,cxc,cyc))
> #Assume zero initial velocities : xp(revious)=xc(urrent)
> fx,xt,xp=zero(xc),zero(xc),copy(xc)
> fy,yt,yp=zero(yc),zero(yc),copy(yc)
> 
> xs,xc,xp,xt,fx,ys,yc,yp,yt,fy
> end
> 
> function update_position!(xt::Array{T,2},xc::Array{T,2},xp::Array{T,2},
> fx::Array{T,2},δt, ms) where {T}
> 
> nsx,nsy=size(xc)
> coef = T(δt^2 / ms)
> t2 = T(2)
> Threads.@threads for j ∈ 1:nsy
> @inbounds for i ∈ 1:nsx
> xt[i,j] = t2*xc[i,j]-xp[i,j]+fx[i,j]*coef
> end
> end
> end
> 
> function update_force!(fx::Array{T,2},xc::Array{T,2},ks) where {T}
> nsx,nsy=size(xc)
> mtks=T(-ks)
> t4=T(4)
> Threads.@threads for j ∈ 2:nsy-1
> @inbounds for i ∈ 2:nsx-1
> fx[i,j] = mtks*(t4*xc[i,j]
> -xc[i,j-1]
> -xc[i-1,j]
> -xc[i+1,j]
> -xc[i,j+1])
> end
> end
> end
> 
> function density_shift!(dc,xc,yc,ls)
> ns=size(xc,1) #assuming square array
> r0,r1,r2=1:ns-2,2:ns-1,3:ns # UnitRange
> ds=1/(4ls^2) #steady density
> @views @. dc[r1,r1] = 1/((xc[r2,r1] - xc[r0,r1])*(yc[r1,r2] - yc[r1,r0])) - ds
> end
> 
> function setup_scene(xaxis,dc,nf)
> # scene = Scene(resolution = (1500, 1500))
> dcnode = Observable(dc)
> 
> xM=maximum(xaxis)
> 
> dM=maximum(dc)
> fig = Figure(resolution = (2000, 1400), fontsize = 22)
> 
> ax = Axis3(fig[1, 1], limits=(0,xM,0,xM,-dM/8,dM/8))
> 
> surface!(ax,xaxis,xaxis,lift(d->d,dcnode),colormap=corporate_gradient(),colorrange = (-0.1,0.1))
> # scale!(ax, 1, 1, 0.1*xM/dM)
> 
> # translate_cam!(scene,(10,10,2))
> # update_cam!(scene, Vec3f0(0, 0, 5), Vec3f0(0.01, 0.01, 0))
> # rotate_cam!(scene,(π/4 ,-π/2 ,0))
> # αstep=-(π/4)/(nf÷2)
> # lift(d->rotate_cam!(scene,(0 ,αstep ,0)),dcnode)
> 
> fig,dcnode
> end
> 
> #Triscale colors ;)
> corporate_gradient() = cgrad([colorant"#177272",colorant"#fafafa",colorant"#ee6f40"])
> 
> function advance_nδtpf(xc, xp, xt, fx, sp, ap)
> (;ms, ks) = sp
> (;δt, nδtperframe) = ap
> 
> for _ ∈ 1:nδtperframe # Anonymous iter variable (_)
> update_force!(fx, xc, ks)
> update_position!(xt, xc, xp, fx, δt, ms)
> xc, xp, xt = xt, xc, xp # Julia's swap
> end
> xc, xp, xt
> end
> 
> function animate_makie(sp,ip,ap,V)
> (;δt,nδt,nδtperframe) = ap
> xs,xc,xp,xt,fx,ys,yc,yp,yt,fy=initialize_arrays(sp,ip,V)
> 
> dc=zero(xc)
> cdc=Array(xc)
> density_shift!(dc,xc,yc,sp.ls)
> @. cdc = dc 
> xaxis=Array(xs[:,1])
> yaxis=Array(ys[1,:])
> nf=nδt÷nδtperframe
> 
> scene,dcnode=setup_scene(xaxis,cdc,nf)
> display(scene)
>     
> t=0.0
> tdynamic=0.0
> for _ ∈ 1:nf
> tdynamic+= @elapsed begin
> xc,xp,xt=advance_nδtpf(xc,xp,xt,fx,sp,ap)
> yc,yp,yt=advance_nδtpf(yc,yp,yt,fy,sp,ap)
> end
> density_shift!(dc,xc,yc,sp.ls)
> dcnode[] = dc
> sleep(0.01)
> t+=nδtperframe*δt
> end
> 
> display_perf(tdynamic,xc,nδt) #Display GFLOPS and GB/s
> end
> 
> function go()
> sp=SpringParam(ls=0.1,ms=1,ks=2.5,ns=1000)
> ip=InitParam(λ=3.0,shift=0.1,pos=sp.ls*sp.ns/2)
> ap=AnimParam(δt=0.2,nδt=5000,nδtperframe=20)
> 
> V=Array{Float64,2}
> 
> @time animate_makie(sp,ip,ap,V)
> end
> go()
> 
> ```

I obtain (M1 Max):

```julia
julia> include("SingleSpring.jl")
21.1 GFLOPS
101.0 GB/s
  8.506742 seconds (1.46 M allocations: 1.164 GiB, 0.86% gc time, 1.18% compilation time)

```

with a final image :

 ![Capture d’écran 2023-01-22 à 16.59.27](https://global.discourse-cdn.com/julialang/original/3X/9/3/93a20b4627a66543159f5abeaa6dec61d20300b8.jpeg)

---

_[View the full topic](https://discourse.julialang.org/t/taking-advantage-of-apple-m1/93310)._
