Chức năng biểu thức chuyển đổi kiểu R()()

Question

Tôi đã cố gắng viết một chương trình trong R để thực hiện phương pháp Newton. Tôi đã hầu như thành công, nhưng có hai tiếng rên rỉ nhỏ đã làm phiền tôi. Dưới đây là mã của tôi:

Newton<-function(f,f.,guess){ 
    #f <- readline(prompt="Function? ") 
    #f. <- readline(prompt="Derivative? ") 
    #guess <- as.numeric(readline(prompt="Guess? ")) 
    a <- rep(NA, length=1000) 
    a[1] <- guess 
    a[2] <- a[1] - f(a[1])/f.(a[1]) 
    for(i in 2:length(a)){ 
     if(a[i] == a[i-1]){ 
      break 
     } 
     else{ 
      a[i+1] <- a[i] - f(a[i])/f.(a[i]) 
     } 
    } 
    a <- a[complete.cases(a)] 
    return(a) 
} 

1. tôi không thể nhận R để nhận ra các chức năng f và f. nếu tôi thử sử dụng readline() để nhắc nhở cho người dùng nhập vào. Tôi nhận được lỗi "Lỗi trong Newton(): không thể tìm thấy chức năng" f. "" Tuy nhiên, nếu tôi nhận xét ra các readlines (như trên), xác định f và f. trước, sau đó tất cả mọi thứ hoạt động tốt.

2. Tôi đã cố gắng làm cho R tính toán đạo hàm của hàm. Vấn đề là đối tượng lớp mà R có thể lấy các dẫn xuất biểu tượng là expression(), nhưng tôi muốn lấy đạo hàm của một số function() và cho nó một số function(). Tóm lại, tôi gặp sự cố với việc chuyển đổi loại giữa expression() và function().

Tôi có một giải pháp xấu xí nhưng hiệu quả cho đi function()-expression(). Với hàm f, D(body(f)[[2]],"x") sẽ cung cấp cho đạo hàm của f. Tuy nhiên, đầu ra này là expression() và tôi không thể chuyển nó trở lại thành function(). Tôi có cần sử dụng eval() hoặc gì đó không? Tôi đã cố gắng hết hạn, nhưng không có kết quả. Ví dụ:

g <- expression(sin(x)) 
g[[1]] 
sin(x) 
f <- function(x){g[[1]]} 
f(0) 
sin(x) 

khi những gì tôi muốn là f (0) = 0 vì tội lỗi (0) = 0.

EDIT: Cảm ơn tất cả! Đây là mã mới của tôi:

Newton<-function(f,f.,guess){ 
    g<-readline(prompt="Function? ") 
    g<-parse(text=g) 
    g.<-D(g,"x") 
    f<-function(x){eval(g[[1]])} 
    f.<-function(x){eval(g.)} 
    guess<-as.numeric(readline(prompt="Guess? ")) 
    a<-rep(NA, length=1000) 
    a[1]<-guess 
    a[2]<-a[1]-f(a[1])/f.(a[1]) 
    for(i in 2:length(a)){ 
     if(a[i]==a[i-1]){break 
     }else{ 
     a[i+1]<-a[i]-f(a[i])/f.(a[i]) 
     } 
    } 
a<-a[complete.cases(a)] 
#a<-a[1:(length(a)-1)] 
return(a) 
} 

Answer 1

1. vấn đề đầu tiên này phát sinh vì readline đọc trong một chuỗi văn bản, trong khi những gì bạn cần là một biểu hiện. Bạn có thể sử dụng parse() để chuyển đổi chuỗi văn bản để một biểu thức:

   f <-readline(prompt="Function? ") 
   sin(x) 
   f 
   # [1] "sin(x)" 
   
   f <- parse(text = f) 
   f 
   # expression(sin(x)) 
   
   g <- D(f, "x") 
   g 
   # cos(x) 
   

2. Để vượt qua trong các giá trị cho các đối số trong lời gọi hàm trong biểu thức (! Whew), bạn có thể eval() nó trong một môi trường chứa cung cấp giá trị. Độc đáo, R sẽ cho phép bạn để cung cấp những giá trị trong một danh sách cung cấp cho envir= đối số của eval():

   > eval(f, envir=list(x=0)) 
   # [1] 0 
   

Answer 2

Josh đã trả lời câu hỏi của bạn

Đối với phần 2 bạn có thể đã sử dụng

g <- expression(sin(x)) 

g[[1]] 
# sin(x) 

f <- function(x){ eval(g[[1]]) } 

f(0) 
# [1] 0 
f(pi/6) 
# [1] 0.5 

Answer 3

BTW, gần đây đã viết một món đồ chơi tính toán mô hình fractal dựa trên sự hội tụ gốc của phương pháp Newton trong mặt phẳng phức tạp, tôi có thể khuyên bạn nên quăng vào tôi mã như sau (nơi danh sách đối số của hàm chính bao gồm "func" và "varname").

func<- gsub(varname, 'zvar', func) 
    funcderiv<- try(D(parse(text=func), 'zvar')) 
    if(class(funcderiv) == 'try-error') stop("Can't calculate derivative") 

Nếu bạn thận trọng hơn, bạn có thể bao gồm một cuộc tranh cãi "funcderiv", và quấn mã của tôi trong

if(missing(funcderiv)){blah blah} 

Ahh, tại sao không: đây là chức năng hoàn chỉnh của tôi cho mọi người sử dụng và thưởng thức :-)

