Adding a constraint for speed of sound

[cpe.sk]

Hello,

I’m trying to simulate waves in a domain with an obstacle. My model works fine without adding the obstacle constraint, however, when I’m trying to add this “c” (speed of soundwave) map for the speed of sound as a constraint in my simulation,

I exported the values of x, y, and c_output to an NPZ file from my external script, and I’m trying to import them into my modulus code as such:

# Read in the C (speed of sound) map generated using finite difference simulator
    obstacle_filename = to_absolute_path(f"FDM_data/obstacle.npz")
    c_map = np.load(obstacle_filename)
    t = np.array(c_map["t_input"]).flatten()[:, None].astype(np.float32)
    y = np.array(c_map["y_input"]).flatten()[:, None].astype(np.float32)
    x = np.array(c_map["x_input"]).flatten()[:, None].astype(np.float32)
    c_obstacle = np.array(c_map["c_output"]).astype(np.float32)
    mesh_x, mesh_y, mesh_t = np.meshgrid(
        x,y,t, indexing="ij"
    )
    wave_speed_invar = {}
    wave_speed_invar["x"] = np.expand_dims(mesh_x.flatten(), axis=-1)
    wave_speed_invar["y"] = np.expand_dims(mesh_y.flatten(), axis=-1)
    wave_speed_invar["t"] = np.expand_dims(mesh_t.flatten(), axis=-1)
    
    wave_speed_outvar = {}
    wave_speed_outvar["c"] = c_obstacle
            
    # add speed of sound constraint
    c_speed= PointwiseConstraint.from_numpy(
        nodes=nodes, invar=wave_speed_invar, outvar=wave_speed_outvar, batch_size=1024
    )
    domain.add_constraint(c_speed, "C_Speed")

I receive this error:

Error executing job with overrides: []
Traceback (most recent call last):
  File "wave_2d_inhomo_RVE.py", line 454, in run
    slv.solve()
  File "C:\ProgramData\Anaconda3\envs\modulus_env\lib\site-packages\modulus-22.9-py3.8.egg\modulus\solver\solver.py", line 159, in solve
    self._train_loop(sigterm_handler)
  File "C:\ProgramData\Anaconda3\envs\modulus_env\lib\site-packages\modulus-22.9-py3.8.egg\modulus\trainer.py", line 519, in _train_loop
    self.load_data(static=True)
  File "C:\ProgramData\Anaconda3\envs\modulus_env\lib\site-packages\modulus-22.9-py3.8.egg\modulus\solver\solver.py", line 61, in load_data
    self.domain.load_data(static)
  File "C:\ProgramData\Anaconda3\envs\modulus_env\lib\site-packages\modulus-22.9-py3.8.egg\modulus\domain\domain.py", line 122, in load_data
    constraint.load_data_static()
  File "C:\ProgramData\Anaconda3\envs\modulus_env\lib\site-packages\modulus-22.9-py3.8.egg\modulus\domain\constraint\continuous.py", line 94, in load_data_static
    self.load_data()
  File "C:\ProgramData\Anaconda3\envs\modulus_env\lib\site-packages\modulus-22.9-py3.8.egg\modulus\domain\constraint\continuous.py", line 81, in load_data
    invar, true_outvar, lambda_weighting = next(self.dataloader)
  File "C:\ProgramData\Anaconda3\envs\modulus_env\lib\site-packages\modulus-22.9-py3.8.egg\modulus\domain\constraint\constraint.py", line 238, in __iter__
    for batch in dataloader:
  File "C:\ProgramData\Anaconda3\envs\modulus_env\lib\site-packages\torch\utils\data\dataloader.py", line 628, in __next__
    data = self._next_data()
  File "C:\ProgramData\Anaconda3\envs\modulus_env\lib\site-packages\torch\utils\data\dataloader.py", line 671, in _next_data
    data = self._dataset_fetcher.fetch(index)  # may raise StopIteration
  File "C:\ProgramData\Anaconda3\envs\modulus_env\lib\site-packages\torch\utils\data\_utils\fetch.py", line 60, in fetch
    data = self.dataset[possibly_batched_index]
  File "C:\ProgramData\Anaconda3\envs\modulus_env\lib\site-packages\modulus-22.9-py3.8.egg\modulus\dataset\continuous.py", line 39, in __getitem__
    outvar = _DictDatasetMixin._idx_var(self.outvar, idx)
  File "C:\ProgramData\Anaconda3\envs\modulus_env\lib\site-packages\modulus-22.9-py3.8.egg\modulus\dataset\dataset.py", line 117, in _idx_var
    idx_var[key] = value[idx]
IndexError: index 553925 is out of bounds for dimension 0 with size 100

