New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Py-side SCF, PCM, EFP, Mol parsing 2 #953

Merged
merged 79 commits into from Aug 2, 2018

Conversation

Projects
None yet
10 participants
@loriab
Member

loriab commented Apr 10, 2018

Description

This replaces #919, which I don't want to overwrite just yet b/c it has history and some review comments. This is up-to-date wrt master (10 Apr) and has had a recent pyvib2 rebased through (9 Apr). All tests pass. This is not ready to look at b/c it still includes #834 . When that makes it into master, I'll rebase it through again and then this PR should only include fully-working py-side scf and mol parsing (create_molecule_from_string is gone in this code). This is my last sequence PR before splitting away qcdb. I suggest it form the start of v1.3 .

Since originally posted (+27k/–10k), this PR has had #834 subtracted from it (now +13k/–7k). It is built atop #965, so I wouldn't review until that gets pulled in and subtracted from this. What does need answering soonish is whether for 1.2/1.3 purposes, the remaining commits of this PR should be separated further.

Now this has had #965 rebased through, it should be stable for a good while. It now has size (+9k/–6k). There doesn't seem much interest in promoting the below to 1.2, so no need for further subdivision.

This has now had v1.2 rebased through so I sincerely hope it's been kneaded for the last time.

Todos: Py-side SCF/PCM/EFP (1.3)

  • Developer Interest
    • SCF iterations moved py-side
    • includes control of MOM, DIIS, damping, soscf, frac, efp
    • SCF finalize, stability, post-iterations, printing moved py-side
    • Rework a good bit of Wfn::HF::common_init() to minimize the convergence helper controllers that are c-side
    • A couple export/def bugs fixed in Wfn::HF and more fns exported and moved to public in class
    • PCMSolver interface reworked for py-side
    • libefp EFP object moved from c-side P::e.EFP to a py-side attached EFP attribute on the psi4.core.Molecule object.
    • Psi4 forgot about the c-side EFP object in bin, lib, globals, exports (still in options)
    • EFP/EFP calcs moved purely py-side
    • SCF/EFP calcs carry out operations on the EFP object py-side
    • efp_torque no longer in P::e but EFP TORQUE in P::e:arrays

