Descriptor Base Classes

The user-facing Descriptors chapter documents the concrete two-body, three-body, and many-body descriptors. The base classes below are the API hooks for adding new descriptors and embedding functions.

class D2_Base : public tadah::models::D_Base

Base class for all two-body type descriptors.

All two-body descriptors must inherit this class.

Subclassed by tadah::models::D2_BP, tadah::models::D2_Blip, tadah::models::D2_Dummy, tadah::models::D2_EAM, tadah::models::D2_LJ, tadah::models::D2_MIE, tadah::models::D2_ZBL, tadah::models::D2_mJoin

Public Functions

virtual void calc_aed(const int Zi, const int Zj, const double rij, const double rij_sq, tadah::core::aed_type &aed, const double scale = 1) const = 0

Calculate AED.

Calculate Atomic Energy Descriptor for the atom local environment.

virtual void calc_dXijdri(const int Zi, const int Zj, const double rij, const double rij_sq, tadah::core::fd_type &fd_ij, const double scale = 1) const = 0

Calculate FD.

Calculate Force Descriptor for the atom local environment.

For two body this is essentially df/dr where f is a descriptor function. Do not include negative sign (as you mifht think E=-df/dr). The negative sign is included elsewhere in the code.

Computes x-direction only.

The resulting vector must be scaled by the unit directional vector delij/rij where delij = r_i - r_j, and rij is |r_i - r_j|

virtual void calc_all(const int Zi, const int Zj, const double rij, const double rij_sq, tadah::core::aed_type &aed, tadah::core::fd_type &fd_ij, const double scale = 1) const = 0

Calculate AED + FD.

inline virtual bool calc_dDdtheta(const std::string&, int, const int, const int, const double, const double, tadah::core::aed_type&, const double) const

Closed-form derivative of AED with respect to a named hyperparameter.

Used by HPO analytical / hybrid gradients (Phase 2) to assemble dPhi/dtheta one pair at a time, mirroring how calc_aed assembles Phi. The default returns false: “this descriptor does not own

`theta_name`”. Concrete descriptors override for the (theta_name, arg_pos) pairs they implement.

Calling convention:

  • theta_name is the OPTIM key (e.g. “CGRID2B”, “SGRID2B”, “RCUT2B”); always upper-case.

  • arg_pos is the 1-based index into the vector parameter (k for “CGRID2B[k]”), or -1 for scalar parameters.

  • dD_buf is the pre-zeroed aed-shaped accumulator. Concrete implementations += their contribution at fidx + c where c is the descriptor-component index.

  • Return true iff this call wrote something into dD_buf

    . The aggregator uses the return value to detect “no descriptor owned

    this theta” and route the gradient via FD fallback (HYBRID) or raise (ANALYTICAL).

  • Cross-cutting parameters (RCUT2B / RCUTMB) flow through the cutoff (fcut->calc_drcut); each D2 descriptor that uses a non-Cut_Dummy cutoff contributes additively.

Default impl returns false (no-op). New descriptors can be added without implementing this method; HYBRID will FD-fallback on their parameters until an analytical version lands.

inline virtual bool calc_dXijdri_dtheta(const std::string&, int, const int, const int, const double, const double, tadah::core::fd_type&, const double) const

Closed-form derivative of FD with respect to a named hyperparameter.

Used by HPO analytical / hybrid gradients (Phase 2) to assemble dPhi/dtheta force rows one pair at a time, mirroring how calc_dXijdri assembles fd_ij. The default returns false: “this

descriptor does not own `theta_name` for the force chain”. Concrete descriptors override for the (theta_name, arg_pos) pairs they implement.

Computes the radial-direction (column 0) component only; the caller scales by the unit directional vector delij/rij to produce the 3-axis derivative, matching the descriptors_calc.hpp axis expansion pattern.

Calling convention parallels calc_dDdtheta:

  • theta_name is the OPTIM key.

  • arg_pos is the 1-based vector index, or -1 for scalar params.

  • dfd_ij is the pre-zeroed fd-shaped accumulator; this method += into column 0 (radial). Columns 1, 2 are left untouched so the caller’s axis-expansion step (mirroring descriptors_calc.hpp:233-238) can run unchanged on the result.

  • Returns true iff something was written into dfd_ij.

Default impl returns false; concrete descriptors override for the theta keys they support.

inline void calc_fd_approx(const int Zi, const int Zj, double r, tadah::core::fd_type &fd, double h = 1e-8) const

Central difference approximation to FD

Calculate Force Descriptor for the atom local environment using central difference approximation.

