Duas linhas de uma matriz são ortogonais se seu produto interno for igual a zero. Chame uma matriz com todas as linhas ortogonais pareadas e uma matriz ortogonal . Uma matriz circulante é aquela em que cada vetor de linha é girado um elemento para a direita em relação ao vetor de linha anterior. Nós estaremos interessados ​​apenas em matrizes onde as entradas são -1ou 1.

Tarefa

Escrever código para contar como muitos diferente n/2por n, matrizes circulantes ortogonais quanto possível em 2 minutos (mesmo para n).

Entrada

O código não tem entrada. Ele pode tentar quaisquer valores iguais de nque goste. Por exemplo, o código pode tentar tudo o nque é multiplicado, 4começando pelo menor e também tentar n = 2.

Resultado

O número de matrizes circulantes ortogonais que você encontrou. Também deve haver uma opção simples no seu código para permitir a saída das próprias matrizes.

Ponto

O número de matrizes circulantes que você encontrou.

Dicas

Ortogonais n/2por nmatrizes circulantes só existem quando né um múltiplo de 4ou né menor que 4.

Um exemplo de matriz circulante ortogonal é:

-1  1 -1 -1
-1 -1  1 -1

Dicas para uma abordagem ingênua

A abordagem mais ingênua é apenas repetir todas as matrizes possíveis. Isso pode ser acelerado usando as duas observações a seguir.

• Para testar a ortogonalidade de uma matriz circulante, precisamos comparar apenas cada linha com a primeira. Isso é implementado no código de exemplo.

• Podemos iterar sobre as palavras de Lyndon e, se encontrarmos uma matriz ortogonal, multiplique pelo número de rotações possíveis. Esta é uma idéia ainda não testada, portanto pode ser um erro.

Código de amostra

Esta é uma resposta python muito simples e ingênua. Eu executei usando timeout 120.

import itertools
def check_orthogonal(row):
    for i in xrange(1,int(n/2)):
        if (sum(row[j]*row[(j+i) % n] for j in xrange(n)) != 0):
            return False
    return True

counter = 0
for n in xrange(4,33,4):
    for row in itertools.product([-1,1],repeat = n):
        if check_orthogonal(row):
            counter +=1
            print "Counter is ", counter, ". n = ", n

Testes de correção

Pois n = 4,8,12,16,20,24,28, o número de matrizes distintas que você deve obter é 12,40,144,128,80,192,560, respectivamente.

Níveis de grandiosidade

A julgar pelo código de exemplo, apresento dois níveis de grandiosidade que qualquer resposta pode aspirar.

• A grandiosidade do nível de prata é alcançada obtendo uma pontuação ou 1156 .

• O nível de ouro de grandiosidade é ficar mais alto que isso.

Línguas e bibliotecas

Você pode usar qualquer idioma ou biblioteca que desejar (que não foi projetada para este desafio). No entanto, para fins de pontuação, executarei seu código na minha máquina, portanto, forneça instruções claras sobre como executá-lo no Ubuntu.

Minha máquina Os horários serão executados na minha máquina. Esta é uma instalação padrão do Ubuntu em um processador de 8 GB AMD FX-8350 de oito núcleos. Isso também significa que eu preciso poder executar seu código.

Respostas principais

• 332 por flawr em oitava

• 404 por RT em Python

• 744 por Sample Solution usando pypy

• 1156 por Thomas Kwa usando Java . Nível de prata impressionante!

• 1588 por Reimer Behrends em OCaml . Incrível nível ouro!

OCaml, 1588 (n = 36)

Esta solução usa a abordagem usual de padrão de bits para representar vetores de -1s e 1s. O produto escalar é calculado como de costume, obtendo-se o xor de vetores de dois bits e subtraindo n / 2. Os vetores são ortogonais se seu xor tiver exatamente n / 2 bits definidos.

As palavras de Lyndon não são, por si só, úteis como uma representação normalizada para isso, pois excluem qualquer padrão que seja uma rotação em si. Eles também são relativamente caros de calcular. Portanto, esse código usa uma forma normal um pouco mais simples, que exige que a sequência consecutiva mais longa de zeros após a rotação (ou um deles, se houver múltiplos) ocupe os bits mais significativos. Daqui resulta que o bit menos significativo é sempre 1.

Observe também que qualquer vetor candidato deve ter pelo menos n / 4 (e no máximo 3n / 4). Portanto, consideramos apenas vetores com n / 4 ... n / 2 bits definidos, pois podemos derivar outros via complemento e rotação (na prática, todos esses vetores parecem ter entre n / 2-2 e n / 2 + 2) , mas isso também parece ser difícil de provar).

