Post Processing#

Reference result#

ReferenceResult#

This module provides classes for handling and post-processing PhaseFieldX simulation results.

class phasefieldx.PostProcessing.ReferenceResult.AllResults(folder_path)[source]#

Bases: object

Class to handle PhaseFieldX results from a given folder.

This class automatically loads and organizes various types of result files including convergence data, energy files, reaction forces, degrees of freedom, input parameters, and ParaView visualization files from a PhaseFieldX simulation output directory.

Parameters:
folder_pathstr

Path to the folder containing PhaseFieldX results. The folder should contain various output files with extensions like .conv, .energy, .reaction, .dof, .input, and potentially ParaView solution folders.

Attributes:
folder_pathstr

Path to the results folder.

inputdict

Dictionary containing input parameters from .input files.

convergence_filesdict

Dictionary of pandas DataFrames containing convergence data from .conv files.

energy_filesdict

Dictionary of pandas DataFrames containing energy data from .energy files.

reaction_filesdict

Dictionary of pandas DataFrames containing reaction force data from .reaction files.

dof_filesdict

Dictionary of pandas DataFrames containing degree of freedom data from .dof files.

file_nameslist

List of all file names found in the results folder.

folder_nameslist

List of all subdirectory names found in the results folder.

paraview_filenameslist, optional

List of ParaView solution file names, if available.

labelstr

Label for plotting and identification purposes. Default is ‘label’.

colorstr

Color specification for plotting. Default is ‘k’ (black).

paraviewParaviewResult, optional

ParaviewResult object for handling visualization data.

Methods

get_paraview(path)

Initialize ParaView result handling for the specified path.

set_label(label)

Set the label for this result set.

set_color(color)

Set the color for plotting this result set.

Notes

The class automatically categorizes files based on their extensions:

  • .conv files: Convergence information

  • .energy files: Energy-related data

  • .reaction files: Reaction force data

  • .dof files: Degrees of freedom data

  • .input files: Input parameters (converted to dictionary format)

ParaView visualization files are detected if a “paraview-solutions_vtu” subdirectory exists in the results folder.

Examples

>>> # Load results from a simulation folder
>>> results = AllResults("./simulation_results")
>>> 
>>> # Set label and color for plotting
>>> results.set_label("My Simulation")
>>> results.set_color("red")
>>> 
>>> # Access energy data
>>> energy_data = results.energy_files['total.energy']
>>> print(energy_data.head())
>>> 
>>> # Access input parameters
>>> parameters = results.input['parameters.input']
>>> print(f"Length scale: {parameters.get('l', 'N/A')}")
get_paraview(path)[source]#

Initialize ParaView result handling for visualization.

Parameters:
pathstr

Path to the ParaView solution files or results folder containing ParaView visualization data.

Notes

This method creates a ParaviewResult object that can be used for post-processing and visualization of the simulation results.

set_color(color)[source]#

Set the color for plotting this result set.

Parameters:
colorstr

