From 753fd7a192d5e4453f03df96b1a63bc9d8efad80 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Stefanos=20Carlstr=C3=B6m?= <stefanos.carlstrom@gmail.com>
Date: Tue, 25 Oct 2022 10:34:29 +0200
Subject: [PATCH 1/2] Calculus of variations for OrbitalRadialOverlap

---
 src/tensors.jl  | 18 ++++++++++++++++++
 test/tensors.jl | 10 +++++++++-
 2 files changed, 27 insertions(+), 1 deletion(-)

diff --git a/src/tensors.jl b/src/tensors.jl
index 4931db4..4b300b8 100644
--- a/src/tensors.jl
+++ b/src/tensors.jl
@@ -350,6 +350,24 @@ EnergyExpressions.isdependent(o::OrbitalRadialOverlap, orb::O) where O = o.b ==
 Base.show(io::IO, o::OrbitalRadialOverlap) =
     write(io, "⟨$(o.a)|$(o.b)⟩ᵣ")
 
+"""
+    diff(ab::OrbitalRadialOverlap, o::O)
+
+Vary the orbital radial overlap ⟨a|b⟩ᵣ with respect to |o⟩.
+"""
+Base.diff(ab::OrbitalRadialOverlap, o::O) where O =
+    o == ab.b ? NBodyEquation(Conjugate(ab.a),
+                              IdentityOperator{1}()) : 0
+
+"""
+    diff(ab::OrbitalRadialOverlap, o::Conjugate{O})
+
+Vary the orbital radial overlap ⟨a|b⟩ᵣ with respect to ⟨o|.
+"""
+Base.diff(ab::OrbitalRadialOverlap, o::Conjugate{O}) where O =
+    conj(o) == ab.a ? NBodyEquation(ab.b,
+                                    IdentityOperator{1}()) : 0
+
 # * Tensor operators
 
 """
diff --git a/test/tensors.jl b/test/tensors.jl
index 84e7a2d..4b46e2f 100644
--- a/test/tensors.jl
+++ b/test/tensors.jl
@@ -2,7 +2,7 @@
     import AngularMomentumAlgebra: components, component,
     RadialOperator, RadialGradientOperator,
     OrbitalRadialOverlap, OrbitalRadialMatrixElement, radial_integral, integrate_spinor
-    import EnergyExpressions: NBodyTerm, NBodyMatrixElement, OrbitalMatrixElement
+    import EnergyExpressions: NBodyTerm, NBodyMatrixElement, OrbitalMatrixElement, NBodyEquation
     @test system(Tensor) == FullSystem()
 
     𝐂⁵ = SphericalTensor(5)
@@ -59,6 +59,14 @@
         @test isdependent(oro, so"1s₀β")
 
         @test string(oro) == "⟨1s₀α|1s₀β⟩ᵣ"
+
+        eq(orb, op, args...) = NBodyEquation(orb, op, args...)
+        I₁ = IdentityOperator{1}()
+
+        @test diff(oro, conj(so"1s₀α")) == eq(so"1s₀β", I₁)
+        @test iszero(diff(oro, conj(so"1s₀β")))
+        @test iszero(diff(oro, so"1s₀α"))
+        @test diff(oro, so"1s₀β") == eq(conj(so"1s₀α"), I₁)
     end
 
     @testset "TensorOperator" begin

From 0343e3bbf8aa784dad8a4778d3dc80a20293d7a9 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Stefanos=20Carlstr=C3=B6m?= <stefanos.carlstrom@gmail.com>
Date: Tue, 25 Oct 2022 10:35:06 +0200
Subject: [PATCH 2/2] Bumped version

---
 Project.toml | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/Project.toml b/Project.toml
index 5fe9b56..4bdbd60 100644
--- a/Project.toml
+++ b/Project.toml
@@ -1,7 +1,7 @@
 name = "AngularMomentumAlgebra"
 uuid = "0bb1b1ac-0216-11e9-361b-bb110c659ce1"
 authors = ["Stefanos Carlström <stefanos.carlstrom@gmail.com>"]
-version = "0.3.0"
+version = "0.3.1"
 
 [deps]
 AtomicLevels = "10933b4c-d60f-11e8-1fc6-bd9035a249a1"