API Usage#

Summary of user API functions#

The user-facing APIs that can be called from within a DFT code are declared in librpa.h (C), librpa.hpp (C++), and librpa_f03.f90 (Fortran 2003).

These APIs fall into three groups:

  • global environment management (initialization and finalization)

  • input parsing

  • computation

They are centered around two management structures:

Note

Some advanced options and diagnostic output paths are marked Experimental in the runtime parameter guide. They are exposed for development and validation workflows, but are not as broadly tested as the default RPA, EXX, and G0W0 paths.

Handler object#

The LibrpaHandler (C) / librpa::Handler (C++) / type(LibrpaHandler) (Fortran) is the central object that manages the LibRPA computation. It encapsulates all input data and computation state.

The typical workflow is:

  1. Create a handler using the MPI communicator

  2. Set input data via setter methods

  3. Call computation functions

  4. Free the handler when done

C:

LibrpaHandler *h = librpa_create_handler(MPI_COMM_WORLD);
// ... set input data and compute ...
librpa_destroy_handler(h);

C++:

librpa::Handler h(MPI_COMM_WORLD);
// ... set input data and compute ...
h.free();

Fortran:

type(LibrpaHandler) :: h
call h%init(MPI_COMM_WORLD)
! ... set input data and compute ...
call h%free()

Environment setup#

Runtime environment for LibRPA should be properly setup to parse the input data and call computation functions. They are managed by the environment setup APIs (C/C++/Fortran):

Any calls to input parsing and computation APIs should be surrounded by calls to librpa_init_global and librpa_finalize_global.

To configure runtime parameters for computation, an LibrpaOptions object should be created, modified and passed to computation functions. Parameters are attributes of LibrpaOptions object, and can be initialized using librpa_init_options. For API calls, set the time-frequency grid type explicitly before computation. For example, to configure a minimax grid, the number of frequency points, and Green’s function threshold:

// C API
LibrpaOptions opts;
librpa_init_options(&opts);
opts.tfgrids_type = LIBRPA_TFGRID_MINIMAX;
opts.nfreq = 16;
opts.gf_threshold = 1.e-3;
librpa_set_output_level(LIBRPA_VERBOSE_INFO);
// Optional: read restart data from a directory distinct from output_dir.
librpa_set_restart_from_dir(&opts, "previous_run");
// C++ API
librpa::Options opts;
opts.tfgrids_type = LIBRPA_TFGRID_MINIMAX;
opts.nfreq = 16;
opts.gf_threshold = 1.e-3;
librpa::set_output_level(LIBRPA_VERBOSE_INFO);
opts.set_restart_from_dir("previous_run");
! Fortran API
type(LibrpaOptions) :: opts
call opts%init()
opts%tfgrids_type = LIBRPA_TFGRID_MINIMAX
opts%nfreq = 16
opts%gf_threshold = 1.0d-3
call librpa_set_output_level(LIBRPA_VERBOSE_INFO)
call opts%set_restart_from_dir("previous_run")

The API uses enum/integer constants from the public headers, while the driver reads string values such as tfgrids_type = minimax and maps them to the same internal options. Driver-only parameters, such as task and input_dir, are not part of LibrpaOptions.

Input parsing#

Input parsing APIs common to all tasks include (C/C++/Fortran):

These functions parse information from the hosting DFT code to LibRPA, such as:

  • dimension of lattice vectors, orbital basis and auxiliary basis

  • k-mesh points, eigenvalues and occupation numbers

  • wave functions of mean-field calculation

  • local RI coefficients

  • Coulomb matrices in auxiliary basis

  • data for band structure calculation

They have to be called after the runtime environment is initialized, see Environment setup above.

set_kgrids_kvec takes the k-grid dimensions, the number of loaded SCF k-points, their Cartesian vectors, and optional weights. If the loaded SCF k-points are symmetry-reduced for the Coulomb calculation, call set_kq_mapping with the 0-based mapping from each loaded SCF k-point to its representative Coulomb q-point.

The wavefunction and auxiliary-basis setters accept optional angular-shell metadata. Symmetry-enabled calculations should also provide the basis convention and real-space symmetry operations through set_basis_convention and set_symmetry_operations; the convention values must match those used by the calling electronic-structure code.

The C and C++ APIs additionally provide packed-complex Coulomb setters:

The C functions take interleaved real/imaginary double buffers; the C++ wrappers take std::complex<double> buffers. The Fortran Coulomb setters already accept native complex arrays and therefore do not have separate _packed names.

Analytic head/wing calculations can receive velocity matrices through librpa_set_velocity_matrix() or librpa_set_velocity_matrix_packed(); the corresponding C++ methods are Handler::set_velocity_matrix() and Handler::set_velocity_matrix_packed(). These setters are not currently exposed by the Fortran wrapper. Matrix elements are ordered as [spin][k-point][Cartesian component][row state][column state].

For spinor workflows, initialize the handler with two spinor components, then use the spinor wavefunction setters. Internally, LibRPA selects the normal or spinor path from the mean-field n_spinor attribute through get_n_spinor(). Driver users can select spinor-format input with use_spinor_wfc; API users should pass nspinor = 2 to set_scf_dimension instead.

Computation#

After data have been correctly set up by input parsing functions, computation functions can be called to get quantities of interest:

General

RPA correlation energy:

Exact exchange:

G0W0 self-energy:

G0W0 quasiparticle energies:

These functions return the correlation self-energy together with quasiparticle energies obtained using opts.option_qpe_solver. The C++ methods return a librpa::G0W0QpeResult. They are not currently exposed by the Fortran wrapper; Fortran users can retrieve correlation self-energies through the get_g0w0_sigc_* methods.

G0W0 spectral function:

The C spectral-function APIs take an optional packed complex self-energy output as the last argument. Pass nullptr when only the spectral function is needed. The C++ get_g0w0_spectral_function_* methods return only spectral-function values; use Handler::get_g0w0_spectral_function_with_sigc_kgrid() or Handler::get_g0w0_spectral_function_with_sigc_band_k() to receive a librpa::G0W0SpectralFunctionResult. The Fortran methods expose the continued self-energy as an optional final argument.

Schematic examples#

C++#

The following C++ code gives an impression of how LibRPA APIs should be called in the host DFT code:

#include <librpa.hpp>
#include <mpi.h>

int main()
{
    MPI_Init(...);
    // If LibRI is enabled, MPI thread support should be set to MPI_THREAD_MULTIPLE
    // MPI_Init_thread(&argc, &argv, MPI_THREAD_MULTIPLE, &provided);

    // Before: setup and SCF calculation of hosting program

    // *** Initialize LibRPA runtime environment ***
    librpa::init_global();

    // *** Set up runtime parameters ***
    librpa::Options opts;
    opts.tfgrids_type = LIBRPA_TFGRID_MINIMAX;
    opts.nfreq = 16;
    opts.gf_threshold = 1e-3;
    librpa::set_output_level(LIBRPA_VERBOSE_INFO);

    // *** Create handler ***
    librpa::Handler h(MPI_COMM_WORLD);

    // *** Input parsing ***
    h.set_scf_dimension(nspins, nkpts, nstates, nbasis);
    h.set_wg_ekb_efermi(nspins, nkpts, nstates, wg, ekb, efermi);
    h.set_latvec_and_G(lat_mat, G_mat);
    h.set_kgrids_kvec(nk1, nk2, nk3, nkpts, kvecs, kweights);
    h.set_kq_mapping(map_q_ks);
    h.set_ao_basis_wfc(nbs_wfc);
    h.set_ao_basis_aux(nbs_aux);

    // Optional symmetry metadata
    h.set_basis_convention(bloch_phase, bloch_ratom, angular_order,
                           negative_m_convention, positive_m_convention);
    h.set_symmetry_operations(n_symops, row_conv, rotmats, translations);

    // Loop over atoms and R cells to set LRI coefficients
    for (int i = 0; i < natoms; ++i) {
        for (int j = 0; j < natoms; ++j) {
            for (auto& R : r_cells) {
                h.set_lri_coeff(LIBRPA_ROUTING_AUTO, i, j, nbasis_i, nbasis_j, naux_mu, R.data(), Cs.data());
            }
        }
    }

    // Loop over k-points to set Coulomb matrices
    for (int ik = 0; ik < nkpts; ++ik) {
        h.set_aux_bare_coulomb_k_atom_pair(ik, i, j, naux_i, naux_j, Vq_real.data(), Vq_imag.data(), vq_threshold);
        // For std::complex<double> input, use
        // h.set_aux_bare_coulomb_k_atom_pair_packed(..., Vq.data(), vq_threshold);
    }

    // *** Compute RPA correlation energy ***
    std::vector<std::complex<double>> rpa_contrib_ibzk(nkpts_ibz);
    double Ec = h.get_rpa_correlation_energy(opts, rpa_contrib_ibzk);

    // *** Compute exact exchange ***
    h.build_exx(opts);
    std::vector<double> vexx = h.get_exx_pot_kgrid(opts, nspins, iks_this, i_state_low, i_state_high);

    // *** Compute G0W0 self-energy ***
    h.build_g0w0_sigma(opts);
    std::vector<std::complex<double>> sigc = h.get_g0w0_sigc_kgrid(opts, nspins, iks_this, i_state_low, i_state_high, vxc, vexx);
    librpa::G0W0QpeResult qpe = h.get_g0w0_qpe_kgrid(
        opts, nspins, iks_this, i_state_low, i_state_high, vxc, vexx);
    librpa::G0W0SpectralFunctionResult spectral =
        h.get_g0w0_spectral_function_with_sigc_kgrid(
            opts, nspins, iks_this, i_state_low, i_state_high, omegas_real, vxc, vexx);

    // *** Clean up ***
    h.free();
    librpa::finalize_global();

    MPI_Finalize();
}

Fortran#

use librpa_f03
use mpi

implicit none

type(LibrpaOptions) :: opts
type(LibrpaHandler) :: h
integer :: ierr
integer :: comm, nkpts_ibz
complex(dp), allocatable :: contrib_ibzk(:)

call MPI_Init(ierr)
comm = MPI_COMM_WORLD

! Initialize LibRPA
call librpa_init_global()

! Set up options
call opts%init()
opts%tfgrids_type = LIBRPA_TFGRID_MINIMAX
opts%nfreq = 16
opts%gf_threshold = 1.0d-3

! Initialize handler
call h%init(comm)

! Input parsing
call h%set_scf_dimension(nspins, nkpts, nstates, nbasis)
call h%set_wg_ekb_efermi(nspins, nkpts, nstates, wg, ekb, efermi)
call h%set_latvec_and_G(lat, recplatt)
call h%set_kgrids_kvec(nk1, nk2, nk3, nkpts, kpoints, kweights)
call h%set_kq_mapping(nkpts, map_q_ks)
call h%set_basis_convention(bloch_phase, bloch_ratom, angular_order, &
                            negative_m_convention, positive_m_convention)
call h%set_symmetry_operations(n_symops, row_conv, rotmats, translations)
! ... more input calls ...

! Compute RPA correlation energy
allocate(contrib_ibzk(nkpts_ibz))
Ec = h%get_rpa_correlation_energy(opts, nkpts_ibz, contrib_ibzk)

! Clean up
call h%free()

! Finalize LibRPA
call librpa_finalize_global()

call MPI_Finalize(ierr)