Source code for qgs.params.params

"""
    The model's parameters module
    =============================

    This module defines the main classes containing the model configuration parameters.
    The parameters are typically specified as :class:`~.params.parameter.Parameter` objects.

    There are seven types of parameters arranged in classes:

    * :class:`ScaleParams` contains the model scale parameters. These parameters are used to scale and
      `nondimentionalize`_ the :class:`~.params.parameter.Parameter` of the other parameters classes according to
      their :attr:`~.params.parameter.Parameter.units` attribute.
    * :class:`AtmosphericParams` contains the atmospheric dynamical parameters.
    * :class:`AtmosphericTemperatureParams` containing the atmosphere's temperature and heat-exchange parameters.
    * :class:`OceanicParams` contains the oceanic dynamical parameters.
    * :class:`OceanicTemperatureParams` contains the ocean's temperature and heat-exchange parameters.
    * :class:`GroundParams` contains the ground dynamical parameters (e.g. orography).
    * :class:`GroundTemperatureParams` contains the ground's temperature and heat-exchange parameters.

    These parameters classes are regrouped into a global structure :class:`QgParams` which also contains

    * spectral modes definition of the model
    * physical constants
    * parameters derived from the ones provided by the user
    * helper functions to initialize and parameterize the model

    This global parameters structure is used by the other modules to construct the model's ordinary differential
    equations.

    Warning
    -------

    If a model's parameter is set to `None`, it is assumed to be disabled.


    ---------------------

    Description of the classes
    --------------------------

    .. _nondimentionalize: https://en.wikipedia.org/wiki/Nondimensionalization
"""

import numpy as np
import pickle
import warnings
from abc import ABC

from qgs.params.parameter import Parameter
from qgs.basis.fourier import contiguous_channel_basis, contiguous_basin_basis
from qgs.basis.fourier import ChannelFourierBasis, BasinFourierBasis

from sympy import simplify

# TODO: - store model version in a variable somewhere
#       - force the user to define the aspect ratio n at parameter object instantiation


