matlab解含参数高次方程,用matlab解一个含参数的一元三次方程，高手帮忙看看这个程序有什么问题，谢谢！ - 程序语言 - 小木虫 - 学术 科研 互动社区...

用matlab解一个含参数的一元三次方程：

C3*P(z)^3+C1*P(z)+D(z)-P_av=0

z是一线性数组，C3，C1常数，P_av是P在不同z的根的平均值。

思路是首先假设一个初始值P_av0, 然后解得P的一组根，取平均值和P_av0比较，如果收敛，则为所求的根，否则把新的平均值代回方程，再求根，直到收敛。

程序如下：

function P_real=P_distribution(a,b,e,q,d,N,m,z,P_r,Dsc)

clc

format long

% define parameters

a=-6.36e7;

b=1.38e9;

e=8.85e-11;

q=1.61e-19;

d=1e-7;

N=2e24;

% coeffecients

c3=b*e;

c1=a*e+1;

% start value

m=0;

P_av0=0.2;          % initial value for the P average

P_real=zeros(1,11);

%start the caculation

while 1

while 1

z=m*1e-8;

Dsc=-q*N*(d/2-z);

c0=Dsc-P_av0;

P_eq=[c3,0,c1,c0];

P_r=roots(P_eq);   %get the roots

for k=1:length(P_r)

P_rr=isreal(P_r(k));   %get the real roots

end

pp=P_r(P_rr);

P_real(1,m+1)=abs(pp);       %assign the real roots to the matrix

m=m+1;           %solve the equation again for another value

if m>10

break

end

end

P_av1=mean(P_real);    %get the average value

if abs(P_av1-P_av0)<1e-11   %less than the tolerance

break

else

P_av0=P_av1;

end

end

程序在matlab中无法运行，不知道问题出在哪里，高手帮忙看看，多谢了！