feat(librpa): export SOC magnetic symmetry metadata

docs/advanced/interface/LibRPA.md

@@ -0,0 +1,37 @@
+# LibRPA reader-v1 SOC symmetry export
+
+For an LCAO RPA export consumed by a reader-v1 LibRPA build, set:
+
+```text
+basis_type                  lcao
+rpa                         1
+out_librpa_reader_version   1
+```
+
+SOC symmetry export additionally requires `nspin 4`, `lspinorb 1`, and
+`symmetry 1`. Magnetic systems must define the intended magnetic moments in
+`STRU`.
+
+## `stru_out.txt` symmetry block
+
+After the existing lattice, atom, and k-point records, ABACUS writes the
+spatial operation table when symmetry is enabled:
+
+```text
+N row
+<9 integer rotation entries and 3 fractional translations, repeated N times>
+```
+
+The same 12-field spatial rows are used for ordinary space groups and magnetic
+groups. A magnetic group's antiunitary spatial operations are written after its
+unitary operations. For `nspin=4`, ABACUS appends:
+
+```text
+spin_symmetry <grey_group> 2
+<antiunitary flag for each spatial operation>
+```
+
+`grey_group` is `1` for a nonmagnetic SOC calculation and `0` for a magnetic
+group. The `2` indicates that LibRPA reconstructs the spin rotation from the
+spatial operation; ABACUS does not export the old `symrot_*` sidecar files as
+part of this reader-v1 interface.

docs/advanced/interface/index.rst

@@ -20,3 +20,4 @@ Interfaces to Other Softwares
    candela
    TB2J
    migration-guide-csr-format
+   LibRPA

source/source_lcao/module_ri/librpa_stru_symmetry.h

@@ -0,0 +1,79 @@
+#ifndef LIBRPA_STRU_SYMMETRY_H
+#define LIBRPA_STRU_SYMMETRY_H
+
+#include "source_base/matrix3.h"
+#include "source_base/vector3.h"
+#include "source_cell/module_symmetry/symmetry.h"
+
+#include <cmath>
+#include <iomanip>
+#include <ostream>
+#include <stdexcept>
+#include <string>
+
+namespace RpaLriDetail
+{
+inline int checked_near_int(const double value, const std::string& context)
+{
+    const double rounded = std::round(value);
+    if (std::abs(value - rounded) > 1e-8)
+    {
+        throw std::runtime_error(context + " is not close to an integer.");
+    }
+    return static_cast<int>(rounded);
+}
+
+// LibRPA derives the SOC spin rotation from each spatial operation. ABACUS
+// therefore writes one spatial table for both ordinary and magnetic groups,
+// followed by the antiunitary flags required to reconstruct that action.
+inline void write_librpa_symmetry_block(std::ostream& output, const ModuleSymmetry::Symmetry& symmetry, const int nspin)
+{
+    if (symmetry.nrotk <= 0)
+    {
+        return;
+    }
+
+    const int n_anti = symmetry.magnetic_nspin4 ? symmetry.nrotk_anti : 0;
+    const auto write_operation
+        = [&output](const ModuleBase::Matrix3& rotation, const ModuleBase::Vector3<double>& translation) {
+              output << std::setw(4) << checked_near_int(rotation.e11, "symmetry rotation e11") << std::setw(4)
+                     << checked_near_int(rotation.e12, "symmetry rotation e12") << std::setw(4)
+                     << checked_near_int(rotation.e13, "symmetry rotation e13") << std::setw(4)
+                     << checked_near_int(rotation.e21, "symmetry rotation e21") << std::setw(4)
+                     << checked_near_int(rotation.e22, "symmetry rotation e22") << std::setw(4)
+                     << checked_near_int(rotation.e23, "symmetry rotation e23") << std::setw(4)
+                     << checked_near_int(rotation.e31, "symmetry rotation e31") << std::setw(4)
+                     << checked_near_int(rotation.e32, "symmetry rotation e32") << std::setw(4)
+                     << checked_near_int(rotation.e33, "symmetry rotation e33") << std::setw(24) << std::scientific
+                     << std::setprecision(15) << translation.x << std::setw(24) << std::scientific
+                     << std::setprecision(15) << translation.y << std::setw(24) << std::scientific
+                     << std::setprecision(15) << translation.z << std::endl;
+          };
+
+    output << (symmetry.nrotk + n_anti) << " row" << std::endl;
+    for (int isym = 0; isym < symmetry.nrotk; ++isym)
+    {
+        write_operation(symmetry.gmatrix[isym], symmetry.gtrans[isym]);
+    }
+    for (int isym = 0; isym < n_anti; ++isym)
+    {
+        write_operation(symmetry.gmatrix_anti[isym], symmetry.gtrans_anti[isym]);
+    }
+
+    if (nspin == 4)
+    {
+        // source=2 means LibRPA reconstructs U_s from det(Q) Q.
+        output << "spin_symmetry " << (symmetry.magnetic_nspin4 ? 0 : 1) << " 2" << std::endl;
+        for (int isym = 0; isym < symmetry.nrotk; ++isym)
+        {
+            output << 0 << std::endl;
+        }
+        for (int isym = 0; isym < n_anti; ++isym)
+        {
+            output << 1 << std::endl;
+        }
+    }
+}
+} // namespace RpaLriDetail
+
+#endif

