.. _oscilacije_inverzni:

Inverzni problem
==================

Kao što smo već u više navrata naglasili, prava snaga NMPFZ dolazi do izražaja kod problema u kojima je potrebno identifikovati nepoznate parametre. Razlog tome je što NMPFZ metoda inverzne probleme tretira na isti način kao direktne. Važno je samo da se problem postavi ne dozvoljavajući višeznačnost, što ćemo demonstrirati na sledećem primeru. Zamislimo eksperiment postavljen kao na :numref:`opruga` u kome nam je poznat samo koeficijent trenja :math:`\mu=0,6`, a ne znamo ni masu kuglice :math:`m` ni koeficijent elastičnosti opruge :math:`k`. Zadatak nam je da identifikujemo ova dva parametra tako što ćemo pustiti kuglicu da osciluje i na osnovu njenog kretanja dati NMPFZ da odredi elemente koji nedostaju.

Recimo da smo vizuelno utvrdili da je u pitanju **podprigušeni slučaj**. Prva strategija merenja koja bi mogla da bude realno izvodljiva je da kuglicu otklonimo za :math:`x_0=2` i samo pustimo (:math:`v_0=0`), a onda štopericom odredimo vremenske trenutke u kojima kuglica prolazi kroz ravnotežni položaj i da te trenutke zabeležimo:

+--------+-------+-------+-------+-------+-------+
| **t**  | 1,18  | 3,27  | 5,37  | 7,46  | 9,55  |
+--------+-------+-------+-------+-------+-------+
| **x**  | 0     | 0     | 0     | 0     | 0     |
+--------+-------+-------+-------+-------+-------+

Dakle, pored graničnih i početnih uslova :math:`(x_0=2,v_0=0)` uvešćemo i granične uslove tipa ``PointSet`` koji definišu vrednost modelovane funkcije u pojedinim tačkama. Programski kod koji implementira ovaj inverzni problem dat je na sledećem listingu.

.. code-block:: python
    :caption: Inverzni problem prigušenih oscilacija u 1D. Nepoznati parametri su ``m`` i ``k``.
    :linenos:

    import deepxde as dde
    from deepxde.backend import tf
    import numpy as np
    import matplotlib.pyplot as plt

    m = dde.Variable(0.5)
    mu = 0.6
    k = dde.Variable(2.)

    x0, v0 = 2, 0

    delta = mu / (2*m)
    w0 = tf.sqrt(k/m)

    # Da li je tacka blizu t=0 (provera pocetnog uslova)
    def boundary_l(t, on_initial):
        return on_initial and np.isclose(t[0], 0)

    # Jednacina ODE
    def Ode(t,x):
            dxdt = dde.grad.jacobian(x,t)
            dxdtt = dde.grad.hessian(x,t)
            return m * dxdtt + mu * dxdt + k*x
        
    # x(0)=x0
    def bc_func1(inputs, outputs, X):
        return outputs - x0

    # x'(0)=v0
    def bc_func2(inputs, outputs, X):
        return dde.grad.jacobian(outputs, inputs, i=0,j=None) - v0
        
    # Resava se na domenu t=(0,10)
    interval = dde.geometry.TimeDomain(0, 10)

    # Pocetni uslovi
    ic1 = dde.icbc.OperatorBC(interval, bc_func1, boundary_l)
    ic2 = dde.icbc.OperatorBC(interval, bc_func2, boundary_l)

    bc_x = np.array([1.18, 3.27, 5.37, 7.46, 9.55]).reshape(6,1)
    bc_y = np.array([0, 0, 0, 0, 0 ]).reshape(6,1)
    ic3 = dde.icbc.PointSetBC(bc_x, bc_y, component=0)

    # Definsanje problema, granicnih uslova, broja kolokacionih tacaka
    data = dde.data.TimePDE(interval, Ode, [ic1, ic2, ic3], 200, 20, solution=func, num_test=100)
        
    layers = [1] + [30] * 3 + [1]
    activation = "tanh"
    init = "Glorot uniform"
    net = dde.nn.FNN(layers, activation, init)

    model = dde.Model(data, net)

    # Callback funkcija koja stampa varijablu na svakih 1000 epoha
    variable1 = dde.callbacks.VariableValue(k, period=1000)
    variable2 = dde.callbacks.VariableValue(m, period=1000)

    model.compile("adam", lr=.001, loss_weights=[0.01, 1, 1, 1], metrics=["l2 relative error"], external_trainable_variables=[k,m])
    losshistory, train_state = model.train(epochs=50000, callbacks=[variable1, variable2])

