diff --git a/test/test_model/test_solid_mechanics_model/test_cohesive/test_cohesive_intrinsic_impl/CMakeLists.txt b/test/test_model/test_solid_mechanics_model/test_cohesive/test_cohesive_intrinsic_impl/CMakeLists.txt index 8403ec4ed..9053cd8c0 100644 --- a/test/test_model/test_solid_mechanics_model/test_cohesive/test_cohesive_intrinsic_impl/CMakeLists.txt +++ b/test/test_model/test_solid_mechanics_model/test_cohesive/test_cohesive_intrinsic_impl/CMakeLists.txt @@ -1,39 +1,41 @@ #=============================================================================== # @file CMakeLists.txt # @author Marco Vocialta # @date Thu Apr 26 15:14:00 2012 # # @section LICENSE # # Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne) # Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides) # # Akantu is free software: you can redistribute it and/or modify it under the # terms of the GNU Lesser General Public License as published by the Free # Software Foundation, either version 3 of the License, or (at your option) any # later version. # # Akantu is distributed in the hope that it will be useful, but WITHOUT ANY # WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR # A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more # details. # # You should have received a copy of the GNU Lesser General Public License # along with Akantu. If not, see . # # @section DESCRIPTION # #=============================================================================== register_test(test_cohesive_intrinsic_impl test_cohesive_intrinsic_impl.cc) #=============================================================================== configure_file( ${CMAKE_CURRENT_SOURCE_DIR}/material.dat ${CMAKE_CURRENT_BINARY_DIR}/material.dat COPYONLY ) + +file(COPY material.dat DESTINATION .) file(COPY simple.msh DESTINATION .) file(MAKE_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}/paraview) file(COPY simple.msh DESTINATION .) \ No newline at end of file diff --git a/test/test_model/test_solid_mechanics_model/test_cohesive/test_cohesive_intrinsic_impl/stiffness_matrix.verified b/test/test_model/test_solid_mechanics_model/test_cohesive/test_cohesive_intrinsic_impl/stiffness_matrix.verified new file mode 100644 index 000000000..81770e0ca --- /dev/null +++ b/test/test_model/test_solid_mechanics_model/test_cohesive/test_cohesive_intrinsic_impl/stiffness_matrix.verified @@ -0,0 +1,308 @@ +%%MatrixMarket matrix coordinate real symmetric +36 36 306 +1 1 4.5 +1 2 -1.85037170770859e-16 +1 5 1 +1 6 8.32667268468867e-17 +1 31 -2.12792746386488e-16 +1 32 0.5 +1 3 -4 +1 4 -2.77555756156289e-16 +1 11 -1.55431223447522e-15 +1 12 -5.55111512312579e-17 +1 9 3.02535774210355e-15 +1 10 -2 +2 2 4.5 +2 5 -0.5 +2 6 0.5 +2 31 0.5 +2 32 -6.47630097698008e-17 +2 3 2 +2 4 -2 +2 11 0 +2 12 -6.66133814775094e-16 +2 9 -2 +2 10 1.99840144432528e-15 +5 5 4.5 +5 6 -1.5 +5 31 0.5 +5 32 0 +5 3 -4 +5 4 2 +5 11 -2 +5 12 2.22044604925031e-16 +5 9 7.7715611723761e-16 +5 10 -2.22044604925031e-16 +6 6 4.5 +6 31 -0.5 +6 32 1 +6 3 -3.19189119579733e-16 +6 4 -2 +6 11 2 +6 12 -4 +6 9 -2.35922392732846e-16 +6 10 3.33066907387547e-16 +31 31 4.5 +31 32 2.32317978986721e-16 +31 3 -1.11022302462516e-16 +31 4 -4.93038065763132e-32 +31 11 -2 +31 12 2 +31 9 1.19348975147204e-15 +31 10 -2 +32 32 5.10406262854645 +32 3 0 +32 4 -2.22044604925031e-16 +32 11 0 +32 12 -4 +32 9 -2 +32 10 9.29811783123569e-16 +3 3 12 +3 4 -2 +3 11 1.22124532708767e-15 +3 12 -2 +3 9 -4 +3 10 2 +4 4 12 +4 11 -2 +4 12 -8.88178419700125e-16 +4 9 2 +4 10 -8 +11 11 12 +11 12 -2 +11 9 -8 +11 10 2 +12 12 12 +12 9 2 +12 10 -4 +9 9 24 +9 10 -4 +10 10 24 +1 13 0.5 +1 14 -0.5 +1 15 -4.16333634234434e-16 +1 16 -1.11022302462516e-16 +1 7 -2 +1 8 2 +2 13 -6.47630097698008e-17 +2 14 1 +2 15 7.40148683083438e-17 +2 16 -1.0547118733939e-15 +2 7 4.44089209850063e-16 +2 8 -4 +31 13 1 +31 14 -1.23259516440783e-32 +31 15 -4 +31 16 4.64635957973443e-16 +31 7 2.77555756156289e-17 +31 8 0 +32 13 -0.5 +32 14 0.197968685726773 +32 15 2 +32 16 -1.39593737145355 +32 7 -5.55111512312578e-17 +32 8 1.38777878078145e-17 +13 13 4.5 +13 14 -1.5 +13 9 0 +13 10 0 +13 15 -4 +13 16 2 +13 7 -2 +13 8 0 +14 14 5.10406262854645 +14 9 -1.11022302462516e-16 +14 10 4.44089209850063e-16 +14 15 -3.53613655510927e-16 +14 16 -1.39593737145355 +14 7 2 +14 8 -4 +9 15 -4 +9 16 2 +9 7 -8 +9 8 2 +10 15 2 +10 16 -8 +10 7 2 +10 8 -4 +15 15 12 +15 16 -2 +15 7 1.77635683940025e-15 +15 8 -2 +16 16 14.4162505141858 +16 7 -2 +16 8 2.22044604925031e-15 +7 7 12 +7 8 -2 +8 8 12 +33 33 4.5 +33 34 1.10000786939369e-16 +33 17 1 +33 18 1.20274161001059e-16 +33 29 -3.05311331771918e-16 +33 30 0.5 +33 35 -4 +33 36 -9.64236319054764e-16 +33 23 -1.22124532708767e-15 +33 24 -1.29526019539602e-16 +33 21 7.07767178198537e-16 +33 22 -2 +34 34 5.10406262854646 +34 17 -0.5 +34 18 0.197968685726773 +34 29 0.5 +34 30 -5.55111512312577e-17 +34 35 2 +34 36 -1.39593737145355 +34 23 5.55111512312578e-17 +34 24 -5.55111512312578e-16 +34 21 -2 +34 22 -2.4980018054066e-16 +17 17 4.5 +17 18 -1.5 +17 29 0.5 +17 30 0 +17 35 -4 +17 36 2 +17 23 -2 +17 24 1.66533453693773e-16 +17 21 3.33066907387547e-16 +17 22 -1.66533453693773e-16 +18 18 5.10406262854645 +18 29 -0.5 +18 30 1 +18 35 2.00957989624968e-16 +18 36 -1.39593737145355 +18 23 2 +18 24 -4 +18 21 -2.22044604925031e-16 +18 22 5.55111512312578e-16 +29 29 4.5 +29 30 0 +29 35 -1.11022302462516e-16 +29 36 1.66533453693773e-16 +29 23 -2 +29 24 2 +29 21 2.88657986402541e-15 +29 22 -2 +30 30 4.5 +30 35 0 +30 36 -2.22044604925031e-16 +30 23 0 +30 24 -4 +30 21 -2 +30 22 1.60982338570648e-15 +35 35 12 +35 36 -2 +35 23 7.7715611723761e-16 +35 24 -2 +35 21 -4 +35 22 2 +36 36 14.4162505141858 +36 23 -2 +36 24 0 +36 21 2 +36 22 -8 +23 23 12 +23 24 -2 +23 21 -8 +23 22 2 +24 24 12 +24 21 2 +24 22 -4 +21 21 24 +21 22 -4 +22 22 24 +33 25 0.5 +33 26 -0.5 +33 27 2.35922392732846e-16 +33 28 -3.88578058618805e-16 +33 19 -2 +33 20 2 +34 25 3.70074341541719e-17 +34 26 1 +34 27 -2.96059473233375e-16 +34 28 -8.32667268468867e-17 +34 19 -2.46519032881566e-32 +34 20 -4 +29 25 1 +29 26 0 +29 27 -4 +29 28 0 +29 19 3.33066907387547e-16 +29 20 0 +30 25 -0.5 +30 26 0.5 +30 27 2 +30 28 -2 +30 19 0 +30 20 1.66533453693773e-16 +25 25 4.5 +25 26 -1.5 +25 21 6.66133814775094e-16 +25 22 -2.22044604925031e-16 +25 27 -4 +25 28 2 +25 19 -2 +25 20 -1.11022302462516e-16 +26 26 4.5 +26 21 -5.55111512312578e-16 +26 22 1.33226762955019e-15 +26 27 5.55111512312578e-16 +26 28 -2 +26 19 2 +26 20 -4 +21 27 -4.00000000000001 +21 28 2 +21 19 -8 +21 20 2 +22 27 2 +22 28 -8.00000000000001 +22 19 2 +22 20 -4 +27 27 12 +27 28 -2 +27 19 8.88178419700125e-16 +27 20 -2 +28 28 12 +28 19 -2 +28 20 1.33226762955019e-15 +19 19 12 +19 20 -2 +20 20 12 +31 17 -1.19130347991335e-31 +31 18 -2.32317978986721e-16 +31 33 5.95651739956675e-32 +31 34 1.23259516440783e-32 +31 35 -1.19130347991335e-31 +31 36 -4.64635957973443e-16 +32 17 -2.32317978986721e-16 +32 18 -0.604062628546454 +32 33 1.23259516440783e-32 +32 34 0.302031314273227 +32 35 -4.64635957973443e-16 +32 36 -0.604062628546454 +13 17 5.95651739956675e-32 +13 18 -1.23259516440783e-32 +13 33 -1.19130347991335e-31 +13 34 2.32317978986721e-16 +13 35 -1.19130347991335e-31 +13 36 4.64635957973443e-16 +14 17 -1.23259516440783e-32 +14 18 0.302031314273227 +14 33 2.32317978986721e-16 +14 34 -0.604062628546455 +14 35 4.64635957973443e-16 +14 36 -0.604062628546455 +15 17 -1.19130347991335e-31 +15 18 -4.64635957973443e-16 +15 33 -1.19130347991335e-31 +15 34 4.64635957973443e-16 +15 35 -4.7652139196534e-31 +15 36 0 +16 17 -4.64635957973443e-16 +16 18 -0.604062628546454 +16 33 4.64635957973443e-16 +16 34 -0.604062628546454 +16 35 0 +16 36 -2.41625051418582 diff --git a/test/test_model/test_solid_mechanics_model/test_cohesive/test_cohesive_intrinsic_impl/test_cohesive_intrinsic_impl.cc b/test/test_model/test_solid_mechanics_model/test_cohesive/test_cohesive_intrinsic_impl/test_cohesive_intrinsic_impl.cc index 9c7484355..17d68da1c 100644 --- a/test/test_model/test_solid_mechanics_model/test_cohesive/test_cohesive_intrinsic_impl/test_cohesive_intrinsic_impl.cc +++ b/test/test_model/test_solid_mechanics_model/test_cohesive/test_cohesive_intrinsic_impl/test_cohesive_intrinsic_impl.cc @@ -1,209 +1,211 @@ /** * @file test_cohesive.cc * @author Marco Vocialta * @date Fri Feb 24 14:32:31 2012 * * @brief Test for cohesive elements * * @section LICENSE * * Copyright (©) 2010-2011 EPFL (Ecole Polytechnique Fédérale de Lausanne) * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides) * * Akantu is free software: you can redistribute it and/or modify it under the * terms of the GNU Lesser General Public License as published by the Free * Software Foundation, either version 3 of the License, or (at your option) any * later version. * * Akantu is distributed in the hope that it will be useful, but WITHOUT ANY * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR * A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more * details. * * You should have received a copy of the GNU Lesser General Public License * along with Akantu. If not, see . * */ /* -------------------------------------------------------------------------- */ #include #include #include /* -------------------------------------------------------------------------- */ #include "aka_common.hh" #include "mesh.hh" #include "mesh_io.hh" #include "mesh_io_msh.hh" #include "mesh_utils.hh" #include "solid_mechanics_model_cohesive.hh" #include "material.hh" #if defined(AKANTU_USE_IOHELPER) # include "io_helper.hh" #endif /* -------------------------------------------------------------------------- */ using namespace akantu; int main(int argc, char *argv[]) { initialize(argc, argv); debug::setDebugLevel(dblError); const UInt spatial_dimension = 2; const ElementType type = _triangle_6; Mesh mesh(spatial_dimension); MeshIOMSH mesh_io; mesh_io.read("simple.msh", mesh); SolidMechanicsModelCohesive model(mesh); /// model initialization model.initFull("material.dat", _static); const Mesh & mesh_facets = model.getMeshFacets(); const ElementType type_facet = mesh.getFacetElementType(type); UInt nb_facet = mesh_facets.getNbElement(type_facet); Vector facet_insertion; Real * bary_facet = new Real[spatial_dimension]; for (UInt f = 0; f < nb_facet; ++f) { mesh_facets.getBarycenter(f, type_facet, bary_facet); - if (bary_facet[1] == 1) facet_insertion.push_back(f); + if (bary_facet[1] == 1){ + facet_insertion.push_back(f);} } delete[] bary_facet; model.insertCohesiveElements(facet_insertion); /// boundary conditions Vector & boundary = model.getBoundary(); UInt nb_nodes = mesh.getNbNodes(); UInt nb_element = mesh.getNbElement(type); Vector &position = const_cast &> (mesh.getNodes()); Vector & displacement = model.getDisplacement(); for (UInt n = 0; n < nb_nodes; ++n) { if (std::abs(position(n,1))< Math::getTolerance()){ boundary(n, 1) = true; displacement(n,1) = 0.0; } if ((std::abs(position(n,0))< Math::getTolerance())&& (position(n,1)< 1.1)){ boundary(n, 0) = true; displacement(n,0) = 0.0; } if ((std::abs(position(n,0)-1)< Math::getTolerance())&&(std::abs(position(n,1)-1)< Math::getTolerance())){ boundary(n, 0) = true; displacement(n,0) = 0.0; } if (std::abs(position(n,1)-2)< Math::getTolerance()){ boundary(n, 1) = true; } } #if defined(AKANTU_USE_IOHELPER) iohelper::ElemType paraview_type = iohelper::TRIANGLE2; /// initialize the paraview output iohelper::DumperParaview dumper; dumper.SetPoints(model.getFEM().getMesh().getNodes().values, spatial_dimension, nb_nodes, "implicit"); dumper.SetConnectivity((int *)model.getFEM().getMesh().getConnectivity(type).values, paraview_type, nb_element, iohelper::C_MODE); dumper.AddNodeDataField(model.getDisplacement().values, spatial_dimension, "displacements"); dumper.AddNodeDataField(model.getVelocity().values, spatial_dimension, "velocity"); dumper.AddNodeDataField(model.getAcceleration().values, spatial_dimension, "acceleration"); dumper.AddNodeDataField(model.getForce().values, spatial_dimension, "applied_force"); dumper.AddNodeDataField(model.getResidual().values, spatial_dimension, "forces"); dumper.AddElemDataField(model.getMaterial(0).getStrain(type).values, spatial_dimension*spatial_dimension, "strain"); dumper.AddElemDataField(model.getMaterial(0).getStress(type).values, spatial_dimension*spatial_dimension, "stress"); dumper.SetEmbeddedValue("displacements", 1); dumper.SetEmbeddedValue("applied_force", 1); dumper.SetEmbeddedValue("forces", 1); dumper.SetPrefix("paraview"); dumper.Init(); dumper.Dump(); #endif const MaterialCohesive & mat_coh = dynamic_cast< const MaterialCohesive &> (model.getMaterial(1)); const Vector & opening = mat_coh.getOpening(_cohesive_2d_6); //const Vector & traction = mat_coh.getTraction(_cohesive_2d_6); - //model.getStiffnessMatrix().saveMatrix("K.mtx"); model.updateResidual(); const Vector & residual = model.getResidual(); UInt max_step = 1000; Real increment = 3./max_step; Real error_tol = 10e-6; std::ofstream fout; fout.open("output"); /// Main loop for ( UInt nstep = 0; nstep < max_step; ++nstep){ Real norm = 10; UInt count = 0; for (UInt n = 0; n < nb_nodes; ++n) { if (std::abs(position(n,1)-2)< Math::getTolerance()){ displacement(n,1) += increment; } } do{ std::cout << "Iter : " << ++count << " - residual norm : " << norm << std::endl; model.assembleStiffnessMatrix(); - model.getStiffnessMatrix().saveMatrix("K.mtx"); + if ((nstep == 0)&&(count == 1)) + model.getStiffnessMatrix().saveMatrix("stiffness_matrix.lastout"); + model.solveStatic(); model.updateResidual(); } while(!model.testConvergenceResidual(1e-5, norm) && (count < 100)) ; std::cout << "Step : " << nstep << " - residual norm : " << norm << std::endl; #if defined(AKANTU_USE_IOHELPER) dumper.Dump(); #endif Real resid = 0; for (UInt n = 0; n < nb_nodes; ++n) { if (std::abs(position(n,1)-2)< Math::getTolerance()){ resid += residual.values[spatial_dimension* n + 1]; } } Real analytical = exp(1) * std::abs(opening.values[3]) * exp (-std::abs(opening.values[3])/0.5)/0.5; - // the residual force is comparing with the theoretical value of the cohesive law + //the residual force is comparing with the theoretical value of the cohesive law error_tol = std::abs((- resid - analytical)/analytical); fout << nstep << " " << -resid << " " << analytical << " " << error_tol << std::endl; if (error_tol > 1e-4) return EXIT_FAILURE; - } + } fout.close(); // finalize(); return EXIT_SUCCESS; }