How to read data from an input and form two-electrons integral?

Question

I am doing a code that can calculate the hartree-fock energy. I want to form the two-electrons integral by reading the data from the input like this one below. The first four number are the indices of array(Matrix) and not all the indices are written because some of them are the same such as 2 1 1 1 == 1 2 1 1 and the float number is the value of index and value of my two-electrons that I want to form an array using these float number

1    1    1    1      0.7746059439198979
2    1    1    1      0.4441350695399573
2    1    2    1      0.2970603935859659
2    2    1    1      0.5696940113278337
2    2    2    1      0.4441350695399575
2    2    2    2      0.7746059439198979

I use my code and I didn't get what I want.

My code

electron2 = []
for line in open('Two_Electron.txt'):
    r,c,d,e,v = line.split()
    r = int(r)-1
    c = int(c)-1
    d = int(d)-1
    e = int(e)-1    
    v = float(v)
    if c == 0: 
        electron2.append( [v] )
    else:
        electron2[-1].append(v)

for i in range(len(electron2)-1):
    for j in range(i+1,len(electron2)):
        for k in range(len(electron2)-1):
            for l in range(k+2,len(electron2)):
                electron2[i].append( electron2[j][i] )
                electron2[k].append( electron2[l][k] )              
                
V_ee = np.array(electron2)

print(V_ee)

The output I got is

[list([0.7746059439198979, 0.4441350695399573, 0.2970603935859659, 0.2970603935859659, 0.2970603935859659, 0.2970603935859659])
 list([0.4441350695399573, 0.5696940113278337])
 list([0.2970603935859659, 0.5696940113278337, 0.4441350695399575, 0.7746059439198979])]

The output I want to get

[[[[0.77460594 0.4441351 ]
   [0.4441351  0.56969403]]

  [[0.4441351  0.29706043]
   [0.29706043 0.4441351 ]]]


 [[[0.4441351  0.29706043]
   [0.29706043 0.4441351 ]]

  [[0.56969403 0.4441351 ]
   [0.4441351  0.77460594]]]]

What did I do wrong?

Answer 1

I'm not sure how your current code is supposed to work. Focusing just on the first half:

electron2 = []
for line in open('Two_Electron.txt'):
    r,c,d,e,v = line.split()
    r = int(r)-1
    c = int(c)-1
    d = int(d)-1
    e = int(e)-1    
    v = float(v)
    if c == 0: 
        electron2.append( [v] )
    else:
        electron2[-1].append(v)

It creates three nested lists:

[list([0.7746059439198979, 0.4441350695399573])
 list([0.4441350695399573, 0.5696940113278337])
 list([0.2970603935859659, 0.5696940113278337, 0.4441350695399575, 0.7746059439198979])]

Since the remaining part of the code just appends to these list (and, at least in this case, just the first one), its difficult to see how you intend to convert these three lists to a 4D array.

To me, the conceptually simplest way to do this would be to just take these values and fill them into the symmetry equivalent elements of an already initialized 4D array.

NBasis=2 #Could be determined with earlier code or just pre-entered
electron2 = np.zeros((NBasis,NBasis,NBasis,NBasis))
for line in open('Two_Electron.txt'):
    r,c,d,e,v = line.split()
    r = int(r)-1
    c = int(c)-1
    d = int(d)-1
    e = int(e)-1    
    val = float(v)

    #Fill in all symmetry equivalent elements
    electron2[r,c,d,e]=val                             
    electron2[c,r,d,e]=val                                    
    electron2[r,c,e,d]=val    
    electron2[c,r,e,d]=val    
    
    #Need to take complex conjugate if not all real integrals
    #Assume integrals are indexed with chemist notation
    electron2[d,e,r,c]=val    
    electron2[e,d,r,c]=val    
    electron2[d,e,c,r]=val    
    electron2[e,d,c,r]=val    

If you have control over the input you are reading in, I would just put NBasis at the beginning of it to avoid having to set it manually or having to go through a more complicated process to determine it automatically.

There are definitely other things that could be done to make this cleaner/faster/etc, but this is probably the most direct way to do this.