source/source_lcao/module_ri/rpa_lri.hpp

@@ -24,6 +24,7 @@
 
 #include "rpa_lri.h"
 #include "exx_lri.h"
+#include "librpa_stru_symmetry.h"
 #include "source_basis/module_ao/elem_basis_idx_orb.h"
 #include "source_base/global_function.h"
 #include "source_estate/elecstate_lcao.h"
@@ -2720,6 +2721,10 @@ void RPA_LRI<T, Tdata>::out_struc(const UnitCell& ucell)
             ofs << (ik + 1) << std::endl;
         }
     }
+    if (ModuleSymmetry::Symmetry::symm_flag == 1 && ucell.symm.nrotk > 0)
+    {
+        RpaLriDetail::write_librpa_symmetry_block(ofs, ucell.symm, this->runtime.input.nspin);
+    }
     ofs.close();
     return;
 }

source/source_lcao/module_ri/test/CMakeLists.txt

@@ -16,4 +16,9 @@ AddTest(
   LIBS parameter
   SOURCES abfs_vec3_order_test.cpp
 )
+AddTest(
+  TARGET MODULE_RI_librpa_stru_symmetry_test
+  LIBS base symmetry device
+  SOURCES librpa_stru_symmetry_test.cpp
+)
 install(DIRECTORY support DESTINATION ${CMAKE_CURRENT_BINARY_DIR})

source/source_lcao/module_ri/test/librpa_stru_symmetry_test.cpp

@@ -0,0 +1,107 @@
+#include "source_lcao/module_ri/librpa_stru_symmetry.h"
+
+#include <cmath>
+#include <cstdlib>
+#include <iostream>
+#include <sstream>
+#include <string>
+#include <vector>
+
+namespace
+{
+void require(const bool condition, const char* message)
+{
+    if (!condition)
+    {
+        std::cerr << "FAILED: " << message << '\n';
+        std::exit(1);
+    }
+}
+
+std::vector<std::string> split_lines(const std::string& text)
+{
+    std::vector<std::string> lines;
+    std::istringstream input(text);
+    std::string line;
+    while (std::getline(input, line))
+    {
+        lines.push_back(line);
+    }
+    return lines;
+}
+
+void require_spatial_row(const std::string& line)
+{
+    std::istringstream input(line);
+    std::vector<std::string> fields;
+    std::string field;
+    while (input >> field)
+    {
+        fields.push_back(field);
+    }
+    require(fields.size() == 12, "a symmetry row must contain 9 integer and 3 translation fields");
+    for (int i = 0; i < 9; ++i)
+    {
+        std::size_t consumed = 0;
+        (void)std::stoi(fields[static_cast<std::size_t>(i)], &consumed);
+        require(consumed == fields[static_cast<std::size_t>(i)].size(), "rotation fields must be integers");
+    }
+    for (int i = 9; i < 12; ++i)
+    {
+        std::size_t consumed = 0;
+        const double value = std::stod(fields[static_cast<std::size_t>(i)], &consumed);
+        require(consumed == fields[static_cast<std::size_t>(i)].size() && std::isfinite(value),
+                "translation fields must be finite floating-point values");
+    }
+}
+
+void test_magnetic_block()
+{
+    ModuleSymmetry::Symmetry symmetry;
+    symmetry.nrotk = 2;
+    symmetry.nrotk_anti = 1;
+    symmetry.magnetic_nspin4 = true;
+    symmetry.gmatrix[0] = ModuleBase::Matrix3();
+    symmetry.gmatrix[1] = ModuleBase::Matrix3(0.0, -1.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 1.0);
+    symmetry.gtrans[1] = ModuleBase::Vector3<double>(0.5, 0.0, 0.0);
+    symmetry.gmatrix_anti[0] = ModuleBase::Matrix3(0.0, 1.0, 0.0, -1.0, 0.0, 0.0, 0.0, 0.0, 1.0);
+    symmetry.gtrans_anti[0] = ModuleBase::Vector3<double>(0.125, 0.25, 0.375);
+
+    std::ostringstream output;
+    RpaLriDetail::write_librpa_symmetry_block(output, symmetry, 4);
+    const auto lines = split_lines(output.str());
+    require(lines.size() == 8, "magnetic SOC block must contain 3 spatial and 3 spin rows");
+    require(lines[0] == "3 row", "unitary and antiunitary operations share one spatial block");
+    require_spatial_row(lines[1]);
+    require_spatial_row(lines[2]);
+    require_spatial_row(lines[3]);
+    require(lines[3].find("1.250000000000000e-01") != std::string::npos,
+            "antiunitary spatial operation must follow unitary operations");
+    require(lines[4] == "spin_symmetry 0 2", "LibRPA must reconstruct the magnetic spin action");
+    require(lines[5] == "0" && lines[6] == "0" && lines[7] == "1", "only the final operation is antiunitary");
+}
+
+void test_nonmagnetic_soc_block()
+{
+    ModuleSymmetry::Symmetry symmetry;
+    symmetry.nrotk = 1;
+    symmetry.magnetic_nspin4 = false;
+    symmetry.gmatrix[0] = ModuleBase::Matrix3();
+
+    std::ostringstream output;
+    RpaLriDetail::write_librpa_symmetry_block(output, symmetry, 4);
+    const auto lines = split_lines(output.str());
+    require(lines.size() == 4, "nonmagnetic SOC block must contain one spatial and one spin row");
+    require(lines[0] == "1 row", "nonmagnetic block must retain its spatial operation");
+    require(lines[2] == "spin_symmetry 1 2", "grey-group SOC metadata must be preserved");
+    require(lines[3] == "0", "the nonmagnetic spatial operation must be unitary");
+}
+} // namespace
+
+int main()
+{
+    test_magnetic_block();
+    test_nonmagnetic_soc_block();
+    std::cout << "LibRPA stru_out symmetry tests passed\n";
+    return 0;
+}