fc: a pointer to initialised cutoff object. r: a separation between two particles fd: appropriate size container to store computation h: central difference parameter

class DM_Base : public tadah::models::D_Base

Base class for all many-body type descriptors.

All many-body descriptors must inherit this class.

Subclassed by tadah::models::DM_Blip, tadah::models::DM_Dummy, tadah::models::DM_EAD, tadah::models::DM_EAM, tadah::models::DM_Poly< ChannelFn, RadialBasis, AggMode >, tadah::models::DM_mEAD< F >, tadah::models::DM_mJoin

Public Functions

virtual void calc_aed(tadah::core::rho_type &rho, tadah::core::aed_type &aed) const = 0

Calculate AED.

Calculate Atomic Energy Descriptor for the atom local environment.

virtual void calc_dXijdri_dXjidri(const int Zi, const int Zj, const double rij, const double rij_sq, const tadah::core::Vec3d &vec_ij, tadah::core::rho_type &rhoi, tadah::core::rho_type &rhoj, tadah::core::fd_type &fd_ij, const double scale = 1) const = 0

Calculate FD.

Calculate Force Descriptor between to atoms.

This method works for half NN lists and linear models.

\[ \mathbf{fd}_{ij} \mathrel{+}= \frac{\partial \mathbf{X}_{ij}}{\partial \mathbf{r}_i} + \frac{\partial \mathbf{X}_{ji}}{\partial \mathbf{r}_i} \]

virtual void calc_dXijdri(const int Zi, const int Zj, const double rij, const double rij_sq, const tadah::core::Vec3d &vec_ij, tadah::core::rho_type &rhoi, tadah::core::fd_type &fd_ij, const double scale = 1) const = 0

Calculate FD.

Calculate Force Descriptor between to atoms.

This method works for full NN lists and all models.

mode=0 indicates only x-direction is computed and the resulting descriptor must be multiplied by delx/r.

\[ \mathbf{fd}_{ij} \mathrel{+}= \frac{\partial \mathbf{X}_{ij}}{\partial \mathbf{r}_i} \]

virtual void init_rhoi(tadah::core::rho_type &rhoi) = 0

Resize arrays for a density and F’.

virtual void calc_rho(const int Zi, const int Zj, const double rij, const double rij_sq, const tadah::core::Vec3d &vec_ij, tadah::core::rho_type &rho, const double scale = 1) const = 0

Calculate density.

inline virtual bool calc_drhoi_dtheta(DMThetaKey, int, const int, const int, const double, const double, const tadah::core::Vec3d&, tadah::core::rho_type&, const double) const

Closed-form per-pair contribution to drho_i/dtheta.

Many-body analog of D2_Base::calc_dDdtheta. The MB design-matrix builder accumulates rho_i in a per-pair pre-pass (calc_rho), then computes aed = F(rho) in a per-atom finalisation step (calc_aed). The dPhi/dtheta walker mirrors that structure with TWO methods:

  1. calc_drhoi_dtheta — per-pair drho_i/dtheta accumulator (analog of calc_rho).

  2. calc_daed_dtheta — per-atom (rho_i, drho_i/dtheta) -> daed/dtheta reducer (analog of calc_aed).

Calling convention mirrors D2_Base::calc_dDdtheta:

  • theta_name is the OPTIM key, upper-case (e.g. “CGRIDMB”).

  • arg_pos is 1-based for vector parameters, -1 for scalars.

  • drhoi is the per-atom accumulator, pre-zeroed by the caller.

  • Returns true iff this call wrote anything into drhoi.

  • Cross-cutting parameters (RCUTMB) flow through fcut->calc_drcut; each DM descriptor that uses a non-Cut_Dummy cutoff contributes additively.

Default returns false (no-op). DM_Dummy and DM descriptors that are numerical-only (DM_EAM, DM_mEAD) keep the default.

inline bool calc_drhoi_dtheta(const std::string &theta_name, int arg_pos, const int Zi, const int Zj, const double rij, const double rij_sq, const tadah::core::Vec3d &vec_ij, tadah::core::rho_type &drhoi, const double scale = 1) const

Back-compat string-keyed shim. Decodes once and dispatches to the enum-keyed virtual above. NOT virtual — derived classes override the enum form. Layer-1 tests still call this signature.

inline virtual bool calc_daed_dtheta(DMThetaKey, int, const tadah::core::rho_type&, const tadah::core::rho_type&, tadah::core::aed_type&, const double) const

Closed-form per-atom contribution to daed/dtheta.

