Метод Ньютона-Рафсона с использованием Фортрана 90 - PullRequest
Новая
       74

Метод Ньютона-Рафсона с использованием Фортрана 90

2 голосов
/ 31 марта 2019

Я разработал код ниже для нахождения корней данного полинома. Он работает нормально, но мне нужно адаптировать его, чтобы найти все корни, а не просто остановиться, когда он сходится Как мне это сделать? Я думал о создании внешнего цикла do для значений x, но я не уверен, что это правильный подход. Любая помощь приветствуется, спасибо заранее. 

PROGRAM nr
integer :: i
real :: x, f, df
write(*,*) "x=?"
read (*,*)  x           
write (*,*) '# Initial value: x=',x

 do i=1,100
     f= x**4 - 26*(x**3) + 131*(x**2) - 226*x + 120
     df = 4*(x**3) - 3.0*26*(x**2) + 2.0*131*x - 226
     write (*,*) i,x,f,df
     x = x-f/df
 end do

 write (*,*) '#x = ',x




END PROGRAM

1 Ответ

3 голосов
/ 01 апреля 2019

Возможный алгоритм поиска всех корней многочлена P состоит в:

  • Начните с некоторого X0 и найдите корень R, используя алгоритм Ньютона.
  • Деление P на (X-R): деление является точным (с точностью до числовой ошибки), поскольку R является корнем. (этот шаг называется дефляция )
  • Перезапуск с начала, если частное имеет степень> 1.

Есть несколько тонкостей:

  • Если ваш многочлен имеет действительные коэффициенты, а X0 является действительным, вы найдете реальный корень, только если он есть. Чтобы найти сложные корни, вам нужно начать с комплексного X0 и, конечно, использовать сложную арифметику.
  • Не все значения X0 гарантируют сходимость, поскольку Ньютон-Рафсон является локальным методом. См. теорему Ньютона-Канторовича и бассейнов притяжения метода Ньютона.
  • Если есть несколько корней, конвергенция там намного медленнее. Есть способы приспособить метод Ньютона, чтобы справиться с этим.
  • Численные ошибки, вносимые на этапе дефляции, складываются и обычно приводят к низкой точности последующих корней. Это особенно проблема для многочленов высокой степени (насколько велика на самом деле зависит). В крайних случаях вычисленные корни могут быть довольно далеки от точных корней. Кроме того, некоторые полиномы являются более «сложными», чем другие: см., Например, полином Уилкинсона .
  • Существуют способы поиска информации о корнях (связанные с абсолютными значениями, окружности, окружающие корни ...), это может помочь в выборе хорошего X0 в алгоритме: см. Геометрические свойства полиномиальных корней . Если вы ищете только реальные корни из реальных многочленов, вы сначала находите оценки, а затем теорема Штурма , чтобы найти интервалы, охватывающие корни. Другой подход в реальном случае состоит в том, чтобы сначала найти корни производной P '(и найти ее корни, сначала найти корни P' 'и т. Д.), Так как корни P' отделяют корни P (за исключением несколько корней).
  • Другое предлагаемое чтение: Корни полиномов . Обратите внимание, что есть намного лучшие алгоритмы, чем Ньютон + дефляция. Существуют также алгоритмы, предназначенные для поиска всех корней.

Однако, если вас интересует только этот конкретный полином, я предлагаю сначала взглянуть на проблему, так как это гораздо проще, чем в общем случае: WolframAlpha . Здесь, начиная с целочисленных значений X0, будет работать достаточно хорошо ...

...