diff --git a/README.md b/README.md index 968aafe6..1058d6b2 100644 --- a/README.md +++ b/README.md @@ -117,7 +117,7 @@ Examples which only include mechanics and fracture are with `examples/mechanics` The first example is an elastic wave propagating through a cube from an initial Gaussian radial displacement profile from [1]. Assuming the build paths above, the example can be run with: ``` -./CabanaPD/build/install/bin/ElasticWave CabanaPD/examples/inputs/elastic_wave.json +./CabanaPD/build/install/bin/ElasticWave CabanaPD/examples/mechanics/inputs/elastic_wave.json ``` The second example is the Kalthoff-Winkler experiment [2], where an impactor @@ -125,26 +125,31 @@ causes crack propagation at an angle from two pre-notches on a steel plate. The example can be run with: ``` -./CabanaPD/build/install/bin/KalthoffWinkler CabanaPD/examples/inputs/kalthoff_winkler.json +./CabanaPD/build/install/bin/KalthoffWinkler CabanaPD/examples/mechanics/inputs/kalthoff_winkler.json ``` -The third example is crack branching in a pre-notched soda-lime glass plate due to traction loading [3]. The example can be -run with: +The third example is crack branching in a pre-notched soda-lime glass plate due to traction loading [3]. The example can be run with: ``` -./CabanaPD/build/install/bin/CrackBranching CabanaPD/examples/inputs/crack_branching.json +./CabanaPD/build/install/bin/CrackBranching CabanaPD/examples/mechanics/inputs/crack_branching.json +``` + +The fourth example is a fragmenting cylinder due to internal pressure [4]. The example can be run with: + +``` +./CabanaPD/build/install/bin/FragmentingCylinder CabanaPD/examples/mechanics/inputs/fragmenting_cylinder.json ``` ### Thermomechanics Examples which demonstrate temperature-dependent mechanics and fracture are with `examples/thermomechanics`. -The first example is thermoelastic deformation in a homogeneous plate due to linear thermal loading [4]. The example can be run with: +The first example is thermoelastic deformation in a homogeneous plate due to linear thermal loading [5]. The example can be run with: ``` ./CabanaPD/build/install/bin/ThermalDeformation CabanaPD/examples/thermomechanics/thermal_deformation.json ``` -The second example is crack initiation and propagation in an alumina ceramic plate due to a thermal shock caused by water quenching [5]. The example can be run with: +The second example is crack initiation and propagation in an alumina ceramic plate due to a thermal shock caused by water quenching [6]. The example can be run with: ``` ./CabanaPD/build/install/bin/ThermalCrack CabanaPD/examples/thermomechanics/thermal_crack.json @@ -164,9 +169,11 @@ Kunze, and L.W. Meyer, eds., Vol 1, DGM Informationsgesellschaft Verlag (1988) [3] F. Bobaru and G. Zhang, Why do cracks branch? A peridynamic investigation of dynamic brittle fracture, International Journal of Fracture 196 (2015): 59–98. -[4] D. He, D. Huang, and D. Jiang, Modeling and studies of fracture in functionally graded materials under thermal shock loading using peridynamics, Theoretical and Applied Fracture Mechanics 111 (2021): 102852. +[4] D.J. Littlewood, M.L. Parks, J.T. Foster, J.A. Mitchell, and P. Diehl, The peridigm meshfree peridynamics code, Journal of Peridynamics and Nonlocal Modeling 6 (2024): 118–148. + +[5] D. He, D. Huang, and D. Jiang, Modeling and studies of fracture in functionally graded materials under thermal shock loading using peridynamics, Theoretical and Applied Fracture Mechanics 111 (2021): 102852. -[5] C.P. Jiang, X.F. Wu, J. Li, F. Song, Y.F. Shao, X.H. Xu, and P. Yan, A study of the mechanism of formation and numerical simulations of crack patterns in ceramics subjected to thermal shock, Acta Materialia 60 (2012): 4540–4550. +[6] C.P. Jiang, X.F. Wu, J. Li, F. Song, Y.F. Shao, X.H. Xu, and P. Yan, A study of the mechanism of formation and numerical simulations of crack patterns in ceramics subjected to thermal shock, Acta Materialia 60 (2012): 4540–4550. ## Contributing diff --git a/examples/mechanics/CMakeLists.txt b/examples/mechanics/CMakeLists.txt index df1b8e15..a10b05f1 100644 --- a/examples/mechanics/CMakeLists.txt +++ b/examples/mechanics/CMakeLists.txt @@ -7,4 +7,7 @@ target_link_libraries(KalthoffWinkler LINK_PUBLIC CabanaPD) add_executable(CrackBranching crack_branching.cpp) target_link_libraries(CrackBranching LINK_PUBLIC CabanaPD) -install(TARGETS ElasticWave KalthoffWinkler CrackBranching DESTINATION ${CMAKE_INSTALL_BINDIR}) +add_executable(FragmentingCylinder fragmenting_cylinder.cpp) +target_link_libraries(FragmentingCylinder LINK_PUBLIC CabanaPD) + +install(TARGETS ElasticWave KalthoffWinkler CrackBranching FragmentingCylinder DESTINATION ${CMAKE_INSTALL_BINDIR}) diff --git a/examples/mechanics/fragmenting_cylinder.cpp b/examples/mechanics/fragmenting_cylinder.cpp new file mode 100644 index 00000000..8c245cc4 --- /dev/null +++ b/examples/mechanics/fragmenting_cylinder.cpp @@ -0,0 +1,165 @@ +/**************************************************************************** + * Copyright (c) 2022-2023 by Oak Ridge National Laboratory * + * All rights reserved. * + * * + * This file is part of CabanaPD. CabanaPD is distributed under a * + * BSD 3-clause license. For the licensing terms see the LICENSE file in * + * the top-level directory. * + * * + * SPDX-License-Identifier: BSD-3-Clause * + ****************************************************************************/ + +#include +#include + +#include "mpi.h" + +#include + +#include + +// Simulate an expanding cylinder resulting in fragmentation +void fragmentingCylinderExample( const std::string filename ) +{ + // ==================================================== + // Use default Kokkos spaces + // ==================================================== + using exec_space = Kokkos::DefaultExecutionSpace; + using memory_space = typename exec_space::memory_space; + + // ==================================================== + // Read inputs + // ==================================================== + CabanaPD::Inputs inputs( filename ); + + // ==================================================== + // Material parameters + // ==================================================== + double rho0 = inputs["density"]; + double K = inputs["bulk_modulus"]; + // double G = inputs["shear_modulus"]; // Only for LPS. + double sc = inputs["critical_stretch"]; + double delta = inputs["horizon"]; + delta += 1e-10; + // For PMB or LPS with influence_type == 1 + double G0 = 9 * K * delta * ( sc * sc ) / 5; + // For LPS with influence_type == 0 (default) + // double G0 = 15 * K * delta * ( sc * sc ) / 8; + + // ==================================================== + // Discretization + // ==================================================== + // FIXME: set halo width based on delta + std::array low_corner = inputs["low_corner"]; + std::array high_corner = inputs["high_corner"]; + std::array num_cells = inputs["num_cells"]; + int m = std::floor( delta / + ( ( high_corner[0] - low_corner[0] ) / num_cells[0] ) ); + int halo_width = m + 1; // Just to be safe. + + // ==================================================== + // Force model + // ==================================================== + using model_type = CabanaPD::ForceModel; + model_type force_model( delta, K, G0 ); + // using model_type = + // CabanaPD::ForceModel; + // model_type force_model( delta, K, G, G0 ); + + // ==================================================== + // Particle generation + // ==================================================== + // Does not set displacements, velocities, etc. + auto particles = std::make_shared< + CabanaPD::Particles>( + exec_space(), low_corner, high_corner, num_cells, halo_width ); + + // ==================================================== + // Custom particle initialization + // ==================================================== + double x_center = 0.5 * ( low_corner[0] + high_corner[0] ); + double y_center = 0.5 * ( low_corner[1] + high_corner[1] ); + double Rout = inputs["cylinder_outer_radius"]; + double Rin = inputs["cylinder_inner_radius"]; + + // Do not create particles outside given cylindrical region + auto init_op = KOKKOS_LAMBDA( const int, const double x[3] ) + { + double rsq = ( x[0] - x_center ) * ( x[0] - x_center ) + + ( x[1] - y_center ) * ( x[1] - y_center ); + if ( rsq < Rin * Rin || rsq > Rout * Rout ) + return false; + return true; + }; + particles->createParticles( exec_space(), init_op ); + + auto rho = particles->sliceDensity(); + auto x = particles->sliceReferencePosition(); + auto v = particles->sliceVelocity(); + auto f = particles->sliceForce(); + + double vrmax = inputs["max_radial_velocity"]; + double vrmin = inputs["min_radial_velocity"]; + double vzmax = inputs["max_vertical_velocity"]; + double zmin = low_corner[2]; + + auto dx = particles->dx; + double height = inputs["system_size"][2]; + double factor = inputs["grid_perturbation_factor"]; + + using pool_type = Kokkos::Random_XorShift64_Pool; + using random_type = Kokkos::Random_XorShift64; + pool_type pool; + int seed = 456854; + pool.init( seed, particles->n_local ); + auto init_functor = KOKKOS_LAMBDA( const int pid ) + { + // Density + rho( pid ) = rho0; + + // Perturb particle positions + auto gen = pool.get_state(); + for ( std::size_t d = 0; d < 3; d++ ) + { + auto rand = + Kokkos::rand::draw( gen, 0.0, 1.0 ); + x( pid, d ) += ( 2.0 * rand - 1.0 ) * factor * dx[d]; + } + pool.free_state( gen ); + + // Velocity + double zfactor = ( ( x( pid, 2 ) - zmin ) / ( 0.5 * height ) ) - 1; + double vr = vrmax - vrmin * zfactor * zfactor; + v( pid, 0 ) = + vr * Kokkos::cos( Kokkos::atan2( x( pid, 1 ), x( pid, 0 ) ) ); + v( pid, 1 ) = + vr * Kokkos::sin( Kokkos::atan2( x( pid, 1 ), x( pid, 0 ) ) ); + v( pid, 2 ) = vzmax * zfactor; + }; + particles->updateParticles( exec_space{}, init_functor ); + + // ==================================================== + // Simulation run + // ==================================================== + + // Define empty pre-notch + CabanaPD::Prenotch<0> prenotch; + + auto cabana_pd = CabanaPD::createSolverFracture( + inputs, particles, force_model, prenotch ); + cabana_pd->init(); + cabana_pd->run(); +} + +// Initialize MPI+Kokkos. +int main( int argc, char* argv[] ) +{ + MPI_Init( &argc, &argv ); + Kokkos::initialize( argc, argv ); + + fragmentingCylinderExample( argv[1] ); + + Kokkos::finalize(); + MPI_Finalize(); +} diff --git a/examples/mechanics/inputs/fragmenting_cylinder.json b/examples/mechanics/inputs/fragmenting_cylinder.json new file mode 100644 index 00000000..a6b7041b --- /dev/null +++ b/examples/mechanics/inputs/fragmenting_cylinder.json @@ -0,0 +1,20 @@ +{ + "num_cells" : {"value": [51, 51, 101]}, + "system_size" : {"value": [0.05, 0.05, 0.1], "unit": "m"}, + "grid_perturbation_factor" : {"value": 0.1}, + "density" : {"value": 7800, "unit": "kg/m^3"}, + "bulk_modulus" : {"value": 130e+9, "unit": "Pa"}, + "shear_modulus" : {"value": 78e+9, "unit": "Pa"}, + "critical_stretch" : {"value": 0.02}, + "horizon" : {"value": 0.00417462, "unit": "m"}, + "cylinder_outer_radius" : {"value": 0.025, "unit": "m"}, + "cylinder_inner_radius" : {"value": 0.02, "unit": "m"}, + "max_radial_velocity" : {"value": 200, "unit": "m/s"}, + "min_radial_velocity" : {"value": 50, "unit": "m/s"}, + "max_vertical_velocity" : {"value": 100, "unit": "m/s"}, + "final_time" : {"value": 2.5e-4, "unit": "s"}, + "timestep" : {"value": 1.7e-07, "unit": "s"}, + "timestep_safety_factor" : {"value": 0.70}, + "output_frequency" : {"value": 50}, + "output_reference" : {"value": false} +} diff --git a/src/CabanaPD_Particles.hpp b/src/CabanaPD_Particles.hpp index c1151770..c9e9816d 100644 --- a/src/CabanaPD_Particles.hpp +++ b/src/CabanaPD_Particles.hpp @@ -137,7 +137,7 @@ class Particles std::shared_ptr< Cabana::Grid::LocalGrid>> local_grid; - double dx[dim]; + Kokkos::Array dx; int halo_width;