NVIDIA Modulus Sym v1.3.0
Sym v1.3.0

Moving Time Window: Taylor Green Vortex Decay

Some of the previous tutorials have shown solving the transient problems with continuous time approach. This tutorial presents a moving time window approach to solve a complex transient Navier-Stokes problem. In this tutorial, you will learn:

  1. How to solve sequences of problems/domains in Modulus Sym.

  2. How to set up periodic boundary conditions.

Note

This tutorial assumes you have completed tutorial on 1D Wave Equation.

As mentioned in tutorial 1D Wave Equation, solving transient simulations with only the continuous time method can be difficult for long time durations. This tutorial shows how this can be overcome by using a moving time window. The example problem is the 3D Taylor-Green vortex decay at Reynolds number 500. The Taylor-Green vortex problem is often used as a benchmark to compare solvers and in this case validation data is generated using a spectral solver. The domain is a cube of length \(2\pi\) with periodic boundary conditions on all sides. The initial conditions you will use are,

(188)\[\begin{split}\begin{aligned} u(x,y,z,0) &= sin(x)cos(y)cos(z) \\ v(x,y,z,0) &= -cos(x)sin(y)cos(z) \\ w(x,y,z,0) &= 0 \\ p(x,y,z,0) &= \frac{1}{16} (cos(2x) + cos(2y))(cos(2z) + 2) \\\end{aligned}\end{split}\]

This problem shows solving the time dependent incompressible Navier-Stokes equations with a density \(1\) and viscosity \(0.002\). Note that because there are periodic boundaries on all sides you do not need boundary conditions Fig. 114.

taylor_green_initial_conditions.png

Fig. 114 Taylor-Green vortex initial conditions.

The moving time window approach works by iteratively solving for small time windows to progress the simulation forward. The time windows use the previous window as new initial conditions. The continuous time method is used for solving inside a particular window. A figure of this method can be found here Fig. 115 for a hypothetical 1D problem. Learning rate decay is restarted after each time window. A similar approach can be found here 1.

moving_time_window_fig.png

Fig. 115 Moving Time Window Method

The case setup for this problem is similar to many of the previous tutorials except for two key differences. This example shows how to set up a sequence of domains to iteratively solve for.

Note

The python script for this problem can be found at examples/taylor_green/.

Sequence of Train Domains

First, construct your geometry similar to the previous problems.

Also, define values for how large the time window will be. In this case solve to \(1\) unit of time and then construct all needed nodes.

Copy
Copied!
            

# make geometry for problem channel_length = (0.0, 2 * np.pi) channel_width = (0.0, 2 * np.pi) channel_height = (0.0, 2 * np.pi) box_bounds = {x: channel_length, y: channel_width, z: channel_height} # define geometry rec = Box( (channel_length[0], channel_width[0], channel_height[0]), (channel_length[1], channel_width[1], channel_height[1]), )

The architecture can be created as below.

Copy
Copied!
            

# make network for current step and previous step flow_net = FullyConnectedArch( input_keys=[Key("x"), Key("y"), Key("z"), Key("t")], output_keys=[Key("u"), Key("v"), Key("w"), Key("p")], periodicity={"x": channel_length, "y": channel_width, "z": channel_height}, layer_size=256, ) time_window_net = MovingTimeWindowArch(flow_net, time_window_size) # make nodes to unroll graph on nodes = ns.make_nodes() + [time_window_net.make_node(name="time_window_network")]

Note that the periodicity is set when creating out network architecture. This will force the network by construction to give periodic solutions on the given bounds. Now two domains are defined, one for the initial conditions and one to be used for all future time windows.

Copy
Copied!
            

