# This will give the wrong result:
2^100
# An exact way to compute the answer is using arbitrary-precision integers
BigInt(2)^100
# Usually faster than arbitrary-precision computations are
# floating point operations (might not be true for this simple case, however)
2.0^100
Compute $$ 15! \qquad 100! \qquad \left(\begin{array}{c} 100 \\ 15 \end{array}\right) $$ with the Julia you know so far.
res15 = 1
for i in 1:15
res15 = res15 * i
end
res100 = BigInt(1)
for i in 1:100
res100 = res100 * i
end
res85 = BigInt(1)
for i in 1:(100 - 15)
res85 = res85 * i
end
res100 / res85 / res15
function acceleration(x)
if x > 1 || x < -1
-4 * sign(x) * (abs(x) - 1)^3
else
0
end
end
x = 0.
v = 1.
for i in 1:5
x, v = euler(acceleration, 0.1, x, v)
end
println(x, " ", v)
#
# A shorter, but equivalent way using the "do" Syntax
#
x = 0.
v = 1.
for i in 1:5
x, v = euler(0.1, x, v) do x
if x > 1 || x < -1
-4 * sign(x) * (abs(x) - 1)^3
else
0
end
end
end
println(x, " ", v)
Which of the following type are subtypes of another?
Try to guess first and then verify by using the operator <:
.
Float64 AbstractFloat Integer
Number AbstractArray Complex
Real Any Nothing
The following type chains (subtype <:
supertype) are true:
Float64 <: AbstractFloat <: Real <: Number <: Any
Integer <: Real <: Number <: Any
Complex <: Number <: Any
Nothing <: Any
AbstractArray <: Any
ones(Int, 3, 5) + ones(Int, 3, 5)
# or
2ones(Int, 3, 5)
reshape(collect(1:16), 4, 4) / 10
We want improve our capabilities to perform simulations to two dimensions using the 2D harmonic potential $$ V\left( \begin{array}{c} x_1 \\ x_2 \end{array}\right) = x_1^2 + x_2^2 $$ and the acceleration map $$ \vec{A}_V = - \nabla V. $$ The adapted forward-Euler scheme is $$\left\{\begin{array}{l} \vec{v}^{(n+1)} = \vec{v}^{(n)} + \vec{A}_V(\vec{x}^{(n)}) \Delta t\\ \vec{x}^{(n+1)} = \vec{x}^{(n)} + \vec{v}^{(n)} \Delta t\\ \end{array}\right. .$$
Change the euler
function we introduced in 02_Functions_Types_Dispatch.ipynb
accordingly (Hint: Suprisingly few changes are needed) and run the dynamics for five steps using
$$x_n = \left( \begin{array}{c} 0 \\ 0 \end{array}\right) \qquad v_n = \left( \begin{array}{c} 1 \\ 0 \end{array}\right) .$$
Check against your previous implementation.
# The previous euler implementation can be used without any additional change:
euler(A, Δt, xn, vn) = xn + vn * Δt, vn + A(xn) * Δt
# Assuming x to be the vector (x_1 \\ x_2), the derivative of V can be written as:
A2d(x) = -2x
# Therefore propagating 5 steps is equal to:
x = [0, 0]
v = [1, 0]
for i in 1:5
x, v = euler(A2d, 0.1, x, v)
end
println(x)
println(v)
Exercise:
A much more stable integrator than the euler
we used so far is the verlocity Verlet:
# One way to code the verlet function is:
function verlet(A, Δt, xn, vn)
An = A(xn)
x_next = xn + vn * Δt + An/2 * Δt^2
v_next = vn + (An + A(x_next)) / 2 * Δt
x_next, v_next
end
Program the total potential function for a matrix $\textbf{x}$. A useful function is norm
from the LinearAlgebra
package.
# One solution:
function Vtot(Vpair, x)
n_particles = size(x, 2)
accu = zero(eltype(x)) # Get a zero of the appropriate type
for i in 1:n_particles, j in i+1:n_particles
accu += Vpair(norm(x[:, i] .- x[:, j]))
end
accu
end