From 02a7639237d8dda1f4ad4c8d274bfdb799efe72b Mon Sep 17 00:00:00 2001
From: Lisa Julia Nebel <lisa_julia.nebel@tu-dresden.de>
Date: Thu, 3 Sep 2020 11:18:08 +0200
Subject: [PATCH] Add Proximal-Newton solver to film-on-substrate
---
src/film-on-substrate.cc | 68 ++++++++++++++++++++++--------------
src/film-on-substrate.parset | 23 ++++++++++--
2 files changed, 62 insertions(+), 29 deletions(-)
diff --git a/src/film-on-substrate.cc b/src/film-on-substrate.cc
index 3d3fa4d8..6348e753 100644
--- a/src/film-on-substrate.cc
+++ b/src/film-on-substrate.cc
@@ -60,6 +60,7 @@
#include <dune/gfe/localgeodesicfeadolcstiffness.hh>
#include <dune/gfe/cosseratvtkwriter.hh>
#include <dune/gfe/geodesicfeassembler.hh>
+#include <dune/gfe/riemannianpnsolver.hh>
#include <dune/gfe/riemanniantrsolver.hh>
#include <dune/gfe/vertexnormals.hh>
#include <dune/gfe/surfacecosseratenergy.hh>
@@ -152,8 +153,9 @@ int main (int argc, char *argv[]) try
int numLevels = parameterSet.get<int>("numLevels");
int numHomotopySteps = parameterSet.get<int>("numHomotopySteps");
const double tolerance = parameterSet.get<double>("tolerance");
- const int maxTrustRegionSteps = parameterSet.get<int>("maxTrustRegionSteps");
+ const int maxSolverSteps = parameterSet.get<int>("maxSolverSteps");
const double initialTrustRegionRadius = parameterSet.get<double>("initialTrustRegionRadius");
+ const double initialRegularization = parameterSet.get<double>("initialRegularization");
const int multigridIterations = parameterSet.get<int>("numIt");
const int nu1 = parameterSet.get<int>("nu1");
const int nu2 = parameterSet.get<int>("nu2");
@@ -577,26 +579,6 @@ int main (int argc, char *argv[]) try
GeodesicFEAssembler<FEBasis,TargetSpace> assembler(gridView, &localGFEADOLCStiffness);
- // /////////////////////////////////////////////////
- // Create a Riemannian trust-region solver
- // /////////////////////////////////////////////////
-
- RiemannianTrustRegionSolver<FEBasis,TargetSpace> solver;
- solver.setup(*grid,
- &assembler,
- x,
- dirichletDofs,
- tolerance,
- maxTrustRegionSteps,
- initialTrustRegionRadius,
- multigridIterations,
- mgTolerance,
- mu, nu1, nu2,
- baseIterations,
- baseTolerance,
- instrumented);
-
- solver.setScaling(parameterSet.get<FieldVector<double,6> >("trustRegionScaling"));
////////////////////////////////////////////////////////
// Set Dirichlet values
@@ -625,16 +607,50 @@ int main (int argc, char *argv[]) try
x[j].r = ddV[j];
x[j].q.set(dOV[j]);
}
+ // /////////////////////////////////////////////////
+ // Create the solver and solve
+ // /////////////////////////////////////////////////
+ if (parameterSet.get<std::string>("solvertype") == "multigrid") {
+ RiemannianTrustRegionSolver<FEBasis,TargetSpace> solver;
+ solver.setup(*grid,
+ &assembler,
+ x,
+ dirichletDofs,
+ tolerance,
+ maxSolverSteps,
+ initialTrustRegionRadius,
+ multigridIterations,
+ mgTolerance,
+ mu, nu1, nu2,
+ baseIterations,
+ baseTolerance,
+ instrumented);
+
+ solver.setScaling(parameterSet.get<FieldVector<double,6> >("trustRegionScaling"));
+ solver.setInitialIterate(x);
+ solver.solve();
+
+ x = solver.getSol();
+ } else { //parameterSet.get<std::string>("solvertype") == "cholmod"
+ RiemannianProximalNewtonSolver<FEBasis,TargetSpace> solver;
+ solver.setup(*grid,
+ &assembler,
+ x,
+ dirichletDofs,
+ tolerance,
+ maxSolverSteps,
+ initialRegularization,
+ instrumented);
+ solver.setInitialIterate(x);
+ solver.solve();
+
+ x = solver.getSol();
+ }
// /////////////////////////////////////////////////////
// Solve!
// /////////////////////////////////////////////////////
- solver.setInitialIterate(x);
- solver.solve();
-
- x = solver.getSol();
-
std::cout << "Overall calculation took " << overallTimer.elapsed() << " sec." << std::endl;
/////////////////////////////////
diff --git a/src/film-on-substrate.parset b/src/film-on-substrate.parset
index 645c408b..21104183 100644
--- a/src/film-on-substrate.parset
+++ b/src/film-on-substrate.parset
@@ -50,14 +50,31 @@ initialDeformation = "x"
# Solver parameters
#############################################
+# Inner solver, cholmod or multigrid
+solvertype = cholmod
+
# Number of homotopy steps for the Dirichlet boundary conditions
numHomotopySteps = 1
-# Tolerance of the trust region solver
+# Tolerance of the solver
tolerance = 1e-3
-# Max number of steps of the trust region solver
-maxTrustRegionSteps = 20
+# Max number of solver steps
+maxSolverSteps = 300
+
+# Measure convergence
+instrumented = 0
+
+#############################################
+# Solver parameters specific for proximal newton solver using cholmod
+#############################################
+
+# initial regularization parameter
+initialRegularization = 10000000000
+
+#############################################
+# Solver parameters specific for trust-region solver using multigrid solver
+#############################################
trustRegionScaling = 1 1 1 0.01 0.01 0.01
--
GitLab