core_mesh.cpp File Reference

basic functions to work with mesh More...

#include "core_mesh.h"

Include dependency graph for core_mesh.cpp:

This graph shows which files directly or indirectly include this file:

Go to the source code of this file.

Macros
#define	ORDER   1

Functions
template<typename T >
static int	sgn (T val)
void	get_k_vec (size_t N, size_t index, int *k_vec)
void	get_k_vec (size_t N, size_t index, Vec_3D< int > &k_vec)
double	get_k_sq (size_t N, size_t index)
template<typename T >
static std::enable_if< std::is_integral< T >::value, T >::type	get_per (T vec, size_t per)
template<typename T >
static std::enable_if<!std::is_integral< T >::value, T >::type	get_per (T vec, size_t per)
template<typename T >
void	get_per (Vec_3D< T > &position, size_t per)
template<typename T >
void	get_per (Vec_3D< T > &position, size_t perx, size_t pery, size_t perz)
void	get_per (std::vector< Particle_v< double >> &particles, const size_t per)
template<typename T >
T	get_distance_1D (const T x_1, const T x_2, const size_t per)
template<typename T >
T	get_distance_1D (const T x_1, const int x_2, const size_t per)
template<typename T >
T	get_sgn_distance_1D (T x_from, T x_to, size_t per)
double	get_distance (const Vec_3D< double > &x_1, const Vec_3D< double > &x_2, int per)
Vec_3D< double >	get_sgn_distance (const Vec_3D< double > &x_from, const Vec_3D< double > &x_to, size_t per)
template<size_t N>
static double	wgh_sch (const Vec_3D< double > &x, const Vec_3D< int > &y, size_t mesh_num)
template<>
double	wgh_sch< 1 > (const Vec_3D< double > &x, const Vec_3D< int > &y, size_t mesh_num)
template<>
double	wgh_sch< 2 > (const Vec_3D< double > &x, const Vec_3D< int > &y, size_t mesh_num)
void	assign_to (Mesh &field, const Vec_3D< double > &position, const double value)
void	assign_to (std::vector< Mesh > &field, const Vec_3D< double > &position, const Vec_3D< double > &value)
void	assign_from (const Mesh &field, const Vec_3D< double > &position, double &value, double mod)
void	assign_from (const std::vector< Mesh > &field, const Vec_3D< double > &position, Vec_3D< double > &value, double mod)
void	fftw_execute_dft_r2c (const FFTW_PLAN_TYPE &p_F, Mesh &rho)
	compute forward (real to complex) FFT on mesh (inplace) More...
void	fftw_execute_dft_c2r (const FFTW_PLAN_TYPE &p_B, Mesh &rho)
	compute backward (complex to real) FFT on mesh (inplace) More...
void	fftw_execute_dft_r2c_triple (const FFTW_PLAN_TYPE &p_F, std::vector< Mesh > &rho)
	compute three forward (real to complex) FFTs on vector of meshes (inplace) More...
void	fftw_execute_dft_c2r_triple (const FFTW_PLAN_TYPE &p_B, std::vector< Mesh > &rho)
	compute three backward (complex to real) FFTs on vector of meshes (inplace) More...
template void	get_per (Vec_3D< int > &, size_t)
template void	get_per (Vec_3D< int > &, size_t, size_t, size_t)
template void	get_per (Vec_3D< size_t > &, size_t)
template void	get_per (Vec_3D< size_t > &, size_t, size_t, size_t)
template void	get_per (Vec_3D< double > &, size_t)
template void	get_per (Vec_3D< double > &, size_t, size_t, size_t)

Detailed Description

basic functions to work with mesh

Author

Michal Vrastil

Date

2018-07-11

Definition in file core_mesh.cpp.

Macro Definition Documentation

#define ORDER   1

Definition at line 12 of file core_mesh.cpp.

Function Documentation

void assign_from	(	const Mesh &	field,
		const Vec_3D< double > &	position,
		double &	value,
		double	mod
	)

Definition at line 224 of file core_mesh.cpp.

References IT< points >::iter(), Mesh::N, and IT< points >::vec.

Referenced by kick_step_no_momentum(), App_Var_Chi::ChiImpl::kick_step_w_chi(), App_Var_Chi::ChiImpl::kick_step_w_chi_no_momentum(), kick_step_w_momentum(), kick_step_w_momentum_pm(), and anonymous_namespace{mod_frozen_potential.cpp}::kick_step_w_pp().

{
    IT<ORDER+1> it(position);
    do{
        #pragma omp atomic
        value += field(it.vec) * mod * wgh_sch<ORDER>(position, it.vec, field.N);
    } while( it.iter() );
}

