Вычислительные методы

Партнерка на США и Канаду по недвижимости, выплаты в крипто

30% recurring commission
Выплаты в USDT
Вывод каждую неделю
Комиссия до 5 лет за каждого referral

ВЫЧИСЛИТЕЛЬНЫЕ МЕТОДЫ

ПРАКТИЧЕСКАЯ РАБОТА № 1

АНАЛИЗ ОШИБОК ОКРУГЛЕНИЯ НА ПРИМЕРЕ ВЫЧИСЛЕНИЯ СИНУСА С ПОМОЩЬЮ РАЗЛОЖЕНИЯ В РЯД ТЕЙЛОРА

1. Ошибки округления при вычислениях на компьютере.

Одним из коренных отличий численных расчетов на компьютере от аналитических вычислений заключается в том, что точность записи чисел на компьютере ограничена. Во-первых, ограничен диапазон их представления, а во-вторых, представление вещественных и комплексных чисел не всегда является точным. Например, все вещественные числа типа real(4) расположены в интервале от –3.402328E+38 до 3.402328E+38. Кроме того, между числами представимыми на компьютере существуют разрывы, то есть, для каждого числа value > 0 существует минимальное положительное число emin, которое в сумме с value даст число большее, чем value. Для положительных чисел e меньших, чем emin

value + e = value.

Чем меньше value, тем меньше emin.

Это свойство машинной арифметики можно трактовать и по другому. Запись числа производится 7 значащими цифрами для обычной точности real(4) и 15 значащими цифрами для двойной точности real(8).

Рассмотрим ставший уже классическим пример вычисления синуса (или косинуса) с помощью разложения в ряд Тейлора:

Считается, что этот ряд годится для любого конечного x, а ошибка, происходящая от ограничения ряда конечным числом членов, не превосходит по модулю первого отброшенного члена.

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

НЕ нашли? Не то? Что вы ищете?

2. Алгоритм вычисления синуса.

Ниже представлен алгоритм вычисления синуса с помощью суммирования ряда Тейлора. Используются следующие переменные:

Degree – входное значение угла в градусах. EPS – точность вычислений, задаваемая пользователем. X – входное значение угла в радианах. PI – число p. S – сумма ряда. TERM – текущий член ряда. ZN – параметр для расчета знаменателя множителя перехода от члена ряда ai к ai+1. SINX – значение sin(x), вычисленное с помощью встроенной программы.

3. Задание.

3.1. Написать на ФОРТРАНЕ программу расчета sin(x) с помощью суммирования ряда Тейлора.

3.2. Проанализировать ошибку в вычислении синуса в зависимости от величины x для программ с обыкновенной точностью и c двойной точностью. Считать точным значение sin(x), полученное при обращении к встроенной в Intel Visual Fortran функции sin(x). Результаты представить в виде графиков.

ПРАКТИЧЕСКАЯ РАБОТА № 2

НАХОЖДЕНИЕ КОРНЕЙ УРАВНЕНИЯ МЕТОДОМ ДЕЛЕНИЯ ОТРЕЗКА ПОПОЛАМ

1. Описание метода и его реализация на ФОРТРАН’е

Подпись: Метод состоит в том, что заданный отрезок вероятного нахождения корня уравнения делится пополам (координата ), а затем производится проверка на совпадение знака значения функции в точках и . Если эти знаки совпадают, то переменной присваивается значение . Если знаки не совпадают, то присваивается значение . Этот процесс циклически продолжается до тех пор, пока значение не станет меньше некоторого заданного , определяющего точность нахождения корня уравнения. Блок-схема подпрограммы представлена на рис. 2. Ее текст на ФОРТРАН’е можно представить в следующем виде.

Рис. 1.

subroutine eqsolv (a, b, fun, eps, xs)

a1=a

b1=b

1 xs=(a1+b1)/2.

if (fun(a1)fun(xs).gt.0.) go to 2

b1=xs

go to 3

2 a1=xs

3 if (abs(a1-b1).gt. eps) go to 1

end subroutine eqsolv

Здесь – левая граница отрезка, – правая граница отрезка, – имя подпрограммы-функции для уравнения , – точность вычисления корня.

2. Пример

Найдем решение уравнения на отрезке [0., 1.5]. Программа для нахождения корня с использованием подпрограммы eqsolv может быть записана в виде:

program example1

external fun1

a=0.

b=1.5

eps=1.e-5

call eqsolv(a, b, fun1, eps, xs)

print , ‘а = ‘, a, ‘b = ‘, b, ‘корень‘, xs

end program example1

function fun1(x)

fun1=cos(x)-x

end function fun1

После работы программы будет найдено значение корня, равное 0.7390881.

3. Использование библиотечной подпрограммы ZBREN

В математической библиотеке IMSL, входящей в версию Fortran PowerStation 4.0, имеется подпрограмма ZBREN вычисления корня внутри отрезка [a, b]. В ней используется метод деления отрезка пополам, метод секущих и обратная квадратичная интерполяция [1]. На каждом шаге производится выбор одного из этих методов для вычисления следующего приближения.

