Overwriting default variables in Input File?

I’ve been working with Abinit’s ujdet function and utility and have run into a problem the documentation does not resolve. I need Abinit’s U determination function to be activated, thus setting macro_uj > 0. But the available macro_uj options, (i.e., macro_uj = 1,2 or 3) the input atvshift is overwritten. I need to preserve the input atvshift and determine U simultaneously.

The variable information for macro_uj says that it "Sets proper input values for the determination of U and J i.e. for pawujat (first atom treated with PAW+U), irdwfk (=1), tolvrs (=10^(-8)), nstep (=255), diemix (=0.45), atvshift (pawujat) pawujv. Do not overwrite these variables manually unless you know what you are doing.” This implies that there is some way to overwrite these variables manually. How do I do this? I assure you, I know what I’m doing.

Hi Lmacenul,

your question is very specific - I have written to @amadon @torrent and their previous student Donat Adams who did the implementation, to see who can help you. Short of hacking the code to remove the “overwrite” I don’t know what else to suggest, but there may be collateral damage.

Best

Matthieu

1 Like

Hi Lmacenul,

If you really “know what you are doing”, the sentence in the manual suggests that you should hack the pawuj_ini function in the src/65_paw/m_paw_uj.F90 file.
But, as Matthieu wrote, there may be collateral damage… and also the result might be unphysical.

This implementation is rather old and it difficult to remember all the choices for the variables.
(perhaps @amadon…)

Regards,
Marc

 do iuj=0,ndtset
   !write(std_out,*)'pawuj_ini iuj ',iuj
   dtpawuj(iuj)%diemix=0.45_dp
   dtpawuj(iuj)%iuj=0
   dtpawuj(iuj)%nat=0
   dtpawuj(iuj)%ndtset=1
   dtpawuj(iuj)%nspden=1
   dtpawuj(iuj)%macro_uj=0
   dtpawuj(iuj)%option=1
   dtpawuj(iuj)%pawujat=1
   dtpawuj(iuj)%pawujga=one
   dtpawuj(iuj)%pawprtvol=1
   dtpawuj(iuj)%ph0phiint=one
   dtpawuj(iuj)%dmatpuopt=3
   dtpawuj(iuj)%pawujrad=3.0_dp
   dtpawuj(iuj)%pawrad=20.0_dp
   !Allocate arrays
   !write(std_out,*)'pawuj_ini before arrays'
   ABI_MALLOC(dtpawuj(iuj)%rprimd,(3,3))
   ABI_MALLOC(dtpawuj(iuj)%scdim,(3))
   ABI_MALLOC(dtpawuj(iuj)%wfchr,(nwfchr))
   dtpawuj(iuj)%rprimd=reshape((/ 1,0,0,0,1,0,0,0,1/),(/ 3,3 /))
   dtpawuj(iuj)%scdim=reshape((/ 250,0,0 /),(/3 /))
   dtpawuj(iuj)%wfchr=(/ (0,im1=1,nwfchr) /)
   if (iuj>0) then
     dtpawuj(iuj)%iuj=-1
   end if
 end do
1 Like

Thanks, Mattieu! I appreciate the message forward.

1 Like

Hello @torrent. Thanks so much for your response.

I ended up having to hack the code, although I avoided intervening with atvshift to avoid collateral as best as I could.

Parsing through the source code, I found that macro_uj is only ever subject to constraints if it is a) greater than zero or b) not equal to zero. To acquiesce to these, I encoded another option for macro_uj = 4, which is both greater than zero and not equal to zero and, for all intents and purposes, should trigger all constraints mutually triggered by the other three macro_uj options. Then I added the following lines:

To 44_abitypes_defs/m_dtset.F90

if (dtsets(idtset)%macro_uj>0) then
     ! Level shift atom with amplitude pawujv
     dtsets(idtset)%atvshift(:,:,pawujat)=dtsets(idtset)%pawujv
     ! Case level shift only on one spin channel
     if ((dtsets(idtset)%macro_uj==2.or.dtsets(idtset)%macro_uj==3).and.dtsets(idtset)%nsppol==2) then
       dtsets(idtset)%atvshift(:,2,pawujat)=0_dp
     end if

!# CHANGE TO CODE LM 2021
     if (dtsets(idtset)%macro_uj==4.and.dtsets(idtset)%nsppol==2) then
       dtsets(idtset)%atvshift(:,2,pawujat)= -dtsets(idtset)%pawujv
     end if
!# END CHANGE LM 2021

   end if ! macro_uj

To 65_paw/m_paw_uj.F90

!Calculation of response matrix

 do jdtset=1,4
   if (nspden==1) then
     chih(jdtset,1:nat_org)=dtpawuj(jdtset)%occ(1,:)
   else if (macro_uj==1.and.nspden==2) then
     chih(jdtset,1:nat_org)=dtpawuj(jdtset)%occ(1,:)+dtpawuj(jdtset)%occ(2,:)
   else if (macro_uj==2.and.nspden==2) then
     chih(jdtset,1:nat_org)=dtpawuj(jdtset)%occ(1,:)
   else if (macro_uj==3.and.nspden==2) then
     chih(jdtset,1:nat_org)=dtpawuj(jdtset)%occ(2,:)
!# CHANGE TO CODE LM 2021
   else if (macro_uj==4.and.nspden==2) then
     chih(jdtset,1:nat_org)=dtpawuj(jdtset)%occ(1,:)+dtpawuj(jdtset)%occ(2,:)
!# END CHANGE LM 2021
   end if
   vsh(jdtset)=dtpawuj(jdtset)%vsh(1,pawujat)
   if (pawprtvol==3) then
     write(message,fmt='(2a,i3,a,f15.12)') ch10,' Potential shift vsh(',jdtset,') =',vsh(jdtset)
     call wrtout(std_out,message,'COLL')
     write(message,fmt='( a,i3,a,120f15.9)') ' Occupations occ(',jdtset,') ',chih(jdtset,1:nat_org)
     call wrtout(std_out,message,'COLL')
   end if
 end do

Selecting macro_uj = 4 now makes atvshift for the second spin channel equal to negative pawujv, keeping atvshift for the first spin channel equal to positive pawujv. From what I could find, there should be no consequences for having a negative atvshift value.

This seemed to do the trick. Are there any other considerations I’ve missed, though? I want to make sure the answers are still physical.

Again, thank you for your response!