question on combining observational information and physics in an ode

[bollegijscoding]

Hi,

Short question:

When considering both observational data and applying information from an ODE with dde.data.PDE, does num_domain relate to points which add to the loss by considering also the solution at these points? Or by considering only the loss related to the ODE?

Question elaborated with example:

I have an ODE which I would like to solve with deepxde considering observations occurring in part of the domain. As a toy problem, I am trying to solve the damped harmonic oscillator case as described nicely by Ben Moseley (e.g. https://benmoseley.blog/blog/ and https://github.com/benmoseley/harmonic-oscillator-pinn).

Hence, I am considering

$\frac{d^{2} y}{dt^{2}} + a \frac{dy}{dt} + b y = 0$, on the domain $y \in [0, 1.25]$ with boundary conditions $y(t=0) = 1$ and $\frac{dy}{dt} (t=0) = 0$. The variables $a$ and $b$ are set to 4 and 400 respectively.

This problem has a closed form solution which is given by (see https://benmoseley.blog/blog/ and references therein):

def func_osc(t):
    d, w0 = a /2 , np.sqrt(b)
    w = np.sqrt(w0**2-d**2)
    phi = np.arctan(-d/w)
    A = 1/(2*np.cos(phi))
    cos = np.cos(phi+w*t)
    sin = np.sin(phi+w*t)
    exp = np.exp(-d*t)
    y  = exp*2*A*cos
    return y

I model the ODE and boundary conditions as follows:

def ode_osc(t, y):
    a = 4
    b = 400
    dy_dt = dde.grad.jacobian(y, t)
    d2y_dt2 = dde.grad.hessian(y, t)
    return d2y_dt2 + a * dy_dt + b * y 

# y(0) = 1
ic1 = dde.icbc.IC(geom, lambda y: 1, lambda _, on_initial: on_initial)

# dy/dt (0) = 0
def boundary_l(t, on_initial):
    return on_initial and np.isclose(t[0], 0)

def bc_func2(inputs, outputs, X):
    return dde.grad.jacobian(outputs, inputs, i=0, j=None)

ic2 = dde.icbc.OperatorBC(geom, bc_func2, boundary_l)

I want to consider observational points within the region of the domain between 0 and 0.35. Hence, I introduce an numpy array to mark these points

time2evaluate = np.array(np.linspace(0, 0.35, 15)).reshape(-1,1)

The neural net is just a FNN and we define the dde.data as follows:

def run_test(
    num_test, loss_weights_0, 
    loss_weights_1, loss_weights_2,  num_domain, lr
):
  geom = dde.geometry.TimeDomain(0, 1.25)
  data = dde.data.PDE(
        geometry=geom, 
        pde=ode_osc, 
        bcs=[ic1, ic2], 
        num_domain=num_domain, 
        num_boundary=2, 
        solution=func_osc, 
        anchors=time2evaluate, 
        num_test=num_test
    )
  
    layer_size = [1] + [32] * 3 + [1]
    activation = "tanh"
    initializer = "Glorot normal"
    net = dde.nn.FNN(layer_size, activation, initializer)
  
    model = dde.Model(data, net)
  
    model.compile(
        "adam", 
        lr=lr, 
        metrics=["mean squared error"], 
        loss_weights=[loss_weights_0, loss_weights_1, loss_weights_2] 
    )
    
    losshistory, train_state = model.train(epochs=20000)
    return losshistory, train_state, dde

I have played around a lot with of the values of the parameters num_test, loss_weights_0, loss_weights_1, loss_weights_2, num_domain, lr; yet the best I have achieved uptil now considers

num_test=25000, 
loss_weights_0=0.0001, 
loss_weights_1=1, 
loss_weights_2=1, 
num_domain=500, 
lr=0.001, 

I only get reasonable results, when I consider a reasonable large number of num_domain (as above 500).

Setting num_domain to a small value (e.g. 10) leads to less appropriate results, see figure below. I am confused however, on the role of num_domain with respect to the loss function. In this setting, does num_domain add to the loss by considering also the solution at these points? Or by considering only the loss related to the ODE?

Many thanks in advance!

Bollegijs

[lululxvi]

When considering both observational data and applying information from an ODE with dde.data.PDE, does num_domain relate to points which add to the loss by considering also the solution at these points? Or by considering only the loss related to the ODE?

Only PDE loss.