[docs]class Params(ABC): """Base class for a model's parameters container. Parameters ---------- dic: dict(float or Parameter), optional A dictionary with the parameters names and values to be assigned. """ _name = "" def __init__(self, dic=None): self.set_params(dic)
[docs] def set_params(self, dic): """Set the specified parameters values. Parameters ---------- dic: dict(float or Parameter) A dictionary with the parameters names and values to be assigned. """ if dic is not None: for key, val in zip(dic.keys(), dic.values()): if key in self.__dict__.keys(): if isinstance(self.__dict__[key], Parameter): if isinstance(val, Parameter): self.__dict__[key] = val else: d = self.__dict__[key].__dict__ self.__dict__[key] = Parameter(val, input_dimensional=d['_input_dimensional'], units=d['_units'], description=d['_description'], scale_object=d['_scale_object'], return_dimensional=d['_return_dimensional']) else: self.__dict__[key] = val
def __str__(self): s = "" for key, val in zip(self.__dict__.keys(), self.__dict__.values()): if 'params' not in key and key[0] != '_': if val is None: pass elif isinstance(val, Parameter): if val.input_dimensional: units = val.units efval = val.dimensional_value else: efval = val.nondimensional_value if val.nondimensional_value == val.dimensional_value: units = "" else: units = "[nondim]" s += "'" + key + "': " + str(efval) + " " + units + " (" + val.description + "),\n" elif isinstance(val, (np.ndarray, list, tuple)) and isinstance(val[0], Parameter): for i, v in enumerate(val): if v.input_dimensional: units = v.units efval = v.dimensional_value else: efval = v.nondimensional_value if v.nondimensional_value == v.dimensional_value: units = "" else: units = "[nondim]" s += "'" + key + "["+str(i+1)+"]': " + str(efval) + " " + units + " (" + v.description + "),\n" else: s += "'"+key+"': "+str(val)+",\n" return s def _list_params(self): return self._name+" Parameters:\n"+self.__str__()
[docs] def print_params(self): """Print the parameters contained in the container.""" print(self._list_params())
@staticmethod def create_params_array(values, input_dimensional=None, units=None, scale_object=None, description=None, return_dimensional=None): if hasattr(values, "__iter__"): ls = len(values) if not isinstance(input_dimensional, list): if input_dimensional is None: input_dimensional = True idx = ls * [input_dimensional] else: idx = input_dimensional if not isinstance(units, list): if units is None: units = "" u = ls * [units] else: u = units if not isinstance(description, list): if description is None: description = "" d = ls * [description] else: d = description if not isinstance(scale_object, list): s = ls * [scale_object] else: s = scale_object if not isinstance(return_dimensional, list): if return_dimensional is None: return_dimensional = False rd = ls * [return_dimensional] else: rd = return_dimensional arr = list() for i, val in enumerate(values): arr.append(Parameter(val, input_dimensional=idx[i], units=u[i], scale_object=s[i], description=d[i], return_dimensional=rd[i])) else: arr = values * [Parameter(0.e0, input_dimensional=input_dimensional, units=units, scale_object=scale_object, description=description, return_dimensional=return_dimensional)] return np.array(arr, dtype=object) def __repr__(self): s = super(Params, self).__repr__()+"\n"+self._list_params() return s
[docs] def load_from_file(self, filename, **kwargs): """Function to load previously saved Params object with the method :meth:`save_to_file`. Parameters ---------- filename: str The file name where the Params object was saved. kwargs: dict Keyword arguments to pass to the :mod:`pickle` module method. """ f = open(filename, 'rb') tmp_dict = pickle.load(f, **kwargs) f.close() self.__dict__.clear() self.__dict__.update(tmp_dict)
[docs] def save_to_file(self, filename, **kwargs): """Function to save the Params object to a file with the :mod:`pickle` module. Parameters ---------- filename: str The file name where to save the Params object. kwargs: dict Keyword arguments to pass to the :mod:`pickle` module method. """ f = open(filename, 'wb') pickle.dump(self.__dict__, f, **kwargs) f.close()
[docs]class ScaleParams(Params): """Class containing the model scales parameters. Parameters ---------- dic: dict(float or Parameter), optional A dictionary with the parameters names and values to be assigned. Attributes ---------- scale: Parameter The characteristic meridional space scale, :math:`L_y = \\pi \\, L`, in meters [:math:`m`]. f0: Parameter Coriolis parameter, in [:math:`s^{-1}`]. n: Parameter Model domain aspect ratio, :math:`n = 2 L_y/L_x` . rra: Parameter Earth radius, in meters [:math:`m`]. phi0_npi: Parameter Latitude expressed in fraction of :math:`\\pi` . deltap: Parameter Difference of pressure between the center of the two atmospheric layers, in [:math:`Pa`]. """ _name = "Scale" def __init__(self, dic=None): Params.__init__(self, dic) # ----------------------------------------------------------- # Scale parameters for the ocean and the atmosphere # ----------------------------------------------------------- self.scale = Parameter(5.e6, units='[m]', description="characteristic space scale (L*pi)", return_dimensional=True) self.f0 = Parameter(1.032e-4, units='[s^-1]', description="Coriolis parameter at the middle of the domain", return_dimensional=True) self.n = Parameter(1.3e0, input_dimensional=False, description="aspect ratio (n = 2 L_y / L_x)") self.rra = Parameter(6370.e3, units='[m]', description="earth radius", return_dimensional=True) self.phi0_npi = Parameter(0.25e0, input_dimensional=False, description="latitude expressed in fraction of pi") self.deltap = Parameter(5.e4, units='[Pa]', description='pressure difference between the two atmospheric layers', return_dimensional=True) self.Ha = Parameter(8500., units='[m]', description="Average height of the 500 hPa pressure level at midlatitude", return_dimensional=True) self.set_params(dic) # ---------------------------------------- # Some derived parameters (Domain, beta) # ---------------------------------------- @property def L(self): """Parameter: Typical length scale :math:`L` of the model, in meters [:math:`m`].""" return Parameter(self.scale / np.pi, units=self.scale.units, description='Typical length scale L', return_dimensional=True) @property def L_y(self): """Parameter: The meridional extent :math:`L_y = \\pi \\, L` of the model's domain, in meters [:math:`m`].""" return Parameter(self.scale, units=self.scale.units, description='The meridional extent of the model domain', return_dimensional=True) @property def L_x(self): """Parameter: The zonal extent :math:`L_x = 2 \\pi \\, L / n` of the model's domain, in meters [:math:`m`].""" return Parameter(2 * self.scale / self.n, units=self.scale.units, description='The zonal extent of the model domain', return_dimensional=True) @property def phi0(self): """Parameter: The reference latitude :math:`\\phi_0` at the center of the domain, expressed in radians [:math:`rad`].""" return Parameter(self.phi0_npi * np.pi, units='[rad]', description="The reference latitude of the center of the domain", return_dimensional=True) @property def beta(self): """Parameter: The meridional gradient of the Coriolis parameter at :math:`\\phi_0`, expressed in [:math:`m^{-1} s^{-1}`]. """ return Parameter(self.L / self.rra * np.cos(self.phi0) / np.sin(self.phi0), input_dimensional=False, units='[m^-1][s^-1]', scale_object=self, description="Meridional gradient of the Coriolis parameter at phi_0")
[docs]class AtmosphericParams(Params): """Class containing the atmospheric parameters. Parameters ---------- scale_params: ScaleParams The scale parameters object of the model. dic: dict(float or Parameter), optional A dictionary with the parameters names and values to be assigned. Attributes ---------- kd: Parameter Atmosphere bottom friction coefficient [:math:`s^{-1}`]. kdp: Parameter Atmosphere internal friction coefficient [:math:`s^{-1}`]. sigma: Parameter Static stability of the atmosphere [:math:`[m^2 s^{-2} Pa^{-2}`]. """ _name = "Atmospheric" def __init__(self, scale_params, dic=None): Params.__init__(self, dic) self._scale_params = scale_params # Parameters for the atmosphere self.kd = Parameter(0.1, input_dimensional=False, scale_object=scale_params, units='[s^-1]', description="atmosphere bottom friction coefficient") self.kdp = Parameter(0.01, input_dimensional=False, scale_object=scale_params, units='[s^-1]', description="atmosphere internal friction coefficient") self.sigma = Parameter(0.2e0, input_dimensional=False, scale_object=scale_params, units='[m^2][s^-2][Pa^-2]', description="static stability of the atmosphere") self.set_params(dic) @property def sig0(self): """Parameter: Static stability of the atmosphere divided by 2.""" return Parameter(self.sigma / 2, input_dimensional=False, scale_object=self._scale_params, units='[m^2][s^-2][Pa^-2]', description="0.5 * static stability of the atmosphere")
[docs]class AtmosphericTemperatureParams(Params): """Class containing the atmospheric temperature parameters. Parameters ---------- scale_params: ScaleParams The scale parameters object of the model. dic: dict(float or Parameter), optional A dictionary with the parameters names and values to be assigned. Attributes ---------- hd: None or Parameter Newtonian cooling coefficient. Newtonian cooling is disabled if `None`. thetas: None or ~numpy.ndarray(float) Coefficients of the Newtonian cooling spectral decomposition (non-dimensional). Newtonian cooling is disabled if `None`. gamma: None or Parameter Specific heat capacity of the atmosphere [:math:`J m^{-2} K^{-1}`]. Heat exchange scheme is disabled if `None`. C: None or ~numpy.ndarray(Parameter) Spectral decomposition of the constant short-wave radiation of the atmosphere [:math:`W m^{-2}`]. Heat exchange scheme is disabled if `None`. eps: None or Parameter Emissivity coefficient for the grey-body atmosphere Heat exchange scheme is disabled if `None`. T0: None or Parameter Stationary solution for the 0-th order atmospheric temperature [:math:`K`]. Heat exchange scheme is disabled if `None`. sc: None or Parameter Ratio of surface to atmosphere temperature Heat exchange scheme is disabled if `None`. hlambda: None or Parameter Sensible + turbulent heat exchange between ocean/ground and atmosphere [:math:`W m^{-2} K^{-1}`]. Heat exchange scheme is disabled if `None`. """ _name = "Atmospheric Temperature" def __init__(self, scale_params, dic=None): Params.__init__(self, dic) self._scale_params = scale_params self.hd = Parameter(0.045, input_dimensional=False, units='[s]', scale_object=scale_params, description="Newtonian cooling coefficient") self.thetas = None # Radiative equilibrium mean temperature decomposition on the model's modes self.gamma = None self.C = None self.eps = None self.T0 = None self.sc = None self.hlambda = None self.set_params(dic)
[docs] def set_insolation(self, value, pos=None): """Function to define the spectral decomposition of the constant short-wave radiation of the atmosphere (insolation) :math:`C_{{\\rm a}, i}` (:attr:`~.AtmosphericTemperatureParams.C`). Parameters ---------- value: float, int or iterable Value to set. If a scalar is given, the `pos` parameter should be provided to indicate which component to set. If an iterable is provided, create a vector of spectral decomposition parameters corresponding to it. pos: int, optional Indicate in which component to set the `value`. """ # TODO: - check for the dimensionality of the arguments if isinstance(value, (float, int)) and pos is not None and self.C is not None: self.C[pos] = Parameter(value, units='[W][m^-2]', scale_object=self._scale_params, description="spectral component "+str(pos+1)+" of the short-wave radiation of the atmosphere", return_dimensional=True) elif hasattr(value, "__iter__"): self._create_insolation(value) else: warnings.warn('A scalar value was provided, but without the `pos` argument indicating in which ' + 'component of the spectral decomposition to put it: Spectral decomposition unchanged !' + 'Please specify it or give a vector as `value`.')
def _create_insolation(self, values): if hasattr(values, "__iter__"): dim = len(values) else: dim = values values = dim * [0.] d = ["spectral component "+str(pos+1)+" of the short-wave radiation of the atmosphere" for pos in range(dim)] self.C = self.create_params_array(values, units='[W][m^-2]', scale_object=self._scale_params, description=d, return_dimensional=True)
[docs] def set_thetas(self, value, pos=None): """Function to define the spectral decomposition of the Newtonian cooling :math:`\\theta^\\star` (:attr:`~.AtmosphericTemperatureParams.thetas`). Parameters ---------- value: float, int or iterable Value to set. If a scalar is given, the `pos` parameter should be provided to indicate which component to set. If an iterable is provided, create a vector of spectral decomposition parameters corresponding to it. pos: int, optional Indicate in which component to set the `value`. """ # TODO: - check for the dimensionality of the arguments if isinstance(value, (float, int)) and pos is not None and self.thetas is not None: self.thetas[pos] = Parameter(value, scale_object=self._scale_params, description="spectral components "+str(pos+1)+" of the temperature profile", return_dimensional=False, input_dimensional=False) elif hasattr(value, "__iter__"): self._create_thetas(value) else: warnings.warn('A scalar value was provided, but without the `pos` argument indicating in which ' + 'component of the spectral decomposition to put it: Spectral decomposition unchanged !' + 'Please specify it or give a vector as `value`.')
def _create_thetas(self, values): if hasattr(values, "__iter__"): dim = len(values) else: dim = values values = dim * [0.] d = ["spectral component "+str(pos+1)+" of the temperature profile" for pos in range(dim)] self.thetas = self.create_params_array(values, scale_object=self._scale_params, description=d, return_dimensional=False, input_dimensional=False)
[docs]class OceanicParams(Params): """Class containing the oceanic parameters Parameters ---------- scale_params: ScaleParams The scale parameters object of the model. dic: dict(float or Parameter), optional A dictionary with the parameters names and values to be assigned. Attributes ---------- gp: Parameter Reduced gravity in [:math:`m \\, s^{-2}`]. r: Parameter Friction coefficient at the bottom of the ocean in [:math:`s^{-1}`]. h: Parameter Depth of the water layer of the ocean, in meters [:math:`m`]. d: Parameter The strength of the ocean-atmosphere mechanical coupling in [:math:`s^{-1}`]. """ _name = "Oceanic" def __init__(self, scale_params, dic=None): Params.__init__(self, dic) self._scale_params = scale_params self.gp = Parameter(3.1e-2, units='[m][s^-2]', return_dimensional=True, scale_object=scale_params, description='reduced gravity') self.r = Parameter(1.e-8, units='[s^-1]', scale_object=scale_params, description="frictional coefficient at the bottom of the ocean") self.h = Parameter(5.e2, units='[m]', return_dimensional=True, scale_object=scale_params, description="depth of the water layer of the ocean") self.d = Parameter(1.e-8, units='[s^-1]', scale_object=scale_params, description="strength of the ocean-atmosphere mechanical coupling") self.set_params(dic)
[docs]class OceanicTemperatureParams(Params): """Class containing the oceanic temperature parameters Parameters ---------- scale_params: ScaleParams The scale parameters object of the model. dic: dict(float or Parameter), optional A dictionary with the parameters names and values to be assigned. Attributes ---------- gamma: None or Parameter Specific heat capacity of the ocean [:math:`J m^{-2} K^{-1}`]. Heat exchange scheme is disabled if `None`. C: None or ~numpy.ndarray(Parameter) Spectral Decomposition of the constant short-wave radiation of the ocean [:math:`W m^{-2}`]. Heat exchange scheme is disabled if `None`. T0: None or Parameter Stationary solution for the 0-th order oceanic temperature [:math:`K`]. Heat exchange scheme is disabled if `None`. """ _name = "Oceanic Temperature" def __init__(self, scale_params, dic=None): Params.__init__(self, dic) self._scale_params = scale_params self.gamma = Parameter(2.e8, units='[J][m^-2][K^-1]', scale_object=scale_params, return_dimensional=True, description='specific heat capacity of the ocean') self.C = None self.T0 = None self.set_params(dic)
[docs] def set_insolation(self, value, pos=None): """Function to define the spectral decomposition of the constant short-wave radiation of the ocean (insolation) :math:`C_{{\\rm o}, i}` (:attr:`~.OceanicTemperatureParams.C`). Parameters ---------- value: float, int or iterable Value to set. If a scalar is given, the `pos` parameter should be provided to indicate which component to set. If an iterable is provided, create a vector of spectral decomposition parameters corresponding to it. pos: int, optional Indicate in which component to set the `value`. """ if isinstance(value, (float, int)) and pos is not None and self.C is not None: self.C[pos] = Parameter(value, units='[W][m^-2]', scale_object=self._scale_params, description="spectral component "+str(pos)+" of the short-wave radiation of the ocean", return_dimensional=True) elif hasattr(value, "__iter__"): self._create_insolation(value) else: warnings.warn('A scalar value was provided, but without the `pos` argument indicating in which ' + 'component of the spectral decomposition to put it: Spectral decomposition unchanged !' + 'Please specify it or give a vector as `value`.')
def _create_insolation(self, values): if hasattr(values, "__iter__"): dim = len(values) else: dim = values values = dim * [0.] d = ["spectral component "+str(pos)+" of the short-wave radiation of the ocean" for pos in range(dim)] self.C = self.create_params_array(values, units='[W][m^-2]', scale_object=self._scale_params, description=d, return_dimensional=True)
[docs]class GroundParams(Params): """Class containing the ground parameters Parameters ---------- scale_params: ScaleParams The scale parameters object of the model. dic: dict(float or Parameter), optional A dictionary with the parameters names and values to be assigned. Attributes ---------- hk: None or ~numpy.ndarray(float) Orography spectral decomposition coefficients (non-dimensional), an array of shape (:attr:`~QgParams.nmod` [0],). Orography is disabled (flat) if `None`. orographic_basis: str String to select which component basis modes to use to develop the orography in series. Can be either 'atmospheric' or 'ground'. Default to 'atmospheric'. """ _name = "Ground" def __init__(self, scale_params, dic=None): Params.__init__(self, dic) self._scale_params = scale_params self.hk = None # spectral orography coefficients self.orographic_basis = "atmospheric" self.set_params(dic)
[docs] def set_orography(self, value, pos=None, basis="atmospheric"): """Function to define the spectral decomposition of the orography profile :math:`h_k` (:attr:`~.GroundParams.hk`). Parameters ---------- value: float, int or iterable Value to set. If a scalar is given, the `pos` parameter should be provided to indicate which component to set. If an iterable is provided, create a vector of spectral decomposition parameters corresponding to it. pos: int, optional Indicate in which component to set the `value`. basis: str, optional Indicate which basis should be used to decompose the orography. Can be either `atmospheric`, `oceanic` or `ground`. Default to `atmospheric`. """ # TODO: - check for the dimensionality of the arguments # - check that inner products are symbolic if basis is not 'atmospheric' self.orographic_basis = basis if isinstance(value, (float, int)) and pos is not None and self.hk is not None: self.hk[pos] = Parameter(value, scale_object=self._scale_params, description="spectral components "+str(pos+1)+" of the orography", return_dimensional=False, input_dimensional=False) elif hasattr(value, "__iter__"): self._create_orography(value) else: warnings.warn('A scalar value was provided, but without the `pos` argument indicating in which ' + 'component of the spectral decomposition to put it: Spectral decomposition unchanged !' + 'Please specify it or give a vector as `value`.')
def _create_orography(self, values): if hasattr(values, "__iter__"): dim = len(values) else: dim = values values = dim * [0.] d = ["spectral component "+str(pos+1)+" of the orography" for pos in range(dim)] self.hk = self.create_params_array(values, scale_object=self._scale_params, description=d, return_dimensional=False, input_dimensional=False)
[docs]class GroundTemperatureParams(Params): """Class containing the ground temperature parameters Parameters ---------- scale_params: ScaleParams The scale parameters object of the model. dic: dict(float or Parameter), optional A dictionary with the parameters names and values to be assigned. Attributes ---------- gamma: None or Parameter Specific heat capacity of the ground [:math:`J m^{-2} K^{-1}`]. Heat exchange scheme is disabled if `None`. C: None or ~numpy.ndarray(Parameter) Spectral decomposition of the constant short-wave radiation of the ground [:math:`W m^{-2}`]. Heat exchange scheme is disabled if `None`. T0: None or Parameter Stationary solution for the 0-th order ground temperature [:math:`K`]. Heat exchange scheme is disabled if `None`. """ _name = "Ground Temperature" def __init__(self, scale_params, dic=None): Params.__init__(self, dic) self._scale_params = scale_params self.gamma = Parameter(2.e8, units='[J][m^-2][K^-1]', scale_object=scale_params, return_dimensional=True, description='specific heat capacity of the ground') self.C = None self.T0 = None self.set_params(dic)
[docs] def set_insolation(self, value, pos=None): """Function to define the decomposition of the constant short-wave radiation of the ground (insolation) :math:`C_{{\\rm g}, i}` (:attr:`~.GroundTemperatureParams.C`). Parameters ---------- value: float, int or iterable Value to set. If a scalar is given, the `pos` parameter should be provided to indicate which component to set. If an iterable is provided, create a vector of spectral decomposition parameters corresponding to it. pos: int, optional Indicate in which component to set the `value`. """ # TODO: - check for the dimensionality of the arguments if isinstance(value, (float, int)) and pos is not None and self.C is not None: self.C[pos] = Parameter(value, units='[W][m^-2]', scale_object=self._scale_params, description="spectral component "+str(pos+1)+" of the short-wave radiation of the ground", return_dimensional=True) elif hasattr(value, "__iter__"): self._create_insolation(value) else: warnings.warn('A scalar value was provided, but without the `pos` argument indicating in which ' + 'component of the spectral decomposition to put it: Spectral decomposition unchanged !' + 'Please specify it or give a vector as `value`.')
def _create_insolation(self, values): if hasattr(values, "__iter__"): dim = len(values) else: dim = values values = dim * [0.] d = ["spectral component "+str(pos+1)+" of the short-wave radiation of the ground" for pos in range(dim)] self.C = self.create_params_array(values, units='[W][m^-2]', scale_object=self._scale_params, description=d, return_dimensional=True)
[docs]class QgParams(Params): """General qgs parameters container. Parameters ---------- dic: dict(float or Parameter), optional A dictionary with the parameters names and values to be assigned. scale_params: None or ScaleParams, optional Scale parameters instance. If `None`, create a new ScaleParams instance. Default to None. Default to `None`. atmospheric_params: bool, None or AtmosphericParams, optional Atmospheric parameters instance. If 'True`, create a new AtmosphericParams instance. If `None`, atmospheric parameters are disabled. Default to `True`. atemperature_params: bool, None or AtmosphericTemperatureParams, optional Atmospheric temperature parameters instance. If 'True`, create a new AtmosphericTemperatureParams instance. If `None`, atmospheric temperature parameters are disabled. Default to `True`. oceanic_params: bool, None or OceanicParams, optional Oceanic parameters instance. If 'True`, create a new OceanicParams instance. If `None`, oceanic parameters are disabled. Default to `None`. otemperature_params: bool, None or OceanicTemperatureParams, optional Oceanic temperature parameters instance. If 'True`, create a new OceanicTemperatureParams instance. If `None`, oceanic temperature parameters are disabled. Default to `None`. ground_params: bool, None or GroundParams, optional Ground parameters instance. If 'True`, create a new GroundParams instance. If `None`, ground parameters are disabled. Default to `True`. gtemperature_params: bool, None or GroundTemperatureParams, optional Ground temperature parameters instance. If 'True`, create a new GroundTemperatureParams instance. If `None`, ground temperature parameters are disabled. Default to `None`. dynamic_T: bool, optional Whether to use a fixed or a dynamical reference temperature if the heat exchange scheme is activated. Default to `False`. T4: bool, optional Use or not the :math:`T^4` forcing for the evolution of the temperature field if the heat exchange is activated. Activate also the dynamical 0-th temperature. Default to `False`. Attributes ---------- scale_params: ScaleParams Scale parameters instance. atmospheric_params: None or AtmosphericParams Atmospheric parameters instance. If `None`, atmospheric parameters are disabled. atemperature_params: None or AtmosphericTemperatureParams Atmospheric temperature parameters instance. If `None`, atmospheric temperature parameters are disabled. oceanic_params: None or OceanicParams Oceanic parameters instance. If `None`, oceanic parameters are disabled. ground_params: None or GroundParams Ground parameters instance If `None`, ground parameters are disabled. gotemperature_params: None, OceanicTemperatureParams or GroundTemperatureParams Ground or Oceanic temperature parameters instance. If `None`, ground and oceanic temperature parameters are disabled. time_unit: float Dimensional unit of time to be used to represent the data. rr: Parameter `Gas constant`_ of `dry air`_ in [:math:`J \\, kg^{-1} \\, K^{-1}`]. sb: float `Stefan-Boltzmann constant`_ in [:math:`J \\, m^{-2} \\, s^{-1} \\, K^{-4}`]. dynamic_T: bool Whether to use a fixed or a dynamical reference temperature if the heat exchange scheme is activated. The atmospheric and possibly oceanic (or ground) basis must be reset if this parameter is changed. T4: bool Use or not the :math:`T^4` forcing for the evolution of the temperature field if the heat exchange is activated. .. _Gas constant: https://en.wikipedia.org/wiki/Gas_constant .. _dry air: https://en.wikipedia.org/wiki/Gas_constant#Specific_gas_constant .. _Stefan-Boltzmann constant: https://en.wikipedia.org/wiki/Stefan%E2%80%93Boltzmann_constant """ _name = "General" def __init__(self, dic=None, scale_params=None, atmospheric_params=True, atemperature_params=True, oceanic_params=None, otemperature_params=None, ground_params=True, gtemperature_params=None, dynamic_T=False, T4=False): Params.__init__(self, dic) # General scale parameters object (Mandatory param block) if scale_params is None: self.scale_params = ScaleParams(dic) else: self.scale_params = scale_params # Atmospheric parameters object if atmospheric_params is True: self.atmospheric_params = AtmosphericParams(self.scale_params, dic=dic) else: self.atmospheric_params = atmospheric_params # Atmospheric temperature parameters object if atmospheric_params is True: self.atemperature_params = AtmosphericTemperatureParams(self.scale_params, dic=dic) else: self.atemperature_params = atemperature_params if oceanic_params is True: self.oceanic_params = OceanicParams(self.scale_params, dic) else: self.oceanic_params = oceanic_params if ground_params is True: self.ground_params = GroundParams(self.scale_params, dic) else: self.ground_params = ground_params if otemperature_params is True: self.gotemperature_params = OceanicTemperatureParams(self.scale_params, dic) else: self.gotemperature_params = otemperature_params if gtemperature_params is True: self.gotemperature_params = GroundTemperatureParams(self.scale_params, dic) else: self.gotemperature_params = gtemperature_params self._atmospheric_basis = None self._oceanic_basis = None self._ground_basis = None self._number_of_atmospheric_modes = 0 self._number_of_oceanic_modes = 0 self._number_of_ground_modes = 0 self._ams = None self._oms = None self._gms = None self.dynamic_T = dynamic_T self.T4 = T4 # Force dynamic temperatures if T4 tendencies are activated if T4: self.dynamic_T = T4 self._atmospheric_latex_var_string = list() self._atmospheric_var_string = list() self._oceanic_latex_var_string = list() self._oceanic_var_string = list() self._ground_latex_var_string = list() self._ground_var_string = list() self._components_units = [r'm$^2$s$^{-1}$', r'K', r'm$^2$s$^{-1}$', r'K'] self.time_unit = 'days' # Physical constants self.rr = Parameter(287.058e0, return_dimensional=True, units='[J][kg^-1][K^-1]', scale_object=self.scale_params, description="gas constant of dry air") self.sb = Parameter(5.67e-8, return_dimensional=True, units='[J][m^-2][s^-1][K^-4]', scale_object=self.scale_params, description="Stefan-Boltzmann constant") self.set_params(dic) # ----------------------------------------------------------- # Derived Quantities (Parameters) # ----------------------------------------------------------- @property def LR(self): """float: Reduced Rossby deformation radius :math:`L_\\text{R} = \\sqrt{g' \\, h } / f_0` .""" op = self.oceanic_params scp = self.scale_params if op is not None: try: return np.sqrt(op.gp * op.h) / scp.f0 except: return None else: return None @property def G(self): """float: The :math:`G = - L^2/L_R^2` parameter.""" scp = self.scale_params if self.LR is not None: try: return -scp.L**2 / self.LR**2 except: return None else: return None @property def Cpgo(self): """float: The :math:`C\'_{{\\rm g/\\rm o},i} = R C_{{\\rm g/\\rm o},i} / (\\gamma_{\\rm g/\\rm o} L^2 f_0^3)` parameter.""" gotp = self.gotemperature_params scp = self.scale_params if gotp is not None: try: return gotp.C / (gotp.gamma * scp.f0) * self.rr / (scp.f0 ** 2 * scp.L ** 2) except: return None else: return None @property def Lpgo(self): """float: The :math:`\\lambda\'_{{\\rm g/\\rm o}} = \\lambda/(\\gamma_{\\rm g/\\rm o} f_0)` parameter.""" atp = self.atemperature_params gotp = self.gotemperature_params scp = self.scale_params if atp is not None and gotp is not None: try: return atp.hlambda / (gotp.gamma * scp.f0) except: return None else: return None @property def Cpa(self): """float: The :math:`C\'_{{\\rm a},i} = R C_{{\\rm a},i} / (2 \\gamma_{\\rm a} L^2 f_0^3)` parameter.""" atp = self.atemperature_params scp = self.scale_params if atp is not None: try: return atp.C / (atp.gamma * scp.f0) * self.rr / (scp.f0 ** 2 * scp.L ** 2) / 2 except: return None else: return None @property def Lpa(self): """float: The :math:`\\lambda\'_{\\rm a} = \\lambda / (\\gamma_{\\rm a} f_0)` parameter.""" atp = self.atemperature_params scp = self.scale_params if atp is not None: try: return atp.hlambda / (atp.gamma * scp.f0) except: return None else: return None @property def sbpgo(self): """float: Long wave radiation lost by ground/ocean to the atmosphere :math:`s_{B,{\\rm g/\\rm o}} = 4\\,\\sigma_B \\, T_{{\\rm a},0}^3 / (\\gamma_{\\rm g/\\rm o} f_0)` in the linearized temperature model equations.""" gotp = self.gotemperature_params scp = self.scale_params if gotp is not None and not self.dynamic_T: try: return 4 * self.sb * gotp.T0 ** 3 / (gotp.gamma * scp.f0) except: return None else: return None @property def sbpa(self): """float: Long wave radiation from atmosphere absorbed by ground/ocean :math:`s_{B,{\\rm a}} = 4\\,\\epsilon_{\\rm a}\\, \\sigma_B \\, T_{{\\rm a},0}^3 / (\\gamma_{\\rm g/\\rm o} f_0)` in the linearized temperature model equations.""" atp = self.atemperature_params gotp = self.gotemperature_params scp = self.scale_params if gotp is not None and atp is not None and not self.dynamic_T: try: return 8 * atp.eps * self.sb * atp.T0 ** 3 / (gotp.gamma * scp.f0) except: return None else: return None @property def LSBpgo(self): """float: Long wave radiation from ground/ocean absorbed by atmosphere :math:`S_{B,{\\rm g/\\rm o}} = 2\\,\\epsilon_{\\rm a}\\, \\sigma_B \\, T_{{\\rm a},0}^3 / (\\gamma_{\\rm a} f_0)` in the linearized temperature model equations.""" atp = self.atemperature_params gotp = self.gotemperature_params scp = self.scale_params if atp is not None and gotp is not None and not self.dynamic_T: try: return 2 * atp.eps * self.sb * gotp.T0 ** 3 / (atp.gamma * scp.f0) except: return None else: return None @property def LSBpa(self): """float: Long wave radiation lost by atmosphere to space & ground/ocean :math:`S_{B,{\\rm a}} = 8\\,\\epsilon_{\\rm a}\\, \\sigma_B \\, T_{{\\rm a},0}^3 / (\\gamma_{\\rm a} f_0)` in the linearized temperature model equations.""" atp = self.atemperature_params scp = self.scale_params if atp is not None and not self.dynamic_T: try: return 8 * atp.eps * self.sb * atp.T0 ** 3 / (atp.gamma * scp.f0) except: return None else: return None @property def T4sbpgo(self): """float: Long wave radiation lost by ground/ocean to the atmosphere :math:`s_{B,{\\rm g/\\rm o}} = \\sigma_B \\, L^6 \\, f_0^5 / (\\gamma_{\\rm g/\\rm o} R^3)` in the :math:`T^4` model equations.""" gotp = self.gotemperature_params scp = self.scale_params if gotp is not None: try: return self.sb * scp.L ** 6 * scp.f0 ** 5 / (gotp.gamma * self.rr ** 3) except: return None else: return None @property def T4sbpa(self): """float: Long wave radiation from atmosphere absorbed by ground/ocean :math:`s_{B,{\\rm a}} = 16 \\,\\epsilon_{\\rm a}\\, \\sigma_B \\, L^6 \\, f_0^5 / (\\gamma_{\\rm g/\\rm o} R^3)` in the :math:`T^4` model equations.""" atp = self.atemperature_params gotp = self.gotemperature_params scp = self.scale_params if gotp is not None and atp is not None: try: return 16 * atp.eps * self.sb * scp.L ** 6 * scp.f0 ** 5 / (gotp.gamma * self.rr ** 3) except: return None else: return None @property def T4LSBpgo(self): """float: Long wave radiation from ground/ocean absorbed by atmosphere :math:`S_{B,{\\rm g/\\rm o}} = \\frac{1}{2} \\, \\epsilon_{\\rm a}\\, \\sigma_B \\, L^6 \\, f_0^5 / (\\gamma_{\\rm a} R^3)` in the :math:`T^4` model equations.""" atp = self.atemperature_params scp = self.scale_params if atp is not None: try: return 0.5 * atp.eps * self.sb * scp.L ** 6 * scp.f0 ** 5 / (atp.gamma * self.rr ** 3) except: return None else: return None @property def T4LSBpa(self): """float: Long wave radiation lost by atmosphere to space & ground/ocean :math:`S_{B,{\\rm a}} = 16 \\,\\epsilon_{\\rm a}\\, \\sigma_B \\, L^6 \\, f_0^5 / (\\gamma_{\\rm a} R^3)` in the :math:`T^4` model equations.""" atp = self.atemperature_params scp = self.scale_params if atp is not None: try: return 16 * atp.eps * self.sb * scp.L ** 6 * scp.f0 ** 5 / (atp.gamma * self.rr ** 3) except: return None else: return None # The following properties might be refactored if the unit system of the model get more widespread across modules. @property def streamfunction_scaling(self): """float: Dimensional scaling of the streamfunction fields.""" return self.scale_params.L**2 * self.scale_params.f0 @property def temperature_scaling(self): """float: Dimensional scaling of the temperature fields.""" return self.streamfunction_scaling * self.scale_params.f0 / self.rr @property def geopotential_scaling(self): """float: Dimensional scaling of the geopotential height.""" return self.scale_params.f0 / 9.81
[docs] def set_params(self, dic): """Set the specified parameters values. Parameters ---------- dic: dict(float or Parameter) A dictionary with the parameters names and values to be assigned. """ if dic is not None: for key, val in zip(dic.keys(), dic.values()): if key in self.__dict__.keys(): self.__dict__[key] = val if 'scale_params' in self.__dict__.keys(): self.scale_params.set_params(dic) if 'atmospheric_params' in self.__dict__.keys(): if self.atmospheric_params is not None: self.atmospheric_params.set_params(dic) if 'atemperature_params' in self.__dict__.keys(): if self.atemperature_params is not None: self.atemperature_params.set_params(dic) if 'oceanic_params' in self.__dict__.keys(): if self.oceanic_params is not None: self.oceanic_params.set_params(dic) if 'ground_params' in self.__dict__.keys(): if self.ground_params is not None: self.ground_params.set_params(dic) if 'otemperature_params' in self.__dict__.keys(): if self.gotemperature_params is not None: self.gotemperature_params.set_params(dic) if 'gtemperature_params' in self.__dict__.keys(): if self.gotemperature_params is not None: self.gotemperature_params.set_params(dic)
[docs] def print_params(self): """Print all the parameters in the container.""" s = self._list_params()+"\n" if 'scale_params' in self.__dict__.keys(): s += self.scale_params._list_params()+"\n" if 'atmospheric_params' in self.__dict__.keys(): if self.atmospheric_params is not None: s += self.atmospheric_params._list_params()+"\n" if 'atemperature_params' in self.__dict__.keys(): if self.atemperature_params is not None: s += self.atemperature_params._list_params()+"\n" if 'oceanic_params' in self.__dict__.keys(): if self.oceanic_params is not None: s += self.oceanic_params._list_params()+"\n" if 'ground_params' in self.__dict__.keys(): if self.ground_params is not None: s += self.ground_params._list_params()+"\n" if 'gotemperature_params' in self.__dict__.keys(): if self.gotemperature_params is not None: s += self.gotemperature_params._list_params() + "\n" print("Qgs v0.2.8 parameters summary") print("=============================\n") print(s)
@property def ndim(self): """int: Total number of variables of the model.""" return self.variables_range[-1] @property def nmod(self): """(int, int): Atmospheric and ground/oceanic number of modes.""" if self._number_of_oceanic_modes != 0: return [self._number_of_atmospheric_modes, self._number_of_oceanic_modes] else: return [self._number_of_atmospheric_modes, self._number_of_ground_modes] @property def var_string(self): """list(str): List of model's variable names.""" ls = list() for var in self._atmospheric_var_string+self._oceanic_var_string+self._ground_var_string: ls.append(var) return ls @property def latex_var_string(self): """list(str): List of model's variable names, ready for use in latex.""" ls = list() for var in self._atmospheric_latex_var_string+self._oceanic_latex_var_string+self._ground_latex_var_string: ls.append(r'{\ '[0:-1] + var + r'}') return ls @property def variables_range(self): """list(int): List of the variables indices upper bound per component.""" natm = self.nmod[0] ngoc = self.nmod[1] vr = list() vr.append(natm) vr.append(vr[-1] + natm) if self.dynamic_T: vr[-1] += 1 if ngoc > 0: vr.append(vr[-1] + ngoc) if self._oceanic_basis is not None: vr.append(vr[-1] + ngoc) if self.dynamic_T: vr[-1] += 1 return vr @property def number_of_variables(self): """list(int): Number of variables per component.""" vr = self.variables_range return [vr[0]] + [vr[i] - vr[i-1] for i in range(1, len(vr))] @property def latex_components_units(self): """list(str): The units of every model's components variables, as a list of latex strings.""" return self._components_units
[docs] def get_variable_units(self, i): """Return the units of a model's variable as a string containing latex symbols. Parameters ---------- i: int The number of the variable. Returns ------- str: The string with the units of the variable. """ if i >= self.ndim: warnings.warn("Variable " + str(i) + " doesn't exist, cannot return its units.") return None else: if i < self.variables_range[0]: return self._components_units[0] if self.variables_range[0] <= i < self.variables_range[1]: return self._components_units[1] if self.oceanic_basis is not None: if self.variables_range[1] <= i < self.variables_range[2]: return self._components_units[2] if self.variables_range[2] <= i < self.variables_range[3]: return self._components_units[3] if self.ground_basis is not None: if self.variables_range[1] <= i < self.variables_range[2]: return self._components_units[3]
@property def dimensional_time(self): """float: Return the conversion factor between the non-dimensional time and the dimensional time unit specified in :attr:`.time_unit`""" c = 24 * 3600 if self.time_unit == 'hours': c = 3600 if self.time_unit == 'days': c = 24 * 3600 if self.time_unit == 'years': c = 24 * 3600 * 365 return 1 / (self.scale_params.f0 * c) # ------------------------------------------------------------------- # Config setters to be used with symbolic inner products # ------------------------------------------------------------------- @property def atmospheric_basis(self): """Basis: The atmospheric basis of functions used to project the PDEs onto.""" return self._atmospheric_basis @atmospheric_basis.setter def atmospheric_basis(self, basis): self._ams = None self._oms = None self._gms = None self._atmospheric_basis = basis self._number_of_atmospheric_modes = len(basis.functions) if self.dynamic_T: self._atmospheric_basis.functions.insert(0, simplify("1")) if self.ground_params is not None and self.ground_params.orographic_basis == "atmospheric": self.ground_params.set_orography(self._number_of_atmospheric_modes * [0.e0]) if self.atemperature_params is not None: self.atemperature_params.set_thetas(self._number_of_atmospheric_modes * [0.e0]) @property def oceanic_basis(self): """Basis: The oceanic basis of functions used to project the PDEs onto.""" return self._oceanic_basis @oceanic_basis.setter def oceanic_basis(self, basis): self._ams = None self._oms = None self._gms = None self._oceanic_basis = basis self._number_of_ground_modes = 0 self._number_of_oceanic_modes = len(basis) if self.dynamic_T: self._oceanic_basis.functions.insert(0, simplify("1")) if self.atemperature_params is not None: # disable the Newtonian cooling self.atemperature_params.thetas = None self.atemperature_params.hd = None self.atemperature_params.gamma = Parameter(1.e7, units='[J][m^-2][K^-1]', scale_object=self.scale_params, description='specific heat capacity of the atmosphere', return_dimensional=True) if self.dynamic_T: self.atemperature_params.set_insolation((self.nmod[0] + 1) * [0.e0]) self.atemperature_params.set_insolation(100.0, 0) self.atemperature_params.set_insolation(100.0, 1) else: self.atemperature_params.set_insolation(self.nmod[0] * [0.e0]) self.atemperature_params.set_insolation(100.0, 0) self.atemperature_params.T0 = Parameter(270.0, units='[K]', scale_object=self.scale_params, return_dimensional=True, description="stationary solution for the 0-th order atmospheric temperature") self.atemperature_params.eps = Parameter(0.76e0, input_dimensional=False, description="emissivity coefficient for the grey-body atmosphere") self.atemperature_params.sc = Parameter(1., input_dimensional=False, description="ratio of surface to atmosphere temperature") self.atemperature_params.hlambda = Parameter(20.00, units='[W][m^-2][K^-1]', scale_object=self.scale_params, return_dimensional=True, description="sensible+turbulent heat exchange between ocean/ground and atmosphere") if self.gotemperature_params is not None: if self.dynamic_T: self.gotemperature_params.set_insolation((self.nmod[0] + 1) * [0.e0]) self.gotemperature_params.set_insolation(350.0, 0) self.gotemperature_params.set_insolation(350.0, 1) else: self.gotemperature_params.set_insolation(self.nmod[0] * [0.e0]) self.gotemperature_params.set_insolation(350.0, 0) self.gotemperature_params.T0 = Parameter(285.0, units='[K]', scale_object=self.scale_params, return_dimensional=True, description="stationary solution for the 0-th order oceanic temperature") # if setting an ocean, then disable the orography if self.ground_params is not None: self.ground_params.hk = None @property def ground_basis(self): """Basis: The ground basis of functions used to project the PDEs onto.""" return self._ground_basis @ground_basis.setter def ground_basis(self, basis): self._ams = None self._oms = None self._gms = None self._ground_basis = basis self._number_of_ground_modes = len(basis) self._number_of_oceanic_modes = 0 if self.dynamic_T: self._ground_basis.functions.insert(0, simplify("1")) if self.atemperature_params is not None: # disable the Newtonian cooling self.atemperature_params.thetas = None self.atemperature_params.hd = None self.atemperature_params.gamma = Parameter(1.e7, units='[J][m^-2][K^-1]', scale_object=self.scale_params, description='specific heat capacity of the atmosphere', return_dimensional=True) if self.dynamic_T: self.atemperature_params.set_insolation((self.nmod[0] + 1) * [0.e0]) self.atemperature_params.set_insolation(100.0, 0) self.atemperature_params.set_insolation(100.0, 1) else: self.atemperature_params.set_insolation(self.nmod[0] * [0.e0]) self.atemperature_params.set_insolation(100.0, 0) self.atemperature_params.T0 = Parameter(270.0, units='[K]', scale_object=self.scale_params, return_dimensional=True, description="stationary solution for the 0-th order atmospheric temperature") self.atemperature_params.eps = Parameter(0.76e0, input_dimensional=False, description="emissivity coefficient for the grey-body atmosphere") self.atemperature_params.sc = Parameter(1., input_dimensional=False, description="ratio of surface to atmosphere temperature") self.atemperature_params.hlambda = Parameter(20.00, units='[W][m^-2][K^-1]', scale_object=self.scale_params, return_dimensional=True, description="sensible+turbulent heat exchange between ocean/ground and atmosphere") if self.gotemperature_params is not None: # if orography is disabled, enable it! if self.ground_params is not None: if self.ground_params.hk is None: if self.ground_params.orographic_basis == 'atmospheric': self.ground_params.set_orography(self._number_of_atmospheric_modes * [0.e0]) else: self.ground_params.set_orography(self._number_of_ground_modes * [0.e0]) self.ground_params.set_orography(0.1, 1) if self.dynamic_T: self.gotemperature_params.set_insolation((self.nmod[0] + 1) * [0.e0]) self.gotemperature_params.set_insolation(350.0, 0) self.gotemperature_params.set_insolation(350.0, 1) else: self.gotemperature_params.set_insolation(self.nmod[0] * [0.e0]) self.gotemperature_params.set_insolation(350.0, 0) self.gotemperature_params.T0 = Parameter(285.0, units='[K]', scale_object=self.scale_params, return_dimensional=True, description="stationary solution for the 0-th order oceanic temperature")
[docs] def set_atmospheric_modes(self, basis, auto=False): """Function to configure the atmospheric modes (basis functions) used to project the PDEs onto. Parameters ---------- basis: Basis Basis object containing the definition of the atmospheric modes. auto: bool, optional Automatically instantiate the parameters container needed to describe the atmospheric models parameters. Default is False. Examples -------- >>> from qgs.params.params import QgParams >>> from qgs.basis.fourier import contiguous_channel_basis >>> q = QgParams() >>> atm_basis = contiguous_channel_basis(2, 2, 1.5) >>> q.set_atmospheric_modes(atm_basis) >>> q.atmospheric_basis [sqrt(2)*cos(y), 2*sin(y)*cos(n*x), 2*sin(y)*sin(n*x), sqrt(2)*cos(2*y), 2*sin(2*y)*cos(n*x), 2*sin(2*y)*sin(n*x), 2*sin(y)*cos(2*n*x), 2*sin(y)*sin(2*n*x), 2*sin(2*y)*cos(2*n*x), 2*sin(2*y)*sin(2*n*x)] """ if auto: if self.atemperature_params is None: self.atemperature_params = AtmosphericTemperatureParams(self.scale_params) if self.atmospheric_params is None: self.atmospheric_params = AtmosphericParams(self.scale_params) self.atmospheric_basis = basis self._atmospheric_latex_var_string = list() self._atmospheric_var_string = list() for i in range(1, self.nmod[0] + 1): self._atmospheric_latex_var_string.append(r'psi_{{\rm a},' + str(i) + "}") self._atmospheric_var_string.append(r'psi_a_' + str(i)) if self.dynamic_T: self._atmospheric_latex_var_string.append(r', T_{{\rm a},0}') self._atmospheric_var_string.append(r'T_a_0') for i in range(1, self.nmod[0] + 1): self._atmospheric_latex_var_string.append(r'theta_{{\rm a},' + str(i) + "}") self._atmospheric_var_string.append(r'theta_a_' + str(i))
[docs] def set_oceanic_modes(self, basis, auto=True): """Function to configure the oceanic modes (basis functions) used to project the PDEs onto. Parameters ---------- basis: Basis Basis object containing the definition of the oceanic modes. auto: bool, optional Automatically instantiate or not the parameters container needed to describe the oceanic models parameters. Default is True. Examples -------- >>> from qgs.params.params import QgParams >>> from qgs.basis.fourier import contiguous_channel_basis, contiguous_basin_basis >>> q = QgParams() >>> atm_basis = contiguous_channel_basis(2, 2, 1.5) >>> oc_basis = contiguous_basin_basis(2, 4, 1.5) >>> q.set_atmospheric_modes(atm_basis) >>> q.set_oceanic_modes(oc_basis) >>> q.oceanic_basis [2*sin(y)*sin(0.5*n*x), 2*sin(2*y)*sin(0.5*n*x), 2*sin(3*y)*sin(0.5*n*x), 2*sin(4*y)*sin(0.5*n*x), 2*sin(y)*sin(1.0*n*x), 2*sin(2*y)*sin(1.0*n*x), 2*sin(3*y)*sin(1.0*n*x), 2*sin(4*y)*sin(1.0*n*x)] """ if self._atmospheric_basis is None: # Presently, the ocean can not yet be set independently of an atmosphere. print('Atmosphere modes not set up. Add an atmosphere before adding an ocean!') print('Oceanic setup aborted.') return if auto: if self.gotemperature_params is None or isinstance(self.gotemperature_params, GroundTemperatureParams): self.gotemperature_params = OceanicTemperatureParams(self.scale_params) if self.oceanic_params is None: self.oceanic_params = OceanicParams(self.scale_params) self.ground_params = None self._ground_basis = None self.oceanic_basis = basis self._oceanic_latex_var_string = list() self._oceanic_var_string = list() self._ground_latex_var_string = list() self._ground_var_string = list() for i in range(1, self.nmod[1] + 1): self._oceanic_latex_var_string.append(r'psi_{\rm o,' + str(i) + "}") self._oceanic_var_string.append(r'psi_o_' + str(i)) if self.dynamic_T: self._oceanic_latex_var_string.append(r', T_{{\rm o},0}') self._oceanic_var_string.append(r'T_o_0') for i in range(1, self.nmod[1] + 1): self._oceanic_latex_var_string.append(r'delta T_{{\rm o},' + str(i) + "}") self._oceanic_var_string.append(r'delta_T_o_' + str(i))
[docs] def set_ground_modes(self, basis=None, auto=True): """Function to configure the ground modes (basis functions) used to project the PDEs onto. Parameters ---------- basis: None or Basis, optional Basis object containing the definition of the ground modes. If `None`, use the basis of the atmosphere. Default to `None`. auto: bool, optional Automatically instantiate or not the parameters container needed to describe the ground models parameters. Default is True. Examples -------- >>> from qgs.params.params import QgParams >>> from qgs.basis.fourier import contiguous_channel_basis, contiguous_basin_basis >>> q = QgParams() >>> atm_basis = contiguous_channel_basis(2, 2, 1.5) >>> q.set_atmospheric_modes(atm_basis) >>> q.set_ground_modes() >>> q.ground_basis [sqrt(2)*cos(y), 2*sin(y)*cos(n*x), 2*sin(y)*sin(n*x), sqrt(2)*cos(2*y), 2*sin(2*y)*cos(n*x), 2*sin(2*y)*sin(n*x), 2*sin(y)*cos(2*n*x), 2*sin(y)*sin(2*n*x), 2*sin(2*y)*cos(2*n*x), 2*sin(2*y)*sin(2*n*x)] """ if self._atmospheric_basis is None: # Presently, the ground can not yet be set independently of an atmosphere. print('Atmosphere modes not set up. Add an atmosphere before adding the ground!') print('Ground setup aborted.') return if auto: if self.gotemperature_params is None or isinstance(self.gotemperature_params, OceanicTemperatureParams): self.gotemperature_params = GroundTemperatureParams(self.scale_params) if self.ground_params is None: self.ground_params = GroundParams(self.scale_params) self.oceanic_params = None self._oceanic_basis = None if basis is not None: self.ground_basis = basis else: self.ground_basis = self._atmospheric_basis self._oceanic_var_string = list() self._oceanic_latex_var_string = list() self._ground_latex_var_string = list() self._ground_var_string = list() if self.dynamic_T: self._oceanic_latex_var_string.append(r', T_{{\rm g},0}') self._oceanic_var_string.append(r'T_g_0') for i in range(1, self.nmod[1] + 1): self._ground_latex_var_string.append(r'delta T_{\rm g,' + str(i) + "}") self._ground_var_string.append(r'delta_T_g_' + str(i))
# ------------------------------------------------------------------- # Specific basis setters # -------------------------------------------------------------------
[docs] def set_atmospheric_channel_fourier_modes(self, nxmax, nymax, auto=False, mode='analytic'): """Function to configure and set the basis for contiguous spectral blocks of atmospheric modes on a channel. Parameters ---------- nxmax: int Maximum x-wavenumber to fill the spectral block up to. nymax: int Maximum :math:`y`-wavenumber to fill the spectral block up to. auto: bool, optional Automatically instantiate the parameters container needed to describe the atmospheric models parameters. Default is `False`. mode: str, optional Mode to set the inner products: Either `analytic` or `symbolic`: `analytic` for inner products computed with formula or `symbolic` using `Sympy`_. Default to `analytic`. Examples -------- >>> from qgs.params.params import QgParams >>> q = QgParams() >>> q.set_atmospheric_channel_fourier_modes(2, 2) >>> q.ablocks array([[1, 1], [1, 2], [2, 1], [2, 2]]) .. _Sympy: https://www.sympy.org/ """ if mode == 'symbolic': basis = contiguous_channel_basis(nxmax, nymax, self.scale_params.n) self.set_atmospheric_modes(basis, auto) else: self._set_atmospheric_analytic_fourier_modes(nxmax, nymax, auto)
[docs] def set_oceanic_basin_fourier_modes(self, nxmax, nymax, auto=True, mode='analytic'): """Function to configure and set the basis for contiguous spectral blocks of oceanic modes on a closed basin. Parameters ---------- nxmax: int Maximum x-wavenumber to fill the spectral block up to. nymax: int Maximum :math:`y`-wavenumber to fill the spectral block up to. auto: bool, optional Automatically instantiate the parameters container needed to describe the atmospheric models parameters. Default is `True`. mode: str, optional Mode to set the inner products: Either `analytic` or `symbolic`. `analytic` for inner products computed with formula or `symbolic` using `Sympy`_. Default to `analytic`. Examples -------- >>> from qgs.params.params import QgParams >>> q = QgParams() >>> q.set_atmospheric_channel_fourier_modes(2, 2) >>> q.set_oceanic_basin_fourier_modes(2, 4) >>> q.oblocks array([[1, 1], [1, 2], [1, 3], [1, 4], [2, 1], [2, 2], [2, 3], [2, 4]]) .. _Sympy: https://www.sympy.org/ """ if mode == 'symbolic': basis = contiguous_basin_basis(nxmax, nymax, self.scale_params.n) self.set_oceanic_modes(basis, auto) else: self._set_oceanic_analytic_fourier_modes(nxmax, nymax, auto)
[docs] def set_ground_channel_fourier_modes(self, nxmax=None, nymax=None, auto=True, mode='analytic'): """Function to configure and set the basis for contiguous spectral blocks of ground modes on a channel. Parameters ---------- nxmax: int Maximum x-wavenumber to fill the spectral block up to. nymax: int Maximum :math:`y`-wavenumber to fill the spectral block up to. auto: bool, optional Automatically instantiate the parameters container needed to describe the atmospheric models parameters. Default is `True`. mode: str, optional Mode to set the inner products: Either `analytic` or `symbolic`. `analytic` for inner products computed with formula or `symbolic` using `Sympy`_. Default to `analytic`. Examples -------- >>> from qgs.params.params import QgParams >>> q = QgParams() >>> q.set_atmospheric_channel_fourier_modes(2,4) >>> q.set_ground_channel_fourier_modes() >>> q.gblocks array([[1, 1], [1, 2], [1, 3], [1, 4], [2, 1], [2, 2], [2, 3], [2, 4]]) .. _Sympy: https://www.sympy.org/ """ if mode == "symbolic": if nxmax is not None and nymax is not None: basis = contiguous_channel_basis(nxmax, nymax, self.scale_params.n) else: basis = None self.set_ground_modes(basis, auto) else: self._set_ground_analytic_fourier_modes(nxmax, nymax, auto)
# ------------------------------------------------------------------- # Model configs setter to be used with analytic inner products # ------------------------------------------------------------------- @property def ablocks(self): """~numpy.ndarray(int): Spectral blocks detailing the model's atmospheric modes :math:`x`- and :math:`y`-wavenumber. Array of shape (:attr:`~QgParams.nmod` [0], 2).""" return self._ams @ablocks.setter def ablocks(self, value): self._ams = value basis = ChannelFourierBasis(self._ams, self.scale_params.n) self._atmospheric_basis = basis namod = 0 for i in range(self.ablocks.shape[0]): if self.ablocks[i, 0] == 1: namod += 3 else: namod += 2 self._number_of_atmospheric_modes = namod if self.ground_params is not None: self.ground_params.orographic_basis = 'atmospheric' self.ground_params.set_orography(namod * [0.e0]) self.ground_params.set_orography(0.1, 1) if self.atemperature_params is not None: self.atemperature_params.set_thetas(namod * [0.e0]) self.atemperature_params.set_thetas(0.1, 0) @property def oblocks(self): """~numpy.ndarray(int): Spectral blocks detailing the model's oceanic modes :math:`x`-and :math:`y`-wavenumber. Array of shape (:attr:`~QgParams.nmod` [1], 2).""" return self._oms @oblocks.setter def oblocks(self, value): self._oms = value self._gms = None basis = BasinFourierBasis(self._oms, self.scale_params.n) self._oceanic_basis = basis self._ground_basis = None if self.atemperature_params is not None: # disable the Newtonian cooling self.atemperature_params.thetas = None # np.zeros(self.nmod[0]) self.atemperature_params.hd = None # Parameter(0.0, input_dimensional=False) self.atemperature_params.gamma = Parameter(1.e7, units='[J][m^-2][K^-1]', scale_object=self.scale_params, description='specific heat capacity of the atmosphere', return_dimensional=True) self.atemperature_params.set_insolation(self.nmod[0] * [0.e0]) self.atemperature_params.set_insolation(100.0, 0) self.atemperature_params.eps = Parameter(0.76e0, input_dimensional=False, description="emissivity coefficient for the grey-body atmosphere") self.atemperature_params.T0 = Parameter(270.0, units='[K]', scale_object=self.scale_params, return_dimensional=True, description="stationary solution for the 0-th order atmospheric temperature") self.atemperature_params.sc = Parameter(1., input_dimensional=False, description="ratio of surface to atmosphere temperature") self.atemperature_params.hlambda = Parameter(20.00, units='[W][m^-2][K^-1]', scale_object=self.scale_params, return_dimensional=True, description="sensible+turbulent heat exchange between ocean/ground and atmosphere") if self.gotemperature_params is not None: self._number_of_ground_modes = 0 self._number_of_oceanic_modes = self.oblocks.shape[0] self.gotemperature_params.set_insolation(self.nmod[0] * [0.e0]) self.gotemperature_params.set_insolation(350.0, 0) self.gotemperature_params.T0 = Parameter(285.0, units='[K]', scale_object=self.scale_params, return_dimensional=True, description="stationary solution for the 0-th order oceanic temperature") # if setting an ocean, then disable the orography if self.ground_params is not None: self.ground_params.hk = None @property def gblocks(self): """~numpy.ndarray(int): Spectral blocks detailing the model's ground modes :math:`x`-and :math:`y`-wavenumber. Array of shape (:attr:`~QgParams.nmod` [1], 2).""" return self._gms @gblocks.setter def gblocks(self, value): self._oms = None self._gms = value basis = ChannelFourierBasis(self._gms, self.scale_params.n) self._oceanic_basis = None self._ground_basis = basis if self.atemperature_params is not None: # disable the Newtonian cooling self.atemperature_params.thetas = None # np.zeros(self.nmod[0]) self.atemperature_params.hd = None # Parameter(0.0, input_dimensional=False) self.atemperature_params.gamma = Parameter(1.e7, units='[J][m^-2][K^-1]', scale_object=self.scale_params, description='specific heat capacity of the atmosphere', return_dimensional=True) self.atemperature_params.set_insolation(self.nmod[0] * [0.e0]) self.atemperature_params.set_insolation(100.0, 0) self.atemperature_params.eps = Parameter(0.76e0, input_dimensional=False, description="emissivity coefficient for the grey-body atmosphere") self.atemperature_params.T0 = Parameter(270.0, units='[K]', scale_object=self.scale_params, return_dimensional=True, description="stationary solution for the 0-th order atmospheric temperature") self.atemperature_params.sc = Parameter(1., input_dimensional=False, description="ratio of surface to atmosphere temperature") self.atemperature_params.hlambda = Parameter(20.00, units='[W][m^-2][K^-1]', scale_object=self.scale_params, return_dimensional=True, description="sensible+turbulent heat exchange between ocean/ground and atmosphere") if self.gotemperature_params is not None: gmod = 0 for i in range(self.gblocks.shape[0]): if self.ablocks[i, 0] == 1: gmod += 3 else: gmod += 2 self._number_of_ground_modes = gmod self._number_of_oceanic_modes = 0 # if orography is disabled, enable it! if self.ground_params is not None: self.ground_params.orographic_basis = 'atmospheric' if self.ground_params.hk is None: self.ground_params.set_orography(self.nmod[0] * [0.e0]) self.ground_params.set_orography(0.1, 1) self.gotemperature_params.set_insolation(self.nmod[0] * [0.e0]) self.gotemperature_params.set_insolation(350.0, 0) self.gotemperature_params.T0 = Parameter(285.0, units='[K]', scale_object=self.scale_params, return_dimensional=True, description="stationary solution for the 0-th order oceanic temperature") def _set_atmospheric_analytic_fourier_modes(self, nxmax, nymax, auto=False): res = np.zeros((nxmax * nymax, 2), dtype=int) i = 0 for nx in range(1, nxmax + 1): for ny in range(1, nymax+1): res[i, 0] = nx res[i, 1] = ny i += 1 if auto: if self.atemperature_params is None: self.atemperature_params = AtmosphericTemperatureParams(self.scale_params) if self.atmospheric_params is None: self.atmospheric_params = AtmosphericParams(self.scale_params) self.ablocks = res self._atmospheric_latex_var_string = list() self._atmospheric_var_string = list() for i in range(self.nmod[0]): self._atmospheric_latex_var_string.append(r'psi_{\rm a,' + str(i + 1) + "}") self._atmospheric_var_string.append(r'psi_a_' + str(i + 1)) for i in range(self.nmod[0]): self._atmospheric_latex_var_string.append(r'theta_{\rm a,' + str(i + 1) + "}") self._atmospheric_var_string.append(r'theta_a_' + str(i + 1)) def _set_oceanic_analytic_fourier_modes(self, nxmax, nymax, auto=True): if self._ams is None: print('Atmosphere modes not set up. Add an atmosphere before adding an ocean!') print('Oceanic setup aborted.') return res = np.zeros((nxmax * nymax, 2), dtype=int) i = 0 for nx in range(1, nxmax + 1): for ny in range(1, nymax+1): res[i, 0] = nx res[i, 1] = ny i += 1 if auto: if self.gotemperature_params is None or isinstance(self.gotemperature_params, GroundTemperatureParams): self.gotemperature_params = OceanicTemperatureParams(self.scale_params) if self.oceanic_params is None: self.oceanic_params = OceanicParams(self.scale_params) self.ground_params = None self.oblocks = res self._oceanic_latex_var_string = list() self._oceanic_var_string = list() self._ground_latex_var_string = list() self._ground_var_string = list() for i in range(self.nmod[1]): self._oceanic_latex_var_string.append(r'psi_{\rm o,' + str(i + 1) + "}") self._oceanic_var_string.append(r'psi_o_' + str(i + 1)) for i in range(self.nmod[1]): self._oceanic_latex_var_string.append(r'delta T_{\rm o,' + str(i + 1) + "}") self._oceanic_var_string.append(r'delta_T_o_' + str(i + 1)) def _set_ground_analytic_fourier_modes(self, nxmax=None, nymax=None, auto=True): if self._ams is None: print('Atmosphere modes not set up. Add an atmosphere before adding the ground!') print('Ground setup aborted.') return if nxmax is None or nymax is None: res = self._ams.copy() else: res = np.zeros((nxmax * nymax, 2), dtype=int) i = 0 for nx in range(1, nxmax + 1): for ny in range(1, nymax+1): res[i, 0] = nx res[i, 1] = ny i += 1 if auto: if self.gotemperature_params is None or isinstance(self.gotemperature_params, OceanicTemperatureParams): self.gotemperature_params = GroundTemperatureParams(self.scale_params) if self.ground_params is None: self.ground_params = GroundParams(self.scale_params) self.oceanic_params = None self.gblocks = res self._oceanic_var_string = list() self._oceanic_latex_var_string = list() self._ground_latex_var_string = list() self._ground_var_string = list() for i in range(self.nmod[1]): self._ground_latex_var_string.append(r'delta T_{\rm g,' + str(i + 1) + "}") self._ground_var_string.append(r'delta_T_g_' + str(i + 1))