Functions module

The functions module contains utility functions used by the model, as well as some high-level functions to handle it.

Sparse matrix operation module

This module supports operations and functions on the sparse tensors defined in the QgsTensor objects.

qgs.functions.sparse_mul.sparse_mul2(coo, value, vec)[source]

Sparse multiplication of a tensor with one vector: \(A_{i,j} = {\displaystyle \sum_{k=0}^{\mathrm{ndim}}} \, \mathcal{T}_{i,j,k} \, a_k\)

Warning

It is a Numba-jitted function, so it cannot take a sparse.COO sparse tensor directly. The tensor coordinates list and values must be provided separately by the user.

In principle, this will be solved later in sparse, see https://github.com/pydata/sparse/issues/378.

Parameters:
  • coo (ndarray(int)) – A 2D array of shape (n_elems, 3), a list of n_elems tensor coordinates corresponding to each value provided.

  • value (ndarray(float)) – A 1D array of shape (n_elems,), a list of value in the tensor

  • vec (ndarray(float)) – The vector \(a_k\) to contract the tensor with. Must be of shape (ndim + 1,).

Returns:

The matrix \(A_{i,j}\), of shape (ndim + 1, ndim + 1).

Return type:

ndarray(float)

qgs.functions.sparse_mul.sparse_mul3(coo, value, vec_a, vec_b)[source]

Sparse multiplication of a tensor with two vectors: \(v_i = {\displaystyle \sum_{j,k=0}^{\mathrm{ndim}}} \, \mathcal{T}_{i,j,k} \, a_j \, b_k\)

Warning

It is a Numba-jitted function, so it cannot take a sparse.COO sparse tensor directly. The tensor coordinates list and values must be provided separately by the user.

In principle, this will be solved later in sparse, see https://github.com/pydata/sparse/issues/378.

Parameters:
  • coo (ndarray(int)) – A 2D array of shape (n_elems, 3), a list of n_elems tensor coordinates corresponding to each value provided.

  • value (ndarray(float)) – A 1D array of shape (n_elems,), a list of value in the tensor

  • vec_a (ndarray(float)) – The vector \(a_j\) to contract the tensor with. Must be of shape (ndim + 1,).

  • vec_b (ndarray(float)) – The vector \(b_k\) to contract the tensor with. Must be of shape (ndim + 1,).

Returns:

The vector \(v_i\), of shape (ndim + 1,).

Return type:

ndarray(float)

qgs.functions.sparse_mul.sparse_mul4(coo, value, vec_a, vec_b, vec_c)[source]

Sparse multiplication of a rank-5 tensor with three vectors: \(A_{i, j} = {\displaystyle \sum_{k,l,m=0}^{\mathrm{ndim}}} \, \mathcal{T}_{i,j,k,l, m} \, a_k \, b_l \, c_m\)

Warning

It is a Numba-jitted function, so it cannot take a sparse.COO sparse tensor directly. The tensor coordinates list and values must be provided separately by the user.

In principle, this will be solved later in sparse, see https://github.com/pydata/sparse/issues/378.

Parameters:
  • coo (ndarray(int)) – A 2D array of shape (n_elems, 5), a list of n_elems tensor coordinates corresponding to each value provided.

  • value (ndarray(float)) – A 1D array of shape (n_elems,), a list of value in the tensor

  • vec_a (ndarray(float)) – The vector \(a_j\) to contract the tensor with. Must be of shape (ndim + 1,).

  • vec_b (ndarray(float)) – The vector \(b_k\) to contract the tensor with. Must be of shape (ndim + 1,).

  • vec_c (ndarray(float)) – The vector \(c_l\) to contract the tensor with. Must be of shape (ndim + 1,).

Returns:

The matrix \(A_{i, j}\), of shape (ndim + 1, ndim + 1).

Return type:

ndarray(float)

qgs.functions.sparse_mul.sparse_mul5(coo, value, vec_a, vec_b, vec_c, vec_d)[source]

Sparse multiplication of a rank-5 tensor with four vectors: \(v_i = {\displaystyle \sum_{j,k,l,m=0}^{\mathrm{ndim}}} \, \mathcal{T}_{i,j,k,l,m} \, a_j \, b_k \, c_l \, d_m\)

