Código golf uma matriz ortogonal aleatória

9

Uma matriz ortogonal é uma matriz quadrada com entradas reais cujas colunas e linhas são vetores unitários ortogonais (ou seja, vetores ortonormais).

Isso significa que M ^ TM = I, onde I é a matriz de identidade e ^ T significa transposição da matriz.

Observe que isso é ortogonal e não "ortogonal especial", portanto o determinante de M pode ser 1 ou -1.

O objetivo desse desafio não é a precisão da máquina; portanto, se M ^ TM = I estiver dentro de quatro casas decimais, tudo ficará bem.

A tarefa é escrever código que pega um número inteiro positivo n > 1e gera uma matriz ortogonal aleatória n por n . A matriz deve ser escolhida aleatoriamente e uniformemente dentre todas as n por n matrizes ortogonais. Nesse contexto, "uniforme" é definido em termos da medida de Haar, que requer essencialmente que a distribuição não mude se multiplicada por qualquer matriz ortogonal livremente escolhida. Isso significa que os valores da matriz serão valores de ponto flutuante no intervalo de -1 a 1.

A entrada e a saída podem ser de qualquer forma que você achar conveniente.

Por favor, mostre um exemplo explícito do seu código em execução.

Você não pode usar nenhuma função de biblioteca existente que crie matrizes ortogonais. Esta regra é um pouco sutil, então vou explicar mais. Essa regra proíbe o uso de qualquer função existente que receba alguma (ou nenhuma) entrada e produza uma matriz de tamanho de pelo menos n por n que é garantidamente ortogonal. Como exemplo extremo, se você quiser a matriz de identidade n por n, precisará criar você mesma.

Você pode usar qualquer biblioteca gerador de números aleatórios padrão para escolher números aleatórios de sua escolha.

Seu código deve ser concluído dentro de alguns segundos por n < 50.


fonte
Então, é proibido usar a matriz de identidade integrada?
JungHwan Min
@JHM Você não pode usá-lo para criar pelo menos uma matriz de identidade n por n.
Que tal diag? Cria uma matriz diagonal que é de fato ortogonal, mas nem sempre ortonormal.
10136 Karl Napf
Este parece ser um exemplo de "faça X sem Y", que - portanto, o consenso - deve ser evitado.
flawr
11
Matrizes diagonais não são ortogonais; portanto diag, tudo bem.
Angs

Respostas:

7

Haskell, 169 150 148 141 132 131 131 bytes

import Numeric.LinearAlgebra
z=(unitary.flatten<$>).randn 1
r 1=asRow<$>z 1
r n=do;m<-r$n-1;(<>diagBlock[m,1]).haussholder 2<$>z n

Estenda recursivamente uma matriz ortogonal de tamanho n-1adicionando 1 ao canto inferior direito e aplique uma reflexão aleatória de Household. randnfornece uma matriz com valores aleatórios de uma distribuição gaussiana e z dfornece um vetor de unidade uniformemente distribuído em ddimensões.

haussholder tau vretorna a matriz I - tau*v*vᵀque não é ortogonal quando vnão é um vetor de unidade.

Uso:

*Main> m <- r 5
*Main> disp 5 m
5x5
-0.24045  -0.17761   0.01603  -0.83299  -0.46531
-0.94274   0.12031   0.00566   0.29741  -0.09098
-0.02069   0.30417  -0.93612  -0.13759   0.10865
 0.02155  -0.83065  -0.35109   0.32365  -0.28556
-0.22919  -0.41411   0.01141  -0.30659   0.82575
*Main> (<1e-14) . maxElement . abs $ tr m <> m - ident 5
True
Angs
fonte
Fazer a 1×1matriz leva muito espaço para o meu gosto, um caso especial apenas para obter zero, a partir de uma variável aleatória gaussiana: / (Sem ele, há uma chance infinitesimal para obter uma coluna de zero)
Angs
Gosto do seu espírito de torná-lo totalmente correto, mas acho que você pode abandonar esse requisito. No meu código, há também uma chance de que 2 linhas sejam linearmente dependentes e ninguém se importe.
10136 Karl Napf
@KarlNapf bem, eu descobri uma maneira de perder dois bytes de parte qualquer maneira, então problema parcialmente resolvido :)
Angs
Ah bem, excluindo os meus comentários ...
Karl Napf
Sempre feliz quando uma resposta Haskell vence!
4

