PAWpySeed
Parallel C/Python package for numerical analysis of PAW DFT wavefunctions
Functions
projector.h File Reference
#include "linalg.h"
Include dependency graph for projector.h:

Go to the source code of this file.

Functions

ppot_tget_projector_list (int num_els, int *labels, int *ls, double *wave_grids, double *projectors, double *aewaves, double *pswaves, double *rmaxs, double grid_encut)
 
real_proj_site_tprojector_values (int num_sites, int *labels, double *coords, double *lattice, double *reclattice, ppot_t *pps, int *fftg)
 
real_proj_site_tsmooth_pw_values (int num_N, int *Nlst, int *labels, double *coords, double *lattice, double *reclattice, ppot_t *pps, int *fftg)
 
void onto_projector_helper (band_t *band, double complex *x, real_proj_site_t *sites, int num_sites, double *lattice, double *reclattice, double *kpt, int num_cart_gridpts, int *fftg, projection_t *projections)
 
void get_aug_freqs_helper (band_t *band, double complex *x, real_proj_site_t *sites, int num_sites, double *lattice, double *reclattice, double *kpt, int num_cart_gridpts, int *fftg, projection_t *projections)
 
void onto_projector (kpoint_t *kpt, int band_num, real_proj_site_t *sites, int num_sites, int *G_bounds, double *lattice, double *reclattice, int num_cart_gridpts, int *fftg)
 
void onto_projector_ncl (kpoint_t *kpt, int band_num, real_proj_site_t *sites, int num_sites, int *G_bounds, double *lattice, double *reclattice, int num_cart_gridpts, int *fftg)
 
void onto_smoothpw (kpoint_t *kpt, int band_num, real_proj_site_t *sites, int num_sites, int *G_bounds, double *lattice, double *reclattice, int num_cart_gridpts, int *fftg)
 
void get_aug_freqs (kpoint_t *kpt, int band_num, real_proj_site_t *sites, int num_sites, int *G_bounds, double *lattice, double *reclattice, int num_cart_gridpts, int *fftg)
 
void add_num_cart_gridpts (ppot_t *pp_ptr, double *lattice, int *fftg)
 
void make_pwave_overlap_matrices (ppot_t *pp_ptr)
 
void setup_projections (pswf_t *wf, ppot_t *pps, int num_elems, int num_sites, int *fftg, int *labels, double *coords)
 
void overlap_setup_real (pswf_t *wf_R, pswf_t *wf_S, int *labels_R, int *labels_S, double *coords_R, double *coords_S, int *N_R, int *N_S, int *N_RS_R, int *N_RS_S, int num_N_R, int num_N_S, int num_N_RS)
 
void overlap_setup_recip (pswf_t *wf_R, pswf_t *wf_S, int *labels_R, int *labels_S, double *coords_R, double *coords_S, int *N_R, int *N_S, int *N_RS_R, int *N_RS_S, int num_N_R, int num_N_S, int num_N_RS)
 
void compensation_terms (double complex *overlap, int BAND_NUM, pswf_t *wf_S, pswf_t *wf_R, int num_M, int num_N_R, int num_N_S, int num_N_RS, int *M_R, int *M_S, int *N_R, int *N_S, int *N_RS_R, int *N_RS_S, int *proj_labels, double *proj_coords, int *ref_labels, double *ref_coords, int *fft_grid)
 
void compensation_terms_recip (double complex *overlap, int BAND_NUM, pswf_t *wf_S, pswf_t *wf_R, int num_M, int num_N_R, int num_N_S, int num_N_RS, int *M_R, int *M_S, int *N_R, int *N_S, int *N_RS_R, int *N_RS_S, int *proj_labels, double *proj_coords, int *ref_labels, double *ref_coords, int *fft_grid)
 
double * besselt (double *r, double *k, double *f, double encut, int N, int l)
 

Detailed Description

This file contains the functions needed to perform overlap operator evaulations between bands of different wavefunctions. It uses a similar algorithm to VASP for the onto_projector and projector_values evaluation.

Function Documentation

void add_num_cart_gridpts ( ppot_t pp_ptr,
double *  lattice,
int *  fftg 
)

Calculates the maximum number of grid points that can be contained within the projector sphere of pp_ptr given the lattice and fftg (FFT grid dimensions, 3D vector of int), and stores this value to pp_ptr.

double* besselt ( double *  r,
double *  k,
double *  f,
double  encut,
int  N,
int  l 
)

DEPRECATED

O(N^2) spherical Bessel transform

void compensation_terms ( double complex *  overlap,
int  BAND_NUM,
pswf_t wf_S,
pswf_t wf_R,
int  num_M,
int  num_N_R,
int  num_N_S,
int  num_N_RS,
int *  M_R,
int *  M_S,
int *  N_R,
int *  N_S,
int *  N_RS_R,
int *  N_RS_S,
int *  proj_labels,
double *  proj_coords,
int *  ref_labels,
double *  ref_coords,
int *  fft_grid 
)

Calculates the components of the overlap operator in the augmentation regions of each ion in the lattice.

