FlowSieve  3.4.0
FlowSieve Coarse-Graining Documentation
functions.hpp File Reference

Collection of all computation-related functions. More...

#include <stdio.h>
#include <stdlib.h>
#include <vector>
#include <string>
#include <map>
#include <mpi.h>
#include "constants.hpp"
Include dependency graph for functions.hpp:
This graph shows which files directly or indirectly include this file:

Go to the source code of this file.

Classes

class  dataset
 Class to store main variables. More...
 
class  Timing_Records
 Class for storing internal timings. More...
 
class  InputParser
 Class to process command-line arguments. More...
 

Functions

void compute_areas (std::vector< double > &areas, const std::vector< double > &longitude, const std::vector< double > &latitude)
 Compute the area of each computational cell. More...
 
size_t Index (const int Itime, const int Idepth, const int Ilat, const int Ilon, const int Ntime, const int Ndepth, const int Nlat, const int Nlon)
 Convenience tool to convert physical index (time, depth, lat, lon) to a logical index. More...
 
void Index1to4 (const size_t index, int &Itime, int &Idepth, int &Ilat, int &Ilon, const int &Ntime, const int &Ndepth, const int &Nlat, const int &Nlon)
 Convenience tool to convert logical index to physical index (time, depth, lat, lon). More...
 
double distance (const double lon1, const double lat1, const double lon2, const double lat2, const double Llon=0, const double Llat=0)
 Compute the distance (in metres) between two points in the domain. More...
 
void compute_local_kernel (std::vector< double > &local_kernel, const double scale, const dataset &source_data, const int Ilat, const int Ilon, const int LAT_lb, const int LAT_ub)
 Compute an array of the kernel values from a given reference point to every other point in the domain. More...
 
void KE_from_vels (std::vector< double > &KE, std::vector< double > *u1, std::vector< double > *u2, std::vector< double > *u3, const std::vector< bool > &mask, const double rho0=constants::rho0)
 Compute (local) KE density from the provided velocities. More...
 
void vel_Spher_to_Cart (std::vector< double > &u_x, std::vector< double > &u_y, std::vector< double > &u_z, const std::vector< double > &u_r, const std::vector< double > &u_lon, const std::vector< double > &u_lat, const dataset &source_data)
 Wrapper that applies vel_Spher_to_Cart_at_point to every point in the domain. More...
 
void vel_Spher_to_Cart_at_point (double &u_x, double &u_y, double &u_z, const double u_r, const double u_lon, const double u_lat, const double lon, const double lat)
 Convert single spherical velocity to Cartesian velocity. More...
 
void vel_Cart_to_Spher (std::vector< double > &u_r, std::vector< double > &u_lon, std::vector< double > &u_lat, const std::vector< double > &u_x, const std::vector< double > &u_y, const std::vector< double > &u_z, const dataset &source_data)
 Wrapper that applies vel_Spher_to_Cart_at_point to every point in the domain. More...
 
void vel_Cart_to_Spher_at_point (double &u_r, double &u_lon, double &u_lat, const double u_x, const double u_y, const double u_z, const double lon, const double lat)
 Convert single Cartesian velocity to spherical velocity. More...
 
void filtering (const dataset &source_data, const std::vector< double > &scales, const MPI_Comm comm=MPI_COMM_WORLD)
 Main filtering driver. More...
 
void filtering_helmholtz (const dataset &source_data, const std::vector< double > &scales, const MPI_Comm comm=MPI_COMM_WORLD)
 Main filtering driver for Helmholtz decomposed data. More...
 
void apply_filter_at_point (std::vector< double *> &coarse_val, const std::vector< const std::vector< double > *> &fields, const dataset &source_data, const int Itime, const int Idepth, const int Ilat, const int Ilon, const int LAT_lb, const int LAT_ub, const double scale, const std::vector< bool > &use_mask, const std::vector< double > &local_kernel, const std::vector< double > *weight=NULL)
 Compute filtered field at a single point. More...
 
double kernel (const double distance, const double scale)
 Primary kernel function coarse-graining procedure (G in publications) More...
 
double kernel_alpha (void)
 Compute alpha value for kernel (for baroclinic transfer) More...
 
void compute_vorticity_at_point (double &vort_r_tmp, double &vort_lon_tmp, double &vort_lat_tmp, double &div_tmp, double &OkuboWeiss_tmp, const dataset &source_data, const std::vector< double > &u_r, const std::vector< double > &u_lon, const std::vector< double > &u_lat, const int Itime, const int Idepth, const int Ilat, const int Ilon)
 Compute the vorticity, divergnce, and OkuboWeiss at a point. More...
 
void compute_vorticity (std::vector< double > &vort_r, std::vector< double > &vort_lon, std::vector< double > &vort_lat, std::vector< double > &vel_div, std::vector< double > &OkuboWeiss, const dataset &source_data, const std::vector< double > &u_r, const std::vector< double > &u_lon, const std::vector< double > &u_lat, const MPI_Comm comm=MPI_COMM_WORLD)
 Wrapper for computing vorticity. More...
 
