# Determine DFT+U values by linear response

I am currently trying to determine how DFT+U values can be determined self consistently. I see that the bare and converged linear response matrix contain the difference in occupancies from the ground state. I am wondering if my understanding of the approach is correct.

1. Calculate ground state and record the final occupancies
2. Calculate perturbed state (bare) by applying U to the ground state charge density then taking the first set of occupancies from the first SCF step as the bare response matrix (bare - ground)
3. Converge the bare state and take the final set of occupancies as the converged response matrix (converged - ground)

Then you end up with two matrices, the bare response matrix and the converged response matrix. The U value is then taken as follows.

U = (Converged-1) - (Bare-1)

My question is maybe naive, but how does this final step of subtracting two matrices give a U value.

I am using this as a resource to understand what is going on in this method.

If it at all matters, I would like to use CASTEP or GPAW as my calculator.

The $$U$$ in this case is actually a matrix, since on Prof. Kulik's page that you linked, she is explaining multi-site $$U$$ calculations. The perturbation of a single site may also affect the response at a different site (e.g. from your link, perturbing Mn can affect O).
When doing the perturbations, you record the changes in all the occupations at different sites. Then when you do the matrix inversion these initial off-diagonal elements will have an affect on the resulting inverted matrix's diagonal elements. The diagonal elements $$U_{ii}$$ are of course the values of $$U$$ for site $$i$$. This type of calculation takes into account the effect of neighboring atoms on a given atom's response to the perturbation, as opposed to the single-site $$U$$ (which you could calculate separately for each site) that does not.