Processing math: 100%

gasci_pchb_doubles_select_particles Module



Contents


Variables

Type Visibility Attributes Name Initial
type(PCHB_ParticleSelection_vals_t), private, parameter :: PCHB_particle_selection_vals = PCHB_ParticleSelection_vals_t()

Abstract Interfaces

abstract interface

  • private pure function GetPgen_t(this, nI, i_sg, I, J)

    Arguments

    Type IntentOptional Attributes Name
    class(ParticleSelector_t), intent(in) :: this
    integer, intent(in) :: nI(nEl)

    The determinant in nI-format and the supergroup index

    integer, intent(in) :: i_sg

    The determinant in nI-format and the supergroup index

    integer, intent(in) :: I

    The particles.

    integer, intent(in) :: J

    The particles.

    Return Value real(kind=dp)

abstract interface

  • private subroutine Finalize_t(this)

    Arguments

    Type IntentOptional Attributes Name
    class(ParticleSelector_t), intent(inout) :: this

abstract interface

  • private subroutine Draw_t(this, nI, ilutI, i_sg, elecs, srcs, p)

    Arguments

    Type IntentOptional Attributes Name
    class(ParticleSelector_t), intent(in) :: this
    integer, intent(in) :: nI(nEl)

    The determinant in nI-format and the supergroup index

    integer(kind=n_int), intent(in) :: ilutI(0:nIfD)

    The determinant in bitmask format

    integer, intent(in) :: i_sg

    The determinant in nI-format and the supergroup index

    integer, intent(out) :: elecs(2)

    The chosen particles I,J and their index in nI. It is guaranteed that scrs(1) < srcs(2).

    integer, intent(out) :: srcs(2)

    The chosen particles I,J and their index in nI. It is guaranteed that scrs(1) < srcs(2).

    real(kind=dp), intent(out) :: p

    The probability of drawing p({I,J})|Di. This is the probability of drawing two particles from a given determinant Di regardless of order.


Derived Types

type, public, extends(EnumBase_t) ::  PCHB_ParticleSelection_t

Components

Type Visibility Attributes Name Initial
integer, public :: val
character(len=9), private :: str

Type-Bound Procedures

generic, public :: operator(==) => eq_EnumBase_t
generic, public :: operator(/=) => neq_EnumBase_t
procedure , public , :: to_str Function

type, public ::  PCHB_ParticleSelection_vals_t

Components

Type Visibility Attributes Name Initial
type(PCHB_ParticleSelection_t), public :: UNIF_UNIF = PCHB_ParticleSelection_t(1, 'UNIF-UNIF')

Both particles are drawn uniformly. We draw from p(I)|Di and then p(J|I)JDi and both probabilites come from the PCHB weighting scheme. We guarantee that I and J are occupied. We draw ˜p(I)|Di uniformly and then p(J|I)JDi The second distribution comes from the PCHB weighting scheme. We guarantee that I and J are occupied. We draw ˜p(I)|Di uniformly and then p(J|I)J. The second distribution comes from the PCHB weighting scheme. We guarantee that I is occupied.

type(PCHB_ParticleSelection_t), public :: FULL_FULL = PCHB_ParticleSelection_t(2, 'FULL-FULL')

Both particles are drawn uniformly. We draw from p(I)|Di and then p(J|I)JDi and both probabilites come from the PCHB weighting scheme. We guarantee that I and J are occupied. We draw ˜p(I)|Di uniformly and then p(J|I)JDi The second distribution comes from the PCHB weighting scheme. We guarantee that I and J are occupied. We draw ˜p(I)|Di uniformly and then p(J|I)J. The second distribution comes from the PCHB weighting scheme. We guarantee that I is occupied.

type(PCHB_ParticleSelection_t), public :: UNIF_FULL = PCHB_ParticleSelection_t(3, 'UNIF-FULL')

Both particles are drawn uniformly. We draw from p(I)|Di and then p(J|I)JDi and both probabilites come from the PCHB weighting scheme. We guarantee that I and J are occupied. We draw ˜p(I)|Di uniformly and then p(J|I)JDi The second distribution comes from the PCHB weighting scheme. We guarantee that I and J are occupied. We draw ˜p(I)|Di uniformly and then p(J|I)J. The second distribution comes from the PCHB weighting scheme. We guarantee that I is occupied.

