From the previous posts, I knew that yambo_nl and ypp_nl can function in calculating two-photon absorption of materials, which is a great function. But I am still confused with the detailed process.
In the equations 15~18 in PHYSICAL REVIEW B 98, 165126 (2018) (https://journals.aps.org/prb/pdf/10.110 ... .98.165126), the polarizations result from the excitation of different amplitudes E, E/2 and E/4 were used, for example, P(E), P(E/2), P(E/4). The calculation of yambo_nl and ypp_nl output a lot of files, for example, o.YPP-X_probe_order_0, o.YPP-X_probe_order_1, ....
Code: Select all
X^3(w: w, w, -w) = 8/3 * ( P_int1 - 6*P_int2 + 8*P_int3 ) / (E**3)I also noticed that there are different decisions for getting P(E).
and interestingly, the examples 1 and 2 use X1_int1, X1_int2, and X1_int3 instead of P_int1, P_int2, P_int3. I have no idea whether this treatment is reasonable.
1.
In the post viewtopic.php?t=2141, the PE(E) (or more exactly the X1_int1, X1_int2 and X1_int3 was used) was read from file o.YPP-X_probe_order_1.
I did three separate calculations for field intensities I, I/4, I/16 with the same frequency range. Then I collected the "o.YPP-X_probe_order_1" for the three intensities and renamed them: o.YPP-X_probe_int_1_order_1 (X1_int1, for I), o.YPP-X_probe_int_2_order_1 (X1_int2, for I/4), o.YPP-X_probe_int_3_order_1 (X1_int3, for I/16).
The X3(w: w, w, -w) for TPA was calculated using:
X3(w: w, w, -w) = 8/3 * ( X1_int1 - 3*X1_int2 + 2*X1_int3 ) / (E**2)
where E=E'/(3*10**4), E' is "Efield Amplitude" in r_nloptics.
2.
In https://github.com/yambo-code/yambo-min ... tion_PP.py, the PE(E) was also read from file **_order_1.
Code: Select all
file_name=file_begin+"_int_"+str(iR+1)+"_order_1"
XHI[iR,:,:]=np.genfromtxt(file_name,comments="#")
X_in[:] = (1j*XHI[:,freq,2*dim+1]+XHI[:,freq,2*(dim+1)])
X_out=Extract_third_order(X_in)*Divide_Efield**2
for example the P_2[:,:,3] corresponding to iX=3 and the iX denotes the string after "_order_"
Code: Select all
for iR in range(0,args.nR):
for iX in range(0,args.nX+1):
file_name=file_begin+"_int_"+str(iR+1)+"_order_"+str(iX)
print("Reading %s " % file_name)
try:
XHI[iR,iX,:,:]=np.genfromtxt(file_name,comments="#")
except:
print("Error reading file "+file_name+" !! ")
sys.exit(1)
P_2[:,0,iX]=1j*XHI[1,iX,:,1]+XHI[1,iX,:,2] # x
P_2[:,1,iX]=1j*XHI[1,iX,:,3]+XHI[1,iX,:,4] # y
P_2[:,2,iX]=1j*XHI[1,iX,:,5]+XHI[1,iX,:,6] # z
P_2[:,:,iX]=P_2[:,:,iX]/Divid_Efield[iX]*Scale_factor[iX]
XHI3 = 8.0/3.0*(P[:,:,3] - 6.0*P_2[:,:,3] + 8.0*P_4[:,:,3])*Divid_Efield[3]
The tutorial on http://www.attaccalite.com/lumen/shg_in_AlAs.html suggest that the **_order_2 file reports χ^2
The data in http://www.attaccalite.com/lumen/thg_in_silicon.html suggest that o-THG6.YPP-X_probe_order_3 reports the X^3.Run ypp and it will produce a file called o.YPP-X_probe_order_2 with the χ^2_{xyz} , columns 6 and 7. You can plot the result:
So does it suggest that the file **_order_1 reports the X^1，**_order_2 reports the X^2，and **_order_3 reports the X^3? That is, the data saved in files **_order_n reports corresponding X^(n) coefficients?
However, from the abovementioned posts and the codes in Github, p(E) were read from files **_order_1 or **_order_3, which suggests that the file **_order_XX records the polarization but not the χ^(x) coefficients.
Thank you very much for your help, and again Yambo code is great!
Best regards,
Zhipeng Huang
--------------------
Tongji University, Shanghai, PR China
