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

Я разработал приведенный ниже код для поиска корней данного многочлена. Он отлично работает, но мне нужно адаптировать его, чтобы найти все корни, а не просто останавливаться, когда он сходится. Как мне это сделать? Я думал о создании внешнего цикла 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

fortran polynomials scientific-computing newtons-method
person Insight    schedule 31.03.2019    source источник
comment
Ньютон-Рафсон не сделает этого за вас, он поможет найти корень, вот и все. В общем случае вам нужно что-то вроде en.wikipedia.org/wiki/Bisection_method, или для конкретный случай многочленов, вы можете проверить соответствующие предложения, en.wikipedia.org/wiki/Descartes %27_rule_of_signs, en.wikipedia.org/wiki/Sturm%27s_theorem , en.wikipedia.org/wiki/Budan%27s_theorem и особенно их сочетание, en.wikipedia.org/wiki/Real-root_isolation   -  person tevemadar    schedule 31.03.2019
comment
Что касается самих методов, вы можете поискать идеи/помощь на math.stackexchange.com, пока речь не идет о Fortran или программирование, но математика.   -  person tevemadar    schedule 31.03.2019
comment
Я позволю себе не согласиться, я думаю, что можно придумать способ сохранить корень с учетом критерия сходимости, возможно, переменной ошибки, а затем повторно применить алгоритм на этот раз для другого начального x. Это будет делаться до тех пор, пока не будут найдены все корни, цикл остановится, когда количество корней станет степенью многочлена. Я просто не уверен, как это сделать на фортране.   -  person Insight    schedule 31.03.2019
comment
Да, именно это и решают предлагаемые алгоритмы. Просто они не называются Ньютон-Рафсон.   -  person tevemadar    schedule 31.03.2019
comment
Дело в том, что мне нужно использовать Ньютон-Рафсон. Вопрос в моем домашнем задании задан конкретно, чтобы использовать этот метод для поиска всех корней. Еще раз, я благодарю вас за ваши предложения, но я предпочитаю ждать чьей-либо поддержки.   -  person Insight    schedule 31.03.2019
comment
Вы можете комбинировать теги с +, поэтому вы можете попробовать выполнить поиск, например math.stackexchange.com/questions/tagged/newton. -raphson+roots или math.stackexchange.com/questions/tagged/polynomials+roots , а то и все три, другие теги и т.д.   -  person tevemadar    schedule 31.03.2019
comment
Если это конкретно корни этого многочлена, когда у вас есть корень, вы можете разделить его, а затем использовать метод Ньютона для многочлена более низкого порядка и повторять, пока не достигнете тривиального случая. Во всем этом есть много числовых тонкостей, которые я не буду притворяться, что полностью понимаю, но я предполагаю, что для вашей домашней работы этого будет достаточно. Кстати, также изучите метод Хорнера для оценки многочленов.   -  person Ian Bush    schedule 31.03.2019
comment
Вы должны искать из нескольких отправных точек. Сделайте цикл и попробуйте с несколькими начальными x. Или изменяйте функцию каждый раз, когда вы находите корень, как предлагает Ян Буш,   -  person Vladimir F    schedule 31.03.2019
comment
Как упоминалось выше, цикл с несколькими начальными значениями. Обратите внимание, что этот метод может дать сбой даже при наличии корней!. Добавьте счетчик = счетчик + 1 в этот цикл, который вы сбрасываете на ноль для каждого начального предположения. Откажитесь от предположения с помощью команды выхода примерно после 20 итераций. (Этот метод очень быстрый).   -  person Dan Sp.    schedule 31.03.2019
comment
Однако с этим подходом есть проблемы. X должен быть целым числом, чтобы цикл работал, и его нельзя переопределить внутри цикла, поэтому мне пришлось бы каким-то образом изменить выражение x-f/df.   -  person Insight    schedule 31.03.2019
comment
Вы можете вычислить x из i, t тривиально, вообще без проблем.   -  person Vladimir F    schedule 31.03.2019
comment
О, вы имеете в виду сделать что-то вроде x = x + i? Я пробовал это, но никакой оценки сходимости не было бы, алгоритм, по сути, просто анализирует значение df для многих точек и надеется найти корень, вся суть метода Ньютона теряется. Если это не то, что вы имели в виду, пожалуйста, уточните.   -  person Insight    schedule 01.04.2019
comment
@Insight Это или просто x = i*что-то. Но не анализировать df по большому количеству точек, а начинать метод Ньютона-Рафсона с гораздо меньшего количества точек! Как сказано в ответе ниже, метод Ньютона-Рафсона является локальным методом. С разными начальными точками вы получите разный корень. Область сходимости к любому заданному корню может быть очень непредсказуемой en.wikipedia.org/wiki/Newton_fractal   -  person Vladimir F    schedule 01.04.2019

Ответы (1)


arrow_upward
3
arrow_downward

Возможный алгоритм нахождения всех корней многочлена P состоит в следующем:

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

Есть некоторые тонкости:

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

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

person Community    schedule 01.04.2019