Warning

It is a Numba-jitted function, so it cannot take a sparse.COO sparse tensor directly. The tensor coordinates list and values must be provided separately by the user.

In principle, this will be solved later in sparse, see https://github.com/pydata/sparse/issues/378.

Parameters:
  • coo (ndarray(int)) – A 2D array of shape (n_elems, 5), a list of n_elems tensor coordinates corresponding to each value provided.

  • value (ndarray(float)) – A 1D array of shape (n_elems,), a list of value in the tensor

  • vec_a (ndarray(float)) – The vector \(a_j\) to contract the tensor with. Must be of shape (ndim + 1,).

  • vec_b (ndarray(float)) – The vector \(b_k\) to contract the tensor with. Must be of shape (ndim + 1,).

  • vec_c (ndarray(float)) – The vector \(c_l\) to contract the tensor with. Must be of shape (ndim + 1,).

  • vec_d (ndarray(float)) – The vector \(d_m\) to contract the tensor with. Must be of shape (ndim + 1,).

Returns:

The vector \(v_i\), of shape (ndim + 1,).

Return type:

ndarray(float)

Tendencies definition module

This module provides functions to create the tendencies functions of the model, based on its parameters.

qgs.functions.tendencies.create_atmo_thermo_tendencies(params, return_atmo_thermo_tensor=False)[source]

Returns a function which return a part of the atmospheric thermodynamic tendencies useful for computing the vertical wind.

Parameters:
  • params (QgParams) – The parameters fully specifying the model configuration.

  • return_atmo_thermo_tensor (bool) – If True, return the tendencies tensor of these tencencies. Default to False.

Returns:

  • f (callable) – The numba-jitted tendencies function.

  • atmo_thermo_tensor (AtmoThermoTensor) – If return_atmo_thermo_tensor is True, the tendencies tensor of the system.

qgs.functions.tendencies.create_tendencies(params, return_inner_products=False, return_qgtensor=False)[source]

Function to handle the inner products and tendencies tensors construction. Returns the tendencies function \(\boldsymbol{f}\) determining the model’s ordinary differential equations:

\[\dot{\boldsymbol{x}} = \boldsymbol{f}(\boldsymbol{x})\]

which is for the model’s integration.

It returns also the linearized tendencies \(\boldsymbol{\mathrm{J}} \equiv \boldsymbol{\mathrm{D}f} = \frac{\partial \boldsymbol{f}}{\partial \boldsymbol{x}}\) (Jacobian matrix) which are used by the tangent linear model:

\[\dot{\boldsymbol{\delta x}} = \boldsymbol{\mathrm{J}}(\boldsymbol{x}) \cdot \boldsymbol{\delta x}\]
Parameters:
  • params (QgParams) – The parameters fully specifying the model configuration.

  • return_inner_products (bool) – If True, return the inner products of the model. Default to False.

  • return_qgtensor (bool) – If True, return the tendencies tensor of the model. Default to False.

Returns:

  • f (callable) – The numba-jitted tendencies function.

  • Df (callable) – The numba-jitted linearized tendencies function.

  • inner_products ((AtmosphericInnerProducts, OceanicInnerProducts)) – If return_inner_products is True, the inner products of the system.

  • qgtensor (QgsTensor) – If return_qgtensor is True, the tendencies tensor of the system.

Symbolic matrix operation module

This module supports operations and functions on the symbolic sparse tensors defined in the SymbolicQgsTensor objects.

qgs.functions.symbolic_mul.symbolic_sparse_mult2(dic, vec_a)[source]

Symbolic multiplication of a tensor with one vector: \(A_{i,j} = {\displaystyle \sum_{k=0}^{\mathrm{ndim}}} \, \mathcal{T}_{i,j,k} \, a_k\)

Parameters:
  • dic (dict(Symbol)) – A dictionary whose keys are the coordinates of the tensor, and the dictionary values are the values of the tensor.

  • vec_a (list(Symbol)) – The list \(a_k\) to contract the tensor with entries of Sympy Symbols. Must be of shape (ndim + 1,).

Returns:

