Привет! Мне нужно поменять xtb в своем коде, где xtb используется для оптимизации геометрии соли и расчета HOMO и LUMO, а поменять нужно на psi4. Метод и базис нужно выбрать самый быстрый и наиболее легкий, чтобы можно было провести расчеты на слабом железе. Вот код в который нужно модифицировать "import pandas as pd
from rdkit import Chem
from rdkit.Chem import AllChem, Descriptors, Descriptors3D, MolSurf, GraphDescriptors, Fragments
from rdkit.Chem.EState import EState_VSA
import subprocess
import tempfile
import os
import re
import shutil
from tqdm import tqdm

# Установка локали для XTB
os.environ['LANG'] = 'C.UTF-8'
os.environ['LC_NUMERIC'] = 'C'

# Конфигурация
input_excel_file = 'BaseOfA.xlsx'
output_excel_file = 'dataset_with_descriptors_xtb_opt.xlsx'
inchi_column_name = 'InChI'
morgan_radius = 2
morgan_nbits = 1024

KNOWN_ANIONS = {
    '[Cl-]': 'Cl-', '[Br-]': 'Br-', '[I-]': 'I-', '[F-]': 'F-',
    'F[P-](F)(F)(F)(F)F': 'PF6-', 'F[B-](F)(F)F': 'BF4-',
    '[O-]S(=O)(=O)[O-]': 'SO4^2-',
    'C(F)(F)(F)S(=O)(=O)[N-]S(=O)(=O)C(F)(F)(F)': 'Tf2N-',
    'BF4': 'BF4-', 'PF6': 'PF6-', '[IH-]': 'I-', 'I': 'I-',
}

EMPTY_ROW_TEMPLATE = {
    'Error': '', 'HOMO_eV_XTB': float('nan'), 'LUMO_eV_XTB': float('nan'), 
    'HOMO_LUMO_gap_eV_XTB': float('nan'),
    'Asphericity': float('nan'), 'Eccentricity': float('nan'), 
    'InertialShapeFactor': float('nan'), 'RadiusOfGyration': float('nan'), 
    'SpherocityIndex': float('nan'), 'NPR1': float('nan'), 'NPR2': float('nan'),
    'NumHeteroatoms_cation': float('nan'), 'fr_NH2_cation': float('nan'), 
    'fr_Ar_N_cation': float('nan'), 'fr_quatN_calc': float('nan'),
    'MolLogP_cation': float('nan'), 'NumQuatNitrogen_cation': float('nan'),
    'AnionType': "N/A", 'MorganFingerprint_cation': "N/A",
    'Chi1n_cation': float('nan'), 'Chi4v_cation': float('nan'), 
    'Chi4n_cation': float('nan'), 'VSA_EState10_cation': float('nan'),
    'PEOE_VSA4_cation': float('nan'), 'PEOE_VSA6_cation': float('nan'),
    'Kappa1_cation': float('nan'), 'Kappa2_cation': float('nan'), 
    'Kappa3_cation': float('nan'), 'MolWt_cation': float('nan')
}

def get_xtb_optimized_mol(mol, charge=0, xtb_path='xtb'):
    """
    Генерирует начальную 3D-структуру в RDKit, затем оптимизирует её через XTB
    """
    temp_dir = tempfile.mkdtemp()
    xyz_path = os.path.join(temp_dir, "input.xyz")
    
    try:
        # Генерация начальной 3D-структуры в RDKit
        mol_3d = Chem.AddHs(mol)
        params = AllChem.ETKDGv3()
        params.randomSeed = 0xf00d
        success = AllChem.EmbedMolecule(mol_3d, params)
        if success != 0:
            # Если RDKit не справился, используем случайные координаты
            params.useRandomCoords = True
            AllChem.EmbedMolecule(mol_3d, params)
            
        # Сохраняем XYZ файл
        Chem.MolToXYZFile(mol_3d, xyz_path)
        
        # Запускаем XTB для оптимизации
        xtb_cmd = [xtb_path, xyz_path, "--gfn", "2", "--opt", "--chrg", str(charge)]
        result = subprocess.run(xtb_cmd, cwd=temp_dir, capture_output=True, text=True, timeout=300)
        
        # Проверяем успешность оптимизации
        if result.returncode != 0:
            shutil.rmtree(temp_dir)
            return None, result.stderr
        
        # Считываем оптимизированную структуру из xtbopt.xyz
        optimized_xyz_path = os.path.join(temp_dir, "xtbopt.xyz")
        if not os.path.exists(optimized_xyz_path):
            shutil.rmtree(temp_dir)
            return None, "Optimized XYZ file not found"
        
        # Читаем XYZ обратно в RDKit
        optimized_mol = Chem.MolFromXYZFile(optimized_xyz_path)
        if optimized_mol is None:
            return None, "Failed to read optimized XYZ"
        
        # Добавляем H, если их нет
        optimized_mol = Chem.AddHs(optimized_mol)
        
        # Удаляем временную папку
        shutil.rmtree(temp_dir)
        return optimized_mol, None
    
    except Exception as e:
        shutil.rmtree(temp_dir)
        return None, str(e)

