Head and wings corrections for k-dTDA #58
Conversation
| if ct.flags.f_contiguous: | ||
| c = np.asfortranarray(ct.transpose(order_ct)) | ||
| elif ct.flags.c_contiguous: | ||
| c = np.ascontiguousarray(ct.transpose(order_ct)) | ||
| else: | ||
| c = ct.transpose(order_ct) |
There was a problem hiding this comment.
can you explain these changes? does it just intend the output to be contiguous?
my original code allows ct to be in any contiguous order, and then doesn't therein constrain the contiguity of c. Sometimes a transpose can just be equivalent to reading the array with the opposite contiguity, if the transpose indices are c that isn't actually needed.
if we want to constrain the contiguity of the output, I think it would be better to just do
ct = ct.reshape(shape_ct, order="A")
c = ct.transpose(order_ct)
c = np.asarray(c, order="A")
There was a problem hiding this comment.
The intention is to try to maintain the contiguous nature of ct, particularly given that you can pass a and b as say c contiguous and then have an output that is not. This is particularly problematic with the ao2mo methods which check the contiguity. Happy though to implement it as mentioned in the comment.
There was a problem hiding this comment.
This doesn't quite work as the contiguity has already been changed with the transpose, so order="A" just keeps it as F contiguous
There was a problem hiding this comment.
that's fine -- the output doesn't even need to be contiguous really. That being said, I can't think of a situation when it won't be, since it's an output from either C or F implementation of DGEMM. at and bt are either, and the output is fine being either as well.
There was a problem hiding this comment.
I can't seem to see it in another PR, but just want to check this isn't change in another PR.
There was a problem hiding this comment.
i've not changed it, just revert this on this PR
momentGW/pbc/tda.py
Outdated
| if self.fc: | ||
| q = np.array([1e-3, 0, 0]).reshape(1, 3) | ||
| self.q_abs = self.kpts.cell.get_abs_kpts(q) | ||
| self.qij = self.build_pert_term(self.q_abs[0]) |
There was a problem hiding this comment.
q_abs is O(1) to compute right? If so, can we make it a property instead.
I'd prefer not storing qij in the class, but if it's needed in several places then maybe it's for the best. Maybe a more descriptive variable name?
There was a problem hiding this comment.
Changed the variable name, still not great but not sure what makes more sense.
momentGW/pbc/tda.py
Outdated
| wing_tmp = wing_tmp.real * 2 | ||
| wing_tmp *= -(np.sqrt(cell_vol / (4 * (np.pi**3))) * q0**2) | ||
|
|
||
| eta[kp, q][x, n] += util.einsum("p,pq->pq", wing_tmp, original) |
There was a problem hiding this comment.
Same comment as the moments_dd. Is it possible to have this function be unchanged, but before we do the convolution, there is something like:
if self.fc:
eta = self._add_eta_fc_correction(moments_dd)
# Construct the self-energy moments
moments_occ, moments_vir = self.convolve(eta)
There was a problem hiding this comment.
Have been able to do this and split moments_dd as well.
momentGW/pbc/tda.py
Outdated
| """ | ||
| kpts = self.kpts | ||
| head = np.zeros((self.nkpts, self.nmom_max + 1), dtype=object) | ||
| HW_const = np.sqrt(4.0 * np.pi) / np.linalg.norm(self.q_abs[0]) |
There was a problem hiding this comment.
no caps in the variable name
momentGW/pbc/tda.py
Outdated
| total_vol = self.kpts.cell.vol * self.nkpts | ||
| q0 = (6 * np.pi**2 / total_vol) ** (1 / 3) | ||
| norm_q_abs = np.linalg.norm(self.q_abs[0]) | ||
| HW_const = np.sqrt(4.0 * np.pi) / norm_q_abs |
| @@ -265,3 +422,8 @@ def kpts(self): | |||
| def nkpts(self): | |||
| """Get the number of k-points.""" | |||
| return self.gw.nkpts | |||
|
|
|||
| @property | |||
…from_master # Conflicts: # momentGW/gw.py # momentGW/pbc/gw.py
Implementation of the head and wings corrections for k-dTDA.
Due to a large number of PRs between when this was started and now this PR, have moved all commits relating to the HW corrections from a previous branch to a new branch from master.