void assign_from	(	const std::vector< Mesh > &	field,
		const Vec_3D< double > &	position,
		Vec_3D< double > &	value,
		double	mod
	)

< reuse the same weigh for every field in std::vector

Definition at line 233 of file core_mesh.cpp.

References IT< points >::iter(), IT< points >::vec, and w.

{
    IT<ORDER+1> it(position);
    FTYPE_t w;
    do{
        w = mod * wgh_sch<ORDER>(position, it.vec, field[0].N); 
        for (unsigned int i = 0; i < 3; i++)
        {
            #pragma omp atomic
            value[i] += field[i](it.vec) * w;
        }
    } while( it.iter() );
}

void assign_to	(	Mesh &	field,
		const Vec_3D< double > &	position,
		const double	value
	)

Definition at line 201 of file core_mesh.cpp.

References IT< points >::iter(), Mesh::N, and IT< points >::vec.

Referenced by get_rho_from_par(), and get_vel_from_par().

{
    IT<ORDER+1> it(position);
    do{
        #pragma omp atomic
        field(it.vec) += value * wgh_sch<ORDER>(position, it.vec, field.N);
    } while( it.iter() );
}

void assign_to	(	std::vector< Mesh > &	field,
		const Vec_3D< double > &	position,
		const Vec_3D< double > &	value
	)

< reuse the same weigh for every field in std::vector

Definition at line 210 of file core_mesh.cpp.

References IT< points >::iter(), IT< points >::vec, and w.

{
    IT<ORDER+1> it(position);
    FTYPE_t w;
    do{
        w = wgh_sch<ORDER>(position, it.vec, field[0].N); 
        for (size_t i = 0; i < 3; i++)
        {
            #pragma omp atomic
            field[i](it.vec) += value[i] * w;
        }
    } while( it.iter() );
}

void fftw_execute_dft_c2r	(	const FFTW_PLAN_TYPE &	p_B,
		Mesh &	rho
	)

compute backward (complex to real) FFT on mesh (inplace)

Parameters

p_B	plan for backward transformation
rho	mesh upon which the transformation is performed

Definition at line 253 of file core_mesh.cpp.

References Mesh::complex(), FFTW_EXEC_C2R, and Mesh_base< T >::real().

Referenced by fftw_execute_dft_c2r_triple(), anonymous_namespace{test_chameleon.cpp}::get_grav_pot(), App_Var_AA::pot_corr(), and anonymous_namespace{chameleon.cpp}::ChiSolver< T >::set_linear_sol_at_level().

{
    FFTW_EXEC_C2R(p_B, rho.complex(), rho.real());
}

void fftw_execute_dft_c2r_triple	(	const FFTW_PLAN_TYPE &	p_B,
		std::vector< Mesh > &	rho
	)

compute three backward (complex to real) FFTs on vector of meshes (inplace)

Parameters

p_B	plan for backward transformations
rho	vector of meshes upon which the transformations are performed

Definition at line 263 of file core_mesh.cpp.

References fftw_execute_dft_c2r(), and get_per().

Referenced by App_Var_AA::AAImpl::aa_convolution(), App_Var_Chi::ChiImpl::get_chi_force(), App_Var_FP_mod::pot_corr(), App_Var< T >::pot_corr(), and App_Var< T >::Impl< T >::set_init_cond().

{
    for (unsigned int i = 0; i < 3; i++) fftw_execute_dft_c2r(p_B, rho[i]);
}

void fftw_execute_dft_r2c	(	const FFTW_PLAN_TYPE &	p_F,
		Mesh &	rho
	)

compute forward (real to complex) FFT on mesh (inplace)

Parameters

p_F	plan for forward transformation
rho	mesh upon which the transformation is performed

Definition at line 247 of file core_mesh.cpp.

References Mesh::complex(), FFTW_EXEC_R2C, Mesh::N, pow(), and Mesh_base< T >::real().

Referenced by App_Var_AA::AAImpl::aa_convolution(), fftw_execute_dft_r2c_triple(), App_Var_Chi::ChiImpl::gen_pow_spec_binned(), gen_rho_dist_k(), App_Var< T >::Impl< T >::get_binned_power_spec(), App_Var_Chi::ChiImpl::get_chi_force(), anonymous_namespace{test_chameleon.cpp}::get_grav_pot(), anonymous_namespace{chameleon.cpp}::ChiSolver< T >::set_linear_sol_at_level(), and App_Var_PM::upd_pos().

{
    FFTW_EXEC_R2C(p_F, rho.real(), rho.complex());
    rho /= pow((FTYPE_t)rho.N, 3); //< normalization
}

