The DFT solver, as written

si_it_check.viv

Generated by RUNNING the independent numpy reference at the deck's settings and pinning what it printed; then the published LDA values; then the figure.

View on GitHub

# Si: what the calculation must reproduce, and the picture of it.

#

# dft/run.sh si it

#

# GENERATED by dft/ref/make_deck.py, which RAN dft/ref/pwdft_ref.py -- an independent numpy implementation that diagonalizes H densely --

# at exactly this deck's settings and pinned what it printed.

# A disagreement here is a bug in one of the two implementations.

GAMMA_P = NSEG

LAST = NPTS - 1.0

test "SCF converged: density, energy, and every occupied band's residual" {

assert_lt(DRHO, 1.0e-8)

assert_lt(EIGRES, EIG_TOL)

assert_lt(ITER, MAXIT)

}

test "every band-path point converged" {

assert_lt(BANDRES, EIG_TOL)

}

test "total energy agrees with the numpy reference" {

assert_near(ETOT, 0.480080446924, 1.0e-7)

}

test "band energies at L, Gamma and X agree with the numpy reference (eV)" {

assert_near(EVB[0.0, 0.0], -9.643898540, 1.0e-5) # L, band 1

assert_near(EVB[0.0, 1.0], -7.020423871, 1.0e-5) # L, band 2

assert_near(EVB[0.0, 2.0], -1.209812691, 1.0e-5) # L, band 3

assert_near(EVB[0.0, 3.0], -1.209812691, 1.0e-5) # L, band 4

assert_near(EVB[0.0, 4.0], 1.425065296, 1.0e-5) # L, band 5

assert_near(EVB[0.0, 5.0], 3.326441180, 1.0e-5) # L, band 6

assert_near(EVB[GAMMA_P, 0.0], -11.983838003, 1.0e-5) # Gamma, band 1

assert_near(EVB[GAMMA_P, 1.0], -0.000000000, 1.0e-5) # Gamma, band 2

assert_near(EVB[GAMMA_P, 2.0], -0.000000000, 1.0e-5) # Gamma, band 3

assert_near(EVB[GAMMA_P, 3.0], 0.000000000, 1.0e-5) # Gamma, band 4

assert_near(EVB[GAMMA_P, 4.0], 2.536605820, 1.0e-5) # Gamma, band 5

assert_near(EVB[GAMMA_P, 5.0], 2.536605820, 1.0e-5) # Gamma, band 6

assert_near(EVB[LAST, 0.0], -7.837469277, 1.0e-5) # X, band 1

assert_near(EVB[LAST, 1.0], -7.837469277, 1.0e-5) # X, band 2

assert_near(EVB[LAST, 2.0], -2.877092576, 1.0e-5) # X, band 3

assert_near(EVB[LAST, 3.0], -2.877092576, 1.0e-5) # X, band 4

assert_near(EVB[LAST, 4.0], 0.614392968, 1.0e-5) # X, band 5

assert_near(EVB[LAST, 5.0], 0.614392968, 1.0e-5) # X, band 6

}

#: REFPATH the reference's lowest 6 bands at all 61 path points, eV

