# HW5 Ex1

unknown
python
a year ago
2.8 kB
2
Indexable
Never
import numpy as np
import torch
import torch.nn as nn
from torch import optim
import matplotlib.pyplot as plt

# υπολογισμός παραγώγου του  u
def deriv_u(u, x):

# υπολογισμός παραγώγου του  w
def deriv_w(w, x):

neurons = 50
netu = nn.Sequential(nn.Linear(1, neurons), nn.Sigmoid(), nn.Linear(neurons,1,bias = False))

neurons = 50
netw = nn.Sequential(nn.Linear(1, neurons), nn.Sigmoid(), nn.Linear(neurons,1,bias = False))

# πεδίο ορισμού του t
x = torch.linspace(0.0, 2*np.pi, 100, requires_grad = True).reshape(-1,1)

loss_fun = loss_u + loss_w

optimizer = optim.Adam(list[netu.parameters()+list(netw.parameters())], lr = 0.01)

# u

epochs = 1000
loss_hist_u = []
u_hist = []
for i in range (epochs):
u = x + x*netu(x)
dudx = deriv_u(u,x)
loss_u = nn.MSELoss()
loss_w = nn.MSELoss()
loss = loss_fun(dudx, w(x))
loss_hist_u.append(loss.item())
loss.backward()
optimizer.step()
u_hist.append(u.detach())

# εύρεση της λύσης του  u  με κλασική μέθοδο

from scipy.integrate import solve_ivp
t = np.linspace(0, 2*np.pi, 100)
def ode1(x,u):
return w(x)

solu = solve_ivp(ode1, [0, 2*np.pi], [1], t_eval=t, method = 'RK45')

# w

epochs = 1000
loss_hist_w = []
w_hist = []
for i in range (epochs):
w = (1-x)+x*netw(x)
dwdx = deriv_w(w,x)
loss_u = nn.MSELoss()
loss_w = nn.MSELoss()
loss = loss_fun(dwdx, -u(x))
loss_hist_w.append(loss.item())
loss.backward()
optimizer.step()
w_hist.append(w.detach())

# εύρεση της λύσης του  w  με κλασική μέθοδο

from scipy.integrate import solve_ivp
t = np.linspace(0, 2*np.pi, 100)
def ode2(x,w):
return -u(x)

solw = solve_ivp(ode2, [0, 2*np.pi], [1], t_eval=t, method = 'RK45')

# plot

plt.figure()
plt.plot(x.detach(), u.detach(), 'b', label = "neural network")
plt.plot(t,solu.y[0], 'r--', label = "numerical solution with RK 45")
plt.legend()
plt.show()

plt.figure()
plt.plot(loss_hist_u, 'b')
plt.show()

plt.figure()
plt.plot(x.detach(), w.detach(), 'b', label = "neural network")
plt.plot(t,solw.y[0], 'r--', label = "numerical solution with RK 45")
plt.legend()
plt.show()

plt.figure()
plt.plot(loss_hist_w, 'b')
plt.show()

# HOMEWORK 5 ΑΣΚΗΣΗ 1