# build Newton-Raphson fractal 
#define: f(z) the convergence per Newton's method is 
# zn+1 = zn - f(zn)/f'(zn) 
#record which root each starting z0 converges to, 
# and to get even nicer coloring, record the number of iterations to get there. 
# Inputs: 
# func: character string, including the variable. E.g., 'x+ 2*x^2' or 'sin(x)' 
# varname: character string indicating the variable name 
# zreal: vector(preferably) of Re(z) 
# zim: vector of Im(z) 
# rootprec: convergence precision for the NewtonRaphson algorithm 
# maxiter: safety switch, maximum iterations, after which throw an error 
# 
nrfrac<-function(func='z^5 - 1 ', varname = 'z', zreal= seq(-5,5,by=.1), zim, rootprec=1.0e-5, maxiter=1e4, drawplot=T, drawiterplot=F, ...) { 
    zreal<-as.vector(zreal) 
    if(missing(zim)) zim <- as.vector(zreal) 
# precalculate F/F' 
    # check for differentiability (in R's capability) 
    # and make sure to get the correct variable name into the function 
    func<- gsub(varname, 'zvar', func) 
    funcderiv<- try(D(parse(text=func), 'zvar')) 
    if(class(funcderiv) == 'try-error') stop("Can't calculate derivative") 
# Interesting "feature" of deparse : default is to limit each string to 60 or64 
# chars. Need to avoid that here. Doubt I'd ever see a derivative w/ more 
# than 500 chars, the max allowed by deparse. To do it right, 
# need sum(nchar(funcderiv)) as width, and even then need to do some sort of 
# paste(deparse(...),collapse='') to get a single string 
    nrfunc <- paste(text='(',func,')/(',deparse(funcderiv, width=500),')', collapse='') 
# first arg to outer() will give rows 
# Stupid Bug: I need to REVERSE zim to get proper axis orientation 
    zstart<- outer(rev(zim*1i), zreal, "+") 
    zindex <- 1:(length(zreal)*length(zim)) 
    zvec <- data.frame(zdata=as.vector(zstart), zindex=zindex,  itermap=rep(0,length(zindex)), badroot=rep(0,length(zindex)), rooterr=rep(0,length(zindex))) 

#initialize data.frame for zout. 
    zout=data.frame(zdata=rep(NA,length(zstart)), zindex=rep(NA,length(zindex)),  itermap=rep(0,length(zindex)), badroot=rep(0,length(zindex)), rooterr=rep(0,length(zindex))) 
    # a value for rounding purposes later on; yes it works for rootprec >1 
    logprec <- -floor(log10(rootprec)) 
    newtparam <- function(zvar) {} 
    body(newtparam)[2] <- parse(text=paste('newz<-', nrfunc, collapse='')) 
    body(newtparam)[3] <- parse(text=paste('return(invisible(newz))')) 
    iter <- 1 
    zold <- zvec # save zvec so I can return original values 
    zoutind <- 1 #initialize location to write solved values 
    while (iter <= maxiter & length(zold$zdata)>0) { 
     zold$rooterr <- newtparam(zold$zdata) 
     zold$zdata <- zold$zdata - zold$rooterr 
     rooterr <- abs(zold$rooterr) 
     zold$badroot[!is.finite(rooterr)] <- 1 
     zold$zdata[!is.finite(rooterr)] <- NA 
# what if solvind = FFFFFFF? -- can't write 'nothing' to zout 
     solvind <- (zold$badroot >0 | rooterr<rootprec) 
      if(sum(solvind)>0) zout[zoutind:(zoutind-1+sum(solvind)),] <- zold[solvind,] 
    #update zout index to next 'empty' row 
     zoutind<-zoutind + sum(solvind) 
# update the iter count for remaining elements: 
     zold$itermap <- iter 
# and reduce the size of the matrix being fed back to loop 
     zold<-zold[!solvind,] 
     iter <- iter +1 
    # just wonder if a gc() call here would make any difference 
# wow -- it sure does 
     gc() 
    } # end of while 
# Now, there may be some nonconverged values, so: 
# badroot[] is set to 2 to distinguish from Inf/NaN locations 
     if(zoutind < length(zindex)) { # there are nonconverged values 
# fill the remaining rows, i.e. zout.index:length(zindex) 
      zout[(zoutind:length(zindex)),] <- zold # all of it 
      zold$badroot[] <- 2 # yes this is safe for length(badroot)==0 
      zold$zdata[]<-NA #keeps nonconverged values from messing up results 
      } 
# be sure to properly re-order everything... 
    zout<-zout[order(zout$zindex),] 
    zout$zdata <- complex(re=round(Re(zout$zdata),logprec), im=round(Im(zout$zdata),logprec)) 
    rootvec <- factor(as.vector(zout$zdata), labels=c(1:length(unique(na.omit(as.vector(zout$zdata)))))) 
    #convert from character, too! 
    rootIDmap<-matrix(as.numeric(rootvec), nr=length(zim)) 
# to colorize very simply: 
    if(drawplot) { 
      colorvec<-rainbow(length(unique(as.vector(rootIDmap)))) 
     imagemat<-rootIDmap 
     imagemat[,]<-colorvec[imagemat] #now has color strings 
     dev.new() 
# all '...' arguments used to set up plot 
     plot(range((zreal)),range((zim)), t='n',xlab='real',ylab='imaginary',...) 
     rasterImage(imagemat, range(zreal)[1], range(zim)[1], range(zreal)[2], range(zim)[2], interp=F)  
     } 

    outs <- list(rootIDmap=rootIDmap, zvec=zvec, zout=zout, nrfunc=nrfunc) 
    return(invisible(outs)) 
}