Matlab中的Newton-Raphson方法 [英] Newton-Raphson Method in Matlab
问题描述
我是Matlab的新手，我需要创建一个函数，该函数以x = a开始逼近，对Newton-Raphson方法进行n次迭代.此开始近似值不算作中间项，另一个要求是需要for循环.我看过其他类似的问题，但就我而言，我不想使用while循环.
I am new to matlab and I need to create a function that does n iterations of the Newton-Raphson method with starting approximation x = a. This starting approximation does not count as an interation and another requirement is that a for loop is required. I have looked at other similar questions posted but in my case I do not want to use a while loop.
这是我的输入应该是的:
This is what my inputs are supposed to be:
mynewton(f,a,n) which takes three inputs:
f: A function handle for a function of x.
a: A real number.
n: A positive integer.
这是到目前为止的代码.
And here is my code so far.
function r=mynewton(f,a,n)
syms x;
z=f(x);
y=a;
for i=1:n
y(i+1)=y(i)-(z(i)/diff(z(i)));
end
r=y
end
当我尝试调用该函数时，出现错误消息:
When I try to call the function I get an error saying:
Error in MuPAD command: DOUBLE cannot convert the input expression into a double array.
If the input expression contains a symbolic variable, use the VPA function instead.
Error in mynewton (line 6)
y(i+1)=y(i)-(z(i)/diff(z(i)));
问题是，我如何使用此VPA功能?授予我的其余代码可能也不是100％正确的，但是任何解决vpa问题或修复代码其他部分的帮助都将不胜感激.
The question is how do I use this VPA function? Granted the rest of my code probably isn't 100% correct either but any help that addresses the vpa problem or fixes the other parts of my code would be greatly appreciated.
谢谢！
推荐答案
使用牛顿-拉夫森(Newton-Raphson)技术，有两点不太正确……但肯定可以解决！解决此问题后，您无需再说到VPA错误.
There are two things that aren't quite correct with your Newton-Raphson technique... but certainly fixable! After we fix this, there is no need for the VPA error that you're speaking of.
第一个是迭代本身.回忆一下牛顿-拉夫森技术的定义:
The first one is the iteration itself. Recall the definition of the Newton-Raphson technique:
(来源: mit.edu )
(source: mit.edu)
对于下一次迭代，请使用上一次迭代的值.您正在做的是使用循环计数器并将其替换为f(x)，这是不正确的.它必须是上一次迭代的值.
For the next iteration, you use the previous iteration's value. What you're doing is using the loop counter and substituting this into your f(x), which is not correct. It must be the previous iteration's value.
如果您看一下如何编写函数，您将以符号方式定义函数，但是您尝试将数字值替换为函数.不幸的是，这在MATLAB中不起作用.如果您实际上要替换值，则必须使用 subs .这将为您替代x的函数或函数正在使用的任何自变量的实际值.完成此操作后，您的值仍然是sym类型.您需要将其强制转换为 double ，以便能够在数字上使用它.
If you take a look at how you're coding your function, you define your function symbolically, yet you are trying to substitute numeric values into your function. This unfortunately doesn't fly with MATLAB. If you actually want to substitute values in, you have to use subs. This will substitute an actual value for you as a function of x or whichever independent variable your function is using. Once you do this, your value is still a sym type. You need to cast this as double in order to be able to use this numerically.
也为了提高效率，不需要将y制成数组.只需将其设置为一个值即可在每次迭代时进行更新.说完所有这些后，您的代码将更新为如下所示.请注意，我在循环之前采用了函数的派生形式，以减少所需的计算量.我还将牛顿-拉夫逊迭代的分子和分母项分解开来，以使情况更清楚，并使之对于subs更可口.事不宜迟:
Also for efficiency, there is no need to make y an array. Just make this a single value that updates itself at each iteration. With all of this being said, your code is updated to look like this. Mind you, I took the derivative of the function before the loop to decrease the amount of computation you need to take. I also split up the numerator and denominator terms for the Newton-Raphson iterations to make things clear, and to make this more palatable for subs. Without further ado:
function r = mynewton(f,a,n)
syms x;
z = f(x);
diffZ = diff(z); %// Edit - Include derivative
y = a; %// Initial root
for idx = 1 : n
numZ = subs(z,x,y); %// Numerator - Substitute f(x) for f(y)
denZ = subs(diffZ,x,y); %// Denominator - Substitute for f'(x) for f'(y)
y = y - double(numZ)/double(denZ); %// Update - Cast to double to get the numerical value
end
r = y; %// Send to output
end
请注意，我在循环中将i替换为idx.原因是因为实际上不建议使用i或j作为循环索引，因为这些字母被保留来表示复数.如果您看 Shai 的这篇文章，您会发现使用这些变量作为循环索引实际上要慢一些:在Matlab中使用i和j作为变量
Take note that I replaced i with idx in the loop. The reason why is because it is actually not recommended to use i or j as loop indices as these letters are reserved to represent complex numbers. If you take a look at this post by Shai, you'll see that it's actually slower to use these variables as loop indices: Using i and j as variables in Matlab
在任何情况下，为了进行验证，假设我们的函数为y = sin(x)，而我的初始根目录为x0 = 2，经过5次迭代，我们可以这样做:
In any case, to test this out, suppose our function was y = sin(x), and my initial root was x0 = 2, with 5 iterations, we do:
f = @(x) sin(x);
r = mynewton(f, 2, 5)
r =
3.1416
这与我们对sin(x)的了解一致，因为sin(x)的截距位于pi的整数倍处. x0 = 2位于pi附近，因此可以正常工作.
This agrees with our knowledge of sin(x), as the intercepts of sin(x) are located at integer multiples of pi. x0 = 2 is located near pi so this does work as we expected.
您的原始代码在每次迭代中都将根的值存储在y中.如果确实要这样做，则必须修改代码，使其看起来像这样.请记住，我已经预先分配了y来提高效率:
Your original code had the values of the root stored at each iteration in y. If you really want to do that, you'll have to modify your code so that it looks something like this. Bear in mind that I pre-allocated y to make things more efficient:
function r = mynewton(f,a,n)
syms x;
z = f(x);
diffZ = diff(z);
y = zeros(1,n+1); %// Pre-allocate output array
y(1) = a; %// First entry is the initial root
for idx = 1 : n
numZ = subs(z,x,y(idx)); %// Remember to use PREVIOUS guess for next guess
denZ = subs(diffZ,x,y(idx));
y(idx+1) = y(idx) - double(numZ)/double(denZ); %// Place next guess in right spot
end
r = y; %// Send to output
end
通过使用与上面完全相同的参数运行此代码，我们得到:
By running this code using the exact same parameters as above, we get:
f = @(x) sin(x);
r = mynewton(f, 2, 5)
r =
2.0000 4.1850 2.4679 3.2662 3.1409 3.1416
r中的每个值告诉您在特定迭代中对根的猜测.数组的第一个元素是初始猜测(当然).下一个值是在牛顿-拉夫森根的每次迭代中的猜测.请注意，数组的最后一个元素是我们的最终迭代，大约等于pi.
Each value in r tells you the guess of the root at that particular iteration. The first element of the array is the initial guess (of course). The next values are the guesses at each iteration of the Newton-Raphson root. Take note that the last element of the array is our final iteration, which is roughly equal to pi.
这篇关于Matlab中的Newton-Raphson方法的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！
