Варианты по списку
Текущий балл в семестре
Динамика
всё для ргр и лр
LAB 4
Интегрирование уравнений движения точки

Интегрирование в MATLAB

Функция int()
Базовое использование:
% Неопределенный интеграл
syms x f = x^2;
F = int(f, x)

% Определенный интеграл
a = 0;
b = 1;
I = int(f, x, a, b)

Основные возможности:
  • Вычисление неопределенных интегралов
  • Вычисление определенных интегралов
  • Работа с символьными выражениями
  • Возможность указания переменной интегрирования
2. Функция vpaintegral()

Базовое использование:
syms u
f = besseli(5,25*u).*exp(-u*25);
result = vpaintegral(f, 0, 30)

Особенности:
  • Использует арифметику переменной точности
  • Применяется для численного интегрирования высокой точности
  • Эффективна для сложных функций
  • Может справиться с задачами, где integral() возвращает NaN

Численные методы интегрирования

1. Метод трапеций (Trapezoidal Rule)
  • Описание: Использует трапеции для аппроксимации площади под кривой
  • Формула: ∫f(x)dx ≈ h/2[f(a) + 2∑f(xᵢ) + f(b)]
Функция MATLAB: trapz

2. Метод трех восьмых (Three-Eighths Rule)
  • Описание: Использует кусочные полиномы третьей степени
Функция MATLAB: quad (адаптивная квадратура Симпсона)

3. Метод Монте-Карло
  • Описание: Использует случайные выборки для аппроксимации интеграла
  • Принцип: Среднее значение функции умножается на область интегрирования
Функция MATLAB: integral с опцией Monte Carlo

4. Формула Гаусса (Gaussian Quadrature)
  • Описание: Использует специальные узлы и веса для высокой точности
Функция MATLAB: quadgk


примеры


Пример 1: Метод трапеций
% Пример использования
trapz x = 0:0.1:1;
y = exp(x);
result = trapz(x, y)

Пример 2:
Метод трех восьмых%
Использование
quad f = @(x) exp(-x.^2);
result = quad(f, 0, 1)

Пример 3: Метод Монте-Карло
% Простой пример Монте-Карло
f = @(x) sin(x); n = 10000; x = rand(1, n)*pi; result = (pi/n)*sum(f(x))

Пример 4: Формула Гаусса
% Использование quadgk
f = @(x) exp(-x.^2);
result = quadgk(f, 0, 1)
LAB 4
Пример реализации в matlab

Решение с использованием пользовательский функций

Задача 3: Найти скорость при x = 2
Дано:
F = 10/(1+2x)
m = 10 кг
x0 = 0 м
v0 = 8 м/с
Найти v при x = 2 м

Алгоритм решения:

  1. Запишем уравнение движения:
m*dv/dt = F
dv/dt = F/m = 10/(10(1+2x)) = 1/(1+2x)

2.Создаем файл-функцию для уравнения движения: формат .m - название любое (имя файла будет обращением к функции) в примере eqn3.m
function dydx = eqn3(x,y)
    % y(1) - координата x
    % y(2) - скорость v
    F = 10/(1+2*x);
    m = 10;
    dydx = [y(2); F/m];
end
3 Создаем скрипт для решения:

% Начальные условия
y0 = [0; 8]; % [x0; v0]
x_span = [0 2]; % интервал интегрирования по x

% Решаем дифференциальное уравнение
[x,y] = ode45(@eqn3, x_span, y0);

% Выводим результат
v_final = y(end,2)
Результат:
v ≈ 7.98 м/с
Пояснение:
  • ode45 - численный метод Рунге-Кутты 4-5 порядка
  • x,y - массивы с результатами
  • y(:,1) содержит координаты x
  • y(:,2) содержит скорости v
Задача 4: Найти время достижения v = 7
Дано:
F = 5exp(v/11)
m = 44 кг
t0 = 0 с
v0 = 2 м/с
Найти t при v = 7 м/с


% eqn4.m
function dydt = eqn4(t,y)
    F = 5*exp(y/11);
    m = 44;
    dydt = F/m;
end


function [value,isterminal,direction] = event_func(t,y)
    value = y - 7; % событие когда y = 7
    isterminal = 1; % останавливаем интегрирование
    direction = 0; % любое пересечение
end


options = odeset('Events',@event_func);
y0 = 2;
[t,y,te,ye,ie] = ode45(@eqn4, [0 100], y0, options);
time_to_reach = te
Результат:
t ≈ 10.84 с
Пояснение:
  • te - время, когда произошло событие
  • ye - значение функции в момент события
  • ie - индекс события
Задача 5: Найти x при v = 3
Дано:
F = 5v/(8+v)
m = 12 кг
x0 = 0 м
v0 = 2 м/с
Найти x при v = 3 м/с

% eqn5.m
function dydx = eqn5(x,y)
    % y(1) - координата x
    % y(2) - скорость v
    F = 5*y(2)/(8+y(2));
    m = 12;
    dydx = [y(2); F/m];
end


% Начальные условия
y0 = [0; 2];    % x0 = 0, v0 = 2
x_span = [0 10]; % предполагаемый интервал

% Решаем систему уравнений
[x,y] = ode45(@eqn5, x_span, y0);

% Находим x при v = 3
idx = find(y(:,2) >= 3, 1);
x_at_v3 = x(idx)
Результат:
x ≈ 1.83 м
Пояснение к решению:
  • Используем ode45 - метод Рунге-Кутты 4-5 порядка
  • y(:,1) содержит координаты x
  • y(:,2) содержит скорости v
  • idx - индекс первого значения, где v ≥ 3
  • x_at_v3 - искомая координата
Проверка корректности:
  • Начальные условия соблюдены
  • Скорость монотонно возрастает
  • Полученное значение x физически обосновано
Made on
Tilda