天道酬勤,学无止境

使用费马定理的矩阵求幂(Matrix exponentiation using fermat's theorem)

问题

就像我们使用费马小定理进行模幂运算一样,我只是想知道有没有这样的快速矩阵求幂方法? 我们可以使用费马定理进行矩阵求幂吗? 如果没有,那么有没有比分治求幂法更快的方法?

回答1

受限制的 HTML

  • 允许的HTML标签:<a href hreflang> <em> <strong> <cite> <blockquote cite> <code> <ul type> <ol start type> <li> <dl> <dt> <dd> <h2 id> <h3 id> <h4 id> <h5 id> <h6 id>
  • 自动断行和分段。
  • 网页和电子邮件地址自动转换为链接。

相关推荐
  • 解释基于费马小定理检查素性的代码(Explain a code to check primality based on Fermat's little theorem)
    问题 我发现了一些声称基于费马小定理检查素性的 Python 代码: def CheckIfProbablyPrime(x): return (2 << x - 2) % x == 1 我的问题: 它是如何工作的? 它与费马小定理有什么关系? 这种方法有多准确? 如果它不准确,使用它有什么好处? 我在这里找到了。 回答1 1.它是如何工作的? 费马小定理说,如果一个数x是素数,那么对于任何整数a : 如果我们将两边除以a ,那么我们可以将方程改写如下: 我打算证明这是如何工作的(你的第一个问题),因为在这个 wiki 页面和一些谷歌搜索下有很多很好的证据(比我能提供的更好)。 2. 代码与定理的关系 因此,您发布的函数检查是否(2 << x - 2) % x == 1 。 首先, (2 << x-2)与编写2**(x-1)或数学形式相同: 那是因为<<是逻辑左移运算符,这里有更好的解释。 位移位和乘以 2 的幂之间的关系特定于数字在计算机上的表示方式(以二进制形式),但这一切都归结为 我可以从两边的指数中减去 1,得到 现在,我们从上面知道对于任何数字 a , 假设a = 2 。 这给了我们 好吧,这与2 << (x-2) ! 那么我们可以这样写: 这导致了最终的关系: 现在,mod 的数学版本看起来有点奇怪,但我们可以编写如下等效代码: (2 << x - 2) % x == 1
  • C++ Fermat‘s little theorem费马小定理寻找模逆实现算法(附完整源码)
    C++ Fermat's little theorem费马小定理寻找模逆实现算法 C++ Fermat's little theorem费马小定理寻找模逆实现算法完整源码(定义,实现,main函数测试) C++ Fermat’s little theorem费马小定理寻找模逆实现算法完整源码(定义,实现,main函数测试) #include <iostream> #include <vector> /** Recursive function to calculate exponent in \f$O(\log n)\f$ using binary * exponent. */ int64_t binExpo( 来源:https://blog.csdn.net/it_xiangqiang/article/details/114967019
  • 费马的素数小定理 Java(Fermat's little theorem for prime numbers Java)
    问题 据我所知,如果 P 是素数,那么a^(p-1)-1 mod p =1我想要的是打印两个整数之间的每个素数我写的是: public static void main(String[] args) { Scanner r = new Scanner(System.in); int x = Integer.parseInt(r.nextLine()); for (int i = 0; i < x; i++) { String s = r.nextLine(); BigInteger n = new BigInteger(s.split(" ")[0]); BigInteger m = new BigInteger(s.split(" ")[1]); for (BigInteger j = n; j.compareTo(m) <= 0; j.add(BigInteger.ONE)) { if (isPrime(j)) { System.out.println(j); } } } } private static boolean isPrime(BigInteger num) { BigInteger a = num.subtract(num.divide(new BigInteger("2"))); a = a.modPow(num.subtract(BigInteger.ONE)
  • RSA算法的数学基础
    前言 最近学校开了一门mathematics for cryptography的课,主要讲了一些密码学中常用到的数学原理。CSDN上关于RSA的文章中,只看到了关于RSA算法的一些讲解,对于其涉及的数学原理并没有多少介绍,因此想写点东西与大家分享,同时也是为了自己复习和巩固。第一次写文章,可能会有一些问题和不足,希望大家能够指正。此外,这篇文章主要参考了我们老师Zülfükar Saygı(TOBB ETÜ - University of Economics & Technology)的PPT课件,《信息安全原理及应用》第三版(熊平,朱天清)以及《信息安全数学基础——算法、应用与实践》第二版(任伟)。 1.RSA简介 RSA是1977年由麻省理工学院的Ron Rivest,Adi Shamir和Leonard Adlenman提出的非对称公开密钥密码算法。其可靠性基于大数的因子分解问题。目前为止,仍未找到快速分解因子的算法,所以还没有任何可靠的攻击RSA算法的方式。只要密钥的长度足够长,用RSA加密的信息实际上是不可能破解的。因此,RSA算法已经成为目前应用最广泛的公钥密码。 RSA算法的加密过程如下: 选择两个足够大素数p、q计算n=p·q Ф=(p-1)(q-1)选择一个数e,满足1<e<Ф,且gcd(e,Ф)=1计算d,使得e·d≡1 mod Ф 其中{n,e}为公钥,{n
  • 高斯整数 / 费马平方和定理 / 拉格朗日的四平方定理
    一、高斯整数 形如 。(其中a,b是整数)的复数被称为高斯整数,高斯整数全体记作Z[i]。注意到若 γ=a+bi 是高斯整数,则它是满足如下方程的代数整数: 通常我们使用希腊字母来表示高斯整数,例如α,β,γ和δ。注意到若 n 是一个整数,则 n=n+0i 也是高斯整数。当我们讨论高斯整数的时候,把通常的整数称为有理整数。 加、减、乘运算 高斯整数在加、减、乘运算下是封闭的,正如下面定理所述。 定理1:设 α=x+iy 和 β=w+iz 是高斯整数,其中 x,y,w 和 z 是有理整数,则 α+β,α-β 和 αβ 都是高斯整数。 虽然高斯整数在加、减和乘运算下封闭,但是他们在除法运算下并不封闭,这一点与有理整数类似。此外,若 α=a+bi 是高斯整数,则a的范数 N(α)= 是非负有理整数。 整除性 我们可以像研究有理整数那样去研究高斯整数。整数的许多基本性质可以直接类推到高斯整数上。要讨论高斯整数的这些性质,我们需要介绍高斯整数类似于通常整数的一些概念。特别地,我们需要说明一个高斯整数整除另一个高斯整数的意义,并给出高斯素数的定义。 定义1:设 α 和 β 是高斯整数,我们称α整除β,是指存在一个高斯整数 γ 使得β=αγ。若α整除β,我们记作α|β ;若α 不整除β ,记作αβ 。 高斯整数的整除也满足有理整数整除的一些相同的性质。例如,若α,β和γ 是高斯整数,α|β,β
  • Matrix exponentiation using fermat's theorem
    Like we use fermat's little theorem for modular exponentiation , I was just wondering that is there any such method for fast matrix exponentiation ? Can we use fermat's theorem for matrix exponentiation ? If no then is there any faster method than divide and conquer method of exponentiation ?
  • 心跳调试(后缀和,费马小定理!!!(难))
    传送门 题目描述 Heartache Debug 给定正整数 n,正整数序列 w1,w2,……,wn。 你的心跳是一个栈,初始时栈中仅有一个元素 1。 需要进行若干次如下的「调试」操作: 示例1 输入 3 1 1 1 输出 500000004 示例2 输入 3 2 1 3 输出 250000002 备注: 1≤n≤106, 1≤wi≤109。 解题关键:我们要先了解逆元(逆元的求解思路大概有三种,仔细琢磨吧,难!!!)的相关知识,运用费马小定理进行除法的取余 借助后缀和的思想求概率问题 #include <bits/stdc++.h> using namespace std; typedef long long ll; const ll mod=1e9+7; ll w[1000010],zhui[1000010]; ll qpow(ll a,ll n){ ll ans=1; a=(a%mod+mod)%mod; for(;n;n>>=1){ if(n&1)ans=(ans*a)%mod; a=(a*a)%mod; } return ans; } int main() { int n; ll ans=1; scanf("%d",&n); for(int i=0;i<n;i++)scanf("%lld",&w[i]); for(int i=n-1;i>=0;i--) zhui[i]
  • HDU 4549M斐波那契数列(矩阵快速幂+费马小定理)
    M斐波那契数列 Time Limit: 3000/1000 MS (Java/Others) Memory Limit: 65535/32768 K (Java/Others) Total Submission(s): 993 Accepted Submission(s): 293 Problem Description M斐波那契数列F[n]是一种整数数列,它的定义如下: F[0] = a F[1] = b F[n] = F[n-1] * F[n-2] ( n > 1 ) 现在给出a, b, n,你能求出F[n]的值吗? Input 输入包含多组测试数据; 每组数据占一行,包含3个整数a, b, n( 0 <= a, b, n <= 10^9 ) Output 对每组测试数据请输出一个整数F[n],由于F[n]可能很大,你只需输出F[n]对1000000007取模后的值即可,每组数据输出一行。 Sample Input 0 1 0 6 10 2 Sample Output 0 60 Source 2013金山西山居创意游戏程序挑战赛——初赛(2) 题目大意: 题目易懂,主要是数据太大。自己了解的降幂公式没有这么吊。这个也不需要讨论。如果mod为质数的话,根据费马小定理可以得到:a^p%(mod)=a^(p%(mod-1)). mod为质数 解题思路:根据题意,慢慢可以推公式
  • 使用 scipy 进行矩阵求幂:expm、expm2 和 expm3(Matrix exponentiation with scipy: expm, expm2 and expm3)
    问题 矩阵求幂可以在 python 中使用scipy.linalg库中的函数执行,即expm, expm2, expm3 。 expm使用 Pade 近似; expm2使用特征值分解方法, expm3使用默认项数为 20 的泰勒级数。 在 SciPy 0.13.0 发行说明中指出: 不推荐使用矩阵指数函数 scipy.linalg.expm2 和 scipy.linalg.expm3。 所有用户都应该改用数值上更强大的 scipy.linalg.expm 函数。 尽管expm2和expm3自发布版本 SciPy 0.13.0 以来expm3被弃用,但我发现在许多情况下,这些实现比expm快。 由此,产生了一些问题: 在什么情况下 expm2 和 expm3 会导致数值不稳定? 在什么情况下(例如稀疏矩阵、对称...)每种算法都更快/更精确? 回答1 这在很大程度上取决于对矩阵求幂的这些不同方式的实现细节。 一般而言,我认为特征分解 ( expm2 ) 不太适合稀疏矩阵,因为它可能会消除稀疏性。 应用于非对称矩阵也将更加困难,因为这将需要使用复杂的算术和更昂贵的算法来计算特征分解。 对于泰勒级数方法 ( expm3 ),如果存在固定数量的独立于矩阵范数的项,这听起来很危险。 当计算标量 x 的 e^x 时,泰勒级数中最大的项围绕着 n 接近 x 的项。 但是,这些(已弃用的
  • Numpy 矩阵求幂给出负值(Numpy matrix exponentiation gives negative value)
    问题 我想在斐波那契问题中使用NumPy ,因为它在矩阵乘法中的效率很高。 您知道有一种方法可以使用矩阵[[1, 1], [1, 0]]找到斐波那契数列。 我写了一些非常简单的代码,但是在增加n ,矩阵开始给出负数。 import numpy def fib(n): return (numpy.matrix("1 1; 1 0")**n).item(1) print fib(90) # Gives -1581614984 这可能是什么原因? 注意: linalg.matrix_power也给出负值。 注2:我尝试了从 0 到 100 的数字。它在 47 之后开始给出负值。这是一个大整数问题,因为 NumPy 是用 C 编码的吗? 如果是这样,我该如何解决这个问题? 编辑:将常规 python list矩阵与linalg.matrix_power一起linalg.matrix_power也会产生负面结果。 另外让我补充一点,并非所有结果在 47 之后都是负面的,它是随机发生的。 Edit2:我尝试使用@AlbertoGarcia-Raboso 建议的方法。 它解决了负数问题,但是出现了另一个问题。 它给出的答案为-5.168070885485832e+19 ,其中我需要-51680708854858323072L 。 所以我尝试使用int() ,它将它转换为L ,但现在由于精度损失
  • 有什么快速的矩阵求幂方法?(Is there any fast method of matrix exponentiation?)
    问题 有没有比简单的分而治之算法更快的矩阵求幂方法来计算M n (其中M为矩阵,n为整数)? 回答1 您可以将矩阵分解为特征值和特征向量。 然后你得到 M = V^-1 * D * V 其中V是特征向量矩阵,D是对角矩阵。 要将其提高到第N次方,您将获得类似以下内容的信息: M^n = (V^-1 * D * V) * (V^-1 * D * V) * ... * (V^-1 * D * V) = V^-1 * D^n * V 因为所有的V和V ^ -1项都抵消了。 由于D是对角线的,因此您只需要将一堆(实数)数提高到n次方即可,而不是完整矩阵。 您可以在n的对数时间内完成该操作。 计算特征值和特征向量为r ^ 3(其中r是M的行数/列数)。 取决于r和n的相对大小,这可能更快,也可能更快。 回答2 使用Euler快速幂算法非常简单。 使用下一个算法。 #define SIZE 10 //It's simple E matrix // 1 0 ... 0 // 0 1 ... 0 // .... // 0 0 ... 1 void one(long a[SIZE][SIZE]) { for (int i = 0; i < SIZE; i++) for (int j = 0; j < SIZE; j++) a[i][j] = (i == j); } //Multiply matrix
  • Python中的矩阵求幂(Matrix exponentiation in Python)
    问题 我试图在 Python 中对复杂矩阵求幂,但遇到了一些麻烦。 我正在使用scipy.linalg.expm函数,当我尝试以下代码时,出现了一条相当奇怪的错误消息: import numpy as np from scipy import linalg hamiltonian = np.mat('[1,0,0,0;0,-1,0,0;0,0,-1,0;0,0,0,1]') # This works t_list = np.linspace(0,1,10) unitary = [linalg.expm(-(1j)*t*hamiltonian) for t in t_list] # This doesn't t_list = np.linspace(0,10,100) unitary = [linalg.expm(-(1j)*t*hamiltonian) for t in t_list] 运行第二个实验时的错误是: This works! Traceback (most recent call last): File "matrix_exp.py", line 11, in <module> unitary_t = [linalg.expm(-1*t*(1j)*hamiltonian) for t in t_list] File "/usr/lib/python2.7/dist
  • Miller-Rabin算法判断是否为素数--python代码实现--利用了快速求幂取余的方法
    判断一个数P是否为素数的一般方法: 方法1:k从2开始,到n-1为止,判断P是否可以整除k 改进1:k到sqrt(P)为止 改进2:k从3开始且只考虑奇数 时间复杂度:O(P) 一个更快的算法:Miller-Rabin算法 实现思路:在2到P-1的范围内,随机选择一个数a,若a正好不满足,就输出P不是质数;否则重复m次,若每次的a都满足上述公式,就输出P是质数。 时间复杂度分析:m是一个不大的数,比如10;关键是计算耗时。 利用到的定理和思想: 1.费马定理:如果P是质数,一定有 因此,如果输出判断为不是质数,结果一定是对的。 反之则不一定成立,所以会出现不是质数却被判断为质数的情况。 2.定理:如果P不是质数,则没有一个a或者至少在[2,P-1]中一半数的a,使得不成立。 3.错误概率: 每次错误的概率不大于1/2,所以判定m次,错误概率不大于。所以m越大的话错误概率就越小,最后是一个可接受的值,甚至可以忽略。 属于Monte-Carlo算法,牺牲了一定的正确率但是换来了时间上的高效。 python代码实现 import random # 快速求幂取模 def quickPowMod(a,n,m): re = 1 base = a % m while(n>0): tem = n&1 if (tem): re = (re*base) % m base = (base*base) %
  • 大值 N 的矩阵求幂算法(Matrix Exponentiation Algorithm for large values of N)
    问题 我想计算非常大的 N 值的斐波那契,即。 10^6 复杂度为 O(logN)。 这是我的代码,但它在 30 秒内给出了 10^6 的结果,这非常耗时。帮我指出错误。我必须以模 10^9+7 给出输出。 static BigInteger mod=new BigInteger("1000000007"); BigInteger fibo(long n){ BigInteger F[][] = {{BigInteger.ONE,BigInteger.ONE},{BigInteger.ONE,BigInteger.ZERO}}; if(n == 0) return BigInteger.ZERO; power(F, n-1); return F[0][0].mod(mod); } void power(BigInteger F[][], long n) { if( n == 0 || n == 1) return; BigInteger M[][] = {{BigInteger.ONE,BigInteger.ONE},{BigInteger.ONE,BigInteger.ZERO}}; power(F, n/2); multiply(F, F); if( n%2 != 0 ) multiply(F, M); } void multiply(BigInteger F[][]
  • 在 Python 中加速矩阵向量乘法和求幂,可能通过调用 C/C++(Speeding up matrix-vector multiplication and exponentiation in Python, possibly by calling C/C++)
    问题 我目前工作的一个机器学习项目中-如果有一个数据矩阵Z和矢量rho -我不得不在计算物流损失函数的值和斜率rho 。 计算涉及基本的矩阵向量乘法和 log/exp 运算,并带有避免数值溢出的技巧(在前一篇文章中进行了描述)。 我目前正在使用 NumPy 在 Python 中执行此操作,如下所示(作为参考,此代码运行时间为 0.2 秒)。 虽然这很有效,但我想加快它的速度,因为我在我的代码中多次调用该函数(它代表了我项目中超过 90% 的计算)。 我正在寻找任何方法来改善此代码的运行时间而无需并行化(即只有 1 个 CPU)。 我很高兴在 Python 中使用任何公开可用的包,或者调用 C 或 C++(因为我听说这将运行时间提高了一个数量级)。 预处理数据矩阵Z也可以。 可能更好的计算被利用一些事情是该载体rho通常是稀疏的(具有50%左右的条目= 0),并且通常有多得多的行比列(在大多数情况下n_cols <= 100 ) import time import numpy as np np.__config__.show() #make sure BLAS/LAPACK is being used np.random.seed(seed = 0) #initialize data matrix X and label vector Y n_rows, n_cols = 1e6
  • 数论(2)逆元
    一、逆元理解 逆元是什么? 逆元就是扩大了概念的倒数。 定义:如果 ax≡1(mod M),就称x为在模M下 a的 逆元! 简单地说,如果一个数x满足 ax%M=1,那么x就称为在模M下 a的 逆元! 为什么说是扩大了概念的倒数呢,可见比起以前的倒数,只加了一个条件,即在后边加了一个 “%M",也可以这样理解,我们以前的倒数,也有这个条件,不过M是1。 二、逆元用处 有时候结果会让取模,除法只能用逆元取,如下: (a/b)%M=(ainv[b])%M=(a%M)(inv[b]%M)%M 三、方法 (一、费马小定理) 费马小定理: 如果p是一个质数,而整数a不是p的倍数,则有a^(p-1)≡1(mod p)。 [1] 局限性:只能用于p为素数,且gcd(a,p)为1 改造一下费马小定理式子得到:a*a^(p-2)≡1(mod p) 这样,便和逆元式子对应起来了,a的逆元便是x,即是a^(p-2) 用快速幂就OK了 int quick(int a,int n,int p) { if(n==0) return 1; if(n%2) return a*quick(a,n-1,p)%p; int tmp=quick(a,n/2,p); return tmp*tmp%p; } (二、扩展欧几里得定理) 局限性:gcd(a,M)=1,但并不要求M为素数,所以这个适用范围更广泛,而且稍快,推荐
  • MSVC 中微基准测试的优化障碍:告诉优化器你破坏内存?(Optimization barrier for microbenchmarks in MSVC: tell the optimizer you clobber memory?)
    问题 Chandler Carruth 在他的 CppCon2015 演讲中介绍了两个函数,可用于对优化器进行一些细粒度抑制。 它们对于编写优化器不会简单地使之变得毫无意义的微基准很有用。 void clobber() { asm volatile("" : : : "memory"); } void escape(void* p) { asm volatile("" : : "g"(p) : "memory"); } 这些使用内联汇编语句来改变优化器的假设。 clobber的汇编语句声明其中的汇编代码可以在内存中的任何位置读写。 实际的汇编代码是空的,但优化器不会查看它,因为它是asm volatile 。 当我们告诉它代码可能会在内存中的任何地方读写时,它就会相信它。 这有效地防止优化器在调用clobber之前重新排序或丢弃内存写入,并在调用clobber † 之后强制读取内存。 escape那个,另外使指针p对汇编块可见。 同样,因为优化器不会查看代码可以为空的实际内联汇编代码,并且优化器仍会假设该块使用指针p指向的地址。 这有效地强制任何p指向在内存中而不是在寄存器中,因为汇编块可能会从该地址执行读取。 (这很重要,因为clobber函数不会强制读取或写入编译器决定放入寄存器的任何内容,因为clobber中的汇编语句没有说明任何特定的内容必须对程序集可见。)
  • 卢卡斯定理证明
    在开始介绍卢卡斯定理之前我们先介绍下面的3个定理 1.乘法逆元 如果ax≡1 (mod p),且gcd(a,p)=1(a与p互质),则称a关于模p的乘法逆元为 x ax≡1 (mod p) 这个等式用中文描述就是 a乘一个数x并模p等于1,即 a%p*x%p=res,res%p=1 2.费马小定理:如果p是一个质数,而整数a不是p的倍数,即这两个数互质,则有 a^(p - 1) ≡ 1 mod(p) gcd(a, p) = 1 3.拓展欧几里得 已知整数a、b,扩展欧几里得算法可以在求得a、b的最大公约数的同时, 能找到整数x、y(其中一个很可能是负数),使它们满足贝祖等式ax + by = gcd(a, b) problem_1:为什么我们可以用费马小定理来求得逆元? 证明: 由费马小定理 a^(p - 1) ≡ 1 mod(p) 变形得 ->> a * a^(p - 2) ≡ 1 mod(p) 因为 a, p 互质,a * a^(p - 2) ≡ 1 mod(p) 且 ax≡1 (mod p) 则 x = a^(p - 2) mod(p) problem_2:为什么我们可以用拓展欧几里得来求得逆元? 证明:法1 我们都知道模后的数就是余数,比如12%5=12-5*2=2,18%4=18-4*4=2。(/是程序运算中的除) 那么ax≡1 (mod p)即ax-yp=1.把y写成
  • 计算并打印第n个质数(Calculating and printing the nth prime number)
    问题 我正在尝试计算素数,这已经完成了。 但是我只想计算和打印第n个质数(用户输入),而计算其余的(不会打印),只会打印第n个质数。 到目前为止,这是我写的内容: import java.util.Scanner; /** * Calculates the nth prime number * @author {Zyst} */ public class Prime { public static void main(String[] args) { Scanner input = new Scanner(System.in); int n, i = 2, x = 2; System.out.printf("This program calculates the nth Prime number\n"); System.out.printf("Please enter the nth prime number you want to find: "); n = input.nextInt(); for(i = 2, x = 2; i <= n; i++) { for(x = 2; x < i; x++) { if(i % x == 0) { break; } } if(x == i) { System.out.printf("\n%d is prime", x); } } } }
  • JavaScript分离轴定理(JavaScript Separating Axis Theorem)
    问题 我正在尝试使用JavaScript中的“分离轴定理”来检测两个正方形发生碰撞(一个旋转,一个不旋转)。 尽我所能,我无法弄清楚JavaScript中的样子,也找不到任何JavaScript示例。 请帮忙,用纯数字或JavaScript代码进行的解释将是最有用的。 更新:在研究了许多几何和数学理论之后,我决定在GitHub存储库中推出简化的SAT实现。 您可以在JavaScript中找到SAT的工作副本:https://github.com/ashblue/canvas-sat 回答1 变换多边形 首先,您必须通过应用angle的旋转来变换凸多边形(在这种情况下为正方形)的所有点,以便它们都在同一空间中。 为了将来支持缩放,平移等,我建议通过矩阵变换来实现。 您将必须编写自己的Matrix类或找到一些已经具有此功能的库(我相信有很多选择)。 然后,您将徒然使用以下代码: var transform = new Matrix(); transform.appendRotation(alpha); points = transform.transformPoints(points); 其中points是Point对象左右的数组。 碰撞算法概述 因此,这就是您遇到任何碰撞问题之前的全部内容。 关于碰撞算法,通常的做法是使用以下步骤尝试分离2个凸多边形(在您的情况下为正方形):