res – The matrix \(A_{i,j}\), of shape (ndim + 1, ndim + 1), contained in a dictionary, where the keys are the tensor coordinates, and the values are the tensor values.

Return type:

dict(Symbol)

qgs.functions.symbolic_mul.symbolic_sparse_mult3(dic, vec_a, vec_b)[source]

Symbolic multiplication of a tensor with two vectors: \(v_i = {\displaystyle \sum_{j,k=0}^{\mathrm{ndim}}} \, \mathcal{T}_{i,j,k} \, a_j \, b_k\)

Parameters:
  • dic (dict(Symbol)) – A dictionary whose keys are the coordinates of the tensor, and the dictionary values are the values of the tensor.

  • vec_a (list(Symbol)) – The list \(a_j\) to contract the tensor with entries of Sympy Symbols. Must be of shape (ndim + 1,).

  • vec_b (list(Symbol)) – The list \(b_k\) to contract the tensor with entries of Sympy Symbols. Must be of shape (ndim + 1,).

Returns:

res – The vector \(v_i\), of shape (ndim + 1,), contained in a dictionary, where the keys are the tensor coordinates, and the values are the tensor values.

Return type:

dict(Symbol)

qgs.functions.symbolic_mul.symbolic_sparse_mult4(dic, vec_a, vec_b, vec_c)[source]

Symbolic multiplication of a rank-5 tensor with three vectors: \(A_{i, j} = {\displaystyle \sum_{k,l,m=0}^{\mathrm{ndim}}} \, \mathcal{T}_{i,j,k,l, m} \, a_k \, b_l \, c_m\)

Parameters:
  • dic (dict(Symbol)) – A dictionary where they keys are a tuple of 5 elements which are the coordinates of the tensor values, which are contained in the dictionary values.

  • vec_a (list(Symbol)) – The list \(a_j\) to contract the tensor with entries of Sympy Symbols. Must be of shape (ndim + 1,).

  • vec_b (list(Symbol)) – The list \(b_k\) to contract the tensor with entries of Sympy Symbols. Must be of shape (ndim + 1,).

  • vec_c (list(Symbol)) – The list \(c_l\) to contract the tensor with entries of Sympy Symbols. Must be of shape (ndim + 1,).

Returns:

res – The matrix \(A_{i, j}\), of shape (ndim + 1, ndim + 1), contained in a dictionary, where the keys are the tensor coordinates, and the values are the tensor values.

Return type:

dict(Symbol)

qgs.functions.symbolic_mul.symbolic_sparse_mult5(dic, vec_a, vec_b, vec_c, vec_d)[source]

Symbolic multiplication of a rank-5 tensor with four vectors: \(v_i = {\displaystyle \sum_{j,k,l,m=0}^{\mathrm{ndim}}} \, \mathcal{T}_{i,j,k,l,m} \, a_j \, b_k \, c_l \, d_m\)

Parameters:
  • dic (dict(Symbol)) – A dictionary whose keys are the coordinates of the tensor, and the dictionary values are the values of the tensor.

  • vec_a (list(Symbol)) – The list \(a_j\) to contract the tensor with entries of Sympy Symbols. Must be of shape (ndim + 1,).

  • vec_b (list(Symbol)) – The list \(b_k\) to contract the tensor with entries of Sympy Symbols. Must be of shape (ndim + 1,).

  • vec_c (list(Symbol)) – The list \(c_l\) to contract the tensor with entries of Sympy Symbols. Must be of shape (ndim + 1,).

  • vec_d (list(Symbol)) – The list \(d_m\) to contract the tensor with entries of Sympy Symbols. Must be of shape (ndim + 1,).

Returns:

res – The vector \(v_i\), of shape (ndim + 1,), contained in a dictionary, where the keys are the tensor coordinates, and the values are the tensor values.

Return type:

dict(Symbol)

qgs.functions.symbolic_mul.symbolic_tensordot(a, b, axes=2)[source]

Compute tensor dot product along specified axes of two sympy symbolic arrays

This is based on Numpy tensordot() .

Parameters:
  • a (DenseNDimArray or SparseNDimArray) – Arrays to take the dot product of.

  • b (DenseNDimArray or SparseNDimArray) – Arrays to take the dot product of.

  • axes (int) – Sum over the last axes axes of a and the first axes axes of b in order. The sizes of the corresponding axes must match.

