Решение разностных уравнений с использованием Z-преобразования

Решить разностные уравнения при помощи Z-преобразований в Symbolic Math Toolbox™ с этим рабочим процессом. Для простых примеров по Z-преобразованию смотрите ztrans и iztrans.

Определение: Z-преобразование

Z-преобразование функционального f (n) задано как

$F (z) = \sum_{n = 0}^{\infty} \frac{f (n)}{z^{n}} .$

Концепция: Использование символьных рабочих процессов

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

Объявить уравнения.
Решить уравнения.
Подстановочные значения.
Постройте график результатов.
Анализ результатов.

Рабочий процесс: Решите задачу «Рост кролика» используя Z-преобразование

Объявить уравнения

Можно использовать Z-преобразование, чтобы решить разностные уравнения, такие как известная задача «Рост кролика». Если пара кроликов созревает за один год, а затем каждый год производит другую пару кроликов, население кроликов p (n) за год n описывается этим разностным уравнением.

p (n +2) = p (n +1) + p (<reservedrangesplaceholder0>).

Объявите уравнение как выражение, предполагающее, что правая сторона 0. Потому что n представляет годы, предположим, что n является положительным целым числом. Это предположение упрощает результаты.

syms p(n) z
assume(n>=0 & in(n,'integer'))
f = p(n+2) - p(n+1) - p(n)

f =
p(n + 2) - p(n + 1) - p(n)

Решение уравнений

Найдите Z-преобразование уравнения.

fZT = ztrans(f,n,z)

fZT =
z*p(0) - z*ztrans(p(n), n, z) - z*p(1) + z^2*ztrans(p(n), n, z) - ...
        z^2*p(0) - ztrans(p(n), n, z)

Функция solve решает только для символьных переменных. Поэтому использовать solve, первая замена ztrans(p(n),n,z) с переменными pZT.

syms pZT
fZT = subs(fZT,ztrans(p(n),n,z),pZT)

fZT =
z*p(0) - pZT - z*p(1) - pZT*z - z^2*p(0) + pZT*z^2

Решение для pZT.

pZT = solve(fZT,pZT)

pZT =
-(z*p(1) - z*p(0) + z^2*p(0))/(- z^2 + z + 1)

Вычислим p (n) путем вычисления обратного Z-преобразования pZT. Упростите результат.

pSol = iztrans(pZT,z,n);
pSol = simplify(pSol)

pSol =
2*(-1)^(n/2)*cos(n*(pi/2 + asinh(1/2)*1i))*p(1) + ...
                 (2^(2 - n)*5^(1/2)*(5^(1/2) + 1)^(n - 1)*(p(0)/2 - p(1)))/5 - ...
                 (2*2^(1 - n)*5^(1/2)*(1 - 5^(1/2))^(n - 1)*(p(0)/2 - p(1)))/5

Подстановочные значения

Чтобы построить график результата, сначала замените значения начальных условий. Позвольте p(0) и p(1) быть 1 и 2, соответственно.

pSol = subs(pSol,[p(0) p(1)],[1 2])

pSol =
4*(-1)^(n/2)*cos(n*(pi/2 + asinh(1/2)*1i)) - (3*2^(2 - n)*5^(1/2)* ...
    (5^(1/2) + 1)^(n - 1))/10 + (3*2^(1 - n)*5^(1/2)*(1 - 5^(1/2))^(n - 1))/5

Графическое изображение результатов

Покажите рост населения кроликов с течением времени путем построения графиков pSol.

nValues = 1:10;
pSolValues = subs(pSol,n,nValues);
pSolValues = double(pSolValues);
pSolValues = real(pSolValues);
stem(nValues,pSolValues)
title('Rabbit Population')
xlabel('Years (n)')
ylabel('Population p(n)')
grid on

Анализ результатов

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

Потому что все функции в pSol может быть выражена в терминах exp, переписать pSol на exp. Упростите результат при помощи simplify с 80 дополнительные этапы упрощения. Теперь можно анализировать pSol.

pSol = rewrite(pSol,'exp');
pSol = simplify(pSol,'Steps',80)

pSol =
(2*2^n)/(- 5^(1/2) - 1)^n - (3*5^(1/2)*(1/2 - 5^(1/2)/2)^n)/10 + ...
    (3*5^(1/2)*(5^(1/2)/2 + 1/2)^n)/10 - (3*(1/2 - 5^(1/2)/2)^n)/2 + ...
    (5^(1/2)/2 + 1/2)^n/2

Визуальный осмотр pSol. Заметьте, что pSol является суммой терминов. Каждый термин является отношением, которое может увеличиваться или уменьшаться как n увеличивается. Для каждого термина можно подтвердить эту гипотезу несколькими способами:

Проверьте, установлен ли предел в n = Inf переходит к 0 или Inf при помощи limit.
Постройте график срока для увеличения n и проверяйте поведение.
Вычислим значение при большом значении n.

Для простоты используйте третий подход. Вычислим условия в n = 100, а затем проверьте подход. Во-первых, найти отдельные условия при помощи children, заменить на n, и преобразовать в double.

pSolTerms = children(pSol);
pSolTermsDbl = subs(pSolTerms,n,100);
pSolTermsDbl = double(pSolTermsDbl)

pSolTermsDbl =
   1.0e+20 *
    0.0000   -0.0000    5.3134   -0.0000    3.9604

Результат показывает, что некоторые условия 0 в то время как другие члены имеют большую величину. Предположим, что крупнозвуковые условия приводят к экспоненциальному поведению. Аппроксимация pSol с этими терминами.

idx = abs(pSolTermsDbl)>1; % use arbitrary cutoff
pApprox = pSolTerms(idx);
pApprox = sum(pApprox)

pApprox =
(3*5^(1/2)*(5^(1/2)/2 + 1/2)^n)/10 + (5^(1/2)/2 + 1/2)^n/2

Проверьте гипотезу, построив график ошибки приближения между pSol и pApprox. Как и ожидалось, ошибка переходит к 0 как n увеличивается. Этот результат демонстрирует, как символические вычисления помогают вам анализировать свою задачу.

Perror = pSol - pApprox;
nValues = 1:30;
Perror = subs(Perror,n,nValues);
stem(nValues,Perror)
xlabel('Years (n)')
ylabel('Error (pSol - pApprox)')
title('Error in Approximation')

Ссылки

[1] Andrews, L.C., Shivamoggi, B.K., Integral Transforms for Engineers and Applied Mathematicians, Macmillan Publishing Company, New York, 1986

[2] Crandall, R.E., Projects in Scientific Computation, Springer-Verlag Publishers, New York, 1994

[3] Strang, G., Введение в прикладную математику, Wellesley-Cambridge Press, Wellesley, MA, 1986

Документация