dvojno integriranje

Question

Zanima me kako se lahko dvojno integrira funkcijo z eno spremenljivko s pomočjo scipy.integrate.dlbquad?

Natančneje, imam izraz za pospešek, neznanka je čas(t), integrirati pa želim od 2 do 0 da dobim izraz za pot.

Answer 1

scipy.integrate.dlbquad je namenjen funkcijam z dvema spremenljivkama f(x, y), torej ni primerna za vaš namen. 

V vašem primeru gre v bistvu za reševanje nedoločenega integrala:

​

Imate dve možnosti reševanja

1. Obravnavamo integral F(t) kot diferencialno enačbo: 

      dF/dt = f(t)

      2. Integriranje F(t) čez celotno območje s spreminjajočo zgornjo mejo.

Obe rešitve si poglejmo na primeru:

xdd = sin(t) # pospesek, F(t)

Izračunati želimo hitrosti in pomike od časa 0 do 10 s. To pomeni, da moramo integrirati enkrat oz. dvakrat. Najprej definirajmo začetne pogoje:

x0 = 0  # začetni pomik
xd0 = -1  # začetna hitrost

1. Poglejmo si reševanje z diferencialnimi enačbami:

a) definirajmo vektor stanj (diferencialno enačbo 2. reda pretvorimo na dve diferencialne enačbe prvega reda):

def resevanje_diferencialne_enacbe(state,t):
    # razpkirajmo vektor stanj
    x = state[0]  # pomiki
    xd = state[1]  # hitrosti, 1. differencialna. enačba

    # izracunajmo pospeske
    xdd = np.sin(t)  # 2. differencialna enačba

    # vrnimo resitev 
    return [xd, xdd]

state0 = [x0, xd0]  # Začetni pogoji (pozor: 

# napačni začetni pogoji lahko močno vplivajo na rešitev)

time = np.arange(0.0, 10, 0.001)  # Definirana časovna os

# reševanje

state = integrate.odeint(resevanje_diferencialne_enacbe, state0, time)  

# Izris

plt.plot(time1, state)
plt.xlabel('cas (s)')
plt.ylabel('stanja')
plt.title('Resevanje diferencialne enacbe')
plt.ylim([-1, 1])
plt.legend(('$x$ (m)', '$\dot{x}$ (m/s)'))

Vidimo, da je rešitev pravilna!

2.) Sedaj si poglejmo še reševanje še z spreminjajočo se zgornjo mejo:

def func(x):
    return np.sin(x)

t = np.linspace(0, 10, 1001)

n = len(t)

hitrost = np.zeros(n)
pomik = np.zeros(n)

# Določitev hitrosti in pomikov

for i in range(n):  
    hitrost[i] = xd0 + integrate.quad(func, 0, t[i])[0]

    # Gaussove kvadrature ne moremo več uporabiti, saj nimamo podane funkcije

    pomik[i] = x0 + integrate.trapz(hitrost[:i+1], t[:i+1])   

# Izris

plt.plot(t, pomik)
plt.plot(t, hitrost)
plt.xlabel('cas (s)')
plt.ylabel('stanja')
plt.title('Resevanje diferencialne enacbe')
plt.ylim([-1, 1])
plt.legend(('$x$ (m)', '$\dot{x}$ (m/s)'))

Vidimo, da je tudi v tem primeru rešitev pravilna.

Comments

Kakor sem razumel moramo prvič integrirati z nedoločenim integralom, kar pa s scipy.integrate.quad ni mogoče. Šele pri drugem integriranju določimo meje. S katero funkcijo je torej to možno?

---

Sem popravil odgovor.

Answer 2

Funkcija scipy.integrate.dlbquad je namenjena za funkcije dveh sprejemljivk (x, y), ki jih integriramo po x in y, zato zahteva še y = gfun(x) in y =hfun(x), ki sta funkciji spodnje in zgornje meje krivulje.

Če imaš pospešek, ki je funkcija ene sprejemljivke rajši uporabi sympy funkcijo integrate(). Prvič narediš nedoločeni integral, da dobiš izraz za hitrost, drugič pa določeni integral v mejah (t, (0,2)) da izračunaš pot.

Comments

naloga zahteva numerični izračun, sympy pa izračuna samo simbolično, zato ne vem če se to sme uporabiti pri tej nalogi

---

Simbolično računa integral več kot 1 uro in ga ne uspe izračunati, kar pomeni možen je samo numeričen izračun. Torej z integrate() si ne moremo veliko pomagati