type(PCHB_ParticleSelection_t), public :: UNIF_FAST = PCHB_ParticleSelection_t(4, 'UNIF-FAST')

Both particles are drawn uniformly. We draw from p(I)|Di and then p(J|I)JDi and both probabilites come from the PCHB weighting scheme. We guarantee that I and J are occupied. We draw ˜p(I)|Di uniformly and then p(J|I)JDi The second distribution comes from the PCHB weighting scheme. We guarantee that I and J are occupied. We draw ˜p(I)|Di uniformly and then p(J|I)J. The second distribution comes from the PCHB weighting scheme. We guarantee that I is occupied.

Type-Bound Procedures

procedure , public , nopass :: from_str => from_keyword Function

type, public ::  ParticleSelector_t

Type-Bound Procedures

procedure (Draw_t) , public :: draw
procedure (GetPgen_t) , public :: get_pgen
procedure (Finalize_t) , public :: finalize

type, public, extends(ParticleSelector_t) ::  UniformParticles_t

Type-Bound Procedures

procedure , public :: draw => draw_UniformParticles_t Subroutine
procedure , public :: get_pgen => get_pgen_UniformParticles_t Function
procedure , public :: finalize => finalize_UniformParticles_t Subroutine

type, private, extends(ParticleSelector_t) ::  PC_Particles_t

Components

Type Visibility Attributes Name Initial
type(AliasSampler_1D_t), private :: I_sampler

The shape is (n_supergroup) -> number_of_spin_orbs

type(AliasSampler_2D_t), private :: J_sampler

The shape is (number_of_spin_orbs, n_supergroup) -> number_of_spin_orbs

class(GASSpec_t), private, allocatable :: GAS_spec
type(SuperGroupIndexer_t), private, pointer :: indexer => null()
logical, public :: use_lookup = .false.

Use a lookup for the supergroup index in global_det_data.

logical, public :: create_lookup = .false.

Create and manage! the supergroup index lookup in global_det_data.

Type-Bound Procedures

procedure (Draw_t) , public :: draw
procedure (GetPgen_t) , public :: get_pgen
procedure , public :: init => init_PC_WeightedParticles_t Subroutine
procedure , public :: finalize => finalize_PC_WeightedParticles_t Subroutine

type, public, extends(PC_Particles_t) ::  PC_FullyWeightedParticles_t

Components

Type Visibility Attributes Name Initial
logical, public :: use_lookup = .false.

Use a lookup for the supergroup index in global_det_data.

logical, public :: create_lookup = .false.

Create and manage! the supergroup index lookup in global_det_data.

Type-Bound Procedures

procedure , public :: init => init_PC_WeightedParticles_t Subroutine
procedure , public :: finalize => finalize_PC_WeightedParticles_t Subroutine
procedure , public :: draw => draw_PC_FullyWeightedParticles_t Subroutine
procedure , public :: get_pgen => get_pgen_PC_FullyWeightedParticles_t Function

type, private, extends(PC_Particles_t) ::  PC_WeightedParticles_t

Components

Type Visibility Attributes Name Initial
logical, public :: use_lookup = .false.

Use a lookup for the supergroup index in global_det_data.

logical, public :: create_lookup = .false.

Create and manage! the supergroup index lookup in global_det_data.

Type-Bound Procedures

procedure , public :: init => init_PC_WeightedParticles_t Subroutine
procedure , public :: finalize => finalize_PC_WeightedParticles_t Subroutine
procedure , public :: draw => draw_PC_WeightedParticles_t Subroutine
procedure , public :: get_pgen => get_pgen_PC_WeightedParticles_t Function

type, public, extends(PC_Particles_t) ::  PC_FastWeightedParticles_t

Components

Type Visibility Attributes Name Initial
logical, public :: use_lookup = .false.

Use a lookup for the supergroup index in global_det_data.

logical, public :: create_lookup = .false.

Create and manage! the supergroup index lookup in global_det_data.

Type-Bound Procedures

procedure , public :: init => init_PC_WeightedParticles_t Subroutine
procedure , public :: finalize => finalize_PC_WeightedParticles_t Subroutine
procedure , public :: draw => draw_PC_FastWeightedParticles_t Subroutine
procedure , public :: get_pgen => get_pgen_PC_FastWeightedParticles_t Function

