Enclosure natural convection problem using Deep XDE

[turagaramki]

hi all

I am working on validation of enclosure natural convection simulation using dde. The code is attached in this post.

With whatever i have tried the results are no so encouraging so looking for community help.

`import deepxde as dde
import numpy as np

geom = dde.geometry.Rectangle([0, 0], [1, 1])

# Parameters
rho = 1.0  # Density
mu = rho * 0.099  # Dynamic viscosity (mu = rho * nu)
alpha = 0.099  # Thermal diffusivity
g = 9.81  # Gravitational acceleration
beta = 0.1  # Thermal expansion coefficient

# Define the PDE function for steady-state
def pde(x, y):
    u, v, p, T = y[:, 0:1], y[:, 1:2], y[:, 2:3], y[:, 3:4]
    u_x = dde.grad.jacobian(y, x, i=0, j=0)
    u_y = dde.grad.jacobian(y, x, i=0, j=1)
    v_x = dde.grad.jacobian(y, x, i=1, j=0)
    v_y = dde.grad.jacobian(y, x, i=1, j=1)
    p_x = dde.grad.jacobian(y, x, i=2, j=0)
    p_y = dde.grad.jacobian(y, x, i=2, j=1)
    T_x = dde.grad.jacobian(y, x, i=3, j=0)
    T_y = dde.grad.jacobian(y, x, i=3, j=1)

    # Second derivatives
    u_xx = dde.grad.hessian(y, x, component=0, i=0, j=0)
    u_yy = dde.grad.hessian(y, x, component=0, i=1, j=1)
    v_xx = dde.grad.hessian(y, x, component=1, i=0, j=0)
    v_yy = dde.grad.hessian(y, x, component=1, i=1, j=1)
    T_xx = dde.grad.hessian(y, x, component=3, i=0, j=0)
    T_yy = dde.grad.hessian(y, x, component=3, i=1, j=1)

    # Continuity equation
    continuity = u_x + v_y

    # Momentum equations
    momentum_u = rho * (u * u_x + v * u_y) + p_x - mu * (u_xx + u_yy)
    momentum_v = rho * (u * v_x + v * v_y) + p_y - mu * (v_xx + v_yy) - rho * g * beta * T

    # Energy equation
    energy = u * T_x + v * T_y - alpha * (T_xx + T_yy)

    return [continuity, momentum_u, momentum_v, energy]

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

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

def boundary_top_bottom(x, on_boundary):
    return on_boundary and (np.isclose(x[1], 0) or np.isclose(x[1], 1))

bc_T_left = dde.DirichletBC(geom, lambda x: 100, boundary_left, component=0)
bc_T_right = dde.DirichletBC(geom, lambda x: 0, boundary_right, component=0)
bc_u = dde.DirichletBC(geom, lambda x: 0, boundary_top_bottom, component=1)
bc_v = dde.DirichletBC(geom, lambda x: 0, boundary_top_bottom, component=2)

data = dde.data.PDE(geom, pde, [bc_T_left, bc_T_right, bc_u, bc_v], num_domain=2540, num_boundary=80)
net = dde.maps.FNN([2] + [64] * 6 + [4], "tanh", "Glorot uniform")
model = dde.Model(data, net)

model.compile("adam", lr=1e-3)
losshistory, train_state = model.train(iterations=10000, display_every=50)

`

This code is running but not giving the expected result.

Request all the community to guide me what i am missing or possible suggestions to get the results validated

[vl-dud]

Maybe you should do non-dimensionalization?

[turagaramki]

yes i did that too i felt it is nothing to do with the definition of PDE

def pde1(x, y):
    T, u, v, p = y[:, 0:1], y[:, 1:2], y[:, 2:3], y[:, 3:4]
    dT_x = dde.grad.jacobian(y, x, i=0, j=0)
    dT_y = dde.grad.jacobian(y, x, i=0, j=1)
    du_x = dde.grad.jacobian(y, x, i=1, j=0)
    du_y = dde.grad.jacobian(y, x, i=1, j=1)
    dv_x = dde.grad.jacobian(y, x, i=2, j=0)
    dv_y = dde.grad.jacobian(y, x, i=2, j=1)
    dp_x = dde.grad.jacobian(y, x, i=3, j=0)
    dp_y = dde.grad.jacobian(y, x, i=3, j=1)
    d2u_xx = dde.grad.hessian(y, x, component=1, i=0, j=0)
    d2u_yy = dde.grad.hessian(y, x, component=1, i=1, j=1)
    d2v_xx = dde.grad.hessian(y, x, component=2, i=0, j=0)
    d2v_yy = dde.grad.hessian(y, x, component=2, i=1, j=1)
    d2T_xx = dde.grad.hessian(y, x, component=0, i=0, j=0)
    d2T_yy = dde.grad.hessian(y, x, component=0, i=1, j=1)

    # Parameters
    Pr = 0.71  # Prandtl number
    Ra = 1e6  # Rayleigh number

    # Continuity equation
    continuity = du_x + dv_y

    # Momentum equations
    momentum_u = u * du_x + v * du_y + dp_x - (d2u_xx + d2u_yy) / Ra
    momentum_v = u * dv_x + v * dv_y + dp_y - (d2v_xx + d2v_yy) / Ra - Pr * T

    # Energy equation
    energy = u * dT_x + v * dT_y - (d2T_xx + d2T_yy) / (Ra * Pr)

    return [continuity, momentum_u, momentum_v, energy]

I am presume some settings in the DDE such as network layers, and optimizers for this problem might help to get the solution.

[turagaramki]

One more thing what i have observed is whatever the network architecture i am using the loss value is not reducing it is increasing.

I need some guidance here

Thank you in advance

[turagaramki]

hello @lululxvi

Can you please help me with this?

[pescap]

Hi @turagaramki, here are some initial observations:

• Your bc_T_left = dde.DirichletBC(geom, lambda x: 100, boundary_left, component=0) with a value of 100 might be problematic. There’s a significant difference in the order of magnitude compared to the 0 boundary conditions. If I were you, I’d prioritize imposing this condition as a hard constraint.
• Having 8 loss terms is quite a lot and could lead to difficulties during training. I’d recommend using hard boundary conditions for all sides to simplify the problem.

[turagaramki]

Hi @pescap

Thank you very much for the reply. I will try with less Delta T. Just to be on the same page on the hard boundary conditions.

bc_u = dde.DirichletBC(geom, lambda x: 0, boundary_top_bottom, component=1)
bc_v = dde.DirichletBC(geom, lambda x: 0, boundary_top_bottom, component=2)

These two lines can be removed as per your suggestion. Is my understanding correct?

[pescap]

You can check the faq adn inspire yourself of this tutorial.

You have to adapt the output_transform to your specific case.