.= operator in Julia

Question

In the code taken from: https://tutorials.sciml.ai/html/models/01-classical_physics.html as given below:

# Simple Harmonic Oscillator Problem
using OrdinaryDiffEq, Plots

# Parameters
ω = 1

# Initial Conditions
x₀ = [0.0]
dx₀ = [π/2]
tspan = (0.0, 2π)

ϕ = atan((dx₀[1]/ω)/x₀[1])
A = √(x₀[1]^2 + dx₀[1]^2)

# Define the problem
function harmonicoscillator(ddu,du,u,ω,t)
    ddu .= -ω^2 * u
end

# Pass to solvers
prob = SecondOrderODEProblem(harmonicoscillator, dx₀, x₀, tspan, ω)
sol = solve(prob, DPRKN6())

# Plot
plot(sol, vars=[2,1], linewidth=2, title ="Simple Harmonic Oscillator", xaxis = "Time", yaxis = "Elongation", label = ["x" "dx"])
plot!(t->A*cos(ω*t-ϕ), lw=3, ls=:dash, label="Analytical Solution x")
plot!(t->-A*ω*sin(ω*t-ϕ), lw=3, ls=:dash, label="Analytical Solution dx")

I don't understand the usage of .= operator in the function harmonicoscillator. Using = gives me the wrong answer. So, I am wondering how is .= different from =? It is not vectorizing ddu because RHS is all scalar.

Answer 1

I don't understand the usage of .= operator in the function harmonicoscillator. [...] It is not vectorizing ddu because RHS is all scalar.

It is; u, du, and ddu are not scalars, they are length-1 vectors.

You can ask Julia what the .= syntax means:

julia> Meta.@lower a .= b
:($(Expr(:thunk, CodeInfo(
    @ none within `top-level scope'
1 ─ %1 = Base.broadcasted(Base.identity, b)
│   %2 = Base.materialize!(a, %1)
└──      return %2
))))

which looks a bit involved, but it is essentially a broadcasted assignment, similar to

for i in eachindex(a, b)
    a[i] = b[i]
end

Using = gives me the wrong answer.

Yes, because the DiffEq library expects the function harmonicoscillator to modify the input. If you use just = you create a new variable local to that function rather than modifying the input vector, and that is not visible from the outside.