Skip to content

Commit

Permalink
Bugfix Williams stress field
Browse files Browse the repository at this point in the history
  • Loading branch information
melc-da committed Apr 18, 2024
2 parents 81d4916 + 758baf6 commit b637bc1
Show file tree
Hide file tree
Showing 18 changed files with 1,078 additions and 288 deletions.
4 changes: 2 additions & 2 deletions CITATION.cff
Original file line number Diff line number Diff line change
Expand Up @@ -18,7 +18,7 @@ authors:
given-names: "Eric"
orcid: "https://orcid.org/0000-0002-3479-9143"
title: "Crack Analysis Tool in Python - CrackPy"
version: 1.2.0
version: 1.2.1
doi: tbd
date-released: 2024-04-03
date-released: 2024-04-18
url: "https://github.com/dlr-wf/crackpy"
2 changes: 1 addition & 1 deletion crackpy/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -5,4 +5,4 @@
import crackpy.structure_elements

# package information
__version__ = "1.2.0"
__version__ = "1.2.1"
4 changes: 3 additions & 1 deletion crackpy/crack_detection/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,5 +2,7 @@
import crackpy.crack_detection.deep_learning
import crackpy.crack_detection.pipeline
import crackpy.crack_detection.utils
import crackpy.crack_detection.model
import crackpy.crack_detection.correction
import crackpy.crack_detection.detection
import crackpy.crack_detection.line_intercept
import crackpy.crack_detection.model
2 changes: 1 addition & 1 deletion crackpy/crack_detection/correction.py
Original file line number Diff line number Diff line change
Expand Up @@ -85,7 +85,7 @@ def correct_crack_tip(
max_iter: maximum number of iterations
step_tol: tolerance for the step :math:`\\Delta x`
damper: damper for the step size
method: 'rethore' or 'symbolic_regression'
method: 'rethore', 'symbolic_regression' or 'custom_function'
verbose: If True, print the current iteration
d_x_str: If method='custom_function', provide a function for the correction in x as a string
d_y_str: If method='custom_function', provide a function for the correction in y as a string
Expand Down
28 changes: 11 additions & 17 deletions crackpy/dic_control/utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -167,35 +167,31 @@ def _check_rbmc(self):

def _check_surface_components(self):
"""Checks if all necessary data are calculated."""
gom_value_elements = self.project.actual_elements.filter('type', 'surface_component')
value_element_name = gom_value_elements[0].get('name')

gom_actual_value_elements = self.project.actual_elements.filter('type', 'inspection_surface_componentt')
actual_value_elements = []
for gom_element in gom_actual_value_elements:
actual_value_elements.append(gom_element.get('name'))
gom_surface_component_elements = self.project.actual_elements.filter('type', 'surface_component')
surface_component_element_name = gom_surface_component_elements[0].get('name')

actual_surf_elements = []
for gom_element in self.project.inspection:
actual_surf_elements.append(gom_element.get('name'))

for result_string in self.result_types.keys():
if value_element_name + result_string not in actual_surf_elements:
if surface_component_element_name + result_string not in actual_surf_elements:
if self.result_types[result_string] == "displacement":
distance_restriction = self.disp_directions[result_string]
_ = self.script.inspection.inspect_dimension(
elements=[self.project.actual_elements[value_element_name]],
elements=[self.project.actual_elements[surface_component_element_name]],
distance_restriction=distance_restriction,
nominal_value=0.0,
nominal_value_source='fixed_value',
type=self.result_types[result_string])
else:
_ = self.script.inspection.inspect_dimension(
elements=[self.project.actual_elements[value_element_name]],
elements=[self.project.actual_elements[surface_component_element_name]],
nominal_value=0.0,
nominal_value_source='fixed_value',
type=self.result_types[result_string])
print(f"Creating surface element '{value_element_name + result_string}' against nominal value = 0.0.")
print(f"Creating surface element '{surface_component_element_name + result_string}' "
f"against nominal value = 0.0.")
print("Recalculating...")
self.script.sys.recalculate_project(with_reports=False)
print("...done.")
Expand Down Expand Up @@ -425,12 +421,13 @@ def gather_stage_process_data(self):
# Crack Tip Detection
def setup_cracktip_detection(self, side: str, interp_size: float, offset: tuple, find_path: bool = False,
crack_path_radius: float or None = 50, export_folder: str or None = None):
"""Function to setup the parameters for crack tip detection method.
"""Function to set up the parameters for crack tip detection method.
Args:
side: can either be 'left' or 'right' defining the specimen side
interp_size: size of the interpolation region in mm
offset: tuple of positions of left edge and vertical offset of interpolation region in mm
find_path: if True, crack path will be detected
crack_path_radius: radius w.r.t. crack tip used to calculate the crack path angle in pixels
export_folder: if given, indicates the export folder relative path
Expand Down Expand Up @@ -461,7 +458,6 @@ def detect_actual_crack_tip(self) -> CrackTipInfo:

result_dic = self.get_result_dict(current_stage_index=self._get_current_stage_indx())


data.set_data_manually(
coor_x=result_dic["x_undef"],
coor_y=result_dic["y_undef"],
Expand Down Expand Up @@ -557,7 +553,6 @@ def calc_actual_sifs(self, ct_info: CrackTipInfo, int_props: IntegralProperties

result_dic = self.get_result_dict(current_stage_index=self._get_current_stage_indx())

print("calc_actual_sifs: eps_x and eps_y in [%], eps_xy in [1]!")
data.set_data_manually(
coor_x=result_dic["x_undef"],
coor_y=result_dic["y_undef"],
Expand Down Expand Up @@ -909,7 +904,7 @@ def __init__(self, side: str, interp_size: float, find_path: bool, angle_det_rad
Args:
side: "left" or "right"
interp_size: size in mm of interpolation for the crack tip detection. Choose such that the crack tip
lies inside this region. Moreover, the area should not be larger than he specimen
lies inside this region. Moreover, the area should not be larger than the specimen
angle_det_radius: radius in mm to use to compute the slope of crack path
offset: distance from global origin (i.e. x=0, y=0) to the left edge of interpolation area in mm
export_folder: path to export folder
Expand All @@ -919,14 +914,13 @@ def __init__(self, side: str, interp_size: float, find_path: bool, angle_det_rad
self.interp_size = interp_size
self.offset = offset
self.export_folder = export_folder
device = torch.device('cuda:0' if torch.cuda.is_available() else 'cpu')

self.det = crack_detection.detection.CrackDetection(
side=side,
detection_window_size=interp_size,
offset=offset,
angle_det_radius=angle_det_radius,
device=device)
device="cpu")

# load crack detection models
self.tip_detector = crack_detection.model.get_model('ParallelNets')
Expand Down
39 changes: 38 additions & 1 deletion crackpy/fracture_analysis/crack_tip.py
Original file line number Diff line number Diff line change
Expand Up @@ -35,7 +35,7 @@ def williams_stress_field(a: list or np.array, b: list or np.array, terms: list
+ (n/2 - 1) * np.sin((n/2 - 3) * phi)))
sigma_xy += n/2 * r**(n/2 - 1) * (a[index] * ((n/2 - 1) * np.sin((n/2 - 3) * phi)
- (n/2 + (-1)**n) * np.sin((n/2 - 1) * phi))
- b[index] * ((n/2 - 1) * np.cos((n/2 - 3) * phi)
+ b[index] * ((n/2 - 1) * np.cos((n/2 - 3) * phi)
- (n/2 - (-1)**n) * np.cos((n/2 - 1) * phi)))
return [sigma_x, sigma_y, sigma_xy]

Expand Down Expand Up @@ -72,6 +72,43 @@ def cjp_displ_field(coeffs: list or np.array, phi: float, r: float, material: Ma
return disp_x, disp_y


def cjp_stress_field(coeffs: list or np.array, phi: float, r: float) -> list:
"""Formula for the stress field around the crack tip in real polar coordinates by means of the **five-parameter CJP model**.
[see formulas 3 in Christopher et al. Extension of the CJP model to mixed mode I and mode II (2013)]
Args:
coeffs: Z = (A_r, B_r, B_i, C, E) as in Christopher et al. '13
phi: angle from polar coordinates [rad]
r: radius from polar coordinates [mm]
Returns:
stresses sigma_x, sigma_y, and sigma_xy
"""
A_r, B_r, B_i, C, E = coeffs
sigma_x = (1 / r ** 0.5) * (
-0.5 * (A_r + 4 * B_r + 8 * E) * np.cos(phi / 2) - 0.5 * B_r * np.cos(5 * phi / 2) + \
0.5 * B_i * (np.sin(5*phi/2) + 7 * np.sin(phi/2)) - \
0.5 * E * (np.log(r) * (np.cos(5 * phi / 2) + 3 * np.cos(phi / 2)) + phi * (
np.sin(5 * phi / 2) + 3 * np.sin(phi / 2)))
) - C

sigma_y = (1 / r ** 0.5) * (
0.5 * (A_r - 4 * B_r - 8 * E) * np.cos(phi / 2) + 0.5 * B_r * np.cos(5 * phi / 2) - \
0.5 * B_i * (np.sin(5*phi/2) - np.sin(phi/2)) + \
0.5 * E * (np.log(r) * (np.cos(5 * phi / 2) - 5 * np.cos(phi / 2)) + phi * (
np.sin(5 * phi / 2) - 5 * np.sin(phi / 2)))
)

sigma_xy = (1 / r ** 0.5) * (
-0.5* (A_r * np.sin(phi / 2) + B_r * np.sin(5 * phi / 2)) +\
0.5 * B_i * (np.cos(5 * phi / 2) + 3 * np.cos(phi / 2)) - \
- E * np.sin(phi) * (np.log(r) * np.cos(3 * phi / 2) + phi * np.sin(3 * phi / 2))
)

return [sigma_x, sigma_y, sigma_xy]


def williams_displ_field(a: list or np.array, b: list or np.array, terms: list or np.array,
phi: float, r: float, material: Material) -> tuple:
"""Formula for the displacement fields around the crack tip in polar coordinates by Williams.
Expand Down
43 changes: 25 additions & 18 deletions crackpy/fracture_analysis/data_processing.py
Original file line number Diff line number Diff line change
Expand Up @@ -241,22 +241,25 @@ def transform_data(self, x_shift: float, y_shift: float, angle: float):
self.eps_xy = np.asarray(eps_xy)

# transformation sig' = R * sig * R^T for each node
sig_x, sig_y, sig_xy = [], [], []
for node_i, _ in enumerate(self.sig_x):
stress_i = np.array([[self.sig_x[node_i], self.sig_xy[node_i]],
[self.sig_xy[node_i], self.sig_y[node_i]]])
stress_i = np.dot(stress_i, trafo_matrix.T)
stress_i = np.dot(trafo_matrix, stress_i)
sig_x.append(stress_i[0, 0])
sig_y.append(stress_i[1, 1])
sig_xy.append(stress_i[0, 1])
self.sig_x = np.asarray(sig_x)
self.sig_y = np.asarray(sig_y)
self.sig_xy = np.asarray(sig_xy)
stress_exists = self.sig_x is not None and self.sig_y is not None and self.sig_xy is not None
if stress_exists:
sig_x, sig_y, sig_xy = [], [], []
for node_i, _ in enumerate(self.sig_x):
stress_i = np.array([[self.sig_x[node_i], self.sig_xy[node_i]],
[self.sig_xy[node_i], self.sig_y[node_i]]])
stress_i = np.dot(stress_i, trafo_matrix.T)
stress_i = np.dot(trafo_matrix, stress_i)
sig_x.append(stress_i[0, 0])
sig_y.append(stress_i[1, 1])
sig_xy.append(stress_i[0, 1])
self.sig_x = np.asarray(sig_x)
self.sig_y = np.asarray(sig_y)
self.sig_xy = np.asarray(sig_xy)

# recalculate Von Mises stress and eq. strain
self.eps_vm = self.calc_eps_vm()
self.sig_vm = self._calc_sig_vm()
if stress_exists:
self.sig_vm = self._calc_sig_vm()

def calc_eps_vm(self):
"""Calculate and return Von Mises equivalent strain.
Expand Down Expand Up @@ -354,12 +357,13 @@ def _renumber_nodes(self, old_node_number: int, mapping: dict):

return new_node_number

def to_vtk(self, output_folder: str = None, alpha: float = 1.0):
def to_vtk(self, output_folder: str = None, metadata: bool = True, alpha: float = 1.0):
"""Returns a vtk file of the DIC data.
Args:
output_folder: path to output folder; If None, no vtk file will be saved. The output file name will be the
same as the input file name with the extension .vtk
metadata: if True, metadata will be added to the vtk file
alpha: alpha value for the Delaunay triangulation (only used if no connectivity data is provided)
Returns:
Expand Down Expand Up @@ -425,13 +429,16 @@ def to_vtk(self, output_folder: str = None, alpha: float = 1.0):
mesh.point_data['eps_2 [%]'] = self.eps_2 * 100.0

# add metadata
self.read_header()
for meta_attribute in self.meta_attributes:
meta_data = getattr(self, meta_attribute)
mesh.field_data[meta_attribute] = [meta_data]
if metadata:
self.read_header()
for meta_attribute in self.meta_attributes:
meta_data = getattr(self, meta_attribute)
mesh.field_data[meta_attribute] = [meta_data]

# write vtk file
if output_folder is not None:
if not os.path.exists(output_folder):
os.makedirs(output_folder)
path = os.path.join(output_folder, self.nodemap_name[:-4] + '.vtk')
mesh.save(path, binary=False)
return mesh
Expand Down
3 changes: 0 additions & 3 deletions crackpy/structure_elements/__init__.py
Original file line number Diff line number Diff line change
@@ -1,5 +1,2 @@
import crackpy.structure_elements.project
import crackpy.structure_elements.specimen
import crackpy.structure_elements.experiment
import crackpy.structure_elements.material
import crackpy.structure_elements.data_files
18 changes: 7 additions & 11 deletions crackpy/structure_elements/data_files.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,13 +2,11 @@


class FileStructure:
"""Basic class for File Structures."""
def __init__(self):
pass


class NodemapStructure(FileStructure):
"""Structure of nodemap files generated by Aramis."""
def __init__(
self,
row_length=10,
Expand Down Expand Up @@ -42,13 +40,12 @@ def __init__(


class DataFile:
"""Data File class for data files which are usually .txt-files with a defined file structure given as the
'structure' argument."""
def __init__(self, name: str, folder: str or os.path or None,
structure: FileStructure or None=None,
def __init__(self, name: str, folder: str or os.path = None,
structure: FileStructure = None,
project_name: str = "example project",
specimen_name: str = "specimen name"):
"""Initializes Data File class for MHW data files.
"""Data File class for data files which are usually .txt-files with a defined file structure given as the
'structure' argument.
Args:
name: data file name as string
Expand All @@ -66,19 +63,18 @@ def __init__(self, name: str, folder: str or os.path or None,


class Nodemap(DataFile):
"""Structure of nodemap files generated by Aramis or Ansys."""
def __init__(self, name: str, folder: str or os.path,
structure: NodemapStructure = NodemapStructure(),
project_name: str = "example project",
specimen_name: str = "specimen name"):
"""Initializes Nodemap File class.
"""Nodemap File class.
Args:
name: data file name as string
folder: path to data file
structure: optional given structure for the data file
project_name: optional given project data data refers to (must be given if shepard export required)
specimen_name: optional given specimen name data refers to (must be given if shepard export required)
project_name: optional given project name
specimen_name: optional given specimen name
"""
super().__init__(name, folder, structure, project_name, specimen_name)
Loading

0 comments on commit b637bc1

Please sign in to comment.