#!/usr/bin/octave -q
% generator różnych prostych rozkładów
% (całkowitych i ułamków zwyczajnych)
% do uruchomienia w środowisku Octave/Matlab
% formatowanie wyników zgodne ze składnią LaTeXa
%
% wersja z 2015.12.11
function a = printarray(array,name,fid,normalize = 0)
  fprintf (fid, "%s", name)
  if normalize!=1
    fprintf (fid, " = ")
  end
  fprintf (fid, "\\left[ \\begin{array}{rr} ")
  printf ("\n%s\n",name)
  if fid == 1
    fids = [fid];
  else
    fids = [fid 1];
  end
  for f = fids
    for i = 1:size(array)(1)
      for j = 1:size(array)(2)
        if (abs(round(array(i,j))-array(i,j))>3*eps) && normalize == 0
          if array(i,j)<0
            sign="-";
            snum = array(i,j)*-1;
          else
            sign = "";
            snum = array(i,j);
          end
          s = rats(snum);
          s = strrep(s,"/","}{");
          if f!=1
            fprintf (f, "%s\\frac{%s} ", sign,strtrim(s) )
          else
            fprintf (f, "%6s ", rats(snum) )
          end
        elseif normalize
          if testint(array(i,j))
            if array(i,j) == 0
              fprintf (f, "%d" , 0 )
            else
              fprintf (f, "\\frac{1}{%d} ", array(i,j))
            end
          else
            fprintf (f,"\\frac{1}{\\sqrt{%d}} ", array(i,j)^2 )
          end
        else
          fprintf (f, "%6d ", array(i,j))
        end
          if j<size(array)(2)
            if f!=1
              fprintf (f," &  ")
          else
            fprintf (f,"  ")
          end
        else
          if f!=1
            fprintf (f," \\\\  ")
          else
            fprintf (f," \n")
          end
          end
      end
    end
    if f!=1
      fprintf (f, "\\end{array}\\right]")
    end
  end
end 
function t = testint(real)
  t = abs(round(real)-real)<3*eps;