tests/08_EXX/63_LIBRPA_SOC_SYMMETRY/INPUT

@@ -0,0 +1,30 @@
+INPUT_PARAMETERS
+suffix            autotest
+pseudo_dir        ../../PP_ORB
+orbital_dir       ../../PP_ORB
+calculation       scf
+ntype             1
+nbands            4
+gamma_only        0
+ecutwfc           20
+lspinorb          1
+nspin             4
+smearing_method   gaussian
+smearing_sigma    1e-4
+basis_type        lcao
+ks_solver         genelpa
+dft_functional    pbe
+rpa               1
+rpa_out_vel       1
+out_mat_xc        1
+out_wfc_lcao      1
+out_librpa_reader_version 1
+symmetry          1
+exx_separate_loop 1
+exx_pca_threshold 1e-3
+shrink_abfs_pca_thr 1e-1
+shrink_lu_inv_thr 1e-3
+exx_singularity_correction massidda
+exx_ccp_rmesh_times 1
+rpa_ccp_rmesh_times 1
+exx_real_number   1

tests/08_EXX/63_LIBRPA_SOC_SYMMETRY/KPT

@@ -0,0 +1,4 @@
+K_POINTS
+0
+Gamma
+1 1 2 0 0 0

tests/08_EXX/63_LIBRPA_SOC_SYMMETRY/README

@@ -0,0 +1,4 @@
+Check reader-v1 RPA producer output, SOC GW-preprocessor inputs, and serialized
+spatial and spin symmetry tables for a compact full-relativistic two-atom H
+calculation. This is a serialization contract, not a symmetry-on/off
+energy-equivalence comparison.

tests/08_EXX/63_LIBRPA_SOC_SYMMETRY/STRU

@@ -0,0 +1,22 @@
+ATOMIC_SPECIES
+H 1.008 H_ONCV_PBE_FR-1.0.upf
+
+NUMERICAL_ORBITAL
+H_gga_6au_60Ry_1s.orb
+
+LATTICE_CONSTANT
+1.8897261254578281
+
+LATTICE_VECTORS
+0.000000 1.500000 1.500000
+1.500000 0.000000 1.500000
+1.500000 1.500000 0.000000
+
+ATOMIC_POSITIONS
+Direct
+
+H
+0.0
+2
+0.000000 0.000000 0.000000 0 0 0
+0.500000 0.500000 0.500000 0 0 0

tests/08_EXX/63_LIBRPA_SOC_SYMMETRY/librpa_producer_manifest.json