REFPATH = [[-9.6438985, -7.0204239, -1.2098127, -1.2098127, 1.4250653, 3.3264412],

[-9.6600675, -6.9973610, -1.2078028, -1.2078028, 1.4273643, 3.3287600],

[-9.7068894, -6.9298644, -1.2017740, -1.2017740, 1.4342582, 3.3356428],

[-9.7799525, -6.8223674, -1.1917287, -1.1917287, 1.4457372, 3.3468662],

[-9.8732350, -6.6809069, -1.1772720, -1.1772720, 1.4526022, 3.3622775],

[-9.9818803, -6.5104914, -1.1592475, -1.1592475, 1.4730729, 3.3809151],

[-10.1004063, -6.3165922, -1.1372310, -1.1372310, 1.4980614, 3.4022360],

[-10.2225011, -6.1029051, -1.1106726, -1.1106726, 1.5276298, 3.4257570],

[-10.3498986, -5.8732232, -1.0807301, -1.0807301, 1.5615012, 3.4493973],

[-10.4778442, -5.6297070, -1.0468818, -1.0468818, 1.5997030, 3.4722168],

[-10.6040379, -5.3737204, -1.0088354, -1.0088354, 1.6457143, 3.4938686],

[-10.7282594, -5.1078286, -0.9673875, -0.9673875, 1.6921885, 3.5097651],

[-10.8489609, -4.8327145, -0.9222543, -0.9222543, 1.7425957, 3.5192401],

[-10.9653662, -4.5494072, -0.8735617, -0.8735617, 1.7967324, 3.5201747],

[-11.0762311, -4.2573999, -0.8225891, -0.8225891, 1.8544715, 3.5076674],

[-11.1823779, -3.9605514, -0.7673193, -0.7673193, 1.9152390, 3.4860873],

[-11.2827651, -3.6580595, -0.7091199, -0.7091199, 1.9787776, 3.4512938],

[-11.3774022, -3.3499481, -0.6485958, -0.6485958, 2.0515500, 3.4021610],

[-11.4659113, -3.0392114, -0.5856112, -0.5856112, 2.1187268, 3.3411511],

[-11.5474428, -2.7250739, -0.5209525, -0.5209525, 2.1869208, 3.2688854],

[-11.6205840, -2.4089919, -0.4545757, -0.4545757, 2.2559731, 3.1887562],

[-11.6900457, -2.0931063, -0.3896333, -0.3896333, 2.3225858, 3.1011558],

[-11.7512136, -1.7786907, -0.3243167, -0.3243167, 2.3871862, 3.0104421],

[-11.8066235, -1.4690188, -0.2608915, -0.2608915, 2.4460965, 2.9193177],

[-11.8527288, -1.1662084, -0.2010273, -0.2010273, 2.5009169, 2.8320152],

[-11.8926942, -0.8762104, -0.1454735, -0.1454735, 2.5443783, 2.7507484],

[-11.9254567, -0.6061277, -0.0964072, -0.0964072, 2.5736503, 2.6789196],

[-11.9509776, -0.3673613, -0.0563012, -0.0563012, 2.5815397, 2.6184257],

[-11.9692267, -0.1739367, -0.0256014, -0.0256014, 2.5707386, 2.5737310],

[-11.9801842, -0.0454434, -0.0064917, -0.0064917, 2.5460056, 2.5460056],

[-11.9838380, -0.0000000, -0.0000000, 0.0000000, 2.5366058, 2.5366058],

[-11.9789660, -0.0349409, -0.0215615, -0.0215615, 2.5237184, 2.5640966],

[-11.9643530, -0.1352701, -0.0835025, -0.0835025, 2.4858121, 2.6438116],

[-11.9400076, -0.2899220, -0.1788905, -0.1788905, 2.4250080, 2.7687762],

[-11.9059442, -0.4858158, -0.2987785, -0.2987785, 2.3444149, 2.9313344],

[-11.8612779, -0.7119196, -0.4352983, -0.4352983, 2.2478183, 3.1225671],

[-11.8078753, -0.9595170, -0.5837214, -0.5837214, 2.1385751, 3.3353605],

[-11.7454194, -1.2223619, -0.7376522, -0.7376522, 2.0172437, 3.5653198],

[-11.6728389, -1.4957007, -0.8945203, -0.8945203, 1.8931580, 3.8092456],

[-11.5907049, -1.7767146, -1.0512696, -1.0512696, 1.7661634, 4.0641912],

[-11.4990700, -2.0632765, -1.2059763, -1.2059763, 1.6386096, 4.1823592],

[-11.4003399, -2.3543119, -1.3573133, -1.3573133, 1.5113288, 4.0167186],

[-11.2899226, -2.6478128, -1.5039317, -1.5039317, 1.3882175, 3.8015106],

[-11.1710743, -2.9461080, -1.6452192, -1.6452192, 1.2664448, 3.5572711],

[-11.0421552, -3.2433290, -1.7803657, -1.7803657, 1.1531599, 3.3052621],

[-10.9041173, -3.5414797, -1.9089265, -1.9089265, 1.0464478, 3.0517637],

[-10.7570636, -3.8400940, -2.0304788, -2.0304788, 0.9471351, 2.8025151],

[-10.6011076, -4.1387422, -2.1446779, -2.1446779, 0.8559305, 2.5611603],

[-10.4363740, -4.4370105, -2.2512402, -2.2512402, 0.7734489, 2.3300747],

[-10.2621747, -4.7386807, -2.3493852, -2.3493852, 0.6964981, 2.1100307],

[-10.0803904, -5.0349985, -2.4400488, -2.4400488, 0.6329780, 1.9036346],

[-9.8903977, -5.3296798, -2.5224803, -2.5224803, 0.5779666, 1.7108042],

[-9.6921356, -5.6222234, -2.5964512, -2.5964512, 0.5351504, 1.5379290],

[-9.4859151, -5.9123342, -2.6620426, -2.6620426, 0.5032167, 1.3730077],

[-9.2719445, -6.1995027, -2.7190517, -2.7190517, 0.4825562, 1.2224466],

[-9.0471092, -6.4808973, -2.7661774, -2.7661774, 0.4745250, 1.0871674],

[-8.8210965, -6.7614116, -2.8063059, -2.8063059, 0.4773023, 0.9623785],

[-8.5853238, -7.0371733, -2.8372365, -2.8372365, 0.4923205, 0.8545645],

[-8.3428136, -7.3082926, -2.8593532, -2.8593532, 0.5198220, 0.7607190],

[-8.0924962, -7.5741203, -2.8725455, -2.8725455, 0.5611004, 0.6830821],

[-7.8374693, -7.8374693, -2.8770926, -2.8770926, 0.6143930, 0.6143930]]

