dipc-cc / hubbard

Python tools for mean-field Hubbard models
https://dipc-cc.github.io/hubbard/
GNU Lesser General Public License v3.0
21 stars 8 forks source link

How to change the size of #133

Closed lidonglin12 closed 1 year ago

lidonglin12 commented 1 year ago

I want to change the size of the carbon atoms and delete the lable number, and I have changed the radius (on lines 187, 190, and 192 ) in plot.py, then, I run the example of anthracenes, but nothing has changed. I would be grateful if you could help me with this problem. Thank you in advance.

sofiasanz commented 1 year ago

To remove the atomic sites labels you have to remove line 35 in anthracenes/compute.py (p.annotate()). This is what is adding the labels on the atoms.

With respect to changing the radius of the carbon atoms, did you update your installation? Otherwise you would be using the last installation that does not have the change you did.

lidonglin12 commented 1 year ago

Thank you for your reply.

So, do you mean that for changing the radius of the carbon atoms, the position I modified is correct (on lines 187, 190, and 192 in plot.py), and it is because of the version that it does not work?

I have updated my installation, but it is still not working.

sofiasanz commented 1 year ago

Yes, the code lines that you have modified look correct.

What I mean is that if you do any changes in your local hubbard repository and you have installed the package by doing python -m pip install . or something similar, you should re-install again so the changes appear in your installed version (the one that python is going to use).

lidonglin12 commented 1 year ago

Yes, the code lines that you have modified look correct.

What I mean is that if you do any changes in your local hubbard repository and you have installed the package by doing python -m pip install . or something similar, you should re-install again so the changes appear in your installed version (the one that python is going to use).

Thank you very much for your reply. I tried to update my installation and my modification worked.

In addition, I have one more question. Can I modify the color range? For example, change the color range from (-0.4, 0.4) to (-0.3, 0.3) to make the spin polarization more pronounced.

sofiasanz commented 1 year ago

In addition, I have one more question. Can I modify the color range? For example, change the color range from (-0.4, 0.4) to (-0.3, 0.3) to make the spin polarization more pronounced.

Yes, of course. When plotting the spin polarization you can use the arguments vmax and vmin (that by default are set to +/-0.4) to control the colorbar scale.

For example, do something like plot.SpinPolarization(H, vmin=-0.3, vmax=0.3), where H is the HubbardHamiltonian object.

lidonglin12 commented 1 year ago

In addition, I have one more question. Can I modify the color range? For example, change the color range from (-0.4, 0.4) to (-0.3, 0.3) to make the spin polarization more pronounced.

Yes, of course. When plotting the spin polarization you can use the arguments vmax and vmin (that by default are set to +/-0.4) to control the colorbar scale.

For example, do something like plot.SpinPolarization(H, vmin=-0.3, vmax=0.3), where H is the HubbardHamiltonian object.

Thank you for your patience. I learned a lot. I tried to calculate the spin polarization of the same structure using two scripts ((examples/molecules/kondo-paper/Type-1/fig3.py) and (examples/molecules/ anthracenes/compute.py)), but I obtained different results. Just for spin polarization. So, in the case of the same U value, what is the difference between these two scripts?

sofiasanz commented 1 year ago

Thank you for your patience. I learned a lot.

No problem, glad to hear that!

I tried to calculate the spin polarization of the same structure using two scripts ((examples/molecules/kondo-paper/Type-1/fig3.py) and (examples/molecules/ anthracenes/compute.py)), but I obtained different results. Just for spin polarization. So, in the case of the same U value, what is the difference between these two scripts?

I don't understand this. How are you using these two different scripts (which calculate the spin polarization for different molecules) for the same molecule, and which molecule are you using?

lidonglin12 commented 1 year ago

Just put the type1.xyz file in the anthracenes folder and convert it to .XV file, then rename it 3-anthracene.XV and run the compute.py. Is it unreasonable to do so?

sofiasanz commented 1 year ago

Hi, sorry for the delay in my answer. If you want to use the script examples/molecules/anthracenes/compute.py for the type1.xyz molecule you should also change the line that sets the initial spin polarization for the molecule, as this is critical for the resulting self-consistent solution. For instance, in the script examples/molecules/anthracenes/compute.py the line that sets the polarization: H.set_polarization([1, 6, 15]) works specifically for the geometry of the 3-anthracenes.

Also, the script examples/molecules/kondo-paper/Type-1/fig3.py uses U=3.5 eV for the Coulomb repulsion parameter, while in examples/molecules/anthracenes/compute.py there is a loop over this parameter.

lidonglin12 commented 1 year ago

polarization

Thank you for your kind reply. So if I want to use the script (examples/molecules/anthracenes/compute.py) for a new molecule, besides U and t being tunable parameters, I also have to set the initial spin polarization, right? On what basis should I set the initial spin polarization? Thank you in advance for your patient explanation.

sofiasanz commented 1 year ago

I also have to set the initial spin polarization, right?

Correct. There can be several different local minima that you can find self-consistently. Finding the self-consistent solution is not a guarantee of finding the groundstate. In fact, you should compare the total energy of different solutions and pick the one with the lowest energy (global minima). The different local minima can be accessed by initializing the calculation with different spin density guesses.

