FlowSieve
3.4.0
FlowSieve Coarse-Graining Documentation
|
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"
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 () |
Collection of all computation-related functions.
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().
[in,out] | coarse_val | where to store filtered value |
[in] | fields | fields to filter |
[in] | source_data | dataset class instance containing data (Psi, Phi, etc) |
[in] | Itime,Idepth,Ilat,Ilon | current position in time dimension |
[in] | LAT_lb,LAT_ub | lower/upper boundd on latitude for kernel |
[in] | scale | filtering scale |
[in] | use_mask | array of booleans indicating whether or not to use mask (i.e. zero out land) or to use the array value |
[in] | local_kernel | pre-computed kernel (NULL indicates not provided) |
[in] | weight | pointer to spatial weight (i.e. rho) (NULL indicates not provided) |
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$
[in,out] | uxux_tmp | where to store filtered (u_x)*(u_x) |
[in,out] | uxuy_tmp | where to store filtered (u_x)*(u_y) |
[in,out] | uxuz_tmp | where to store filtered (u_x)*(u_z) |
[in,out] | uyuy_tmp | where to store filtered (u_y)*(u_y) |
[in,out] | uyuz_tmp | where to store filtered (u_y)*(u_z) |
[in,out] | uzuz_tmp | where to store filtered (u_z)*(u_z) |
[in,out] | vort_ux_tmp | where to store filtered (vort_r)*(u_x) |
[in,out] | vort_uy_tmp | where to store filtered (vort_r)*(u_y) |
[in,out] | vort_uz_tmp | where to store filtered (vort_r)*(u_z) |
[in] | u_x,u_y,u_z | fields to filter |
[in] | vort_r | vorticity field to filter |
[in] | source_data | dataset class instance containing data (Psi, Phi, etc) |
[in] | Itime,Idepth,Ilat,Ilon | current position in time dimension |
[in] | LAT_lb,LAT_ub | lower/upper boundd on latitude for kernel |
[in] | scale | filtering scale |
[in] | local_kernel | pre-computed kernel (NULL indicates not provided) |
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).
[in,out] | areas | array in which areas will be stored |
[in] | longitude | longitude vector (1D) |
[in] | latitude | latitude vector (1D) |
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:
NOT implemented:
* 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) ) * ) *
[in,out] | div_J | where to store the computed values (array) |
[in] | u_x,u_y,u_z | coarse Cartesian velocity components |
[in] | uxux,uxuy,uxuz,uyuy,uyuz,uzuz | coarse velocity products (e.g. bar(u*v) ) |
[in] | coarse_p | coarse pressure |
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.
[in,out] | local_kernel | where to store the local kernel |
[in] | scale | Filtering scale |
[in] | source_data | dataset class instance containing data (Psi, Phi, etc) |
[in] | Ilat,Ilon | reference coordinate (kernel centre) |
[in] | LAT_lb,LAT_ub | upper and lower latitudinal bounds for kernel |
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
This computation is applied to the Cartesian velocity components
[in,out] | energy_transfer | where to store the computed values (array) |
[in] | source_data | dataset class instance containing data (Psi, Phi, etc) |
[in] | ux,uy,uz | coarse Cartesian velocity components |
[in] | uxux,uxuy,uxuz,uyuy,uyuz,uzuz | coarse velocity products (e.g. bar(u*v) ) |
[in] | comm | MPI communicator object |
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
This computation is applied to the Cartesian velocity components
[in,out] | energy_transfer | where to store the computed values (array) |
[in] | source_data | dataset class instance containing data (Psi, Phi, etc) |
[in] | ulon,ulat | coarse Spherical velocity components |
[in] | ulon_ulon,ulon_ulat,ulat_ulat | coarse velocity products (e.g. bar(u*v) ) |
[in] | comm | MPI communicator object |
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
This computation is applied to the Cartesian velocity components
[in,out] | energy_transfer | where to store the computed values (array) |
[in] | source_data | dataset class instance containing data (Psi, Phi, etc) |
[in] | ux,uy,uz | coarse Cartesian velocity components |
[in] | uxux,uxuy,uxuz,uyuy,uyuz,uzuz | coarse velocity products (e.g. bar(u*v) ) |
[in] | comm | MPI communicator object |
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.
[in,out] | means | where to store horizontal spatial averages |
[in] | field | array of vector to be averaged |
[in] | areas | 2D array of cell areas |
[in] | Ntime,Ndepth,Nlat,Nlon | (MPI-local) dimension sizes |
[in] | mask | 2D array to distinguish land from water |
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.
[in,out] | vort_r,vort_lon,vort_lat | where to store computed vorticity components (array) |
[in,out] | vel_div | where to store computed velocity divergence (array) |
[in,out] | OkuboWeiss | where to store computed OkuboWeiss (array) |
[in] | source_data | dataset class instance containing data (Psi, Phi, etc) |
[in] | u_r,u_lon,u_lat | velocity components |
[in] | mask | 2D array to distinguish land from water |
[in] | comm | MPI communicator object |
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
[in,out] | vort_r_tmp,vort_lon_tmp,vort_lat_tmp | where to store vorticity components |
[in,out] | div_tmp | where to store velocity divergence |
[in,out] | OkuboWeiss_tmp | where to store OkuboWeiss result |
[in] | source_data | dataset class instance containing data (Psi, Phi, etc) |
[in] | u_r,u_lon,u_lat | velocity components |
[in] | Itime,Idepth,Ilat,Ilon | Current index in time and space |
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 where
This computation is applied to the Cartesian velocity components
[in,out] | enstrophy_transfer | where to store the computed values (array) |
[in] | source_data | dataset class instance containing data (Psi, Phi, etc) |
[in] | ux,uy,uz | coarse Cartesian velocity components |
[in] | coarse_vort_r | coarse radial vorticity |
[in] | vort_ux,vort_uy,vort_uz | coarse vort-velocity products (e.g. bar(omega * u_x) ) |
[in] | comm | MPI communicator object |
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.
[in,out] | longitude,latitude | Coordinate grids to be converted, if appropriate |
double depotential_temperature | ( | const double | p, |
const double | theta | ||
) |
Convert potential temperature to actual temperature following equation 1.154a from Vallis (p.36)
[in] | p | Pressure (Pascales) |
[in] | theta | potential temperature (Celcius) |
[in] | p | pressure at point |
[in] | theta | potential temperature at point (Celsius) |
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.
[in] | lon1,lat1 | coordinates for the first position |
[in] | lon2,lat2 | coordinates for the second position |
[in] | Llon,Llat | physical length of the dimensions |
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
[in] | T | temperature (degrees Celcius) |
[in] | S | salinity (PSU) |
[in] | p | pressure (bars) |
[in] | T | temperature (degrees C) |
[in] | S | salinity (PSU) |
[in] | p | pressure (bars) |
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.
[in] | source_data | dataset class instance containing data (velocities, etc) |
[in] | scales | scales at which to filter the data |
[in] | comm | MPI communicator (default MPI_COMM_WORLD) |
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.
[in] | source_data | dataset class instance containing data (Psi, Phi, etc) |
[in] | scales | scales at which to filter the data |
[in] | comm | MPI communicator (default MPI_COMM_WORLD) |
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.
[in,out] | LAT_lb,LAT_ub | logical index bounds for the filtering domain |
[in] | latitude | 1D latitude vector |
[in] | Ilat | logical index of current position |
[in] | scale | filtering scale (metres) |
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.
[in,out] | LON_lb,LON_ub | logical index bounds for the filtering domain |
[in] | longitude | 1D longitude vector |
[in] | Ilon | logical index of current position |
[in] | centre_lat | latitude for the centre of the filtering kernel |
[in] | curr_lat | current latitude |
[in] | scale | filtering scale (metres) |
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.
[in] | Nlat,Nlon | size of the latitude and longitude grids |
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
[in] | Itime,Idepth,Ilat,Ilon | 4-indices to be converted to a 1-index |
[in] | Ntime,Ndepth,Nlat,Nlon | dimension sizes |
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
[in] | index | logical index to be converted to a 4-index |
[in,out] | Itime,Idepth,Ilat,Ilon | 4-indices to be returned |
[in] | Ntime,Ndepth,Nlat,Nlon | dimension sizes |
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)
[in,out] | KE | Where to storer computed KE (array) |
[in] | u1,u2,u3 | Velocities (pointers to arrays) |
[in] | mask | 2D array to differentiate land from water |
[in] | rho0 | constant density |
double kernel | ( | const double | dist, |
const double | scale | ||
) |
Primary kernel function coarse-graining procedure (G in publications)
[in] | distance | distance for evaluating the kernel |
[in] | scale | filter scale (in metres) |
double kernel_alpha | ( | void | ) |
Compute alpha value for kernel (for baroclinic transfer)
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
[in] | latitude | Latitude grid |
[in,out] | mask | Mask array to differentiate land/water |
[in] | Ntime,Ndepth,Nlat,Nlon | Dimension sizes |
void print_compile_info | ( | const std::vector< double > * | scales | ) |
Print summary of compile-time variables.
This is triggered by passing –version to the executable.
[in] | scales | Scales used for filtering |
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.
[in,out] | field_to_roll | field to be rolled |
[in] | dimension | dimension to roll (currently must be "lon") |
[in] | roll_count | ammount by which to roll (currently must be 1) |
[in] | Ntime,Ndepth,Nlat,Nlon | (MPI-local) dimension sizes |
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.
[in,out] | u_r,u_lon,u_lat | Computed Spherical velocities |
[in] | u_x,u_y,u_z | Cartesian velocities to convert |
[in] | source_data | dataset class instance containing data (Psi, Phi, etc) |
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
Note: we are still using a Spherical coordinate system, we are only converting the velocity fields.
[in,out] | u_r,u_lon,u_lat | Computed Spherical velocities |
[in] | u_x,u_y,u_z | Cartesian velocities to be converted |
[in] | lon,lat | coordinates of the location of conversion |
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.
[in,out] | u_x,u_y,u_z | Computed Cartesian velocities |
[in] | u_r,u_lon,u_lat | Spherical velocities to convert |
[in] | source_data | dataset class instance containing data (Psi, Phi, etc) |
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
Note: we are still using a Spherical coordinate system, we are only converting the velocity fields.
[in,out] | u_x,u_y,u_z | Computed Cartesian velocities |
[in] | u_r,u_lon,u_lat | Spherical velocities to be converted |
[in] | lon,lat | coordinates of the location of conversion |