MthSc 360: Project 2


Here are the rules for Project 2.

The goal of this project is to write a program which will enumerate all of the 4-by-4 magic squares constructed from the numbers 1 through 16. You are to write three variations of this search program, each of which should incorporate techniques that we have studied in class to make the search execute significantly faster.

Start with program Magic4a.m which has been sent to you by email, and which is listed below.

  1. The first version, Magic4b.m, should simply introduce pruning. Do the tests as early as possible. For example, as soon as we have selected k1 thru k8 we can check whether the first two rows have the same sum. The lines
    for k8 = 1:n
       if (~used(k8))
          used(k8) = 1;
    
    can simply be changed to
    for k8 = 1:n
       S1 = k5 + k6 + k7 + k8;
       if (~used(k8)) & (S0 == S1)
          used(k8) = 1;
    
    This assumes that we have already calculated S0 = k1 + k2 + k3 + k4. Five of the remaining 8 tests can be done earlier in the program

  2. The second version, Magic4c.m, should incorporate the consequences of our mathematical insight. Where ever possible calculate the value that must go in a specific square. For example, the lines (from Magic4b.m)
    for k8 = 1:n
       S1 = k5 + k6 + k7 + k8;
       if (~used(k8)) & (S0 == S1)
          used(k8) = 1;
    
    can simply be changed to
    k8 = S0 - k5 - k6 - k7;
    if (1 <= k8) & (k8 <= n)
       if (~used(k8))
          used(k8) = 1;
    
    Note that as we add direct calculations we can eliminate the tests, since we know mathematically that the conditions are satisfied, and that the tests would be redundant.

    In your writeup, you must include a mathematical justification for each calculation and any tests that are eliminated. The mathematics here is very elementary.

  3. The third version, Magic4d.m, should exploit the symmetries of the problem by only checking one representative of each equivalence class. For example, we know that we can rotate the smallest corner to the k1 position. This means that the search on k1 can be changed from
    for k1 = 1:n
    
    to
    for k1 = 1:(n-3)
    
    and that the searches on k4 (as well as k13 and k16) can be changed from
    for k4 = 1:n
    
    to
    for k4 = k1+1:n
    

    When exploiting symmetries, it will probably be useful to rearrange the order in which you search the squares. It will also be advantageous to exploit all the symmetries that you can justify.

    In your writeup, you must include a description of the symmetries that you are exploiting and an explanation of the way you are choosing the representative from each equivalence class.

  4. When you get each version running correctly run it with "PrintOut = 1" and then stop it after it has generated about 10 solutions. Print the command window showing the solutions. This will take quite a while for Magic4b.m. Then, for Magic4c.m and Magic4d.m, set "PrintOut = 0" and run the program to completion. Print the command window showing the total number of magic squares and the time it took your program to find them. Don't bother trying to run Magic4b.m to completetion. However, you should try to time how long it took to get the first few answers and use that time to estimate how long Magic4b.m would take to run to completetion.

  5. Turn in your writeup, a printed copy of each version of your program, and the printouts of the command windows for each version.


%  This program generates all the 4-by-4 magic squares.  
%  It does so by the brute force method of generating all 
%  16! squares containing each of the permutations of the
%  numbers from 1 to 16.  It then checks each completed
%  square to see if it is magic.  
%
%  This version does not use pruning; it does note use
%  the consequences from linear algebra; and it does not
%  exploit any symmetries.
%
%  This program will not succeed in determining all the
%  4-by-4 magic squares in your lifetime.  However, it is
%  a useful skeleton for smarter programs

%  Label the square with the following variables
%    k1   k2   k3   k4
%    k5   k6   k7   k8
%    k9   k10  k11  k12
%    k13  k14  k15  k16


% Set up a counter for the number of solutions.
count = 0

% Set PrintOut to 1 to print out the solutions.
% Set PrintOut to 0 to suppress printing the solutions.
PrintOut = 1

% Start the clock
Tstart = clock

n = 16;
used = zeros(1,n);

for k1 = 1:n
  used(k1) = 1;
  for k2 = 1:n
    if (~used(k2))
      used(k2) = 1;
      for k3 = 1:n
        if (~used(k3))
          used(k3) = 1;
          for k4 = 1:n
            if(~used(k4))
              used(k4) = 1;
              for k5 = 1:n
                if (~used(k5))
                  used(k5) = 1;
                  for k6 = 1:n
                    if (~used(k6))
                      used(k6) = 1;
                      for k7 = 1:n
                        if(~used(k7))
                          used(k7) = 1;
                          for k8 = 1:n
                            if (~used(k8))
                              used(k8) = 1;
                              for k9 = 1:n
                                if (~used(k9))
                                  used(k9) = 1;
                                  for k10 = 1:n
                                    if(~used(k10))
                                      used(k10) = 1;
                                      for k11 = 1:n
                                        if (~used(k11))
                                          used(k11) = 1;
                                          for k12 = 1:n
                                            if (~used(k12))
                                              used(k12) = 1;
                                              for k13 = 1:n
                                                if(~used(k13))
                                                  used(k13) = 1;
                                                  for k14 = 1:n
                                                    if (~used(k14))
                                                      used(k14) = 1;
                                                      for k15 = 1:n
                                                        if (~used(k15))
                                                          used(k15) = 1;
                                                          for k16 = 1:n
                                                            if(~used(k16))
                                                              used(k16) = 1;
                                                              S0 = k1  + k2  + k3  + k4;
                                                              S1 = k5  + k6  + k7  + k8;
                                                              S2 = k9  + k10 + k11 + k12;
                                                              S3 = k13 + k14 + k15 + k16;
                                                              S4 = k1  + k5  + k9  + k13;
                                                              S5 = k2  + k6  + k10 + k14;
                                                              S6 = k3  + k7  + k11 + k15;
                                                              S7 = k4  + k8  + k12 + k16;
                                                              S8 = k1  + k6  + k11 + k16;
                                                              S9 = k4  + k7  + k10 + k13;
                                                              if (S0==S1)&(S0==S2)&(S0==S3)&(S0==S4)&(S0==S5)&(S0==S6)&(S0==S7)&(S0==S8)&(S0==S9)
                                                                if PrintOut
                                                                  M = [k1 k2 k3 k4];
                                                                  M = [M; k5 k6 k7 k8];
                                                                  M = [M; k9 k10 k11 k12];
                                                                  M = [M; k13 k14 k15 k16]
                                                                end
                                                                count = count + 1;
                                                              end
                                                              used(k16) = 0;
                                                            end
                                                          end
                                                          used(k15) = 0;
                                                        end
                                                      end
                                                      used(k14) = 0;
                                                    end
                                                  end
                                                  used(k13) = 0;
                                                end
                                              end
                                              used(k12) = 0;
                                            end
                                          end
                                          used(k11) = 0;
                                        end
                                      end
                                      used(k10) = 0;
                                    end
                                  end
                                  used(k9) = 0;
                                end
                              end
                              used(k8) = 0;
                            end
                          end
                          used(k7) = 0;
                        end
                      end
                      used(k6) = 0;
                    end
                  end
                  used(k5) = 0;
                end
              end
              used(k4) = 0;
            end
          end
          used(k3) = 0;
        end
      end
      used(k2) = 0;
    end
  end
  used(k1) = 0;
end

% Stop the clock.
Tstop = clock
fprintf('Found %g solutions in %g seconds\n', count, etime(Tstop,Tstart))


Professor: Daniel D. Warner, Mathematical Sciences, Clemson University
Send questions and comments to: warner@math.clemson.edu
Last Updated: September 2, 1996