On what basis should I set the initial spin polarization?

Unfortunately there is not a rule to know which spin polarization one should use. However there are a few things that you can take into account. For instance, if there are zigzag sites in your molecule, those probably host localized spins and is a good site to place an unpaired electron. On the other side, you can look into the low energy molecular orbitals (HOMO, HOMO-1, LUMO, LUMO+1, etc) and see where the molecular orbitals show localized features. There is another way to check where the radicals should be, using chemical arguments for these sp2 carbon structures. The idea behind this is by considering the Clar's sextets rule, which states that the gain of resonant rings can justify the appearance of radical sites. You can find a nice explanation to this in the supplemental material of Ref. Nature volume 10, Article number: 200 (2019) (section "Understanding the appearance of spins in the junctions from Clar’s theory").

With respect to the code, to setup the initial spin polarization you can make use of the following function before converging: H.set_spin_polarization(up, dn) Where up and dn are 1D array_like elements (e.g. lists, 1D-numpy arrays, etc) with the atomic indices where you want to place the up- and down-electrons, respectively.

lidonglin12 commented 1 year ago

The different local minima can be accessed by initializing the calculation with different spin density guesses.

You mean I should set different spin densities to optimize my structure and find the one with the lowest energy, right?

if there are zigzag sites in your molecule, those probably host localized spins and is a good site to place an unpaired electron.

After I get the the one with the lowest energy, I need to set the sites that probably host localized spins, like this, H.set_polarization([1, 6, 15]), the numbers represent the this three sites (carbon 1, 6, 15) that probably host localized spins, right?

H.set_spin_polarization(up, dn)

How does setting the spin up and down differ from setting the initial spin like this, H.set_polarization([1, 6, 15])?

Finally, can I verify the correctness of the calculation only by comparing it with the experiment?

Sorry for so many questions. Thank you very much.

sofiasanz commented 1 year ago

You mean I should set different spin densities to optimize my structure and find the one with the lowest energy, right?

Correct. By starting the calculation with different spin polarizations you can end up in different self-consistent solutions. You can check your initial spin configuration, e.g., by plotting the spin polarization after using the function H.set_polarization before convergence.

After I get the the one with the lowest energy, I need to set the sites that probably host localized spins, like this, H.set_polarization([1, 6, 15]), the numbers represent the this three sites (carbon 1, 6, 15) that probably host localized spins, right?

Yes, although by setting the initial spin density like that you are placing only localized spins as up-electrons at sites 1, 6, 15. The second argument of that function allows you to select atomic sites to place also down-electrons. For instance by doing H.set_polarization([1, 6], [15]), you would place two electrons up at sites 1,6 and one down at site 15.

Finally, can I verify the correctness of the calculation only by comparing it with the experiment?

No, of course you can (and should) compare to other methodologies. For instance you could run a DFT calculation, e.g., with SIESTA, and compare your results.

lidonglin12 commented 1 year ago

you would place two electrons up at sites 1,6 and one down at site 15.

For a new molecule, the spin distribution is unknown. I want to predict its spin distribution, and first I have to set the sites that probably host localized spins, and then how do I know the spin should be up or down?

In addition, I used the graphics visualization software (VESTA) to open the clar-goblet.xyz file in the example, if I understood correctly, the sites host the spin should be carbon 7 and 29, but in the script, the setting is H.set_polarization(up=[6], dn=[28]), I don't understand why? image

sofiasanz commented 1 year ago

Dear @lidonglin12,

For a new molecule, the spin distribution is unknown. I want to predict its spin distribution, and first I have to set the sites that probably host localized spins, and then how do I know the spin should be up or down?

You can maybe look at the low energy molecular orbitals (HOMO, HOMO-1, LUMO, LUMO+1, etc) and see where those orbitals localize, maybe this can give you an idea. The label up or down is not important, you can interchange them and obtain the same solution but with opposite names.

In the case of an even number of electrons in your system, I think that you should know at least whether your system is expected to have an antiferromagnetic-like solution or ferromagnetic-like solution. In the first case, you should start your spin density by placing one up electron at one radical site and one down electron at the other radical site. In the other case you can assume, for example, that you have larger number of up electrons than down, and then place only up electrons at the radical sites. To try to predict the total spin of the ground state, if your system is bipartite, then you can use the Lieb's theorem which states that the total spin of the ground state is S=0.5*|N_a - N_b|, where N_a and N_b are the number of a and b sublattices. If S=0 then the system is predicted to have an antiferromagnetic-like solution.

If you have an odd number of electrons then you can assume you have one more up electron and place it at the radical site.

In addition, I used the graphics visualization software (VESTA) to open the clar-goblet.xyz file in the example, if I understood correctly, the sites host the spin should be carbon 7 and 29, but in the script, the setting is H.set_polarization(up=[6], dn=[28]), I don't understand why?

Because python starts counting at 0 and the visualization software starts counting at 1.

sofiasanz commented 1 year ago

Dear @lidonglin12, can we close this issue?

lidonglin12 commented 1 year ago

Dear @lidonglin12, can we close this issue?

Yes, it can be closed. Thank you very much for your kind help.