如何使用 fsolve 函数在 R 中找到多个非线性方程的解

Question

我想用“fsolve”函数求解非线性方程，方程和代码如下，但我只能用“fsolve”函数求一组非线性方程组的解，例如，我在A和B系数中有三个数（A_coeff和B_coeff），按照我的想法是每个数经过公式计算出一组解，那么三个，应该是三组解决方案，我该怎么做才能实现它们

A_coeff<-c(177506.9,177639.3,178039.4)
B_coeff<-c(0.0003485474,0.0005155126,0.0004671370)
C_coeff<-5.511464
D_coeff<-23.39138
E_coeff<-5.0866e+17
F_coeff<-0.9732414
library('pracma')
Para_fun <- function(temp1) {
  new <- sqrt((4*temp1-1)/3)
  return(new)
}
Para_fun2<- function(temp1) {
  new2 <- ceiling(temp1/C_coeff)
  return(new2)
}
F_try<- function(x){
  s_actual <- x[1]
  K_actual <- x[2]
  n_tube <- x[3]
  c( A_coeff/K_actual-s_actual,
     (B_coeff+F_coeff/(E_coeff/Para_fun(n_tube)^(2/3))^0.25)^-1-K_actual,
     Para_fun2(s_actual)*D_coeff-n_tube)}
x0_xinitial_value<- c(20,2000,20)
X_result<- fsolve(F_try, x0_xinitial_value)
X_result$x

Answer 1

解决问题的最简单方法是用循环求解每对 A_coeff 和 B_coeff 的方程组。

将函数 F_try 重新定义为（我重写了代码以使其更易于阅读并减少混淆）

F_try<- function(x,k){
  s_actual <- x[1]
  K_actual <- x[2]
  n_tube <- x[3]
  y <- numeric(length(x))
  y[1] <- A_coeff[k]/K_actual-s_actual
  y[2] <- (B_coeff[k]+F_coeff/(E_coeff/Para_fun(n_tube)^(2/3))^0.25)^-1-K_actual
  y[3] <- Para_fun2(s_actual)*D_coeff-n_tube

y }

参数k是系数A_coeff和B_coeff向量的索引。

如果你这样尝试

X_result <- matrix(NA,nrow=3,ncol=3)
xstart <- x0_xinitial_value
for( k in 1:3){
    z <- fsolve(F_try, xstart,k=k)
    X_result[k,] <- z$x
}
X_result

您将收到一条错误消息

Error in if (norm(s, "F") < tol || norm(as.matrix(ynew), "F") < tol) break :

有消息

  missing value where TRUE/FALSE needed
Calls: fsolve -> broyden
In addition: Warning message:
In sqrt((4 * temp1 - 1)/3) : NaNs produced
Execution halted

目前还不清楚问题出在哪里以及为什么会出现错误。

还有另一个包 nleqslv 可以更深入地了解问题所在。

你可以这样使用它

library(nleqslv)
X_result <- matrix(NA,nrow=3,ncol=3)
xstart <- x0_xinitial_value
for( k in 1:3){
    z <- nleqslv(xstart,F_try,k=k)
    X_result[k,] <- z$x
}
X_result

检查X_result表明第三种解决方案很可能是错误的。长话短说，对于 k=3 和您提供的起始值，算法似乎无法找到解决方案。

一个解决方案是使每个 k 的起始值等于前一个 k 的解决方案。像这样

X_result <- matrix(NA,nrow=3,ncol=3)
xstart <- x0_xinitial_value
for( k in 1:3){
    z <- nleqslv(xstart,F_try,k=k)
    X_result[k,] <- z$x
    xstart <- z$x
}
X_result

导致

          [,1]     [,2]     [,3]
[1,]  72.60480 2444.837 327.4793
[2,] 102.59563 1731.451 444.4362
[3,]  94.16426 1890.732 421.0448

建议为该矩阵的每一行检查 nleqslv 的退出代码以确保找到解决方案。

如何使用 fsolve 函数在 R 中找到多个非线性方程的解

How to find the solutions of multiple nonlinear equations in R using the fsolve function

equation

r