![]() Главная страница Случайная страница КАТЕГОРИИ: АвтомобилиАстрономияБиологияГеографияДом и садДругие языкиДругоеИнформатикаИсторияКультураЛитератураЛогикаМатематикаМедицинаМеталлургияМеханикаОбразованиеОхрана трудаПедагогикаПолитикаПравоПсихологияРелигияРиторикаСоциологияСпортСтроительствоТехнологияТуризмФизикаФилософияФинансыХимияЧерчениеЭкологияЭкономикаЭлектроника |
Метод квадратичной интерполяции.
Исходные функции называются интерполируемыми. Полиномы, полученные после интерполяции, называются интерполяционными полиномами. Чаще всего узлы интерполяции задаются на равномерной сетке. Функция и полином при любых методах интерполяции обязательно в точности совпадают в узлах интерполяции. Если значения функции определяются внутри узлов интерполяции, то это процедура называется интерполяцией, если функция определяется вне узлов рассматриваемого участка интерполяции, то такая процедура называется экстраполяцией. Если для интерполяции используются: 1) две точки, то это линейная интерполяция 2) три – квадратичная и т. д. Пусть исходная функция обозначена Проведем интерполяцию ее полиномами Лагранжа. Линейная интерполяция: 1) Если x=x0, то F(x0)=f(x0)=f0 . Если x=x1, то F(x1)=f(x1)=f1 . Квадратичная интерполяция: 2) Легко проверить, что F(x) проходит точно через узловые точки, задавая последовательно вместо Преобразуем знаменатели: Определим производную от
полагаем
Собираем слагаемые с hm:
Подставим h в «узлы» и определим hm:
Новое значение
function[p, ym, dp, dy, P]=kwadrmin(fungold, a, b, del, eps) %Вход: %fungold функция вводится как строка 'fingold' %[a, b]-начальный интервал поиска минимума %del-допустимое отклонение для абсцисс %eps-допустимое отклонение для ординат %Выход: %p, ym-абсцисса и ордината минимума %dp, dy-грани ошибок по абсциссе и ординате %P-вектор итераций %Обращение из МАTLAB [p, ym, dp, dy, P]=kwadrmin('fungold', 0.5, 1.5, 0.001, 0.0001) p0=a; maxk=30; big=1e6; maxj=20; err=1; k=1; P(k)=p0; cond=0; h=1; if(abs(p0)> 1e4), h=abs(p0)/1e4; end while(k< maxk& err> eps& cond~=5) f1=(feval('fungold', p0+0.00001)-feval('fungold', p0-0.00001))/0.00002; if(f1> 0), h=-abs(h); end p1=p0+h; p2=p0+2*h; pmin=p0; y0=feval('fungold', p0); y1=feval('fungold', p1); y2=feval('fungold', p2); yymin=y0; cond=0; j=0; % while(j< maxj& abs(h)> del& cond==0) if(y0< =y1), p2=p1; y2=y1; h=h/2; p1=p0+h; y1=feval('fungold', p1); else if(y2< y1), y1=y2; p1=p2; h=2*h; p2=p0+2*h; y2=feval('fungold', p2); else cond=-1; end end j=j+1; if(abs(h)> big|abs(p0)> big), cond=5; end end if(cond==5), pmin=p1; ymin=feval('fungold', p1); else % d=4*y1-2*y0-2*y2; if(d< 0), hmin=h*(4*y1-3*y0-y2)/d; else hmin=h/3; cond=4; end pmin=p0+hmin; ymin=feval('fungold', pmin); h=abs(h); h0=abs(hmin); h1=abs(hmin-h); h2=abs(hmin-2*h); % if(h0< h), h=h0; end if(h1< h), h=h1; end if(h2< h), h=h2; end if(h==0), h=hmin; end if(h< del), cond=1; end if(abs(h)> big|abs(pmin)> big), cond=5; end % e0=abs(y0-ymin); e1=abs(y1-ymin); e2=abs(y2-ymin); if(e0~=0& e0< err), err=e0; end if(e1~=0& e1< err), err=e1; end if(e2~=0& e2< err), err=e2; end if(e0~=0& e1==0& e2==0), err=0; end if(err< eps), cond=2; end p0=pmin; k=k+1; P(k)=p0; end if(cond==2& h< del), cond=3; end end p=p0; dp=h; ym=feval('fungold', p); dy=err;
|