Calculate streamlines using ODE

13 views (last 30 days)
Fabio Freschi
Fabio Freschi on 8 Mar 2023
Commented: Fabio Freschi on 24 Aug 2023
Hi Matlab Gurus,
I have read many posts in this forum about my problem but I'm still not capable to solve it. I'm trying with a direct support
Problem
I want to calculate streamlines starting from a 3d vector field function, using the differential equation
the field is not known a-priori in the region, so I prefer not to use the built-in stream3 function and interpolate the data.
The vector field is solenoidal, so the streamlines are closed.
Here is a MWE, with the attached fieldFunction (it calculates the magnetic field due to a circular coil, but it's only an example)
clear variables, close all
% coil geometrical params
a = 0.1;
d = 0;
% starting point
Q0 = [.05 .05 0];
% differential equation
streamLineFunct = @(t,x)fieldFunction(a,d,x(:).').';
% ODE options
options = odeset('OutputFcn',@odeplot,'OutputSel',[1 2 3]);
sol = ode15s(streamLineFunct,[0 0.5],Q0,options);
Q = sol.y.';
% plot
figure, axis equal, hold on, view([-1 1 1])
% plot coil
plot(a*exp(1j*(0:pi/100:2*pi)),'LineWidth',2);
% plot streamline
plot3(Q(:,1),Q(:,2),Q(:,3));
Because the integration variable 's' is non-physical, it is difficult to define a proper tspan for the ode solver. It's possible to see, in fact, that the streamlines make about 2.5 complete circles. A good choice for this streamline of this example is around 0.21, but this guess can be done only after the computation.
Question:
Version 1: how to get the correct streamline without unnecessary computation, stopping the ode solver once the position returns to the starting point?
Version 2: is Event Location the right option to pass to the ode solver? If so, how can I define a suitable event that avoids the repetition of unnecessary calculation and stop the ode once the streamline returns at the original position?
Clarifications
  1. I know that playing with AbsTol and RelTol the streamline can be more accurate, avoiding the visible misalignment.
  2. I already tried to define an Event but due to the initial small integration steps, I can't find a good setup for this property
  3. The real-problem vector function is very accurate but can be time consuming, so the use of a large tspan is not practical.
Thank you in advance for your help

Answers (1)

Ishu
Ishu on 24 Aug 2023
Hi Fabio,
If you try to output the “Q” values then you will notice that the positions do not exactly return to their starting point, there is a slight change in the position. Therefore, there is no repetition in the results after a certain time span.
You can use Events” option to stop the solver when some specific event is detected or the conditions you provided are met.
As in your case your condition to stop the solver is once the position returns to the starting point but this condition is not valid. Using “Events” option will not serve the purpose in this case.
For more information on ode15s solver you can refer below documentation:
Hope it helps.
  1 Comment
Fabio Freschi
Fabio Freschi on 24 Aug 2023
Ishu, thanks for your reply. However it is not solving my question. Check “Clarifications” section

Sign in to comment.

Products


Release

R2022b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!