Construímos essas formas normais a partir do bit menos significativo, observando a restrição de que quaisquer execuções remanescentes de zeros (chamadas de "lacunas" no código) devem seguir nosso requisito de forma normal. Em particular, desde que seja necessário colocar pelo menos mais um bit, deve haver espaço para o intervalo atual e outro que seja pelo menos tão grande quanto o intervalo atual ou qualquer outro intervalo observado até agora.

Também observamos que a lista de resultados é pequena. Portanto, não tentamos evitar duplicatas durante o processo de descoberta, mas simplesmente registramos os resultados em conjuntos por trabalhador e calculamos a união desses conjuntos no final.

Vale ressaltar que o custo de tempo de execução do algoritmo ainda aumenta exponencialmente e a uma taxa comparável à da versão de força bruta; o que isso nos compra é essencialmente uma redução por um fator constante e tem o custo de um algoritmo mais difícil de paralelizar do que a versão de força bruta.

Saída para n até 40:

 4: 12
 8: 40
12: 144
16: 128
20: 80
24: 192
28: 560
32: 0
36: 432
40: 640

O programa está escrito em OCaml, para ser compilado com:

ocamlopt -inline 100 -nodynlink -o orthcirc unix.cmxa bigarray.cmxa orthcirc.ml

Execute ./orthcirc -helppara ver quais opções o programa suporta.

Nas arquiteturas que o suportam, -fno-PICpode oferecer alguns pequenos ganhos adicionais de desempenho.

Isso foi escrito para o OCaml 4.02.3, mas também pode funcionar com versões mais antigas (desde que não sejam muito antigas).

ATUALIZAÇÃO: Esta nova versão oferece melhor paralelização. Observe que ele usa p * (n/4 + 1)threads de trabalho por instância do problema, e alguns deles ainda serão consideravelmente mais curtos que outros. O valor de pdeve ser uma potência de 2. A aceleração em 4-8 núcleos é mínima (talvez em torno de 10%), mas é melhor para um grande número de núcleos n.

let max_n = ref 40
let min_n = ref 4
let seq_mode = ref false
let show_res = ref false
let fanout = ref 8

let bitcount16 n =
  let b2 n = match n land 3 with 0 -> 0 | 1 | 2 -> 1 | _ -> 2 in
  let b4 n = (b2 n) + (b2 (n lsr 2)) in
  let b8 n = (b4 n) + (b4 (n lsr 4)) in
  (b8 n) + (b8 (n lsr 8))

let bitcount_data =
  let open Bigarray in
  let tmp = Array1.create int8_signed c_layout 65536 in
    for i = 0 to 65535 do
      Array1.set tmp i (bitcount16 i)
    done;
    tmp

let bitcount n =
  let open Bigarray in
  let bc n = Array1.unsafe_get bitcount_data (n land 65535) in
  (bc n) + (bc (n lsr 16)) + (bc (n lsr 32)) + (bc (n lsr 48))

module IntSet = Set.Make (struct
  type t = int
  let compare = Pervasives.compare
end)

let worker_results = ref IntSet.empty

let test_row vec row mask n =
  bitcount ((vec lxor (vec lsr row) lxor (vec lsl (n-row))) land mask) * 2 = n

let record vec len n =
  let m = (1 lsl n) - 1 in
  let rec test_orth_circ ?(row=2) vec m n =
    if 2 * row >= n then true
    else if not (test_row vec row m n) then false
    else test_orth_circ ~row:(row+1) vec m n
  in if test_row vec 1 m n &&
        test_orth_circ vec m n then
  begin
    for i = 0 to n - 1 do
      let v = ((vec lsr i) lor (vec lsl (n - i))) land m in
      worker_results := IntSet.add v !worker_results;
      worker_results := IntSet.add (v lxor m) !worker_results
    done
  end

