Current location - Training Enrollment Network - Mathematics courses - Solve the uniform distribution of points on the sphere in matlab.
Solve the uniform distribution of points on the sphere in matlab.
I went back and thought about it today. There is another way to generate evenly distributed points on a sphere.

I tried it, and it was faster than the previous genetic algorithm. It took two minutes to calculate 500 points.

It is also possible that it takes time to calculate convex hull between algorithms.

The new method is based on a physical model, which makes n points on the sphere repel each other.

The repulsive force decreases with the increase of the distance between points, for example, it can be a force that decreases with the square of the distance.

This problem can be compared to the motion of n dots with equal charges on a sphere.

When the movement finally reaches balance and the balls stop moving, the distribution of the balls will be uniform due to the repulsive force.

At this time, the problem becomes a numerical problem of multivariable differential equations, and the final equilibrium state is what we want.

With the simplest Euler difference formula, six state variables of X, Y, Z, vx, vy and vz are set for each point to describe the motion state.

At each step, the motion state of a point is calculated, and the initial state is randomly assigned, with a speed of 0.

Function? []=main()

n = 12; ? Percentage of points

a=rand(N, 1)* 2 * pi; % is evenly distributed according to the random surface, and Mr enters the initial state.

b=asin(rand(N, 1)* 2- 1);

r0=[cos(a)。 *cos(b),sin(a)。 *cos(b),sin(b)];

V0 = zero (size (r0)););

g = 1e-2; % repulsion constant, this value is better.

For what? Ii= 1:200% simulates 200 steps, which generally converges, but in fact you can quit under it.

[rn,vn]=countnext(r0,v0,G); % update status

r0 = rnv0 = vn

end

' plot3(rn(:, 1),rn(:,2),rn(:,3),'.'); Keep? Open; % drawing results

[xx, yy, ZZ]= sphere (50);

h2=surf(xx,yy,ZZ); ? % Draw a unit ball for reference.

set(h2,' edgecolor ',' none ',' facecolor ',' r ',' facealpha ',0.7);

Axis? Equality;

Axis ([- 1? 1? - 1? 1? - 1? 1]);

Keep? Close;

end

Function? 【rn? vn]=countnext(r,v,G)? % Function to update status

%r stores the x, y and z data of each point, and v stores the speed data of each point.

num=size(r, 1);

Dd = zero (3, num, num); ? Percentage of vector difference between points

For what? M= 1: number-1

For what? N=m+ 1: quantity

dd(:,m,n)=(r(m,)-r(n,)';

dd(:,n,m)=-dd(:,m,n);

end

end

L=sqrt(sum(dd。 ^2, 1)); Percentage of distance between points

L(L & lt; 1e-2)= 1e-2; ? % distance too small limit

F = sum (dd. /repmat(L.^3,[3? 1? 1]),3)'; % calculated resultant force

Fr=r.*repmat(dot(F,r,2),[ 1? 3]); ? % Calculate the radial component of resultant force.

Fv = F-Fr; ? % tangential component

rn = r+v; % update coordinates

rn=rn。 /repmat(sqrt(sum(rn。 ^2,2)),[ 1? 3]);

VN = v+ G * Fv; % and new speed

end

Another advantage of this method is that you can see the results of each step.

The following is a process of finding the uniform distribution of 12 points by an algorithm.

The final icosahedron is very close to the regular icosahedron.