Functions

private pure function to_str(options) result(res)

Arguments

Type IntentOptional Attributes Name
class(PCHB_ParticleSelection_t), intent(in) :: options

Return Value character(len=9)

private pure function from_keyword(w) result(res)

Parse a given keyword into the possible particle selection schemes

Arguments

Type IntentOptional Attributes Name
character(len=*), intent(in) :: w

Return Value type(PCHB_ParticleSelection_t)

private pure function get_pgen_UniformParticles_t(this, nI, i_sg, I, J) result(p)

Calculates p({I,J})|Di

Read more…

Arguments

Type IntentOptional Attributes Name
class(UniformParticles_t), intent(in) :: this
integer, intent(in) :: nI(nEl)

The determinant in nI-format and the supergroup index

integer, intent(in) :: i_sg

The determinant in nI-format and the supergroup index

integer, intent(in) :: I

The particles.

integer, intent(in) :: J

The particles.

Return Value real(kind=dp)

private pure function get_pgen_PC_FullyWeightedParticles_t(this, nI, i_sg, I, J) result(p)

Calculates p({I,J})|Di

Read more…

Arguments

Type IntentOptional Attributes Name
class(PC_FullyWeightedParticles_t), intent(in) :: this
integer, intent(in) :: nI(nEl)

The determinant in nI-format and the supergroup index

integer, intent(in) :: i_sg

The determinant in nI-format and the supergroup index

integer, intent(in) :: I

The particles.

integer, intent(in) :: J

The particles.

Return Value real(kind=dp)

private pure function get_pgen_PC_WeightedParticles_t(this, nI, i_sg, I, J) result(p)

Calculates p({I,J})|Di

Read more…

Arguments

Type IntentOptional Attributes Name
class(PC_WeightedParticles_t), intent(in) :: this
integer, intent(in) :: nI(nEl)

The determinant in nI-format and the supergroup index

integer, intent(in) :: i_sg

The determinant in nI-format and the supergroup index

integer, intent(in) :: I

The particles.

integer, intent(in) :: J

The particles.

Return Value real(kind=dp)

private pure function get_pgen_PC_FastWeightedParticles_t(this, nI, i_sg, I, J) result(p)

Calculates p({I,J})|Di

Read more…

Arguments

Type IntentOptional Attributes Name
class(PC_FastWeightedParticles_t), intent(in) :: this
integer, intent(in) :: nI(nEl)

The determinant in nI-format and the supergroup index

integer, intent(in) :: i_sg

The determinant in nI-format and the supergroup index

integer, intent(in) :: I

The particles.

integer, intent(in) :: J

The particles.

Return Value real(kind=dp)

public function get_PCHB_weight(exc) result(res)

If there are three-body excitations, the double excitations actually become determinant-dependent. This function returns a fake determinant independent value in all cases, but evaluating exc for the reference determinant.

Arguments

Type IntentOptional Attributes Name
type(Excite_2_t), intent(in) :: exc

Return Value real(kind=dp)


Subroutines

public subroutine allocate_and_init(PCHB_particle_selection, GAS_spec, IJ_weights, rank_with_info, use_lookup, particle_selector)

Arguments

Type IntentOptional Attributes Name
type(PCHB_ParticleSelection_t), intent(in) :: PCHB_particle_selection
class(GASSpec_t), intent(in) :: GAS_spec
real(kind=dp), intent(in) :: IJ_weights(:,:,:)
integer, intent(in) :: rank_with_info

The intra-node rank that contains the weights to be used all other arr of all other ranks are ignored (and can be allocated with size 0).

logical, intent(in) :: use_lookup
class(ParticleSelector_t), intent(inout), allocatable :: particle_selector

private subroutine draw_UniformParticles_t(this, nI, ilutI, i_sg, elecs, srcs, p)

Arguments

Type IntentOptional Attributes Name
class(UniformParticles_t), intent(in) :: this
integer, intent(in) :: nI(nEl)

The determinant in nI-format and the supergroup index

integer(kind=n_int), intent(in) :: ilutI(0:nIfD)

The determinant in bitmask format

integer, intent(in) :: i_sg

The determinant in nI-format and the supergroup index

