double gyre - compute absolute and relative dispersion

30 views (last 30 days)
Hello everyone! I wrote a routine to integrate numerical particle trajectories for a double gyre flow and now I am supposed to compute the absolute and relative dispersion based on the computed trajectories. Can anyobdy help me in building the for cycle? I believe the cycle has three indices ( i = length of the x-domain, j = length of the y-domain and, k = timesteps) but i have no idea how to build it. you can find attached the scripts that I am using. If you run" integrateDGtemp.m" you can quickly understand what it is about.
doublegyreVEC.m contains the double gyre vector field
rk4singlestep.m contains the Runge Kutta 4 scheme for solving the ODE
integrateDGtemp.m is the script used to plot the initial particle positions and then the double gyre field.
Also, attached here you can see the formulas that I am supposed to use for computing absolute and relative dispersion.
Thank you a lot if you can help.

Accepted Answer

Jiri Hajek
Jiri Hajek on 9 Nov 2022
Hi, seems like you are almost there... But you have to make up your mind on the meaning of your dispersion formulae. In statistics, dispersion is not a concrete function, but rather a quality of a data set, which may be expressed by any of several statistical parameters. (I guess that the i, j, k indices all describe the same set of M particles and consequently the expressions in brackes have the meaning of eulerian distances between pairs of particles.)
To evaluate your formulae, you'll need just a few basic MATLAB functions:
  • diff to compute the differences of positions (the bracketed terms)
  • .* dot product of same-size vectors
  • sum of a vector
Using these, you will be able to calculate single element of the A2 and R2 matricex. Then, to complete these matrices at a time instant, you would need two more for loops.
Jiri Hajek
Jiri Hajek on 10 Nov 2022
Yes, your script would benefit from some enhancement and clear data structure. But you see, the scope and detail of your code is entirely up to you; here we can just discuss issues that cause confusion. To your question:
It seems not necessary to store particle positions for every timestep, to calculate dispersion. You may need storage for plotting or reporting, but let us focus on the algorithm, where you need just the current particle coordinates, say
currentX1 % 1xM array
currentX2 % 1xM array
and initial coordinates
initialX1 % 1xM array
initialX2 % 1xM array
So, your domain has two dimensions, x1 and x2. Then, your dispersion is a 2x2 matrix. You need to calculate each element of the dispersion matrix separately. Now, within your timestep, you need to add a loop to go through the M particles.
A2 = zeros(2); % initialize the dispersion matrix
for m = 1:M
A2(1,1) = A2(1,1) + (currentX1(m)-initialX1(m))^2/M;
A2(2,2) = A2(2,2) + (currentX2(m)-initialX2(m))^2/M;
A2(1,2) = A2(1,2) + (currentX1(m)-initialX1(m))*(currentX2(m)-initialX2(m))/M;
A2(2,1) = A2(1,2); % symmetry
% Here, you have the instantaneous absolute dispersion in A2
Similarly for the relative dispersion.

Sign in to comment.

More Answers (0)


Find more on Programming in Help Center and File Exchange




Community Treasure Hunt

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

Start Hunting!