Обращение к подпрограмме имеет вид:

CALL ZBREN(f, errabs, errrel, a, b, maxfn).

Параметры подпрограммы ZBREN:

f – внешняя функция f(x), корень которой нужно найти. Ее аргумент нельзя изменять внутри функции f ;

errabs – значение абсолютного отклонения для первого критерия остановки процесса нахождения корня . Значение errabs на входе можно положить равным нулю;

errrel – относительное отклонение для второго критерия остановки. Если относительное отклонение меньше или равно errrel, то значение корня полагается равным a;

a, b – две начальные точки, в которых f(a) и f(b) имеют разные знаки. Их значения изменяются в процессе вычислений. При выходе из подпрограммы значение b равно наилучшему приближению корня.

maxfn – на входе равно максимальному числу вычислений значения f(x), необходимому для сходимости процесса. На выходе из подпрограммы maxfn равно реальному числу вычислений значения f(x).

4. Задание

Написать программу решения уранения, рассмотренного в пункте 2, с помощью подпрограммы ZBREN.

ПРАКТИЧЕСКАЯ РАБОТА № 3

ЧИСЛЕННОЕ РЕШЕНИЕ СИСТЕМЫ НЕЛИНЕЙНЫХ УРАВНЕНИЙ

1. Метод Ньютона

Подпись: Рассмотрим снова нелинейное уравнение . Можно найти решение этого уравнения, проводя касательную к кривой в точке до пресечения с осью (рис. 2). В найденной точке персечения снова проводим касательную и повторяем этот процесс до тех пор, пока значение на очередном шаге не станет близким к корню уравнения .

Такой вычислительный процесс называется методом Ньютона. Его можно представить в виде итерационной формулы

. (2.1)

Этот метод можно обобщить на случай системы двух нелинейных уравнений с двумя неизвестными

Якобиан системы в точке равен

Будем предполагать, что он отличен от нуля. Тогда итерационные формулы принимают вид

Итерационный процесс продолжается до тех пор, пока два последовательных приближения не будут близки друг к другу с заданной точностью.

Аналогичным образом можно построить алгоритм решения системы нелинейных уравнений.

2. Использование библиотечной подпрограммы NEQNF

В подпрограмме NEQNF реализована модификация метода Ньютона с испольованием конечно-разностной аппроксимации якобиана.

Обращение к подпрограмме имеет вид:

CALL NEQNF(fcn, errrel, n, itmax, xguess, x, fnorm).

Параметры подпрограммы NEQNF:

fcn – подпрограмма вычисления значения функций в левых частях системы уравнений. Обращение к ней должно иметь вид CALL fcn(x, f, n), где x – точка, в которой вычисляются значения функций (не должна изменяться внутри fcn), f – вычисленные значения, n – размерность векторов x и f. Подпрограмма fcn должна быть описана оператором EXTERNAL в основной программе;

errel – относительное отклонение двух последовательных приближений решения системы, при котором останавливается итерационная процедура;

n – число уравнений и число неизвестных решаемой системы;

itmax – максимально допустимое число итераций;

xguess – вектор размерности n, содержащий начальное значение корней системы;

x – вектор размерности n. После окончания работы подпрограммы NEQNF содержит значения корней решаемой системы уравнений;

fnorm – скаляр, равный значению в точке .

3. Задание

Написать программу решения системы двух нелинейных уравнений.

ПРАКТИЧЕСКАЯ РАБОТА № 4

ВЫЧИСЛЕНИЕ ОРЕДЕЛЕННЫХ ИНТЕГРАЛОВ

1. Метод Симпсона

В этом методе интервал интегрирования разбивается на равных отрезков длиной . Для вычисления интеграла через каждые три точки разбиения проводится парабола. Складывая площади криволинейных трапеций высотой , легко получить следующее выражение для приближенного значения интеграла от функции

где - значения подынтегральной функции в точках разбиения.

Подпрограмму вычисления определенного интеграла по методу Симпсона можно можно представить на ФОРТРАН’е в следующем виде

subroutine simp (n, a, b, fun, s)

b1=b

dx=(b-a)/n

dx2=2.*dx

x1=a-dx

fa=fun(a)

s=0.

do i=2,n,2

x1=x1+dx2

x2=x1+dx

fb=fun(x2)

s=s+2.*fun(x1)+fb

end do

s=(s+(fa-fb)/2.)*2*dx/3.

end subroutine simp

Здесь - число отрезков разбиения, – левая граница интервала, – правая граница интервала, – имя подпрограммы-функции для вычисления подынтегральной функции , - значение интеграла на выходе из подпрограммы.

2. Пример

Вычислим интеграл . Программа для вычисления интеграла с использованием подпрограммы simp может быть записана в виде:

program example3

external fun

a=1.

b=3.