void apply_filter_at_point_for_quadratics (double &uxux_tmp, double &uxuy_tmp, double &uxuz_tmp, double &uyuy_tmp, double &uyuz_tmp, double &uzuz_tmp, double &vort_ux_tmp, double &vort_uy_tmp, double &vort_uz_tmp, const std::vector< double > &u_x, const std::vector< double > &u_y, const std::vector< double > &u_z, const std::vector< double > &vort_r, const dataset &source_data, const int Itime, const int Idepth, const int Ilat, const int Ilon, const int LAT_lb, const int LAT_ub, const double scale, const std::vector< double > &local_kernel)
 Compute filtered quadratic velocities at a single point. More...
 
void compute_Pi (std::vector< double > &energy_transfer, const dataset &source_data, const std::vector< double > &ux, const std::vector< double > &uy, const std::vector< double > &uz, const std::vector< double > &uxux, const std::vector< double > &uxuy, const std::vector< double > &uxuz, const std::vector< double > &uyuy, const std::vector< double > &uyuz, const std::vector< double > &uzuz, const MPI_Comm comm=MPI_COMM_WORLD)
 Compute the energy transfer through the current filter scale. More...
 
void compute_Pi_shift_deriv (std::vector< double > &energy_transfer, const dataset &source_data, const std::vector< double > &ux, const std::vector< double > &uy, const std::vector< double > &uz, const std::vector< double > &uxux, const std::vector< double > &uxuy, const std::vector< double > &uxuz, const std::vector< double > &uyuy, const std::vector< double > &uyuz, const std::vector< double > &uzuz, const MPI_Comm comm=MPI_COMM_WORLD)
 Compute the energy transfer through the current filter scale. More...
 
void compute_Pi_Helmholtz (std::vector< double > &energy_transfer, const dataset &source_data, const std::vector< double > &ulon, const std::vector< double > &ulat, const std::vector< double > &ulon_ulon, const std::vector< double > &ulon_ulat, const std::vector< double > &ulat_ulat, const MPI_Comm comm=MPI_COMM_WORLD)
 Compute the energy transfer through the current filter scale. More...
 
void compute_Z (std::vector< double > &enstrophy_transfer, const dataset &source_data, const std::vector< double > &ux, const std::vector< double > &uy, const std::vector< double > &uz, const std::vector< double > &coarse_vort_r, const std::vector< double > &vort_ux, const std::vector< double > &vort_uy, const std::vector< double > &vort_uz, const MPI_Comm comm=MPI_COMM_WORLD)
 Compute the enstrophy transfer through the current filter scale. More...
 
void compute_Lambda_rotational (std::vector< double > &Lambda_rot, const std::vector< double > &coarse_vort_r, const std::vector< double > &coarse_vort_lon, const std::vector< double > &coarse_vort_lat, const std::vector< double > &coarse_rho, const std::vector< double > &coarse_p, const int Ntime, const int Ndepth, const int Nlat, const int Nlon, const std::vector< double > &longitude, const std::vector< double > &latitude, const std::vector< bool > &mask, const double scale_factor)
 
void compute_Lambda_full (std::vector< double > &Lambda, const std::vector< double > &coarse_u_r, const std::vector< double > &coarse_u_lon, const std::vector< double > &coarse_u_lat, const std::vector< double > &tilde_u_r, const std::vector< double > &tilde_u_lon, const std::vector< double > &tilde_u_lat, const std::vector< double > &coarse_p, const int Ntime, const int Ndepth, const int Nlat, const int Nlon, const std::vector< double > &longitude, const std::vector< double > &latitude, const std::vector< bool > &mask)
 
void compute_Lambda_nonlin_model (std::vector< double > &Lambda_nonlin, const std::vector< double > &coarse_u_r, const std::vector< double > &coarse_u_lon, const std::vector< double > &coarse_u_lat, const std::vector< double > &coarse_rho, const std::vector< double > &coarse_p, const int Ntime, const int Ndepth, const int Nlat, const int Nlon, const std::vector< double > &longitude, const std::vector< double > &latitude, const std::vector< bool > &mask, const double scale_factor)
 
double depotential_temperature (const double p, const double theta)
 Convert potential temperature to actual temperature following equation 1.154a from Vallis (p.36) More...
 
double equation_of_state (const double T, const double S, const double p)
 UNESCO 1981 equation of state. More...
 
void compute_div_transport (std::vector< double > &div_J, const dataset &source_data, const std::vector< double > &u_x, const std::vector< double > &u_y, const std::vector< double > &u_z, const std::vector< double > &uxux, const std::vector< double > &uxuy, const std::vector< double > &uxuz, const std::vector< double > &uyuy, const std::vector< double > &uyuz, const std::vector< double > &uzuz, const std::vector< double > &coarse_p, const MPI_Comm comm=MPI_COMM_WORLD)
 Compute KE transport caused by div(J) More...
 