let show vec n =
  for i = 0 to n / 2 - 1 do
    let vec' = (vec lsr i) lor (vec lsl (n - i)) in
    for j = 0 to n-1 do
      match (vec' lsr (n-j)) land 1 with
      | 0 -> Printf.printf "  1"
      | _ -> Printf.printf " -1"
    done; Printf.printf "\n"
  done; Printf.printf "\n"; flush stdout

let rec build_normalized ~prefix ~plen ~gap ~maxgap ~maxlen ~bits ~fn =
  if bits = 0 then
    fn prefix plen maxlen
  else begin
    let room = maxlen - gap - plen - bits in
    if room >= gap && room >= maxgap then begin
      build_normalized
        ~prefix:(prefix lor (1 lsl (plen + gap)))
        ~plen:(plen + gap + 1)
        ~gap:0
        ~maxgap:(if gap > maxgap then gap else maxgap)
        ~maxlen
        ~bits:(bits - 1)
        ~fn;
      if room > gap + 1 && room > maxgap then
        build_normalized ~prefix ~plen ~gap:(gap + 1) ~maxgap ~maxlen ~bits ~fn
    end
  end

let rec log2 = function
| 0 -> -1
| n -> 1 + (log2 (n lsr 1))

let rec test_gap n pat =
  if n land pat = 0 then true
  else if pat land 1 = 0 then test_gap n (pat lsr 1)
  else false

let rec test_gaps n maxlen len =
  let fill k = (1 lsl k) -1 in
  if len = 0 then []
  else if test_gap n ((fill maxlen) lxor (fill (maxlen-len))) then
    len :: (test_gaps n maxlen (len-1))
  else test_gaps n maxlen (len-1)

let rec longest_gap n len =
  List.fold_left max 0 (test_gaps n len len)

let start_search low lowbits maxlen bits fn =
  let bits = bits - (bitcount low) in
  let plen = log2 low + 1 in
  let gap = lowbits - plen in
  let maxgap = longest_gap low lowbits in
  worker_results := IntSet.empty;
  if bits >= 0 then
    build_normalized ~prefix:low ~plen ~gap ~maxgap ~maxlen ~bits ~fn;
  !worker_results

let spawn f x =
  let open Unix in
  let safe_fork () = try fork() with _ -> -1 in
  let input, output = pipe () in
  let pid = if !seq_mode then -1 else safe_fork() in
  match pid with
  | -1 -> (* seq_mode selected or fork() failed *)
     close input; close output; (fun () -> f x)
  | 0 -> (* child process *)
    close input;
    let to_parent = out_channel_of_descr output in
    Marshal.to_channel to_parent (f x) [];
    close_out to_parent; exit 0
  | pid -> (* parent process *)
    close output;
    let from_child = in_channel_of_descr input in
    (fun () -> 
      ignore (waitpid [] pid);
      let result = Marshal.from_channel from_child in
      close_in from_child; result)

let worker1 (n, k) =
  start_search 1 1 n k record

let worker2 (n, k, p) =
  start_search (p * 2 + 1) (log2 !fanout + 1) n k record

let spawn_workers n =
  let queue = Queue.create () in
  if n = 4 || n = 8 then begin
    for i = n / 4 to n / 2 do
      Queue.add (spawn worker1 (n, i)) queue
    done
  end else begin
    for i = n / 2 downto n / 4 do
      for p = 0 to !fanout - 1 do
        Queue.add (spawn worker2 (n, i, p)) queue
      done
    done
  end;
  Queue.fold (fun acc w -> IntSet.union acc (w())) IntSet.empty queue

let main () =
  if !max_n > 60 then begin
    print_endline "error: cannot handle n > 60";
    exit 1
  end;
  min_n := max !min_n 4;
  if bitcount !fanout <> 1 then begin
    print_endline "error: number of threads must be a power of 2";
    exit 1;
  end;
  for n = !min_n to !max_n do
    if n mod 4 = 0 then
    let result = spawn_workers n in
    Printf.printf "%2d: %d\n" n (IntSet.cardinal result);
    if !show_res then
      IntSet.iter (fun v -> show v n) result;
    flush stdout
  done

let () =
  let args =[("-m", Arg.Set_int min_n, "min size of the n by n/2 matrix");
             ("-n", Arg.Set_int max_n, "max size of the n by n/2 matrix");
             ("-p", Arg.Set_int fanout, "parallel fanout");
             ("-seq", Arg.Set seq_mode, "run in single-threaded mode");
             ("-show", Arg.Set show_res, "display list of results") ] in
  let usage = ("Usage: " ^
               (Filename.basename Sys.argv.(0)) ^
               " [-n size] [-seq] [-show]") in
  let error _ = Arg.usage args usage; exit 1 in
  Arg.parse args error usage;
  main ()

Java, 1156 matrizes

Isso usa bitmasking bastante ingênuo e leva menos de 15 segundos para n = 28 na minha máquina.

Matrizes circulantes são determinadas por suas primeiras linhas. Portanto, eu represento as primeiras linhas das matrizes como vetores de bits: 1 e 0 representam -1 e 1. Duas linhas são ortogonais quando o número de bits definidos quando eles são xor'd juntos é n / 2.

import java.util.Arrays;

class Main {

    static void dispmatrix(long y,int N)
    {
        int[][] arr = new int[N/2][N];
        for(int row=0; row < N/2; row++)
        {
            for(int col=0; col < N; col++)
            {
                arr[row][col] = (int) ((y >>> (N+row-col)) & 1L);
            }
        }
        System.out.println(Arrays.deepToString(arr));
    }

  public static void main(String[] args) {

    int count = 0;
    boolean good;
    boolean display = false;
    long y;
    for(int N=4; N <= 28 ;N += 4)
    {
    long mask = (1L << N) - 1;
    for(long x=0; x < (1L<<N); x++)
    {
        good = true;
        y = x + (x << N);

        for(int row = 1; row < N/2; row++)
        {
            good &= N/2 == Long.bitCount((y ^ (y >>> row)) & mask);
        }

        if(good)
        {
            if(display) dispmatrix(y,N);
            count++;
        }

    }
    System.out.println(count);
    }
  }
}

Não consigo fazer o Eclipse funcionar no momento, então isso foi testado no repl.it.

Aqui está o número de primeiras linhas ortogonais às primeiras r linhas depois para n = 28:

[268435456, 80233200, 23557248, 7060320, 2083424, 640304, 177408, 53088, 14896, 4144, 2128, 1008, 1008, 560]

Otimizações:

• Como a ortogonalidade não muda com uma rotação cíclica das duas linhas, basta verificar se a primeira linha é ortogonal ao restante.
• Em vez de fazer manualmente um deslocamento cíclico de N bits N / 2 vezes, armazeno os bits a serem deslocados na parte superior do longe, em seguida, uso um único andcom os N bits inferiores para extrair os necessários.

Possíveis otimizações adicionais:

• Gere palavras de Lyndon. De acordo com meus cálculos, isso só faz sentido se as palavras de Lyndon puderem ser geradas em menos de ~ 1000 ciclos cada.
• Se as palavras de Lyndon demorarem muito, ainda podemos verificar apenas os vetores de bits que começam 00e usar suas rotações (e rotações do NOT) quando encontrarmos uma matriz ortogonal. Podemos fazer isso porque, 0101...01e 1010...10não são possíveis, as primeiras linhas e todas as outras contêm a 00ou a 11.
• Ramifique talvez na metade quando a matriz é provavelmente ortogonal. No entanto, não sei quanto custará a ramificação, então precisarei testar.
• Se o procedimento acima funcionar, comece com uma linha diferente da linha 1?
• Talvez haja alguma maneira matemática de eliminar algumas possibilidades. Não sei o que seria isso.
• Escreva em C ++, é claro.

Python (matrizes 404 no i5-5300U)

Principalmente postar isso como um ponto de partida para outras pessoas melhorarem, isso pode ser muito limpo, paralelizado etc.

import numpy
import itertools
import time

def findCirculantOrthogonalRows(n, t1, timeLimit):
  validRows = []
  testMatrix = numpy.zeros((n//2, n))
  identityMatrixScaled = n*numpy.identity(n//2)
  outOfTime = False
  for startingRowTuple in itertools.product([1, -1], repeat=n):
     for offset in range(n//2):
       for index in range(n):
         testMatrix[offset][index] = startingRowTuple[(index-offset) % n]

     if(numpy.array_equal(identityMatrixScaled, testMatrix.dot(testMatrix.transpose()))):
      validRows.append(startingRowTuple)

    if(time.clock() - t1 >= timeLimit):
      outOfTime = True
      break

  return (validRows, outOfTime)

n = 4
validFirstRows = []
t1 = time.clock()
timeLimit = 120
fullOutput = True

while(True):
  print('calling with', n)
  (moreRows, outOfTime) = findCirculantOrthogonalRows(n, t1, timeLimit)

  if(len(moreRows) > 0):
    validFirstRows.extend(moreRows)

  if(outOfTime == True):
    break

  n += 4

print('Found', len(validFirstRows), 'circulant orthogonal matrices in', timeLimit, 'seconds')

if(fullOutput):
  counter = 1
  for r in validFirstRows:
    n = len(r)
    matrix = numpy.zeros((n//2, n))
    for offset in range(n//2):
      for index in range(n):
        matrix[offset][index] = r[(index-offset) % n]
    print('matrix #', counter, ':\n', matrix)
    counter += 1
    input()

Matlab / Octave, matrizes 381/328

Também apenas a abordagem ingênua, tentando todas as combinações possíveis.

counter = 0;
%ok: 2,4,8
%none: 
tic
for n=[2,4,8,12,16,20];
    m=n/2;
    N=2^n-1;
    for k=1:N

        k/N; %remove ; in order to see progress
        v=(dec2bin(k,n)-'0')*2-1;               %first row

        z=toeplitz([v(1) fliplr(v(m+2:n))], v); %create circulante matrix
        w = z * z.';
        if norm(w - diag(diag(w))) < eps
            counter = counter+1;
            %n    %remove % in order to see the matrices
            %z
        end
        if toc>120;
            break;
        end
    end
end
counter