高位数的阶乘计算(简单的matlab,Python以及C语言实现)

1 篇文章 0 订阅

今天有同学问我1000万的阶乘用计算机怎么算,所以我花了大概20分钟想了一个思路,分享给大家,主要面向新手,大佬轻喷。
大家都知道,计算机计算的数字精度是有限制的,就整数而言对于不同的CPU,不同的操作系统,不同的语言,不同的编译器等都可能导致计算位数精度的不同。对于一般的编程语言,在常用的平台上处理精度至多不会超过2的63。话不多说,我们来看算法本身:
先说思路,逻辑就是小学教的乘法计算方法。具体操作方法是,第一步把要计算的数的每一位都提出来然后依次放入数组,然后把该数组直接乘以下一个数,对得到的数组进行进位运算,再乘下一个。

function b=jiechen(a)
%将要计算的数的每一位数字依次放入一个一维数组
a_1=a;
len=floor(log10(a))+1;
for i=1:len
    b(i)=mod(a_1,10);
    a_1=floor(a_1/10);
end
%开始计算
for j=a-1:-1:1
    b=b*j;
    len=size(b,2);
    %进位运算
    for i=1:len-1
        if b(i)>=10
            b(i+1)=b(i+1)+floor(b(i)/10);
            b(i)=mod(b(i),10);
        end
    end
    while b(end)>=10
        c=b(end);
        b(end)=mod(b(end),10);
        b=[b,floor(c/10)];
    end
end
%逆序输出
b=b(end:-1:1);

100的阶乘计算结果如下:
在这里插入图片描述
使用matlab自带的factorial,prod,gamma函数计算的高位数阶乘并不准确。如果为了简单验证可以使用Windows计算器
在这里插入图片描述
按照这个思路,你也可以轻松将matlab代码转换成任何你想要的编程语言代码。

用Python也非常简单,下面给出代码:

import math
import numpy as np

a=10000
a_1=a
b=np.array([a_1%10])
len=math.floor(np.log10(a_1))+1
for i in np.arange(1,len):
    a_1=math.floor(a_1/10)
    b=np.append(b,a_1%10)
print(b)
for j in range(a-1,1,-1):
    b=b*j
    len=np.size(b)
    #进位运算
    for i in range(0,len-1):
        if b[i]>=10:
            b[i+1]=b[i+1]+int(b[i]/10)
            b[i]=b[i]%10

    while b[-1]>=10:
        c=b[-1]
        b[-1]=b[-1]%10
        b=np.append(b,int(c/10))
    
#逆序输出
print(list(reversed(b)),np.size(b))

当然实际上用Python不必这么麻烦,Python可以精确处理长整型的各种运算,而且计算速度非常快。
例如:

>>> math.factorial(100)
93326215443944152681699238856266700490715968264381621468592963895217599993229915608941463976156518286253697920827223758251185210916864000000000000000000000000

使用C语言也很简单

#include<stdio.h>
#include<math.h>
#define N 10000

int main()
{
    float a, a_1,len,c ;
    int i, j;
    int b[N]= { 0 };
    printf("请输入要计算阶乘的数字:");
    scanf_s("%f", &a);
    a_1 = a;
    len = floor(log10(a)) + 1;
    for (i = 0;i<len;i++)
    {
        b[i] = (int)a_1%10;
        printf("%d\n", b[i]);
        a_1 = floor(a_1 / 10);
    }
    /*开始计算*/
    for (j = a - 1; j > 0; j--)
    {
        i = N;
        while (b[i - 1] == 0)
        {
            i--;
            if (b[i] != 0)break;
        }
        len = i;
        for (i = 0; i < len; i++)
        {
            b[i] = b[i] * j;
        }
       /*进位运算*/
        for (i = 0; i < N-1 ; i++)
        {
            if (b[i] >= 10)
            {
                b[i + 1] = b[i + 1] + floor(b[i] / 10);
                b[i] = b[i]%10;
            }
        }
    }
    i = N;
    while (b[i - 1] == 0)
    {
        i--;
        if (b[i] != 0)break;
    }
    len = i;
    for (i=len-1;i>=0; i--)
    {
        printf("%d", b[i]);
    }    
}
  • 0
    点赞
  • 4
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值