0

在代码中取自：https ://tutorials.sciml.ai/html/models/01-classical_physics.html ，如下所示：

# 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")

我不明白 .= 运算符在函数谐波振荡器中的用法。使用 = 给了我错误的答案。所以，我想知道与 有什么.=不同=？它不是矢量化ddu的，因为 RHS 都是标量。

4

1 回答 1

2

我不明白 .= 运算符在函数谐波振荡器中的用法。[...] 它不是向量化 ddu 因为 RHS 都是标量。

这是; u, du, 和ddu不是标量，它们是长度为 1 的向量。

你可以问 Julia.=语法是什么意思：

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
))))

这看起来有点复杂，但它本质上是一个广播分配，类似于

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

使用 = 给了我错误的答案。

是的，因为 DiffEq 库希望函数harmonicoscillator修改输入。如果您只使用=该函数，则创建该函数的局部新变量，而不是修改输入向量，并且从外部看不到。

于 2021-04-24T15:17:21.693 回答