WIP: Return bond_atoms from crack tip finder and use for kick

libAtoms · Mar 12, 2024 · 7456502 · 7456502
1 parent b28dafc
commit 7456502
Show file tree

Hide file tree

Showing 2 changed files with 6 additions and 22 deletions.
diff --git a/matscipy/fracture_mechanics/thin_strip_utils.py b/matscipy/fracture_mechanics/thin_strip_utils.py
@@ -694,7 +694,7 @@ def find_strip_crack_tip(self,final_crack_state,bondlength,bulk_nn,step_tolerant
 
         print(f'Found crack tip at position {tip_pos}')
 
-        return tip_pos[0], tip_pos_y
+        return tip_pos[0], tip_pos_y[0], bond_atoms
 
 def write_potential_and_buffer(atoms,lammps_filename):
     #get the potential and buffer array from atoms

diff --git a/scripts/fracture_mechanics/thin_strip_lammps.py b/scripts/fracture_mechanics/thin_strip_lammps.py
@@ -286,26 +286,10 @@
             plt.close()
 
             # ------------- find crack tip ----------- #
-            # tmp_bondlength = bondlength
-            # found_tip=False
-            # for j in range(10):
-            #     try:
-            #         # print(tmp_bondlength,bulk_nn)
-            #         bond_atoms = find_tip_coordination(final_crack_state,bondlength=tmp_bondlength,bulk_nn=bulk_nn,calculate_midpoint=True)
-            #         tip_pos = (final_crack_state.get_positions()[bond_atoms,:][:,0])
-            #         #print(tip_pos[0],tip_pos[1])
-            #         assert np.abs(tip_pos[0]-tip_pos[1]) < 1
-            #         found_tip=True
-            #         break
-            #     except AssertionError:
-            #         tmp_bondlength += 0.01
-            #         #keep trying till a crack tip is found
-            # if not found_tip:
-            #     raise RuntimeError('Lost crack tip!')
-            # tip_pos = tip_pos[0]
-            # print(f'Found crack tip at position {tip_pos}')
-            tip_pos,tip_pos_y = tsb.find_strip_crack_tip(final_crack_state,bondlength,bulk_nn,step_tolerant=step_tolerant)
+
+            tip_pos,tip_pos_y,bond_atoms = tsb.find_strip_crack_tip(final_crack_state,bondlength,bulk_nn,step_tolerant=step_tolerant)
             print('y pos', tip_pos_y)
+
             # ------------check for an arrested crack ------------ #
             #only start checking after the first initial_damping steps
             if (i >= initial_damping_time) or (not initial_damp):
@@ -340,8 +324,8 @@
                 simple_strip = tsb.build_thin_strip(strip_width,strip_height,strip_thickness,vacuum)
                 x_pos = simple_strip.get_positions()[:,0]
                 y_pos = simple_strip.get_positions()[:,1]
-                tip_x_pos = tip_pos
-                tip_y_pos = tip_pos_y
+                tip_x_pos = x_pos[bond_atoms[0]]
+                tip_y_pos = [y_pos[bond_atoms[0]], y_pos[bond_atoms[1]]]
                 full_mask = (x_pos>(tip_x_pos-10)) & (x_pos<(tip_x_pos+40))
 
                 #mask is full mask and atoms which have a y position equal to that of the crack tip within numerical resolution