Page 1 of 1

How to output the microscopic current density j(r,t) for constructing the gauge-invariant polarization density P(r,t)?

Posted: Tue May 12, 2026 7:47 am
by chichenyulong
Dear YAMBO Developers and Users,

I am working on the second-order nonlinear optical response (SHG) in crystalline solids. My goal is to obtain the gauge-invariant, spatially-resolved induced polarization density P(r,t).

After carefully studying the theoretical foundations of the YAMBO code, I understand that:

yambo_nl implements the modern theory of polarization via the Berry phase (Souza et al., PRB 69, 085106) to compute the macroscopic polarization P(t).

This elegant formulation directly yields the cell-averaged P(t) as a geometric phase of the occupied manifold, intentionally bypassing the ill-defined microscopic P(r,t) in periodic systems.

Therefore, the o.polarization_F1 output contains only the macroscopic, integrated signal—the spatially-resolved P(r,t) is not directly available.

Given this, I believe the only physically rigorous route to obtain a gauge-invariant P(r,t) is to time-integrate the microscopic current density j(r,t):
text

P(r,t) = ∫ j(r,t′) dt′

The underlying theory is consistent: the current density computed from the evolving Bloch states satisfies the continuity equation ∇·j = -∂ρ/∂t and is a well-defined, gauge-invariant vector field at each point in the unit cell.

My primary questions are:

Microscopic current density output: Is there any existing mechanism (keyword, verbosity flag, or hidden feature) to output the spatially-resolved current density j(r, t) on the FFT grid during a yambo_nl propagation? I see EvalCurrent for the macroscopic j(t), but I need the grid-resolved quantity.

Code modification guidance: If this feature is not yet available, could you please point me to:

The specific subroutine and source file where j(r) is computed on the FFT grid just before being spatially integrated?

Any advice on the recommended way to write this data to disk (e.g., in cube or NETCDF format)?

Reference approach: Is there any established workflow within the YAMBO community for constructing a spatially-resolved polarization density from the real-time dynamics?

For completeness, I should mention that I also explored the route of obtaining P(r,t) from ρ(r,t) via the continuity equation. However, the non-uniqueness of that inversion convinced me that the current-density route is the only physically sound approach.

Thank you sincerely for your time and for making YAMBO such a powerful tool. Any guidance would be immensely appreciated.

Re: How to output the microscopic current density j(r,t) for constructing the gauge-invariant polarization density P(r,t

Posted: Sat May 16, 2026 11:15 pm
by Davide Sangalli
Dear chichenyulong,
in general, the history of microscopic quantities is not printed, also because the corresponding databases would be rather large. However, we do have options to print something which would do.

First of all, please notice that both yambo_rt and yambo_nl work in the KS basis. Accordingly, what you can get are quantities in such basis set.
From them, you can re-construct quantities in real space by loading the equilibrium KS wave-functions.

With yambo_rt you can add the input flag

Code: Select all

SaveGHistory
to get a netcdf4 database with rho_{nmk}(t). You can select how often this is printed to disk.
This is for example discussed here:
https://wiki.yambo-code.eu/wiki/index.p ... ix_to_disk

For example, ypp can load the history of the density matrix to get the time dependent density n(r,t)

However, as you probably know, yambo_rt does not contain a correct coupling with external fields beyond the linear regime. A similar option has been recently coded in yambo_nl, but (i) it is not tested-nor documented, and (ii) you need to switch to the lumen fork of the yambo code: https://gitlab.com/lumen-code/lumen

There (you can also use release 2.0) you can use

Code: Select all

SaveVbhistory
Which contains the history of the time evolution of the occupied KS states, projected onto the equilibrium KS states.
This is used if yambo_nl runs in "pump and probe mode", which should be automatically activated if you set in input

Code: Select all

% NLEnRange
 -1.00 | -1.00 |
%
NLAngSteps=0
Best,
D.

P.S.: please add your signature in your user settings. It is a rule of the forum