Оценить неопределенный интеграл численно в Matlab / Mathematica, что он не может сделать символически

Question

Я пытаюсь вычислить интеграл от функции в Matlab и Mathematica, которую программное обеспечение не может сделать символически.

Вот мой код MatLab до сих пор, но я понимаю, что он может быть не очень полезным, как есть.

f = @(t) asin(0.5*sin(t));
a = @(t) sin(t);
F = int(f,t)   % Matlab can't do this
F = 
int(asin(sin(t)/2), t)
A = int(a,t)   % This works
A =
-cos(t)

dt = 1/(N-1); % some small number
for i=1:N
    F(i) = integral(f,(i-1)*dt,i*dt);
    A(i) = integral(a,(i-1)*dt,i*dt);
end

Оба вычисления в цикле for дают грубое приближение f или a, а не их интегралов после умножения на dt.

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

Я пытаюсь получить или символическое уравнение для интеграла, или приближение интеграла отфункция в каждом месте.

Так что мой вопрос тогда будет , если у меня есть функция f, что MatLab и Mathematica не могут легко взять интеграл от

1. Могу ли я аппроксимировать интеграл напрямую с помощью интегрального калькулятора, кроме стандартных? * (int, integral, trapz)

или

Могу ли я сначала аппроксимировать функцию конечными разностями, а затем символически вычислить интеграл?

Answer 1

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

Для двух функций f и g см.ниже

T = 1;  % Period
NT = 1;  % Number of periods
dt = 0.01; % time interval
time = 0:dt:NT*T;  % time

syms t
x = K*sin(2*pi*t+B);   % edit as appropriate

% f = A/tanh(K)*tanh(K*sin(2*pi*t+p))
% g = A/asin(K)*asin(K*sin(2*pi*t+p))

найдено формул здесь

f = A1/tanh(K1)*(2^(2*1)-1)*2^(2*1)*bernoulli(2*1)/factorial(2*1)*x^(2*1-1);
% |K1|<pi/2
g = A2/asin(K2)*factorial(2*0)/(2^(2*0)*factorial(0)^2*(2*0+1))*x^(2*0+1);
% |K2|<1

в принятом ответе нет таких ограничений

N = 60;
for k=2:N
    a1 = (2^(2*k)-1)*2^(2*k)*bernoulli(2*k)/factorial(2*k);
    f = f + A1/tanh(K1)*a1*x^(2*k-1);

    a2 = factorial(2*k)/(2^(2*k)*factorial(k)^2*(2*k+1));
    g = g + A2/asin(K2)*a*x^(2*k+1);
end

MATLAB может рассчитать sin^n(t) для n, являющегося целым числом.

F = int(f,t);
phi = double(subs(F,t,time));

G = int(g,t);
psi = double(subs(G,t,time));

Answer 2

Ваш код почти в порядке, просто

for i=1:N
    F(i) = integral(f,0,i*dt);
end

Вы также можете сделать

F(1)=integral(f,0,dt)
for i=2:N
    F(i) = F(i-1)+integral(f,(i-1)*dt,i*dt);
end

Второй вариант, безусловно, более эффективен

Поскольку примитив действительноF (x) = int (f (x), 0, x) (0 определяет некоторую постоянную), и для достаточно малого dx вы показали, что f (x) = int (f (x), x, x + dx)/ дх я.Вы доказали, что встроенная функция MATLAB выполняет свою работу.

Например, давайте возьмем = , функция, приведенная выше, вычислит , если вы хотите вычислить , просто замените 0 выше на постоянную a, которая вам нравится.

сейчас и так выдолжен получить F, содержащий дискретность