void compensation_terms_recip ( double complex *  overlap,
int  BAND_NUM,
pswf_t wf_S,
pswf_t wf_R,
int  num_M,
int  num_N_R,
int  num_N_S,
int  num_N_RS,
int *  M_R,
int *  M_S,
int *  N_R,
int *  N_S,
int *  N_RS_R,
int *  N_RS_S,
int *  proj_labels,
double *  proj_coords,
int *  ref_labels,
double *  ref_coords,
int *  fft_grid 
)
void get_aug_freqs ( kpoint_t kpt,
int  band_num,
real_proj_site_t sites,
int  num_sites,
int *  G_bounds,
double *  lattice,
double *  reclattice,
int  num_cart_gridpts,
int *  fftg 
)
void get_aug_freqs_helper ( band_t band,
double complex *  x,
real_proj_site_t sites,
int  num_sites,
double *  lattice,
double *  reclattice,
double *  kpt,
int  num_cart_gridpts,
int *  fftg,
projection_t projections 
)
ppot_t* get_projector_list ( int  num_els,
int *  labels,
int *  ls,
double *  wave_grids,
double *  projectors,
double *  aewaves,
double *  pswaves,
double *  rmaxs,
double  grid_encut 
)

Returns a point to a list of ppot_t objects, one for each element in a POTCAR file. Called as a helper function by Wavefunction.make_c_projectors

void make_pwave_overlap_matrices ( ppot_t pp_ptr)

Calculates <phi_i|phi_j>, <phit_i|phit_j>, and <(phi_i-phit_i)|(phi_j-phit_j)> for all onsite i and j for the element represented by pp_ptr.

void onto_projector ( kpoint_t kpt,
int  band_num,
real_proj_site_t sites,
int  num_sites,
int *  G_bounds,
double *  lattice,
double *  reclattice,
int  num_cart_gridpts,
int *  fftg 
)

Calculates <p_i|psit_nk> for all i={R,epsilon,l,m} in the structure for one band

void onto_projector_helper ( band_t band,
double complex *  x,
real_proj_site_t sites,
int  num_sites,
double *  lattice,
double *  reclattice,
double *  kpt,
int  num_cart_gridpts,
int *  fftg,
projection_t projections 
)

Helper function for onto_projector, which performs the FFT of the wavefunction into real space and calculates from <p_i|psit_nk> from the grid points found in projector_values.

void onto_projector_ncl ( kpoint_t kpt,
int  band_num,
real_proj_site_t sites,
int  num_sites,
int *  G_bounds,
double *  lattice,
double *  reclattice,
int  num_cart_gridpts,
int *  fftg 
)
void onto_smoothpw ( kpoint_t kpt,
int  band_num,
real_proj_site_t sites,
int  num_sites,
int *  G_bounds,
double *  lattice,
double *  reclattice,
int  num_cart_gridpts,
int *  fftg 
)

Calculates <(phi_i-phit_i)|psit_nk> for all i={R,epsilon,l,m} in the structure for one band

void overlap_setup_real ( pswf_t wf_R,
pswf_t wf_S,
int *  labels_R,
int *  labels_S,
double *  coords_R,
double *  coords_S,
int *  N_R,
int *  N_S,
int *  N_RS_R,
int *  N_RS_S,
int  num_N_R,
int  num_N_S,
int  num_N_RS 
)

Much more efficient version of overlap_setup. Calculates three overlap terms for when bands have different structures: <(phi1_i-phit1_i)|psit2_n2k> <(phi2_i-phit2_i)|psit1_n1k> <(phi1_i-phit1_i)|(phi2_i-phit2_i)>

void overlap_setup_recip ( pswf_t wf_R,
pswf_t wf_S,
int *  labels_R,
int *  labels_S,
double *  coords_R,
double *  coords_S,
int *  N_R,
int *  N_S,
int *  N_RS_R,
int *  N_RS_S,
int  num_N_R,
int  num_N_S,
int  num_N_RS 
)
real_proj_site_t* projector_values ( int  num_sites,
int *  labels,
double *  coords,
double *  lattice,
double *  reclattice,
ppot_t pps,
int *  fftg 
)

Finds the coordinates on the FFT grid that fall within each projection sphere and stores the values of the projectors at that those points, as well as the coordinates of those points, in real_proj_site_t objects, one for each site in the structure.

void setup_projections ( pswf_t wf,
ppot_t pps,
int  num_elems,
int  num_sites,
int *  fftg,
int *  labels,
double *  coords 
)

Evaluates <p_i|psit_nk> for all bands and kpoints of wf.

real_proj_site_t* smooth_pw_values ( int  num_N,
int *  Nlst,
int *  labels,
double *  coords,
double *  lattice,
double *  reclattice,
ppot_t pps,
int *  fftg 
)

Finds the coordinates on the FFT grid that fall with each augmentation sphere (NOT projector sphere) and stores the values of the difference between the partial waves at those points, as well as the coordinates of those points, in real_proj_site_t objects, one for each site in the structure.