Feature: DFT+U for noncollinear spin calculation with PW base

Author: dyzheng

docs/advanced/input_files/input-main.md

@@ -1235,6 +1235,12 @@ Note: In new angle mixing, you should set `mixing_beta_mag >> mixing_beta`. The
 - **Description**: To determine the number of old iterations' `drho` used in slope calculations.
 - **Default**: `mixing_ndim`
 
+### sc_os_ndim
+
+- **Type**: int
+- **Description**: To determine the number of old iterations to judge oscillation, it occured,  more accurate lambda with DeltaSpin method would be calculated, only for PW base.
+- **Default**: 5
+
 ### chg_extrap
 
 - **Type**: String

source/CMakeLists.txt

@@ -47,6 +47,7 @@ list(APPEND device_srcs
 
   module_hamilt_pw/hamilt_pwdft/kernels/force_op.cpp
   module_hamilt_pw/hamilt_pwdft/kernels/stress_op.cpp
+  module_hamilt_pw/hamilt_pwdft/kernels/onsite_op.cpp
   module_hamilt_pw/hamilt_pwdft/kernels/wf_op.cpp
   module_hamilt_pw/hamilt_pwdft/kernels/vnl_op.cpp
   module_base/kernels/math_op.cpp
@@ -60,6 +61,7 @@ if(USE_CUDA)
     module_hamilt_pw/hamilt_pwdft/kernels/cuda/ekinetic_op.cu
     module_hamilt_pw/hamilt_pwdft/kernels/cuda/meta_op.cu
     module_hamilt_pw/hamilt_stodft/kernels/cuda/hpsi_norm_op.cu
+    module_hamilt_pw/hamilt_pwdft/kernels/cuda/onsite_op.cu
     module_basis/module_pw/kernels/cuda/pw_op.cu
     module_hsolver/kernels/cuda/dngvd_op.cu
     module_hsolver/kernels/cuda/math_kernel_op.cu
@@ -83,6 +85,7 @@ if(USE_ROCM)
     module_hamilt_pw/hamilt_pwdft/kernels/rocm/veff_op.hip.cu
     module_hamilt_pw/hamilt_pwdft/kernels/rocm/ekinetic_op.hip.cu
     module_hamilt_pw/hamilt_pwdft/kernels/rocm/meta_op.hip.cu
+    module_hamilt_pw/hamilt_pwdft/kernels/rocm/onsite_op.hip.cu
     module_hamilt_pw/hamilt_stodft/kernels/rocm/hpsi_norm_op.hip.cu
     module_basis/module_pw/kernels/rocm/pw_op.hip.cu
     module_hsolver/kernels/rocm/dngvd_op.hip.cu

source/Makefile.Objects

@@ -634,10 +634,13 @@ OBJS_SRCPW=H_Ewald_pw.o\
     forces_nl.o\
     forces_cc.o\
     forces_scc.o\
+    forces_onsite.o\
+    onsite_proj_pw.o\
     fs_nonlocal_tools.o\
     fs_kin_tools.o\
     force_op.o\
     stress_op.o\
+    onsite_op.o\
     wf_op.o\
     vnl_op.o\
     global.o\
@@ -662,6 +665,7 @@ OBJS_SRCPW=H_Ewald_pw.o\
     stress_func_loc.o\
     stress_func_nl.o\
     stress_func_us.o\
+    stress_func_onsite.o\
     stress_pw.o\
     of_stress_pw.o\
     symmetry_rho.o\
@@ -672,7 +676,9 @@ OBJS_SRCPW=H_Ewald_pw.o\
     elecond.o\
     sto_tool.o\
     sto_elecond.o\
-    sto_dos.o
+    sto_dos.o\
+    onsite_projector.o\
+    onsite_proj_tools.o
 
 OBJS_VDW=vdw.o\
     vdwd2_parameters.o\
@@ -690,7 +696,8 @@ OBJS_DFTU=dftu.o\
       dftu_io.o\
       dftu_tools.o\
       dftu_occup.o\
-      dftu_hamilt.o
+      dftu_hamilt.o\
+      dftu_pw.o
 
 OBJS_DELTASPIN=basic_funcs.o\
       cal_mw_from_lambda.o\

source/module_base/CMakeLists.txt

@@ -54,6 +54,7 @@ add_library(
     spherical_bessel_transformer.cpp
     cubic_spline.cpp
     parallel_2d.cpp
+    projgen.cpp
     module_mixing/mixing_data.cpp
     module_mixing/mixing.cpp
     module_mixing/plain_mixing.cpp

source/module_basis/module_nao/projgen.cpp renamed to source/module_base/projgen.cpp

@@ -14,7 +14,6 @@ if(ENABLE_LCAO)
     two_center_table.cpp
     two_center_integrator.cpp
     two_center_bundle.cpp
-    projgen.cpp
   )
 
   if(ENABLE_COVERAGE)

source/module_basis/module_nao/projgen.h renamed to source/module_base/projgen.h

@@ -7,7 +7,7 @@
 // FIXME: should update with pyabacus
 // #include "module_io/orb_io.h"
 
-#include "projgen.h"
+#include "module_base/projgen.h"
 
 #include <fstream>
 #include <iostream>

source/module_basis/module_nao/CMakeLists.txt

@@ -14,7 +14,6 @@ AddTest(
     ../atomic_radials.cpp
     ../radial_set.cpp
     ../numerical_radial.cpp
-    ../projgen.cpp
     ../../module_ao/ORB_atomic_lm.cpp
     ../../module_ao/ORB_atomic.cpp
     ../../../module_io/orb_io.cpp
@@ -84,7 +83,6 @@ AddTest(
     ../pswfc_radials.cpp
     ../radial_set.cpp
     ../numerical_radial.cpp
-    ../projgen.cpp
     ../sphbes_radials.cpp
     ../../module_ao/ORB_atomic_lm.cpp
     ../../module_ao/ORB_atomic.cpp
@@ -104,7 +102,6 @@ AddTest(
     ../pswfc_radials.cpp
     ../sphbes_radials.cpp
     ../radial_set.cpp
-    ../projgen.cpp
     ../numerical_radial.cpp
     ../two_center_bundle.cpp
     ../two_center_integrator.cpp
@@ -131,7 +128,6 @@ AddTest(
     ../real_gaunt_table.cpp
     ../radial_collection.cpp
     ../atomic_radials.cpp
-    ../projgen.cpp
     ../beta_radials.cpp
     ../hydrogen_radials.cpp
     ../pswfc_radials.cpp
@@ -158,7 +154,6 @@ AddTest(
     ../pswfc_radials.cpp
     ../sphbes_radials.cpp
     ../radial_set.cpp
-    ../projgen.cpp
     ../numerical_radial.cpp
     ../../../module_io/orb_io.cpp
   LIBS parameter ${math_libs} device base container orb 
@@ -179,7 +174,6 @@ AddTest(
     ../pswfc_radials.cpp
     ../sphbes_radials.cpp
     ../radial_set.cpp
-    ../projgen.cpp
     ../numerical_radial.cpp
     ../../../module_io/orb_io.cpp
   LIBS parameter ${math_libs} device base container orb 

source/module_basis/module_nao/atomic_radials.cpp

@@ -1,4 +1,4 @@
-#include "module_basis/module_nao/projgen.h"
+#include "module_base/projgen.h"
 #include "gtest/gtest.h"
 
 #include "module_base/math_integral.h"