void compute_spatial_average (std::vector< double > &means, const std::vector< double > &field, const std::vector< double > &areas, const int Ntime, const int Ndepth, const int Nlat, const int Nlon, const std::vector< bool > &mask)
 Given a field, compute the horizontal average. More...
 
void get_lat_bounds (int &LAT_lb, int &LAT_ub, const std::vector< double > &latitude, const int Ilat, const double scale)
 Get latitude bounds for coarse-graining. More...
 
void get_lon_bounds (int &LON_lb, int &LON_ub, const std::vector< double > &longitude, const int Ilon, const double centre_lat, const double curr_lat, const double scale)
 Get longitude integratation bounds for a specific latitude. More...
 
void print_compile_info (const std::vector< double > *scales=NULL)
 Print summary of compile-time variables. More...
 
int get_omp_chunksize (const int Nlat, const int Nlon)
 Compute openmp chunk size for OMP for loops. More...
 
void convert_coordinates (std::vector< double > &longitude, std::vector< double > &latitude)
 Apply degree->radian conversion, if necessary. Applied in-place. More...
 
void mask_out_pole (const std::vector< double > &latitude, std::vector< bool > &mask, const int Ntime, const int Ndepth, const int Nlat, const int Nlon)
 If the grid includes the poles (lat = 90 degrees), then mask it out. More...
 
void roll_field (std::vector< double > &field_to_roll, const std::string dimension, const int roll_count, const int Ntime, const int Ndepth, const int Nlat, const int Nlon)
 Roll field along dimension. More...
 
void extend_latitude_to_poles (const std::vector< double > &original_latitude, std::vector< double > &extended_latitude, int &orig_Ilat_start, const bool IS_DEGREES=false, const MPI_Comm comm=MPI_COMM_WORLD)
 
void extend_field_to_poles (std::vector< double > &array_to_extend, const dataset &source_data, const std::vector< double > &extended_latitude, const int Ilat_start)
 
void extend_mask_to_poles (std::vector< bool > &mask_to_extend, const dataset &source_data, const std::vector< double > &extended_latitude, const int Ilat_start, const bool extend_val=constants::FILTER_OVER_LAND)
 
bool string_to_bool (std::string str)
 Convert string to boolean.
 
void print_header_info ()
 

Detailed Description

Collection of all computation-related functions.

Function Documentation

◆ apply_filter_at_point()

void apply_filter_at_point ( std::vector< double *> &  coarse_vals,
const std::vector< const std::vector< double > *> &  fields,
const dataset source_data,
const int  Itime,
const int  Idepth,
const int  Ilat,
const int  Ilon,
const int  LAT_lb,
const int  LAT_ub,
const double  scale,
const std::vector< bool > &  use_mask,
const std::vector< double > &  local_kernel,
const std::vector< double > *  weight 
)

Compute filtered field at a single point.

Computes the integral of the provided field with the kernel().

Parameters
[in,out]coarse_valwhere to store filtered value
[in]fieldsfields to filter
[in]source_datadataset class instance containing data (Psi, Phi, etc)
[in]Itime,Idepth,Ilat,Iloncurrent position in time dimension
[in]LAT_lb,LAT_ublower/upper boundd on latitude for kernel
[in]scalefiltering scale
[in]use_maskarray of booleans indicating whether or not to use mask (i.e. zero out land) or to use the array value
[in]local_kernelpre-computed kernel (NULL indicates not provided)
[in]weightpointer to spatial weight (i.e. rho) (NULL indicates not provided)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ apply_filter_at_point_for_quadratics()

void apply_filter_at_point_for_quadratics ( double &  uxux_tmp,
double &  uxuy_tmp,
double &  uxuz_tmp,
double &  uyuy_tmp,
double &  uyuz_tmp,
double &  uzuz_tmp,
double &  vort_ux_tmp,
double &  vort_uy_tmp,
double &  vort_uz_tmp,
const std::vector< double > &  u_x,
const std::vector< double > &  u_y,
const std::vector< double > &  u_z,
const std::vector< double > &  vort_r,
const dataset source_data,
const int  Itime,
const int  Idepth,
const int  Ilat,
const int  Ilon,
const int  LAT_lb,
const int  LAT_ub,
const double  scale,
const std::vector< double > &  local_kernel 
)

Compute filtered quadratic velocities at a single point.

Computes the integral of each quadratic Cartesian velocity with the kernel().

In particular, the quadratic terms being filtered are: $u_xu_x$, $u_xu_y$, $u_xu_z$, $u_yu_y$, $u_yu_z$, $u_zu_z$

