.. _akustika_rupa: Rešavanje na domenu oblika kvadrata sa šupljinom ================================================ U odnosu na prethodni primer :ref:`akustika_ravan` dodajemo šupljinu u sredini kvadratnog domena i propisujemo odgovarajuće Nojmanove granične uslove na obodu šupljine. Da se podsetimo, domen problema :math:`\Omega` je kvadrat stranice :math:`L, \, L=1`, iz koga isključujemo krug poluprečnika :math:`R=\frac{1}{4}`. Za talasni broj :math:`k_0=2 \pi n` i :math:`n=1`, rešavamo Helmholcovu jednačinu: .. math:: -\frac{\partial^2 u}{\partial x^2} - \frac{\partial^2 u}{\partial y^2} - k_0^2 u = f, \quad \text{ u } \Omega, gde je član koji specificira izvor :math:`f(x,y)=k_0^2 \sin(k_0x) \sin(k_0 y)`. Oblik domena može se videti na :numref:`akustika-rupa`. .. _akustika-rupa: .. figure:: akustika2.png :width: 60% Postavka problema i granični uslovi Postoji analitičko rešenje ovog problema i ono glasi: .. math:: u(x,y) = \sin(k_0 x) \sin(k_0 y). Ovde ćemo specificirati Dirihleove granične uslove prema analitičkom rešenju na spoljnoj granici domena :math:`\Omega` koju označavamo sa :math:`\Gamma_{outer}`: .. math:: u(x,y) \left| \right. \Gamma_{outer} = \sin(k_0 x) \sin(k_0 y), \quad (x,y) \in \Gamma_{outer} Na sličan način, možemo da definišemo granični uslov na unutrašnjoj granici, ovog puta Nojmanov: .. math:: :label: eq:njuman (\nabla u \left| \right. \Gamma_{inner} \cdot n) (x,y) = \left[ k_0 \cos(k_0 x) \sin(k_0 y), k_0 \sin(k_0 x) \cos(k_0 y) \right] \cdot n = g(x,y), gde je :math:`n` vektor normale. Konciznije napisano, Nojmanov granični uslov na unutrašnjoj granici glasi: .. math:: \nabla u \left| \right. \Gamma_{inner} \cdot n = g. Implementacija ---------------- Na sledećem listingu su dati glavni detalji implementacije. Namerno su izostavljeni delovi koji su irelevantni za samo rešavanje, kao što je crtanje dijagrama. Celokupna skripta se, kao i ostale, nalazi u repozitorijumu sa primerima. .. code-block:: python :caption: Rešenje problema prostiranja stojećeg talasa u 2D domenu sa šupljinom :linenos: import deepxde as dde import matplotlib.pyplot as plt import numpy as np from deepxde.backend import tf sin = tf.sin # Opsti parametri n = 2 length = 1 R = 1 / 4 precision_train = 15 precision_test = 30 weight_inner = 10 weight_outer = 100 iterations = 5000 learning_rate = 1e-3 num_dense_layers = 3 num_dense_nodes = 350 activation = "sin" k0 = 2 * np.pi * n wave_len = 1 / n # Parcijalna diferencijalna jednacina def pde(x, y): dy_xx = dde.grad.hessian(y, x, i=0, j=0) dy_yy = dde.grad.hessian(y, x, i=1, j=1) f = k0**2 * sin(k0 * x[:, 0:1]) * sin(k0 * x[:, 1:2]) return -dy_xx - dy_yy - k0**2 * y - f # Egzaktno resenje def func(x): return np.sin(k0 * x[:, 0:1]) * np.sin(k0 * x[:, 1:2]) # Da li je tacka na granici? def boundary(_, on_boundary): return on_boundary # Njumanovi granicni uslovi prema egzaktnom resenju def neumann(x): grad = np.array([ k0 * np.cos(k0 * x[:, 0:1]) * np.sin(k0 * x[:, 1:2]), k0 * np.sin(k0 * x[:, 0:1]) * np.cos(k0 * x[:, 1:2]),]) normal = -inner.boundary_normal(x) normal = np.array([normal]).T result = np.sum(grad * normal, axis=0) return result # Geometrija outer = dde.geometry.Rectangle([-length / 2, -length / 2], [length / 2, length / 2]) inner = dde.geometry.Disk([0, 0], R) # Da li je tacka na spoljnoj granici? def boundary_outer(x, on_boundary): return on_boundary and outer.on_boundary(x) # Da li je tacka na unutrasnjoj granici? def boundary_inner(x, on_boundary): return on_boundary and inner.on_boundary(x) # Iskljuci krug iz kvadrata geom = outer - inner hx_train = wave_len / precision_train nx_train = int(1 / hx_train) hx_test = wave_len / precision_test nx_test = int(1 / hx_test) # Na unutrasnjoj granici Njuman, na spoljnoj Dirihleovi bc_inner = dde.icbc.NeumannBC(geom, neumann, boundary_inner) bc_outer = dde.icbc.DirichletBC(geom, func, boundary_outer) data = dde.data.PDE( geom, pde, [bc_inner, bc_outer], num_domain=nx_train**2, num_boundary=16 * nx_train, solution=func, num_test=nx_test**2, ) net = dde.nn.FNN( [2] + [num_dense_nodes] * num_dense_layers + [1], activation, "Glorot uniform" ) model = dde.Model(data, net) loss_weights = [1, weight_inner, weight_outer] model.compile("adam", lr=learning_rate, metrics=["l2 relative error"], loss_weights=loss_weights) losshistory, train_state = model.train(iterations=iterations) Koristićemo *Tensorflow* kao *backend* u svim našim primerima, ali treba imati u vidu da okvir *DeepXDE* podržava i *PyTorch* i još neke. Nakon standardne specifikacije opštih parametara i hiper-parametara, kao u primeru :ref:`akustika_ravan`, uz jedinu modifikaciju dodavanja nešto više neurona po sloju (350), definišemo Nojmanov granični uslov prema jednačini :math:numref:`eq:njuman`: .. code-block:: python def neumann(x): grad = np.array([ k0 * np.cos(k0 * x[:, 0:1]) * np.sin(k0 * x[:, 1:2]), k0 * np.sin(k0 * x[:, 0:1]) * np.cos(k0 * x[:, 1:2]),]) normal = -inner.boundary_normal(x) normal = np.array([normal]).T result = np.sum(grad * normal, axis=0) return result Kao što se vidi iz koda, postoje uslužne funkcije koje računaju normale na pravilne granice u kolokacionim tačkama. Ponderske težine graničnih uslova u obuci ``weight_inner`` i ``weight_outer`` takođe spadaju u neku vrstu hiper-parametara, pa i njima treba posvetiti pažnju uz nekoliko manuelnih proba. Dalje, sledi specifikacija geometrije problema kao razlike kvadrata i diska: .. code-block:: python outer = dde.geometry.Rectangle([-length / 2, -length / 2], [length / 2, length / 2]) inner = dde.geometry.Disk([0, 0], R) geom = outer - inner Ostatak skripte je sličan primeru bez šupljine :ref:`akustika_ravan`, pa ga nećemo dodatno pojašnjavati. Dovoljno je reći da pažnju treba obratiti da bude dovoljno kolokacionih tačaka na spoljnoj i na unutrašnjoj granici. Rezultati ------------ Dobijeni rezultati su prikazani u formi konturnog grafika na :numref:`rezultati-rupa`. Oko unutrašnje granice prikazani su pravci vektora normala. .. _rezultati-rupa: .. figure:: rezultati-rupa.png :width: 90% Rezultati primera kvadratnog domena sa šupljinom Mera relativne greške modela iznosi 0,048. Uz obraćanje posebne pažnje na forsiranje graničnih uslova, zatim arhitekturu NMPFZ i najzad tip aktivacione funkcije, uspeli smo da dobijemo prilično dobro rešenje. Čitalac može samostalno da proba kako bi promena frekvencije (a samim tim i talasne dužine), gustine kolokacionih tačaka, arhitekture, uticala na proces obuke modela.