Companion to calc_drhoi_dtheta. Called by the dPhi/dtheta walker once per atom AFTER:

  • calc_rho has filled rhoi[rfidx..] = rho_i,

  • calc_aed has filled rhoi[drfidx..] = F(rho_i),

  • calc_drhoi_dtheta has filleddrhoi[rfidx..] = drho_i/dtheta`.

Chain rule for F(rho_i, theta):

daed/dtheta = (dF/drho) (drho/dtheta) + (dF/dtheta)|_{rho fixed}

The first term reads F’(rho) from rhoi[drfidx..] and dotproducts with drhoi[rfidx..]. The second term is non-zero only for cross-cutting parameters that enter the per-L multinomial prefactor f[L,o] = L!/(lx! ly! lz!) / (4 rc)^L — namely RCUTMB.

Returns true iff this call wrote anything into daed_out.

Default returns false.

inline bool calc_daed_dtheta(const std::string &theta_name, int arg_pos, const tadah::core::rho_type &rhoi, const tadah::core::rho_type &drhoi, tadah::core::aed_type &daed_out, const double scale = 1) const

Back-compat string-keyed shim — decodes once, calls enum virtual.

inline virtual bool calc_dXijdri_dtheta(DMThetaKey, int, const int, const int, const double, const double, const tadah::core::Vec3d&, const tadah::core::rho_type&, const tadah::core::rho_type&, tadah::core::fd_type&, const double) const

Closed-form per-pair contribution to d(fd_ij)/dtheta — the many-body sibling of D2_Base::calc_dXijdri_dtheta. Walker calls this once per (atom i, neighbour j) pair AFTER:

  • calc_rho has filled rhoi[rfidx..] = rho_i,

  • calc_aed has filled rhoi[drfidx..] = f[L,o] F(rho_i),

  • calc_drhoi_dtheta has filleddrhoi_dtheta[rfidx..] = drho_i/dtheta`.

Unlike the D2 sibling (which returns a single radial column to be axis-expanded by the walker), the MB descriptor’s calc_dXijdri writes all three axis components directly because the per-orbital polynomial poly = x^lx * y^ly * z^lz couples to direction. This override mirrors that signature — it writes the full axis-expanded d(fd_ij)/dtheta into fd_ij_out(row, 0..2).

Chain rule sketch (for DM_mEAD<F>, expanding path):

fd_ij(p, k) = scale * w_j * f[L,o] * F’(rho_i, p) * term_k, term_k = g * (dpoly/dr_k) + dg * poly * (v_k / r)

with g, dg the radial Gaussian and its r-derivative (both cutoff-weighted). Differentiating w.r.t. theta sees three groups:

  1. Radial theta (CGRIDMB[k], SGRIDMB[k], RCUTMB) — touch (g, dg) directly AND rho via calc_drhoi_dtheta. The F’(rho) chain needs F’’(rho) via F_Base::calc_d2F.

  2. Embedding theta (CEMBFUNC[i], SEMBFUNC[i]) — touch F’ directly via F_Base::calc_dFp_dtheta; rho is fixed.

  3. Cross-cutting RCUTENV (F_mEnv only) — same direct F’ branch.

Calling convention mirrors calc_drhoi_dtheta / calc_daed_dtheta:

  • theta_name is the OPTIM key, upper-case.

  • arg_pos is 1-based for vector parameters, -1 for scalars.

  • fd_ij_out is the per-pair, axis-expanded buffer, pre-zeroed by the caller.

  • Returns true iff this call wrote anything into fd_ij_out.

Default returns false (descriptor doesn’t own this theta on the force side; orchestrator routes to needs_fd).

inline bool calc_dXijdri_dtheta(const std::string &theta_name, int arg_pos, const int Zi, const int Zj, const double rij, const double rij_sq, const tadah::core::Vec3d &vec_ij, const tadah::core::rho_type &rhoi, const tadah::core::rho_type &drhoi_dtheta, tadah::core::fd_type &fd_ij_out, const double scale = 1) const

Back-compat string-keyed shim — decodes once, calls enum virtual.

inline void calc_dXijdri(const int Zi, const int Zj, tadah::core::Vec3d &del, tadah::core::rho_type &rhoi, tadah::core::fd_type &fdij, double scale = 1) const

This is just a convenient wrapper.

inline void calc_rho_mb(const int Zi, const int Zj, tadah::core::Vec3d &del, tadah::core::rho_type &rhoi, double scale = 1) const

This is just a convenient wrapper.

inline virtual size_t rhoi_size() const