integer, intent(out) :: elecs(2)

The chosen particles I,J and their index in nI. It is guaranteed that scrs(1) < srcs(2).

integer, intent(out) :: srcs(2)

The chosen particles I,J and their index in nI. It is guaranteed that scrs(1) < srcs(2).

real(kind=dp), intent(out) :: p

The probability of drawing p({I,J})|Di. This is the probability of drawing two particles from a given determinant Di regardless of order.

private subroutine finalize_UniformParticles_t(this)

Arguments

Type IntentOptional Attributes Name
class(UniformParticles_t), intent(inout) :: this

private subroutine init_PC_WeightedParticles_t(this, GAS_spec, weights, rank_with_info, use_lookup, create_lookup)

Arguments

Type IntentOptional Attributes Name
class(PC_Particles_t), intent(inout) :: this
class(GASSpec_t), intent(in) :: GAS_spec
real(kind=dp), intent(in) :: weights(:,:,:)

w(J, I, i_sg), The weight of picking J after picking I under supergroup i_sg

integer, intent(in) :: rank_with_info

The intra-node rank that contains the weights to be used all other arr of all other ranks are ignored (and can be allocated with size 0).

logical, intent(in) :: use_lookup
logical, intent(in) :: create_lookup

private subroutine draw_PC_FullyWeightedParticles_t(this, nI, ilutI, i_sg, elecs, srcs, p)

Arguments

Type IntentOptional Attributes Name
class(PC_FullyWeightedParticles_t), intent(in) :: this
integer, intent(in) :: nI(nEl)

The determinant in nI-format and the supergroup index

integer(kind=n_int), intent(in) :: ilutI(0:nIfD)

The determinant in bitmask format

integer, intent(in) :: i_sg

The determinant in nI-format and the supergroup index

integer, intent(out) :: elecs(2)

The chosen particles I,J and their index in nI. It is guaranteed that scrs(1) < srcs(2).

integer, intent(out) :: srcs(2)

The chosen particles I,J and their index in nI. It is guaranteed that scrs(1) < srcs(2).

real(kind=dp), intent(out) :: p

The probability of drawing p({I,J})|Di. This is the probability of drawing two particles from a given determinant Di regardless of order.

private subroutine draw_PC_WeightedParticles_t(this, nI, ilutI, i_sg, elecs, srcs, p)

Arguments

Type IntentOptional Attributes Name
class(PC_WeightedParticles_t), intent(in) :: this
integer, intent(in) :: nI(nEl)

The determinant in nI-format and the supergroup index

integer(kind=n_int), intent(in) :: ilutI(0:nIfD)

The determinant in bitmask format

integer, intent(in) :: i_sg

The determinant in nI-format and the supergroup index

integer, intent(out) :: elecs(2)

The chosen particles I,J and their index in nI. It is guaranteed that scrs(1) < srcs(2).

integer, intent(out) :: srcs(2)

The chosen particles I,J and their index in nI. It is guaranteed that scrs(1) < srcs(2).

real(kind=dp), intent(out) :: p

The probability of drawing p({I,J})|Di. This is the probability of drawing two particles from a given determinant Di regardless of order.

private subroutine draw_PC_FastWeightedParticles_t(this, nI, ilutI, i_sg, elecs, srcs, p)

Arguments

Type IntentOptional Attributes Name
class(PC_FastWeightedParticles_t), intent(in) :: this
integer, intent(in) :: nI(nEl)

The determinant in nI-format and the supergroup index

integer(kind=n_int), intent(in) :: ilutI(0:nIfD)

The determinant in bitmask format

integer, intent(in) :: i_sg

The determinant in nI-format and the supergroup index

integer, intent(out) :: elecs(2)

The chosen particles I,J and their index in nI. It is guaranteed that scrs(1) < srcs(2).

integer, intent(out) :: srcs(2)

The chosen particles I,J and their index in nI. It is guaranteed that scrs(1) < srcs(2).

real(kind=dp), intent(out) :: p

The probability of drawing p({I,J})|Di. This is the probability of drawing two particles from a given determinant Di regardless of order.

private subroutine finalize_PC_WeightedParticles_t(this)

Arguments

Type IntentOptional Attributes Name
class(PC_Particles_t), intent(inout) :: this