Re: [Wien] Optical properties with SO coupling
Dear Gavin, I will describe my observation: I have calculated optical (epzz) and magneto-optical (K, for example K=epxy for M001) spectra of permittivity elements for bcc Fe. The electronic structure calculations are spin polarized, with spin-orbit, run by commands: runsp_lapw -p runsp_lapw -p -so -cc 0.01 -ec 0.001 -s lapw1 x lapw2 -p -fermi -up -so x optic -p -up -so x joint -p -up For w2k version 16.1, the calculated spectra corresponds to the experimental spectra (for both epzz and K). For w2k version 17.1, the calculated spectra are half-value for both epzz and K, compared to the experiment. Figures comparing spectra are here: http://alma.karlov.mff.cuni.cz/hamrle/w2kfig/Fe_Imepzz_compare.pdf http://alma.karlov.mff.cuni.cz/hamrle/w2kfig/Fe_ReK_compare.pdf In this example, I used permittivity spectra read directly from case.jointup files (I do not use output of kram). In the figures: - solid (noisy) line is output from case.jointup - the symbols are smeared spectra - black '+' are the experimental spectra - blue 'o' and red 'x' are spectra calculated by w2k version 17 - green '+' and yellow '*' are spectra calculated by w2k version 16 - y-axis denotes permittivity*E (in eV). That is why I have concluded that joint function in w2k version 17 has a bug in calculation of the optical permittivity. But I have not tested non-magnetic cases, I did it only for bcc Fe (sp+so). Hoping it helps. If I can help more, please let me know.. With my regards Jaro On 04/10/17 16:40, Gavin Abo wrote: Dear Jaro, I thought the spin-polarized SO optic normalization was broken in older versions of WIEN2k and was fixed in 17.1: https://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/msg16011.html Is it still broken? Kind Regards, Gavin On 10/3/2017 4:30 PM, Jaroslav Hamrle wrote: Hallo, to calculate optical properties of Ni, after calculating electronic structure being spin-polarized and being with spin-orbit, do: 1) create both case.inop (your file looks correct) and case.injoint Example of case.injoint is: example of case.injoint === 1 : LOWER,UPPER and (optional) UPPER-VAL BANDINDEX 0. 0.00100 1. : EMIN DE EMAX FOR ENERGYGRID IN ryd eV : output units eV / ryd / cm-1 4 : SWITCH 9 : NUMBER OF COLUMNS 0.1 0.1 0.3 : BROADENING (FOR DRUDE MODEL - switch 6,7 - ONLY) SWITCH: 0...JOINTDOS FOR EACH BAND COMBINATION 1...JOINTDOS AS SUM OVER ALL BAND COMBINATIONS 2...DOS FOR EACH BAND 3...DOS AS SUM OVER ALL BANDS 4...Im(EPSILON) 5...Im(EPSILON) for each band combination 6...INTRABAND contributions 7...INTRABAND contributions including band analysis = end example case.injoint Now, you have to decide if you want to calculate optics at finer k-mesh than electronic structure, or the same mesh. In case electronic structure is calculated with k-mesh 30x30x30, it is good enough for Imxy and MOKE. 2a) when keeping the same k-mesh for optical calculations as for electronic structure, do: x lapw2 -p -fermi -up -so x optic -p -up -so (your command in your email is opticc, i.e. complex variant of command optic; opticc should be used when the structure lacks point symmetry, which Ni does not) x joint -up x kram -up 2b) when you want mesh for optical calculations to be finer, do: x kgen -so (to generate finer mesh) in third line in case.in2, change value of TETRA to be 101 x lapw1 -p -up x lapw1 -p -dn x lapwso -up -p x lapw2 -p -fermi -up -so x optic -so -up -p x joint -up -p x kram -up 3) When using w2k version 17.1, there is a bug in the function joint when electronic structure is spin-polarized case with so. In this case, all optical constant outgoing function joint have half values for w2k ver 17.1 compared to previous w2k versions. So either use w2k version 16.1 or smaller, or with w2k version 17.1, simply multiply all optical constants by factor 2. Hoping it helps Best regards Jaro ___ Wien mailing list Wien@zeus.theochem.tuwien.ac.at http://zeus.theochem.tuwien.ac.at/mailman/listinfo/wien SEARCH the MAILING-LIST at: http://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/index.html -- -- Mgr. Jaroslav Hamrle, Ph.D. Institute of Physics, room F232 Faculty of Mathematics and Physics Charles University Ke Karlovu 5 121 16 Prague Czech Republic tel: +420-95155 1340 email: ham...@karlov.mff.cuni.cz -- ___ Wien mailing list Wien@zeus.theochem.tuwien.ac.at http://zeus.theochem.tuwien.ac.at/mailman/listinfo/wien SEARCH the MAILING-LIST at: http://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/index.html
Re: [Wien] Optical properties with SO coupling
Dear Lokanath, I dont know why you have just two outgoing columns in the case.outputjointup. I always calculate full permittivity tensor from both case.inop and the case.injoint. And it always works fine for me. In the case of the case.inop, calculation of full optical tensor means calculation of 10 columns (9 optical constant + one energy), with case.injoint being as === start of an example case.inop 99 1 number of k-points, first k-point -5.0 3.0 Emin, Emax for matrix elements, NBvalMAX 9 number of choices (columns in *outmat): 2: hex or tetrag. case 1 Re xx 2 3 Re zz 4 5 6 7 8 9 OFF ON/OFF writes MME to unit 4 Choices: 1..Re 2..Re 3..Re 4..Re 5..Re 6..Re 7..Im 8..Im 9..Im = end of an example case.inop === Maybe problem is, that you calculated three optical constants in the case.inop, but asked for nine optical constants in the case.injoint. With my regards Jaroslav On 04/10/17 03:18, lokanath patra wrote: Dear Karel and Hamrle, Thank you for your reply. I have followed the steps already you have mentioned. I used switch 6 and obtained the plasma frequency. Then used it in case.inkram as my system is a metal. Then again I used switch 4 in case.injoint file and did the calculations before proceeding for x kram -up. My question is as follows: (1)When I proceed with x joint -up (with switch 4 in case.injoint) , the case.outputjointup file coming as follows # Energy [eV] Im_eps_xx Im_eps_zz 0.0 0.E+00 0.E+00 0.01361 0.E+00 0.E+00 0.02721 0.E+00 0.E+00 0.04082 0.33667974E+02 0.97123076E+01 0.05442 0.45512899E+02 0.14941513E+02 0.06803 0.56560294E+02 0.25979180E+02 0.08163 0.72797967E+02 0.51973076E+02 0.09524 0.14197996E+03 0.13384727E+03 0.10885 0.16008265E+03 0.16068757E+03 0.12245 0.17348056E+03 0.18197670E+03 0.13606 0.17206392E+03 0.18134116E+03 0.14966 0.16432457E+03 0.17616364E+03 0.16327 0.15422491E+03 0.16631451E+03 . Why there is no column for xy component when my case.symmatup has all the three components? How to calculate epsilon or absorptivity in xy direction also? Best Regards, Lokanath On Wed, Oct 4, 2017 at 4:00 AM, Jaroslav Hamrle <ham...@karlov.mff.cuni.cz <mailto:ham...@karlov.mff.cuni.cz>> wrote: Hallo, to calculate optical properties of Ni, after calculating electronic structure being spin-polarized and being with spin-orbit, do: 1) create both case.inop (your file looks correct) and case.injoint Example of case.injoint is: example of case.injoint === 1 : LOWER,UPPER and (optional) UPPER-VAL BANDINDEX 0. 0.00100 1. : EMIN DE EMAX FOR ENERGYGRID IN ryd eV : output units eV / ryd / cm-1 4 : SWITCH 9 : NUMBER OF COLUMNS 0.1 0.1 0.3 : BROADENING (FOR DRUDE MODEL - switch 6,7 - ONLY) SWITCH: 0...JOINTDOS FOR EACH BAND COMBINATION 1...JOINTDOS AS SUM OVER ALL BAND COMBINATIONS 2...DOS FOR EACH BAND 3...DOS AS SUM OVER ALL BANDS 4...Im(EPSILON) 5...Im(EPSILON) for each band combination 6...INTRABAND contributions 7...INTRABAND contributions including band analysis = end example case.injoint Now, you have to decide if you want to calculate optics at finer k-mesh than electronic structure, or the same mesh. In case electronic structure is calculated with k-mesh 30x30x30, it is good enough for Imxy and MOKE. 2a) when keeping the same k-mesh for optical calculations as for electronic structure, do: x lapw2 -p -fermi -up -so x optic -p -up -so (your command in your email is opticc, i.e. complex variant of command optic; opticc should be used when the structure lacks point symmetry, which Ni does not) x joint -up x kram -up 2b) when you want mesh for optical calculations to be finer, do: x kgen -so (to generate finer mesh) in third line in case.in2, change value of TETRA to be 101 x lapw1 -p -up x lapw1 -p -dn x lapwso -up -p x lapw2 -p -fermi -up -so x optic -so -up -p x joint -up -p x kram -up 3) When using w2k version 17.1, there is a bug in the function joint when electronic structure is spin-polarized case with so. In this case, all optical constant outgoing function joint have half values for w2k ver 17.1 compared to previous w2k versions. So
Re: [Wien] Optical properties with SO coupling
Hallo, to calculate optical properties of Ni, after calculating electronic structure being spin-polarized and being with spin-orbit, do: 1) create both case.inop (your file looks correct) and case.injoint Example of case.injoint is: example of case.injoint === 1 : LOWER,UPPER and (optional) UPPER-VAL BANDINDEX 0. 0.00100 1. : EMIN DE EMAX FOR ENERGYGRID IN ryd eV : output units eV / ryd / cm-1 4 : SWITCH 9 : NUMBER OF COLUMNS 0.1 0.1 0.3 : BROADENING (FOR DRUDE MODEL - switch 6,7 - ONLY) SWITCH: 0...JOINTDOS FOR EACH BAND COMBINATION 1...JOINTDOS AS SUM OVER ALL BAND COMBINATIONS 2...DOS FOR EACH BAND 3...DOS AS SUM OVER ALL BANDS 4...Im(EPSILON) 5...Im(EPSILON) for each band combination 6...INTRABAND contributions 7...INTRABAND contributions including band analysis = end example case.injoint Now, you have to decide if you want to calculate optics at finer k-mesh than electronic structure, or the same mesh. In case electronic structure is calculated with k-mesh 30x30x30, it is good enough for Imxy and MOKE. 2a) when keeping the same k-mesh for optical calculations as for electronic structure, do: x lapw2 -p -fermi -up -so x optic -p -up -so (your command in your email is opticc, i.e. complex variant of command optic; opticc should be used when the structure lacks point symmetry, which Ni does not) x joint -up x kram -up 2b) when you want mesh for optical calculations to be finer, do: x kgen -so (to generate finer mesh) in third line in case.in2, change value of TETRA to be 101 x lapw1 -p -up x lapw1 -p -dn x lapwso -up -p x lapw2 -p -fermi -up -so x optic -so -up -p x joint -up -p x kram -up 3) When using w2k version 17.1, there is a bug in the function joint when electronic structure is spin-polarized case with so. In this case, all optical constant outgoing function joint have half values for w2k ver 17.1 compared to previous w2k versions. So either use w2k version 16.1 or smaller, or with w2k version 17.1, simply multiply all optical constants by factor 2. Hoping it helps Best regards Jaro On 03/10/17 21:04, lokanath patra wrote: Dear Wien2K experts, I am trying to calculate the optical properties of Ni with SO coupling. I have prepared my case.inop file which looks like 9 1 number of k-points, first k-point -5.0 3.0 Emin, Emax for matrix elements, NBvalMAX 3 number of choices (columns in *outmat): 2: hex or tetrag. case 1 Re xx 3 Re zz 7 Im xy OFF ON/OFF writes MME to unit 4 Then I proceeded with x opticc -so -up and got the three components i.e. Re xx, Re zz and Im xy in case.symmatup file. From here, how should I proceed with the calculations so that I can get all these components in the case.absorpup file after x kram -up ? As I want to proceed with MOKE calculations, I want the xy components also. Please help me. Best regards, Lokanath -- Lokanath Patra Ph.D Scholar Dept. of Physics School of Applied and Basic Sciences Central University of Tamil Nadu Thiruvarur Tamil Nadu - 610101 Ph no - +91-8675834507 <http://www.avg.com/email-signature?utm_medium=email_source=link_campaign=sig-email_content=webmail> Virus-free. www.avg.com <http://www.avg.com/email-signature?utm_medium=email_source=link_campaign=sig-email_content=webmail> ___ Wien mailing list Wien@zeus.theochem.tuwien.ac.at http://zeus.theochem.tuwien.ac.at/mailman/listinfo/wien SEARCH the MAILING-LIST at: http://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/index.html -- -- Mgr. Jaroslav Hamrle, Ph.D. Institute of Physics, room F232 Faculty of Mathematics and Physics Charles University Ke Karlovu 5 121 16 Prague Czech Republic tel: +420-95155 1340 email: ham...@karlov.mff.cuni.cz -- ___ Wien mailing list Wien@zeus.theochem.tuwien.ac.at http://zeus.theochem.tuwien.ac.at/mailman/listinfo/wien SEARCH the MAILING-LIST at: http://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/index.html
Re: [Wien] Optical properties with SO coupling
Dear Gerhard, dear Gavin, as written in last email of Gavin, I think: In case of spin-polarized case without so, the optical transitions for up, down spins are written in case.jointup, case.jointdn, each one containing only optical transitions for corresponding up/down electrons. Program kram then sums up contributions from optical transitions from both electrons spins. In case of spin-polarized case with so, up and down spins are not described separately and written only to file 'up', like case.jointup. That is why in case with so coupling, all optical transitions should be contained with file case.jointup. That is why I think the current version of program joint should be corrected in the way that case.jointup should contain all optical transitions. However, now, output in case.jointup is somewhat artificially divided by too. Of course, it can be corrected in code kram, but I think it is not good idea. It should be corrected in the joint program. Then, x kram -up is just fine. With my regards Jaroslav On 05/10/17 15:14, Fecher, Gerhard wrote: Hi Jaroslav, if you check only case.jointup it has possibly only half the value because the other half is supposed to be in case.jointdn (with SO they should be the same) Did you try to copy case.jointup to case.jointdn (or run in addition everything for dn) and then addjoint then the factor 2 is included in the case.joint The problem with spinpolarisation and SO is that case.jointup is the only necessary and produced, however, kram expects that case.joint exists that's why one has to do the copy by hand (not rename, then the factor 2 will be missing, again !) (one might also run optic and joint both in addition for dn, before addjoint) Indeed, it would be nice if this behaviopur could be changed (maybe by some switch(es) to kram : e.g.: -so -up) Just for curiosity, I wonder whether and how crossterms are respected, the selcction rules on the total angular momentum j' = j, j+-1 together with those on the magnetic quantum number mj allow spin flip transitions even though the dipole operator does not act on the spin ! Ciao Gerhard DEEP THOUGHT in D. Adams; Hitchhikers Guide to the Galaxy: "I think the problem, to be quite honest with you, is that you have never actually known what the question is." Dr. Gerhard H. Fecher Institut of Inorganic and Analytical Chemistry Johannes Gutenberg - University 55099 Mainz and Max Planck Institute for Chemical Physics of Solids 01187 Dresden Von: Wien [wien-boun...@zeus.theochem.tuwien.ac.at] im Auftrag von Jaroslav Hamrle [ham...@karlov.mff.cuni.cz] Gesendet: Donnerstag, 5. Oktober 2017 09:57 An: wien@zeus.theochem.tuwien.ac.at Betreff: Re: [Wien] Optical properties with SO coupling Dear Gavin, I will describe my observation: I have calculated optical (epzz) and magneto-optical (K, for example K=epxy for M001) spectra of permittivity elements for bcc Fe. The electronic structure calculations are spin polarized, with spin-orbit, run by commands: runsp_lapw -p runsp_lapw -p -so -cc 0.01 -ec 0.001 -s lapw1 x lapw2 -p -fermi -up -so x optic -p -up -so x joint -p -up For w2k version 16.1, the calculated spectra corresponds to the experimental spectra (for both epzz and K). For w2k version 17.1, the calculated spectra are half-value for both epzz and K, compared to the experiment. Figures comparing spectra are here: http://alma.karlov.mff.cuni.cz/hamrle/w2kfig/Fe_Imepzz_compare.pdf http://alma.karlov.mff.cuni.cz/hamrle/w2kfig/Fe_ReK_compare.pdf In this example, I used permittivity spectra read directly from case.jointup files (I do not use output of kram). In the figures: - solid (noisy) line is output from case.jointup - the symbols are smeared spectra - black '+' are the experimental spectra - blue 'o' and red 'x' are spectra calculated by w2k version 17 - green '+' and yellow '*' are spectra calculated by w2k version 16 - y-axis denotes permittivity*E (in eV). That is why I have concluded that joint function in w2k version 17 has a bug in calculation of the optical permittivity. But I have not tested non-magnetic cases, I did it only for bcc Fe (sp+so). Hoping it helps. If I can help more, please let me know.. With my regards Jaro On 04/10/17 16:40, Gavin Abo wrote: Dear Jaro, I thought the spin-polarized SO optic normalization was broken in older versions of WIEN2k and was fixed in 17.1: https://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/msg16011.html Is it still broken? Kind Regards, Gavin On 10/3/2017 4:30 PM, Jaroslav Hamrle wrote: Hallo, to calculate optical properties of Ni, after calculating electronic structure being spin-polarized and being with spin-orbit, do: 1) create both case.inop (your file looks correct) and case.injoint Example of case.injoint is: example of case.injoint === 1
Re: [Wien] different MLD for bcc structure for magnetic equivalent directions M001, M010 and M100
else has thought", Albert Szent-Gyorgi www.numis.northwestern.edu <http://www.numis.northwestern.edu> ; Corrosion in 4D: MURI4D.numis.northwestern.edu <http://MURI4D.numis.northwestern.edu> Partner of the CFW 100% program for gender equity, www.cfw.org/100-percent <http://www.cfw.org/100-percent> Co-Editor, Acta Cryst A ___ Wien mailing list Wien@zeus.theochem.tuwien.ac.at http://zeus.theochem.tuwien.ac.at/mailman/listinfo/wien SEARCH the MAILING-LIST at: http://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/index.html -- -- Mgr. Jaroslav Hamrle, Ph.D. Institute of Physics, room F232 Faculty of Mathematics and Physics Charles University Ke Karlovu 5 121 16 Prague Czech Republic tel: +420-95155 1340 email: ham...@karlov.mff.cuni.cz -- ___ Wien mailing list Wien@zeus.theochem.tuwien.ac.at http://zeus.theochem.tuwien.ac.at/mailman/listinfo/wien SEARCH the MAILING-LIST at: http://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/index.html
Re: [Wien] different MLD for bcc structure for magnetic equivalent directions M001, M010 and M100
Dear Laurence, thank you for your detailed answer. I have tried all your suggestions, - I changed case.in0 with increased oversampling by factor two (new parameters LUSE 26 and IFFTfactor 4) start of case.in0 --- TOT XC_LDA (XC_PBE,XC_PBESOL,XC_WC,XC_MBJ,XC_REVTPSSS) NR2V IFFT 26 (R2V) 24 24 24 4.00 1 min IFFT-parameters, enhancement factor, iprint end of case.in0 --- - I also tried to impose strong convergence criteria to be -cc 0.0001 -ec 0.0001 However, in both cases, the ghost MLD remains practically identical as when using my default values (default Fe + convergence -cc 0.1 -ec 0.01) Also, final :PUP 'Current' parameters remained practically the same for all the calculations (by about 2 digits), being like (for M001) PW CHANGE H K L Current Change Residue :PUP001: 0 0 0 2.10719649E-02 5.090E-10 -7.530E-09 :PUP002: 0 -1 -1 3.67084665E-04 -7.004E-10 -3.572E-10 :PUP003: 1 -1 0 1.82124988E-04 -3.669E-10 -2.211E-10 :PUP004: 0 0 -2 -1.87471938E-03 6.192E-11 7.254E-10 :PUP005: 0 -2 0 -3.75090680E-03 1.125E-10 1.378E-09 :PUP006: 1 -1 -2 -3.46372731E-03 1.026E-10 1.378E-09 :PUP007: 1 -2 -1 -6.92804324E-03 2.418E-10 2.776E-09 :PUP008: 0 -2 -2 -7.14014083E-04 8.677E-11 2.032E-10 :PUP009: 2 -2 0 -3.57036516E-04 4.298E-11 1.121E-10 :PUP010: 0 -1 -3 2.62159615E-04 9.199E-11 -1.588E-10 :PUP011: 0 -3 -1 2.62289930E-04 9.844E-11 -1.555E-10 :PUP012: 1 -3 0 2.62219244E-04 9.133E-11 -1.551E-10 Unfortunately, all reflections seems to be allowed for bcc (H+K+L is even), forbidden reflections of bcc are (H+K+L=odd), so I can not see how they get close to zero ;-)) But the idea is excellent. Thank you again and with my best regards Jaroslav On 27/11/17 15:27, Laurence Marks wrote: Let me clarify slightly my comment about symmetry -- as I realized the explanation (I think) and can also suggest something that might help. First, concerning symmetry the explanation is I believe simple. If the problem has a real symmetry operation such as inversion which is being removed, then the Jacobian at the solution has zero's for charge disturbances that break this symmetry. Because of this noise due to numerical accuracy has a large effect, and almost certainly one has to tighten the convergence criteria particularly -cc. You can monitor this by looking at the :PUPXXX values in case.scfm and look how well the forbidden reflections have converged to zero. Second, do not be surprised about numerical issues. While the calculations are done in double precision, there are many large sums and in some cases double sums, and also numerical integrations/differentiation. Any large sum or numerical integration/differentiation in general reduces the numerical accuracy. Hence even though double precision has an accuracy of 1D-15 the sum may only be accurate to 1D-10 or even 1D-7. Also, the Intel ifort compiler will reduce the numerical accuracy for speed if one is not careful. One thing that may help is to increase the oversampling in case.in0 for VXC, both that of the PW's and of the CLMs. A standard test is to use LDA and see if the problem goes away, since oversampling is much less relevant for this. Of course your problem may have nothing to do with any of this ___ Wien mailing list Wien@zeus.theochem.tuwien.ac.at http://zeus.theochem.tuwien.ac.at/mailman/listinfo/wien SEARCH the MAILING-LIST at: http://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/index.html
Re: [Wien] different MLD for bcc structure for magnetic equivalent directions M001, M010 and M100
Hi Gerhard, I know that due to SO, the electronic structure calculated for 100, 010 and 001 magnetization directions are different. The problem I have is following: I have three calculated electronic structures of bcc Fe, with magnetizations along 001, 010 and 100. Then, for any cubic structure, the permittivity tensor elements (ep_ij) with the same relations with respect to the magnetization should be equal in all three calculated structures. For example, symmetry clearly states that diagonal permittivity elements parallel to magnetization direction must equal ep_xx (for M=100) = ep_yy (for M=010) = ep_zz (for M=001). My problem is, that for calculated bcc Fe they do not equal. More specifically, they do not equal solely for bcc structures in wien2k, with disagreement upto 1%. For simple cubic and fcc structures they do equal, with tiny disagreement upto 0.01% Any help how to overcome this would be very helpful Thank you and with my best regards Jaroslav On 25/11/17 14:13, Fecher, Gerhard wrote: Hi Jaroslav, with SO, 001 is not equivalent to 001 or 010, if the magnetisation is along 001 this you see easily from the changed symmetry after initializing SO (symmetso) regards from Dresden Ciao Gerhard DEEP THOUGHT in D. Adams; Hitchhikers Guide to the Galaxy: "I think the problem, to be quite honest with you, is that you have never actually known what the question is." Dr. Gerhard H. Fecher Institut of Inorganic and Analytical Chemistry Johannes Gutenberg - University 55099 Mainz and Max Planck Institute for Chemical Physics of Solids 01187 Dresden Von: Wien [wien-boun...@zeus.theochem.tuwien.ac.at] im Auftrag von Jaroslav Hamrle [ham...@karlov.mff.cuni.cz] Gesendet: Freitag, 24. November 2017 16:36 An: wien@zeus.theochem.tuwien.ac.at Betreff: [Wien] different MLD for bcc structure for magnetic equivalent directions M001, M010 and M100 Dear colleagues, We have found non-physical asymmetry related with equivalent magnetization directions, when calculating electronic structure for bcc Fe: We want to calculate magnetic linear dichroism, MLD, defined as a difference between diagonal permittivity element being parallel, perpendicular to direction of magnetization, respectively. MLD=epzz - (epxx+epyy)/2 for M001 MLD=epyy - (epxx+epzz)/2 for M010 MLD=epxx - (epyy+epzz)/2 for M100 Obviously, MLD calculated for different equivalent magnetization directions should be identical. But they are not, MLD calculated for 001 is different to MLD calculated for 010 and 100 (MLD for 010 and 100 are identical). In most cases, we used k-mesh 30x30x30, exgange LDA (choice 5), with convergence criteria runsp_lapw -so -cc 0.01 -ec 0.001 -s lapw1 and the convergence was reached. * We tested this asymmetry also for fcc structures (Ni, Co, Co2MnSi). We also tested simple cubic structure (bcc Fe, defined as a simple cubic structure with two Fe atoms). In all those cases, the asymmetry disappears. On the other hand, it also appeared also in bcc Ni. Hence, the asymmetry seems to be specifically related with bcc structure. * this asymmetry can be observed already in energy levels (files case.energysoup). Hence, we think, the asymmetry is not a feature of optics. Namely, there is a very good agreement for energies for M010 and M100 (in example below difference is below 2e-7Ry), but much bigger difference between energies for M001 and (M010,M100) (in example below max. difference is 18e-6 Ry for band 5). Therefore it seems that this problem arises in either lapw0 or lapw1 for bcc structure. To demonstrate the difference, we show energy levels for the first k-point (in vicinity of the Gamma point shifted in 111 direction from the Gamma point): Fe30M001: 0.E-01 0.E-01 0.E-01 155 18 8.0 1 -3.4390104377017581 2 -3.4064979309023942 3 -3.3508627657180750 4 -3.2276472567243979 5 -3.1955089683446780 6 -3.1702455400854954 7 -7.1658179115217727E-002 8 -4.3723732810772589E-002 9 0.37296762299903474 10 0.37521967189559313 Fe30M010: 0.E-01 0.E-01 0.E-01 155 18 8.0 1 -3.4390110394480322 2 -3.4064968725403300 3 -3.3508644682352022 4 -3.2276486274720977 5 -3.1954902103327028 6 -3.1702472318057655 7 -7.1659013996950252E-002 8 -4.3723316415832839E-002 9 0.37296632778787425 10 0.37521816821120640 Fe30M100: 0.E-01 0.E-01 0.E-01 155 18 8.0 1 -3.4390109925234049 2 -3.4064968346346225 3 -3.3508643700301919 4 -3.2276485335559135 5 -3.195490
[Wien] problem with local rotation matrices in bcc Fe
ferent magnetization directions, where local rotation matrix rotates magnetization to local x/y direction compared to those rotating to local z direction? Numerically, the results are very similar, but for example in MLD, we can see some difference, but it can be due to numerical issues. For example, LM for 110 magnetization is lm: 0 0 2 0 2 2 4 0 4 2 4 4 6 0 6 2 6 4 6 6 which seems to me to be build for preferred z-direction. On the other hand, band energies equals within 1e-6 Ry (see below) between [110] and [011], having magnetization in the local coordinates in x,z directions, respectively. 3) do you think/suggest, is it a good practise to replace local rotation matrix by a matrix, which always rotates magnetization to [0 0 1] direction in the local coordinates? I understand there is also structural symmetry, but I'm not sure how limiting it is for simple cases such as bcc Fe. Thank you for any help and with our best regards Jaroslav Hamrle Ondrej Stejskal === magnetization [1 1 1] == --> from case.outsymso: magnetization direction: [1 1 1] theta , phi [deg] of the magnetization direction: 54.7356 45. local rotation matrix R: 0.4082 -0.7071 0.5774 0.4082 0.7071 0.5774 -0.8165 0 0.5774 --> R rotates magnetization in local coordinates to: [0 0 1] ---> calculated using [0 0 1]*inv(R) (should be equal to original magnetization direction): direction [0 0 1]* inv(R): [1 1 1] /sqrt(3) theta, phi [deg] of direction [0 0 1]* inv(R): 54.7356 45. local rotational matrix inv(Rz(phi)*Ry(theta)) 0.4082 -0.7071 0.5774 0.4082 0.7071 0.5774 -0.8165 0 0.5774 === magnetization [1 1 -1] == --> from case.outsymso: magnetization direction: [1 1 -1] theta , phi [deg] of the magnetization direction: 125.2644 45. local rotation matrix R: 0.4082 -0.7071 -0.5774 0.4082 0.7071 -0.5774 0.8165 0 0.5774 --> R rotates magnetization in local coordinates to: [0 0 -1] ---> calculated using [0 0 1]*inv(R) (should be equal to original magnetization direction): direction [0 0 1]* inv(R): [-1 -1 1] /sqrt(3) theta, phi [deg] of direction [0 0 1]* inv(R): 54.7356 -135. local rotational matrix inv(Rz(phi)*Ry(theta)) -0.4082 -0.7071 0.5774 -0.4082 0.7071 0.5774 -0.8165 0. -0.5774 === magnetization [1 -1 1] == --> from case.outsymso: magnetization direction: [1 -1 1] theta , phi [deg] of the magnetization direction: 54.7356 -45. local rotation matrix R: 0.4082 0.7071 0.5774 -0.4082 0.7071 -0.5774 -0.8165 0 0.5774 --> R rotates magnetization in local coordinates to: [0 0 1] ---> calculated using [0 0 1]*inv(R) (should be equal to original magnetization direction): direction [0 0 1]* inv(R): [1 -1 1] /sqrt(3) theta, phi [deg] of direction [0 0 1]* inv(R): 54.7356 -45. local rotational matrix inv(Rz(phi)*Ry(theta)) 0.4082 0.7071 0.5774 -0.4082 0.7071 -0.5774 -0.8165 0 0.5774 === magnetization [-1 1 1] == --> from case.outsymso: magnetization direction: [-1 1 1] theta , phi [deg] of the magnetization direction: 54.7356 135. local rotation matrix R: 0.4082 0.7071 0.5774 -0.4082 0.7071 -0.5774 -0.8165 0 0.5774 --> R rotates magnetization in local coordinates to: [-0.94281 0 -0.3] ---> calculated using [0 0 1]*inv(R) (should be equal to original magnetization direction): direction [0 0 1]* inv(R): [1 -1 1] /sqrt(3) theta, phi [deg] of direction [0 0 1]* inv(R): 54.7356 -45. local rotational matrix inv(Rz(phi)*Ry(theta)) -0.4082 -0.7071 -0.5774 0.4082 -0.7071 0.5774 -0.8165 0 0.5774 === magnetization [0 0 1] == --> from case.outsymso: magnetization direction: [0 0 1] theta , phi [deg] of the magnetization direction: 0 0 local rotation matrix R: 1 0 0 0 1 0 0 0 1 --> R rotates magnetization in local coordinates to: [0 0 1] ---> calculated using [0 0 1]*inv(R) (should be equal to original magnetization direction): direction [0 0 1]* inv(R): [0 0 1] theta, phi [deg] of direction [0 0 1]* inv(R): 0 0 local rotational matrix inv(Rz(phi)*Ry(theta)) 1 0 0 0 1 0 0 0 1 === magnetization [0 0 -1] == --> from case.outsymso: magnetization direction: [0 0 -1] theta , phi [deg] of the magnetization direction: 180 0 local rotation matrix R: 1 0 0 0 1 0 0 0 1 --> R rotates magnetization
Re: [Wien] Magnetocrystalline anisotropy
Dear Xavier, your problem somewhat resembles me my problem I had when calculating magnetic linear dochroism (MLD) on bcc Fe. The similarity is that we both want to see small changes in electronic structure when rotating magnetic field direction. What help me: 1) run fine convergence criteria, such as runsp_lapw -p -cc 0.01 -ec 0.01 2) as suggested Prof. Blaha, it was important to increase k-mesh (in my case up to 100), and apply fine BZ integration (TEMP or TEMPS) with small value as 0.001, not default TETRA. for example change case.in2 by using command sed '3s/^/TEMP 0.001 /' $file.in2 > $file.in2_TEMPnew more here" https://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/msg16815.html https://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/msg16844.html 3) for some magnetization direction, I had problem with either wrong, either suspicious values of local rotation matrix, https://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/msg16894.html Problem was that for some external magnetization directions, the local direction of magnetization was not [001]. In one case (external M along [-111]), the local magnetization direction was [-0.94281 0 -0.3] which I think is not correct, and MLD was wrong too. In some case, local magnetization direction was along x or along y, which I dont know if it is correct. On one hand, eigenenergies agreed perfectly, but anyway I saw small change in MLD in those cases. But as a blind suggestion for you, try if local rotation matrices are correct. Namely try if mag_glob*R = mag_loc where mag_glob is your (external i.e. in global coordinates) magnetization direction, R is local rotation matrix for each atom (can be found in case.struct or case.outsymso) and mag_loc is local magnetization direction, which in my (maybe naive and wrong) understanding should be [001]. Hoping it helps With my best regards Jaroslav On 09/01/18 09:44, Xavier Rocquefelte wrote: Dear Colleagues I recently obtained a surprising result concerning the calculation of the magnetocrystalline anisotropy energy (MAE) of SeCuO3. This compound has a monoclinic symmetry (SG. P21/n) and is known to be antiferromagnetically ordered at low temperature. Here I provide the results obtained for two magnetic orders, named FM and AFM1 (see attached document) : https://filesender.renater.fr/?s=download=1da93a22-9592-3a7e-ba2e-1533fcae45d2 These calculations have been done using WIEN2k_17, GGA = PBE, RKMAX = 6, kmesh = 5 4 4 and in P1 symmetry. The results are the same using RKMAX = 7. The AFM1 order is the more stable one, as expected. However, as shown in the document the MAE of AFM1 order is not symmetric, which is not expected. In contrast the MAE for FM order is symmetric. Based on the recent discussion "zigzag potential", it seems to me that the AFM1 MAE should be symmetric, because the magnetic moment is a pseudo-vector. Is it possible that the present problem is related to the fact that in the present implementation of the spin-orbit coupling we neglect the off-diagonal terms? Do you have any idea about the problem we are facing? Does someone observe such unusual MAE for other systems? Best Regards Xavier ___ Wien mailing list Wien@zeus.theochem.tuwien.ac.at http://zeus.theochem.tuwien.ac.at/mailman/listinfo/wien SEARCH the MAILING-LIST at: http://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/index.html -- ------ Mgr. Jaroslav Hamrle, Ph.D. Institute of Physics, room F232 Faculty of Mathematics and Physics Charles University Ke Karlovu 5 121 16 Prague Czech Republic tel: +420-95155 1340 email: ham...@karlov.mff.cuni.cz -- ___ Wien mailing list Wien@zeus.theochem.tuwien.ac.at http://zeus.theochem.tuwien.ac.at/mailman/listinfo/wien SEARCH the MAILING-LIST at: http://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/index.html
[Wien] wrong sign of the off-diagonal magneteoptic permittivity spectra (file case.jointup)
Dear collegues, based on our measurements of MOKE spectra, we think that the magnetooptic off-diagonal permittivity in the file case.jointup has opposite sign, as compared to the experiment. For example, in the case of [001] magnetization, column Re_eps_xy in the file case.jointup has opposite sign, as compared to the experimental spectra. Namely, we made careful sign tests on bcc Fe. Then, according to our experiment, Re_eps_xy should be negative for energy about 2 eV, but the permittivity calculated by wien2k, Re_eps_xy, has positive sign. This statement is further supported for example by spectra in article Oppeneer et al, PRB45 10924 (1992), https://journals.aps.org/prb/abstract/10.1103/PhysRevB.45.10924 figure 3b, where after conversion of the off-diagonal conductivity to the permitivity (ep=i*sigma/(omega*ep0)) it provides negative real part of the off-diagonal permittivity at about 2eV as well. Furthermore, we have written our own integration over BZ of the transition matrix elements, and transition matrix elements (files case.symmatup, case.mommat2up) seems to us to be sign-correct. We think, wrong sign in Re_eps_xy simply appears in conversion from conductivity (calculated by Kubo formula, OPT1.f) to the outgoing permittivity. We think this step appears on lines 965 and 985 of file joint.f, w2k version 18.2: VV=EPSF*DENSTY(IB,J,icol)/ENG2(j) where division by the energy makes basically conversion from the conductivity to the permittivity. However, as the conversion itself handles complex numbers: eps = i/omega * sig the conversion between real and imaginary parts have opposite signs eps_Im = sigma_Re/omega (this is correct in code) eps_Re= - sigma_Im/omega (this is the expression for the off-diagonal permittivity induced by magnetization, and in this case, we think minus sign is missing in the code, namely for Re_eps_xy, Re_eps_xz and Re_eps_yz. ) Probably, those arguments will apply also for X-ray off-diagonal spectra elements (line 906 in joint.f). Hoping my comment will be useful to improve w2k code. With my best regards Jaroslav Hamrle & Ondrej Stejskal -- -- Mgr. Jaroslav Hamrle, Ph.D. Institute of Physics, room F232 Faculty of Mathematics and Physics Charles University Ke Karlovu 5 121 16 Prague Czech Republic tel: +420-95155 1340 email: ham...@karlov.mff.cuni.cz -- ___ Wien mailing list Wien@zeus.theochem.tuwien.ac.at http://zeus.theochem.tuwien.ac.at/mailman/listinfo/wien SEARCH the MAILING-LIST at: http://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/index.html
