北太天元用根式逼近浮点数

% 北太天元用一个根式逼近浮点数

% ( m sqrt(p) + k sqrt(q) )/ n

% 例如 sin(15度) = ( sqrt(6)-sqrt(2) ) /4

% x = sin(deg2rad(15)) 

% genshiEx1(x);

function [str, m,p,k,q,n] = genshiEx1(x)

   if(nargin ~= 1)

      help genshiEx1;

      return;

   end

   if( abs(x) >= 1)

      fprintf("仅仅对|x|<1 来做逼近，因此下面做的是 rem(%f,1)\n", x);

      x = rem(x,1);

   end   

   % 设定忍量

   tol = 1e-7;

   %暴力搜索

   % 设定暴力搜索的范围

   N = 10;

   for n = 1:N

      for m = 0:n

         for k = 0:n

            for p = 0:3*n

               for q = 0:3*n

                  x0 = (m*sqrt(p) +k*sqrt(q))/n;

                  x1 = (-m*sqrt(p)+k*sqrt(q))/n;

                  %x2 = (m*sqrt(p)-k*sqrt(q))/n;

                  %x3 = (-m*sqrt(p)-k*sqrt(q))/n;

                  A = [x0, x1, -x1, -x0];

                  [r,c] = find( abs(A-x) < tol );

                  if length(c) > 0

                     col = c(1);

                     if(col == 1)

                        fprintf("(%d sqrt(%d)+%d sqrt(%d))/%d\n" , m, p, k, q, n);

                        str = [ '(', num2str(m,0), ' x sqrt(', num2str(p,0), ')', '+', num2str(k,0), ...

                              ' x sqrt(', num2str(q), '))/', num2str(n,0)];      

                     elseif (col == 2)

                        fprintf("(-%d sqrt(%d)+%d sqrt(%d))/%d\n" , m, p, k, q, n);

                        str = [ '(-', num2str(m,0), ' x sqrt(', num2str(p,0), ')', '+', num2str(k,0), ...

                              ' x sqrt(', num2str(q), '))/', num2str(n,0)];      

                     elseif (col == 3)

                        fprintf("(%d sqrt(%d)-%d sqrt(%d))/%d\n" , m, p, k, q, n);

                        str = [ '(', num2str(m,0), ' x sqrt(', num2str(p,0), ')', '-', num2str(k,0), ...

                              ' x sqrt(', num2str(q), '))/', num2str(n,0)];      

                     elseif (col == 4)

                        fprintf("(-%d sqrt(%d)-%d sqrt(%d))/%d\n" , m, p, k, q, n);

                        str = [ '(-', num2str(m,0), ' x sqrt(', num2str(p,0), ')', '-', num2str(k,0), ...

                              ' x sqrt(', num2str(q), '))/', num2str(n,0)];      

                     end

                     return;

                  end

               end

            end

         end

      end

   end

   fprintf("根式逼近没有成功, 你可以试试对x^2 来逼近\n");

   str = num2str(x);

   m =0;

   p = 0;

   k = 0;

   q = 0;

   n = 0;

end