Parameters
[in,out]uxux_tmpwhere to store filtered (u_x)*(u_x)
[in,out]uxuy_tmpwhere to store filtered (u_x)*(u_y)
[in,out]uxuz_tmpwhere to store filtered (u_x)*(u_z)
[in,out]uyuy_tmpwhere to store filtered (u_y)*(u_y)
[in,out]uyuz_tmpwhere to store filtered (u_y)*(u_z)
[in,out]uzuz_tmpwhere to store filtered (u_z)*(u_z)
[in,out]vort_ux_tmpwhere to store filtered (vort_r)*(u_x)
[in,out]vort_uy_tmpwhere to store filtered (vort_r)*(u_y)
[in,out]vort_uz_tmpwhere to store filtered (vort_r)*(u_z)
[in]u_x,u_y,u_zfields to filter
[in]vort_rvorticity field to filter
[in]source_datadataset class instance containing data (Psi, Phi, etc)
[in]Itime,Idepth,Ilat,Iloncurrent position in time dimension
[in]LAT_lb,LAT_ublower/upper boundd on latitude for kernel
[in]scalefiltering scale
[in]local_kernelpre-computed kernel (NULL indicates not provided)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ compute_areas()

void compute_areas ( std::vector< double > &  areas,
const std::vector< double > &  longitude,
const std::vector< double > &  latitude 
)

Compute the area of each computational cell.

Works for both spherical and Cartesian coordinated (based on constants.hpp) as well as uniform and non-uniform grids (also based on constants.hpp).

Parameters
[in,out]areasarray in which areas will be stored
[in]longitudelongitude vector (1D)
[in]latitudelatitude vector (1D)

◆ compute_div_transport()

void compute_div_transport ( std::vector< double > &  div_J,
const dataset source_data,
const std::vector< double > &  u_x,
const std::vector< double > &  u_y,
const std::vector< double > &  u_z,
const std::vector< double > &  uxux,
const std::vector< double > &  uxuy,
const std::vector< double > &  uxuz,
const std::vector< double > &  uyuy,
const std::vector< double > &  uyuz,
const std::vector< double > &  uzuz,
const std::vector< double > &  coarse_p,
const MPI_Comm  comm 
)

Compute KE transport caused by div(J)

Currently implements:

  • advection by coarse-scale velocity
  • pressure-induced transport
  • advection by fine-scale velocity

NOT implemented:

  • Diffusion
*  J_transport =   0.5 * rho0 * | u_l |^2 * u_l
*                + P_l * u_l
*                - nu * 0.5 * rho * grad( | u_l |^2 )
*                + rho0 * u_l * tau(u_l, u_l)
*
*
*  (Spherical)
*     div(J) = ( 
*                (1 / (r * cos(lat)) ) d/dlon (J_lon),
*                (1 / (r * cos(lat)) ) d/dlat (J_lat * cos(lat)),
*                (1 /  r^2           ) d/dr   (J_r * r^2)
*              )
*
*  (Cartesian)
*     div(J) = ( 
*                d/dx (J_x),
*                d/dy (J_y),
*                d/dz (J_z) 
*              )
*
*
*  Term 1: 0.5 * rho0 * | u_l |^2 * u_l
*      This is advection of large-scale KE by the large-scale 
*      velocity. 
*
*      (index form: 0.5 * rho0 * [ (u_i*u_i) * u_j ]    )
*      (   of grad: 0.5 * rho0 * [ (u_i*u_i) * u_j ],j  )
*               = 0.5 * rho0 * (u_i*u_i),j * u_j 
*               =       rho0 * u_i * u_i,j * u_j
*
*
*
*  Term 2: P_l * u_l
*      Transport caused by pressure 
*
*      (index form:  p * u_j     )
*      (   of grad: (p * u_j),j  )
*
*
*
*  Term 3: - nu * 0.5 * rho * grad( | u_l |^2 )
*      This is diffusion
*      NOT YET IMPLEMENTED
*
*
*
*  Term 4: rho0 * u_l * tau(u, u)
*      This is advection of the large-scale KE
*      by the small-scale flow.
*
*      (index form:  rho0 *   u_i * tau_ij     )
*      (   of grad:  rho0 * [ u_i * tau_ij ],j )
*        = rho0 * ( u_i,j * tau_ij  + u_i * tau_ij,j )
*        = rho0 * 
*            (  
*               u_i,j * ( bar(u_i*u_j)   - bar(u_i  )*bar(u_j) )
*             + u_i   * ( bar(u_i*u_j),j - bar(u_i,j)*bar(u_j) )
*            )
*  
Parameters
[in,out]div_Jwhere to store the computed values (array)
[in]u_x,u_y,u_zcoarse Cartesian velocity components
[in]uxux,uxuy,uxuz,uyuy,uyuz,uzuzcoarse velocity products (e.g. bar(u*v) )
[in]coarse_pcoarse pressure
Here is the call graph for this function:
Here is the caller graph for this function:

◆ compute_local_kernel()

void compute_local_kernel ( std::vector< double > &  local_kernel,
const double  scale,
const dataset source_data,
const int  Ilat,
const int  Ilon,
const int  LAT_lb,
const int  LAT_ub 
)

Compute an array of the kernel values from a given reference point to every other point in the domain.

(ref_ilat, ref_ilon) is the reference point from which the kernel values are computed.

LAT_lb and LAT_ub are the (pre-computed) latitudinal bounds for the kernel.

