ptr1 => ptr2 T ptr1, ptr2; ptr1 = ptr2; ptr1 = ptr2 ptr1 = ptr2;
subroutine init_fcimc_fn_pointers() character(*), parameter :: this_routine = "init_fcimc_fn_pointers" ! Almost all excitation generators in NECI are Full CI generators. gen_all_excits => gen_all_excits_default ! Select the excitation generator. if (tSD_GAS()) then call class_managed(generate_excitation, gen_all_excits) else if (t_3_body_excits .and. .not. (t_mol_3_body .or. t_ueg_3_body)) then if (t_uniform_excits) then generate_excitation => gen_excit_uniform_k_space_hub_transcorr else generate_excitation => gen_excit_k_space_hub_transcorr end if else if (t_ueg_3_body) then if (tTrcorrExgen) then generate_two_body_excitation => gen_ueg_excit else if (TLatticeGens) then generate_two_body_excitation => gen_rand_excit endif generate_excitation => gen_excit_mol_tc elseif(t_impurity_excitgen) then generate_excitation => gen_excit_impurity_model elseif ((t_back_spawn_option .or. t_back_spawn_flex_option)) then if (tHUB .and. tLatticeGens) then ! for now the hubbard + back-spawn still uses the old ! genrand excit gen generate_excitation => gen_excit_back_spawn_hubbard else if (tUEGNewGenerator .and. tLatticeGens) then generate_excitation => gen_excit_back_spawn_ueg_new else if (tUEG .and. tLatticeGens) then generate_excitation => gen_excit_back_spawn_ueg else generate_excitation => gen_excit_back_spawn end if else if (tUEGNewGenerator) then generate_excitation => gen_ueg_excit else if (tPickVirtUniform) then ! pick-uniform-random-mag is on if (tReltvy) then generate_excitation => gen_rand_excit_Ex_Mag else call stop_all(this_routine, "Excitation generator has not been set!") end if else if (tGenHelWeighted) then generate_excitation => gen_excit_hel_weighted else if (tGen_4ind_2) then generate_excitation => gen_excit_4ind_weighted2 else if (tGen_4ind_weighted) then generate_excitation => gen_excit_4ind_weighted else if (tGen_4ind_reverse) then generate_excitation => gen_excit_4ind_reverse else if (tGUGA) then if (tgen_guga_crude) then if (t_k_space_hubbard) then generate_excitation => gen_excit_k_space_hub else if (t_new_real_space_hubbard) then generate_excitation => gen_excit_rs_hubbard else generate_excitation => gen_excit_4ind_weighted2 end if else generate_excitation => generate_excitation_guga end if else if (t_pcpp_excitgen) then generate_excitation => gen_rand_excit_pcpp else if (t_fci_pchb_excitgen) then call class_managed(generate_excitation, gen_all_excits) else if (t_k_space_hubbard) then if (t_3_body_excits) then if (t_uniform_excits) then generate_excitation => gen_excit_uniform_k_space_hub_transcorr else if (t_mixed_excits) then generate_excitation => gen_excit_mixed_k_space_hub_transcorr else generate_excitation => gen_excit_k_space_hub_transcorr end if else if (t_uniform_excits) then generate_excitation => gen_excit_uniform_k_space_hub else generate_excitation => gen_excit_k_space_hub end if end if else if (t_new_real_space_hubbard) then if (t_trans_corr_hop) then if (t_hole_focus_excits) then generate_excitation => gen_excit_rs_hubbard_transcorr_hole_focus else if (t_uniform_excits) then generate_excitation => gen_excit_rs_hubbard_transcorr_uniform else generate_excitation => gen_excit_rs_hubbard_transcorr end if else if (t_spin_dependent_transcorr) then generate_excitation => gen_excit_rs_hubbard_spin_dependent_transcorr else generate_excitation => gen_excit_rs_hubbard end if else if (t_tJ_model) then generate_excitation => gen_excit_tJ_model else if (t_heisenberg_model) then generate_excitation => gen_excit_heisenberg_model else generate_excitation => gen_rand_excit end if ! yes, fortran pointers work this way ! Pointer assignment with => ! Fortran !> ptr1 => ptr2 ! C !> T *ptr1, *ptr2; !> ptr1 = ptr2; ! Copy/value assignment with = ! Fortran !> ptr1 = ptr2 ! C !> *ptr1 = *ptr2; ! if we are using the 3-body excitation generator, embed the chosen excitgen ! in the three-body one if (t_mol_3_body) then generate_two_body_excitation => generate_excitation generate_excitation => gen_excit_mol_tc end if ! Do the same for HPHF if (tHPHF) then exc_generator_for_HPHF => generate_excitation generate_excitation => gen_hphf_excit end if ! In the main loop, we only need to find out if a determinant is ! connected to the reference det or not (so no ex. level above 2 is ! required). Except in some cases where we need to know the maximum ! excitation level if (tTruncSpace .or. tHistSpawn .or. tCalcFCIMCPsi) then max_calc_ex_level = nel else if (t_3_body_excits) then max_calc_ex_level = 3 else max_calc_ex_level = 2 end if end if ! How many children should we spawn given an excitation? if (t_real_time_fciqmc) then attempt_create => attempt_create_realtime else if (tTruncCas .or. tTruncSpace .or. & tPartFreezeCore .or. tPartFreezeVirt .or. tFixLz .or. & (tUEG .and. .not. tLatticeGens) .or. tTruncNOpen .or. t_trunc_nopen_diff) then if (tHPHF .or. tSemiStochastic) then attempt_create => attempt_create_trunc_spawn else attempt_create => att_create_trunc_spawn_enc end if else attempt_create => attempt_create_normal end if ! In attempt create, how should we evaluate the off diagonal matrix ! elements between a parent and its (potentially) spawned offspring? if (tHPHF) then if (tGenMatHEL) then get_spawn_helement => hphf_spawn_sign else get_spawn_helement => hphf_off_diag_helement_spawn end if ! new guga addition: do not need to recalculate Helement else if (tGUGA) then ! use hphf_routine also, since it does exactly what needed get_spawn_helement => hphf_spawn_sign else get_spawn_helement => get_helement_det_only end if ! When calling routines to generate all possible connections, this ! routine is called to generate the corresponding Hamiltonian matrix ! elements. if (tHPHF) then get_conn_helement => hphf_off_diag_helement_spawn else get_conn_helement => get_helement_det_only end if ! Once we have generated the children, do we need to encode them? if (.not. (tHPHF .or. tGen_4ind_weighted .or. tGUGA)) then encode_child => FindExcitBitDet else encode_child => null_encode_child end if ! What message should we display for a particle bloom? if (tAddToInitiator) then bloom_warn_string = & '("Bloom of more than n_add on ", a, " excit: A max of ", f10.2, " particles created. ", i8, " blooms occurred.")' else ! Use this variable to store the bloom cutoff level. InitiatorWalkNo = 3.0_dp bloom_warn_string = & '("Bloom of more than 3 on ", a, " excit: A max of ", f10.2, " particles created. ", i8, " blooms occurred.")' end if bloom_max = 0 if (tPreCond) then attempt_die => attempt_die_precond else attempt_die => attempt_die_normal end if extract_bit_rep_avsign => extract_bit_rep_avsign_no_rdm fill_rdm_diag_currdet => fill_rdm_diag_currdet_norm select case (sfTag) case (0) scaleFunction => powerScaleFunction case (1) scaleFunction => expScaleFunction case (2) scaleFunction => negScaleFunction case (3) scaleFunction => expCOScaleFunction case default call stop_all(this_routine, "Invalid scale function specified") end select if (tLinearAdaptiveShift) then shiftFactorFunction => linearShiftFactorFunction else if (tAutoAdaptiveShift) then shiftFactorFunction => autoShiftFactorFunction else shiftFactorFunction => constShiftFactorFunction end if ! select the procedure that returns all connected determinants. if (t_k_space_hubbard) then gen_all_excits => gen_all_excits_k_space_hubbard else if (t_new_real_space_hubbard) then gen_all_excits => gen_all_excits_r_space_hubbard end if end subroutine init_fcimc_fn_pointers