あるケミストの独り言(winchemwinの日記)

ケミスト(化学者)の視点で、面白そうな情報(シミュレーション関係など)を発信

Pyscfによる量子化学計算-Tkinterアプリ_6


 引き続きpyscfの計算実施用アプリについて紹介したいと思います。
 今回はマリケンの電荷計算の部分になります。
以下のコードはこれまでと同様に計算用の関数設定です。Tkinterアプリからの入力値(method, function, baseset等)を受けて、計算の実行を行っています。methodの違い(HF, DFT, MP2)はif文で分岐させています。Muliken 電荷は、以前の記事でも紹介したようにmulliken_pop()の関数を用いることで計算可能です。計算結果はmul_popに格納されます。MP2レベルでの計算(そもそもMullken電荷はあまり高レベルでの計算には適してはいませんが)ではエラーが発生しましたので、計算できない旨のアラートを表示させるようにしています。

def Mul_Charge():
    # Calculation of Mulliken Charges
    
    # Display 'start'
    Label15=ttk.Label(text=u'Calculations(Mulliken analysis) were started', font=("Times","12"))
    Label15.place(x=20, y=570)
    
    # Molecular orbital analysis (method in HF and MP2, function in DFT)
    if meth=='HF':   
        mol_mul=gto.M(atom=comp_geo, basis=base,verbose=4,output=svname+'-hfmul.xyz', max_memory = mem)
        mol_mul.spin=mul
        mol_mul.charge=chr
        lib.num_threads(thre)
        mf_hf_mul=scf.RHF(mol_mul)
        mf_hf_mul.chkfile=svname+'-hfmul.chk'
        mf_hf_mul.kernel()    
        
        mul_pop=mf_hf_mul.mulliken_pop(mol_mul)
        
        
    elif meth=='DFT':
        mol_mul=gto.M(atom=comp_geo, basis=base,verbose=4,output=svname+'-dftmul.xyz', max_memory = mem)
        mol_mul.spin=mul
        mol_mul.charge=chr
        lib.num_threads(thre)
        mf_dft_mul=dft.RKS(mol_mul)
        mf_dft_mul.xc=func
        mf_dft_mul=mf_dft_mul.newton()
        mf_dft_mul.chkfile=svname+'-dftmul.chk'
        mf_dft_mul.kernel()
        
        mul_pop=mf_dft_mul.mulliken_pop(mol_mul)
        
        
    elif meth=='MP2':
        # MP2 ではMulliken Charge 計算でエラー。計算不可の表示
    
        Label15.place_forget()
        Label11_1=ttk.Label(text=u'Mulliken Charge calculations by MP2 can not be performed', font=("Times","12"))
        Label11_1.place(x=20, y=570)
        
        return
        
    # Display 'finish'
    Label15.place_forget()
    Label16=ttk.Label(text=u'Calculation was finished', font=("Times","12"))
    Label16.place(x=20, y=570)   

以下は電荷データを表示させるコードです。電荷に関するデータは第二要素に格納されているため、そのデータを「Charge」として取り出しています。またそれらデータを各原子に対応させたるためfor文以下の処理を行うと共に、pandasのデータフレームとして表示させる形にしています。表示についてはこれまでと同様にTkinterのToplevelのサブウィンドウを活用する形です。また、以前にも紹介させていただいた化学あたらしいカタチさんの記事を参考に化学構造上に電荷マッピングした図を作成し、別の画面で表示させるようにもしています(sub_window4=tk.Toplevel()以下)。

# Calculation of Mulliken charges
   
    charge=mul_pop[1]

    atom_symbol=[]
    mulliken_data=[]

    for n in range(len(mol_mul._atom)):
        print (mol_mul.atom_symbol(n), round(charge[n],4))
        atom_symbol.append(mol_mul.atom_symbol(n))
        mulliken_data.append(round(charge[n],4))

    Mullikendf=pd.DataFrame(({'Atom': atom_symbol, 'Mulliken Charge':mulliken_data}))

    # Display results of Mullken data
    sub_window3=tk.Toplevel()
    sub_window3.title('Calculations results (Mulliken Data)')
    sub_window3.geometry('720x540')
    
    LabelS_9=ttk.Label(sub_window3, text='Mulliken Charges', font=("Arial","14", 'bold'))
    LabelS_9.place(x=10, y=10, width=300)
    
    LabelS_10=ttk.Label(sub_window3, text=Mullikendf, font=("Arial","12"))
    LabelS_10.place(x=10, y=50)
    
    # 電荷の画像描画(化学あたらしいカタチさんの記事参考)
    sub_window4=tk.Toplevel()
    sub_window4.title('Calculations results (Mulliken Charges)')
    sub_window4.geometry('720x540')
    
    mol_file1=Chem.MolFromXYZFile(comp_geo)
    rdDetermineBonds.DetermineBonds(mol_file1)
    print (Chem.MolToMolBlock(mol_file1))
    mol_file1H=Chem.AddHs(mol_file1)

    fig=SimilarityMaps.GetSimilarityMapFromWeights(mol_file1H, mul_pop[1], colorMap='RdBu')
    fig.savefig('fig1.png', bbox_inches='tight')
    canvas1 = tk.Canvas(sub_window4, bg="white", height=500, width=700)
    canvas1.pack()
    
    figure1=tk.PhotoImage(file='fig1.png', height=500, width=700)
    canvas1.create_image(0, 0, image=figure1, anchor=tk.NW)
    
    sub_window4.mainloop()
  

以上、今回はマリケンの電荷計算の関数のコードについて紹介させていただきました。
次回は双極子モーメント計算の関数のコードについて紹介したいと思います。