This is my code ，the results are incorrect; they seem to be unaffected by the terrain. `device = torch.device("cuda:0" if torch.cuda.is_available() else "cpu") g = 9.81 xmin, xmax, tmax = 0., 30., 40. input_dim = 2 output_dim = 2

def Slope(x): if isinstance(x, np.ndarray): # 如果输入是 numpy 数组，转换为 torch.Tensor x = torch.tensor(x, dtype=torch.float32) y = torch.where(x < 7, 0., torch.where(x < 13, 0.05 (x - 7), torch.where(x < 15, 0.3, torch.where(x < 18, 0.3 - 0.1 (x - 15), torch.where(x < 20, 0., 0.04 * (x - 20)))))) return y

Slope_x 函数：支持 PyTorch 张量

def Slope_x(x): if isinstance(x, np.ndarray): # 如果输入是 numpy 数组，转换为 torch.Tensor x = torch.tensor(x, dtype=torch.float32) y = torch.where(x < 7, 0., torch.where(x < 13, 0.05, torch.where(x < 15, 0., torch.where(x < 18, -0.1, torch.where(x < 20, 0., 0.04))))) return y

still water height

def sth(x): if isinstance(x, np.ndarray): # 如果输入是 numpy 数组，转换为 torch.Tensor x = torch.tensor(x, dtype=torch.float32) y = torch.where(x < 7, 0.4, torch.where(x < 13, 0.4-0.05 (x - 7), torch.where(x < 15, 0.1, torch.where(x < 18, 0.1 + 0.1 (x - 15), torch.where(x < 20, 0.4, 0.4-0.04 * (x - 20)))))) return y

Define bottom topography

def get_bottom(bottom_name): assert (isinstance(bottom_name, str))

bottom5 = lambda x: Slope(x)
bottom5_x = lambda x: Slope_x(x)

if bottom_name == 'b5':
    return {'b': bottom5, 'b_x': bottom5_x}
else:
    raise ValueError('Invalid bottom name, {} bottom not implemented'.format(bottom_name))

def _swe(bottom_name, g, source_term=True): def swe(X, U): bottom = get_bottom(bottom_name) bx = bottom['b_x'] x, t = X[:, 0:1], X[:, 1:2] h = U[:, 0:1] u = U[:, 1:2]

    U1 = h
    U2 = h * u

    F1 = h * u
    F2 = h * u * u + 0.5 * g * h * h

    F1_x = dde.grad.jacobian(F1, X, i=0, j=0)
    F2_x = dde.grad.jacobian(F2, X, i=0, j=0)
    U1_t = dde.grad.jacobian(U1, X, i=0, j=1)
    U2_t = dde.grad.jacobian(U2, X, i=0, j=1)

    if source_term:
        h = 2 + torch.cos(x) * torch.cos(t)
        v = torch.sin(x) * torch.sin(t) / h
        S = torch.sin(x) * torch.cos(t) * (1 + v ** 2 - g * h) + 2 * v * torch.cos(x) * torch.sin(t) + g * h * bx(x)
    else:
        S = torch.zeros_like(x)
    b_x = bx(x)
    equaz_1 = U1_t + F1_x
    equaz_2 = U2_t + F2_x
    #+ g * h * b_x - S
    return [equaz_1, equaz_2]

return swe

def oninitial(, on_initial): return on_initial

def boundary(_, on_boundary): return on_boundary

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

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

def func_ic_h(X): if isinstance(X, np.ndarray): # 检查输入是否为 numpy.ndarray X = torch.tensor(X, dtype=torch.float32) x = X[:, 0:1] h_ic = 0.4 + 0.1 torch.sin(2 torch.pi / 2.02 * x) return h_ic

def _transform_output(U0): def transform_output(x, U): t = x[:, 1:2] h = U[:, 0:1] u = U[:, 1:2] u_new = u * t + U0(x) h_new = h return torch.cat([h_new,u_new], axis=1) return transform_output

def func_ic_u(x): return 0.0

def func_bc_h1(X): T = 2.02 # 周期 (秒) A = 0.01 # 振幅 (cm) if isinstance(X, np.ndarray): # 检查输入是否为 numpy.ndarray X = torch.tensor(X, dtype=torch.float32) t = X[:, 1:2]

x_l = torch.ones_like(t) * xmin

h_lbc = 0.4 + A * torch.sin(2 * torch.pi / T * t)
return h_lbc

def func_bc_h2(X): if isinstance(X, np.ndarray): # 检查输入是否为 numpy.ndarray X = torch.tensor(X, dtype=torch.float32) t = X[:, 1:2] H2 = 0.43 + 0.01 torch.sin(2 torch.pi / 2.02 * t) return 0.42

def func_bc_u1(X): if isinstance(X, np.ndarray): # 检查输入是否为 numpy.ndarray X = torch.tensor(X, dtype=torch.float32) t = X[:, 1:2] v = 0.0495torch.cos(3.1105t+90) return v

def func_bc_h2(X):

if isinstance(X, np.ndarray): # 检查输入是否为 numpy.ndarray

X = torch.tensor(X, dtype=torch.float32)

H = ref_value[30, :]

print(H)

return H

model part

geom = dde.geometry.Interval(xmin, xmax) timedomain = dde.geometry.TimeDomain(0.0, tmax) geomtime = dde.geometry.GeometryXTime(geom, timedomain)

IC_h = dde.IC(geomtime, func_ic_h, on_initial, component=0)

IC_u = dde.IC(geomtime, func_ic_u, on_initial, component=1)

BC_h1 = dde.DirichletBC(geomtime, func_bc_h1, boundary_l, component=0) BC_h2 = dde.DirichletBC(geomtime, func_bc_h2, boundary_r, component=0)

BC_h3 = dde.icbc.PointSetBC(12.5, ref_value[:, 37])

BC_u1 = dde.DirichletBC(geomtime, func_bc_u1, boundary_l, component=1)

BC = [IC_h, IC_u, BC_h1, BC_h2, BC_u1] swe = _swe('b5', g, source_term=False) data = dde.data.TimePDE(geomtime, swe, BC, num_domain=100, num_boundary=2, num_initial=100, num_test=100) net = dde.maps.FNN(layer_sizes=[input_dim] + [60] * 5 + [output_dim], activation="tanh", kernel_initializer="Glorot normal") model = dde.Model(data, net) resampler = dde.callbacks.PDEPointResampler(period=100) model.compile('adam', lr=0.0005, decay = ("inverse time", 1000, 0.3), loss_weights=0.5) transform_output = _transform_output(func_ic_u) net.apply_output_transform(lambda x, U: transform_output(x, U)) losshistory, train_state = model.train(epochs=20_000)`

lululxvi / deepxde

1D shallow water #1844