Issue reading a2f file written from the code in superconductivity calculation

[jarinfrench]

Dear EPW Community,

I am encountering a list-directed I/O syntax error while running the EPW module of Quantum Espresso. Here are the details of the issue:

Error Message:

Code: Select all

forrtl: severe (59): list-directed I/O syntax error, unit 65, file .../zrn.a2f
Image              PC                Routine            Line        Source
epw.x              0000000000D9076E  for_read_seq_lis_     Unknown  Unknown
epw.x              0000000000D8FC4C  for_read_seq_lis      Unknown  Unknown
epw.x              00000000006D3E89  Unknown               Unknown  Unknown
epw.x              00000000005FF195  Unknown               Unknown  Unknown
epw.x              000000000041A549  Unknown               Unknown  Unknown
epw.x              000000000041973D  Unknown               Unknown  Unknown
libc-2.28.so       0000149C79C487E5  __libc_start_main     Unknown  Unknown
epw.x              000000000041965E  Unknown               Unknown  Unknown

Input File:

Code: Select all

&inputepw
  prefix      = 'zrn',
  amass(1)    = 91.2240,
  amass(2)    = 14.0067,
  outdir      = './',
  dvscf_dir   = './save',

  elph        = .true.,
  ep_coupling = .true.,
  epwwrite    = .true.,
  epwread     = .false.,

  vme         = 'wannier',
  scdm_proj   = .true.,
  nbndsub     =  16,
  bands_skipped = 'exclude_bands = 1-5',

  wannierize  = .true.,
  num_iter    = 300,
  !dis_win_max = 21,
  !nosym = .true.,
  !dis_froz_min= 11.1,
  !dis_froz_max= 19.8,
  !proj(1)     = 'Th:d',
  !proj(2)     =  'N:p',

  auto_projections = .true.,
  wdata(1) = 'bands_plot = .true.',
  wdata(2) = 'begin kpoint_path',
  wdata(3) = 'G 0.00 0.00 0.00 X 0.00 0.50 0.50',
  wdata(4) = 'X 0.00 0.50 0.50 W 0.25 0.50 0.75',
  wdata(5) = 'W 0.25 0.50 0.75 L 0.50 0.50 0.50',
  wdata(6) = 'L 0.50 0.50 0.50 K 0.375 0.375 0.75',
  wdata(7) = 'K 0.375 0.375 0.75 G 0.00 0.00 0.00',
  wdata(8) = 'G 0.00 0.00 0.00 L 0.50 0.50 0.50',
  wdata(9) = 'end kpoint_path',
  wdata(10) = 'bands_plot_format = gnuplot',

  fsthick  = 2,
  degaussw = 0.1,
  degaussq = 0.05,

  ephwrite = .true.,
  eliashberg = .true.,

  liso       = .true.,
  limag      = .true.,
  lpade      = .true.,
  lacon      = .true.,

  nsiter     = 500,
  npade      = 40,
  broyden_beta = -0.8,
  conv_thr_iaxis = 1.0d-3,
  conv_thr_racon = 1.0d-3,

  wscut  = 0.1,
  muc = 0.1,
  nstemp = 21,
  temps = 0.25, 10,

  mp_mesh_k = .true.,

  nkf1         = 18,
  nkf2         = 18,
  nkf3         = 18,

  nqf1         = 18,
  nqf2         = 18,
  nqf3         = 18,

  nk1         = 8,
  nk2         = 8,
  nk3         = 8,

  nq1         = 8,
  nq2         = 8,
  nq3         = 8,
/

The zrn.a2f file is written by the EPW code and appears to be formatted correctly:

Code: Select all

w[meV] a2f and integrated 2*a2f/w for   10 smearing values
   0.1421312   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000
   0.2842625   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000
   0.4263937   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000
...
  70.7813582   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000 125.9278102 166.8877107 167.9911906 168.0345032 168.0873925 168.1522559 168.2292301 168.3184796 168.4201989 168.5346146
  70.9234894   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000 125.9278102 166.8877107 167.9911906 168.0345032 168.0873925 168.1522559 168.2292301 168.3184796 168.4201989 168.5346146
  71.0656207   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000   0.0000000 125.9278102 166.8877107 167.9911906 168.0345032 168.0873925 168.1522559 168.2292301 168.3184796 168.4201989 168.5346146
 Integrated el-ph coupling
  #          125.9278102 166.8877107 167.9911906 168.0345032 168.0873925 168.1522559 168.2292301 168.3184796 168.4201989 168.5346146
 Phonon smearing (meV)
  #            0.0500000   0.1000000   0.1500000   0.2000000   0.2500000   0.3000000   0.3500000   0.4000000   0.4500000   0.5000000
Electron smearing (eV)   0.1000000
Fermi window (eV)   0.4000000
Summed el-ph coupling  167.9408767

Description: The error occurs when the EPW module attempts to read the zrn.a2f file. I have checked the format of the input file and the zrn.a2f file, and they appear to be correctly structured. The input file is also formatted correctly as far as I am able to tell.

System Information:

Quantum Espresso Version: 7.2
Compiler: Intel OneAPI 2023.2.1
Operating System: Rocky Linux 8.9
I would greatly appreciate any guidance or suggestions on how to resolve this issue. Thank you for your assistance.

Best regards,

Jarin French
Idaho National Laboratory

[hmori]

Dear Jarin,

I'm not sure of the exact reason why you received the error message when reading the a2f file. However, it appears that your zrn.a2f file may not have been calculated using the specified input parameters. The zrn.a2f file indicates a Fermi window of 0.4 eV, while you specified fsthick as 2 eV. Could you please double-check that both your input file and zrn.a2f are located in the same directory?

