Working version with random pre-notch positions and orientations on x…

…y-plane

examples/random_cracks.cpp

@@ -48,8 +48,8 @@ int main( int argc, char* argv[] )
         std::array<double, 3> high_corner = { high_x, high_y, high_z };
 
         // Time
-        // double t_final = 43e-6;
-        double t_final = 50e-8;
+        double t_final = 43e-6;
+        // double t_final = 50e-8;
         double dt = 5e-8;
         double output_frequency = 5;
         bool output_reference = true;
@@ -72,20 +72,18 @@ int main( int argc, char* argv[] )
         // Multiple random pre-notches
 
         // Number of pre-notches
-        constexpr int Npn = 10;
+        constexpr int Npn = 200;
 
         // Minimum and maximum pre-notch length
-        double minl = 0.005;
-        double maxl = 0.01;
+        double minl = 0.001;
+        double maxl = 0.0025;
 
         // Initialize pre-notch arrays
         Kokkos::Array<Kokkos::Array<double, 3>, Npn> notch_positions;
         Kokkos::Array<Kokkos::Array<double, 3>, Npn> notch_v1;
         Kokkos::Array<Kokkos::Array<double, 3>, Npn> notch_v2;
 
-        // Fixed v1 and v2 vectors (temporatily)
-        double L_prenotch = height / 2.0;
-        Kokkos::Array<double, 3> v1 = { L_prenotch, 0, 0 };
+        // Fixed pre-notch v2 vector
         Kokkos::Array<double, 3> v2 = { 0, 0, thickness };
 
         // Loop over pre-notches
@@ -94,8 +92,6 @@ int main( int argc, char* argv[] )
             // Random numbers for pre-notch position
             double random_number_x = std::rand() / double( RAND_MAX );
             double random_number_y = std::rand() / double( RAND_MAX );
-            // std::cout << n  << " : rand 1: " << random_number_x << " : rand
-            // 2: " << random_number_y << " Pi: " << M_PI << std::endl;
 
             // Coordinates of one endpoint of the pre-notch (random)
             // Note: the addition and subtraction of "maxl" ensures the prenotch
@@ -107,37 +103,26 @@ int main( int argc, char* argv[] )
                 ( low_y + maxl ) +
                 ( ( high_y - maxl ) - ( low_y + maxl ) ) * random_number_y;
             Kokkos::Array<double, 3> p0 = { Xc1, Yc1, low_z };
-            // std::cout << n  << " " << Xc1 << " " << Yc1 << " " << low_z <<
-            // std::endl;
 
             // Assign pre-notch position
             notch_positions[n] = p0;
 
-            //  Pre-notch length on x-y plane
+            //  Pre-notch length on xy-plane
             double random_number_l = std::rand() / double( RAND_MAX );
             double l = minl + ( maxl - minl ) * random_number_l;
 
-            // Random orientation
+            // Random orientation on xy-plane
             double random_number_theta = std::rand() / double( RAND_MAX );
             double theta = M_PI * random_number_theta;
 
-            // std::cout << n  << " : rand 3: " << random_number_l << " : rand
-            // 4: " << random_number_theta << std::endl; std::cout << "l = " <<
-            // " " << l << " ; theta = " << 180*theta/M_PI << " deg." <<
-            // std::endl; std::cout << "l*cos(theta) = " << " " << l * cos(
-            // theta ) << " ; l * sin( theta) = " << l * sin( theta) << " " <<
-            // std::endl;
-
-            // Assign pre-notch length
+            // Assign pre-notch v1 vector
             Kokkos::Array<double, 3> v1_temp = { l * cos( theta ),
                                                  l * sin( theta ), 0 };
             notch_v1[n] = v1_temp;
+            notch_v2[n] = v2;
         }
 
-        // This line needs to be extended to have a list of v1, so that each
-        // pre-notch has a different orientation
-        CabanaPD::Prenotch<Npn> prenotch( v1, v2, notch_positions );
-        // CabanaPD::Prenotch<Npn> prenotch( notch_v1, v2, notch_positions );
+        CabanaPD::Prenotch<Npn> prenotch( notch_v1, notch_v2, notch_positions );
 
         // Choose force model type.
         using model_type =