PS: 希望你这个不是作业，作业要自己做XD primefactor <- function(num){ p <- which(num %% 1:num == 0) isPrime <- function(num){ return(sum(num %% 1:num == 0) == 2) } p <- p[sapply(p, isPrime)] e <- vector('numeric', length(p)) for (i in 1:length(p)) { repeat{ if (num %% p[i]^e[i] != 0) break e[i] = e[i] + 1 } } return(list('pfactor'=p, 'multi'=e - 1)) } st = proc.time() primefactor(12345678) # $pfactor # [1] 2 3 47 14593 # # $multi # [1] 1 2 1 1 proc.time() - st # user system elapsed # 20.53 5.54 26.29 env: [email protected] machine with RRO-3.2.2 in windows 7 SP1 在p <- which(num %% 1:num == 0)跟isPrime这两步时要小心记忆体会被用完 如果担心的话可以拆分成repeat去做 例如： num = 12345678 p = vector('numeric', 1000) k = 1 step_p = 100000 repeat{ tmp = which(num %% k:(k+step_p) == 0) if (length(tmp) > 0) p[(sum(p!=0)+1):(sum(p!=0)+length(tmp))] = tmp + k - 1 if (k > num) break k = k + step_p } p = p[p!=0] # 1 2 3 6 9 18 47 # 94 141 282 423 846 14593 29186 # 43779 87558 131337 262674 685871 1371742 2057613 # 4115226 6172839 12345678 给一个更快的方法(直接用向量化去做)： primefactorFast <- function(num){ p <- which(num %% 1:num == 0) isPrime <- function(num){ return(sum(num %% 1:num == 0) == 2) } p <- p[sapply(p, isPrime)] e <- vector('numeric', length(p)) e_basic <- 1 step_e <- 3 subset_e <- 1:length(p) repeat{ m = num %% sweep(matrix(replicate(step_e, p[subset_e]), length(p[subset_e])), 2, e_basic:(e_basic + step_e - 1), '^') == 0 e[subset_e] = e[subset_e] + rowSums(m) e_basic = e_basic + step_e subset_e = which(apply(m,1,all)) if (is.null(subset_e) || length(subset_e) == 0) break } return(list('pfactor'=p, 'multi'=e)) } st = proc.time() primefactorFast(12345678) # $pfactor # [1] 2 3 47 14593 # # $multi # [1] 1 2 1 1 proc.time() - st # user system elapsed # 0.06 0.00 0.06 PS: 更动step_e可以得到更好的效能 final version: http://pastebin.com/zkD3N2kp test 2^11*3^6*47 (10^7.8) 只需要6.5秒 ※ 引述《haolihy (好厉害)》之铭言： : [问题类型]: : 效能谘询(我想让R 跑更快) : [软体熟悉度]: : 入门(写过其他程式，只是对语法不熟悉) : [问题叙述]: : 质因数分解 : [程式范例]: : 想写一个程式对正整数n作质因数分解 : 目标是n<10^8 可以1min跑完 : 我的程式大概要6~7min 改进空间应该很大 : 而太慢的主因应该是列出<n质数表的程式太慢 : 请大家赐教 : 质数表的程式: : primeY <- function(num){ : sqnum <- floor(sqrt(num)) : x <- c(2:num) : y <- c() : repeat{ : if(x[1] > sqnum) break : y <- c(y, x[1]) : x <- x[x %% x[1] != 0] : } : return(c(y,x)) : } : 质因数分解的程式: : primefactor <- function(num){ : p <- primeY(num)[which(num %% primeY(num) == 0)] : q <- c() : e <- c() : for(i in 1:length(p)){ : e[i] <- 1 : q[i] <- num / p[i] : repeat{ : if(q[i] %% p[i] != 0) break : e[i] <- e[i] + 1 : q[i] <- q[i] / p[i] : } : } : return(list('pfactor'=p, : 'multi'=e)) : } : [环境叙述]: : R 3.2.2 --

※ 发信站: 批踢踢实业坊(ptt.cc), 来自: 140.109.73.190 ※ 文章网址: https://webptt.com/cn.aspx?n=bbs/R_Language/M.1446631193.A.0DC.html ※ 编辑: celestialgod (140.109.73.190), 11/04/2015 18:43:36 ※ 编辑: celestialgod (180.218.154.163), 11/04/2015 19:51:09