Parameters
[in,out]local_kernelwhere to store the local kernel
[in]scaleFiltering scale
[in]source_datadataset class instance containing data (Psi, Phi, etc)
[in]Ilat,Ilonreference coordinate (kernel centre)
[in]LAT_lb,LAT_ubupper and lower latitudinal bounds for kernel
Here is the call graph for this function:
Here is the caller graph for this function:

◆ compute_Pi()

void compute_Pi ( std::vector< double > &  energy_transfer,
const dataset source_data,
const std::vector< double > &  ux,
const std::vector< double > &  uy,
const std::vector< double > &  uz,
const std::vector< double > &  uxux,
const std::vector< double > &  uxuy,
const std::vector< double > &  uxuz,
const std::vector< double > &  uyuy,
const std::vector< double > &  uyuz,
const std::vector< double > &  uzuz,
const MPI_Comm  comm 
)

Compute the energy transfer through the current filter scale.

In particular, computes $ \rho_0 * ( u_i \tau_{ij,j} - (u_i \tau_{ij})_{,j} ) $

This computation is applied to the Cartesian velocity components

Parameters
[in,out]energy_transferwhere to store the computed values (array)
[in]source_datadataset class instance containing data (Psi, Phi, etc)
[in]ux,uy,uzcoarse Cartesian velocity components
[in]uxux,uxuy,uxuz,uyuy,uyuz,uzuzcoarse velocity products (e.g. bar(u*v) )
[in]commMPI communicator object
Here is the call graph for this function:
Here is the caller graph for this function:

◆ compute_Pi_Helmholtz()

void compute_Pi_Helmholtz ( std::vector< double > &  energy_transfer,
const dataset source_data,
const std::vector< double > &  ulon,
const std::vector< double > &  ulat,
const std::vector< double > &  ulon_ulon,
const std::vector< double > &  ulon_ulat,
const std::vector< double > &  ulat_ulat,
const MPI_Comm  comm 
)

Compute the energy transfer through the current filter scale.

In particular, computes $ \rho_0 * ( u_i \tau_{ij,j} - (u_i \tau_{ij})_{,j} ) $

This computation is applied to the Cartesian velocity components

Parameters
[in,out]energy_transferwhere to store the computed values (array)
[in]source_datadataset class instance containing data (Psi, Phi, etc)
[in]ulon,ulatcoarse Spherical velocity components
[in]ulon_ulon,ulon_ulat,ulat_ulatcoarse velocity products (e.g. bar(u*v) )
[in]commMPI communicator object
Here is the call graph for this function:

◆ compute_Pi_shift_deriv()

void compute_Pi_shift_deriv ( std::vector< double > &  energy_transfer,
const dataset source_data,
const std::vector< double > &  ux,
const std::vector< double > &  uy,
const std::vector< double > &  uz,
const std::vector< double > &  uxux,
const std::vector< double > &  uxuy,
const std::vector< double > &  uxuz,
const std::vector< double > &  uyuy,
const std::vector< double > &  uyuz,
const std::vector< double > &  uzuz,
const MPI_Comm  comm 
)

Compute the energy transfer through the current filter scale.

In particular, computes $ \rho_0 * ( u_i \tau_{ij,j} - (u_i \tau_{ij})_{,j} ) $

This computation is applied to the Cartesian velocity components

Parameters
[in,out]energy_transferwhere to store the computed values (array)
[in]source_datadataset class instance containing data (Psi, Phi, etc)
[in]ux,uy,uzcoarse Cartesian velocity components
[in]uxux,uxuy,uxuz,uyuy,uyuz,uzuzcoarse velocity products (e.g. bar(u*v) )
[in]commMPI communicator object
Here is the call graph for this function:

◆ compute_spatial_average()

void compute_spatial_average ( std::vector< double > &  means,
const std::vector< double > &  field,
const std::vector< double > &  areas,
const int  Ntime,
const int  Ndepth,
const int  Nlat,
const int  Nlon,
const std::vector< bool > &  mask 
)

Given a field, compute the horizontal average.

The result, means, is a function of time and depth.

Parameters
[in,out]meanswhere to store horizontal spatial averages
[in]fieldarray of vector to be averaged
[in]areas2D array of cell areas
[in]Ntime,Ndepth,Nlat,Nlon(MPI-local) dimension sizes
[in]mask2D array to distinguish land from water
Here is the call graph for this function:
Here is the caller graph for this function:

◆ compute_vorticity()

void compute_vorticity ( std::vector< double > &  vort_r,
std::vector< double > &  vort_lon,
std::vector< double > &  vort_lat,
std::vector< double > &  vel_div,
std::vector< double > &  OkuboWeiss,
const dataset source_data,
const std::vector< double > &  u_r,
const std::vector< double > &  u_lon,
const std::vector< double > &  u_lat,
const MPI_Comm  comm 
)

Wrapper for computing vorticity.

This wrapper applies compute_vorticity_at_point() at each point in the grid. If the compiler flag COMP_VORT is set to false, then this procedure is skipped.

