去网上抄了个画图的函数,但是报错了,求各位大神解答
#sol_equ用来解六次方程
sol_equ<-function(coe1,coe2,coe3,coe4,coe5,coe6,coe7,b1,b2,e){
#定义了一个六次函数cubic
cubic<-function(x){
return(coe1*x^6+coe2*x^5+coe3*x^4+coe4*x^3+coe5*x^2+coe6*x+coe7)
}
k=0
while(2>1){
if(k>15){#设置最大循环次数是15
return(NA)
}
k=k+1
y1<-cubic(b1)
y2<-cubic(b2)
b0<-(b1+b2)/2
if(cubic(b0)==0){
return(b0)
}
if(y1==0){
return(b1)
}
if(y2==0){
return(b2)
}
y0<-cubic(b0)
if(y1*y0<0){
b2<-b0
}
if(y2*y0<0){
b1<-b0
}
if(abs(b1-b2)<e){
result<-((b1+b2)/2)
break
}
}
return(result)
}
#设置y,z 坐标值
y<-seq(-1.5,1.5,0.01)
z<-seq(-1.5,1.5,0.01)
x<-NULL
#下面是对于每个(y,z)(共301×301)求得x(x>0)
for(i in 1:301){
m<-NULL
for(j in 1:301){
a<-z[j]^2+9*y^2/4-1
b<–z[j]^3
c<–9*y^2*z[j]^3/80
r<-sol_equ(1,0,3*a,0,3*a^2+b,0,a^3+c,0,2,0.01)
m<-c(m,r)
}
x<-cbind(x,m)
}
x
library(rgl)
#当x>0时作图
persp3d(y,z,t(x),theta=60,phi=30,col=’red’)
#当x<0时作图
persp3d(y,z,-t(x),theta=60,phi=30,col=’red’,add=T)
报错如图
相关分类