% Problem set 4 Q1 (2nd part)
% Maximum number of breeding triangles.

% The question tells us that all the adjacency matrices are of this form:
% A = [ 0  x1   0  x2  x3 ; ...
%       x1  0   x4 x5  x6 ; ...
%        0  x4  0  x7  x8 ; ...
%       x2  x5  x7  0   0 ; ...
%       x3  x6  x8  0   0 ];
% where the sums of the rows are [ 22 ; 46 ; 12 ; 16 ; 20 ].

% If x is a valid vector to use in A, then it must satisfy
% B*x = [22 ; 46 ; 12 ; 16 ; 20 ], where
B = [ 1 1 1 0 0 0 0 0 ;
    1 0 0 1 1 1 0 0 ;
    0 0 0 1 0 0 1 1 ;
    0 1 0 0 1 0 1 0 ;
    0 0 1 0 0 1 0 1 ];

% To find the maximum number of triangles, let's loop over every possible
% adjacency matrix:
max_triangles = 0;    % keep max triangles
C = zeros( 5 );       % keep corresponding adjacency graph
for x1 = 0:22         % pick breeding number between Wildtype and Big Eyes
    for x2 = 0:16     % pick breeding number between Wildtype and Small ears
        for x5 = 0:16 % pick breeding number between Big eyes and Small ears
            % x1, x2, x5 fix the graph, solve for the other breeding numbers:
            y = B(:,[3 4 6 7 8]) \ ( [22 ; 46 ; 12 ; 16 ; 20 ] - B(:,[1:2 5]) * [ x1 ; x2 ; x5 ] );
            x = [ x1 ; x2 ; y(1:2) ; x5 ; y(3:end) ];
            
            if ( any( x < 0 ) )
                continue
            end
            % Adjacency matrix:
            A = [ 0  x(1)   0  x(2)  x(3) ; ...
                x(1)  0   x(4) x(5)  x(6) ; ...
                0  x(4)  0  x(7)  x(8) ; ...
                x(2)  x(5)  x(7)  0   0 ; ...
                x(3)  x(6)  x(8)  0   0 ];
            
            % Eigenvalues:
            e = eig( A );
            
            % Keep current maximum:
            if ( sum(e.^3)/6 > max_triangles )
                C = A;
                max_triangles = sum(e.^3)/6;
            end
        end
    end
end

% Here are the maximum triangles and one possible adjacency matrix that
% achieves that many triangles:
adjacency_matrix = round( C )
max_triangles = round( max_triangles )