Todos: Py-side Mol Parsing (could be split out for 1.2)

  • User Interest
    • molecule {...} will take xyz (element symbol or atomic number) or psi4 formats.
    • overall chgmult now specifiable through 1 3\n-- before any fragments.
    • can specify atomic number, mass number, mass, ghosting, extra-label in mol spec.
  • Developer Interest
    • *.Molecule forgot charge_specified and multiplicity_specified attributes. They were only used once in deciding mult defaulting in the SCF and that logic is now handled at mol parsing time.
    • *.Molecule forgot create_molecule_from_string() and good riddance. Also *.Molecule init_with_xyz (all the format variations can be handled by plain from_string. There's an optional dtype if you want to be specific.
    • Remove singlet/doublet defaulting logic in Wavefunction::common_init(), as this is taken care of by mol init
    • In keywords (read_options), EFP_POL, QMEFP_POL, EFP_POL_DAMPING, all POL --> IND.
    • Main molutil.set_geometry() that all psithon and psiapi calls use to process user molecule switch from c-side parsing to qcdb.molparse.from_string() parsing. This fn handles QM mol, EFP frag, and PubChem contents.
    • libefp accessed through PylibEFP and -DENABLE_libefp=ON activates both
    • psi4.core.Molecule.molecule_from_string(string) constructor calls py-side parsing. Can override string contents with fix_com, fix_orientation, fix_symmetry args. string can have QM mol, EFP frag, and pubchem contents.

Questions

  • @robertodr, feel free to elaborate PCMSolver changes
  • Everyone the previous PR #965 gets both Mol classes to a point where they can be instantiated by py-side parsing (or arrays, or dict). Do we want the second set of changes above where inputparser is calling py-side parsing into 1.2? It will require PylibEFP. SCF/EFP will remain c-side.
  • @robertodr #953 (comment) full test suite (incl. pcm) passes with form_F commented out

Checklist

Status

  • Ready for review
  • Ready for merge

@dgasmith dgasmith added this to the Psi4 1.3 milestone Apr 13, 2018

@loriab

This comment has been minimized.

Member

loriab commented Apr 15, 2018

Say, @psi-rking, does the 6 signify anything special here? Looks like another was confused, too (https://github.com/psi4/psi4/blob/master/psi4/src/psi4/optking/geom_gradients_io.cc#L376). If it was 3 * nat, it'd look a lot like Cfour's FCMFINAL file.

@JonathonMisiewicz

This comment has been minimized.

Contributor

JonathonMisiewicz commented Apr 16, 2018

I know that 6*natom is also found in the header of the second derivative output file from WDA's INTDER, but I'm still not sure why. In addition to my findif work, I'm using INTDER for my work at UGA, so I'm doubly interested in learning what it means.

Our anharm program calls this number NSTORE and NNSTOR, but it never does anything with them.

@psi-rking

This comment has been minimized.

Contributor

psi-rking commented Apr 16, 2018

I haven't the slightest idea! But I was mimicking the old file15 output:
http://vergil.chemistry.gatech.edu/psi/devel/progman/node47.html

The description here is plain as day "sixtimesn" but there is no explanation - a memory issue on an ancient architecture? An equally interesting question would be "who knows the answer to this question?"

@loriab

This comment has been minimized.

Member

loriab commented Apr 16, 2018

That's really great – it acknowledges our curiosity but does not satisfy it. Sounds like an office-by-office survey at CCQC would be required.

"The cartesian Hessian matrix is found in file15. The first line of this file gives the number of atoms (n) and, in case you are curious, six times the number of atoms (sixtimesn)."

@JonathonMisiewicz

This comment has been minimized.

Contributor

JonathonMisiewicz commented Apr 16, 2018

I'll look into this. There aren't many people who might know here, which (sadly) cuts the detective work short.

@psi-rking

This comment has been minimized.

Contributor

psi-rking commented Apr 16, 2018

@loriab

This comment has been minimized.

Member

loriab commented Apr 16, 2018

Origin aside, do I gather correctly that optking (or its successors; optking being the only visible consumer of the .hess file) isn't dependent on the 6 and if I could absorb the format into FCMFINAL, it'd be ok with your code, @psi-rking?

@lothian

This comment has been minimized.

Contributor

lothian commented Apr 16, 2018

@loriab loriab removed this from the Psi4 1.3 milestone Apr 16, 2018

@loriab loriab referenced this pull request Apr 16, 2018

Merged

Mol & Vib polishing #965

10 of 11 tasks complete
@psi-rking

This comment has been minimized.

Contributor

psi-rking commented Apr 16, 2018

@loriab

This comment has been minimized.

Member

loriab commented Apr 16, 2018

@psi-rking, I should explicitly generate both and do a diff, but the 3*nat vs. 6*nat is the sole difference I see between Psi and Cfour hessian formats. I was hoping to completely collapse the format names, but I guess I'll keep the 6 one generate-able (hondo or intder name? I'll have to confirm hondo involvement) in case it turns up to matter. Only one line of non-shared code. So the current plan is for both to be generate-able but Psi printing to switch to 3 for Cfour consistency.

@JonathonMisiewicz

This comment has been minimized.

Contributor

JonathonMisiewicz commented Apr 18, 2018

Dr. Allen said the convention predates him, and he doesn't know why it exists either. It looks like this bit of quantum chemistry lore has been lost.

@CDSherrill

This comment has been minimized.

Member

CDSherrill commented Apr 18, 2018

@loriab

This comment has been minimized.

Member

loriab commented Apr 18, 2018

fwiw, I searched in the HONDO manual and couldn't find the 6. But if anyone runs across a HSSFIL spec, take a look. I put a bounty comment on it in the code. Thanks for your investigations, @JonathonMisiewicz

@loriab

This comment has been minimized.

Member

loriab commented Jul 18, 2018

scfitertopy is ready for review. again.

# 0 X Y Z XX YY ZZ XY XZ YZ
prefacs = np.array([ 1, 1, 1, 1, 1/3, 1/3, 1/3, 2/3, 2/3, 2/3,
1/15, 1/15, 1/15, 3/15, 3/15, 3/15, 3/15, 3/15, 3/15, 6/15])

This comment has been minimized.

@andysim

andysim Jul 20, 2018

Member

I know we're dropping py2 support, but I don't think this division works correctly there. If we're dropping it now, that's totally find with me, else we should probably have 1.0/3, etc.

This comment has been minimized.

@loriab

loriab Jul 20, 2018

Member

it's got from __future__ import division so I think it's ok.

This comment has been minimized.

@andysim

andysim Jul 20, 2018

Member

Oh, good call. As you were!

outfile->Printf("\n");
outfile->Printf(" ---------------------------------------------------------\n");
outfile->Printf(" SCF\n");
outfile->Printf(" by Justin Turney, Rob Parrish, Andy Simmonett\n");

This comment has been minimized.

@andysim

andysim Jul 20, 2018

Member

This overhaul is more than enough justification to add your name to the list of authors of this chunk of code, Lori. Especially with the number of times you've had to rebase this branch! If we're going to continue the Psi3 convention of tagging the code segments, it may as well be accurate

This comment has been minimized.

@dgasmith

dgasmith Jul 20, 2018

Member

+1 here

This comment has been minimized.

@robertodr

robertodr Jul 23, 2018

Contributor

👍

This comment has been minimized.

@loriab

loriab Aug 2, 2018

Member

Thanks very much. I'll consider my banner-place mostly earned, and if I ever add truly new functionality, I will append. :-)

@andysim

I think this is going to be extremely helpful going forwards, and it seems like a number of people are already developing from this branch. I vote we get this one in sooner rather than later, to allow people to play with convergence schemes and embedding methods. LGTM!

@dgasmith

dgasmith approved these changes Jul 20, 2018 edited

This is going to be a great change. Overall LGTM, not showstoppers found- just a few nitpicks.

geom = pubchemre.sub(process_pubchem_command, geom)
geom = filter_comments(geom)
molecule = core.Molecule.create_molecule_from_string(geom)
molrec = qcdb.molparse.from_string(geom,

This comment has been minimized.

@dgasmith

dgasmith Jul 20, 2018

Member

Doesn't have to be in this PR, but I have been thinking it would be nice to give either a string or a JSON representation here. I think we have the tech to do both easily now with a isinstance switch.

This comment has been minimized.

@loriab

loriab Jul 20, 2018

Member

Certainly do-able for the Mol constructor. Just have to check if dict arg is a valid schmea repr here and call molparse.from_schema, which will have to be created by moving Molecule.from_schmea deeper into molparse.

When that's in place, can decide if psi4.geometry should expand its takings from {string} to {string, json}. I've no objection.

from psi4.driver.p4util.exceptions import ConvergenceError, ValidationError
from psi4 import core
from .efp import get_qm_atoms_opts, modify_Fock_permanent, modify_Fock_induced

This comment has been minimized.

@dgasmith

dgasmith Jul 20, 2018

Member

This could be a potential point to start adding a __all__ = [] to the top of these files. Looking at it, I dont think we actually want these functions getting out as they are "private" and attaching them to the class should still work fine. Just noticed that psi4.sys/psi4.os works and had me thinking more in these directions lately.

This comment has been minimized.

@loriab

loriab Jul 20, 2018

Member

__all__ only controls from psi4 import *, I believe. psi4.os will still be accessible and in dir(psi4). that's why I started deleting some temps https://github.com/psi4/psi4/blob/master/psi4/__init__.py#L73 . We can delete more. But since __all__ doesn't actually restrict access, for us it only serves as guide-to-eye, and I think public/_private can do for that.

Another tidbit I found is that program chokes if you from psi4 import * and __all__ = ['DO_NOT_WILD_IMPORT'] is present b/c that member is not in module. We could disable wild imports!

This comment has been minimized.

@dgasmith

dgasmith Jul 20, 2018

Member

Yea __all__ will only effect from psi4.submodule import *, at the same time we do not say from psi4.driver import os. os comes from many places, I was mostly thinking that this particular would be a good spot to start the __all__ revolution.

SCFE_old = 0.0
SCFE = 0.0
Drms = 0.0
while True:

This comment has been minimized.

@dgasmith

dgasmith Jul 20, 2018

Member

Stylistic, but would recommend:

for iter in range(core.get_option('SCF', 'MAXITER')):
    ...
else:
    raise ConvergenceError("""SCF iterations""", self.iteration_)

This comment has been minimized.

@loriab

loriab Jul 20, 2018

Member

In general I entirely agree with for/else over "while True"/break. But in this case where self.iterations_ is the true counter and we may not enter the loop at =0 (thus throwing off where to give up upon exceeding MAXITER), I wonder if having the shadow iter isn't more confusing.

SCFE += self.compute_E()
if efp_enabled:
global efp_Dt

This comment has been minimized.

@dgasmith

dgasmith Jul 20, 2018

Member

Can we do this in a non-global way? A user input could mess with this I believe.

This comment has been minimized.

@loriab

loriab Jul 20, 2018

Member

honestly, no. I'd like to, but I don't think it's possible (open to suggestions). There's a callback fn with a very specific signature that the C library is expecting. No way to get the density into it except by global (self is out b/c that's part of the signature; no trouble using member data if this is in a C++ class).

This comment has been minimized.

@dgasmith

dgasmith Jul 20, 2018

Member

Hmm, do you have the function to point me to?

# EFP: Add efp contribution to energy
efp_Dt = self.Da().clone()
efp_Dt.add(self.Db())
efp_wfn_dependent_energy = self.molecule().EFP.get_wavefunction_dependent_energy()

This comment has been minimized.

@dgasmith

dgasmith Jul 20, 2018

Member

SCFE += self.molecule...

outfile->Printf("\n");
outfile->Printf(" ---------------------------------------------------------\n");
outfile->Printf(" SCF\n");
outfile->Printf(" by Justin Turney, Rob Parrish, Andy Simmonett\n");

This comment has been minimized.

@dgasmith

dgasmith Jul 20, 2018

Member

+1 here

@@ -792,7 +783,7 @@ std::vector<SharedMatrix> RHF::cphf_solve(std::vector<SharedMatrix> x_vec, doubl
return ret_vec;
}
int RHF::soscf_update() {
int RHF::soscf_update(float soscf_conv, int soscf_min_iter, int soscf_max_iter, int soscf_print) {

This comment has been minimized.

@dgasmith

dgasmith Jul 20, 2018

Member

float -> double?

This comment has been minimized.

@loriab

loriab Jul 20, 2018

Member

I don't know what I was thinking

C_DGEMM('N','T',nso,nso,na,1.0,Ca[0],nmo,Ca[0],nmo,0.0,Da[0],nso);
C_DGEMM('N','T',nso,nso,nb,1.0,Ca[0],nmo,Ca[0],nmo,0.0,Db[0],nso);
if (na == 0) memset(static_cast<void*>(Da[0]), '\0', sizeof(double) * nso * nso);
if (nb == 0) memset(static_cast<void*>(Db[0]), '\0', sizeof(double) * nso * nso);

This comment has been minimized.

@dgasmith

dgasmith Jul 20, 2018

Member

Nitpick, but while I see it. Can we do Da_->zero(); Db_->zero(); above the nirrep_ loop and remove the memset?

}
int ROHF::soscf_update()
{
int ROHF::soscf_update(float soscf_conv, int soscf_min_iter, int soscf_max_iter, int soscf_print) {

This comment has been minimized.

@dgasmith

dgasmith Jul 20, 2018

Member

Ditto from above.

This comment has been minimized.

@loriab

loriab Jul 20, 2018

Member

There aren't two memset pairs in rohf, but I went ahead and applied pattern to cuhf, rhf, rohf, uhf

C_DGEMM('N','T',nso,nso,na,1.0,Ca[0],nmo,Ca[0],nmo,0.0,Da[0],nso);
C_DGEMM('N','T',nso,nso,nb,1.0,Cb[0],nmo,Cb[0],nmo,0.0,Db[0],nso);
if (na == 0) ::memset(static_cast<void*>(Da[0]), '\0', sizeof(double) * nso * nso);
if (nb == 0) ::memset(static_cast<void*>(Db[0]), '\0', sizeof(double) * nso * nso);

This comment has been minimized.

@dgasmith

dgasmith Jul 20, 2018

Member

Ditto from above.

This was referenced Jul 20, 2018

@PeterKraus PeterKraus referenced this pull request Jul 27, 2018

Merged

Do not freeze electrons on ghost atoms. #1109

5 of 5 tasks complete
@loriab

This comment has been minimized.

Member

loriab commented Aug 2, 2018

Ok, this is rebased and passing 435/437 full tests. Will know the full truth in half and hour. So potentially ready for merge.

EDIT: full tests pass

@jturney jturney merged commit fdbb220 into psi4:master Aug 2, 2018

1 check passed

continuous-integration/travis-ci/pr The Travis CI build passed
Details

@JonathonMisiewicz JonathonMisiewicz referenced this pull request Aug 3, 2018

Merged

FINDIF Overhaul #1024

24 of 24 tasks complete
@tyzhang1993

This comment has been minimized.

Contributor

tyzhang1993 commented on psi4/driver/molutil.py in 7b11d4f Aug 15, 2018

@loriab This raise from syntax is not supported in python2. Could you please consider changing it to a python2 friendly version? Thanks!

This comment has been minimized.

Member

loriab replied Aug 15, 2018

It's true that we haven't (much) explicitly stripped out py2 compatibility, but we have stopped supporting py2 past the v1.2 release. What are you working with that's py2-only?

This comment has been minimized.

Contributor

tyzhang1993 replied Aug 15, 2018

Thank you for the specification! I don't think we have code rely on py2. We will just need to upgrade the python on our cluster to a compatible version.

This comment has been minimized.

Member

loriab replied Aug 15, 2018

sounds good. I recommend Linux-source-3.6-nightly from http://vergil.chemistry.gatech.edu/nu-psicode/install-v1.2.1.html

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment