Arbitrary External potentials

Jens · January 22, 2016, 6:18pm

Hello,

I would like to implement external potentials in the code, which cannot be expressed as point charges and thus have to integrated on a grid. I already had a good look at the code, but I am really unsure where the best place/object to integrate would be. As of now, I am thinking of using the Vxc grid for the purpose, the external potential would be read in from a file as well.

So the workflow would be
Run psi4
print vxc grid coordinates to
stop
Run psi4 again
read in potential
run DFT or whatever

As I am not extremely familiar with Psi4, a word of wisdom or advice would be very beneficial.

Otherwise great code and great project.

Cheers

jgonthier · January 26, 2016, 3:57pm

I am not familiar with the real-space integration grids part of the code, but I think you should be able to do all that in a single Psi4 run… Do you need to compute the external potential on the grid points with an external software ?

Jens · January 28, 2016, 10:25pm

Hello, yes that is exactly the case. Do you by accident know, if effective core potentials will be implemented any time soon?

Cheers Jens

sherrill · January 29, 2016, 3:39am

@crawdad might know about complicated external potentials from his work on PCM.

For simpler external potentials, Andy Simmonett has revamped our abilities and this will be included and documented in the next release.

Unfortunately no current plans to implement effective core potentials. We will do it eventually, but are very busy with other things right now.

David

susilehtola · January 2, 2018, 2:04pm

I realize this is an old discussion, but … ECPs have been already implemented. For an arbitrary external potential, the code should probably go somewhere around

github.com

psi4/psi4/blob/master/psi4/src/psi4/libscf_solver/hf.cc#L878


// FROZEN_UOCC takes precedence, FREEZE_UOCC directive has second priority
if (options_["FROZEN_UOCC"].has_changed()) {
    if (options_["FROZEN_UOCC"].size() != epsilon_a_->nirrep()) {
        throw PSIEXCEPTION("The FROZEN_UOCC array has the wrong dimensions");
    }
    for (int h = 0; h < epsilon_a_->nirrep(); h++) {
        frzvpi_[h] = options_["FROZEN_UOCC"][h].to_integer();
    }
} else {
    int nfzv = options_.get_int("NUM_FROZEN_UOCC");
    // Print out orbital energies.
    std::vector<std::pair<double, int> > pairs;
    for (int h = 0; h < epsilon_a_->nirrep(); ++h) {
        for (int i = 0; i < epsilon_a_->dimpi()[h]; ++i) pairs.push_back(std::make_pair(epsilon_a_->get(h, i), h));
        frzvpi_[h] = 0;
    }
    sort(pairs.begin(), pairs.end(), std::greater<std::pair<double, int> >());


    for (int i = 0; i < nfzv; ++i) frzvpi_[pairs[i].second]++;
}
}