Python 2D DDD-Crack Prototype

This directory contains a small Python prototype for fixed-crack strong coupling between two-dimensional discrete dislocation dynamics and a regularized crack COD solve.

The code is intentionally a research scaffold. It keeps the same broad workflow as the existing MATLAB 2dDDD code and the fixed-crack ParaDiS branch:

1. compute DDD forces and glide-projected velocities;
2. solve a fixed crack COD problem from the current DDD geometry and applied stress;
3. feed crack feedback forces back into DDD motion;
4. write accepted-state DDD and crack diagnostics.

Requirements

• Python 3.10+
• NumPy
• Matplotlib only for plotting commands

No YAML dependency is required. Cases use JSON plus CSV.

Run The Example

From the repository root:

cd py2d_crack_ddd
python -m crackddd2d run cases/fixed_crack_strong/config.json
python -m crackddd2d run cases/fixed_crack_strong/no_crack.json
python -m crackddd2d run cases/fixed_crack_strong/lagged_crack.json

run uses output.root from the JSON file. Passing --out OUT overrides that path for temporary or comparison runs.

Plot a result:

python -m crackddd2d plot-results results/fixed_crack_strong

Create a 2D motion GIF from accepted-state outputs:

python -m crackddd2d animate cases/fixed_crack_strong/config.json results/fixed_crack_strong --out results/fixed_crack_strong/motion.gif

Create a second animation with a sigma_xy stress background:

python -m crackddd2d animate-stress cases/fixed_crack_strong/config.json results/fixed_crack_strong --out results/fixed_crack_strong/stress_xy.gif --grid 80

Scan finite-boundary correction parameters for one case:

python -m crackddd2d scan-boundary cases/fixed_crack_strong/config.json --out results/fixed_crack_strong/boundary_scan.csv --samples 8,12,16,20 --offsets 0.4,0.8,1.6

Run a systematic crack-mesh convergence check for one JSON case:

python -m crackddd2d convergence cases/classic_benchmarks/internal_crack_mode_i_uniform_tension/config.json --out results/internal_mode_i_convergence --n-elements 21,41,81

Li 2025 Fig. 7 style cases use the same run pipeline. They are static-like coupled cases with zero mobility and one accepted time step, and their JSON files enable final-state stress maps:

python -m crackddd2d run cases/li2025_fig7/fig7_phi0_rd30.json
python -m crackddd2d run cases/li2025_fig7/fig7_phi60_rd40.json

For static post-processing or debugging of any single JSON case, field-map and path-sample remain available. Their default field mode is an unbounded superposition:

python -m crackddd2d field-map cases/li2025_fig7/fig7_phi0_rd30.json --out results/li2025_phi0_field --grid 320
python -m crackddd2d path-sample cases/li2025_fig7/fig7_phi0_rd30.json --out results/li2025_phi0_paths --samples 96

Those secondary outputs use total = dislocation + crack + external, with a zero boundary_correction column in the reports. To include the current finite-box traction-free correction explicitly in these static commands, pass --field-mode finite_traction_free. The primary run command does not need that flag: it derives the stress-map field mode from the JSON boundary_correction and periodicity settings. field-map writes field_map.png and field_map_report.json; path-sample writes path_sample.png, path_samples.csv, and path_sample_report.json.

The animations draw the simulation box, boundary condition label, applied stress, standard 2D dislocation glyphs (⊥ for positive b_x, ⊤ for negative b_x), and crack elements colored by the configured animation COD component. Set output.animations.color_crack_by to opening, sliding, magnitude, or none; opening is the default. Set output.animations.crack_color_range: [min, max] to fix the animation crack colorbar range. In animate-stress, the background color is the post-processed total sigma_xy field from DDD dislocations, a regularized mixed-mode crack-COD stress contribution, and the applied stress_xy. The plotting code requests Times New Roman for all text. If the font is not installed on the machine, Matplotlib will fall back to an available font.

Outputs

Each run writes:

• history.csv: step, time, stress, strain, speed, energy, crack residual, total COD norm.
• Coupling diagnostics in history.csv: sliding/opening COD norms, tangential/normal crack traction and residual norms, COD max, crack feedback force norm, optional finite-boundary correction force/source norms, boundary traction residuals, fallback flag, boundary residual norm, total force norm, feedback/total force ratio, coupled iteration count, absolute and relative position/COD/force/traction residuals, total coupled residuals, convergence flag, failure flag, relaxation used, retry count, crack-tip SIF estimates, COD fit error, and nearest dislocation-tip distance.
• dislocations/dislocations_step_*.csv: accepted DDD state.
• crack/crack_step_*.csv: crack element centers, sliding/opening COD, tangential/normal crack-face traction, and tangential/normal residual. When a crack is enabled, crack/crack_step_000000.csv is written after the initial COD solve so animations and diagnostics start from a consistent accepted crack state.
• summary.json: final diagnostics and final positions.
• stress_maps/*.png and stress_maps/stress_map_report.json: optional JSON-driven final-state maps from output.stress_maps. The available field decompositions are dislocation, crack, external, boundary_correction, and total; components are sigma11, sigma22, and sigma12.
• motion.gif: optional 2D motion animation generated by animate or by output.animations.motion.
• stress_xy.gif: optional stress-background animation generated by animate-stress or by output.animations.stress.

Useful output options:

• root: default result directory used by run CONFIG.
• progress: set enabled: true and frequency: N to print a simple timestep progress line every N steps.
• animations: set enabled: true plus motion and/or stress to generate GIFs automatically at the end of run; fps, dpi, stress_grid, motion_file, and stress_file control the output. color_crack_by selects the animated crack overlay (opening, sliding, magnitude, or none), and crack_color_range: [min, max] fixes its colorbar range.
• stress_maps: final static stress-map controls, including grid, fields, components, crack-element overlay, crack coloring mode, and dislocation overlay. Set stress_color_range: [min, max] to override the stress colorbar range, and crack_color_range: [min, max] to override the crack COD colorbar range. Leave either value unset to use the automatic range.

Useful crack options for the DDBEM/HBIE backend:

• element_order: constant or linear_tip_enriched. The latter uses node-based straight-crack displacement discontinuities with constrained crack tips and a distributed equivalent-dislocation density for the traction matrix and field reconstruction.
• tip_enrichment.fit_elements: number of near-tip elements used by the multi-window COD/SIF extrapolation diagnostic.
• near_singular: controls source-line/tip sampling diagnostics for stress reconstruction; tolerance, max_subdivision, and mask_singular_points are reported in stress-field metrics.

Useful convergence options:

• n_elements: crack mesh sizes used by python -m crackddd2d convergence.
• surface_panels: optional free-surface panel sizes paired with each crack mesh size for edge/free-surface cases.
• metrics: reported diagnostic groups; currently sif, cod, traction_residual, and field_decomposition.

Run Tests

From the repository root:

python -m unittest discover -s py2d_crack_ddd/tests -t py2d_crack_ddd

Current Scope

Implemented:

• infinite-boundary 2D DDD force kernel ported from 2dDDD/forceInf.m;
• x-periodic and xy-periodic DDD force kernels ported from 2dDDD/forceXPBC.m and forceCut.m, with independent x/y image offsets for the xy-periodic image cut-off kernel;
• glide-plane mobility projection;
• fixed line crack with mixed-mode constant elements and sliding/opening COD unknowns and a straight-crack finite-element DDBEM/HBIE traction matrix for the main simulation and benchmark paths;
• no-crack, lagged-crack, and strong-crack coupling modes;
• strong-crack residual diagnostics with position, COD, force, and crack traction components, including relative normalized residuals;
• straight-crack Mode-I COD/SIF benchmark diagnostics for the calibrated ddbemt2d_hbie backend;
• optional straight-crack node-based linear_tip_enriched DDBEM/HBIE path with distributed equivalent-dislocation density, near-source sampling diagnostics, multi-window SIF extrapolation, and JSON/CLI convergence reports;
• optional finite rectangular free-boundary correction for non-periodic cases using physical-boundary traction collocation and offset equivalent dislocation sources, with per-edge residual diagnostics, fixed or scaled regularization, condition-number fallback, parameter scanning, and optional PK force feedback into DDD motion;
• optional strong-coupling retry with reduced timestep when convergence fails;
• JSON/CSV case loading and deterministic text outputs;
• post-processing GIF animation with box, load, boundary condition, dislocation, and COD-colored crack overlays;
• post-processing sigma_xy stress-background GIF from accepted dislocation and crack states.

Not implemented yet:

• fully general high-order HBIE/DDBEM singular and near-singular quadrature for arbitrary curved crack geometry;
• crack growth or remeshing;
• three-dimensional ParaDiS topology/restart/MPI behavior.

Physics Notes

The ddbemt2d_hbie backend is now the primary straight-crack mixed-mode kernel. It assembles a constant-element finite displacement-discontinuity traction matrix from the same 2D plane-strain edge-dislocation Green function used for field reconstruction. Non-self interactions use Gaussian line integration; self collocation uses a principal-value numerical rule with a finite local self term. This is a strict straight-crack constant-element DDBEM/HBIE implementation, but it is still not a high-order or arbitrary curved crack singular quadrature implementation.

The stress-background animation includes the mixed-mode crack-COD stress term from the same finite crack-element Green-function integration used by the solver, with sliding and opening components along the local crack tangent and normal.

The same mixed-mode crack-COD stress path is also used for the current crack feedback force: [sliding, opening] COD -> sigma_crack(x) -> Peach-Koehler force. This keeps the main simulation, stress animation, and static post-processing on a shared crack representation. Remaining calibration work is focused on higher-order singular/near-singular treatment and finite-boundary coupling, not on replacing the crack field with an analytic reference solution.

For non-periodic finite boxes, the optional boundary_correction block reduces the total traction residual on the four outer edges. The residual is evaluated from dislocation + crack + external + boundary_correction; therefore a true four-side traction-free calculation attempts to drive the total edge traction toward zero, not only the traction produced by internal dislocation/crack sources. Set enabled: true and apply_to_forces: true to include this correction in DDD motion; periodic and semi-periodic cells disable it automatically. The preferred finite-boundary mode is now strict_ddbem; the older collocation_proxy remains available as a fast offset-source diagnostic.

loading.mode controls whether loading.applied_stress contributes to this finite-boundary residual. The default uniform_stress mode includes it. Use none when the finite-boundary correction should ignore the external stress field. A nonzero uniform applied stress together with four traction-free outer edges is a modeling choice: the boundary correction will cancel the applied stress traction on the box edges, so the result should not be interpreted as a prescribed-traction boundary-value problem.

Useful boundary_correction controls are samples_per_edge, source_offset, regularization, regularization_mode (fixed or scaled), max_condition_number, mode, quadrature_order, and near_singular_subdivision. The default mode is collocation_proxy, which uses the fast offset equivalent-source correction. Set mode: "strict_ddbem" to use real four-side constant boundary panels and Gaussian line integration for the boundary influence matrix and domain stress reconstruction. If the influence matrix condition number exceeds max_condition_number, or if the correction would increase the traction residual, the solver falls back to zero boundary correction for that evaluation and records the fallback in diagnostics.

The first strict mode is a constant-element 2D DDBEM implementation with principal-value-style self evaluation from symmetric Gauss points. It removes the artificial offset contour used by collocation_proxy and reports tangential/normal residual components and maximum corrected traction, but it is still not a final high-order HBIE implementation with analytic singular and near-singular quadrature for arbitrary crack/boundary geometry.

Strong coupling is monitored with a weighted residual made from DDD position increment, COD increment, total force increment, and crack traction residual. Both absolute and normalized relative residuals are written. If retry_on_failure is enabled in the JSON coupling block, an unconverged strong-coupling step is retried with dt *= dt_reduction up to max_retries. The reported time is the accumulated accepted timestep, so reduced retry steps are reflected in history.csv. If a step still fails and accept_unconverged is false, the run stops with an explicit error instead of silently accepting the step.

The xy-periodic DDD kernel uses a local 3 x 3 image cut-off. This follows the spirit of the MATLAB forceCut.m case but uses independent x and y offsets for the image grid.