Parameters
[in,out]vort_r,vort_lon,vort_latwhere to store computed vorticity components (array)
[in,out]vel_divwhere to store computed velocity divergence (array)
[in,out]OkuboWeisswhere to store computed OkuboWeiss (array)
[in]source_datadataset class instance containing data (Psi, Phi, etc)
[in]u_r,u_lon,u_latvelocity components
[in]mask2D array to distinguish land from water
[in]commMPI communicator object
Here is the call graph for this function:
Here is the caller graph for this function:

◆ compute_vorticity_at_point()

void compute_vorticity_at_point ( double &  vort_r_tmp,
double &  vort_lon_tmp,
double &  vort_lat_tmp,
double &  div_tmp,
double &  OkuboWeiss_tmp,
const dataset source_data,
const std::vector< double > &  u_r,
const std::vector< double > &  u_lon,
const std::vector< double > &  u_lat,
const int  Itime,
const int  Idepth,
const int  Ilat,
const int  Ilon 
)

Compute the vorticity, divergnce, and OkuboWeiss at a point.

Since we're computing derivatives anyways, might as well get a few things out of it.

At the moment, only computes the vort_r component of vorticity

Parameters
[in,out]vort_r_tmp,vort_lon_tmp,vort_lat_tmpwhere to store vorticity components
[in,out]div_tmpwhere to store velocity divergence
[in,out]OkuboWeiss_tmpwhere to store OkuboWeiss result
[in]source_datadataset class instance containing data (Psi, Phi, etc)
[in]u_r,u_lon,u_latvelocity components
[in]Itime,Idepth,Ilat,IlonCurrent index in time and space
Here is the call graph for this function:
Here is the caller graph for this function:

◆ compute_Z()

void compute_Z ( std::vector< double > &  enstrophy_transfer,
const dataset source_data,
const std::vector< double > &  ux,
const std::vector< double > &  uy,
const std::vector< double > &  uz,
const std::vector< double > &  coarse_vort_r,
const std::vector< double > &  vort_ux,
const std::vector< double > &  vort_uy,
const std::vector< double > &  vort_uz,
const MPI_Comm  comm 
)

Compute the enstrophy transfer through the current filter scale.

In particular, computes $ \rho_0 * ( \omega \tau_{j,j} - (\omega \tau_{j})_{,j} ) $ where $ \tau_{j} = \overline{\omega u_j} - \overline{\omega}\, \overline{u_j} $

This computation is applied to the Cartesian velocity components

Parameters
[in,out]enstrophy_transferwhere to store the computed values (array)
[in]source_datadataset class instance containing data (Psi, Phi, etc)
[in]ux,uy,uzcoarse Cartesian velocity components
[in]coarse_vort_rcoarse radial vorticity
[in]vort_ux,vort_uy,vort_uzcoarse vort-velocity products (e.g. bar(omega * u_x) )
[in]commMPI communicator object
Here is the call graph for this function:
Here is the caller graph for this function:

◆ convert_coordinates()

void convert_coordinates ( std::vector< double > &  longitude,
std::vector< double > &  latitude 
)

Apply degree->radian conversion, if necessary. Applied in-place.

This is just to provide a clean wrapper to improve some case file readability. It applies a Cartesian check (via constants.hpp), so it is safe to include by default.

Parameters
[in,out]longitude,latitudeCoordinate grids to be converted, if appropriate

◆ depotential_temperature()

double depotential_temperature ( const double  p,
const double  theta 
)

Convert potential temperature to actual temperature following equation 1.154a from Vallis (p.36)

Parameters
[in]pPressure (Pascales)
[in]thetapotential temperature (Celcius)
Returns
non-potential temperature (Celcius)
Parameters
[in]ppressure at point
[in]thetapotential temperature at point (Celsius)

◆ distance()

double distance ( const double  lon1,
const double  lat1,
const double  lon2,
const double  lat2,
const double  Llon,
const double  Llat 
)

Compute the distance (in metres) between two points in the domain.

