From 1a99f2ceda091a2d74fad227be9f6adee5fb35af Mon Sep 17 00:00:00 2001
From: Oliver Sander <sander@igpm.rwth-aachen.de>
Date: Mon, 4 May 2015 06:10:30 +0200
Subject: [PATCH] Use the BSplineBasis
Unfortunately, this bases needs to be switched on and off using preprocessor
commands. Hopefully I'll find a way to make that prettier in the future.
---
src/harmonicmaps.cc | 27 +++++++++++++++++++++------
1 file changed, 21 insertions(+), 6 deletions(-)
diff --git a/src/harmonicmaps.cc b/src/harmonicmaps.cc
index 4b35c648..0632cc3a 100644
--- a/src/harmonicmaps.cc
+++ b/src/harmonicmaps.cc
@@ -21,6 +21,7 @@
#include <dune/functions/gridfunctions/discretescalarglobalbasisfunction.hh>
#include <dune/functions/functionspacebases/pq1nodalbasis.hh>
#include <dune/functions/functionspacebases/pqknodalbasis.hh>
+#include <dune/functions/functionspacebases/bsplinebasis.hh>
#include <dune/fufem/boundarypatch.hh>
#include <dune/fufem/functiontools/basisinterpolator.hh>
@@ -41,6 +42,7 @@
#include <dune/gfe/geodesicfeassembler.hh>
#include <dune/gfe/riemanniantrsolver.hh>
#include <dune/gfe/embeddedglobalgfefunction.hh>
+#include <dune/gfe/bsplineinterpolate.hh>
// grid dimension
const int dim = 2;
@@ -56,6 +58,8 @@ typedef UnitVector<double,3> TargetSpace;
// Tangent vector of the image space
const int blocksize = TargetSpace::TangentVector::dimension;
+//#define LAGRANGE
+
using namespace Dune;
@@ -85,8 +89,10 @@ int main (int argc, char *argv[]) try
ParameterTreeParser::readOptions(argc, argv, parameterSet);
- // read solver settings
+ // read problem settings
const int numLevels = parameterSet.get<int>("numLevels");
+ const int order = parameterSet.get<int>("order");
+ // read solver settings
const double tolerance = parameterSet.get<double>("tolerance");
const int maxTrustRegionSteps = parameterSet.get<int>("maxTrustRegionSteps");
const double initialTrustRegionRadius = parameterSet.get<double>("initialTrustRegionRadius");
@@ -106,13 +112,14 @@ int main (int argc, char *argv[]) try
shared_ptr<GridType> grid;
FieldVector<double,dim> lower(0), upper(1);
+ array<unsigned int,dim> elements;
if (parameterSet.get<bool>("structuredGrid")) {
lower = parameterSet.get<FieldVector<double,dim> >("lower");
upper = parameterSet.get<FieldVector<double,dim> >("upper");
- array<unsigned int,dim> elements = parameterSet.get<array<unsigned int,dim> >("elements");
+ elements = parameterSet.get<array<unsigned int,dim> >("elements");
grid = StructuredGridFactory<GridType>::createCubeGrid(lower, upper, elements);
} else {
@@ -128,11 +135,15 @@ int main (int argc, char *argv[]) try
//////////////////////////////////////////////////////////////////////////////////
// Construct the scalar function space basis corresponding to the GFE space
//////////////////////////////////////////////////////////////////////////////////
-
- typedef Dune::Functions::PQKNodalBasis<typename GridType::LeafGridView, 3> FEBasis;
-
+#ifdef LAGRANGE
+ typedef Dune::Functions::PQKNodalBasis<typename GridType::LeafGridView, 1> FEBasis;
FEBasis feBasis(grid->leafGridView());
-
+#else
+ typedef Dune::Functions::BSplineBasis<typename GridType::LeafGridView> FEBasis;
+ for (auto& i : elements)
+ i = i<<(numLevels-1);
+ FEBasis feBasis(grid->leafGridView(), lower, upper, elements, order);
+#endif
typedef DuneFunctionsBasis<FEBasis> FufemFEBasis;
FufemFEBasis fufemFeBasis(feBasis);
@@ -160,7 +171,11 @@ int main (int argc, char *argv[]) try
auto pythonInitialIterate = module.get("fdf").toC<std::shared_ptr<FBase>>();
std::vector<TargetSpace::CoordinateType> v;
+#ifdef LAGRANGE
::Functions::interpolate(fufemFeBasis, v, *pythonInitialIterate);
+#else
+ Dune::Functions::interpolate(feBasis, v, *pythonInitialIterate, lower, upper, elements, order);
+#endif
for (size_t i=0; i<x.size(); i++)
x[i] = v[i];
--
GitLab