Changes to fragmenting cylinder example

ORNL · Jul 24, 2024 · 7ec9ba1 · 7ec9ba1
1 parent 3948cdd
commit 7ec9ba1
Show file tree

Hide file tree

Showing 2 changed files with 34 additions and 46 deletions.
diff --git a/examples/mechanics/fragmenting_cylinder.cpp b/examples/mechanics/fragmenting_cylinder.cpp
@@ -38,9 +38,13 @@ void fragmentingCylinderExample( const std::string filename )
     double rho0 = inputs["density"];
     double K = inputs["bulk_modulus"];
     double G = inputs["shear_modulus"]; // Only for LPS.
-    double G0 = inputs["fracture_energy"];
+    double sc = inputs["critical_stretch"];
     double delta = inputs["horizon"];
     delta += 1e-10;
+    double G0 = 9 * K * delta * ( sc * sc ) /
+                5; // For PMB or LPS with influence_type == 1
+    // double G0 = 15*K*delta*(sc*sc)/8;   // For LPS with influence_type == 0
+    // (default)
 
     // ====================================================
     //                  Discretization
@@ -104,46 +108,30 @@ void fragmentingCylinderExample( const std::string filename )
     double dy = particles->dx[1];
     double dz = particles->dx[2];
 
-    double height = inputs["system_size"][0];
+    double height = inputs["system_size"][2];
     auto init_functor = KOKKOS_LAMBDA( const int pid )
     {
         // Density
         rho( pid ) = rho0;
+
         // Perturb particle positions
-        double factor = 0.5;
-        x( pid, 0 ) =
-            x( pid, 0 ) *
-            ( 1 + ( -1 + 2 * ( (double)std::rand() / ( RAND_MAX ) ) ) *
-                      ( factor * dx ) );
-        x( pid, 1 ) =
-            x( pid, 1 ) *
-            ( 1 + ( -1 + 2 * ( (double)std::rand() / ( RAND_MAX ) ) ) *
-                      ( factor * dy ) );
-        x( pid, 2 ) =
-            x( pid, 2 ) *
-            ( 1 + ( -1 + 2 * ( (double)std::rand() / ( RAND_MAX ) ) ) *
-                      ( factor * dz ) );
+        double random_number_x = ( (double)rand() / ( RAND_MAX ) );
+        double random_number_y = ( (double)rand() / ( RAND_MAX ) );
+        double random_number_z = ( (double)rand() / ( RAND_MAX ) );
+
+        double factor = 0.001; // inputs["grid_perturbation_factor"];
+        x( pid, 0 ) = x( pid, 0 ) + ( 2 * random_number_x - 1.0 ) * factor * dx;
+        x( pid, 1 ) = x( pid, 1 ) + ( 2 * random_number_y - 1.0 ) * factor * dy;
+        x( pid, 2 ) = x( pid, 2 ) + ( 2 * random_number_z - 1.0 ) * factor * dz;
+
         // Velocity
-        double zfactor = ( x( pid, 2 ) - zmin ) / ( 0.5 * height ) - 1;
+        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;
-
-        // Perturb nodes
-        double random_number_x = ( (double)rand() / ( RAND_MAX ) );
-        double random_number_y = ( (double)rand() / ( RAND_MAX ) );
-        double random_number_z = ( (double)rand() / ( RAND_MAX ) );
-
-        double magnitude = 0.001;
-        x( pid, 0 ) =
-            x( pid, 0 ) + ( 2 * random_number_x - 1.0 ) * magnitude * dx;
-        x( pid, 1 ) =
-            x( pid, 1 ) + ( 2 * random_number_y - 1.0 ) * magnitude * dy;
-        x( pid, 2 ) =
-            x( pid, 2 ) + ( 2 * random_number_z - 1.0 ) * magnitude * dz;
     };
     particles->updateParticles( exec_space{}, init_functor );
 

diff --git a/examples/mechanics/inputs/fragmenting_cylinder.json b/examples/mechanics/inputs/fragmenting_cylinder.json
@@ -1,19 +1,19 @@
 {
-    "num_cells"              : {"value": [51, 51, 101]},
-    "system_size"            : {"value": [0.05, 0.05, 0.1], "unit": "m"},
-    "density"                : {"value": 7800, "unit": "kg/m^3"},
-    "bulk_modulus"           : {"value": 130e+9, "unit": "Pa"},
-    "shear_modulus"          : {"value": 78e+9, "unit": "Pa"},
-    "fracture_energy"        : {"value": 3.74e+6, "unit": "J/m^2"},
-    "horizon"                : {"value": 4.0e-3, "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": 5e-4, "unit": "s"},
-    "timestep"               : {"value": 2.0e-7, "unit": "s"},
-    "timestep_safety_factor" : {"value": 0.85},
-    "output_frequency"       : {"value": 50},
-    "output_reference"       : {"value": false}
+    "num_cells"                : {"value": [51, 51, 101]},
+    "system_size"              : {"value": [0.05, 0.05, 0.1], "unit": "m"},
+    "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}
 }