Returns:

output – The tensor dot product of the input.

Return type:

Sympy tensor

Symbolic tendencies module

This module provides functions to create a symbolic representation of the tendencies functions of the model in various languages and for various external software.

qgs.functions.symbolic_tendencies.create_auto_file(equations, params, continuation_variables, auto_main_template=None, auto_c_template=None, initialize_params=False, initialize_solution=False)[source]

Creates the AUTO configuration file and the model file. Saves files to specified folder.

Parameters:
  • equations (dict) – Dictionary of the substituted symbolic model equations

  • params (QgParams) – The parameters fully specifying the model configuration.

  • continuation_variables (list(Parameter, ScalingParameter, ParametersArray)) – The variables to not substitute by their numerical value and to leave in the equations. There must be at least one variable in this list and, due to AUTO constraints, its maximum length is 10.

  • auto_main_template (str, optional) – The template to be used to generate the main AUTO file. If not provided, use the default template.

  • auto_c_template (str, optional) – The template to be used to generate the AUTO config file. If not provided, use the default template.

  • initialize_params (bool, optional) – Add lines in the AUTO STPNT function to initialize the parameters. Default to False.

  • initialize_solution (bool, optional) – Add lines in the AUTO STPNT function to initialize the solution. Default to False.

Returns:

  • auto_file (str) – The auto model file as a string

  • auto_config (str) – Auto configuration file as a string

qgs.functions.symbolic_tendencies.create_symbolic_tendencies(params, continuation_variables, atm_ip=None, ocn_ip=None, gnd_ip=None, language='python', return_inner_products=False, return_jacobian=False, return_symbolic_eqs=False, return_symbolic_qgtensor=False)[source]

Function to output the raw symbolic tendencies of the qgs model.

Parameters:
  • params (QgParams) – The parameters fully specifying the model configuration.

  • continuation_variables (list(Parameter, ScalingParameter or ParametersArray) or None) – The variables to not substitute by their numerical value and to leave in the equations. If None, no variables are substituted. If an empty list is provided, then all variables are substituted, providing fully numerical tendencies.

  • atm_ip (AtmosphericSymbolicInnerProducts, optional) – Allows for stored inner products to be input.

  • ocn_ip (OceanSymbolicInnerProducts, optional) – Allows for stored inner products to be input.

  • gnd_ip (GroundSymbolicInnerProducts, optional) – Allows for stored inner products to be input.

  • language (str, optional) – Options for the output language syntax: ‘python’, ‘julia’, ‘fortran’, ‘auto’, ‘mathematica’. Default to ‘python’.

  • return_inner_products (bool, optional) – If True, return the inner products of the model. Default to False.

  • return_jacobian (bool, optional) – If True, return the Jacobian of the model. Default to False.

  • return_symbolic_eqs (bool, optional) – If True, return the substituted symbolic equations.

  • return_symbolic_qgtensor (bool, optional) – If True, return the symbolic tendencies tensor of the model. Default to False.

Returns:

  • funcs (str) – The substituted functions in the language syntax specified, as a string.

  • dict_eq_simplified (dict(~sympy.core.expr.Expr)) – Dictionary of the substituted Jacobian matrix.

  • inner_products (tuple(AtmosphericSymbolicInnerProducts, OceanicSymbolicInnerProducts, GroundSymbolicInnerProducts)) – If return_inner_products is True, the inner products of the system.

  • eq_simplified (dict(~sympy.core.expr.Expr)) – If return_symbolic_eqs is True, dictionary of the model tendencies symbolic functions.

  • agotensor (SymbolicQgsTensor) – If return_symbolic_qgtensor is True, the symbolic tendencies tensor of the system.

qgs.functions.symbolic_tendencies.equation_as_function(equations, params, continuation_variables, language='python')[source]

Converts the symbolic equations to a function in string format in the language syntax specified, or a lambdified python function.