In spherical coordinates,computes the distance between two points on a spherical shell along a great circle. (see https://en.wikipedia.org/wiki/Great-circle_distance) It should avoid floating point issues for the grid scales that we're considering.

The last two arguments (Llon and Llat) give the physical length of the two dimensions. This is used in the case of periodic Cartesian grids. They are otherwise unused.

Parameters
[in]lon1,lat1coordinates for the first position
[in]lon2,lat2coordinates for the second position
[in]Llon,Llatphysical length of the dimensions
Returns
returns the distance (in metres) between two points.
Here is the caller graph for this function:

◆ equation_of_state()

double equation_of_state ( const double  T,
const double  S,
const double  p 
)

UNESCO 1981 equation of state.

Valid for 0 < T < 40 (Celcius), 0 < S < 42 (PSU)

https://link.springer.com/content/pdf/bbm%3A978-3-319-18908-6%2F1.pdf

Parameters
[in]Ttemperature (degrees Celcius)
[in]Ssalinity (PSU)
[in]ppressure (bars)
Returns
density
Parameters
[in]Ttemperature (degrees C)
[in]Ssalinity (PSU)
[in]ppressure (bars)

◆ filtering()

void filtering ( const dataset source_data,
const std::vector< double > &  scales,
const MPI_Comm  comm 
)

Main filtering driver.

This function is the main filtering driver. It sets up the appropriate loop sequences, calls the other funcations (velocity conversions), and calls the IO functionality.

Parameters
[in]source_datadataset class instance containing data (velocities, etc)
[in]scalesscales at which to filter the data
[in]commMPI communicator (default MPI_COMM_WORLD)
Here is the call graph for this function:

◆ filtering_helmholtz()

void filtering_helmholtz ( const dataset source_data,
const std::vector< double > &  scales,
const MPI_Comm  comm 
)

Main filtering driver for Helmholtz decomposed data.

This function is the main filtering driver. It sets up the appropriate loop sequences, calls the other funcations (velocity conversions), and calls the IO functionality.

Parameters
[in]source_datadataset class instance containing data (Psi, Phi, etc)
[in]scalesscales at which to filter the data
[in]commMPI communicator (default MPI_COMM_WORLD)
Here is the call graph for this function:

◆ get_lat_bounds()

void get_lat_bounds ( int &  LAT_lb,
int &  LAT_ub,
const std::vector< double > &  latitude,
const int  Ilat,
const double  scale 
)

Get latitude bounds for coarse-graining.

See the methods documentation (Methods) for more details.

Parameters
[in,out]LAT_lb,LAT_ublogical index bounds for the filtering domain
[in]latitude1D latitude vector
[in]Ilatlogical index of current position
[in]scalefiltering scale (metres)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ get_lon_bounds()

void get_lon_bounds ( int &  LON_lb,
int &  LON_ub,
const std::vector< double > &  longitude,
const int  Ilon,
const double  centre_lat,
const double  curr_lat,
const double  scale 
)

Get longitude integratation bounds for a specific latitude.

See the methods documentation (Methods) for more details.

Parameters
[in,out]LON_lb,LON_ublogical index bounds for the filtering domain
[in]longitude1D longitude vector
[in]Ilonlogical index of current position
[in]centre_latlatitude for the centre of the filtering kernel
[in]curr_latcurrent latitude
[in]scalefiltering scale (metres)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ get_omp_chunksize()

int get_omp_chunksize ( const int  Nlat,
const int  Nlon 
)

Compute openmp chunk size for OMP for loops.

The point of using a function for this is that it standardizes the chunk size across the different functions.

This assumes that the loop being split is a Lat-Lon loop.

Parameters
[in]Nlat,Nlonsize of the latitude and longitude grids
Returns
chunksize that can be passed to OpenMP pragma commands

◆ Index()

size_t Index ( const int  Itime,
const int  Idepth,
const int  Ilat,
const int  Ilon,
const int  Ntime,
const int  Ndepth,
const int  Nlat,
const int  Nlon 
)

Convenience tool to convert physical index (time, depth, lat, lon) to a logical index.

Index is a function to convert a four-point (physical) index (Itime, Idepth, Ilat, Ilon) into a one-point (logical) index to access the double arrays.

Assumes standard CF ordering: time-depth-lat-lon

Parameters
[in]Itime,Idepth,Ilat,Ilon4-indices to be converted to a 1-index
[in]Ntime,Ndepth,Nlat,Nlondimension sizes
Returns
The effective 1-index that corresponds to the 4-index tuplet
Here is the caller graph for this function:

◆ Index1to4()

void Index1to4 ( const size_t  index,
int &  Itime,
int &  Idepth,
int &  Ilat,
int &  Ilon,
const int &  Ntime,
const int &  Ndepth,
const int &  Nlat,
const int &  Nlon 
)

Convenience tool to convert logical index to physical index (time, depth, lat, lon).

Index is a function to convert a one-point (logical) index into a four-point (physical) index (Itime, Idepth, Ilat, Ilon) to access the double arrays.

Assumes standard CF ordering: time-depth-lat-lon

Parameters
[in]indexlogical index to be converted to a 4-index
[in,out]Itime,Idepth,Ilat,Ilon4-indices to be returned
[in]Ntime,Ndepth,Nlat,Nlondimension sizes
Here is the caller graph for this function:

◆ KE_from_vels()

void KE_from_vels ( std::vector< double > &  KE,
std::vector< double > *  u1,
std::vector< double > *  u2,
std::vector< double > *  u3,
const std::vector< bool > &  mask,
const double  rho0 
)

Compute (local) KE density from the provided velocities.

KE is simply computed pointwise as 0.5 * rho0 * (u1^2 + u2^2 + u3^2)

Parameters
[in,out]KEWhere to storer computed KE (array)
[in]u1,u2,u3Velocities (pointers to arrays)
[in]mask2D array to differentiate land from water
[in]rho0constant density
Here is the caller graph for this function:

◆ kernel()

double kernel ( const double  dist,
const double  scale 
)

Primary kernel function coarse-graining procedure (G in publications)

Parameters
[in]distancedistance for evaluating the kernel
[in]scalefilter scale (in metres)
Returns
The kernel value for a given distance and filter scale
Here is the caller graph for this function:

◆ kernel_alpha()

double kernel_alpha ( void  )

Compute alpha value for kernel (for baroclinic transfer)

Returns
Returns a 'size' coefficient for the kernel.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ mask_out_pole()

void mask_out_pole ( const std::vector< double > &  latitude,
std::vector< bool > &  mask,
const int  Ntime,
const int  Ndepth,
const int  Nlat,
const int  Nlon 
)

If the grid includes the poles (lat = 90 degrees), then mask it out.

Modifies mask in-place

Parameters
[in]latitudeLatitude grid
[in,out]maskMask array to differentiate land/water
[in]Ntime,Ndepth,Nlat,NlonDimension sizes
Here is the call graph for this function:

◆ print_compile_info()

void print_compile_info ( const std::vector< double > *  scales)

Print summary of compile-time variables.

This is triggered by passing –version to the executable.

Parameters
[in]scalesScales used for filtering

◆ roll_field()

void roll_field ( std::vector< double > &  field_to_roll,
const std::string  dimension,
const int  roll_count,
const int  Ntime,
const int  Ndepth,
const int  Nlat,
const int  Nlon 
)

Roll field along dimension.

Currently hard-coded to roll along lon dimension only.

Parameters
[in,out]field_to_rollfield to be rolled
[in]dimensiondimension to roll (currently must be "lon")
[in]roll_countammount by which to roll (currently must be 1)
[in]Ntime,Ndepth,Nlat,Nlon(MPI-local) dimension sizes
Here is the call graph for this function:

◆ vel_Cart_to_Spher()

void vel_Cart_to_Spher ( std::vector< double > &  u_r,
std::vector< double > &  u_lon,
std::vector< double > &  u_lat,
const std::vector< double > &  u_x,
const std::vector< double > &  u_y,
const std::vector< double > &  u_z,
const dataset source_data 
)

Wrapper that applies vel_Spher_to_Cart_at_point to every point in the domain.

Parameters
[in,out]u_r,u_lon,u_latComputed Spherical velocities
[in]u_x,u_y,u_zCartesian velocities to convert
[in]source_datadataset class instance containing data (Psi, Phi, etc)
Here is the call graph for this function:

◆ vel_Cart_to_Spher_at_point()

void vel_Cart_to_Spher_at_point ( double &  u_r,
double &  u_lon,
double &  u_lat,
const double  u_x,
const double  u_y,
const double  u_z,
const double  lon,
const double  lat 
)

Convert single Cartesian velocity to spherical velocity.

Convert Cartesian velocities to spherical velocities. (u_x, u_y, u_z) -> (u_r, u_lon, u_lat)

Note: we are using linear velocities (m/s), not angular (rad/s), so

\begin{eqnarray*} u_\lambda = r\cos(\phi)\cdot\hat{u}_\lambda \\ u_\phi = r\cdot\hat{u}_\phi \end{eqnarray*}

Note: we are still using a Spherical coordinate system, we are only converting the velocity fields.

Parameters
[in,out]u_r,u_lon,u_latComputed Spherical velocities
[in]u_x,u_y,u_zCartesian velocities to be converted
[in]lon,latcoordinates of the location of conversion
Here is the caller graph for this function:

◆ vel_Spher_to_Cart()

void vel_Spher_to_Cart ( std::vector< double > &  u_x,
std::vector< double > &  u_y,
std::vector< double > &  u_z,
const std::vector< double > &  u_r,
const std::vector< double > &  u_lon,
const std::vector< double > &  u_lat,
const dataset source_data 
)

Wrapper that applies vel_Spher_to_Cart_at_point to every point in the domain.

Parameters
[in,out]u_x,u_y,u_zComputed Cartesian velocities
[in]u_r,u_lon,u_latSpherical velocities to convert
[in]source_datadataset class instance containing data (Psi, Phi, etc)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ vel_Spher_to_Cart_at_point()

void vel_Spher_to_Cart_at_point ( double &  u_x,
double &  u_y,
double &  u_z,
const double  u_r,
const double  u_lon,
const double  u_lat,
const double  lon,
const double  lat 
)

Convert single spherical velocity to Cartesian velocity.

Convert Spherical velocities to Cartesian velocities. (u_r, u_lon, u_lat) -> (u_x, u_y, u_z)

Note: we are using linear velocities (m/s), not angular (rad/s), so

\begin{eqnarray*} u_\lambda = r\cos(\phi)\cdot\hat{u}_\lambda \\ u_\phi = r\cdot\hat{u}_\phi \end{eqnarray*}

Note: we are still using a Spherical coordinate system, we are only converting the velocity fields.

Parameters
[in,out]u_x,u_y,u_zComputed Cartesian velocities
[in]u_r,u_lon,u_latSpherical velocities to be converted
[in]lon,latcoordinates of the location of conversion
Here is the caller graph for this function: