PINN to obtain temperature distribution using Navier-Stokes equations

[073kukhta]

Dear LuLu,
At first, thanks for the great library.

I am solving Navier-Stokes equations with DeepXDE to get temperature distribution. Input values: x, y, Pr (Prandtl numbers as Z-coordinate), output values: T (Temperature). 'x' and 'y' in my case mean 'y' and 'z' respectively.

Where mean u'y * c' approximated as (nut/Prt * dT/dy) and mean u'z * c' approximated as (nut/Prt * dT/dz).
With boundary conditions:

The problem is that there should be different temperature distributions for different Prandtl numbers. But for Pr = 0.025 and Pr = 1.0 I get similar linear distributions. The distributions that should be (upper for Pr = 0.025, lower for Pr = 1):

But I get something like that:


Here is my code (constant values ​​loaded from the files above):
` from deepxde.backend import tf
import matplotlib.pyplot as plt
import numpy as np
from fluidfoam import *
from deepxde.backend import tf
import deepxde as dde

def pde(x, y):
dnut_dx = tf.convert_to_tensor(dnut[1])
dnut_dy = tf.convert_to_tensor(dnut[2])
nu = tf.convert_to_tensor([9.6899224806201e-05]*9216)
Pr = x[:,2:3]
#Pr = tf.convert_to_tensor([1.0]*9216)
Prt = tf.convert_to_tensor([0.85]*9216)
Ux = U[1]
Uy = U[2]

dy = tf.gradients(y, x)[0]
dy_x = dy[:,0:1]
dy_y = dy[:,1:2]
#dy_x = dde.grad.jacobian(y, x, i = 0, j = 0)
#dy_y = dde.grad.jacobian(y, x, i = 0, j = 1)

ddy = tf.gradients(dy, x)[0]
dy_xx = ddy[:,0:1]
dy_yy = ddy[:,1:2]

#dy_xx = dde.grad.hessian(y, x, i=0, j=0)
#dy_yy = dde.grad.hessian(y, x, i=1, j=1)

return - Ux*dy_x - Uy*dy_y + (nu/Pr + nut/Prt)*(dy_xx + dy_yy) + (1/Prt)*(dy_x*dnut_dx + dy_y*dnut_dy)  

def boundary(_, on_boundary):
return on_boundary

def boundary_bottom(x, on_boundary):
return on_boundary and dde.utils.isclose(x[1], 0)

def boundary_top(x, on_boundary):
return on_boundary and dde.utils.isclose(x[1], 1)

def boundary_left(x, on_boundary):
return on_boundary and dde.utils.isclose(x[0], 0)

def boundary_right(x, on_boundary):
return on_boundary and dde.utils.isclose(x[0], 1)

def func(x):
return 1 - x[:,1:2]

#geom = dde.geometry.Rectangle([0, 0], [1, 1])
geom = dde.geometry.Cuboid([0, 0, 0.025], [1, 1, 1])
bc = dde.icbc.DirichletBC(geom, func, boundary)
bc_bot = dde.icbc.DirichletBC(geom, lambda x:1, boundary_bottom)
bc_top = dde.icbc.DirichletBC(geom, lambda x:0, boundary_top)
bc_left = dde.icbc.DirichletBC(geom, func, boundary_left)
bc_right = dde.icbc.DirichletBC(geom, func, boundary_right)

data = dde.data.PDE(geom, pde, [bc_top, bc_bot], num_domain=1500, num_boundary=500, num_test=1500) #, anchors = X_train)
net = dde.nn.FNN([3] + [128] * 6 + [1], "tanh", "Glorot uniform")
model = dde.Model(data, net)

model.compile("adam", lr=0.0001)#, loss_weights=[1, e+4])

losshistory, train_state = model.train(iterations=5000)
#model.compile("L-BFGS", lr=0.001)
#losshistory, train_state = model.train(iterations=10000)
dde.saveplot(losshistory, train_state, issave=True, isplot=True)`

Any idea how to fix it?
Thanks in advance for your answer.

With best regards,
G.K.