Nishant Kumar
Nishant Kumar

Reputation: 1

Why does my PINN model produce a negative solution?

I'm trying to create a PINN model, for the following Poisson Equation:

Here's the custom loss function:

def customloss(y_t,y_p):  
  with tf.GradientTape(persistent=True) as g:
    with tf.GradientTape(persistent=True) as tape:
      op=model(ip)
    gradi=tape.batch_jacobian(op, ip)
  gradi2=g.batch_jacobian(gradi, ip)
  x=ip[:,0]
  y=ip[:,1]
  uxx=gradi2[:,0,0,0]
  uyy=gradi2[:,0,1,1]
  gradi
  l1=tf.reduce_mean((uxx+uyy+1)**2)
  l2=tf.reduce_mean((op[10000:11000])**2)
  return l1+l2

(I've given 10000:11000 of the input array as points on the boundary: l2 is loss on boundary conditions)

Here's the network:

model=Sequential()
model.add(Dense(10, input_shape=(2,), activation='relu'))
for i in range(3):
  model.add(Dense(50, activation='tanh'))
model.add(Dense(1))

I'm getting something of this sort: enter image description here

Note that most of the solution has negative value over the field.

When I should be getting something along these lines: enter image description here

Why is this? Is the negative solution a valid one, and how do I reach the correct solution?

I have tried deeper networks with 8 and 12 layers, of varying widths. None of these seem to help.

Also tried adding something in the loss function to penalized negative solutions, but it gets messy, and each time the network converges to a different solution.

Upvotes: 0

Views: 277

Answers (1)

L Maxime
L Maxime

Reputation: 124

I dont think the batch_jacobian is the function you want to use.

Here is an example working of the PDE part on TF2 using eager model to help you:

# loss function over the domain that computes the residuals
def lossDom(X,model):
    
    # 'recording' the operations 
    # to unroll the computational graph for gradient calculation
    with tf.GradientTape(persistent=True) as tape:
        x , y = tf.Variable(X[:,0]) , tf.Variable(X[:,1])
        
        # making predictions 
        u = model(tf.stack([x,y],axis=1))
        
        # first derivatives
        u_x = tape.gradient(u,x)
        u_y = tape.gradient(u,y)
    
        # second derivatives
        u_xx = tape.gradient(u_x,x)
        u_yy = tape.gradient(u_y,y)
    
    # values of source term
    fvals = f(X)
    return u_xx + u_yy - fvals  

X represents a tensor of size (N,2) containing the points respecting the pde. The Dirichlet boundary conditions needs to be added. Some will use the boundary points in the pde loss. You can try both approches to see which one gives the best results.

Upvotes: 0

Related Questions