dabench.data ============ .. py:module:: dabench.data .. autoapi-nested-parse:: Data generators and downloaders Classes ------- .. autoapisummary:: dabench.data.Barotropic dabench.data.Data dabench.data.ENSOIndices dabench.data.GCP dabench.data.Lorenz63 dabench.data.Lorenz96 dabench.data.PyQG dabench.data.PyQGJax dabench.data.QGS dabench.data.SQGTurb Package Contents ---------------- .. py:class:: Barotropic(beta = 0.0, rek = 0.0, rd = 0.0, H = 1.0, nx = 256, ny = None, L = 2 * np.pi, x0 = None, delta_t = 0.001, taveint = 1, ntd = 1, time_dim = None, store_as_jax = False, **kwargs) Barotropic model data generator based on pyqg This data class is a wrapper of a "optional" pyqg package. See https://pyqg.readthedocs.io .. rubric:: Notes DEPRECATED Uses default attribute values from pyqg.BTModel: https://pyqg.readthedocs.io/en/latest/api.html#pyqg.BTModel Those values originally come from Mcwilliams 1984: J. C. Mcwilliams (1984). The emergence of isolated coherent vortices in turbulent flow. Journal of Fluid Mechanics, 146, pp 21-43 doi:10.1017/S0022112084001750. :param system_dim: system dimension :param beta: Gradient of coriolis parameter. Units: meters^-1 * seconds^-1. Default is 0. :param rek: Linear drag in lower layer. Units: seconds^-1. Default is 0. :param rd: Deformation radius. Units: meters. Default is 0. :param H: Layer thickness. Units: meters. Default is 1. :param nx: Number of grid points in the x direction. Default is 256. :param ny: Number of grid points in the y direction. Default: nx. :param L: Domain length in x direction. Units: meters. Default is 2*pi. :param x0: the initial conditions. Can also be provided when initializing model object. If provided by both, the generate() arg takes precedence. :param W: Domain width in y direction. Units: meters. Default: L. :param filterfac: amplitdue of the spectral spherical filter. Default is 23.6. :type filterfac: float :param delta_t: Numerical timestep. Units: seconds. :param taveint: Time interval for accumulation of diagnostic averages. For performance purposes, averaging does not have to occur every timestep. Units: seconds. Default is 1 (i.e. every 1000 timesteps when delta_t = 0.001) :param ntd: Number of threads to use. Should not exceed the number of cores on your machine. :param store_as_jax: Store values as jax array instead of numpy array. Default is False (store as numpy). .. py:method:: __advance__() Advances the QG model according to set attributes :returns: Array of absolute potential vorticity (relative potential vorticity + background vorticity). .. py:method:: forecast(n_steps=None, t_final=None, x0=None) Alias for self.generate(), except returns values as output .. py:method:: generate(n_steps = None, t_final = 40, x0 = None) Generates values and times, saves them to the data object .. rubric:: Notes Either provide n_steps or t_final in order to indicate the length of the forecast. These are used to set the values, times, and time_dim attributes. :param n_steps: Number of timesteps. Default is None, which sets :param n_steps to t_final/delta_t: :param t_final: Final time of trajectory. Default is 40, which results in n_steps = 40000 :param x0: the initial conditions. Can also be provided when initializing model object. If provided by both, the generate() arg takes precedence. .. py:class:: Data(system_dim = 3, original_dim = None, random_seed = 37, delta_t = 0.01, store_as_jax = False, x0 = None, **kwargs) Base for all data generator objects. :param system_dim: system dimension :param original_dim: dimensions in original space, e.g. could be 3x3 for a 2d system with system_dim = 9. Defaults to (system_dim), i.e. 1d. :param random_seed: random seed, defaults to 37 :param delta_t: the timestep of the data (assumed uniform) :param store_as_jax: Store values as jax array instead of numpy array. Default is False (store as numpy). .. py:method:: calc_lyapunov_exponents_final(total_time = None, rescale_time = 1, convergence = 0.01, x0 = None) Computes the final Lyapunov Exponents .. rubric:: Notes See self.calc_lyapunov_exponents_series for full info :param total_time: Time to integrate over to compute LEs. Usually there's a tradeoff between accuracy and computation time (more total_time leads to higher accuracy but more computation time). Default depends on model type and are based roughly on how long it takes for satisfactory convergence: For Lorenz63: n_steps=15000 (total_time=150 for delta_t=0.01) For Lorenz96: n_steps=50000 (total_time=500 for delta_t=0.01) :param rescale_time: Time for when the algorithm rescales the propagator to reduce the exponential growth in errors. Default is 1 (i.e. 100 timesteps when delta_t = 0.01). :param convergence: Prints warning if LE convergence is below this number. Default is 0.01. :param x0: initial condition to start computing LE. Needs to be on the attractor (i.e., remove transients). Default is None, which will fallback to use the x0 set during model object initialization. :returns: Lyapunov exponents array of size (system_dim) .. py:method:: calc_lyapunov_exponents_series(total_time = None, rescale_time = 1, convergence = 0.01, x0 = None) Computes the spectrum of Lyapunov Exponents. .. rubric:: Notes Lyapunov exponents help describe the degree of "chaos" in the model. Make sure to plot the output to check that the algorithm converges. There are three ways to make the estimate more accurate: 1. Decrease the delta_t of the model 2. Increase total_time 3. Decrease rescale time (try increasing total_time first) Algorithm: Eckmann 85, https://www.ihes.fr/~ruelle/PUBLICATIONS/%5B81%5D.pdf pg 651 This method computes the full time series of Lyapunov Exponents, which is useful for plotting for debugging. To get only the final Lyapunov Exponent, use self.calc_lyapunov_exponents. :param total_time: Time to integrate over to compute LEs. Usually there's a tradeoff between accuracy and computation time (more total_time leads to higher accuracy but more computation time). Default depends on model type and are based roughly on how long it takes for satisfactory convergence: For Lorenz63: n_steps=15000 (total_time=150 for delta_t=0.01) For Lorenz96: n_steps=50000 (total_time=500 for delta_t=0.01) :param rescale_time: Time for when the algorithm rescales the propagator to reduce the exponential growth in errors. Default is 1 (i.e. 100 timesteps when delta_t = 0.01). :param convergence: Prints warning if LE convergence is below this number. Default is 0.01. :param x0: initial condition to start computing LE. Needs to be on the attractor (i.e., remove transients). Default is None, which will fallback to use the x0 set during model object initialization. :returns: Lyapunov exponents for all timesteps, array of size (total_time/rescale_time - 1, system_dim) .. py:method:: generate(n_steps = None, t_final = None, x0 = None, M0 = None, return_tlm = False, stride = None, **kwargs) Generates a dataset and returns xarray state vector. .. rubric:: Notes Either provide n_steps or t_final in order to indicate the length of the forecast. :param n_steps: Number of timesteps. One of n_steps OR t_final must be specified. :param t_final: Final time of trajectory. One of n_steps OR t_final must be specified. :param x0: initial conditions state vector of shape (system_dim) :param M0: the initial condition of the TLM matrix computation shape (system_dim, system_dim). :param return_tlm: specifies whether to compute and return the integrated Jacobian as a Tangent Linear Model for each timestep. :param stride: specify how many steps to skip in the output data versus the model timestep (delta_t). :param \*\*kwargs: arguments to the integrate function (permits changes in convergence tolerance, etc.). :returns: Xarray Dataset of output vector, and if return_tlm=True Xarray DataArray of TLMs corresponding to the system trajectory. .. py:method:: load_netcdf(filepath = None, include_vars = None, exclude_vars = None, years_select = None, dates_select = None) Loads values from netCDF file, saves them in values attribute :param filepath: Path to netCDF file to load. If not given, defaults to loading ERA5 ECMWF SLP data over Japan from 2018 to 2021. :param include_vars: Data variables to load from NetCDF. If None (default), loads all variables. Can be used to exclude bad variables. :param exclude_vars: Data variabes to exclude from NetCDF loading. If None (default), loads all vars (or only those specified in include_vars). It's recommended to only specify include_vars OR exclude_vars (unless you want to do extra typing). :param years_select: Years to load (ints). If None, loads all timesteps. :param dates_select: Dates to load. Elements must be datetime date or datetime objects, depending on type of time indices in NetCDF. If both years_select and dates_select are specified, time_stamps overwrites "years" argument. If None, loads all timesteps. .. py:method:: rhs_aux(x, t) The auxiliary model used to compute the TLM. :param x: State vector with size (system_dim) :param t: Array of times with size (time_dim) :returns: State vector [size: (system_dim,)] :rtype: dxaux (ndarray) .. py:method:: save_netcdf(ds, filename) Saves values in values attribute to netCDF file :param ds: Xarray dataset :param filepath: Path to netCDF file to save .. py:class:: ENSOIndices(file_dict = None, var_types = None, system_dim = None, store_as_jax = False, **kwargs) Gets ENSO indices from CPC website .. rubric:: Notes Source: https://www.cpc.ncep.noaa.gov/data/indices/ :param system_dim: system dimension :param store_as_jax: Store values as jax array instead of numpy array. Default is False (store as numpy). :param file_dict: Lists of files to get. Dict keys are type of data: 'wnd': Wind 'slp': Sea level pressure 'soi': Southern Oscillation Index 'soi3m' Southerm Oscillation Index 3-month running mean 'sst': SST indices 'desst': Detrended Nino3.4 index 'rsst': Regional SST index (North and South Atlantic, and Tropics) 'olr': Outgoing long-wave radiance 'cpolr': Central Pacific OLR Dict values are individual files from the website, see full list at https://www.cpc.ncep.noaa.gov/data/indices/ Default is {'wnd': ['zwnd200'], 'slp': ['darwin']} :param var_types: List of variables within file to get. Dict keys are type of data (see list in file_dict description). Dict values are type of variable: 'ori' = Original 'ano' = Anomaly 'std' = Standardized Anomaly sst (sea surface temp) has nino3, nino34, and nino12 options, see for more info: https://climatedataguide.ucar.edu/climate-data/nino-sst-indices-nino-12-3-34-4-oni-and-tni rsst (regional sea surface temp) has prefixes: 'sa_' is South Atlantic, 'na_' is North Atlantic, 'tr_' is Tropics. Default depends on data type (first string in each value list) 'wnd'='ori', 'lp'='ori', 'soi'='ano', 'soi3m'='ori', 'eqsoi'='std', 'sst'='nino12', 'desst'='ori', 'rsst'='na', 'olr'='ori', 'cpolr'='ano' .. py:method:: generate() Alias for _load_gcp_era5 .. py:class:: GCP(variables = ['2m_temperature'], date_start = '2020-01-01', date_end = '2020-12-31', min_lat = 19.8554808619, max_lat = 23.1886107447, min_lon = -84.9749110583, max_lon = -74.1780248685, store_as_jax = False, **kwargs) Loads ERA5 data from Google Cloud Platform .. rubric:: Notes Source: https://cloud.google.com/storage/docs/public-datasets/era5 Data is hourly :param variables: Names of ERA5 variables to load. For description of variables, see: https://github.com/google-research/arco-era5?tab=readme-ov-file#full_37-1h-0p25deg-chunk-1zarr-v3 Default is ['2m_temperature'] (Air temperature at 2 metres) :param date_start: Start of time range to download, in 'yyyy-mm-dd' format. Can also just specify year ('yyyy') or year and month ('yyyy-mm'). Default is '2020-06-01'. :param date_end: End of time range to download, in 'yyyy-mm-dd' format. Can also just specify year ('yyyy') or year and month ('yyyy-mm'). Default is '2020-9-30'. :param min_lat: Minimum latitude for bounding box. If None, loads global data (which can be VERY large). Bounding box default covers Cuba. :param max_lat: Max latitude for bounding box (see min_lat for info). :param min_lon: Min latitude for bounding box (see min_lat for info). :param max_lon: Max latitude for bounding box (see min_lat for info). :param store_as_jax: Store values as jax array instead of numpy array. Default is False (store as numpy). .. py:method:: generate() Alias for _load_gcp_era5 .. py:method:: load() Alias for _load_gcp_era5 .. py:class:: Lorenz63(sigma = 10.0, rho = 28.0, beta = 8.0 / 3.0, delta_t = 0.01, x0 = jnp.array([-10.0, -15.0, 21.3]), system_dim = 3, values = None, store_as_jax = False, **kwargs) Lorenz 63 model data generator. :param sigma: Lorenz 63 param. Default is 10., the original value used in Lorenz, 1963. https://doi.org/10.1175/1520-0469(1963)020<0130:DNF>2.0.CO;2 :param rho: Lorenz 63 param. Default is 28., the value used in Lorenz, 1963 (see DOI above) :param beta: Lorenz 63 param. Default is 8./3., the value used in Lorenz, 1963 (see DOI above) :param delta_t: length of one time step :param x0: Initial state, array of floats of size (system_dim). Default is jnp.array([-10.0, -15.0, 21.3]), which is the system state after a 6000 step spinup with delta_t=0.01 and initial conditions [0., 1., 0.], a spinup which replicates the simulation described in Lorenz, 1963. :param system_dim: system dimension. Must be 3 for Lorenz63. :param store_as_jax: Store values as jax array instead of numpy array. Default is False (store as numpy). .. py:method:: Jacobian(x) Jacobian of the L63 system :param x: state vector with shape (system_dim) :returns: ndarray of Jacobian matrix with shape (system_dim, system_dim) .. py:method:: rhs(x, t = None) vector field of Lorenz 63 :param x: state vector with shape (system_dim) :param t: times vector. Needed as argument slot for odeint but unused :returns: vector field of Lorenz 63 with shape (system_dim) .. py:class:: Lorenz96(forcing_term = 8.0, delta_t = 0.05, x0 = None, system_dim = 36, values = None, store_as_jax = False, **kwargs) Lorenz 96 model data generator. .. rubric:: Notes Default values come from Lorenz, 1996: eapsweb.mit.edu/sites/default/files/Predicability_a_Problem_2006.pdf :param forcing_term: Forcing constant for Lorenz96, prevents energy from decaying to 0. Default is 8.0. :param x0: Initial state vector, array of floats of size (system_dim). For system_dim of 5, 6, or 36, defaults are the final state of a 14400 timestep spinup starting with an initial state (x0) of all 0s except the first element, which is set to 0.01. This setup is taken from Lorenz, 1996. Defaults are provided from spinups usin delta_t = 0.05 (Lorenz, 1996 original value) and delta_t = 0.01 (more frequently used today) For all other system_dim settings, default is all 0s except the first element, which is set to 0.01. :param system_dim: System dimension, must be between 4 and 40. Default is 36. :param delta_t: Length of one time step. Default is 0.05 from Lorenz, 1996, but on modern computers 0.01 is often used. :param store_as_jax: Store values as jax array instead of numpy array. Default is False (store as numpy). .. py:method:: Jacobian(x) Computes the Jacobian of the Lorenz96 system :param x: state vector with shape (system_dim) :returns: ndarray of Jacobian matrix with shape (system_dim, system_dim) .. py:method:: rhs(x, t = None) Computes vector field of Lorenz 96 :param x: state vector with shape (system_dim) :param t: times vector. Required as argument slot for some numerical integrators but unused :returns: vector field of Lorenz 96 with shape (system_dim) .. py:class:: PyQG(beta=1.5e-11, rd=15000.0, delta=0.25, H1=500, H2=None, U1=0.025, U2=0.0, x0=None, twrite=10000, nx=64, ny=None, delta_t=7200, ntd=1, time_dim=None, values=None, times=None, store_as_jax=False, **kwargs) PyQG quasi-geotropic model data generator. The PyQG class is simply a wrapper of a "optional" pyqg package. See https://pyqg.readthedocs.io .. rubric:: Notes DEPRECATED Uses default attribute values from pyqg.QGModel: https://pyqg.readthedocs.io/en/latest/api.html#pyqg.QGModel :param beta: Gradient of coriolis parameter. Units: meters^-1 * seconds^-1 :type beta: float :param rek: Linear drag in lower layer. Units: seconds^-1 :type rek: float :param rd: Deformation radius. Units: meters. :type rd: float :param delta: Layer thickness ratio (H1/H2) :type delta: float :param U1: Upper layer flow. Units: m/s :type U1: float :param U2: Lower layer flow. Units: m/s :type U2: float :param H1: Layer thickness (sets both H1 and H2). :type H1: float :param nx: Number of grid points in the x direction. :type nx: int :param ny: Number of grid points in the y direction (default: nx). :type ny: int :param L: Domain length in x direction. Units: meters. :type L: float :param W: Domain width in y direction. Units: meters (default: L). :type W: float :param filterfac: amplitude of the spectral spherical filter (originally 18.4, later changed to 23.6). :type filterfac: float :param delta_t: Numerical timestep. Units: seconds. :type delta_t: float :param twrite: Interval for cfl writeout. Units: number of timesteps. :type twrite: int :param tmax: Total time of integration (overwritten by t_final). Units: seconds. :type tmax: float :param ntd: Number of threads to use. Should not exceed the number of cores on your machine. :type ntd: int :param store_as_jax: Store values as jax array instead of numpy array. Default is False (store as numpy). :type store_as_jax: bool .. py:method:: __advance__() Advances the QG model according to set attributes :returns: Array of absolute potential vorticity (relative potential vorticity + background vorticity). .. py:method:: forecast(n_steps=None, t_final=None, x0=None) Alias for self.generate(), except returns values as output .. py:method:: generate(n_steps=None, t_final=None, x0=None) Generates values and times, saves them to the data object .. rubric:: Notes Either provide n_steps or t_final in order to indicate the length of the forecast. These are used to set the values, times, and time_dim attributes. :param n_steps: Number of timesteps. One of n_steps OR t_final must be specified. :type n_steps: int :param t_final: Final time of trajectory. One of n_steps OR t_final must be specified. :type t_final: float :param x0: the initial conditions. Can also be provided when initializing model object. If provided by both, the generate() arg takes precedence. :type x0: ndarray, optional .. py:class:: PyQGJax(beta = 1.5e-11, rd = 15000.0, delta = 0.25, H1 = 500, H2 = None, U1 = 0.025, U2 = 0.0, x0 = None, nx = 64, ny = None, delta_t = 7200, random_seed = 37, store_as_jax = False, **kwargs) PyQGJax quasi-geotropic model data generator. The PyQGJax class is simply a wrapper of the "optional" pyqg-jax package. See https://pyqg-jax.readthedocs.io .. rubric:: Notes Uses default attribute values from pyqg_jax.QGModel: https://pyqg.readthedocs.io/en/latest/api.html#pyqg.QGModel :param beta: Gradient of coriolis parameter. Units: meters^-1 * seconds^-1 :param rd: Deformation radius. Units: meters. :param delta: Layer thickness ratio (H1/H2) :param H1: Layer thickness (sets both H1 and H2 if H2 not specified). :param H2: Layer 2 thickness. :param U1: Upper layer flow. Units: m/s :param U2: Lower layer flow. Units: m/s :param x0: the initial conditions. Can also be provided when initializing model object. If provided by both, the generate() arg takes precedence. :param nx: Number of grid points in the x direction. :param ny: Number of grid points in the y direction (default: nx). :param delta_t: Numerical timestep. Units: seconds. :param store_as_jax: Store values as jax array instead of numpy array. Default is False (store as numpy). .. py:method:: forecast(n_steps=None, t_final=None, x0=None) Alias for self.generate(), except returns values as output .. py:method:: generate(n_steps = None, t_final = 40, x0 = None) Generates values and times, saves them to the data object .. rubric:: Notes Either provide n_steps or t_final in order to indicate the length of the forecast. :param n_steps: Number of timesteps. One of n_steps OR t_final must be specified. :param t_final: Final time of trajectory. One of n_steps OR t_final must be specified. :param x0: the initial conditions. Can also be provided when initializing model object. If provided by both, the generate() arg takes precedence. .. py:class:: QGS(model_params = None, x0 = None, delta_t = 0.1, system_dim = None, store_as_jax = False, random_seed = 37, **kwargs) QGS quasi-geostrophic model data generator. The QGS class is simply a wrapper of an *optional* qgs package. See https://qgs.readthedocs.io/ :param model_params: qgs parameter object. See: https://qgs.readthedocs.io/en/latest/files/technical/configuration.html#qgs.params.params.QgParams If None, will use defaults specified by: De Cruz, et al. (2016). Geosci. Model Dev., 9, 2793-2808. :param x0: Initial state vector, array of floats. Default is: :param delta_t: Numerical timestep. Units: seconds. :param store_as_jax: Store values as jax array instead of numpy array. Default is False (store as numpy). .. py:method:: Jacobian(x, t = 0) Jacobian of the qgs system :param x: State vector of size (system_dim) :param t: times vector. Required as argument slot for some numerical integrators but unused. :returns: Jacobian matrix of size (system_dim, system_dim) .. py:method:: calc_lyapunov_exponents_series(total_time = None, rescale_time = 1, convergence = 0.01, x0 = None) Computes the spectrum of Lyapunov Exponents. .. rubric:: Notes Lyapunov exponents help describe the degree of "chaos" in the model. Make sure to plot the output to check that the algorithm converges. There are three ways to make the estimate more accurate: 1. Decrease the delta_t of the model 2. Increase total_time 3. Decrease rescale time (try increasing total_time first) Algorithm: Eckmann 85, https://www.ihes.fr/~ruelle/PUBLICATIONS/%5B81%5D.pdf pg 651 This method computes the full time series of Lyapunov Exponents, which is useful for plotting for debugging. To get only the final Lyapunov Exponent, use self.calc_lyapunov_exponents. :param total_time: Time to integrate over to compute LEs. Usually there's a tradeoff between accuracy and computation time (more total_time leads to higher accuracy but more computation time). Default depends on model type and are based roughly on how long it takes for satisfactory convergence: For Lorenz63: n_steps=15000 (total_time=150 for delta_t=0.01) For Lorenz96: n_steps=50000 (total_time=500 for delta_t=0.01) :param rescale_time: Time for when the algorithm rescales the propagator to reduce the exponential growth in errors. Default is 1 (i.e. 100 timesteps when delta_t = 0.01). :param convergence: Prints warning if LE convergence is below this number. Default is 0.01. :param x0: initial condition to start computing LE. Needs to be on the attractor (i.e., remove transients). Default is None, which will fallback to use the x0 set during model object initialization. :returns: Lyapunov exponents for all timesteps, array of size (total_time/rescale_time - 1, system_dim) .. py:method:: generate(n_steps = None, t_final = None, x0 = None, M0 = None, return_tlm = False, stride = None, **kwargs) Generates a dataset and assigns values and times to the data object. .. rubric:: Notes Either provide n_steps or t_final in order to indicate the length of the forecast. :param n_steps: Number of timesteps. One of n_steps OR t_final must be specified. :type n_steps: int :param t_final: Final time of trajectory. One of n_steps OR t_final must be specified. :type t_final: float :param M0: the initial condition of the TLM matrix computation shape (system_dim, system_dim). :type M0: ndarray :param return_tlm: specifies whether to compute and return the integrated Jacobian as a Tangent Linear Model for each timestep. :type return_tlm: bool :param x0: initial conditions state vector of shape (system_dim) :type x0: ndarray :param stride: specify how many steps to skip in the output data versus the model timestep (delta_t). :type stride: int :param \*\*kwargs: arguments to the integrate function (permits changes in convergence tolerance, etc.). :returns: Xarray Dataset of output vector, and if return_tlm=True Xarray DataArray of TLMs corresponding to the system trajectory. .. py:method:: rhs(x, t = 0) Vector field (tendencies) of qgs system :param x: State vector of size (system_dim) :param t: times vector. Required as argument slot for some numerical integrators but unused. :returns: Vector field of qgs .. py:method:: rhs_aux(x, t) The auxiliary model used to compute the TLM. :param x: State vector with size (system_dim) :param t: Array of times with size (time_dim) :returns: State vector of size (system_dim,) .. py:class:: SQGTurb(pv = None, f = 0.0001, nsq = 0.0001, L = 20000000.0, H = 10000.0, U = 30.0, r = 0.0, tdiab = 10.0 * 86400, diff_order = 8, diff_efold = 86400.0 / 3, symmetric = True, dealias = True, precision = 'single', tstart = 0, delta_t = 900, store_as_jax = False, **kwargs) SQGTurb model data generator. :param pv: Potential vorticity array. If None (default), loads data from 57600 step spinup with initial conditions taken from Jeff Whitaker's original implementation: https://github.com/jswhit/sqgturb. 57600 steps matches the "nature run" spin up in that repository. :param system_dim: The dimension of the system state :param delta_t: model time step (seconds) :param x0: Initial state, array of floats of size (system_dim). :param f: coriolis :param nqr: Brunt-Vaisalla (buoyancy) freq squared :param L: size of square domain :param H: height of upper boundary :param U: basic state velocity at z = H :param r: Ekman damping (at z=0) :param tdiab: thermal relaxation damping :param diff_order: hyperdiffusion order :param diff_efold: hyperdiff time scale :param symmetric: symmetric jet, or jet with U=0 at sf :param dealias: if True, dealiasing applied using 2/3 rule :param precision: 'single' or 'double'. Default is 'single' :param tstart: initialize time counter :param delta_t: the timestep of the data (assumed uniform) :param store_as_jax: Store values as jax array instead of numpy array. Default is False (store as numpy). .. py:method:: fft2(pv) Alias method for FFT of PV .. py:method:: fft2_2dto1d(pv) Runs FFT then maps from 2D to 1D .. py:method:: ifft2(pvspec) Alias method for inverse FFT of PV Spectral .. py:method:: ifft2_2dto1d(pvspec) Runs inverse FFT then maps from 2D to 1D .. py:method:: integrate(f, x0, t_final, delta_t = None, include_x0 = True, t = None, **kwargs) Advances pv forward number of timesteps given by self.n_steps. .. note:: If pv not specified, use pvspec instance variable. :param f: right hand side (rhs) of the ODE. Not used, but needed to function with generate() from _data.Data(). :type f: function :param x0: potential vorticity (pvspec) initial condition in spectral space :type x0: ndarray .. py:method:: map1dto2d(x) Maps 1D state vector to 2D PV .. py:method:: map1dto2d_fft2(x) Maps for 1D to 2D then runs FFT .. py:method:: map1dto2d_ifft2(x) Maps for 1D to 2D then runs inverse FFT .. py:method:: map2dto1d(pv) Maps 2D PV to 1D system state .. py:method:: rhs(pvspec) Computes pv deriv on z=0, inverts pv to get streamfunction.