Parameters:
  • equations (dict(Expr)) – Dictionary of the substituted symbolic model equations.

  • params (QgParams) – The parameters fully specifying the model configuration.

  • continuation_variables (list(Parameter, ScalingParameter, ParametersArray) or None) – The variables to not substitute by their numerical value and to leave in the equations. If None, no variables are substituted. If an empty list is provided, then all variables are substituted, providing fully numerical tendencies.

  • language (str, optional) –

    Language syntax that the equations are returned in. Options are:

    • python

    • fortran

    • julia

    • auto

    • mathematica

    Default to python.

Returns:

f_output – Output is a function as a string in the specified language syntax.

Return type:

str

qgs.functions.symbolic_tendencies.equations_to_string(equations)[source]

Converts the symbolic equations, held in a dict format, to a dict of strings.

Parameters:

equations (dict(Expr)) – Dictionary of the substituted symbolic model equations.

Returns:

Dictionary of the substituted symbolic model equations.

Return type:

dict(str)

qgs.functions.symbolic_tendencies.format_equations(equations, params, save_loc=None, language='python', print_equations=False)[source]

Function formats the equations, in the programming language specified, and saves the equations to the specified location. The variables in the equation are substituted if the model variable is input.

Parameters:
  • equations (dict(Expr)) – Dictionary of symbolic model equations.

  • params (QgParams) – The parameters fully specifying the model configuration.

  • save_loc (str, optional) – Location to save the outputs as a .txt file.

  • language (str, optional) –

    Language syntax that the equations are returned in. Options are:

    • python

    • fortran

    • julia

    • auto

    • mathematica

    Default to python.

  • print_equations (bool, optional) – If True, equations are printed by the function, if False, equation strings are returned by the function. Defaults to False

Returns:

equation_dict – Dictionary of symbolic model equations, that have been substituted with numerical values.

Return type:

dict(Expr)

qgs.functions.symbolic_tendencies.jacobian_as_function(equations, params, continuation_variables, language='python')[source]

Converts the symbolic equations of the jacobain to a function in string format in the language syntax specified, or a lambdified python function.

Parameters:
  • equations (dict(Expr)) – Dictionary of the substituted symbolic model equations.

  • params (QgParams) – The parameters fully specifying the model configuration.

  • continuation_variables (list(Parameter, ScalingParameter, ParametersArray) or None) – The variables to not substitute by their numerical value and to leave in the equations. If None, no variables are substituted. If an empty list is provided, then all variables are substituted, providing fully numerical tendencies.

  • language (str, optional) –

    Language syntax that the equations are returned in. Options are:

    • python

    • fortran

    • julia

    • auto

    • mathematica

    Default to python.

Returns:

f_output – Output is a function as a string in the specified language syntax

Return type:

str

qgs.functions.symbolic_tendencies.translate_equations(equations, language='python')[source]

Function to output the model equations as a string in the specified language syntax.

Parameters:
  • equations (dict(str), list, str) – Dictionary, list, or string of the symbolic model equations.

  • language (str, optional) –

    Language syntax that the equations are returned in. Options are:

    • python

    • fortran

    • julia

    • auto

    • mathematica

    Default to python.

Returns:

str_eq – Dictionary of strings of the model equations.

Return type:

dict(str)

Utility functions module

This module has some useful functions for the model.

qgs.functions.util.add_to_dict(dic, loc, value)[source]

Adds value to dictionary dic, with the dictionary key of loc. If the dictionary did not have a key of loc before, a new key is made.

Parameters:
  • dic (dict) – Dictionary to add the value to.

  • loc – Item of the dictionary to add the value to.

  • value – Value to add.

qgs.functions.util.normalize_matrix_columns(a)[source]

Numba-jitted function to normalize the columns of a 2D array to one.

Parameters:

a (ndarray) – The 2D array to column-normalize.

Returns:

The normalized array.

Return type:

ndarray

qgs.functions.util.reverse(a)[source]

Numba-jitted function to reverse a 1D array.

Parameters:

a (ndarray) – The 1D array to reverse.

Returns:

The reversed array.

Return type:

ndarray

qgs.functions.util.solve_triangular_matrix(a, b)[source]

Solve a triangular linear matrix equation \(ax = b\).

Parameters:
  • a (ndarray) – The 2D array \(a\).

  • b (ndarray) – The 2D array \(b\).

Returns:

The 2D array of solution \(x\).

Return type:

ndarray