@@ -0,0 +1,20 @@
+{
+  "output_dir": "OUT.librpa",
+  "reference_dir": "producer_reference",
+  "required": [
+    {"pattern": "stru_out.txt", "kind": "stru", "symmetry_rows": 1, "spin_symmetry_source": 2, "spin_symmetry_grey_group": 1, "reference": true, "abs_tol": 1e-7},
+    {"pattern": "band_out.txt", "kind": "band", "reference": true, "abs_tol": 1e-7},
+    {"pattern": "basis_*", "kind": "text", "reference": true, "abs_tol": 1e-7},
+    {"pattern": "bz_sampling_out", "kind": "text", "reference": true, "abs_tol": 1e-7},
+    {"pattern": "v1_Cs_shrinked_data_*.txt", "kind": "lri_coeff_v1", "reference": true, "abs_tol": 1e-7},
+    {"pattern": "v1_Cs_data_*.txt", "kind": "lri_coeff_v1", "reference": true, "abs_tol": 1e-7},
+    {"pattern": "v1_shrink_sinvS_*.txt", "kind": "shrink_sinvs_v1", "reference": true, "abs_tol": 1e-7},
+    {"pattern": "v1_coulomb_cut_iq_*.dat", "kind": "coulomb_v1", "reference": true, "abs_tol": 1e-7},
+    {"pattern": "v1_coulomb_full_iq_*.dat", "kind": "coulomb_v1", "reference": true, "abs_tol": 1e-7},
+    {"pattern": "KS_eigenvector_*.dat", "kind": "gauge"},
+    {"pattern": "velocity_matrix", "kind": "gauge"},
+    {"base_dir": ".", "pattern": "OUT.autotest/vxc_out.dat", "kind": "text", "reference": true, "abs_tol": 1e-7},
+    {"base_dir": ".", "pattern": "OUT.autotest/KPT.info", "kind": "text", "reference": true, "abs_tol": 1e-7},
+    {"base_dir": ".", "pattern": "OUT.autotest/wfk*s4_nao.txt", "kind": "wfc_nao", "expected_count": 2}
+  ]
+}

tests/08_EXX/63_LIBRPA_SOC_SYMMETRY/producer_reference/OUT.autotest/KPT.info

@@ -0,0 +1,10 @@
+                               nkstot now = 2
+ K-POINTS DIRECT COORDINATES
+ KPOINTS    DIRECT_X    DIRECT_Y    DIRECT_Z  WEIGHT
+       1  0.00000000  0.00000000  0.00000000  0.5000
+       2  0.00000000  0.00000000  0.50000000  0.5000
+                                   nkstot = 2                                                            ibzkpt
+K-POINTS REDUCTION ACCORDING TO SYMMETRY
+     KPT    DIRECT_X    DIRECT_Y    DIRECT_Z     IBZ    DIRECT_X    DIRECT_Y    DIRECT_Z
+       1  0.00000000  0.00000000  0.00000000       1  0.00000000  0.00000000  0.00000000
+       2  0.00000000  0.00000000  0.50000000       2  0.00000000  0.00000000  0.50000000

tests/08_EXX/63_LIBRPA_SOC_SYMMETRY/producer_reference/OUT.autotest/vxc_out.dat

@@ -0,0 +1,11 @@
+2
+1
+4
+-4.0274539288505579e-01	-1.0959264372970836e+01
+-4.0274539288505579e-01	-1.0959264372970836e+01
+-5.2769536070469791e-01	-1.4359327427498375e+01
+-5.2769536070469791e-01	-1.4359327427498375e+01
+-4.4986613169200201e-01	-1.2241485456459218e+01
+-4.4986613169200168e-01	-1.2241485456459207e+01
+-4.4986613138922499e-01	-1.2241485448220232e+01
+-4.4986613138922515e-01	-1.2241485448220237e+01

tests/08_EXX/63_LIBRPA_SOC_SYMMETRY/producer_reference/band_out.txt

@@ -0,0 +1,15 @@
+2
+1
+4
+4
+0.023189820575926
+     1     1
+    1   1.000000000000000       -0.443892166842433      -12.078925533247524
+    2   1.000000000000000       -0.443892166842433      -12.078925533247524
+    3   0.000000000000000        0.783096439472231       21.309147320668909
+    4   0.000000000000000        0.783096439472231       21.309147320668909
+     2     1
+    1   0.500000003527986        0.023189820575614        0.631027390851974
+    2   0.500000003527986        0.023189820575614        0.631027390851974
+    3   0.499999996427854        0.023189820576243        0.631027390869097
+    4   0.499999996427854        0.023189820576243        0.631027390869097

tests/08_EXX/63_LIBRPA_SOC_SYMMETRY/producer_reference/basis_aux_out

@@ -0,0 +1,4 @@
+         1         2    abacus
+         1         1
+         1         1
+         0

tests/08_EXX/63_LIBRPA_SOC_SYMMETRY/producer_reference/basis_aux_shrink_out

@@ -0,0 +1,4 @@
+         1         2    abacus
+         1         1
+         1         1
+         0