Как визуализировать решения уравнения u(x,t)?

Question

[Глеб]

я уже строил графики в случае x and f(x) , но появились проблемы и куча ошибок, при построении x и u(x,t).
Хочу визуализировать сеточную функцию, как возможно сделать? просто plot(x, res) не работает.

from sympy.interactive import printing
printing.init_printing(use_latex=True)
from sympy import*
import sympy as sp
import matplotlib.pyplot as plt

x, t = sp.symbols('x t')
a, b, c, T = sp.symbols('a b c T')

u = sp.Function('u')(x, t)
u0 = sp.Function('u')(x, 0)
fi = sp.Function('fi')(x, t)
psi = sp.Function('psi')(x)

a = 0.2; b = 2.2; c = -0.4; T = 1.9
fi = log(17*t) - 2*x
psi = exp(-x/2)

diffeq = Eq(u.diff(t) - c*u.diff(x), fi)
diffeq1 = Eq(u0, psi)

display(diffeq)
display(diffeq1)

pdsolve(diffeq)
res = pdsolve(diffeq)

Comments

[Алан Гибизов]

Пока что наблюдаю, что вы определяете fi и psi как некие объекты-функции из sympy, а в следующих строках переопределяете их как просто уравнения, т.е. предыдущее действие

fi = sp.Function('fi')(x, t)
psi = sp.Function('psi')(x)

мне кажется, теряет смысл.

[Глеб]

Алан Гибизов, учту, спасибо. Ну, я в интернете не нашел , как решать неоднородные дифуры в частных в .py, подумал, что так можно.

[Алан Гибизов]

Андрей, признаться, очень отдаленно понимаю, что за матмагию вы вытворяете...
Правильно я понимаю, что Вы хотите получить в итоге что-то вроде вот такого графика?

[Глеб]

Алан Гибизов,
скорее вот такой:

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

import numpy as np
import matplotlib.pyplot as plt

def solve(n,m):
    h = 0.001
    c = -0.4
    T = 1.9

    xmin = 0.2
    xmax = 2.2
    
    r = 0.1
    tau = r*h
    r = tau*c/h
    t = []

    for i in range(n): 
      t.append(i*tay)

    X = np.linspace(xmin,xmax,m)

    U = np.zeros((n, m))
    
    for i in range(0, n):
        U[i, 0] = np.exp(-X[i]/2)
    
    for j in range(1, m):
        for i in range(1, n):
            U[i, j] = (1-c*r)*U[i,j-1] + c*r*U[i-1,j] + tau*(np.log(1+t[i])-2*X[j])

    plt.plot(X, U[:,i])
    plt.show()
    
solve(100,100)

[Алан Гибизов]

Андрей, тьфу, нашел опечатку - не tay а tau и всё заработало.

[Глеб]

Алан Гибизов, ну если что, я вот этот фрагмент реализую на [a,b] и константой еще c.

[Алан Гибизов]

Андрей, мне, конечно, очень трудно с непривычки продираться сквозь ясные любому математику эпсилон-символы :) но у меня закралось следующее соображение:
для решения этого уравнения надо подобрать какое-то одно решение, и все остальные решения по сетке строятся от этого решения. Я не уверен, это у меня из беглого просмотра приведенных вами сканов и вот этого куска кода, нарытого в инете, навеяло:

from sympy import *
x, b1, b2 = symbols("x b1 b2")
f = x/(x+exp(b1-b2*x))
res = {b1:29.3930964972769,b2:0.327159886574049}
plot(f.subs(res), (x, 0, 100))

Возможно, я вообще не то несу, тут конечно нужен математик-аналитик, желательно дружный с пайтоном...

[Глеб]

Алан Гибизов, Идея про 1 слой и построение по нему верная. Но нам и задан слой U = Exp(-x/2) в момент времени t=0.
Далее мы по формуле U[i;j] как раз строим остальные слои, рекуррентно как раз.
Я и пробую этого добиться, но подводят конечно знания python. Все можно решить одними циклами for , но получается он их очень долго считать будет.

[Алан Гибизов]

Поигрался немного с этим самым sympy,
вот чего нарыл:
1.

вариант кода

from sympy.interactive import printing
printing.init_printing(use_latex=True)
from sympy import*
import sympy as sp
import matplotlib.pyplot as plt

x, t = sp.symbols('x t')
a, b, c, T = sp.symbols('a b c T')

u = sp.Function('u')(x, t)
u0 = sp.Function('u')(x, 0)
phi = sp.Function('phi')(x, t)
psi = sp.Function('psi')(x)

a = 0.2; b = 2.2; c = -0.4; T = 1.9
#phi = log(17*t) - 2*x
#psi = exp(-x/2)

diffeq = Eq(u.diff(t) - c*u.diff(x), phi)
diffeq1 = Eq(u0, psi)

display(diffeq)
display(diffeq1)

pdsolve(diffeq)
res = pdsolve(diffeq)

2. после твоего кода или после вышеуказанного кода можно сделать print(res) и будет в разных случаях разный вывод.

[Алан Гибизов]

кроме того, у res есть метод args, который, в свою очередь, выдает кортеж из двух аргументов. Я не знаю, можно ли их как-то еще использовать, но видимо можно.
т.е.
print(res.args)
печатает кортеж аргументов. К ним можно обращаться по индексам 0 и 1:
print(res.args[1])

[Алан Гибизов]

Андрей, кстати, обрати внимание, что написание fi не дает правильной буквы, а написание phi - дает правильную букву эпсилон-записей.