Conversation
|
@mrcinv rebased dn03 tudi na new-main. Sedaj ok? |
| Because no simple parametric form is available, we approximate this curve numerically by integrating the associated ODE and then use the trajectory to estimate its length. | ||
|
|
||
| # Integrating with Explicit Euler | ||
| We integrate the ODE using the explicit Euler method. Given a current point $u_k$ on the curve and step size $h$, the next point is computed as |
There was a problem hiding this comment.
@kurtaga Explicit Euler method is order 1 and it needs orders of magnitude more calculations than higher order methods to achieve sufficient accuracy. The method used should be at least of order 4.
| ``` | ||
|
|
||
| # Curve Length | ||
| The algorithm produces a trajectory that traces the implicit intersection curve. The total arc length is approximated by summing distances between consecutive points. |
There was a problem hiding this comment.
@kurtaga Given that you normalized the cross product, the parameter from the ODE will be the natural parameter on the curve. So there is no need to perform an additional sum.
| @@ -0,0 +1,31 @@ | |||
|
|
|||
| # ODE right-hand side | |||
| function rhs(u, ∇F1, ∇F2) | |||
|
|
||
|
|
||
| # explicit Euler step | ||
| function explicit_euler_step(u, h, ∇F1, ∇F2) |
There was a problem hiding this comment.
Create a generic method for euler step that takes one function and not the gradients. Generic methods are easier to test and design. There is less room for bugs.
| function explicit_euler_step(u, h, ∇F1, ∇F2) | |
| function explicit_euler_step(u, h, rhs) |
Add a docstring.
| end | ||
|
|
||
| # integrate | ||
| function integrate_explicit_euler(u0, h, ∇F1, ∇F2; max_steps=250000, min_steps=100) |
| for k in 1:max_steps | ||
| u = explicit_euler_step(u, h, ∇F1, ∇F2) | ||
| push!(traj, u) | ||
| if k > min_steps && norm3(u - u0) < tol |
There was a problem hiding this comment.
You are effectively calculating period of the solution. Checking for the period this way is limited by the step size. So the accuracy is limited to the step size. To get higher accuracy, you should formulate the problem as a nonlinear equation ($y(t_0) = y(t_0 + p)$) and use methods to solve nonlinear equations (Newton, bisection, etc.)
| for k in 1:max_steps | ||
| u = explicit_euler_step(u, h, ∇F1, ∇F2) | ||
| push!(traj, u) | ||
| if k > min_steps && norm3(u - u0) < tol |
There was a problem hiding this comment.
To make the method integrate_explicit_euler generic, you can pass the stop condition as an argument.
function integrate_explicit_euler(u0, h, ∇F1, ∇F2; stop_condition=(u, t)->false, max_steps=250000, min_steps=100)
...
if stop_condition(u, t)
break
end
...
endand then the stop condition can be passed as an argument
integrate_explicit_euler(u0, h, ∇F1, ∇F2; stop_condition=(u, t)->abs(t)>1 && norm3(u -u0) < tol)| end | ||
|
|
||
| @testset "Integration" begin | ||
| # For this toy case, tangent always points +z, so trajectory is a straight line |
There was a problem hiding this comment.
It would be easier to test generic method for integration, as you can pass arguments with known solutions.
2 participants
No description provided.