def get_homo_lumo_from_xtb_log(log_content):
    """
    Извлекает HOMO/LUMO из лога XTB
    """
    homo_ev = lumo_ev = float('nan')
    
    # Ищем HOMO/LUMO в конце лога
    for line in reversed(log_content.split('\n')):
        if "HOMO:" in line and "LUMO:" in line:
            match_homo = re.search(r"HOMO:\s*(-?\d+\.\d+)", line)
            match_lumo = re.search(r"LUMO:\s*(-?\d+\.\d+)", line)
            if match_homo and match_lumo:
                homo_ev = float(match_homo.group(1))
                lumo_ev = float(match_lumo.group(1))
                break
    
    # Если не найдены, ищем в секции SUMMARY
    if pd.isna(homo_ev) or pd.isna(lumo_ev):
        for line in log_content.split('\n'):
            if "SUMMARY" in line:
                for subline in log_content.split('\n')[log_content.index(line):]:
                    if "(HOMO)" in subline:
                        homo_ev = float(re.search(r"(-?\d+\.\d+)\s*eV", subline).group(1))
                    if "(LUMO)" in subline:
                        lumo_ev = float(re.search(r"(-?\d+\.\d+)\s*eV", subline).group(1))
                    if not pd.isna(homo_ev) and not pd.isna(lumo_ev):
                        break
    return homo_ev, lumo_ev

def process_inchi(inchi_str):
    results = EMPTY_ROW_TEMPLATE.copy()
    
    # Парсим InChI
    mol = Chem.MolFromInchi(inchi_str, sanitize=True, removeHs=False)
    if mol is None:
        results['Error'] = "Invalid InChI"
        return results
    
    # Разделяем на фрагменты
    frags = Chem.GetMolFrags(mol, asMols=True, sanitizeFrags=True)
    cation = None
    anions = []
    
    for frag in frags:
        charge = Chem.GetFormalCharge(frag)
        if charge > 0:
            cation = frag
        elif charge < 0:
            anions.append({
                'smiles': Chem.MolToSmiles(frag),
                'mol': frag,
                'charge': charge
            })
    
    if cation is None:
        results['Error'] = "No cation fragment found"
        return results
    
    # Оптимизация структуры через XTB
    optimized_mol, error = get_xtb_optimized_mol(cation, charge=Chem.GetFormalCharge(cation))
    if optimized_mol is None:
        results['Error'] = f"XTB optimization failed: {error}"
        return results
    
    # Считываем HOMO/LUMO из лога XTB
    with open(os.path.join(tempfile.tempdir, "xtb.out"), 'r') as f:
        log_content = f.read()
    homo_ev, lumo_ev = get_homo_lumo_from_xtb_log(log_content)
    
    # Обновляем результаты
    results.update({
        'HOMO_eV_XTB': homo_ev,
        'LUMO_eV_XTB': lumo_ev,
        'HOMO_LUMO_gap_eV_XTB': lumo_ev - homo_ev if not pd.isna(lumo_ev) else float('nan')
    })
    
    # Рассчитываем RDKit дескрипторы на оптимизированной структуре
    rdkit_desc = get_rdkit_descriptors(optimized_mol, cation, anions)
    results.update(rdkit_desc)
    
    return results