void fftw_execute_dft_r2c_triple	(	const FFTW_PLAN_TYPE &	p_F,
		std::vector< Mesh > &	rho
	)

compute three forward (real to complex) FFTs on vector of meshes (inplace)

Parameters

p_F	plan for forward transformations
rho	vector of meshes upon which the transformations are performed

Definition at line 258 of file core_mesh.cpp.

References fftw_execute_dft_r2c().

Referenced by App_Var< T >::Impl< T >::print_vel_pwr().

{
    for (unsigned int i = 0; i < 3; i++) fftw_execute_dft_r2c(p_F, rho[i]);
}

double get_distance	(	const Vec_3D< double > &	x_1,
		const Vec_3D< double > &	x_2,
		int	per
	)

Definition at line 115 of file core_mesh.cpp.

References get_distance_1D(), pow2(), and sqrt().

{
    FTYPE_t dst = 0;
    for (unsigned int i = 0; i < 3; i++) dst+= pow2(get_distance_1D(x_1[i], x_2[i], per));
    return sqrt(dst);
}

template<typename T >

T get_distance_1D ( const T  x_1, const T  x_2, const size_t  per  )

inline

Definition at line 90 of file core_mesh.cpp.

References growth_allz::T.

Referenced by get_distance(), get_distance_1D(), wgh_sch< 1 >(), and wgh_sch< 2 >().

{   
    T d = std::abs(x_2 - x_1);
    if (d <= per / 2) return d; // most probable, first condition
    else if (d <= per) return per - d;
    else d = fmod(d, T(per)); // fmod unlikely to evaluate, speed up code
    if (d <= per / 2) return d;
    else return per - d;
}

template<typename T >

T get_distance_1D ( const T  x_1, const int  x_2, const size_t  per  )

inline

Definition at line 101 of file core_mesh.cpp.

References get_distance_1D(), and growth_allz::T.

{ return get_distance_1D(x_1, (T)x_2, per); }

double get_k_sq	(	size_t	N,
		size_t	index
	)

Definition at line 38 of file core_mesh.cpp.

References get_k_vec(), and pow2().

Referenced by gen_pot_k(), gen_rho_w_pow_k(), and anonymous_namespace{chameleon.cpp}::ChiSolver< T >::get_chi_k().

{
    int k_vec[3];
    FTYPE_t tmp = 0;
    get_k_vec(N, index, k_vec);
    for (unsigned int i = 0; i<3; i++) tmp += pow2(k_vec[i]);
    return tmp;
}

void get_k_vec	(	size_t	N,
		size_t	index,
		int *	k_vec
	)

Definition at line 20 of file core_mesh.cpp.

Referenced by gen_displ_k_S2(), get_k_sq(), pwr_spec_k(), pwr_spec_k_init(), and vel_pwr_spec_k().

{
    k_vec[0] = index / ((N/2 + 1)*N);
    k_vec[1] = (index / (N/2 + 1)) % N;
    k_vec[2] = index % (N/2 + 1);
    
    for (unsigned int i = 0; i<2; i++) if (k_vec[i] > N/2) k_vec[i] -= N; // k_vec[2] is ALWAYS less or equal than N/2 (real FFTW)
}

void get_k_vec	(	size_t	N,
		size_t	index,
		Vec_3D< int > &	k_vec
	)

Definition at line 29 of file core_mesh.cpp.

{
    k_vec[0] = index / ((N/2 + 1)*N);
    k_vec[1] = (index / (N/2 + 1)) % N;
    k_vec[2] = index % (N/2 + 1);
    
    for (unsigned int i = 0; i<2; i++) if (k_vec[i] > N/2) k_vec[i] -= N; // k_vec[2] is ALWAYS less or equal than N/2 (real FFTW)
}

template<typename T >

static std::enable_if<std::is_integral<T>::value, T>::type get_per ( T  vec, size_t  per  )

static

Definition at line 48 of file core_mesh.cpp.

References Catch::Generators::value().

Referenced by fftw_execute_dft_c2r_triple(), get_per(), set_pert_pos(), stream_kick_stream(), and TEST_CASE().

{
    static_assert(std::is_integral<T>::value, "Integral required.");
    if ((vec >= (int)per) || (vec < 0) ){
        vec %= (int)per;
        return (vec < 0) ? vec + (int)per : vec;
    }
    else return vec;
}

template<typename T >

static std::enable_if<!std::is_integral<T>::value, T>::type get_per ( T  vec, size_t  per  )

static

Definition at line 59 of file core_mesh.cpp.

References Catch::Generators::value().

{
    static_assert(!std::is_integral<T>::value, "Non-integral required.");
    return ((vec >= per) || (vec < 0) ) ? vec - per * floor( vec / per ) : vec;
}

