الگوریتم مورد استفاده در تابع roots

[ehsan_faal]

سلام.
من باید ریشه های یک چند جمله ای رو توی زبان C++‎‎‎ بدست بیارم. میخوام که الگوریتمی که توی تابع roots متلب استفاده میشه رو استفاده کنم.
تا جایی که تحقیق کردم متلب مقادیر ویژه ماتریش Comapnion رو حساب میکنه. من هم همون الگوریتم رو برای ++C نوشتم منتها ریشه ها خیلی با جوابهای متلب فاصله داره، کسی میدونه ایراد کارم کجاست؟


vector<complex<double>> roots(
const std::vector<double> &coeffs)
{
int matsz = coeffs.size() - 1;
vector<complex<double>> vret;
vector<double> C(coeffs.size());
transform(coeffs.rbegin(),coeffs.rend(),C.begin(),[](auto& item){
return item;
});

MatrixXd companion_mat=MatrixXd::Zero(matsz,matsz);

for(int n = 0; n < matsz; n++){
for(int m = 0; m < matsz; m++){

if(n == m + 1)
companion_mat(n,m) = 1.0;

if(m == matsz-1)
companion_mat(n,m) =
-C[n] / C.back();

}
}

MatrixXcd eig = companion_mat.eigenvalues();

for(int i = 0; i < matsz; i++)
vret.push_back( eig(i) );

return vret;

}

جوابهای ++C:
Root Finding Section:
Interested Polynomial Coefficients:
4.16257 , 0 , -11.1309 , 0 , 10.9619 , 0 , -4.98668 , 0 , 1.28174 , 0
(0.939933,0.13167)
(0.939933,-0.13167)
(-0.939931,0.131673)
(-0.939931,-0.131673)
(-1.2409,0.678747)
(-1.2409,-0.678747)
(1.2409,0.678747)
(1.2409,-0.678747)

جوابهای متلب:
0.000000000000000 + 0.000000000000000i
-1.043428789214894 + 0.146171858075669i
-1.043428789214894 - 0.146171858075669i
-0.620284451427052 + 0.339282162499915i
-0.620284451427052 - 0.339282162499915i
1.043428789214895 + 0.146171858075669i
1.043428789214895 - 0.146171858075669i
0.620284451427052 + 0.339282162499915i
0.620284451427052 - 0.339282162499915i

(مطمئن نبودم که این تاپیک رو کجا باید ایجاد کنم ولی به نظرم چون میخوام راجع به الگوریتم استفاده شده توی متلب سوال کنم اینجا جای مناسبتریه)
با تشکر

[Desaghi]

بردار ضرائب را normalize کنید.

[rahnema1]

سلام
این سورس octave هست. به همین صورت اعمال کنید

function r = roots (v)

if (nargin != 1 || min (size (v)) > 1)
print_usage ();
elseif (any (isnan(v) | isinf(v)))
error ("roots: inputs must not contain Inf or NaN");
endif

n = numel (v);
v = v(:);

## If v = [ 0 ... 0 v(k+1) ... v(k+l) 0 ... 0 ], we can remove the
## leading k zeros and n - k - l roots of the polynomial are zero.

if (isempty (v))
f = v;
else
f = find (v ./ max (abs (v)));
endif
m = numel (f);

if (m > 0 && n > 1)
v = v(f(1):f(m));
l = max (size (v));
if (l > 1)
A = diag (ones (1, l-2), -1);
A(1,:) = -v(2:l) ./ v(1);
r = eig (A);
if (f(m) < n)
tmp = zeros (n - f(m), 1);
r = [r; tmp];
endif
else
r = zeros (n - f(m), 1);
endif
else
r = [];
endif

endfunction

[sim-power]

پیدا کردن حداقل و حداکثر از منحنی (http://sim-power.ir/%D9%BE%DB%8C%D8%AF%D8%A7-%DA%A9%D8%B1%D8%AF%D9%86-%D8%AD%D8%AF%D8%A7%D9%82%D9%84-%D9%88-%D8%AD%D8%AF%D8%A7%DA%A9%D8%AB%D8%B1-%D8%A7%D8%B2-%D9%85%D9%86%D8%AD%D9%86%DB%8C)
اگر ما حداقل و حداکثر را از یک نمودار جستجو کنیم، ما به ط.ر اساسی بیش ترین و کمترین نقاط روی نمودار از یک تابع در محل خاص یا محدوده خاصی از مقادیر از متغیر نمادین را جستجو می کنیم.
برای تابع y=f(x) نقاط روی نمودار که در نمودار شیب صفر دارند نقاط ایستا نامیده می شوند. به عبارت دیگر نقاط ایستا f’(x)=0 هستند.
برای پیدا کردن نقاط ایستا از تابع دیفرانسیل ما ، نیاز به مجموعه مشتق معادله صفر و حل معادلات داریم.
مثال
اجازه دهید نقاط ایستا از تابع f(x)=2x3+3x2-12x+17 را در مراحل زیر بگیریم:
1.اول اجازه دهید تابع را وارد کنیم و نمودار آن را رسم کنیم:
syms x
y =2*x^3+3*x^2-l2*x +17;% defining the function
ezplot(y)