关于将一个数分解成四个数平方和的算法matlab

偶尔敲代码的创作猿

已于 2022-05-29 04:41:25 修改

阅读量1.2k

点赞数 1

文章标签： matlab 算法矩阵

于 2022-05-29 04:38:50 首次发布

本文链接：https://blog.csdn.net/weixin_55339858/article/details/125026334

版权

怎么把一个大数分解成四个小数的平方和呢？

理论基础

拉格朗日四平方数和定理

每个正整数都可以表示为至多四个正整数的平方和。或者是，每个正整数都可以表示为四个整数的平方和形式。

2 = 0^2 + 0^2 + 1^2 + 1^2

10 = 0^2 + 0^2 +1^2 + 3^2

1234 = 30^2 + 18^2 + 3^2 + 1^2

在历史上也出现过二平方和定理和三平方和定理。证明方式是二次剩余法，或二次互反律，这里不加以证明。

在已知了拉格朗日四平方和定理的基础上，可以通过枚举，累加的方式通过四层循环来寻找满足条件的四个数。

但是在数字比较大的时候，四层循环是很费劲的，十分耗时。那么我们可以用分解质因数的方式，先把大数分解成若干个质因数的乘积，求出质因数的四平方和，然后再通过质因数的四数求出大数大的四数。

高斯恒等式

如果有正整数x和y,

x = x1^2+x2^2+x3^2+x4^2
y = y1^2+y2^2+y3^2+y4^2

那么
x*y = z1^2+z2^2+z3^2+z4^2
z1 = x1*y1+x2*y1+x3*y1+x4*y1
z1 = x1*y2-x2*y1+x3*y4-x4*y3
z1 = x1*y3-x2*y4-x3*y1+x4*y2
z1 = x1*y4+x2*y3-x3*y2-x4*y1

矩阵形式

1 1 1 1
2 -1 4 -3
3 -4 -1 2
4 3 -2 -1
学过高数的同学可以用矩阵乘法推一下。

在有了高斯恒等式之后，我们就可以按照步骤来做了。

操作步骤

第一步分解质因数

分解质因数

分解质因数的方式已经在上一篇文章，讲欧拉定理的时候讲过了。感兴趣的同学可以去看一看。还是一样的方法，先高斯筛质数，然后通过能不能整除，来判断质因数。

%计算一个数的质因数
function e = EULER_E(n)
isnot_prime = zeros(1 , n);%判断是否是合数
prime = [];%质数集
primen = [];%质因数集
if(n == 1)
    e = 1;
    return;
end
for i = 2:n
    if (isnot_prime(i) == 0)
        prime(end + 1) = i;
    end
    for j = 1 : length(prime)
        if((i*(prime(j))>n))
            break;
        end
        isnot_prime(i*(prime(j))) = 1;
        if(mod(i , prime(j)) == 0)
            break;
        end
    end
end
for k = 1:length(prime)
    if(mod(n , prime(k)) == 0)
        primen(end + 1) = prime(k);
    end
end
e = primen;

下一步，求解质因数的四平方数。

求解四平方数

用四层循环嵌套的方式来求。通过暴力破解的累加法。同时为了节省时间，可以只循环到sqrt(n)，即根号n,同时也可以用字典排列的方式来节省时间，比如第一层循环为i = 1：n,第二层可以是j = 1:i;依此类推。

代码如下：

%计算一个数的四个平方和数（拉格朗日四平方数定理）
function f = FOUR_number_Square_and(n)
for i = 0 : floor(sqrt(n))
    for j = 0 : floor(sqrt(n))
        for k = 0 : floor(sqrt(n))
            for l = 0 : floor(sqrt(n))
                if((i^2 + j^2 + k^2 + l^2) == n)
                    f = [i , j , k , l];
                end
            end
        end
    end
end

第三步，应用高斯恒等式。

应用高斯恒等式

这一步内容比较多。首先要判断n = q1^a1+q2^a2+...+qn^an ，各个质因数的次数。

然后，用高斯恒等式，分别先将每个质因数都乘一次，然后再去乘以多余的次数。

代码如下：

%拉格朗日四平方数定理实验
clc;clear;
w = input('n = ');FNSA = [];FNSA_X = [];n = w;prime = [];
primen = EULER_E(w);count = zeros(1 ,length(primen));
for i = 1 : length(primen)
    k = 0;
    while (mod(n , primen(i)) == 0)
        n = n/primen(i);
        k = k + 1;
    end
    count(i) = k;
end
for i = 1 : length(primen)
    temp = FOUR_number_Square_and(primen(i));
    for j = 1 : 4
        FNSA(end+1) = temp(j);
    end
