Как найти максимальное значение функции в Sympy?

Question

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

Первый

tau,t,t_r,omega,p0=symbols('tau,t,t_r,omega,p0',positive=True) 
h=expand(sin(omega*(t-tau))) 
f=simplify(integrate(p0*tau/t_r*h,(tau,0,t_r))+integrate(p0*h,(tau,t_r,t))) 

Конечной целью является получение максимального абсолютного значения f (Переменная t). Прямой путь:

df=diff(f,t) 
sln=solve(simplify(df),t) 
simplify(f.subs(t,sln[1])) 

Вот результат, я пробовал много способов, но я не могу упростить дальнейшее продвижение.

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

df=expand_trig(diff(expand(f)**2,t)) 
sln=solve(df,t) 
simplify(f.subs(t,sln[2])) 

Похоже, что ответ почти такой же, как раз в другой форме.

Ожидаемый ответ является sinc функция плюс константа следующим образом:

Таким образом, вопрос заключается в том, чтобы получить окончательное представление.

Второй может быть немного сложнее. Вопрос можно уменьшить, чтобы найти максимальное значение f=sin(pi*t/t_r)-T/2/t_r*sin(2*pi/T*t), в котором t_r и T являются двумя параметрами. Максимум находится в разном пике при изменении отношения t_r и T. И я не нашел способ решить эту проблему в Sympy. Любое предложение? Ответ может быть представлен на следующем рисунке.

Answer 1

Проблема заключается в log(exp(I*omega*t_r/2)) срок. SymPy не уменьшает это до I*omega*t_r/2. SymPy не упрощает это, потому что в общем случае log(exp(x)) != x, а скорее log(exp(x)) = x + 2*pi*I*n для некоторого целого n. Но в этом случае, если вы замените log(exp(I*omega*t_r/2)) на omega*t_r/2 или omega*t_r/2 + 2*pi*I*n, он будет таким же, потому что он просто добавит 2*pi*I*n внутри sin.

Я не мог понять, какие функции, которые заставляют это упрощение, но самый простой способ это просто сделать замену:

In [18]: print(simplify(f.subs(t,sln[1]).subs(log(exp(I*omega*t_r/2)), I*omega*t_r/2))) 
p0*(omega*t_r - 2*sin(omega*t_r/2))/(omega**2*t_r) 

Это выглядит как ответ вы ищете, за исключением абсолютного значения (Я не уверен, откуда они должны прийти).