Color specification for matplotlib plotting. Can be a single character (‘r’, ‘g’, ‘b’, ‘k’, etc.), color name (‘red’, ‘blue’, etc.), or hex color code (‘#FF0000’, etc.).

Examples

>>> results = AllResults("./simulation_results")
>>> results.set_color("red")
>>> # or
>>> results.set_color("#FF0000")
>>> # or
>>> results.set_color("r")
set_label(label)[source]#

Set the label for this result set.

Parameters:
labelstr

Label to identify this result set in plots and analysis. This label will be used in plot legends and other visual representations.

Examples

>>> results = AllResults("./simulation_results")
>>> results.set_label("High Resolution Mesh")

Paraview#

Paraview Result Processing#

This module provides a class for loading and processing results from Paraview solution files. It is designed to read VTU files generated by a simulation and organize the data for further analysis or visualization.

class phasefieldx.PostProcessing.paraview.ParaviewResult(folder_path, steps)[source]#

Bases: object

A class to load and process Paraview solution files.

Parameters:
folder_pathstr

Path to the folder containing the Paraview solution files.

stepslist of int

List of time steps to load the corresponding solution files.

Attributes:
folder_pathstr

Path to the folder containing the Paraview solution files.

paraview_folderstr

Path to the folder where Paraview solution files are stored.

stepslist of int

List of time steps to load the corresponding solution files.

datanumpy.ndarray

Array of PyVista datasets loaded from the solution files.

filenumpy.ndarray

Array of file paths to the loaded solution files.

Measure crack length 2D#

Measure crack length in 2D phase-field simulations#

This script provides a post-processing workflow for measuring crack length in 2D phase-field fracture simulations, following the procedure presented in Castillón et al.[1]. It is designed to work with simulation outputs in VTU format (e.g., from dolfinx or FEniCSx), where the phase-field variable (‘phi’) indicates the presence of a crack.

Workflow overview#

  1. VTU to PNG conversion: For each .vtu file in the simulation output folder, the mesh is visualized,

    and the crack region (where phi > 0.95) is highlighted in red. The resulting image is saved as a PNG for further processing.

  2. Crack skeletonization and length measurement: Each PNG image is processed to extract the red crack

    region, skeletonize it (reducing it to a 1-pixel-wide centerline), and count the number of pixels in the skeleton. The pixel count is converted to a physical length using the known specimen size and image resolution.

  3. Saving results: The measured crack lengths for all time steps are saved to a text file

    (crack_lengths.txt) with a column header ‘a’ for easy loading into pandas or other analysis tools.

  4. Visualization: Optionally, a GIF is generated showing the crack evolution, with the measured crack

    length overlaid on each frame.

Requirements#

  • pyvista

  • matplotlib

  • scikit-image

  • numpy

  • imageio

Usage example#

Set your simulation folder and specimen size, then call:

folder_simulation = “path/to/your/simulation/folder” measure_crack(folder_simulation, physical_horizontal=1.0, physical_vertical=1.0)

This will generate PNG images, measure crack lengths, and create a GIF in the specified folder.

Functions#

  • get_black_bbox(image): Finds the bounding box of the specimen in an image.

  • generate_crack_images(input_folder, output_folder): Converts VTU files to PNG images with cracks highlighted.

  • measure_crack_lengths(output_folder, dx, dy): Measures crack lengths from PNG images and saves results.

  • generate_crack_gif(output_folder, crack_lengths, gif_name): Creates a GIF showing crack evolution.

  • measure_crack(folder_simulation, physical_horizontal, physical_vertical): Main workflow function.

phasefieldx.PostProcessing.measure_crack_length_2d.generate_crack_gif(output_folder, crack_lengths, gif_name='crack_evolution.gif')[source]#

Generate a GIF showing crack evolution, overlaying the measured crack length on each frame.

Parameters:
output_folderstr

Path to the folder containing PNG images.

crack_lengthsarray-like

Crack lengths to overlay on each frame.

gif_namestr, optional

Name of the output GIF file (default: “crack_evolution.gif”).

phasefieldx.PostProcessing.measure_crack_length_2d.generate_crack_images(input_folder, output_folder)[source]#

Generate PNG images from .vtu files, highlighting the crack region in red.

For each .vtu file in the input folder, this function: - Loads the mesh and phase-field data. - Thresholds the ‘phi’ field to extract the crack region. - Plots the mesh in black and the crack in red. - Saves the result as a PNG image in the output folder.

Parameters:
input_folderstr

Path to the folder containing .vtu files.

output_folderstr

Path to the folder where PNG images will be saved.

phasefieldx.PostProcessing.measure_crack_length_2d.get_black_bbox(image)[source]#

Find the bounding box of the black (background/specimen) region in an image.

Parameters:
imagenp.ndarray

Input image (grayscale or RGB).

Returns:
min_row, max_row, min_col, max_colint

Indices defining the bounding box of the black region.

phasefieldx.PostProcessing.measure_crack_length_2d.measure_crack(folder_simulation, physical_horizontal, physical_vertical)[source]#

Main workflow to measure crack lengths and generate images and GIF.

This function: - Generates PNG images from VTU files, highlighting cracks. - Detects the specimen bounding box and computes pixel size. - Measures crack lengths from images and saves results. - Generates a GIF showing crack evolution.

Parameters:
folder_simulationstr

Path to the simulation folder containing ‘paraview-solutions_vtu’.

physical_horizontalfloat

Physical width of the specimen (used to compute pixel size).

physical_verticalfloat

Physical height of the specimen (used to compute pixel size).

phasefieldx.PostProcessing.measure_crack_length_2d.measure_crack_lengths(output_folder, dx, dy)[source]#

Measure crack lengths from PNG images and save results to a text file.

For each PNG image in the output folder: - Converts the image to HSV and thresholds for red pixels (crack). - Skeletonizes the crack region to obtain a 1-pixel wide centerline. - Counts the number of skeleton pixels and converts to physical length.

Parameters:
output_folderstr

Path to the folder containing PNG images.

dxfloat

Physical size of a pixel in the horizontal direction.

dyfloat

Physical size of a pixel in the vertical direction (not used, but available).

Returns:
crack_lengthsnp.ndarray

Array of measured crack lengths for each image.