summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
authorJaron Kent-Dobias <jaron@kent-dobias.com>2021-11-08 23:05:15 +0100
committerJaron Kent-Dobias <jaron@kent-dobias.com>2021-11-08 23:05:15 +0100
commit6b42ce01cd289567a6109d831fe86a0667d18f37 (patch)
tree22a48f0136b060769bd475642cce0d0c1219a9d6
parent897df7986e9cafcfa0c5f60a424d5a452dceca07 (diff)
downloadcode-6b42ce01cd289567a6109d831fe86a0667d18f37.tar.gz
code-6b42ce01cd289567a6109d831fe86a0667d18f37.tar.bz2
code-6b42ce01cd289567a6109d831fe86a0667d18f37.zip
Working minimization.
-rw-r--r--dynamics.hpp8
-rw-r--r--langevin.cpp7
2 files changed, 9 insertions, 6 deletions
diff --git a/dynamics.hpp b/dynamics.hpp
index 639803b..21afd30 100644
--- a/dynamics.hpp
+++ b/dynamics.hpp
@@ -55,13 +55,13 @@ Vector<Real> findMinimum(const Tensor<Real, p>& J, const Vector<Real>& z0, Real
Matrix<Real> M = ddH - (dH * z.transpose() + z.dot(dH) * Matrix<Real>::Identity(z.size(), z.size()) + (ddH * z) * z.transpose()) / (Real)z.size() + 2.0 * z * z.transpose();
while (g.norm() / z.size() > ε && λ < 1e8) {
- Vector<Real> dz = (M + λ * Matrix<Real>::Identity(z.size(), z.size())).partialPivLu().solve(g);
+ Vector<Real> dz = (M + λ * (Matrix<Real>)abs(M.diagonal().array()).matrix().asDiagonal()).partialPivLu().solve(g);
dz -= z.dot(dz) * z / (Real)z.size();
Vector<Real> zNew = normalize(z - dz);
auto [HNew, dHNew, ddHNew] = hamGradHess(J, zNew);
- if (HNew <= H * 1.01) {
+ if (HNew * 1.0001 <= H) {
z = zNew;
H = HNew;
dH = dHNew;
@@ -70,12 +70,10 @@ Vector<Real> findMinimum(const Tensor<Real, p>& J, const Vector<Real>& z0, Real
g = dH - z.dot(dH) * z / (Real)z.size();
M = ddH - (dH * z.transpose() + z.dot(dH) * Matrix<Real>::Identity(z.size(), z.size()) + (ddH * z) * z.transpose()) / (Real)z.size() + 2.0 * z * z.transpose();
- λ /= 1.001;
+ λ /= 2;
} else {
λ *= 1.5;
}
-
- std::cout << "error : " << H << " " << g.norm() / z.size() << " " << λ << std::endl;
}
return z;
diff --git a/langevin.cpp b/langevin.cpp
index a775774..1181173 100644
--- a/langevin.cpp
+++ b/langevin.cpp
@@ -146,8 +146,13 @@ int main(int argc, char* argv[]) {
Vector<Real> zMin = randomMinimum(ReJ, Red, r, ε);
auto [Hr, dHr, ddHr] = hamGradHess(ReJ, zMin);
- Eigen::EigenSolver<Matrix<Real>> eigenS(ddHr - ((ddHr * zMin) * zMin.transpose()) / (Real)zMin.size());
+ Eigen::EigenSolver<Matrix<Real>> eigenS(ddHr - (dHr * zMin.transpose() + zMin.dot(dHr) * Matrix<Real>::Identity(zMin.size(), zMin.size()) + (ddHr * zMin) * zMin.transpose()) / (Real)zMin.size() + 2.0 * zMin * zMin.transpose());
std::cout << eigenS.eigenvalues().transpose() << std::endl;
+ for (unsigned i = 0; i < N; i++) {
+ Vector<Real> zNew = normalize(zMin + 0.01 * eigenS.eigenvectors().col(i).real());
+ std::cout << getHamiltonian(ReJ, zNew) - Hr << " " << real(eigenS.eigenvectors().col(i).dot(zMin)) << std::endl;
+ }
+ std::cout << std::endl;
getchar();
complex_normal_distribution<Real> d(0, 1, 0);