template<typename T >

void get_per	(	Vec_3D< T > &	position,
		size_t	per
	)

Definition at line 66 of file core_mesh.cpp.

References get_per(), and growth_allz::T.

Referenced by Mesh_base< size_t >::operator()(), and Mesh::operator()().

{
    for (T& pos : position) pos = get_per(pos, per);
}

template<typename T >

void get_per	(	Vec_3D< T > &	position,
		size_t	perx,
		size_t	pery,
		size_t	perz
	)

Definition at line 72 of file core_mesh.cpp.

References get_per().

{
    position[0] = get_per(position[0], perx);
    position[1] = get_per(position[1], pery);
    position[2] = get_per(position[2], perz);
}

void get_per	(	std::vector< Particle_v< double >> &	particles,
		const size_t	per
	)

Definition at line 79 of file core_mesh.cpp.

References get_per().

{
    const size_t Np = particles.size();
    #pragma omp parallel for
    for (size_t i = 0; i < Np; i++)
    {
        get_per(particles[i].position, per);
    }
}

template void get_per	(	Vec_3D< int > &	,
		size_t
	)

template void get_per	(	Vec_3D< int > &	,
		size_t	,
		size_t	,
		size_t
	)

template void get_per	(	Vec_3D< size_t > &	,
		size_t
	)

template void get_per	(	Vec_3D< size_t > &	,
		size_t	,
		size_t	,
		size_t
	)

template void get_per	(	Vec_3D< double > &	,
		size_t
	)

template void get_per	(	Vec_3D< double > &	,
		size_t	,
		size_t	,
		size_t
	)

Vec_3D<double> get_sgn_distance	(	const Vec_3D< double > &	x_from,
		const Vec_3D< double > &	x_to,
		size_t	per
	)

Definition at line 122 of file core_mesh.cpp.

References get_sgn_distance_1D(), wgh_sch(), and x.

Referenced by anonymous_namespace{mod_frozen_potential.cpp}::force_short().

{
    Vec_3D<FTYPE_t> dst;
    for (unsigned int i = 0; i < 3; i++) dst[i] = get_sgn_distance_1D(x_from[i], x_to[i], per);
    return dst;
}

template<typename T >

T get_sgn_distance_1D ( T  x_from, T  x_to, size_t  per  )

inline

Definition at line 104 of file core_mesh.cpp.

References sgn(), and growth_allz::T.

Referenced by get_sgn_distance().

{   // return signed (oriented) distance, e.g. x_to - x_from (no periodicity)
    // periodicity makes this a little bit trickier
    T d = x_from - x_to;
    if (std::abs(d) <= per / 2) return d; // most probable, first condition
    else if (std::abs(d) <= per) return d-sgn(d)*per;
    else d = fmod(d, T(per)); // fmod unlikely to evaluate, speed up code
    if (std::abs(d) <= per / 2) return d;
    else return d-sgn(d)*per;
}

template<typename T >

static int sgn ( T  val )

static

Definition at line 15 of file core_mesh.cpp.

References growth_allz::T.

Referenced by get_sgn_distance_1D(), and TEST_CASE().

{
    return (T(0) < val) - (val < T(0));
}

template<size_t N>

static double wgh_sch ( const Vec_3D< double > &  x, const Vec_3D< int > &  y, size_t  mesh_num  )

static

Referenced by get_sgn_distance().

template<>

double wgh_sch< 1 >	(	const Vec_3D< double > &	x,
		const Vec_3D< int > &	y,
		size_t	mesh_num
	)

Definition at line 140 of file core_mesh.cpp.

References get_distance_1D(), and w.

{ // CIC: Cloud in cells
    FTYPE_t dx, w = 1;
    for (unsigned int i = 0; i < 3; i++)
    {
        dx = get_distance_1D(x[i], y[i], mesh_num);
        if (dx > 1) return 0;
        else w *= 1 - dx;
    }
    return w;
}

template<>

double wgh_sch< 2 >	(	const Vec_3D< double > &	x,
		const Vec_3D< int > &	y,
		size_t	mesh_num
	)

Definition at line 152 of file core_mesh.cpp.

References get_distance_1D(), pow2(), and w.

{ // TSC: Triangular shaped clouds
    FTYPE_t dx, w = 1;
    for (unsigned int i = 0; i < 3; i++)
    {
        dx = get_distance_1D(x[i], y[i], mesh_num);
        if (dx > 1.5) return 0;
        else if (dx > 0.5) w *= pow2(1.5 - dx) / 2;
        else w *= 3 / FTYPE_t(4) - dx*dx;
    }
    return w;
}