For clarification, x is a dimension from 0 to 5 with 100 steps and y, t are the same. Also, the values of speed don’t change over time so the values I’m trying to import are constant throughout the time steps.

How can I proceed with training my model while having these constant values of c (speed of soundwave) in my domain? Please advise.

[ngeneva]

Hi @cpe.sk

From a first pass I think your code looks good. I dont see anything obvious. So let me just confirm some items here. First without this line of code:

c_speed= PointwiseConstraint.from_numpy(
        nodes=nodes, invar=wave_speed_invar, outvar=wave_speed_outvar, batch_size=1024
    )
    domain.add_constraint(c_speed, "C_Speed")

your training file runs fine. Please confirm this.

The second thing I would verify is that all you numpy arrays are of the same dimension [N,1]:

 wave_speed_invar["x"].shape
 wave_speed_invar["y"].shape
 wave_speed_invar["t"].shape
wave_speed_outvar["c"].shape

Let me know, and we can hopefully dig in more from here.

[cpe.sk]

Yes, it works fine without adding c_speed.

I made some changes to this block:

    # Read in the C (speed of sound)/obstacle map generated using finite difference simulator
    obstacle_filename = to_absolute_path(f"FDM_data/obstacle.npz")
    c_map = np.load(obstacle_filename)
    t_obs = np.array(c_map["t_input"]).flatten()[:, None].astype(np.float32)
    y_obs = np.array(c_map["y_input"]).flatten()[:, None].astype(np.float32)
    x_obs = np.array(c_map["x_input"]).flatten()[:, None].astype(np.float32)
    c_obs = np.array(c_map["c_output"]).astype(np.float32)
    mesh_x_obs, mesh_y_obs, mesh_t_obs = np.meshgrid(
        x_obs, y_obs, t_obs, indexing="ij"
    )
    wave_speed_invar = {}
    # wave_speed_invar["x"] = np.expand_dims(mesh_x.flatten(), axis=-1)
    # wave_speed_invar["y"] = np.expand_dims(mesh_y.flatten(), axis=-1)
    # wave_speed_invar["t"] = np.expand_dims(mesh_t.flatten(), axis=-1)
    wave_speed_invar["x"] = x_obs
    wave_speed_invar["y"] = y_obs
    wave_speed_invar["t"] = t_obs
    wave_speed_outvar = {}
    wave_speed_outvar["c"] = c_obs
            
    print ("x", wave_speed_invar["x"].shape)
    print ("y",wave_speed_invar["y"].shape)
    print ("t",wave_speed_invar["t"].shape)
    print ("c",wave_speed_outvar["c"].shape)
    # add speed of sound constraint
    c_speed= PointwiseConstraint.from_numpy(
        nodes=nodes, invar=wave_speed_invar, outvar=wave_speed_outvar, batch_size=1024
    )
    domain.add_constraint(c_speed, "C_Speed")

The error is no longer there, but the code gets stuck at this output:

[07:57:30] - JIT using the NVFuser TorchScript backend
[07:57:30] - JitManager: {'_enabled': True, '_arch_mode': <JitArchMode.ONLY_ACTIVATION: 1>, '_use_nvfuser': True, '_autograd_nodes': False}
[07:57:30] - GraphManager: {'_func_arch': False, '_debug': False, '_func_arch_allow_partial_hessian': True}
x (100, 1)
y (100, 1)
t (100, 1)
c (100, 100)
[07:57:57] - Installed PyTorch version 1.13.0 is not TorchScript supported in Modulus. Version 1.13.0a0+d321be6 is officially supported.
[07:57:57] - attempting to restore from: outputs\wave_2d_inhomo_RVE
[07:57:57] - optimizer checkpoint not found
[07:57:57] - model source_network.0.pth not found
[07:57:57] - model speed_network.0.pth not found
[07:57:57] - model wave_network.0.pth not found

The sizes of the numpy arrays are

x (100, 1)
y (100, 1)
t (100, 1)
c (100, 100)

Any suggestions on how to get around this?
Thanks for your patience and support @ngeneva !

[ngeneva]

Hi @cpe.sk

I believe this issue here is that your x,y,t,c arrays all need to be the same size. For your output it also needs to be of the shape [N,1]. The last dimension is the dimension of the variable at a given point (it should also match the output dimension of your speed_net which is 1). In other words each element of your c array is a independent training point.

This also means your input variables need to expand as well (the mesh grid approach you were using kinda works but remember your c_obs is not 3D [only 100x100]). So given the dimensions here you probably want:

x (100*100, 1)
y (100*100, 1)
t (100*100, 1)
c (100*100, 1)

[cpe.sk]

I removed the time dimension and flattened all the numpy arrays. It’s running now. Thank you @ngeneva ! Hopefully, it will produce the targeted results.