def get_rdkit_descriptors(mol_3d, cation_mol, anions_info):
    desc = {}
    
    # Геометрические дескрипторы
    if mol_3d and mol_3d.GetNumConformers() > 0:
        try:
            desc['Asphericity'] = Descriptors3D.Asphericity(mol_3d)
            desc['Eccentricity'] = Descriptors3D.Eccentricity(mol_3d)
            desc['InertialShapeFactor'] = Descriptors3D.InertialShapeFactor(mol_3d)
            desc['RadiusOfGyration'] = Descriptors3D.RadiusOfGyration(mol_3d)
            desc['SpherocityIndex'] = Descriptors3D.SpherocityIndex(mol_3d)
            desc['NPR1'] = Descriptors3D.NPR1(mol_3d)
            desc['NPR2'] = Descriptors3D.NPR2(mol_3d)
        except Exception as e:
            desc.update({
                'Asphericity': float('nan'), 'Eccentricity': float('nan'), 
                'InertialShapeFactor': float('nan'), 'RadiusOfGyration': float('nan'), 
                'SpherocityIndex': float('nan'), 'NPR1': float('nan'), 'NPR2': float('nan')
            })
    
    # Остальные дескрипторы
    try:
        desc['NumHeteroatoms_cation'] = Descriptors.NumHeteroatoms(cation_mol)
        desc['fr_NH2_cation'] = Fragments.fr_NH2(cation_mol)
        desc['fr_Ar_N_cation'] = Descriptors.NumAromaticHeterocycles(cation_mol)
        
        # Четвертичные азоты
        desc['NumQuatNitrogen_cation'] = sum(
            1 for atom in cation_mol.GetAtoms() 
            if atom.GetAtomicNum() == 7 and atom.GetFormalCharge() == 1 and atom.GetDegree() == 4
        )
        desc['fr_quatN_calc'] = desc['NumQuatNitrogen_cation']
        
        desc['MolLogP_cation'] = Descriptors.MolLogP(cation_mol)
    except Exception as e:
        return {"Error": f"RDKit descriptors failed: {str(e)[:50]}", **{k: float('nan') for k in EMPTY_ROW_TEMPLATE}}
    
    # Анионы
    anion_types = []
    for anion in anions_info:
        smiles = anion['smiles']
        if smiles in KNOWN_ANIONS:
            anion_types.append(KNOWN_ANIONS[smiles])
        else:
            try:
                mol = anion['mol']
                formula = Chem.rdMolDescriptors.CalcMolFormula(mol)
                if formula in KNOWN_ANIONS:
                    anion_types.append(KNOWN_ANIONS[formula])
            except:
                anion_types.append("Unknown")
    desc['AnionType'] = ";".join(anion_types) if anion_types else "N/A"
    
    # Morgan Fingerprint
    try:
        mfpgen = AllChem.GetMorganGenerator(radius=2, fpSize=1024)
        fp = mfpgen.GetFingerprintAsBitVect(cation_mol)
        desc['MorganFingerprint_cation'] = fp.ToBitString()
    except:
        desc['MorganFingerprint_cation'] = "Error_MFP"
    
    # Остальные дескрипторы
    try:
        desc['Chi1n_cation'] = GraphDescriptors.Chi1n(cation_mol)
        desc['Chi4v_cation'] = GraphDescriptors.Chi4v(cation_mol)
        desc['Chi4n_cation'] = GraphDescriptors.Chi4n(cation_mol)
        desc['VSA_EState10_cation'] = EState_VSA.EState_VSA_(cation_mol)[9]
        desc['PEOE_VSA4_cation'] = MolSurf.PEOE_VSA_(cation_mol)[3]
        desc['PEOE_VSA6_cation'] = MolSurf.PEOE_VSA_(cation_mol)[5]
        desc['Kappa1_cation'] = Descriptors.Kappa1(cation_mol)
        desc['Kappa2_cation'] = Descriptors.Kappa2(cation_mol)
        desc['Kappa3_cation'] = Descriptors.Kappa3(cation_mol)
        desc['MolWt_cation'] = Descriptors.MolWt(cation_mol)
    except Exception as e:
        desc.update({
            'Chi1n_cation': float('nan'), 'Chi4v_cation': float('nan'), 
            'Chi4n_cation': float('nan'), 'VSA_EState10_cation': float('nan'),
            'PEOE_VSA4_cation': float('nan'), 'PEOE_VSA6_cation': float('nan'),
            'Kappa1_cation': float('nan'), 'Kappa2_cation': float('nan'),
            'Kappa3_cation': float('nan'), 'MolWt_cation': float('nan')
        })
    
    return desc

if __name__ == "__main__":
    # Чтение Excel
    df = pd.read_excel(input_excel_file)
    
    all_descriptors = []
    for idx, row in tqdm(df.iterrows(), total=len(df), desc="Processing"):
        inchi = row[inchi_column_name]
        if not isinstance(inchi, str) or not inchi.startswith("InChI="):
            all_descriptors.append(EMPTY_ROW_TEMPLATE.copy())
            continue
        
        try:
            results = process_inchi(inchi)
            all_descriptors.append(results)
        except Exception as e:
            error_row = EMPTY_ROW_TEMPLATE.copy()
            error_row['Error'] = f"Critical error: {str(e)[:100]}"
            all_descriptors.append(error_row)
    
    # Создаем DataFrame и сохраняем
    descriptors_df = pd.DataFrame(all_descriptors)
    df_final = pd.concat([df, descriptors_df], axis=1)
    df_final.to_excel(output_excel_file, index=False)"