end
row = 3 ;
col = 3 ;
count = 12 ;
FORCE_SYMMETRIC        = 0;
FORCE_WITH_PERMUTATION = 0;
nr_total = 0 ;
for dimension = [2 3]
  row = col =dimension;
  if dimension == 2
    how_many_zeros = 0;
  end
  if dimension == 3
    how_many_zeros = 2;
  end
  if dimension == 4
    how_many_zeros = 20;
  end
  fileout = ["rozklady", int2str(dimension), ".dat"];
  fid = fopen (fileout,"w");
  fprintf (fid,"<ul>\n")
  array_of_A=[];
  for nr = 1:count
    nr_total++;
    if nr <=count/2
      FORCE_SYMMETRIC = 0;
    else
      FORCE_SYMMETRIC = 1;
    end
    if mod(abs(nr)-count/2,count/2) < count/4
      FORCE_INTEGER_EIG = 0;
    else
      FORCE_INTEGER_EIG = 1;
    end
    if (dimension>2)
      FORCE_INTEGER_EIG =0
    end
    ii=0;
    do
    ii++;
    A = randi([-9 9],row,col);
    if FORCE_SYMMETRIC == 1
      for i = 1:size(A)(1)
        for j = i+1:size(A)(2)
          A(i,j) = A(j,i);
      end
    end
  else
    if (A == A')
      continue;
    end
  end
  if sum(sum(A == 0)) > how_many_zeros
    continue;
    if any(all(all(A == array_of_A)))
      continue;
    end
  end
  if FORCE_INTEGER_EIG == 1
    [v lambda]=eig(A);
    if !all(all(testint(48*lambda)))
      continue
    end
  end
  [Q R]=qr(A);
  flag = 0;
  for n = 1:col
    if R(n,n) <0
      flag = 1 ;
      break
    end
  end
  if flag == 1
    continue
  end
  if all(all(testint(R)))
    [L, U, P] = lu(A);
    if (diag(U)!=0)
      if all(diag(L)==1)
        if all(diag(P)==1) && FORCE_WITH_PERMUTATION == 0
          break;
        end
        if !all(diag(P)==1) && FORCE_WITH_PERMUTATION == 1
          break
        end
      end
    end
end
until 0
array_of_A(:,:,nr) = A;
D     = diag([diag(U)']);
Uprim = inv(D)*U;
if FORCE_SYMMETRIC == 1
   if Uprim == L'
   end
end
fprintf (fid,"\t<li>\n\t\t<ul id=\"rozwin\">\n")
fprintf (fid,"\t\t\t<li>\n\t\t\t<a>Przykład %d — $A_{%dx%d}$", nr_total,row,col)
if FORCE_SYMMETRIC
  fprintf (fid, " [m. symetryczna]")
end
if FORCE_WITH_PERMUTATION
  fprintf (fid, " NIEZBĘDNA PERMUTACJA")
end
fprintf (fid,"</a>\n\t\t\t<div class=\"rozwin\">\n")
fprintf(fid,"$$")
printarray(A,"A",fid)
fprintf(fid,"$$")
fprintf (fid,"\n\t\t\t</div>\n\t\t\t</li>\n")
fprintf (fid,"\t\t\t | <li>\n\t\t\t<a>$LU$</a>\n")
fprintf (fid,"\t\t\t<div class=\"rozwin\">")
fprintf(fid,"$$")
printarray(L,"L",fid)
fprintf(fid,";\\quad ")
printarray(U,"U",fid)
fprintf(fid,"$$")
fprintf (fid,"\n\t\t\t</div>\n\t\t\t</li>")
if (Uprim == L')
  fprintf (fid,"\t\t\t | <li>\n\t\t\t<a>$LDU'=LDL^T$</a>\n")
else
  fprintf (fid,"\t\t\t | <li>\n\t\t\t<a>$LDU'$</a>\n")
end
fprintf (fid,"\t\t\t<div class=\"rozwin\">")
fprintf(fid,"$$")
printarray(L,"L",fid)
fprintf(fid,";\\quad ")
printarray(D,"D",fid)
fprintf(fid,";\\quad ")
if (Uprim == L')
  printarray(Uprim,"U'=L^T",fid)
else
  printarray(Uprim,"U'",fid)
end
fprintf(fid,"$$")
fprintf (fid,"\n\t\t\t</div>\n\t\t\t</li>")
fprintf (fid,"\t\t\t | <li>\n\t\t\t<a>$QR$</a>\n")
fprintf (fid,"\t\t\t<div class=\"rozwin\">")
fprintf(fid,"$$")
printarray(Q,"Q",fid)
fprintf(fid,";\\quad ")
printarray(R,"R",fid)
fprintf(fid,"$$")
fprintf (fid,"\n\t\t\t</div>\n\t\t\t</li>")
if FORCE_INTEGER_EIG == 1
  for i = 1 : col
    [~,inx] = min(abs(v(:,i)));
    v(:,i) = v(:,i)/v(inx,i);
  end
fprintf (fid,"\t\t\t | <li>\n\t\t\t<a>$S \\Lambda S^{-1}$</a>\n")
fprintf (fid,"\t\t\t<div class=\"rozwin\">")
fprintf(fid,"$$")
printarray(v,"S",fid)
fprintf(fid,";\\quad ")
printarray(lambda,"\\Lambda",fid)
fprintf(fid,";\\quad ")
printarray(inv(v),"S^{-1}",fid)
fprintf(fid,"$$")
fprintf (fid,"\n\t\t\t</div>\n\t\t\t</li>")
if FORCE_SYMMETRIC ==1
fprintf (fid,"\t\t\t | <li>\n\t\t\t<a>$Q \\Lambda Q^{T}$</a>\n")
fprintf (fid,"\t\t\t<div class=\"rozwin\">")
fprintf(fid,"$$")
M = sqrt(diag([ diag(v'*v)]));
printarray(M,"Q = S\\cdot ",fid,1)
fprintf(fid,"$$")
fprintf (fid,"\n\t\t\t</div>\n\t\t\t</li>")
end
end
  fprintf (fid,"\n\t\t\t</div>\n\t\t\t</li>")
  fprintf (fid,"\n\t\t</ul>\n")
  fprintf (fid,"\t</li>\n")
end
  if fid != 1
    fprintf (fid,"</ul>\n")
    fclose(fid);
  end
end