If they are indeed in the same directory, try moving the zrn.a2f file out of the working directory and rerun the calculation with changing the following input flags so that you do not have to start the calculation from scratch.

Code: Select all

  elph        = .false.,
  ep_coupling = .false.,
  epwwrite    = .false.,
  epwread     = .true.,
  
  wannierize  = .false.,
  ephwrite = .false.,
 

The EPW should recalculate the a2f if the a2f file is not found.

Best regards,
Hitoshi

[hpaudya1]

Hi,

In addition to the Hitoshi's comment, please consider checking the value of el-ph coupling strength. It looks unphysical.

Happy EPWing,
Hari

[Shashi]

Hi Jarin,

Can you also try without using auto_projections? Check if the a2f is correctly written without printing **** in the file. Also, verify the number of lines in the a2f file. The default nqstep = 500 and if you change nqstep while reading the a2f file it will crash.

Regards,
Shashi

[jarinfrench]

Thanks everyone for your input! I started my calculations over from scratch in a new directory, and found that the issue with the incorrectly generated a2f file was due to a previous EPW run that caused some mismatch (as suggested by Hitoshi). Hari's suggestion to check the unphysical values of the el-ph coupling strength also helped in determining the issue.

To Shashi, I'm still not very familiar with the software and am unsure how I would try these calculations without using auto_projections.

Now that the file is correctly generated, there are no problems with reading the file. However, the calculation still fails during calculation of Pade approximant:

Code: Select all

...
     Finish reading a2f file

     Electron-phonon coupling strength =    0.6736182

     Estimated Allen-Dynes Tc =    11.575226 K for muc =    0.10000

     Estimated w_log in Allen-Dynes Tc =    31.576969 meV

     Estimated BCS superconducting gap =     1.755557 meV

     Estimated Tc from machine learning model =    11.507883 K


     temp(  1) =      0.25000 K

     Solve isotropic Eliashberg equations on imaginary-axis

     Total number of frequency points nsiw(     1) =    739
     Cutoff frequency wscut =     0.1001
     linear mixing factor =      0.80000
     mixing factor = 0.2 is used for the first three iterations.

     Actual number of frequency points (     1) =    739 for uniform sampling
        iter      ethr          znormi      deltai [meV]
          1   8.441228E-01   1.580458E+00   2.105738E+00
          2   2.090802E-02   1.580176E+00   2.126163E+00
          3   1.745944E-02   1.579925E+00   2.143931E+00
          4   1.455885E-02   1.579704E+00   2.159500E+00
          5   4.943736E-02   1.578913E+00   2.212644E+00
          6   1.657813E-02   1.578520E+00   2.243430E+00
          7   7.874799E-03   1.578297E+00   2.263692E+00
          8   5.677449E-03   1.578152E+00   2.277981E+00
          9   4.545666E-03   1.578052E+00   2.288391E+00
         10   3.530687E-03   1.577979E+00   2.296079E+00
         11   2.677884E-03   1.577925E+00   2.301790E+00
         12   2.009068E-03   1.577885E+00   2.306040E+00
         13   1.500255E-03   1.577855E+00   2.309204E+00
         14   1.118051E-03   1.577833E+00   2.311561E+00
         15   8.324800E-04   1.577817E+00   2.313315E+00
     Convergence was reached in nsiter =     15

     Temp (itemp =   1) =    0.250 K  Free energy =    -0.002049 meV

     iaxis_imag   :      0.01s CPU      0.85s WALL (       1 calls)


     Pade approximant of isotropic Eliashberg equations from imaginary-axis to real-axis
     Cutoff frequency wscut =     0.1000

        pade    Re[znorm]   Re[delta] [meV]
        296            NaN            NaN
     Convergence was reached for N =    296 Pade approximants

     raxis_pade   :      0.00s CPU      1.07s WALL (       1 calls)

This leads to nsiter iterations for the analytic continuation of isotropic Eliashberg equations from imaginary-axis to real-axis that consist entirely of NaNs. Where should I be looking for a way to address this?

As a side note, this NaN error occurs regardless of linear mixing or Broyden mixing.

[hmori]

Dear Jarin,

The Padé approximation is a relatively robust method, but it can fail to calculate Padé coefficients when using many Matsubara frequencies. This is particularly true because the gap function at higher Matsubara frequencies often contains numerical noise, which can cause failures. Given that there are 739 Matsubara frequencies in this case, such a failure is understandable. You can adjust the number of Matsubara frequencies by changing wcut, but it's important to note that the number of Matsubara frequencies required for accurately calculating the gap function differs from that needed for the Padé approximation, so adjusting wcut for the Padé approximation is not advisable. Here, you can use the input parameter npade. The default value of npade is 90, meaning that 90% of the smallest Matsubara frequencies are used for the Padé approximation. Since you set npade to 40, it appears that you used 296 Matsubara frequencies for the Padé approximation. You could try reducing npade further, but applying the Padé approximation to the gap function at 0.25 K is probably challenging. Even if the Padé approximation cannot be applied, you can still compute the gap function and transition temperature, so it might be best to forgo the Padé approximation at 0.25 K. The Allen-Dynes prediction for Tc is about 11 K, so the gap function can still be obtained at slightly higher temperatures. Since fewer Matsubara frequencies are involved in the ME equation at higher temperatures, the Padé approximation is expected to work out. A good approach would be to set lpade to false, obtain the temperature dependence of the gap function, and then apply the Padé approximation only at some temperatures of interest.

Best regards,
Hitoshi