We present a symbolic algorithm designed to investigate the perturbative expansions of the two-body interacting Hamiltonian. This method is part of the recent developments in symbolic integration methods within many-body perturbation theory (MBPT). The direct application of this technique to the one-dimensional Hubbard model yields the coefficients of the grand potential up to sixth order, expressed in terms of both the interacting potential (U expansions) and high-temperature expansions (β expansions). A key component of this process is the ability to merge the β expansions up to any desired order.

How to run the code:

1. Copy all files in the same directory.
2. Open the Mathematica file "High temperature series expansion of 1D Hubbard model.nb".
3. Execute the code by pressing "Shift+Enter".