From 58d4d5d27d8a784fdd04ded5c1a16824b3d3e1ef Mon Sep 17 00:00:00 2001 From: Tarun Mitruka Date: Fri, 1 Sep 2023 16:34:33 +0200 Subject: [PATCH] fix iks003 and iks007 --- src/iks003_incompressibleLinearElasticity.cpp | 5 ++--- src/iks004_kirchhoffPlate.cpp | 2 -- src/iks007_vonMisesTruss.cpp | 13 +++++++------ 3 files changed, 9 insertions(+), 11 deletions(-) diff --git a/src/iks003_incompressibleLinearElasticity.cpp b/src/iks003_incompressibleLinearElasticity.cpp index 0cc5a68..032f86c 100644 --- a/src/iks003_incompressibleLinearElasticity.cpp +++ b/src/iks003_incompressibleLinearElasticity.cpp @@ -192,7 +192,6 @@ int main(int argc, char **argv) { /// Create assembler auto sparseFlatAssembler = SparseFlatAssembler(fes, dirichletValues); - auto denseFlatAssembler = DenseFlatAssembler(fes, dirichletValues); /// Create function for external forces and stiffness matrix double lambda = 0; @@ -204,7 +203,7 @@ int main(int argc, char **argv) { auto fextFunction = [&](auto &&lambdaLocal, auto &&dLocal) -> auto & { req.insertGlobalSolution(Ikarus::FESolutions::displacement, dLocal) .insertParameter(Ikarus::FEParameter::loadfactor, lambdaLocal); - return denseFlatAssembler.getReducedVector(req); + return sparseFlatAssembler.getReducedVector(req); }; auto KFunction = [&](auto &&lambdaLocal, auto &&dLocal) -> auto & { req.insertGlobalSolution(Ikarus::FESolutions::displacement, dLocal) @@ -218,7 +217,7 @@ int main(int argc, char **argv) { ld.compute(K); if (ld.info() != Eigen::Success) DUNE_THROW(Dune::MathError, "Failed Compute"); - d -= denseFlatAssembler.createFullVector(ld.solve(R)); + d -= sparseFlatAssembler.createFullVector(ld.solve(R)); if (ld.info() != Eigen::Success) DUNE_THROW(Dune::MathError, "Failed Solve"); /// Postprocess diff --git a/src/iks004_kirchhoffPlate.cpp b/src/iks004_kirchhoffPlate.cpp index d207762..4fc8e97 100644 --- a/src/iks004_kirchhoffPlate.cpp +++ b/src/iks004_kirchhoffPlate.cpp @@ -285,8 +285,6 @@ int main(int argc, char **argv) { const auto &rule = Dune::QuadratureRules::rule( ele.type(), 2 * localView.tree().finiteElement().localBasis().order()); for (auto gp : rule) { - const auto gpGlobalPos = geo.global(gp.position()); - const auto w_ex = localwAna(gp.position()); const auto w_fe = localw(gp.position()); l2_error += Dune::power(w_ex - w_fe, 2) * ele.geometry().integrationElement(gp.position()) * gp.weight(); diff --git a/src/iks007_vonMisesTruss.cpp b/src/iks007_vonMisesTruss.cpp index 9cb6691..de1f992 100644 --- a/src/iks007_vonMisesTruss.cpp +++ b/src/iks007_vonMisesTruss.cpp @@ -88,7 +88,7 @@ class Truss : public PowerBasisFE { const ScalarType Egl = 0.5 * (lsquared - LRefsquared) / LRefsquared; - return 0.5 * EA / sqrt(LRefsquared) * Egl * Egl; + return 0.5 * EA * sqrt(LRefsquared) * Egl * Egl; } private: @@ -100,7 +100,7 @@ int main(int argc, char **argv) { /// Construct grid Dune::GridFactory> gridFactory; const double h = 1.0; - const double L = 1.0; + const double L = 2.0; gridFactory.insertVertex({0, 0}); gridFactory.insertVertex({L, h}); gridFactory.insertVertex({2 * L, 0}); @@ -200,17 +200,18 @@ int main(int argc, char **argv) { auto analyticalLoadDisplacementCurve = [&](auto &w) { const double Ltruss = std::sqrt(h * h + L * L); - return EA * Dune::power(h, 3) / Dune::power(Ltruss, 3) + return 2.0 * EA * Dune::power(h, 3) / Dune::power(Ltruss, 3) * (w / h - 1.5 * Dune::power(w / h, 2) + 0.5 * Dune::power(w / h, 3)); }; std::vector x = linspace(0.0, dVec.maxCoeff()); std::vector y1 = transform(x, [&](auto x) { return analyticalLoadDisplacementCurve(x); }); auto p = plot(x, y1, dVec, lambdaVec); - p[0]->line_width(2); + p[0]->line_width(3); p[1]->line_width(2); p[1]->marker(line_spec::marker_style::asterisk); + // save("vonMisesTruss.png"); // f->draw(); - // using namespace std::chrono_literals; - // std::this_thread::sleep_for(5s); + // using namespace std::chrono_literals; + // std::this_thread::sleep_for(5s); }