test "every band at every point of the path agrees with the numpy reference (eV)" {

assert_lt(max(abs(EVB[:, 0:shape(REFPATH)[1]] - REFPATH)), 1.0e-5)

}

# --- the physics: published LDA values, loose on purpose -----------------------

NV = NELEC / 2.0

GAP = 1.0e9

KGAP = 0.0

for p to NPTS {

if EVB[p, NV] < GAP {

GAP = EVB[p, NV]

KGAP = p

}

show GAP

test "the valence-band top at Gamma is a TRIPLET: the density has the crystal's full symmetry" {

assert_near(EVB[GAMMA_P, NV - 1.0] - EVB[GAMMA_P, NV - 3.0], 0.0, 1.0e-4)

}

test "valence band width is about 12 eV" {

assert_near(EVB[GAMMA_P, 0.0], -11.95, 0.15)

}

test "silicon is an INDIRECT-gap semiconductor, minimum near 0.85 of Gamma-X" {

assert_lt(GAP, EVB[GAMMA_P, NV])

assert_near((KGAP - GAMMA_P) / NSEG, 0.85, 0.1)

}

test "the LDA gap is about 0.5 eV -- less than half the measured 1.17 eV" {

assert_near(GAP, 0.5, 0.1)

assert_lt(GAP, 0.6 * 1.17)

}

test "zone-boundary valence states: X4 and L3'" {

assert_near(EVB[LAST, NV - 1.0], -2.87, 0.12)

assert_near(EVB[0.0, NV - 1.0], -1.2, 0.1)

}

figure("Silicon: LDA band structure", "wavevector along the path", "energy from the valence-band maximum (eV)")

for b to NBANDS {

plot(PATHX, EVB[:, b], "band")

}

xtick("L", PATHX[0])

xtick("Γ", PATHX[GAMMA_P])

xtick("X", PATHX[LAST])