From 514229f1fd5a16c263dc8f29d4ef42b8f4f8f237 Mon Sep 17 00:00:00 2001 From: Oliver Sander <sander@igpm.rwth-aachen.de> Date: Wed, 14 Jan 2015 15:22:13 +0000 Subject: [PATCH] Ignore more files [[Imported from SVN: r9997]] --- cosserat-continuum-twisted-strip.parset | 2 ++ ...continuum-wong-pellegrino-wrinkling.parset | 6 +++-- cosserat-continuum-wriggers-l-shape.parset | 13 +++++++---- dune/gfe/trustregionmmgbasesolver.hh | 17 +++++++++++--- finite-strain-elasticity.cc | 9 ++++---- harmonicmaps.cc | 2 ++ mixed-cosserat-continuum.cc | 6 +++-- test/Makefile.am | 12 ++++++++-- test/adolctest.cc | 23 +++++++++++++++++++ test/targetspacetest.cc | 23 +++++++++++-------- 10 files changed, 86 insertions(+), 27 deletions(-) diff --git a/cosserat-continuum-twisted-strip.parset b/cosserat-continuum-twisted-strip.parset index d2cdfd39..6accef9b 100644 --- a/cosserat-continuum-twisted-strip.parset +++ b/cosserat-continuum-twisted-strip.parset @@ -84,6 +84,8 @@ kappa = 1 # Boundary values ############################################# +problem = twisted-strip + ### Python predicate specifying all Dirichlet grid vertices # x is the vertex coordinate dirichletVerticesPredicate = "x[0] < 0.001 or x[0] > 0.0999" diff --git a/cosserat-continuum-wong-pellegrino-wrinkling.parset b/cosserat-continuum-wong-pellegrino-wrinkling.parset index e01e58cf..790a91cd 100644 --- a/cosserat-continuum-wong-pellegrino-wrinkling.parset +++ b/cosserat-continuum-wong-pellegrino-wrinkling.parset @@ -21,7 +21,7 @@ numHomotopySteps = 1 tolerance = 1e-8 # Max number of steps of the trust region solver -maxTrustRegionSteps = 200 +maxTrustRegionSteps = 500 # Initial trust-region radius initialTrustRegionRadius = 0.001 @@ -71,7 +71,7 @@ lambda = 2.1796e+09 mu_c = 0 # Length scale parameter -L_c = 2.5e-5 +L_c = 2.5e-8 # Curvature exponent q = 2 @@ -92,4 +92,6 @@ problem = wong-pellegrino dirichletVerticesPredicate = "x[1] < 0.0001 or x[1] > 0.128 - 0.0001" # Initial deformation +#startFromFile = true +initialIterateFilename = cosserat_iterate_2.vtu initialDeformation = "[x[0] + 0.003*x[1] / 0.128, x[1], 0.002*math.cos(1e4*x[0])]" diff --git a/cosserat-continuum-wriggers-l-shape.parset b/cosserat-continuum-wriggers-l-shape.parset index d9be4481..73364963 100644 --- a/cosserat-continuum-wriggers-l-shape.parset +++ b/cosserat-continuum-wriggers-l-shape.parset @@ -17,10 +17,10 @@ numLevels = 1 numHomotopySteps = 1 # Tolerance of the trust region solver -tolerance = 1e-8 +tolerance = 1e-6 # Max number of steps of the trust region solver -maxTrustRegionSteps = 200 +maxTrustRegionSteps = 500 # Initial trust-region radius initialTrustRegionRadius = 0.1 @@ -76,7 +76,7 @@ L_c = 0.6e-3 q = 2 # Shear correction factor -kappa = 1 +kappa = 0.01 [] @@ -84,6 +84,8 @@ kappa = 1 # Boundary values ############################################# +problem = wriggers-l-shape + ### Python predicate specifying all Dirichlet grid vertices # x is the vertex coordinate dirichletVerticesPredicate = "x[0] < 0.001" @@ -93,7 +95,8 @@ dirichletVerticesPredicate = "x[0] < 0.001" neumannVerticesPredicate = "x[1] < -0.239" ### Neumann values, if needed -neumannValues = -5e6 0 0 +#neumannValues = -2e7 0 0 # Initial deformation -initialDeformation = "[x[0], x[1], -0.1*x[0]*(0.24+x[1])]" +#initialDeformation = "[x[0], x[1], -0.3*x[0]*(0.24+x[1])]" +initialDeformation = "[x[0], x[1], 0 if (x[0] < 0.225 or x[1] <-0.015) else 4*(x[0]-0.225)*(x[1]+0.015)]" diff --git a/dune/gfe/trustregionmmgbasesolver.hh b/dune/gfe/trustregionmmgbasesolver.hh index 48db1958..7b4d9aa2 100644 --- a/dune/gfe/trustregionmmgbasesolver.hh +++ b/dune/gfe/trustregionmmgbasesolver.hh @@ -3,6 +3,7 @@ #ifndef DUNE_GFE_TRUSTREGIONMMGBASESOLVER_HH #define DUNE_GFE_TRUSTREGIONMMGBASESOLVER_HH +#include <dune/istl/solver.hh> #include <dune/istl/umfpack.hh> #include <dune/solvers/solvers/iterativesolver.hh> @@ -101,6 +102,10 @@ void TrustRegionMMGBaseSolver<MatrixType, VectorType>::solve() // Modify Dirichlet rows /////////////////////////////////////////// + for (size_t j=0; j<modifiedStiffness.N(); j++) + for (int k=0; k<blocksize; k++) + modifiedStiffness[j][j][k][k] += 0.0; + for (size_t j=0; j<modifiedStiffness.N(); j++) { auto cIt = modifiedStiffness[j].begin(); @@ -128,6 +133,7 @@ void TrustRegionMMGBaseSolver<MatrixType, VectorType>::solve() ///////////////////////////////////////////////////////////////// Dune::InverseOperatorResult statistics; Dune::UMFPack<MatrixType> solver(modifiedStiffness); + solver.setOption(UMFPACK_PRL, 0); // no output at all solver.apply(*x_, modifiedRhs, statistics); // Model increase? Then the matrix is not positive definite -- fall back to ipopt solver @@ -138,17 +144,21 @@ void TrustRegionMMGBaseSolver<MatrixType, VectorType>::solve() matrix_->umv(*x_, tmp); double modelDecrease = (modifiedRhs*(*x_)) - 0.5 * ((*x_)*tmp); - if (std::isnan(modelDecrease) or modelDecrease < 0) + if (std::isnan(modelDecrease) or modelDecrease < 0 or true) { std::cout << "Model increase: " << -modelDecrease << ", falling back to slower solver" << std::endl; + std::cout << "Total VARIABLES: " << x_->size() << " x " << blocksize << " = " << x_->size()*blocksize << std::endl; + std::cout << "Dirichlet DOFS: " << this->ignoreNodes_->count() << std::endl; + QuadraticIPOptSolver<MatrixType, VectorType> baseSolver; - baseSolver.verbosity_ = NumProc::QUIET; + baseSolver.verbosity_ = NumProc::REDUCED; baseSolver.tolerance_ = 1e-8; *x_ = 0; baseSolver.setProblem(*matrix_, *x_, *rhs_); baseSolver.obstacles_ = obstacles_; + baseSolver.ignoreNodes_ = this->ignoreNodes_; baseSolver.solve(); } @@ -173,12 +183,13 @@ void TrustRegionMMGBaseSolver<MatrixType, VectorType>::solve() std::cout << "Inadmissible step, falling back to slower solver" << std::endl; QuadraticIPOptSolver<MatrixType, VectorType> baseSolver; - baseSolver.verbosity_ = NumProc::QUIET; + baseSolver.verbosity_ = NumProc::REDUCED; baseSolver.tolerance_ = 1e-8; *x_ = 0; baseSolver.setProblem(*matrix_, *x_, *rhs_); baseSolver.obstacles_ = obstacles_; + baseSolver.ignoreNodes_ = this->ignoreNodes_; baseSolver.solve(); } diff --git a/finite-strain-elasticity.cc b/finite-strain-elasticity.cc index 4f995703..51244941 100644 --- a/finite-strain-elasticity.cc +++ b/finite-strain-elasticity.cc @@ -134,7 +134,8 @@ int main (int argc, char *argv[]) try typedef GridType::LeafGridView GridView; GridView gridView = grid->leafGridView(); - typedef P1NodalBasis<GridView,double> FEBasis; +// typedef P1NodalBasis<GridView,double> FEBasis; + typedef P2NodalBasis<GridView,double> FEBasis; FEBasis feBasis(gridView); // ///////////////////////////////////////// @@ -271,14 +272,14 @@ int main (int argc, char *argv[]) try std::vector<FieldVector<double,3> > pointLoads(x.size()); std::fill(pointLoads.begin(), pointLoads.end(), 0); - pointLoads[1372] = parameterSet.get<FieldVector<double,3> >("neumannValues"); - pointLoads[1372] *= 0.5; +// pointLoads[1372] = parameterSet.get<FieldVector<double,3> >("neumannValues"); +// pointLoads[1372] *= 0.5; // ///////////////////////////////////////////////// // Create a Riemannian trust-region solver // ///////////////////////////////////////////////// - TrustRegionSolver<GridType,SolutionType> solver; + TrustRegionSolver<FEBasis,SolutionType> solver; solver.setup(*grid, &assembler, x, diff --git a/harmonicmaps.cc b/harmonicmaps.cc index d37a1c3a..d8de45aa 100644 --- a/harmonicmaps.cc +++ b/harmonicmaps.cc @@ -33,6 +33,7 @@ #include <dune/gfe/unitvector.hh> #include <dune/gfe/realtuple.hh> #include <dune/gfe/harmonicenergystiffness.hh> +#include <dune/gfe/chiralskyrmionenergy.hh> #include <dune/gfe/geodesicfeassembler.hh> #include <dune/gfe/riemanniantrsolver.hh> @@ -229,6 +230,7 @@ int main (int argc, char *argv[]) try typedef P1NodalBasis<typename GridType::LeafGridView,double> FEBasis; FEBasis feBasis(grid.leafGridView()); + GFE::ChiralSkyrmionEnergy<GridType::LeafGridView, FEBasis::LocalFiniteElement, double> chiralSkyrmionEnergy; HarmonicEnergyLocalStiffness<GridType::LeafGridView, FEBasis::LocalFiniteElement, TargetSpace> harmonicEnergyLocalStiffness; GeodesicFEAssembler<FEBasis,TargetSpace> assembler(grid.leafGridView(), diff --git a/mixed-cosserat-continuum.cc b/mixed-cosserat-continuum.cc index 376cb2be..993b0582 100644 --- a/mixed-cosserat-continuum.cc +++ b/mixed-cosserat-continuum.cc @@ -243,12 +243,14 @@ int main (int argc, char *argv[]) try typedef GridType::LeafGridView GridView; GridView gridView = grid->leafGridView(); - typedef P2NodalBasis<GridView,double> DeformationFEBasis; - typedef P1NodalBasis<GridView,double> OrientationFEBasis; + typedef P3NodalBasis<GridView,double> DeformationFEBasis; + typedef P2NodalBasis<GridView,double> OrientationFEBasis; DeformationFEBasis deformationFEBasis(gridView); OrientationFEBasis orientationFEBasis(gridView); + std::cout << "Deformation: " << deformationFEBasis.size() << ", orientation: " << orientationFEBasis.size() << std::endl; + // ///////////////////////////////////////// // Read Dirichlet values // ///////////////////////////////////////// diff --git a/test/Makefile.am b/test/Makefile.am index adf1c5c5..d56b4d36 100644 --- a/test/Makefile.am +++ b/test/Makefile.am @@ -19,6 +19,7 @@ check_PROGRAMS = adolctest \ globalgfetestfunctionbasistest \ harmonicenergytest \ interpolationtest \ + interillustration \ localgeodesicfefunctiontest \ localgeodesicfestiffnesstest \ localgfetestfunctiontest \ @@ -30,11 +31,12 @@ check_PROGRAMS = adolctest \ rotationtest \ svdtest \ targetspacetest \ - true-adolctest + true-adolctest \ + vtkreadertest adolctest_SOURCES = adolctest.cc adolctest_LDFLAGS = $(ADOLC_LDFLAGS) $(PYTHON_LDFLAGS) $(IPOPT_LDFLAGS) -adolctest_LDADD = $(ADOLC_LIBS) $(PYTHON_LIBS) $(IPOPT_LIBS) +adolctest_LDADD = $(ADOLC_LIBS) $(PYTHON_LIBS) $(IPOPT_LIBS) -lmpfr -lgmpxx -lgmp true_adolctest_SOURCES = true-adolctest.cc true_adolctest_LDFLAGS = $(ADOLC_LDFLAGS) $(PYTHON_LDFLAGS) $(IPOPT_LDFLAGS) @@ -64,6 +66,8 @@ harmonicenergytest_SOURCES = harmonicenergytest.cc interpolationtest_SOURCES = interpolationtest.cc +interillustration_SOURCES = interillustration.cc + cosseratenergytest_SOURCES = cosseratenergytest.cc averagedistanceassemblertest_SOURCES = averagedistanceassemblertest.cc @@ -76,6 +80,10 @@ targetspacetest_SOURCES = targetspacetest.cc svdtest_SOURCES = svdtest.cc +vtkreadertest_SOURCES = vtkreadertest.cc +vtkreadertest_CXXFLAGS = -DHAVE_TINYXML2 +vtkreadertest_LDADD = -ltinyxml2 + # don't follow the full GNU-standard # we need automake 1.5 AUTOMAKE_OPTIONS = foreign 1.5 diff --git a/test/adolctest.cc b/test/adolctest.cc index 2e60666a..bfbcc17f 100644 --- a/test/adolctest.cc +++ b/test/adolctest.cc @@ -30,6 +30,7 @@ typedef double FDType; #include <dune/istl/io.hh> #include <dune/fufem/functionspacebases/p2nodalbasis.hh> +#include <dune/fufem/functiontools/basisinterpolator.hh> #include <dune/gfe/rigidbodymotion.hh> @@ -47,7 +48,19 @@ typedef RigidBodyMotion<double,3> TargetSpace; using namespace Dune; +class Identity +: public VirtualFunction<FieldVector<double,2>, FieldVector<double,3> > +{ +public: + void evaluate(const FieldVector<double,2>& in, + FieldVector<double,3>& out) const + { + out[0] = in[0]; + out[1] = in[1]; + out[2] = 0.0; + } +}; /** \brief Assembles energy gradient and Hessian with ADOL-C */ @@ -446,6 +459,7 @@ int main (int argc, char *argv[]) try //////////////////////////////////////////7 // Read initial iterate from file //////////////////////////////////////////7 +#if 0 Dune::BlockVector<FieldVector<double,7> > xEmbedded(x.size()); std::ifstream file("dangerous_iterate", std::ios::in|std::ios::binary); @@ -458,6 +472,15 @@ int main (int argc, char *argv[]) try for (int ii=0; ii<x.size(); ii++) x[ii] = xEmbedded[ii]; +#else + Identity identity; + + std::vector<FieldVector<double,3> > v; + Functions::interpolate(feBasis, v, identity); + + for (size_t i=0; i<x.size(); i++) + x[i].r = v[i]; +#endif // //////////////////////////////////////////////////////////// // Create an assembler for the energy functional diff --git a/test/targetspacetest.cc b/test/targetspacetest.cc index 67aa9588..820f3e81 100644 --- a/test/targetspacetest.cc +++ b/test/targetspacetest.cc @@ -345,7 +345,7 @@ void test() // Test each element in the list for (int i=0; i<nTestPoints; i++) { - testOrthonormalFrame<TargetSpace>(testPoints[i]); + //testOrthonormalFrame<TargetSpace>(testPoints[i]); for (int j=0; j<nTestPoints; j++) { @@ -354,8 +354,13 @@ void test() testPointPair[1] = testPoints[j]; if (diameter(testPointPair) > TargetSpace::convexityRadius) continue; + + TargetSpace p = testPointPair[0]; + TargetSpace q = testPointPair[1]; + std::cout << "p: " << testPointPair[0] << ", q: " << testPointPair[1] << std::endl; + std::cout << TargetSpace::exp(p, TargetSpace::log(p,q)) << std::endl; - testDerivativesOfSquaredDistance<TargetSpace>(testPoints[i], testPoints[j]); + //testDerivativesOfSquaredDistance<TargetSpace>(testPoints[i], testPoints[j]); } @@ -366,18 +371,18 @@ void test() int main() try { - test<RealTuple<double,1> >(); - test<RealTuple<double,3> >(); +// test<RealTuple<double,1> >(); +// test<RealTuple<double,3> >(); test<UnitVector<double,2> >(); test<UnitVector<double,3> >(); test<UnitVector<double,4> >(); - test<Rotation<double,3> >(); - - test<RigidBodyMotion<double,3> >(); - - test<HyperbolicHalfspacePoint<double,2> >(); +// test<Rotation<double,3> >(); +// +// test<RigidBodyMotion<double,3> >(); +// +// test<HyperbolicHalfspacePoint<double,2> >(); } catch (Exception e) { -- GitLab