Analytical intersection volume between two spheres

Compute the overlap volume between 2 spheres defined in an array

Compute the overlap volume between 2 spheres defined in an array.
Computation is vectorized, and intersection volume are computed an
analytical way.

Input: spheres data presented in an array G of four columns.
G contains parameters of the n spheres
. G(1:n,1) - x-coordinate of the center of spheres,
. G(1:n,2) - y-coordinate of the center of spheres,
. G(1:n,3) - z-coordinate of the center of spheres,
. G(1:n,4) - radii of the spheres
Each row of the array contains the information for one sphere.

Input can also be provided in three different vectors. These
vectors can be row or column vectors. The 1st one corresponds to
x-coordinate of the center of spheres, the 2nd one to the
y-coordinate, the 3rd one to the z-coordinate and the 4th one to
the radii of the spheres.
An optional binary argument can be provided to display or not the
result.

Output: Square matrix M(n,n) containing intersection volumes between
spheres
M(i,j) contains the intersection volume between spheres i & j
By definition, M(i,i) corresponds to the volume of sphere i.

Examples:

x = [0,1,5,3,-5];
y = [0,4,3,7,0];
z = [0,4,3,7,0];
r = [1,5,3,2,2];
Display_solution = true;
disp('First call')
M = volume_intersect_sphere_analytical(x,y,z,r,Display_solution);
disp('Second call')
M = volume_intersect_sphere_analytical([x',y',z',r'],false);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
For 2 spheres i & j, three cases are possible depending on the distance
d(i,j) of the centers of the spheres i and j.
Case 1: spheres i & j do not overlap, there is no overlap volume
M(i,j)=0;
Condition: d(i,j)>= ri+rj
M(i,j) = 0;
Case 2: spheres i & j fully overlap, the overlap volume has to be
computed.
Condition: d(i,j)<= abs(ri-rj)
M(i,j) = 4/3*pi*min(ri,rj).^3
Case 3: spheres i & j partially overlap, the overlap volume has to be
computed decomposing the overlap volume.