{'fraction_non_alphanumeric': 0.03725035529957364, 'fraction_numerical': 0.010422120826788341, 'mean_word_length': 4.503293084522503, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 7, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Meta-Learning Dynamics Forecasting Using Task Inference

Current deep learning models for dynamics forecasting struggle with generalization. They can only forecast in a specific domain and fail when applied to systems with different parameters, external forces, or boundary conditions. We propose a model-based meta-learning method called DyAd which can generalize across heterogeneous domains by partitioning them into different tasks. DyAd has two parts: an encoder which infers the time-invariant hidden features of the task with weak supervision, and a forecaster which learns the shared dynamics of the entire domain. The encoder adapts and controls the forecaster during inference using adaptive instance normalization and adaptive padding. Theoretically, we prove that the generalization error of such procedure is related to the task relatedness in the source domain, as well as the domain differences between source and target. Experimentally, we demonstrate that our model outperforms state-of-the-art approaches on both turbulent flow and real-world ocean data forecasting tasks.


Introduction
Modeling dynamical systems with deep learning has shown great success in a wide range of systems from climate science, Internet of Things (IoT) to infectious diseases [23,64,8,33,42,22]. However, the main limitation of previous works is limited generalizability. Most approaches only focus on a specific system and train on past data in order to predict the future. Thus, a new model must be trained to predict a system with different dynamics. Consider, for example, learning fluid dynamics; shown in Fig. 1 are two fluid flows with different degrees of turbulence. Although the flows are governed by the same equations, the difference in buoyant forces would require two separate deep learning models to forecast. Therefore, it is imperative to develop generalizable deep learning models for dynamical systems that can learn and predict well over a large heterogeneous domain. Meta-learning [63,7,13], or learning to learn, improves generalization by learning multiple tasks from the environment. Recent developments in meta-learning have been successfully applied to few-shot classification [46,62], active learning [77], and reinforcement learning [18,20]. However, meta-learning in the context of forecasting high-dimensional physical dynamics has not been studied before. The challenges with meta-learning dynamical systems are unique in that (1) we need to efficiently infer the latent representation of the dynamical system given observed time series data, (2) we need to account for changes in unknown initial and boundary conditions, and (3) we need to model the temporal dynamics across heterogeneous domains.

arXiv:2102.10271v5 [cs.LG] 11 Oct 2022
Our approach is inspired by the fact that similar dynamical systems may share time-invariant hidden features. Even the slightest change in these features may lead to vastly different phenomena. For example, in climate science, fluids are governed by a set of differential equations called Navier-Stokes equations. Some features such as kinematic viscosity and external forces (e.g. gravity), are time-invariant and determine the flow characteristics. By inferring this latent representation, we can model diverse system behaviors from smoothly flowing water to atmospheric turbulence.
Inspired by neural style transfer [26], we propose a model-based meta-learning method, called DyAd, which can rapidly adapt to systems with varying dynamics. DyAd has two parts, an encoder g and a forecaster f . The encoder maps different dynamical systems to time-invariant hidden features representing constants of motion, boundary conditions, and external forces which characterize the system. The forecaster f then takes the learnt hidden features from the encoder and the past system states to forecast the future system state. Controlled by the time-invariant hidden features, the forecaster has the flexibility to adapt to a wide range of systems with heterogeneous dynamics.
Unlike gradient-based meta-learning techniques such as MAML [13], DyAd automatically adapts during inference using an encoder and does not require retraining. Similar to model-based metalearning methods such as MetaNets [46], we employ a two-part design with an adaptable learner which receives task-specific weights. However, for time-series forecasting, since input and output come from the same domain, a support set of labeled data is unnecessary to define the task. The encoder can infer the task directly from the query input.

Our contributions include:
• A novel model-based meta-learning method (DyAd) for dynamic forecasting in a large heterogeneous domain.
• An encoder capable of extracting the time-invariant hidden features of a dynamical system using time-shift invariant model structure and weak supervision.
• A new adaptive padding layer (AdaPad), designed for adapting to boundary conditions.
• Theoretical guarantees for DyAd on the generalization error of task inference in the source domain as well as domain adaptation to the target domain.
• Improved generalization performance on heterogeneous domains such as fluid flow and sea temperature forecasting, even to new tasks outside the training distribution.

Meta-learning in dynamics forecasting
Let x ∈ R d be a d-dimensional state of a dynamical system governed by parameters ψ. The problem of dynamics forecasting is that given a sequence of past states (x 1 , . . . , x t ), we want to learn a map f such that: Here l is the length of the input series, and h is the forecasting horizon in the output. Existing approaches for dynamics forecasting only predict future data for a specific system as a single task.
Here a task refers to forecasting for a specific system with a given set of parameters. The resulting models often generalize poorly to different system dynamics. Thus a new model must be trained to predict for each specific system.
To perform meta-learning, we identify each forecasting task by some parameters c ⊂ ψ, such as constants of motion, external forces, and boundary conditions. We learn multiple tasks simultaneously and infer the task from data. Here we use c for a subset of system parameters ψ, because we usually do not have the full knowledge of the system dynamics. In the turbulent flow example, the state x t is the velocity field at time t. Parameters c can represent Reynolds number, mean vorticity, mean magnitude, or a vector of all three.
Let µ be the data distribution over X × Y representing the function f : X → Y where X = R d×t and Y = R d×k . Our main assumption is that the domain X can be partitioned into separate tasks X = ∪ c∈C X c , where X c is the domain for task c and C is the set of all tasks. The data in the same task share the same set of parameters. Let µ c be the conditional distribution over X c × Y for task c. During training, the model is presented with data drawn from a subset of tasks {(x, y) : (x, y) ∼ µ c , c ∼ C}. Our goal is to learn the function f : X → Y over the whole domain X which can thus generalize across all tasks c ∈ C. To do so, we need to learn the map g : X → C taking x ∈ X c to c in order to infer the task with minimal supervision.

DyAd: Dynamic Adaptation Network
We propose a model-based meta-learning approach for dynamics forecasting. Given multiple forecasting tasks, we propose to learn the function f in two stages. That is, by first inferring the task c from the input x, and then adapting to a specialized forecaster f c : X c → Y for each task. An alternative is to use a single deep neural network to directly model f in one step over the whole domain. But this requires the training set to have good and uniform coverage of the different tasks. If the data distribution µ is highly heterogeneous or the training set is not sampled i.i.d. from the whole domain, then a single model may struggle with generalization.
We hypothesize that by partitioning the domain into different tasks, the model would learn to pick up task-specific features without requiring uniform coverage of the training data. Furthermore, by separating task inference and forecasting into two stages, we allow the forecaster to rapidly adapt to new tasks that never appeared in the training set.
As shown in Fig. 2, our model consists of two parts: an encoder g and a forecaster f . We introduce z c as a time-invariant hidden feature for task c. We assume that c depends linearly on the hidden feature for simplicity and easy interpretation. We design the encoder to infer the hidden feature z c given the input x. We then use z c to adapt the forecaster f to the specific task, i.e., model y = f c (x) as y = f (x, z c ). As the system dynamics are encoded in the input sequence x, we can feed the same input sequence x to a forecaster and generate predictionsŷ = f c (x).

Encoder Network
The encoder maps the input x to the hidden features z c that are time-invariant. To enforce this inductive bias, we encode time-invariance both in the architecture and in the training objective.
Time-Invariant Encoder. The encoder is implemented using 4 Conv 3D layers, each followed by BatchNorm, LeakyReLU, and max-pooling. Note that theoretically, max-pooling is not perfectly shift invariant since 2x2x2 max-pooling is equivariant to shifts of size 2 and only approximately invariant to shifts of size 1. But standard convolutional architectures often include max-pooling layers to boost performance. We convolve both across spatial and temporal dimensions.
After that, we use a global mean-pooling layer and a fully connected layer to estimate the hidden featureẑ c . The task parameter depends linearly on the hidden feature. We use a fully connected layer to compute the parameter estimateĉ.
Since convolutions are equivariant to shift (up to boundary frames) and mean pooling is invariant to shift, the encoder is shift-invariant. In practice, shifting the time sequence one frame forward will add one new frame at the beginning and drop one frame at the end. This creates some change in output value of the encoder. Thus, practically, the encoder is only approximately shift-invariant.
Encoder Training. The encoder network g is trained first. To combat the loss of shift invariance from the change from the boundary frames, we train the encoder using a time-invariant loss. Given two training samples (x (i) , y (i) ) and (x (j) , y (j) ) and their task parameters c, we have loss is an affine transformation of z c . The first term ĉ − c 2 uses weak supervision of the task parameters whenever they are available. Such weak supervision helps guide the learning of hidden feature z c for each task. While not all parameters of the dynamical system are known, we can compute approximate values in the datum c (i) based on our domain knowledge. For example, instead of the Reynolds number of the fluid flow, we can use the average vorticity as a surrogate for task parameters.
is the time-shift invariance loss, which penalizes the changes in latent variables between samples from different time steps. Since the time-shift invariance of convolution is only approximate, this loss term drives the time-shift error even lower. The third term | ẑ (i) c − m| 2 (m is a positive value) prevents the encoder from generating smallẑ (i) c due to time-shift invariance loss. It also helps the encoder to learn more interesting z, even in the absence of weak supervision.
Hidden Features. The encoder learns time-invariant hidden features. These hidden features resemble the time-invariant dimensionless parameters [30] in physical modeling, such as Reynolds number in fluid mechanics.
The hidden features may also be viewed as partial disentanglement of the system state. As suggested by [37,47], our disentanglement method is guided by inductive bias and training objectives. Unlike complete disentanglement, as in e.g. [39], in which the latent representation is factored into timeinvariant and time-varying components, we focus only on time-shift-invariance.
Nonetheless, the hidden features can control the forecaster which is useful for generalization.

Forecaster Network.
The forecaster incorporates the hidden feature z c from the encoder and adapts to the specific forecasting task f c = f (·, z c ). In what follows, we use z for z c . We use two specialized layers, adaptive instance normalization (AdaIN) and adaptive padding (AdaPad). AdaIN has been used in neural style transfer [26,21] to control generative networks. Here, AdaIN may adapt for specific coefficients and external forces. We also introduce a new layer, AdaPad(x, z), which is designed for encoding the boundary conditions of dynamical systems. The backbone of the forecaster network can be any sequence prediction model. We use a design that is similar to ResNet for spatiotemporal sequences. AdaIN. We employ AdaIN to adapt the forecaster network. Denote the channels of input x by x i and let µ(x i ) and σ(x i ) be the mean and standard deviation of channel i. For each AdaIN layer, a particular style is computed s = (µ i , σ i ) i = Az + b, where the linear map A and bias b are learned weights. Adaptive instance normalization is then defined as In essence, the channels are renormalized to the style s.
For dynamics forecasting, the hidden feature z encodes data analogous to the coefficients of a differential equation and external forces on the system. In numerical simulation of a differential equation these coefficients enter as scalings of different terms in the equation and the external forces are added to the combined force equation [24]. Thus in our context AdaIN, which scales channels and adds a global vector, is well-suited to injecting this information.
AdaPad. To complement AdaIN, we introduce AdaPad, which encodes the boundary conditions of each specific dynamical system. Generally when predicting dynamical systems, error is introduced along the boundaries since it is unknown how the dynamics interact with the boundary of the domain, and there may be unknown inflows or outflows. In our method, the inferred hidden feature z may contain the boundary information. AdaPad uses the hidden features to compute the boundary conditions via a linear layer. Then it applies the boundary conditions as padding immediately outside the spatial domain in each layer, as shown in Fig. 3.

Forecaster Training
The forecaster is trained separately after the encoder. The kernels of the convolutions and the mappings of the AdaIN and AdaPad layers are all trained simultaneously as the forecaster network is trained. Denote the true state as y and the predicted state asŷ, we compute the loss per time step ŷ − y 2 for each example. We accumulate the loss over different time steps and generate multi-step forecasts in an autoregressive fashion. In practice, we observe separate training achieves better performances than training end-to-end, see experiments for details.

Theoretical Analysis
The high-level idea of our method is to learn a good representation of the dynamics that generalizes well across a heterogeneous domain, and then adapt this representation to make predictions on new tasks. Our model achieves this by learning on multiple tasks simultaneously and then adapting to new tasks with domain transfer. We prove that learning the tasks simultaneously as opposed to independently results in better generalization (Proposition B.3 in Appendix). We also provide a theoretical decomposition for the generalization error.
Specifically, our hypothesis space has the form {x → f θ (x, g φ (x))} where φ and θ are the weights of the encoder and forecaster respectively. Let X be the error over the entire domain X , that is, for all c. Let enc (g φ ) = E x∼X (L 1 (g(x), g φ (x)) be the encoder error where g : X → C is the ground truth. Let G = {g φ : X → C} be the encoder hypothesis space. Denote the empirical risk of g φ byˆ enc (g φ ). W 1 denotes the Wasserstein distance between tasks. R(G) and R(F) represent the Rademacher complexities of encoder and forecaster. The following Proposition decomposes the generalization error in terms of forecaster error, encoder error, and distance between tasks. Proposition 2.1. Assume the forecaster c → f θ (·, c) is Lipschitz continuous with Lipschitz constant γ uniformly in θ and l ≤ 1/2. Let λ c = min f ∈F ( c (f ) + 1/K K k=1 c k (f )). For large n and probability ≥ 1 − δ, This result helps to quantify the trade-offs with respect to different terms in the error bound in our two-part architecture and the error of the model can be controlled by minimizing the empirical errors of the encoder and forecaster. See the Appendix Proposition B.6 for full proofs.

Related Work
Learning Dynamical Systems. Deep learning models are gaining popularity for learning dynamical systems [60,8,29,5,73,64,51]. An emerging topic is physics-informed deep learning [53,6,11,69,4,3,68,12,10] which integrates inductive biases from physical systems to improve learning. For example, [45,6] incorporated Koopman theory into the architecture. [43] used deep neural networks to solve PDEs with physical laws enforced in the loss functions. [17] and [9] build models upon Hamiltonian and Lagrangian mechanics that respect conservation laws. [68] proposed a hybrid approach by marrying two well-established turbulent flow simulation techniques with deep learning to produce a better prediction of turbulence. [10] propose a physics-informed GAN architecture that uses physics knowledge to inform the learning of both the generator and discriminator models. However, these approaches deal with a specific system dynamics instead of a large heterogeneous domain.
Multi-task learning and Meta-learning. Multi-task learning [65] focuses on learning shared representations from multiple related tasks. Architecture-based MTL methods can be categorized into encoder-focused [35] and decoder-focused [74]. There are also optimization-based MTL methods, such as task balancing methods [27]. But MTL assumes tasks are known a priori instead of inferring the task from data. On the other hand, the aim of meta-learning [63] is to leverage the shared representation to fast adapt to unseen tasks. Based on how the meta-level knowledge is extracted and used, meta-learning methods are classified into model-based [46,1,50,59,78], metric-based [67,61] and gradient-based [13,57,16,76]. Most meta-learning approaches are not designed for forecasting with a few exceptions. [72] proposed to train a domain classifier with weak supervision to help domain adaptation but focused on low-dimensioal time series classification problems. [50] designed a residual architecture for time series forecasting with a meta-learning parallel. [1] proposed a modular meta-learning approach for continuous control. But forecasting physical dynamics poses unique challenges to meta-learning as it requires encoding physical knowledge into our model.
Style Transfer. Our approach is inspired by neural style transfer techniques. In style transfer, a generative network is controlled by an external style vector through adaptive instance normalization between convolutional layers. Our hidden representation bears affinity with the "style" vector in style transfer techniques. Rather than aesthetic style in images, our hidden representation encodes time-invariant features. Style transfer initially appear in non-photorealistic rendering [31]. Recently, neural style transfer [25] has been applied to image synthesis [15,36], videos generation [56], and language translation [52]. For dynamical systems, [58] adapts texture synthesis to transfer the style of turbulence for animation. [28] studies unsupervised generative modeling of turbulent flows but for super-resolution reconstruction rather than forecasting.
Video Prediction. Our work is also related to video prediction. Conditioning on the historic observed frames, video prediction models are trained to predict future frames, e.g., [41,14,75,66,49,14,70,71,32,40,34]. There is also conditional video prediction [48] which achieves controlled synthesis. Many of these models are trained on natural videos from unknown physical processes. Our work is substantially different because we do not attempt to predict object or camera motions. However, our method can be potentially combined with video prediction models to improve generalization.

Datasets
We experiment on three datasets: synthetic turbulent flows, real-world sea surface temperature and ocean currents data. These are difficult to forecast using numerical methods due to unknown external forces and complex dynamics not fully captured by simplified mathematical models.
Turbulent Flow with Varying Buoyancy. We generate a synthetic dataset of turbulent flows with a numerical simulator, PhiFlow 1 . It contains 64×64 velocity fields of turbulent flows in which we vary the buoyant force acting on the fluid from 1 to 25. Each buoyant force corresponds to a forecasting task and there are 25 tasks in total. We use the mean vorticity of each task as partial supervision c as we can directly calculate it from the data. Vorticity can characterize formation and circular motion of turbulent flows.
Sea Surface Temperature. We evaluate on a real-world sea surface temperature data generated by the NEMO ocean engine [38] 2 . We select an area from Pacific ocean range from 01/01/2018 to 12/31/2020. The corresponding latitude and longitude are (-150∼-120, -20∼-50). This area is then divided into 25 64×64 subregions, each is a task since the mean temperature varies a lot along longitude and latitude. For the encoder training, we use season as an additional supervision signal besides the mean temperature of each subregion. In other words, the encoder should be able to infer the mean temperature of the subregion as well as to classify four seasons given the temperature series.
Ocean Currents. We also experiment with the velocity fields of ocean currents from the same region and use the same task division as the sea surface temperature data set. Similar to the turbulent flow data set, we use the mean vorticity of each subregion as the weak-supervision signal.

Baselines
We include several SoTA baselines from meta-learning and dynamics forecasting.
• U-net [55]: Originally developed for biomedical image segmentation, adapted for dynamics forecasting [11]. • ResNet/Unet-c: Above ResNet and Unet with an additional final layer that generates task parameter c and trained with weak-supervision and forecasting loss altogether.  Note that both ResNet-c and Unet-c have access to task parameters c. Modular-attn has a convolutional encoder f that takes the same input x as each module M to generate attention weights, m l=1 . Modular-wt also has the same encoder but to generate weights for combining the convolution parameters of all modules. MetaNet requires samples from test tasks as a support set and MAML needs adaptation retraining on test tasks, while other models do not need any information from the test domains. Thus, we use additional samples of up to 20% of the test set from test tasks. MetaNet uses these as a support set. MAML is retrained on these samples for 10 epoch for adaptation. To demonstrate the generalizability of DyAd, we experimented with ResNet and U-net as our forecaster.

Experiments Setup
For all datasets, we use a sliding-window approach to generate samples of sequences. We evaluate on two scenarios of generalization. For test-future, we train and test on the same task but different time steps. For test-domain, we train and test on different tasks with an 80-20 split. All models are trained to make next step prediction given the history. We forecast in an autoregressive manner to generate multi-step ahead predictions. All results are averaged over 3 runs with random initialization.
Apart from the root mean square error (RMSE), we also report the energy spectrum error (ESE) for ocean current prediction, which quantifies the physical consistency. ESE indicates whether the predictions preserve the correct statistical distribution and obey the energy conservation law, which is a critical metric for physical consistency. See details about energy spectrum and complete experiments details in Appendix A.1.  Table 1 shows the RMSE of multi-step predictions on Turbulent Flows (20 steps), Sea Surface Temperature (10 steps), and Ocean Currents (10 step) in two testing scenarios. We observe that DyAd makes the most accurate predictions in both scenarios across all datasets. Comparing ResNet/Unet-c with DyAd, we observe the clear advantage of task inference with separate training. VarSepNet achieves competitive performances on Ocean Currents (second best) through spatiotemporal disentanglement but cannot adapt to future data. Table 1 also reports ESEs on real-world Ocean Currents. DyAd not only has small RMSE but also obtains the smallest ESE, suggesting it captures the statistical distribution of ocean currents well. Figure 4 shows the target and the predicted velocity norm fields ( √ u 2 + v 2 ) by Unet-c, Modular-wt and DyAd at time step 1, 5, 10 for Turbulent Flows with buoyancy factors 9 and 21 respectively. We can see that DyAd can generate realistic flows with the corresponding characteristics while the baselines have trouble understanding and disentangling the buoyancy factor. Figure 5 shows DyAd, ResNet, U-net, PredRNN predictions on an ocean current sample in the future test set, and we see the shape of predictions by DyAd is closest to the target. These results demonstrate that DyAd not only forecasts well but also accurately captures the physical characteristics of the system.

Controllable Forecast.
DyAd infers the hidden features from data, which allows direct control of the latent space in the forecaster. We vary the encoder input while keeping the forecaster input fixed. Figure 6 right shows the forecasts from DyAd when the encoder is fed with flows having different buoyancy factors c = 5, 15, 25. As expected, with higher buoyancy factors, the predictions from the forecaster become more turbulent. This demonstrates that the encoder can successfully disentangle the latent representation of difference tasks, and control the predictions of the forecaster.

Ablation Study.
We performed several ablation studies of DyAd on the turbulence dataset to understand the contribution of each model component.
Ablation of the Model and Training Strategy. We first performed ablation study of the model architecture, shown in Table 2. We first remove the encoder from DyAd while keeping the same forecaster network (No_enc). The resulting model degrades but still outperforms ResNet. This demonstrates the effectiveness of AdaIN and AdaPad for forecasting. We also tested DyAd with AdaIN only (No_AdaPad), and the performance without AdaPad was slightly worse.
Another notable feature of our model is the ability to infer tasks with weakly supervised signals c. It is important to have a c that is related to the task domain. As an ablative study, we fed the encoder in DyAd with a random c, leading to Wrong_enc. We can see that having the totally wrong supervision may slightly hurt the forecasting performance. We also trained the encoder and the forecaster in DyAd altogether (End2End) but observed worse performance. This validates our hypothesis about the significance of domain partitioning and separate training strategy.
Ablation of Encoder Training Losses. We performed an ablation study of three training loss terms for the encoder to show the necessity of each loss term, shown in Table 3 We tried training the encoder only with the supervision of task parameters (Sup_c). We also tested training with the first two terms but without the magnitude loss (Sup_c+time_inv). We can see that missing any of three loss terms would make the encoder fail to learn the correct task-specific and time-invariant features.
Alternatives to AdaIN. We also tested 5 different alternatives to AdaIN for injecting the hidden feature z c into the forecaster, and reported the results in Table 4. We tested concatenating z c to the input of the forecaster, using z c as the kernels of the forecaster, concatenating z c to the hidden states, adding z c to the hidden states, using z c as the biases in the convolutional layers in the forecaster. AdaIN worked better than any alternative we tested.  Table 3: Ablation study of the encoder training loss terms. We tried training the encoder only with the weak supervision (Sup_c), without the magnitude loss (Sup_c+time_inv).  Table 4: We tested concatenating z c to the input of the forecaster (ConI), using z c as the kernels of the forecaster (KGen), concatenating z c to the hidden states (ConH), adding z c to the hidden states(Sum), using z c as the biases in the convolutional layers (Bias). AdaIN worked better than any alternative we tested.

Conclusion
We propose a model-based meta-learning method, DyAd to forecast physical dynamics. DyAd uses an encoder to infer the parameters of the task and a prediction network to adapt and forecast giving the inferred task. Our model can also leverage any weak supervision signals that can help distinguish different tasks, allowing the incorporation of additional domain knowledge. On challenging turbulent flow prediction and real-world ocean temperature and currents forecasting tasks, we observe superior performance of our model across heterogeneous dynamics. Future work would consider non-grid data such as flows on a graph or a sphere.

A.1 Experiments Details
We use a sliding window approach to generate samples of sequences. For test-future, we train and test on the same tasks but different time steps. For test-domain, we train and test on different tasks with an 80%-20% split. All models are trained to make next step prediction given the previous steps as input. We forecast in an autoregressive manner to generate multi-step ahead predictions. All results are averaged over 3 runs with random initialization.
We tuned learning rate (1e-2∼1e-5), batch size (16∼64), the number of unrolled prediction steps for backpropogation (2∼6), and hidden size (64∼512). We fixed the number of historic input frames as 20. The models that use ResNet, including MAML, MetaNets and DyAd+ResNet, have same number of residual blocks for fair comparison. When we trained the encoder on turbulence and sea surface temperature, we used α = 1 and β = 1 in the Equ. 1. For ocean currents, we used α = 0.2 and β = 0.2. We performed all experiments on 4 V100 GPUs. Table 5 in Appendix A displays the number of parameters of each fine-tuned model on the turbulence dataset.
We compare our model with a series of baselines on the multi-step forecasting with different dynamics. We consider two testing scenarios: (1) dynamics with different initial conditions (test-future) and (2) dynamics with different parameters such as external force (test-domain). The first scenario evaluates the models' ability to extrapolate into the future for the same task. The second scenario is to estimate the capability of the models to generalize across different tasks.
Apart from the root mean square error (RMSE), we also report the energy spectrum error (ESE) for ocean current prediction which quantifies the physical consistency. ESE indicates whether the predictions preserve the correct statistical distribution and obey the energy conservation law, which is a critical metric for physical consistency. The turbulence kinetic energy spectrum E(k) is related to the mean turbulence kinetic energy as where the k is the wavenumber and t is the time step. The spectrum can describe the transfer of energy from large scales of motion to the small scales and provides a representation of the dependence of energy on frequency. Thus, the ESE can indicate whether the predictions preserve the correct statistical distribution and obey the energy conservation law. A trivial example that can illustrate why we need ESE is that if a model simply outputs moving averages of input frames, the accumulated RMSE of predictions might not be high but the ESE would be really big because all the small or even medium eddies are smoothed out.

B Theoretical Analysis
The high-level idea of our method is to learn a good representation of the dynamics that generalizes well across a heterogeneous domain, and then adapt this representation to make predictions on new tasks. Our model achieves this by learning on multiple tasks simultaneously and then adapting to new tasks with domain transfer. We prove that learning the tasks simultaneously as opposed  to independently results in tighter generalization bounds. As these are upper bounds, they do not necessarily imply better generalization, although empirically, this is what we observe. Although we are not able to provide a proof of our empirical finding that an encoder-forecaster architecture outperforms a monolithic architecture, we do provide an analysis of the error, decomposing it into error due to the forecaster, encoder, and domain adaptation. This decomposition helps to quantify the trade-offs involved in our two-part architecture.
Suppose we have K tasks {c k } K k=1 ∼ C, each of which is sampled from a continuous parameter space C. Here c k are the parameters of the task, which can be inferred by the encoder. For each task c k , we have a collection of dataμ c k of size n, sampled from µ k , a shorthand for µ c k .

B.1 Multi-Task Learning Error
We want to bound the true loss using the empirical lossˆ and Rademacher complexity of the hypothesis class F. We can use the classic results from [2]. Define empirical Rademacher complexity for all tasks asR where {σ The true Rademacher complexity is then defined R(F) = E X (R X (F)).
The following theorem restates the main result from [2]. We simplify the statement by using Rademacher complexity rather than the set cover number argument used in the original proof. The assumption of bounded loss is reasonable given a prescribed data domain and data preprocessing. Theorem B.1. ( [2]) Assume the loss is bounded l ≤ 1/2. Given n samples each from K different forecasting tasks µ 1 , · · · , µ k , then with probability at least 1 − δ, the following inequality holds for each f ∈ F in the hypothesis class: k } as independent random variables. For a function φ that satisfies |φ(x (1) , · · · , x (i) , · · · x (n) ) − φ(x (1) , · · · ,x (i) , · · · x (n) )| ≤ c i by McDiarmid's inequality, we have Applying this inequality to the max difference Q(X) = sup f ∈F [ (f )−ˆ (f, X)], then with probability at least 1 − δ, we have where C is a constant depending on the bounds c i . If the loss l ≤ 1/2, then |Q| ≤ 1/2 and so we can take c i = 1 leading to C = 1/ √ 2. A standard computation (see [44] Theorem 3.3) using the law of total expectation shows E X [Q(X)] ≤ 2R(F), which finishes the proof.
The following inequality compares the performance for multi-task learning to learning individual tasks. Let R k (F) be the Rademacher complexity for F over µ k . Lemma B.2. The Rademacher complexity for multi-task learning is bounded R( Proof. We compute the empirical Rademacher complexity, Proposition B.3. Assume n = n k for all tasks k and the loss l is bounded l ≤ 1/2, then the generalization bound given by considering each task individually is which is weaker than the bound of Theorem B.1.
Proof. Applying the classical analog of Theorem B.1 to a single task, we find with probability greater than 1 − δ, Averaging over all tasks yields Since the loss l is bounded l ≤ 1/2, we can take C = C k = 1/ √ 2, giving the generalization upper bound for the joint error of Equation 3. By Lemma B.2 and the fact 1/ √ 2nK ≤ 1/ √ 2n, the bound in Theorem B.1 is strictly tighter.
The upper bound in Theorem B.1 is strictly tighter than that of Proposition B.3 as the first termsˆ (f ) are equal, the second term is smaller by Lemma B.2 and the third is smaller since 1/ √ 2nK ≤ 1/ √ 2n. This helps explain why our multitask learning forecaster has better generalization than learning each task independently. The shared data tightens the generalization bound. Ultimately, though a tighter upper bound suggests lower error, it does not strictly imply it. We further verify this theory in our experiments by comparison to baselines which learn each task independently.

B.2 Domain Adaptation Error
Since we test on c ∼ C outside the training set {c k }, we incur error due to domain adaptation from the source domains µ c1 , . . . , µ c K to target domain µ c with µ being the true distribution. Denote the corresponding empirical distributions of n samples per task byμ c = 1 For different c and c , the domains µ c and µ c may have largely disjoint support, leading to very high KL divergence. However, if c and c are close, samples x c ∼ µ c and x c ∼ µ c may be close in the domain X with respect to the metric · X . For example, if the external forces c and c are close, given x c ∼ µ c there is likely x c ∼ µ c such that the distance between the velocity fields x c − x c is small. This implies the distributions µ c and µ c may be be close in Wasserstein distance W 1 (µ c , µ c ). The bound from [54] applies well to our setting as such: Theorem B.4. Let λ c = min f ∈F c (f ) + 1/K K k=1 c k (f ) . There is N = N (dim(X )) such that for n > N , for any hypothesis f , with probability at least 1 − δ, Proof. We apply [54] Theorem 2 to target domain µ T = µ c and joint source domain µ S = 1/K K k=1 µ c k with empirical samplesμ T =μ c andμ S = 1/K K k=1μ c k .

B.3 Encoder versus Prediction Network Error
Our goal is to learn a joint hypothesis h over the entire domain X in two steps, first inferring the task c and then inferring x t+1 conditioned on c. Error from DyAd may result from either the encoder g φ or the prediction network f θ . Our hypothesis space has the form {x → f θ (x, g φ (x))} where φ and θ are the weights of the encoder and prediction network respectively. Let X be the error over the entire domain X , that is, for all c. Let enc (g φ ) = E x∼X (L 1 (g(x), g φ (x)) be the encoder error where g : X → C is the ground truth. We state a result that decomposes the final error into that attributable to the encoder and that to the forecaster. Proposition B.5. Assume c → f θ (·, c) is Lipschitz continuous with Lipschitz constant γ uniformly in θ. Then we bound X (f θ (·, g φ (·))) ≤ γ enc (g φ ) + E c∼C [ c (f θ (x, c))] (4) where the first term is the error due to the encoder incorrectly identifying the task and the second term is the error due the prediction network alone.
The hypothesis in the second term consists of the prediction network combined with the ground truth task label x → f θ (x, g(x)). In practice, the Lipschitz constant γ of the encoder may be indirectly minimized through weight decay.
Proof. By the triangle inequality and linearity of expectation, By Lipschitz continuity, which, since g(x) = c and by linearity of expectation, by definition of enc and c ,= γ enc (g φ ) + E c∼C [ c (f θ (x, c))].
By combining Theorem B.1, Proposition B.5, and Theorem B.4, we can bound the generalization error in terms of the empirical error of the prediction network on the source domains, the Wasserstein distance between the source and target domains, and the empirical error of the encoder. Let G = {g φ : X → C} be the task encoder hypothesis space. Denote the empirical risk of the encoder g φ with respect to X byˆ enc (g φ ).
In order to apply Theorem B.1 to the risk c and relate it to the empirical risk, we first relate the error on the target domain back to the source domain of our empirical samples. By Theorem B.4, c k (f θ (·, c k )) + W 1 μ c , 1 K K k=1μ c k + 2 log(1/δ) 1/n + 1/(nK) + λ c .
Using linearity of the expectation over c, removing it where there is no dependence on c, and rearranging terms.
Although this result does not settle either the question of end-to-end versus pre-training or encoderforecaster versus monolithic model, it quantifies the trade-offs these choices depend on. We empirically consider both questions in the experiments section.