итерация с фиксированной точкой - PullRequest
1 голос
/ 07 октября 2011

я пытаюсь найти такое p, что для данной функции f (p) мы имеем равенство

p = f (p);Вот код

#include <iostream>
#include<math.h>
using namespace std;
float fixed(float x){
    return (float)(pow(x,3)-4*pow(x,2)-10);
}

int main(){
    float p=0.0;
    float p0=1.5;
    float tol=(float).001;
    int N=25;
    int i=1;
    while(i<N){
        p=(float)fixed(p0);
        if((p-p0)<tol){
          cout<<p<<endl;
          break;
        }
        i=i+1;
        p0=p;
        if(i>N){
          cout<<"solution not found ";
          break;
        }
    }
    return 0;
}

Я пробовал другую начальную точку, также разные допуски, но результат очень глупо -16 или -15. Что-то не так? Код правильный? Пожалуйста, помогите

Ответы [ 3 ]

3 голосов
/ 07 октября 2011

Я думаю, что у вас просто нет ситуации, в которой применим итерационный алгоритм. См. Здесь для некоторых условий.Ваша функция не имеет привлекательной фиксированной точки около 1,5, и алгоритм расходится.

Но почему цифры: ваша функция f(x) = x^3 - 4x - 10, поэтому решение f(x) = x равносильно нахождению нулей f(x) - xи есть только один реальный ноль около 5,35.Как бы то ни было, f'(x) в этот момент очень велик, поэтому даже там итерационный алгоритм не может быть использован.

Алгоритм числового поиска корня может быть более подходящим подходом.

3 голосов
/ 07 октября 2011

Я не уверен, чего вы пытаетесь достичь, но, вероятно, используете:

if(fabs(p-p0)<tol){

вместо:

if((p-p0)<tol){

приблизит вас к тому месту, куда вы хотите пойти.

Причина в том, что все отрицательные значения меньше, чем .001, поэтому, если результатом функции является отрицательное значение (как в этом случае), вы немедленно остановитесь.

Кстати: этот чек:

if(i>N){

никогда не будет правдой. Вы, вероятно, хотели использовать == или >= вместо >.

1 голос
/ 07 октября 2011
fabs(p - p0) <= tol * fabs(p);

Правильная формула для равных значений.В случае специальных диапазонов вы также должны заботиться о NaN и Inf.

#include <limits>
#include <cfloat>

inline bool Equal(const double& x, const double& y)
{
  const int classX (_fpclass(x));

  if(classX != _fpclass(y)) {
    return false;
  }

  switch(classX) {
      case _FPCLASS_SNAN:  
      case _FPCLASS_QNAN: 
      case _FPCLASS_NINF: 
      case _FPCLASS_NZ  : 
      case _FPCLASS_PZ  : 
      case _FPCLASS_PINF: 
        return true;
        break;
      case _FPCLASS_NN  : 
      case _FPCLASS_ND  : 
      case _FPCLASS_PD  : 
      case _FPCLASS_PN  : 
      default:
        break;
  }

  return fabs(x - y) <= std::numeric_limits<double>::epsilon() * fabs(x);
}
...