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

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?

• Can you clarify what the input format means, i.e. what do those numbers indicate? Commented Dec 5, 2021 at 12:15
• @SRMaiti I edited the post again to explain what is the format of input Commented Dec 5, 2021 at 12:29
• It may help if you can try to explain the logic behind your code. I could answer giving a fairly different way of putting these integrals together, but that may not address whatever misconception you have that led to your version of the code.
– Tyberius
Commented Dec 5, 2021 at 21:40
• @Tyberius I want to read the two-electron integral ( the data in the input ) in shape of array (2*2*2*2) and depend on the data in the input the size can be different such as (4*4*4*4), then I will use the two electron integral, kinetic integral and nuclear electron integral to calculate the Fock matrix then starting calculate the electronic energy and so on. Is it clear now, or I can explain more ? I tried to do a loop where I can make a range depend on the size of data, so it will be 2 and I can form an array 2*2 for each loop. Commented Dec 6, 2021 at 10:31

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.

• Thank you so much for your help and you time :) Commented Dec 6, 2021 at 15:58
• @Abd-ElazeemMohamed Happy to help. I might clear some of the earlier comments to avoid confusion for future readers.
– Tyberius
Commented Dec 6, 2021 at 16:00