Python 2 + NumPy, 163 bytes

Agradecemos ao xnor por apontar para usar valores aleatórios distribuídos normais em vez de valores uniformes.

from numpy import*
n=input()
Q=random.randn(n,n)
for i in range(n):
 for j in range(i):u=Q[:,j];Q[:,i]-=u*dot(u,Q[:,i])/dot(u,u)
Q/=(Q**2).sum(axis=0)**0.5
print Q

Usa a ortogonalização de Gram Schmidt em uma matriz com valores aleatórios gaussianos para ter todas as direções.

O código de demonstração é seguido por

print dot(Q.transpose(),Q)

n = 3:

[[-0.2555327   0.89398324  0.36809917]
 [-0.55727299  0.17492767 -0.81169398]
 [ 0.79003155  0.41254608 -0.45349298]]
[[  1.00000000e+00   0.00000000e+00   0.00000000e+00]
 [  0.00000000e+00   1.00000000e+00  -5.55111512e-17]
 [  0.00000000e+00  -5.55111512e-17   1.00000000e+00]]

n = 5:

[[-0.63470728  0.41984536  0.41569193  0.25708079  0.42659843]
 [-0.36418389  0.06244462 -0.82734663 -0.24066123  0.3479231 ]
 [ 0.07863783  0.7048799   0.08914089 -0.64230492 -0.27651168]
 [ 0.67691426  0.33798442 -0.05984083  0.17555011  0.62702062]
 [-0.01095148 -0.45688226  0.36217501 -0.65773717  0.47681205]]
[[  1.00000000e+00   1.73472348e-16   5.37764278e-17   4.68375339e-17
   -2.23779328e-16]
 [  1.73472348e-16   1.00000000e+00   1.38777878e-16   3.33066907e-16
   -6.38378239e-16]
 [  5.37764278e-17   1.38777878e-16   1.00000000e+00   1.38777878e-16
    1.11022302e-16]
 [  4.68375339e-17   3.33066907e-16   1.38777878e-16   1.00000000e+00
    5.55111512e-16]
 [ -2.23779328e-16  -6.38378239e-16   1.11022302e-16   5.55111512e-16
    1.00000000e+00]]

É concluído em um piscar de olhos por n = 50 e alguns segundos em n = 500.

Karl Napf
fonte
Eu não acho isso uniforme. Sua distribuição inicial com um cubo, que tem mais coisas para as diagonais. Gaussianos aleatórios funcionariam porque geram uma distribuição esférica simétrica.
Xnor
@xnor corrigido. Felizmente, isso custou exatamente 1 byte.
Karl Napf
@xnor Ainda mais sorte, este salvou os bytes para-0.5
Karl Napf
Quase, você precisa que a média do normal seja 0, mas isso não é mais do que n.
Xnor
-1

Mathematica, 69 bytes, provavelmente não concorrente

#&@@QRDecomposition@Array[RandomVariate@NormalDistribution[]&,{#,#}]&

QRDecompositionretorna um par de matrizes, a primeira das quais é garantida como ortogonal (e a segunda não é ortogonal, mas triangular superior). Pode-se argumentar que isso tecnicamente obedece à letra da restrição no post: não gera uma matriz ortogonal, mas um par de matrizes.

Mathematica, 63 bytes, definitivamente não concorrente

Orthogonalize@Array[RandomVariate@NormalDistribution[]&,{#,#}]&

Orthogonalizeé inequivocamente proibido pelo OP. Ainda assim, o Mathematica é bem legal né?

Greg Martin
fonte
You may not use any existing library function which creates orthogonal **matrices**.
10136 Karl Napf