HDU 1402 大数乘法 FFT
发布时间:2020-12-14 03:38:41 所属栏目:大数据 来源:网络整理
导读:大整数乘法,给定两个长度不超过10000的整数,返回乘法的结果。 char*multi(char* number_a,char* number_b) 有疑问欢迎提问,本人学通信的,,手上有 数字信号处理 书,,可以把图搬出来解答的 #include "stdafx.h"#include iostream#include cstring#inclu
大整数乘法,给定两个长度不超过10000的整数,返回乘法的结果。 char*multi(char* number_a,char* number_b) 有疑问欢迎提问,本人学通信的,,手上有 数字信号处理 书,,可以把图搬出来解答的 #include "stdafx.h" #include <iostream> #include <cstring> #include <cmath> #define pi acos(-1.0) using namespace std; struct complex { double r,i; complex(double real = 0,double image = 0) { r = real; i = image; } complex operator +(const complex &a) { return complex(r + a.r,i + a.i); } complex operator -(const complex &a) { return complex(r - a.r,i - a.i); } complex operator *(const complex &a) { return complex(r*a.r - i*a.i,r*a.i + i*a.r); } }; void transAddr(complex *x,int len) {// 变址 int i,j,k; complex temp; for (i = 1,j = 0; i < len - 1; i++) {// 0和N-1码位反转结果不变无需交换,略去对其操作。 // 每次循环把j的开头遇到的1全变为0,随后遇到的第一个0变为1,就能得i的二进制码位反转结果。例如1100-->0010 k = len >> 1; while (j & k) {// 从最高位依次往低位检测,若为1,则通过相与变0 j &= ~k; k >>= 1; } j |= k;// 第一次遇到的0变为1 // 交换互为下标码位反转的元素 if (i < j) { temp = x[i]; x[i] = x[j]; x[j] = temp; } } } void fft(complex *x,int len,int factor) {// factor == 1时为DFT,factor==-1为IDFT int i,k; complex g,h; // 蝶形运算中的两个零时存储量G(k)、H(k) transAddr(x,len); // 先变址 for (i = 2; i <= len; i <<= 1)// 记多级蝶形运算的每列中各组序列长度为i,从长度为2的序列开始,后一列序列长度为前一列2倍,直到i为原始长度len为止 { complex Wn(cos(-factor * 2 * pi / i),sin(-factor * 2 * pi / i));// 为旋转增量Wn == e^(-2pi/N)按照欧拉公式展开 for (j = 0; j<len; j += i) {// 某列中各组起始下标 complex w(1,0); // 初始化蝴蝶因子。注:Wn^0 == 1 for (k = j; k < j + i / 2; k++) {// 蝶形运算 g = x[k]; h = w * x[k + i / 2]; x[k] = g + h; x[k + i / 2] = g - h; w = w * Wn; // 旋转蝴蝶因子 } } } if (factor == -1) // IDFT跟DFT区别为定义式中蝴蝶因子指数为正 Wn^(-1) == e^(2pi/N),且大小要乘1/N。注:此处虚部不会再用到,故虚部运算略去。 for (i = 0; i < len; i++) x[i].r /= len; } char* multi(char* number_a,char* number_b) { // 检查空指针 if (!number_a || !number_b) { cout << "error:输入了空指针n"; system("pause"); return NULL; } int len_a = strlen(number_a); int len_b = strlen(number_b); // 检查字符串长度不为空 if (!len_a || !len_b) { cout << "error:输入了空字符串n"; system("pause"); return NULL; } // 乘法结果的长度最多为2N,且FFT要求len为2的幂 int len_final = 1; while (len_final < 2 * len_a || len_final < 2 * len_b) len_final <<= 1; complex *cmplx_a = (complex*)malloc(len_final * sizeof(complex)); complex *cmplx_b = (complex*)malloc(len_final * sizeof(complex)); int *product = (int*)malloc(len_final * sizeof(int)); int i,temp; for (i = 0; i<len_a; i++) { temp = number_a[len_a - 1 - i] - '0'; if (0 <= temp || temp <= 9) {// 检查输入是否全为数字 cmplx_a[i] = complex(temp,0); } else { free(cmplx_a); free(cmplx_b); free(product); cout << "error:输入字符串必须全为数字n"; system("pause"); return NULL; } } for (; i < len_final; i++) { cmplx_a[i] = complex(0,0); } for (i = 0; i<len_b; i++) { temp = number_b[len_b - 1 - i] - '0'; if (0 <= temp || temp <= 9) {// 检查输入是否全为数字 cmplx_b[i] = complex(temp,0); } else { free(cmplx_a); free(cmplx_b); free(product); cout << "error:输入字符串必须全为数字n"; system("pause"); return NULL; } } for (; i < len_final; i++) { cmplx_b[i] = complex(0,0); } fft(cmplx_a,len_final,1);//通过DFT化卷积为相乘 fft(cmplx_b,1); for (i = 0; i < len_final; i++) cmplx_a[i] = cmplx_a[i] * cmplx_b[i]; fft(cmplx_a,-1);//逆变换得结果 // 计算中产生些许误差,故需要四舍五入 for (i = 0; i < len_final; i++) { product[i] = floor(cmplx_a[i].r + 0.5); } // 每个位上的书不超过9,超出则进位 for (i = 0; i < len_final - 1; i++) { product[i + 1] += product[i] / 10; product[i] %= 10; } // 找到最高位在i处 for (i = len_final - 1; i >= 0; i--) { if (product[i] != 0) break; } // 有i+1位数,再加上 ,共需i+2个单元 char *number_c = (char*)malloc((i + 2) * sizeof(char)); char * chrptr = number_c; for (; i >= 0; i--) { *chrptr++ = product[i] + '0'; } *chrptr = ' '; free(cmplx_a); free(cmplx_b); free(product); return number_c; } int _tmain(int argc,_TCHAR* argv[]) { char *number_a = "9999"; char *number_b = "9999"; char *number_c = multi(number_a,number_b); cout << number_c <<endl; system("pause"); return 0; } (编辑:李大同) 【声明】本站内容均来自网络,其相关言论仅代表作者个人观点,不代表本站立场。若无意侵犯到您的权利,请及时与联系站长删除相关内容! |