# make initial condition domain ic_domain = Domain("initial_conditions") # make moving window domain window_domain = Domain("window") # make initial condition ic = PointwiseInteriorConstraint( nodes=nodes, geometry=rec, outvar={ "u": sin(x) * cos(y) * cos(z), "v": -cos(x) * sin(y) * cos(z), "w": 0, "p": 1.0 / 16 * (cos(2 * x) + cos(2 * y)) * (cos(2 * z) + 2), }, batch_size=cfg.batch_size.initial_condition, bounds=box_bounds, lambda_weighting={"u": 100, "v": 100, "w": 100, "p": 100}, parameterization={t_symbol: 0}, ) ic_domain.add_constraint(ic, name="ic") # make constraint for matching previous windows initial condition ic = PointwiseInteriorConstraint( nodes=nodes, geometry=rec, outvar={"u_prev_step_diff": 0, "v_prev_step_diff": 0, "w_prev_step_diff": 0}, batch_size=cfg.batch_size.interior, bounds=box_bounds, lambda_weighting={ "u_prev_step_diff": 100, "v_prev_step_diff": 100, "w_prev_step_diff": 100, }, parameterization={t_symbol: 0}, ) window_domain.add_constraint(ic, name="ic") # make interior constraint interior = PointwiseInteriorConstraint( nodes=nodes, geometry=rec, outvar={"continuity": 0, "momentum_x": 0, "momentum_y": 0, "momentum_z": 0}, bounds=box_bounds, batch_size=4094, parameterization=time_range, ) ic_domain.add_constraint(interior, name="interior") window_domain.add_constraint(interior, name="interior") # add inference data for time slices for i, specific_time in enumerate(np.linspace(0, time_window_size, 10)): vtk_obj = VTKUniformGrid( bounds=[(0, 2 * np.pi), (0, 2 * np.pi), (0, 2 * np.pi)], npoints=[128, 128, 128], export_map={"u": ["u", "v", "w"], "p": ["p"]}, ) grid_inference = PointVTKInferencer( vtk_obj=vtk_obj, nodes=nodes, input_vtk_map={"x": "x", "y": "y", "z": "z"}, output_names=["u", "v", "w", "p"], requires_grad=False, invar={"t": np.full([128**3, 1], specific_time)}, batch_size=100000, ) ic_domain.add_inferencer(grid_inference, name="time_slice_" + str(i).zfill(4)) window_domain.add_inferencer( grid_inference, name="time_slice_" + str(i).zfill(4) )

In the moving time window domain there is a new initial condition that will come from the previous time window. Now that the domain files have been constructed for the initial conditions and future time windows You can put it all together and make your sequential solver.

Sequence Solver

Copy
Copied!
            

# make solver slv = SequentialSolver( cfg, [(1, ic_domain), (nr_time_windows, window_domain)], custom_update_operation=time_window_net.move_window, ) # start solver slv.solve()

Instead of using the normal Solver class use the SequentialSolver class. This class expects a list where each element is a tuple of the number of times solving as well as the given domain. In this case, solve the initial condition domain once and then solve the time window domain multiple times. A custom_update_operation is provided which is called after solving each iteration and in this case is used to update network parameters.

The iterative domains are solved here in a very general way. In subsequent chapters this structure is used to implement a complex iterative algorithm to solve conjugate heat transfer problems.

After solving this problem you can visualize the results. Look in the network directory to see,

current_step.txt initial_conditions window_0000 window_0001...

If you look in any of these directories you see the typical files storing network checkpoint and results. Plotting at time 15 gives the snapshot Fig. 114. To validate the simulation, the results are compared against a spectral solver’s results. Comparing the point-wise error of transient simulations can be misleading as this is a chaotic dynamical system and any small errors will quickly cause large differences. Instead, you can look at the average turbulent kinetic energy decay (TKE) of the simulation. A figure of this can be found here below.

taylor_green_re_500.png

Fig. 116 Taylor-Green vortex at time \(15.0\).

tke_plot.png

Fig. 117 Taylor-Green Turbulent kinetic energy decay.

References

[1]

Colby L Wight and Jia Zhao. Solving allen-cahn and cahn-hilliard equations using the adaptive physics informed neural networks. arXiv preprint arXiv:2007.04542, 2020.

Previous STL Geometry: Blood Flow in Intracranial Aneurysm
Next Electromagnetics: Frequency Domain Maxwell’s Equation
© Copyright 2023, NVIDIA Modulus Team. Last updated on Jan 25, 2024.