n=100.

call simp (n, a, b, fun, s)

print , ‘а = ‘, a, ‘b = ‘, b, ‘интеграл‘, s

end program example3

function fun(x)

fun=exp(x)*cos(x)

end function fun1

После работы программы будет найдено значение интеграла, равное *******.

3. Метод Гаусса

Если подынтегральная функция на одних участках интервала изменяется медленно, а на других – быстро, то деление на отрезки равной длины является невыгодным. Более рационально использовать конечно-разностное представление интеграла в виде

где - некоторые неравноотстоящие узловые точки, - соответствующие весовые множители. На этом и основан метод Гаусса. Можно показать, что для наилучшего приближения к значению интеграла значения должны совпадать с корнями полинома Лежандра , а значения нужно находить из решения некоторой системы линейных уравнений с коэффициентами, зависящими от [1].

4. Использование библиотечной подпрограммы QDAGS

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

Обращение к подпрограмме имеет вид:

CALL QDAGS(f, a, b, errabs, errrel, result, errest).

Параметры подпрограммы QDAGS:

f – подпрограмма-фунуция для вычисления значений подынтегральной функции. Подпрограмма f должна быть описана оператором EXTERNAL в основной программе;

a – нижний предел интегрирования;

b – верхний предел интегрирования;

errabs – желаемая абсолютная погрешность;

errrel – желаемая относительная погрешность;

result – вычисленное значение интеграла;

errest – оценка значения абсолютной погрешности.

5. Задание

Написать программу вычисления интеграла, рассмотренного в пункте 2, с помощью подпрограммы QDAGS.

ПРАКТИЧЕСКАЯ РАБОТА № 5

ИНТЕГРИРОВАНИЕ ОБЫКНОВЕННЫХ ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙ МЕТОДОМ ЭЙЛЕРА

1. Метод Эйлера

Рассмотрим обыкновенное дифференциальное уравнение

(1)

с начальным условием . Два значения решения в точках и связаны соотношением

. (2)

Заметим, что выражение (2) является точным. В численных методах рассчитывается приближенное значение интеграла. В методе Эйлера [1] подынтегральная функция в выражении (2) заменяется ее значением в точке

. (3)

Из выражения (3) следует, что метод Эйлера относится к явным методам и имеет первый порядок точности по шагу . Для уравнений с затуханием ()он является устойчивым при выполнении условия [1]

Метод Эйлера применим и к решению системы дифференциальных уравнений первого порядка вида (1). В этом случае представляет собой вектор некоторой размерности .

ПРАКТИЧЕСКАЯ РАБОТА № 6

ИНТЕГРИРОВАНИЕ ОБЫКНОВЕННЫХ ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙ МЕТОДОМ

«С ПЕРЕШАГИВАНИЕМ»

В методе «с перешагиванием» испльзуется двухшаговая процедура:

(1)

Начальный шаг интегрирования для определения проводится по методу Эйлера

Метод «с перешагиванием» является центрированным методом второго порядка точности по шагу интегрирования.

Методы Рунге-Кутты, основанные на центрировании по времени интеграла с помощью некоторой многошаговой процедуры.

Приведем пример двухшагового метода:

2. Промежуточный шаг:

(1)

2. Основной шаг:

Метод, основанный на усреднении значений производных:

Четырехшаговый метод Рунге-Кутты

где

Задание.

Используя метод «с перешагиванием» решить уравнение движения электрона в поле продольной электрической волны, распространяющейся вдоль x.

(1)

с начальными условиями , ,

где v – скорость электрона, e и m – заряд и масса электрона, соотметственно, E0 – амплитуда напряженности поля волны, w - круговая частота волны, - волновой вектор. Волна распространяется с постоянной скоростью .

Рассмотреть следующие случаи:

2. Начальная скорость электрона много больше скорости волны.

3. Начальная скорость электрона много меньше скорости волны.

4. Начальная скорость электрона близка к скорости распространения волны.

....................................................................................

Для решения системы (1) численными методами следует его представить в безразмерном виде. Введем безразмерные переменные:

, , , .

Тогда система (1) приобретает вид

(2)

Наибольший интерес представляет решение системы (2) в случае 3. Если скорость электрона немногим меньше скорости волны (электрон отстает от волны), но находится в ускоряющем поле волны, то очевидно, он будет приобретать энергию и «догонять» волну.

На рис. 1 представлены: разность между координатой электрона и координатой фронта волны (слева) и электрическое поле, в котором находится электрон (справа).

Уравнения для фазы и энергии в случае электронного циклотронного резонанса

В однородном магнитном поле, направленном вдоль оси Z, с индукцией и однородном правовращающемся с частотой электрическом поле E движется электрон. Движение в плоскости XY. .

Изменение полной энергии электрона, обусловлено работой сил электрического поля

где - релятивистский фактор, - скорость электрона, - фаза электрона (угол между вектором скорости электрона и вектором напряженности электрического поля).