diff --git a/dune/tectonic/globalruinanonlinearity.hh b/dune/tectonic/globalruinanonlinearity.hh index 101fd0396d4bbe05aa8bcc8ebec4d6ded94548e4..2cac0eb69c19147be3ee63367799dd7f9734671f 100644 --- a/dune/tectonic/globalruinanonlinearity.hh +++ b/dune/tectonic/globalruinanonlinearity.hh @@ -18,26 +18,25 @@ template <class MatrixType, class VectorType> class GlobalRuinaNonlinearity : public GlobalNonlinearity<MatrixType, VectorType> { private: - typedef shared_ptr<BlockVector<FieldVector<double, 1>> const> dataptr; + typedef BlockVector<FieldVector<double, 1>> const &dataref; public: using GlobalNonlinearity<MatrixType, VectorType>::dim; - GlobalRuinaNonlinearity(dataptr nodalIntegrals, dataptr a, dataptr mu, - dataptr eta, dataptr normalStress, dataptr b, - dataptr state, dataptr L, double h) - // TODO: can we get the size from another place? - : restrictions(nodalIntegrals->size()) { + GlobalRuinaNonlinearity(dataref nodalIntegrals, dataref a, dataref mu, + dataref eta, dataref normalStress, dataref b, + dataref state, dataref L, double h) + : restrictions(nodalIntegrals.size()) { auto trivialNonlinearity = make_shared<LocalNonlinearity<dim> const>( make_shared<TrivialFunction const>()); for (size_t i = 0; i < restrictions.size(); ++i) { restrictions[i] = - (*nodalIntegrals)[i] == 0 + nodalIntegrals[i] == 0 ? trivialNonlinearity : make_shared<LocalNonlinearity<dim> const>( make_shared<RuinaFunction const>( - (*nodalIntegrals)[i], (*a)[i], (*mu)[i], (*eta)[i], - (*normalStress)[i], (*b)[i], (*state)[i], (*L)[i], h)); + nodalIntegrals[i], a[i], mu[i], eta[i], normalStress[i], + b[i], state[i], L[i], h)); } } diff --git a/src/assemblers.cc b/src/assemblers.cc index b6f8c7264359ff65aa919c0f044cbb45dec7accf..1c5d5ed059269c7c926258f6048aab2734631f35 100644 --- a/src/assemblers.cc +++ b/src/assemblers.cc @@ -56,33 +56,31 @@ template <class MatrixType, class VectorType> Dune::shared_ptr<Dune::GlobalNonlinearity<MatrixType, VectorType> const> assemble_nonlinearity( Dune::ParameterTree const &parset, - Dune::shared_ptr<Dune::BlockVector<Dune::FieldVector<double, 1>>> - nodalIntegrals, - Dune::shared_ptr<Dune::BlockVector<Dune::FieldVector<double, 1>>> state, - double h) { - auto const size = (*nodalIntegrals).size(); + Dune::BlockVector<Dune::FieldVector<double, 1>> const &nodalIntegrals, + Dune::BlockVector<Dune::FieldVector<double, 1>> const &state, double h) { + auto const size = nodalIntegrals.size(); typedef Dune::BlockVector<Dune::FieldVector<double, 1>> SingletonVectorType; // {{{ Assemble terms for the nonlinearity - auto mu = Dune::make_shared<SingletonVectorType>(size); - *mu = parset.get<double>("mu"); + SingletonVectorType mu(size); + mu = parset.get<double>("mu"); - auto normalStress = Dune::make_shared<SingletonVectorType>(size); - *normalStress = parset.get<double>("normalstress"); + SingletonVectorType normalStress(size); + normalStress = parset.get<double>("normalstress"); switch (parset.get<Config::model>("model")) { case Config::Exponential: { - auto a = Dune::make_shared<SingletonVectorType>(size); - *a = parset.get<double>("ruina.a"); + SingletonVectorType a(size); + a = parset.get<double>("ruina.a"); - auto eta = Dune::make_shared<SingletonVectorType>(size); - *eta = parset.get<double>("eta"); + SingletonVectorType eta(size); + eta = parset.get<double>("eta"); - auto b = Dune::make_shared<SingletonVectorType>(size); - *b = parset.get<double>("ruina.b"); + SingletonVectorType b(size); + b = parset.get<double>("ruina.b"); - auto L = Dune::make_shared<SingletonVectorType>(size); - *L = parset.get<double>("ruina.L"); + SingletonVectorType L(size); + L = parset.get<double>("ruina.L"); return Dune::make_shared< Dune::GlobalRuinaNonlinearity<MatrixType, VectorType> const>( diff --git a/src/assemblers.hh b/src/assemblers.hh index d98ac6f8e17487da6e6bc93c60f67a779e25c29d..58a9d780902a763c8becce9eeb4767d7fb25538d 100644 --- a/src/assemblers.hh +++ b/src/assemblers.hh @@ -26,9 +26,7 @@ template <class MatrixType, class VectorType> Dune::shared_ptr<Dune::GlobalNonlinearity<MatrixType, VectorType> const> assemble_nonlinearity( Dune::ParameterTree const &parset, - Dune::shared_ptr<Dune::BlockVector<Dune::FieldVector<double, 1>>> - nodalIntegrals, - Dune::shared_ptr<Dune::BlockVector<Dune::FieldVector<double, 1>>> state, - double h); + Dune::BlockVector<Dune::FieldVector<double, 1>> const &nodalIntegrals, + Dune::BlockVector<Dune::FieldVector<double, 1>> const &state, double h); #endif diff --git a/src/assemblers_tmpl.cc b/src/assemblers_tmpl.cc index 9f8b4bc214393dbd3e0f252cd4e996aa44593cba..388e60cf6139689f2ce4600c981b6c94a4bd28a8 100644 --- a/src/assemblers_tmpl.cc +++ b/src/assemblers_tmpl.cc @@ -33,7 +33,5 @@ template Dune::shared_ptr< Dune::GlobalNonlinearity<MatrixType, VectorType> const> assemble_nonlinearity<MatrixType, VectorType>( Dune::ParameterTree const &parset, - Dune::shared_ptr<Dune::BlockVector<Dune::FieldVector<double, 1>>> - nodalIntegrals, - Dune::shared_ptr<Dune::BlockVector<Dune::FieldVector<double, 1>>> state, - double h); + Dune::BlockVector<Dune::FieldVector<double, 1>> const &nodalIntegrals, + Dune::BlockVector<Dune::FieldVector<double, 1>> const &state, double h); diff --git a/src/one-body-sample.cc b/src/one-body-sample.cc index 15d301843849386426b03f475bd91b658e5465e0..d04c8cc4f17064a16e762894cd967686dc369103 100644 --- a/src/one-body-sample.cc +++ b/src/one-body-sample.cc @@ -198,8 +198,7 @@ int main(int argc, char *argv[]) { VectorType u_diff(finestSize); u_diff = 0.0; // Has to be zero! - auto alpha = Dune::make_shared<SingletonVectorType>(finestSize); - *alpha = alpha_old; + SingletonVectorType alpha(alpha_old); SingletonVectorType vonMisesStress; VectorType rhs; @@ -268,7 +267,7 @@ int main(int argc, char *argv[]) { ++state_fpi) { auto myGlobalNonlinearity = assemble_nonlinearity<MatrixType, VectorType>( - parset.sub("boundary.friction"), nodalIntegrals, alpha, h); + parset.sub("boundary.friction"), *nodalIntegrals, alpha, h); MyConvexProblemType const myConvexProblem(stiffnessMatrix, *myGlobalNonlinearity, rhs); @@ -295,11 +294,10 @@ int main(int argc, char *argv[]) { switch (parset.get<Config::state_model>( "boundary.friction.state.model")) { case Config::Dieterich: - (*alpha)[i] = - state_update_dieterich(h, unorm / L, alpha_old[i]); + alpha[i] = state_update_dieterich(h, unorm / L, alpha_old[i]); break; case Config::Ruina: - (*alpha)[i] = state_update_ruina(h, unorm / L, alpha_old[i]); + alpha[i] = state_update_ruina(h, unorm / L, alpha_old[i]); break; } } @@ -318,7 +316,7 @@ int main(int argc, char *argv[]) { if (parset.get<bool>("writeEvolution")) { double out; neumannFunction.evaluate(time, out); - octave_writer << (*alpha)[first_frictional_node][0] << " " + octave_writer << alpha[first_frictional_node][0] << " " << u[first_frictional_node][0] * 1e6 << " " << out << std::endl; } @@ -329,7 +327,7 @@ int main(int argc, char *argv[]) { if (parset.get<bool>("printVelocitySteppingComparison")) { double const v = u_diff[first_frictional_node].two_norm() / h / L; - double const euler = (*alpha)[first_frictional_node]; + double const euler = alpha[first_frictional_node]; double direct; if (run < 120) { direct = euler; @@ -352,7 +350,7 @@ int main(int argc, char *argv[]) { // Record the coefficient of friction at a fixed node if (parset.get<bool>("printCoefficient")) { double const V = u_diff[first_frictional_node].two_norm(); - double const state = (*alpha)[first_frictional_node]; + double const state = alpha[first_frictional_node]; coefficient_writer << (mu + a * std::log(V * eta) + b * (state - std::log(eta * L))) << std::endl; @@ -360,7 +358,7 @@ int main(int argc, char *argv[]) { } u += u_diff; - alpha_old = *alpha; + alpha_old = alpha; // Compute von Mises stress and write everything to a file if (parset.get<bool>("writeVTK")) { @@ -372,7 +370,7 @@ int main(int argc, char *argv[]) { .assemble(localStressAssembler, vonMisesStress, true); writeVtk<P1Basis, P0Basis, VectorType, SingletonVectorType, GridView>( - p1Basis, u, *alpha, p0Basis, vonMisesStress, leafView, + p1Basis, u, alpha, p0Basis, vonMisesStress, leafView, (boost::format("obs%d") % run).str()); } }