end
for i = 1 : length(primen)
    if(i == 1)
        temp = primen(1);
        for k = 1:4
            FNSA_X(end+1) = FNSA(k);
            tempn = FNSA_X;
        end
        continue;
    end
    temp = primen(i)*temp;
    FNSA_X(1) = abs(FNSA(4*i-3)*tempn(1)+FNSA(4*i-2)*tempn(2)+FNSA(4*i-1)*tempn(3)+FNSA(4*i)*tempn(4));
    FNSA_X(2) = abs(FNSA(4*i-3)*tempn(2)-FNSA(4*i-2)*tempn(1)+FNSA(4*i-1)*tempn(4)-FNSA(4*i)*tempn(3));
    FNSA_X(3) = abs(FNSA(4*i-3)*tempn(3)-FNSA(4*i-2)*tempn(4)-FNSA(4*i-1)*tempn(1)+FNSA(4*i)*tempn(2));
    FNSA_X(4) = abs(FNSA(4*i-3)*tempn(4)+FNSA(4*i-2)*tempn(3)-FNSA(4*i-1)*tempn(2)-FNSA(4*i)*tempn(1));
    tempn = FNSA_X;
end
for j = 1 : length(count)
    while(count(j)>1)
        count(j) = count(j) - 1;
        FNSA_X(1) = abs(FNSA(4*j-3)*tempn(1)+FNSA(4*j-2)*tempn(2)+FNSA(4*j-1)*tempn(3)+FNSA(4*j)*tempn(4));
        FNSA_X(2) = abs(FNSA(4*j-3)*tempn(2)-FNSA(4*j-2)*tempn(1)+FNSA(4*j-1)*tempn(4)-FNSA(4*j)*tempn(3));
        FNSA_X(3) = abs(FNSA(4*j-3)*tempn(3)-FNSA(4*j-2)*tempn(4)-FNSA(4*j-1)*tempn(1)+FNSA(4*j)*tempn(2));
        FNSA_X(4) = abs(FNSA(4*j-3)*tempn(4)+FNSA(4*j-2)*tempn(3)-FNSA(4*j-1)*tempn(2)-FNSA(4*j)*tempn(1));
        tempn = FNSA_X;
    end
end
temp = 0;
for j = 1:4
    temp = temp + FNSA_X(j)^2;
end
if(temp == w)
    fprintf('r %d %d %d %d',FNSA_X(1),FNSA_X(2),FNSA_X(3),FNSA_X(4));
else
    fprintf('w --%d --%d %d %d &d',temp ,FNSA_X(1),FNSA_X(2),FNSA_X(3),FNSA_X(4));
end

小结

对于高斯恒等式在代码里通过矩阵表示会简洁一些。

完整代码就是上面三个代码。两个函数，一个主程序。

另外在编辑的过程中，写高斯恒等式来费劲啦！！！小编就写了一个程序来输出高斯恒等式，为了减少没及时保存而重打一次的工作量。算是没用的小代码吧。

高斯恒等式输出代码

%打印高斯恒等式
clc;clear;
A = [1,1,1,1;2,-1,4,-3;3,-4,-1,2;4,3,-2,-1];
disp(A);
fprintf('x = ');
for i = 1 : 4
    fprintf('x%d^2',i);
    if(i ~= 4)
        fprintf('+');
    else
        fprintf('\n');
    end
end
fprintf('y = ');
for j = 1 : 4
    fprintf('y%d^2',j);
    if(j ~= 4)
        fprintf('+');
    else
        fprintf('\n');
    end
end
 fprintf('x*y = ');
 for k = 1  : 4
    fprintf('z%d^2',k);
    if(k ~= 4)
        fprintf('+');
    end
 end
for i = 1 : 4
    for j = 1 : 4
        if (j == 1)
            fprintf('\n z%d = ',j);
        end
        if(j ~= 1)
            if(A(i , j)>0)
                fprintf('+');
            else
                fprintf('-');
            end
        end
        fprintf('x%d*y%d',j,abs(A(i , j)));
    end
end

输出结果

     1     1     1     1
     2    -1     4    -3
     3    -4    -1     2
     4     3    -2    -1

x = x1^2+x2^2+x3^2+x4^2
y = y1^2+y2^2+y3^2+y4^2
x*y = z1^2+z2^2+z3^2+z4^2
z1 = x1*y1+x2*y1+x3*y1+x4*y1
z1 = x1*y2-x2*y1+x3*y4-x4*y3
z1 = x1*y3-x2*y4-x3*y1+x4*y2
z1 = x1*y4+x2*y3-x3*y2-x4*y1

对了，补充一下源程序的运行结果：

运行结果

n = 12345678
r 2422 348 2507 271

希望各位有所收获。瑞斯拜。如果有不懂的欢迎留言评论。

偶尔敲代码的创作猿

关注

1
点赞
踩
7

收藏

觉得还不错? 一键收藏
打赏
0
评论
关于将一个数分解成四个数平方和的算法matlab

怎么把一个大数分解成四个小数的平方和呢？拉格朗日四平方数和定理：每个正整数都可以表示为至多四个正整数的平方和。或者是，每个正整数都可以表示为四个整数的平方和形式。2 = 0^2 + 0^2 + 1^2 + 1^210 = 0^2 + 0^2 +1^2 + 3^21234 = 30^2 + 18^2 + 3^2 + 1^2在历史上也出现过二平方和定理和三平方和定理。证明方式是二次剩余法，或二次互反律，这里不加以证明。在已知了拉格朗日四平方和定理的基础上，可以通过枚举，累加的方式通过
复制链接

扫一扫