Objasnićemo samo delove koji se razlikuju u odnosu na direktni problem rešen u sekciji :ref:`oscilacije_implementacija`. Za početak, tu su nepoznati parametri koje deklarišemo na sledeći način:

.. code-block:: python

    m = dde.Variable(0.5)
    k = dde.Variable(2.)

U zagradi se daju početne vrednosti parametra. Sledeće linije definišu pomenuti dodatni ``PointSet`` granični uslov (uslove) koji važe u pojedinim tačkama *unutar domena*:

.. code-block:: python

    bc_x = np.array([1.18, 3.27, 5.37, 7.46, 9.55]).reshape(6,1)
    bc_y = np.array([0, 0, 0, 0, 0 ]).reshape(6,1)
    ic3 = dde.icbc.PointSetBC(bc_x, bc_y, component=0)

Kako bismo obezbedili praćenje vrednosti nepoznatih parametara tokom obuke, potrebno je da postavimo tzv. ``callback`` funkcije, koje će se pozivati na svakih 1000 epoha:

.. code-block:: python

    variable1 = dde.callbacks.VariableValue(k, period=1000)
    variable2 = dde.callbacks.VariableValue(m, period=1000)

Prilikom postavljanja modela postavljamo eksterne varijable za treniranje ``k`` i ``m``, dok se pri pozivu treninga navode ``callback`` funkcije:

.. code-block:: python

    model.compile("adam", lr=.001, loss_weights=[0.01, 1, 1, 1], metrics=["l2 relative error"], external_trainable_variables=[k,m])
    losshistory, train_state = model.train(epochs=50000, callbacks=[variable1, variable2])

Nakon završenog obučavanja, dobija se, očigledno pogrešno, rešenje koje je prikazano na :numref:`inverzni-pogresno`. U odnosu na analitičko rešenje koje je postavljeno koristeći vrednosti parametara ``m=1`` i ``k=2,25``, dobijene vrednosti ``m=0,287`` i ``k=1,22`` se puno razlikuju.

.. _inverzni-pogresno:

.. figure:: resenje1-inv-pogresno.png
    :width: 80%

    NMPFZ rešenje inverznog problema sa nepoznatim parametrima. Pronađene vrednosti su ``m=0,287`` i ``k=1,22``. 

Zašto smo dobili ovako loše rešenje? Odgovor se krije u loše postavljenim graničnim uslovima koji dovode do nejednoznačnosti inverznog problema. Naime, postoji više scenarija, tj. parova ``(m,k)`` koji zadovoljavaju granične uslove postavljene samo u tačkama prolaska tega kroz ravnotežni položaj. Očigledno je da moramo da dodamo još neku tačku van ravnotežnog položaja, kako bismo obezbedili jednoznačno rešenje. Zamislimo da smo izmerili i trenutak kada je teg bio na najvećoj negativnoj udaljenosti u odnosu na ravnotežni položaj i kolika je ta udaljenost bila. Dodajmo sada i tu tačku u postavljeni ``PointSet``, koji sada izgleda ovako:

+--------+-------+-------+-------+-------+-------+-------+
| **t**  | 1,18  | 3,27  | 5,37  | 7,46  | 9,55  | 2,12  |
+--------+-------+-------+-------+-------+-------+-------+
| **x**  | 0     | 0     | 0     | 0     | 0     |-1,67  |
+--------+-------+-------+-------+-------+-------+-------+

Pogledom na :numref:`inverzni-ispravno` odmah se vidi da je poklapanje sa analitičkim rešenjem u ovako postavljenom problemu skoro pa idealno.

.. _inverzni-ispravno:

.. figure:: resenje1-inv.png
    :width: 80%

    NMPFZ rešenje inverznog problema sa nepoznatim parametrima. Pronađene vrednosti su ``m=0,98`` i ``k=2,26``. 

Ispravnost rešenja dodatno potvrđuju parametri ``m=0,98`` i ``k=2,26``, čije su vrednosti veoma bliske onima koje su date u analitičkoj postavci ``m=1`` i ``k=2,25``. Na ovaj način smo pokazali da pristup rešavanju inverznog problema, iako metodološki sličan, ima specifičnosti o kojima treba voditi računa. Kod direktnih problema rešenje je uvek jednoznačno, dok kod inverznih moraju da se obezbede odgovarajući uslovi koji u dovoljnoj meri determinišu rešenje.