Return size of the density array.

inline virtual size_t rhoip_size() const

Return size of the derivative of the embedding energy array.

virtual std::string label() const = 0

Return label of this descriptor.

class F_Base

Subclassed by tadah::models::F_Blip, tadah::models::F_RLR, tadah::models::F_SQ, tadah::models::F_SQRT, tadah::models::F_mEnv< F >

Public Functions

virtual size_t output_dim() const = 0

Returns the number of output channels per orbital per L.

For most embedding functions this equals the radial grid size (ng). For expanding functions (e.g. F_Blip) this equals the embedding function grid size.

virtual double calc_F(double rho, size_t p) const = 0

Calculates the embedding function value.

Parameters:

  • rho – Electron density.

  • p – Index for model parameters.

Returns:

Function value F(rho).

virtual double calc_dF(double rho, size_t p) const = 0

Calculates the first derivative of the embedding function.

Parameters:

  • rho – Electron density.

  • p – Index for model parameters.

Returns:

Derivative dF/d(rho).

inline virtual bool calc_dF_dtheta(double, size_t, const std::string&, int, double&) const

Closed-form derivative of calc_F w.r.t. a Context-tunable embedding-function parameter, with rho held fixed.

Called by DM_mEAD<F>::calc_daed_dtheta on the embedding side (theta_name in {CEMBFUNC, SEMBFUNC, RCUTENV, …} — anything F owns rather than DM_mEAD’s radial layer). The radial side (CGRIDMB / SGRIDMB / RCUTMB) is handled by DM_mEAD directly.

Returns true iff this F owns the (theta_name, arg_pos) pair AND the derivative is non-zero at output channel p. Component-local theta (e.g. F_Blip’s CEMBFUNC[i] only touches p == i-1) must return false when p doesn’t match — the orchestrator counts true returns to decide whether the theta is FD-fallback-able.

Defaults to false (no ownership). F_Blip and F_mEnv override.

Parameters:

  • rho – Electron density at this orbital.

  • p – Output-channel index (0..output_dim()-1).

  • theta_name – OPTIM key, upper-case.

  • arg_pos – 1-based index for vector parameters, -1 for scalars.

  • out – Filled with d calc_F / d theta on success.

inline virtual bool calc_d2F(double, size_t, double&) const

Second derivative of F: F’’(rho, p) = d^2 F / d rho^2.

Used by DM_mEAD<F>::calc_dXijdri_dtheta for the radial-theta chain d(F’(rho))/d theta_radial = F’’(rho) * drho/d theta_radial — i.e. how the rho-dependence of theta propagates through F’(rho) into the force descriptor.

Returns true iff this F has an analytical second derivative AND filled out. Default returns false; DM_mEAD falls back by reporting the radial-theta force derivative as unowned and letting the orchestrator route to FD.

F_Blip, F_mEnv override; F_RLR / F_SQRT / F_SQ leave the default (force-side radial closed forms not yet wired for those F’s).

inline virtual bool calc_dFp_dtheta(double, size_t, const std::string&, int, double&) const

Closed-form derivative of F’(rho) (i.e. dF/drho) wrt a Context-tunable embedding-function parameter, with rho held fixed.

Companion to calc_dF_dtheta but for F’ rather than F: used by DM_mEAD<F>::calc_dXijdri_dtheta to chain d(fd_ij)/d theta_emb = scale * w * f[L,o] * (∂F’/∂theta) * term_k for embedding-side theta (CEMBFUNC, SEMBFUNC, RCUTENV) — rho is independent of these, so the only contribution is through F’.

Returns true iff this F owns the (theta_name, arg_pos) pair AND the derivative is non-zero at output channel p. Component-local theta (CEMBFUNC[i] only touches p == i-1) returns false when p doesn’t match. Defaults to false; F_Blip and F_mEnv override.

Public Static Attributes

static constexpr bool rotationally_invariant = false

True only for \(\mathcal{F}(\rho)=\rho^2\) (F_SQ).

DM_mEAD applies \(\mathcal{F}\) to individual orbital densities \(\rho^{\eta,r_s,l_x,l_y,l_z}\) before the rotational-invariance restoring quadratic sum. This is invariant only when \(\mathcal{F}=\rho^2\); any other choice breaks SO(3) invariance for \(L>0\).

static constexpr bool expands_dimensions = false

Whether this embedding function expands dimensions.

When true, the embedding function takes a single rho value (from radial grid of size 1) and produces multiple output channels (one per embedding grid point). DM_mEAD uses this to adjust its inner loops accordingly.