- disp_xThe x displacement variable
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:The x displacement variable
- disp_yThe y displacement variable
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:The y displacement variable
- displacementsThe string of displacements suitable for the problem statement
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:The string of displacements suitable for the problem statement
- penalty0The penalty factor
Default:0
C++ Type:Real
Unit:(no unit assumed)
Controllable:No
Description:The penalty factor
- primary_boundaryThe name of the primary boundary sideset.
C++ Type:BoundaryName
Controllable:No
Description:The name of the primary boundary sideset.
- primary_subdomainThe name of the primary subdomain.
C++ Type:SubdomainName
Controllable:No
Description:The name of the primary subdomain.
- secondary_boundaryThe name of the secondary boundary sideset.
C++ Type:BoundaryName
Controllable:No
Description:The name of the secondary boundary sideset.
- secondary_subdomainThe name of the secondary subdomain.
C++ Type:SubdomainName
Controllable:No
Description:The name of the secondary subdomain.
- secondary_variablePrimal variable on secondary surface.
C++ Type:VariableName
Unit:(no unit assumed)
Controllable:No
Description:Primal variable on secondary surface.
BilinearMixedModeCohesiveZoneModel
Description
The BilinearMixedModeCohesiveZoneModel object computes a global traction vector from a local traction vector that is calculated from a bilinear mixed mode cohesive zone model. Theoretical equations and their application to particle fuel is described and discussed in Jiang et al. (2021).
This object implements the bilinear mixed mode traction separation law using the nodal-based mortar quantities that are employed in the MOOSE formulation. Such quantities are identified with a tilde at a secondary surface node . The definitions of mortar surfaces and mortar quantities are identical to those presented in previous technical reports (Martin Recuero and Yushu (2022) and Recuero et al. (2022)), which are in turn based on existing discretization approaches (Wohlmuth (2011) and Puso et al. (2008)).
Quadratic failure criterion
The softening of the interface in the bilinear mixed mode traction separation law is given by a quadratic failure criterion. The onset of softening of the interface takes place when the traction state reaches the following equality: (1) where , , and are the (positive) normal traction, the shear traction along direction one, and the shear traction along direction two, respectively. These three traction components are computed at each node using interpolated, weighted quantities. , , and are the interlaminar tensile and shear strengths, which are model parameters.
Mode mixity
The displacement on the interface can be split into normal and tangential parts. Therefore, the mixed mode relative displacement can be expressed as: (2) where subscripts , and and denote the normal direction and the two tangential directions, respectively.
The formulation assumes that the interface traction generated by the displacements is proportional; in essence, a traction component is equal to a penalty factor times the displacement jump component (, where is the penalty factor). If , the mixity ratio is obtained as: (3)
Assuming , the uni-mode relative displacements at the beginning of softening (denoted by superscript 0'') are given by: (4)
The mixed mode relative displacement corresponding to the initiation of softening, , is defined by the following multi-part definition: (5)
Mode mixity Delamination propagation prediction
The new CZM mortar implementation also accounts for two types of delamination propagation models. (1) A power law-based criterion and (2) The Benzeggagh-Kenane (BK) criterion (see, e.g.,Benzeggagh and Kenane (1996) and Bui (2011)).
Power law-based delamination criterion
The power law criterion may be described as: (6) where and are the critical energy release rates for mixed modes; whereas and are the critical energy release rates for pure modes. The power exponent is can be obtained experimentally.
During the delamination or debonding process, a relative interface displacement may be reached such that cohesive fores completely disappear, i.e. the cohesive model ceases to influence boundary tractions. Such a mixed-mode displacement value is defined as follows (note that we carry the node-wise definition, ): (7)
Benzeggagh-Kenane (BK) delamination criterion
The BK criterion can be written as (Benzeggagh and Kenane (1996)): (8) where , is the total critical strain energy release rate, and is a model parameter.
The mixed-mode total decohesion displacements are, for the BK criterion: (9)
Constitutive model
The CZM constitutive models relate the weighted interface jump components to the interface traction in the local frame of reference. The traction, in the local frame of reference, can be expressed as: (10) where is the component of the local traction and is an effective stiffness. The constitutive components are computed as: (11) where represents the CZM weighted interface jump with shear components, i.e. it applies shear stiffness, and normal stiffness only if a gap exists. applies a normal traction that prevents interpenetration. A few notes about the previous equation:
The action of can potentially be replaced with (mortar) mechanical contact enforcement.
The local traction is computed nodally with weighted quantities, but is interpolated and applied on the primary and secondary surfaces using Lagrange interpolation functions.
The constitutive model can be expanded as needed in the user object.
The damage, , is computed as follows (we omit the nodal subscript):
Auxiliary output conventions
The local-frame displacement jump and cohesive traction components can be output with CohesiveZoneMortarUserObjectAux. The tangential displacement jump effective value is the tangential magnitude, sqrt(t1^2 + t2^2). The cohesive traction tangential magnitude uses the same tangential norm, while the cohesive traction effective value is the full local traction norm, sqrt(tn^2 + t1^2 + t2^2). In two dimensions, the second tangential component, t2, is zero.
The traction components in the local coordinate system follow the convention used internally in the model, which is that opening is positive and tensile traction is negative.
Viscous regularization
Sudden degradation of stiffness may be give rise to the appearance of large displacements. To alleviate this problematic numerical behavior in quasti-static simulations, we include in the constitutive equations a viscous damage instead:
Use standard mortar contact to enforce interpenetratation
The standard cohesive zone model uses a penalty method to prevent penetration in compression. Alternatively, users may choose to use Mortar contact to enforce the normal constraint. To do this, set set_compressive_traction_to_zero = true to remove the compressive traction from the cohesive zone model, and then add the appropriate mortar contact blocks (e.g., NormalMortarMechanicalContact) to handle normal enforcement.
Couple cohesive and frictional behaviors
Under combined compressive and shear loading, the coupling between cohesive and frictional behavior, controlled by the interface damage variable, induces a progressive increase in frictional stress that accompanies the softening of the cohesive response Venzal et al. (2020). Because the BilinearMixedModeCohesiveZoneModel object inherits from PenaltyWeightedGapUserObject, it can compute frictional forces based on the damage variable. These frictional froces can then be used by TangentialMortarMechanicalContact.
Example of usage:
[UserObjects<<<{"href": "../../syntax/UserObjects/index.html"}>>>]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel<<<{"description": "Computes the bilinear mixed mode cohesive zone model.", "href": "BilinearMixedModeCohesiveZoneModel.html"}>>>
primary_boundary<<<{"description": "The name of the primary boundary sideset."}>>> = 'top_bottom'
secondary_boundary<<<{"description": "The name of the secondary boundary sideset."}>>> = 'bottom_top'
primary_subdomain<<<{"description": "The name of the primary subdomain."}>>> = 10000
secondary_subdomain<<<{"description": "The name of the secondary subdomain."}>>> = 10001
disp_x<<<{"description": "The x displacement variable"}>>> = disp_x
disp_y<<<{"description": "The y displacement variable"}>>> = disp_y
friction_coefficient<<<{"description": "The friction coefficient ruling Coulomb friction equations."}>>> = 0.1 # with 2.0 works
secondary_variable<<<{"description": "Primal variable on secondary surface."}>>> = disp_x
penalty<<<{"description": "The penalty factor"}>>> = 0e6
penalty_friction<<<{"description": "The penalty factor for frictional interaction. If not provided, the normal penalty factor is also used for the frictional problem."}>>> = 1e4
use_physical_gap<<<{"description": "Whether to compute the normal pressure as (penalty / A) * (physical gap) instead of penalty * (weighted gap), where A is the lower-dimensional area assigned to the node on the contacting boundary. This parameter is defaulted to 'true' in the contact action."}>>> = true
correct_edge_dropping<<<{"description": "Whether to enable correct edge dropping treatment for mortar constraints. When disabled any Lagrange Multiplier degree of freedom on a secondary element without full primary contributions will be set (strongly) to 0."}>>> = true
normal_strength<<<{"description": "Tensile strength in normal direction."}>>> = N
shear_strength<<<{"description": "Tensile strength in shear direction."}>>> = 1e3
viscosity<<<{"description": "Viscosity for damage model."}>>> = 1e-3
penalty_stiffness<<<{"description": "Penalty stiffness for CZM."}>>> = 1e6
power_law_parameter<<<{"description": "The power law parameter."}>>> = 2.2
GI_c<<<{"description": "Critical energy release rate in normal direction."}>>> = 1e3
GII_c<<<{"description": "Critical energy release rate in shear direction."}>>> = 1e2
displacements<<<{"description": "The string of displacements suitable for the problem statement"}>>> = 'disp_x disp_y'
[]
[]Input Parameters
- adaptivity_penalty_normalSIMPLEThe augmented Lagrange update strategy used on the normal penalty coefficient.
Default:SIMPLE
C++ Type:MooseEnum
Options:SIMPLE, BUSSETTA
Controllable:No
Description:The augmented Lagrange update strategy used on the normal penalty coefficient.
- aux_lmAuxiliary variable that is utilized together with the penalty approach to interpolate the resulting contact tractions using dual bases.
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:Auxiliary variable that is utilized together with the penalty approach to interpolate the resulting contact tractions using dual bases.
- correct_edge_droppingFalseWhether to enable correct edge dropping treatment for mortar constraints. When disabled any Lagrange Multiplier degree of freedom on a secondary element without full primary contributions will be set (strongly) to 0.
Default:False
C++ Type:bool
Controllable:No
Description:Whether to enable correct edge dropping treatment for mortar constraints. When disabled any Lagrange Multiplier degree of freedom on a secondary element without full primary contributions will be set (strongly) to 0.
- debug_meshFalseWhether this constraint is going to enable mortar segment mesh debug information. An exodusfile will be generated if the user sets this flag to true
Default:False
C++ Type:bool
Controllable:No
Description:Whether this constraint is going to enable mortar segment mesh debug information. An exodusfile will be generated if the user sets this flag to true
- disp_zThe z displacement variable
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:The z displacement variable
- friction_coefficient0The friction coefficient ruling Coulomb friction equations.
Default:0
C++ Type:Real
Unit:(no unit assumed)
Controllable:No
Description:The friction coefficient ruling Coulomb friction equations.
- ghost_higher_d_neighborsFalseWhether we should ghost higher-dimensional neighbors. This is necessary when we are doing second order mortar with finite volume primal variables, because in order for the method to be second order we must use cell gradients, which couples in the neighbor cells.
Default:False
C++ Type:bool
Controllable:No
Description:Whether we should ghost higher-dimensional neighbors. This is necessary when we are doing second order mortar with finite volume primal variables, because in order for the method to be second order we must use cell gradients, which couples in the neighbor cells.
- interpolate_normalsFalseWhether to interpolate the nodal normals (e.g. classic idea of evaluating field at quadrature points). If this is set to false, then non-interpolated nodal normals will be used, and then the _normals member should be indexed with _i instead of _qp
Default:False
C++ Type:bool
Controllable:No
Description:Whether to interpolate the nodal normals (e.g. classic idea of evaluating field at quadrature points). If this is set to false, then non-interpolated nodal normals will be used, and then the _normals member should be indexed with _i instead of _qp
- minimum_projection_angle40Parameter to control which angle (in degrees) is admissible for the creation of mortar segments. If set to a value close to zero, very oblique projections are allowed, which can result in mortar segments solving physics not meaningfully, and overprojection of primary nodes onto the mortar segment mesh in extreme cases. In 3D with mortar_3d_subpatch_plane = GEOMETRIC_NORMAL, this parameter also controls which primary and secondary subpatch normal pairings are admissible before polygon clipping.
Default:40
C++ Type:Real
Unit:(no unit assumed)
Controllable:No
Description:Parameter to control which angle (in degrees) is admissible for the creation of mortar segments. If set to a value close to zero, very oblique projections are allowed, which can result in mortar segments solving physics not meaningfully, and overprojection of primary nodes onto the mortar segment mesh in extreme cases. In 3D with mortar_3d_subpatch_plane = GEOMETRIC_NORMAL, this parameter also controls which primary and secondary subpatch normal pairings are admissible before polygon clipping.
- mortar_3d_qp_mappingNORMAL_PROJECTIONMethod used to map quadrature points on 3D mortar segments to primary and secondary parent-face reference coordinates. This parameter has no effect for 2D mortar.
Default:NORMAL_PROJECTION
C++ Type:MooseEnum
Options:NORMAL_PROJECTION, REFERENCE_INTERPOLATION
Controllable:No
Description:Method used to map quadrature points on 3D mortar segments to primary and secondary parent-face reference coordinates. This parameter has no effect for 2D mortar.
- mortar_3d_subpatch_planeGEOMETRIC_NORMALMethod used to construct the local 3D mortar subpatch planes used for projection and clipping. GEOMETRIC_NORMAL uses the geometric normal of each linear or bilinear secondary subpatch. AVERAGED_NODAL_NORMAL uses the averaged nodal normal on each secondary subpatch.
Default:GEOMETRIC_NORMAL
C++ Type:MooseEnum
Options:GEOMETRIC_NORMAL, AVERAGED_NODAL_NORMAL
Controllable:No
Description:Method used to construct the local 3D mortar subpatch planes used for projection and clipping. GEOMETRIC_NORMAL uses the geometric normal of each linear or bilinear secondary subpatch. AVERAGED_NODAL_NORMAL uses the averaged nodal normal on each secondary subpatch.
- penalty_frictionThe penalty factor for frictional interaction. If not provided, the normal penalty factor is also used for the frictional problem.
C++ Type:Real
Unit:(no unit assumed)
Controllable:No
Description:The penalty factor for frictional interaction. If not provided, the normal penalty factor is also used for the frictional problem.
- periodicFalseWhether this constraint is going to be used to enforce a periodic condition. This has the effect of changing the normals vector for projection from outward to inward facing
Default:False
C++ Type:bool
Controllable:No
Description:Whether this constraint is going to be used to enforce a periodic condition. This has the effect of changing the normals vector for projection from outward to inward facing
- primary_variablePrimal variable on primary surface. If this parameter is not provided then the primary variable will be initialized to the secondary variable
C++ Type:VariableName
Unit:(no unit assumed)
Controllable:No
Description:Primal variable on primary surface. If this parameter is not provided then the primary variable will be initialized to the secondary variable
- set_compressive_traction_to_zeroFalseZero compressive traction (set to true, allowing the use of standard zero-penetration mortar contact constraints in the normal direction).
Default:False
C++ Type:bool
Controllable:No
Description:Zero compressive traction (set to true, allowing the use of standard zero-penetration mortar contact constraints in the normal direction).
- triangulate_trianglesFalseWhether a clipped 3D mortar polygon that is already a triangle should still be subdivided during triangulation. When enabled, already-triangular polygons are subdivided with the centroid-based path because the vertex-fan, ear-clipping, and Delaunay backends cannot refine a triangle any further on their own.
Default:False
C++ Type:bool
Controllable:No
Description:Whether a clipped 3D mortar polygon that is already a triangle should still be subdivided during triangulation. When enabled, already-triangular polygons are subdivided with the centroid-based path because the vertex-fan, ear-clipping, and Delaunay backends cannot refine a triangle any further on their own.
- triangulationcentroidStrategy used to triangulate clipped 3D mortar polygons into mortar segments. The default is 'centroid' to preserve the legacy 3D mortar segmentation behavior.
Default:centroid
C++ Type:MooseEnum
Options:vertex, centroid, ear_clipping, delaunay
Controllable:No
Description:Strategy used to triangulate clipped 3D mortar polygons into mortar segments. The default is 'centroid' to preserve the legacy 3D mortar segmentation behavior.
- use_physical_gapFalseWhether to compute the normal pressure as (penalty / A) * (physical gap) instead of penalty * (weighted gap), where A is the lower-dimensional area assigned to the node on the contacting boundary. This parameter is defaulted to 'true' in the contact action.
Default:False
C++ Type:bool
Controllable:No
Description:Whether to compute the normal pressure as (penalty / A) * (physical gap) instead of penalty * (weighted gap), where A is the lower-dimensional area assigned to the node on the contacting boundary. This parameter is defaulted to 'true' in the contact action.
Optional Parameters
- GII_cCritical energy release rate in shear direction.
C++ Type:MaterialPropertyName
Unit:(no unit assumed)
Controllable:No
Description:Critical energy release rate in shear direction.
- GI_cCritical energy release rate in normal direction.
C++ Type:MaterialPropertyName
Unit:(no unit assumed)
Controllable:No
Description:Critical energy release rate in normal direction.
- mixed_mode_criterionBKOption for mixed mode propagation criterion.
Default:BK
C++ Type:MooseEnum
Options:POWER_LAW, BK
Controllable:No
Description:Option for mixed mode propagation criterion.
- normal_strengthTensile strength in normal direction.
C++ Type:MaterialPropertyName
Unit:(no unit assumed)
Controllable:No
Description:Tensile strength in normal direction.
- penalty_stiffnessPenalty stiffness for CZM.
C++ Type:Real
Unit:(no unit assumed)
Range:penalty_stiffness > 0.0
Controllable:No
Description:Penalty stiffness for CZM.
- power_law_parameterThe power law parameter.
C++ Type:Real
Unit:(no unit assumed)
Controllable:No
Description:The power law parameter.
- regularization_alpha1e-10Regularization parameter for the Macaulay bracket.
Default:1e-10
C++ Type:Real
Unit:(no unit assumed)
Controllable:No
Description:Regularization parameter for the Macaulay bracket.
- shear_strengthTensile strength in shear direction.
C++ Type:MaterialPropertyName
Unit:(no unit assumed)
Controllable:No
Description:Tensile strength in shear direction.
- viscosity0Viscosity for damage model.
Default:0
C++ Type:Real
Unit:(no unit assumed)
Controllable:No
Description:Viscosity for damage model.
Bilinear Mixed Mode Traction Parameters
- allow_duplicate_execution_on_initialFalseIn the case where this UserObject is depended upon by an initial condition, allow it to be executed twice during the initial setup (once before the IC and again after mesh adaptivity (if applicable).
Default:False
C++ Type:bool
Controllable:No
Description:In the case where this UserObject is depended upon by an initial condition, allow it to be executed twice during the initial setup (once before the IC and again after mesh adaptivity (if applicable).
- execution_order_group0Execution order groups are executed in increasing order (e.g., the lowest number is executed first). Note that negative group numbers may be used to execute groups before the default (0) group. Please refer to the user object documentation for ordering of user object execution within a group.
Default:0
C++ Type:int
Controllable:No
Description:Execution order groups are executed in increasing order (e.g., the lowest number is executed first). Note that negative group numbers may be used to execute groups before the default (0) group. Please refer to the user object documentation for ordering of user object execution within a group.
- force_postauxFalseForces the UserObject to be executed in POSTAUX
Default:False
C++ Type:bool
Controllable:No
Description:Forces the UserObject to be executed in POSTAUX
- force_preauxFalseForces the UserObject to be executed in PREAUX
Default:False
C++ Type:bool
Controllable:No
Description:Forces the UserObject to be executed in PREAUX
- force_preicFalseForces the UserObject to be executed in PREIC during initial setup
Default:False
C++ Type:bool
Controllable:No
Description:Forces the UserObject to be executed in PREIC during initial setup
Execution Scheduling Parameters
- control_tagsAdds user-defined labels for accessing object parameters via control logic.
C++ Type:std::vector<std::string>
Controllable:No
Description:Adds user-defined labels for accessing object parameters via control logic.
- enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Controllable:Yes
Description:Set the enabled status of the MooseObject.
- use_displaced_meshTrueWhether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.
Default:True
C++ Type:bool
Controllable:No
Description:Whether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.
Advanced Parameters
- prop_getter_suffixAn optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
C++ Type:MaterialPropertyName
Unit:(no unit assumed)
Controllable:No
Description:An optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
- use_interpolated_stateFalseFor the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.
Default:False
C++ Type:bool
Controllable:No
Description:For the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.
Material Property Retrieval Parameters
Input Files
- (modules/combined/test/tests/czm_pressure_sidesets/czm_pressure_sideset_cavity.i)
- (modules/contact/test/tests/cohesive_zone_model/traction_curve_modeI_3d_inclined.i)
- (modules/contact/test/tests/cohesive_zone_model/traction_curve_modeI_3d.i)
- (modules/contact/examples/mortar_czm/progressive_debonding.i)
- (modules/contact/test/tests/cohesive_zone_model/compression_shear_toggle.i)
- (modules/combined/test/tests/czm_pressure_sidesets/czm_pressure_sideset.i)
- (modules/contact/test/tests/cohesive_zone_model/traction_curve_modeI_viscous.i)
- (modules/contact/test/tests/cohesive_zone_model/traction_curve_modeI.i)
- (modules/contact/test/tests/cohesive_zone_model/traction_curve_modeIII_3d.i)
- (modules/contact/test/tests/cohesive_zone_model/mortar_czm.i)
- (modules/contact/test/tests/cohesive_zone_model/mortar_czm_analysis.i)
- (modules/contact/test/tests/cohesive_zone_model/robustness_pure_shear_limit.i)
- (modules/contact/test/tests/cohesive_zone_model/bilinear_mixed_compare.i)
- (modules/contact/test/tests/cohesive_zone_model/traction_curve_mixed.i)
- (modules/contact/test/tests/cohesive_zone_model/traction_cycle_modeI.i)
- (modules/contact/test/tests/cohesive_zone_model/robustness_brittle_interface.i)
- (modules/contact/test/tests/cohesive_zone_model/traction_curve_modeII.i)
- (modules/contact/examples/mortar_czm/debond_recontact.i)
- (modules/contact/test/tests/cohesive_zone_model/cohesive_friction.i)
- (modules/contact/test/tests/cohesive_zone_model/bilinear_mixed_mortar_only_czm.i)
References
- ML Benzeggagh and MJCS Kenane.
Measurement of mixed-mode delamination fracture toughness of unidirectional glass/epoxy composites with mixed-mode bending apparatus.
Composites science and technology, 56(4):439–449, 1996.[Export]
BibTeX
@article{benzeggagh1996measurement, author = "Benzeggagh, ML and Kenane, MJCS", title = "Measurement of mixed-mode delamination fracture toughness of unidirectional glass/epoxy composites with mixed-mode bending apparatus", journal = "Composites science and technology", volume = "56", number = "4", pages = "439--449", year = "1996" }RIS
TY - JOUR AU - Benzeggagh, ML AU - Kenane, MJCS TI - Measurement of mixed-mode delamination fracture toughness of unidirectional glass/epoxy composites with mixed-mode bending apparatus JO - Composites science and technology PY - 1996 VL - 56 IS - 4 SP - 439 EP - 449 ER -Plain Text
ML Benzeggagh and MJCS Kenane. Measurement of mixed-mode delamination fracture toughness of unidirectional glass/epoxy composites with mixed-mode bending apparatus. Composites science and technology, 56(4):439–449, 1996. - Q.V. Bui.
A modified benzeggagh-kenane fracture criterion for mixed-mode delamination.
Journal of Composite Materials, 45(4):389–413, 2011.[Export]
BibTeX
@article{bui2011modified, author = "Bui, Q.V.", title = "A modified Benzeggagh-Kenane fracture criterion for mixed-mode delamination", journal = "Journal of Composite Materials", volume = "45", number = "4", pages = "389--413", year = "2011", publisher = "SAGE Publications Sage UK: London, England" }RIS
TY - JOUR AU - Bui, Q.V. TI - A modified Benzeggagh-Kenane fracture criterion for mixed-mode delamination JO - Journal of Composite Materials PY - 2011 VL - 45 IS - 4 PB - SAGE Publications Sage UK: London, England SP - 389 EP - 413 ER -Plain Text
Q.V. Bui. A modified benzeggagh-kenane fracture criterion for mixed-mode delamination. Journal of Composite Materials, 45(4):389–413, 2011. - Wen Jiang, Gyanender Singh, Jason D Hales, Aysenur Toptan, Benjamin W Spencer, and Stephen R Novascone.
Efficient failure probability calculations and modeling interface debonding in TRISO particles with BISON.
Technical Report, Idaho National Lab.(INL), Idaho Falls, ID (United States), 2021.[Export]
BibTeX
@techreport{jiang2021efficient, author = "Jiang, Wen and Singh, Gyanender and Hales, Jason D and Toptan, Aysenur and Spencer, Benjamin W and Novascone, Stephen R", title = "Efficient Failure Probability Calculations and Modeling Interface Debonding in {TRISO} Particles with {BISON}", year = "2021", institution = "Idaho National Lab.(INL), Idaho Falls, ID (United States)" }RIS
TY - RPRT AU - Jiang, Wen AU - Singh, Gyanender AU - Hales, Jason D AU - Toptan, Aysenur AU - Spencer, Benjamin W AU - Novascone, Stephen R TI - Efficient Failure Probability Calculations and Modeling Interface Debonding in TRISO Particles with BISON PY - 2021 ER -Plain Text
Wen Jiang, Gyanender Singh, Jason D Hales, Aysenur Toptan, Benjamin W Spencer, and Stephen R Novascone. Efficient failure probability calculations and modeling interface debonding in TRISO particles with BISON. Technical Report, Idaho National Lab.(INL), Idaho Falls, ID (United States), 2021. - Antonio Martin Recuero and Dewen Yushu.
Implement and test 3D mortar contact in bison.
Technical Report, Idaho National Laboratory (INL), Idaho Falls, ID (United States), 2022.[Export]
BibTeX
@techreport{martin2022implement, author = "Martin Recuero, Antonio and Yushu, Dewen", title = "Implement and test {3D} mortar contact in BISON", year = "2022", institution = "Idaho National Laboratory (INL), Idaho Falls, ID (United States)" }RIS
TY - RPRT AU - Martin Recuero, Antonio AU - Yushu, Dewen TI - Implement and test 3D mortar contact in BISON PY - 2022 ER -Plain Text
Antonio Martin Recuero and Dewen Yushu. Implement and test 3D mortar contact in bison. Technical Report, Idaho National Laboratory (INL), Idaho Falls, ID (United States), 2022. - Michael A Puso, TA Laursen, and Jerome Solberg.
A segment-to-segment mortar contact method for quadratic elements and large deformations.
Computer Methods in Applied Mechanics and Engineering, 197(6-8):555–566, 2008.[Export]
BibTeX
@article{puso2008segment, author = "Puso, Michael A and Laursen, TA and Solberg, Jerome", title = "A segment-to-segment mortar contact method for quadratic elements and large deformations", journal = "Computer Methods in Applied Mechanics and Engineering", volume = "197", number = "6-8", pages = "555--566", year = "2008", publisher = "Elsevier" }RIS
TY - JOUR AU - Puso, Michael A AU - Laursen, TA AU - Solberg, Jerome TI - A segment-to-segment mortar contact method for quadratic elements and large deformations JO - Computer Methods in Applied Mechanics and Engineering PY - 2008 VL - 197 IS - 6-8 PB - Elsevier SP - 555 EP - 566 ER -Plain Text
Michael A Puso, TA Laursen, and Jerome Solberg. A segment-to-segment mortar contact method for quadratic elements and large deformations. Computer Methods in Applied Mechanics and Engineering, 197(6-8):555–566, 2008. - Antonio Recuero, Alexander Lindsay, Dewen Yushu, John W Peterson, and Benjamin Spencer.
A mortar thermomechanical contact computational framework for nuclear fuel performance simulation.
Nuclear Engineering and Design, 394:111808, 2022.[Export]
BibTeX
@article{recuero2022mortar, author = "Recuero, Antonio and Lindsay, Alexander and Yushu, Dewen and Peterson, John W and Spencer, Benjamin", title = "A mortar thermomechanical contact computational framework for nuclear fuel performance simulation", journal = "Nuclear Engineering and Design", volume = "394", pages = "111808", year = "2022", publisher = "Elsevier" }RIS
TY - JOUR AU - Recuero, Antonio AU - Lindsay, Alexander AU - Yushu, Dewen AU - Peterson, John W AU - Spencer, Benjamin TI - A mortar thermomechanical contact computational framework for nuclear fuel performance simulation JO - Nuclear Engineering and Design PY - 2022 VL - 394 PB - Elsevier SP - 111808 ER -Plain Text
Antonio Recuero, Alexander Lindsay, Dewen Yushu, John W Peterson, and Benjamin Spencer. A mortar thermomechanical contact computational framework for nuclear fuel performance simulation. Nuclear Engineering and Design, 394:111808, 2022. - V. Venzal, S. Morel, T. Parent, and F. Dubois.
Frictional cohesive zone model for quasi-brittle fracture: mixed-mode and coupling between cohesive and frictional behaviors.
International Journal of Solids and Structures, 198:17–30, 2020.
URL: https://www.sciencedirect.com/science/article/pii/S0020768320301402, doi:https://doi.org/10.1016/j.ijsolstr.2020.04.023.[Export]
BibTeX
@article{VENZAL202017, author = "Venzal, V. and Morel, S. and Parent, T. and Dubois, F.", title = "Frictional cohesive zone model for quasi-brittle fracture: Mixed-mode and coupling between cohesive and frictional behaviors", journal = "International Journal of Solids and Structures", volume = "198", pages = "17-30", year = "2020", issn = "0020-7683", doi = "https://doi.org/10.1016/j.ijsolstr.2020.04.023", url = "https://www.sciencedirect.com/science/article/pii/S0020768320301402" }RIS
TY - JOUR AU - Venzal, V. AU - Morel, S. AU - Parent, T. AU - Dubois, F. TI - Frictional cohesive zone model for quasi-brittle fracture: Mixed-mode and coupling between cohesive and frictional behaviors JO - International Journal of Solids and Structures PY - 2020 VL - 198 DO - https://doi.org/10.1016/j.ijsolstr.2020.04.023 UR - https://www.sciencedirect.com/science/article/pii/S0020768320301402 SP - 17 EP - 30 ER -Plain Text
V. Venzal, S. Morel, T. Parent, and F. Dubois. Frictional cohesive zone model for quasi-brittle fracture: mixed-mode and coupling between cohesive and frictional behaviors. International Journal of Solids and Structures, 198:17-30, 2020. URL: https://www.sciencedirect.com/science/article/pii/S0020768320301402, doi:https://doi.org/10.1016/j.ijsolstr.2020.04.023. - Barbara Wohlmuth.
Variationally consistent discretization schemes and numerical algorithms for contact problems.
Acta Numerica, 20:569–734, 2011.[Export]
BibTeX
@article{wohlmuth2011variationally, author = "Wohlmuth, Barbara", title = "Variationally consistent discretization schemes and numerical algorithms for contact problems", journal = "Acta Numerica", volume = "20", pages = "569--734", year = "2011", publisher = "Cambridge University Press" }RIS
TY - JOUR AU - Wohlmuth, Barbara TI - Variationally consistent discretization schemes and numerical algorithms for contact problems JO - Acta Numerica PY - 2011 VL - 20 PB - Cambridge University Press SP - 569 EP - 734 ER -Plain Text
Barbara Wohlmuth. Variationally consistent discretization schemes and numerical algorithms for contact problems. Acta Numerica, 20:569–734, 2011.
(modules/contact/test/tests/cohesive_zone_model/mortar_czm.i)
[Mesh]
[msh]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymax = 1
nx = 5
ny = 5
boundary_name_prefix = bottom
[]
[msh_id]
type = SubdomainIDGenerator
input = msh
subdomain_id = 1
[]
[msh_two]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = 1
ymax = 2
nx = 5
ny = 5
boundary_name_prefix = top
boundary_id_offset = 10
[]
[msh_two_id]
type = SubdomainIDGenerator
input = msh_two
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'msh_id msh_two_id'
[]
[top_node]
type = ExtraNodesetGenerator
coord = '0 2 0'
input = combined
new_boundary = top_node
[]
[bottom_node]
type = ExtraNodesetGenerator
coord = '0 0 0'
input = top_node
new_boundary = bottom_node
[]
# Build subdomains
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'bottom_top'
input = bottom_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = 'top_bottom'
new_block_name = 'primary_lower'
input = secondary
[]
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
generate_output = 'stress_yy'
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = TaylorExpansion
block = '1 2'
[]
[]
[]
[]
[BCs]
[fix_x]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_node
variable = disp_x
[]
[fix_top]
type = DirichletBC
preset = true
boundary = top_top
variable = disp_x
value = 0
[]
[top]
type = FunctionDirichletBC
boundary = top_top
variable = disp_y
function = 'if(t<=0.3,t,if(t<=0.6,0.3-(t-0.3),0.6-t))'
preset = true
[]
[bottom]
type = DirichletBC
boundary = bottom_bottom
variable = disp_y
value = 0
preset = true
[]
[]
[AuxVariables]
[]
[AuxKernels]
[]
[Materials]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric9
C_ijkl = '1.684e5 0.176e5 0.176e5 1.684e5 0.176e5 1.684e5 0.754e5 0.754e5 0.754e5'
block = '1 2'
[]
[normal_strength]
type = GenericFunctionMaterial
prop_names = 'N'
prop_values = 'if(x<0.5,1,100)*1e4'
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'NEWTON'
line_search = none
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
automatic_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 30
nl_rel_tol = 1e-11
nl_abs_tol = 1e-11
start_time = 0.0
dt = 0.01
end_time = 0.05
dtmin = 0.01
[]
[Outputs]
exodus = true
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
disp_x = disp_x
disp_y = disp_y
friction_coefficient = 0.1 # with 2.0 works
secondary_variable = disp_x
penalty = 0e6
penalty_friction = 1e4
use_physical_gap = true
correct_edge_dropping = true
normal_strength = N
shear_strength = 1e3
viscosity = 1e-3
penalty_stiffness = 1e6
power_law_parameter = 2.2
GI_c = 1e3
GII_c = 1e2
displacements = 'disp_x disp_y'
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
correct_edge_dropping = true
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
correct_edge_dropping = true
[]
[c_x]
type = MortarGenericTraction
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
correct_edge_dropping = true
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
correct_edge_dropping = true
[]
[]
(modules/combined/test/tests/czm_pressure_sidesets/czm_pressure_sideset_cavity.i)
# Pressurization of the sidesets built by MortarCZMDebondedBoundaryUpdater. The
# czm_pressure_sideset.i test is repeated here while the upper block also
# slides tangentially, so the faces paired in the undisplaced configuration are no
# longer aligned. Once the whole interface has failed the blocks carry no load, so the void between
# the two sidesets must equal the integral of the prescribed opening, opening * t * (1 - x / 2),
# over x in [0, 1], which is 3 / 4 * opening * t.
#
# The void is a cavity, connected to an unmeshed plenum of volume Vp, whose gas pressure follows
# the ideal gas law, p (V + Vp) = n R T, and acts on both sidesets, which is only possible if the
# pressure boundary condition picks up the faces added during the run. Gas is injected at a constant
# rate from full separation at t = 0.8 until t = 1, so the pressure first rises, and then falls as
# the peel keeps enlarging the void. Each block is in uniaxial compression, sigma_yy = -p, so that
# - the vertical reactions on lower_bottom and upper_top are +p and -p times the interface length;
# - the lower block displaces by u_y = -p / E * (y + 0.5), and the upper block by
# p / E * (0.5 - y) in addition to the imposed peel;
# - the void area is V = V0 + p / E, with V0 = 3 / 4 * opening * t, as each block of height 0.5
# shortens by p / (2 E).
# These fields are linear, so the elements reproduce them exactly. Substituting V into the gas law
# gives p^2 / E + W p - n R T = 0 with W = V0 + Vp, whose positive root is
# p = E / 2 * (sqrt(W^2 + 4 n R T / E) - W).
opening = 0.1
slide = 0.02
youngs_modulus = 1e9
gas_constant = 8.314
gas_temperature = 300
# Plenum volume, of the order of the void area, which keeps the gas volume positive while the void
# is still closed
plenum_volume = 0.05
# Injection rate chosen so that the pressure peaks near 1e7 at t = 1; p / E = 0.01 is large enough
# to see in the displacements
gas_rate = 2700
[Mesh]
# The interface lies on y = 0, so the void integrand -y n_y vanishes on the vertical closing
# segment at the crack front, which belongs to neither sideset.
[lower]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = -0.5
ymax = 0
nx = 4
ny = 1
boundary_name_prefix = lower
[]
[lower_id]
type = SubdomainIDGenerator
input = lower
subdomain_id = 1
[]
[upper]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = 0
ymax = 0.5
nx = 3
ny = 1
boundary_name_prefix = upper
boundary_id_offset = 10
[]
[upper_id]
type = SubdomainIDGenerator
input = upper
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'lower_id upper_id'
[]
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = secondary_lower
sidesets = lower_top
input = combined
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
new_block_name = primary_lower
sidesets = upper_bottom
input = secondary
[]
add_sideset_names = 'czm_pressure_primary czm_pressure_secondary'
add_sideset_ids = '1000 1001'
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Problem]
extra_tag_vectors = 'ref'
[]
[Physics/SolidMechanics/QuasiStatic/all]
strain = SMALL
add_variables = true
use_automatic_differentiation = true
extra_vector_tags = 'ref'
block = '1 2'
[]
[AuxVariables]
[damage]
[]
[react_y]
block = '1 2'
[]
[]
[AuxKernels]
[damage]
type = CohesiveZoneMortarUserObjectAux
variable = damage
user_object = czm_uo
cohesive_zone_quantity = cohesive_damage
[]
[react_y]
type = ReactionForceAux
vector_tag = 'ref'
v = disp_y
variable = react_y
block = '1 2'
[]
[]
[Functions]
[gas_moles]
type = ParsedFunction
expression = '${gas_rate} * max(min(t, 1) - 0.8, 0)'
[]
[expected_pressure]
type = ParsedFunction
expression = 'W := 0.75 * ${opening} * t + ${plenum_volume};
n := ${gas_rate} * max(min(t, 1) - 0.8, 0);
${youngs_modulus} / 2 * (sqrt(W^2 + 4 * n * ${gas_constant} * ${gas_temperature} / ${youngs_modulus}) - W)'
[]
[]
[BCs]
[fix_x]
type = DirichletBC
boundary = lower_bottom
variable = disp_x
value = 0
[]
[fix_y]
type = DirichletBC
boundary = lower_bottom
variable = disp_y
value = 0
[]
[slide_x]
type = FunctionDirichletBC
boundary = upper_top
variable = disp_x
function = '${slide} * t'
[]
[peel_y]
type = FunctionDirichletBC
boundary = upper_top
variable = disp_y
function = '${opening} * t * (1 - x / 2)'
[]
[CavityPressure]
[void]
boundary = 'czm_pressure_primary czm_pressure_secondary'
R = ${gas_constant}
temperature = gas_temperature
volume = void_area
material_input = gas_moles
additional_volumes = plenum_volume
temperature_of_additional_volumes = gas_temperature
output = cavity_pressure
use_automatic_differentiation = true
# Act on the undisplaced faces so that the load stays purely normal to the interface and the
# solution above is exact
use_displaced_mesh = false
[]
[]
[]
[Materials]
[elasticity_tensor]
type = ADComputeIsotropicElasticityTensor
youngs_modulus = ${youngs_modulus}
poissons_ratio = 0
block = '1 2'
[]
[stress]
type = ADComputeLinearElasticStress
block = '1 2'
[]
[normal_strength]
type = GenericConstantMaterial
prop_names = normal_strength
prop_values = 1e4
[]
[shear_strength]
type = GenericConstantMaterial
prop_names = shear_strength
prop_values = 1e4
[]
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = upper_bottom
secondary_boundary = lower_top
primary_subdomain = primary_lower
secondary_subdomain = secondary_lower
disp_x = disp_x
disp_y = disp_y
secondary_variable = disp_x
use_physical_gap = true
correct_edge_dropping = true
normal_strength = normal_strength
shear_strength = shear_strength
penalty_stiffness = 1e12
power_law_parameter = 2
# Fully separated (damage = 1) once the opening exceeds 2 GI_c / normal_strength = 4e-2, which
# is much larger than the 1e-5 elastic stretch of the blocks at peak traction, so the response
# does not snap back
GI_c = 200
GII_c = 200
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = upper_bottom
secondary_boundary = lower_top
primary_subdomain = primary_lower
secondary_subdomain = secondary_lower
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = upper_bottom
secondary_boundary = lower_top
primary_subdomain = primary_lower
secondary_subdomain = secondary_lower
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[c_x]
type = MortarGenericTraction
primary_boundary = upper_bottom
secondary_boundary = lower_top
primary_subdomain = primary_lower
secondary_subdomain = secondary_lower
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_y]
type = MortarGenericTraction
primary_boundary = upper_bottom
secondary_boundary = lower_top
primary_subdomain = primary_lower
secondary_subdomain = secondary_lower
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[]
[MeshModifiers]
[pressure_sideset_updater]
type = MortarCZMDebondedBoundaryUpdater
block = secondary_lower
damage_variable = damage
cohesive_zone_userobject = czm_uo
primary_boundary = upper_bottom
secondary_boundary = lower_top
primary_subdomain = primary_lower
secondary_subdomain = secondary_lower
debonded_boundary_primary = czm_pressure_primary
debonded_boundary_secondary = czm_pressure_secondary
execute_on = TIMESTEP_END
# Run before the default group so that the postprocessors see faces promoted this step
execution_order_group = -1
[]
[]
[Postprocessors]
# Undisplaced length and first moment of each sideset
[primary_length]
type = AreaPostprocessor
boundary = czm_pressure_primary
[]
[secondary_length]
type = AreaPostprocessor
boundary = czm_pressure_secondary
[]
[primary_moment]
type = FunctionSideIntegral
boundary = czm_pressure_primary
function = x
[]
[secondary_moment]
type = FunctionSideIntegral
boundary = czm_pressure_secondary
function = x
[]
# Area between the two sidesets in the current configuration, which is the cavity volume
[void_area]
type = InternalVolume
boundary = 'czm_pressure_primary czm_pressure_secondary'
component = 1
execute_on = 'initial linear timestep_end'
[]
[expected_void_area_separated]
type = ParsedPostprocessor
expression = '0.75 * ${opening} * t + expected_pressure / ${youngs_modulus}'
pp_names = 'expected_pressure'
use_t = true
[]
# Cavity inputs
[gas_moles]
type = FunctionValuePostprocessor
function = gas_moles
execute_on = 'initial linear timestep_end'
[]
[gas_temperature]
type = ConstantPostprocessor
value = ${gas_temperature}
execute_on = 'initial linear timestep_end'
[]
[plenum_volume]
type = ConstantPostprocessor
value = ${plenum_volume}
execute_on = 'initial linear timestep_end'
[]
# Pressure checks, valid after full separation
[lower_reaction_y]
type = NodalSum
variable = react_y
boundary = lower_bottom
[]
[upper_reaction_y]
type = NodalSum
variable = react_y
boundary = upper_top
[]
[expected_pressure]
type = FunctionValuePostprocessor
function = expected_pressure
[]
# The lower block moves only because of the pressure
[lower_disp_y]
type = PointValue
variable = disp_y
point = '0.5 -0.25 0'
[]
[expected_lower_disp_y]
type = ParsedPostprocessor
expression = '-expected_pressure / ${youngs_modulus} * 0.25'
pp_names = 'expected_pressure'
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type -pc_factor_mat_solver_type'
petsc_options_value = 'lu superlu_dist'
nl_rel_tol = 1e-10
nl_abs_tol = 1e-10
dt = 0.1
end_time = 1.2
[]
[Outputs]
csv = true
[]
(modules/contact/test/tests/cohesive_zone_model/traction_curve_modeI_3d_inclined.i)
# Purpose: Verify mode-I traction-separation with a 45-degree inclined 3D interface.
# Assertion: A rotated 3D interface produces local-frame normal, tangential, effective traction, and damage outputs matching the CSV gold curve.
# Why this exists: Confirms local-frame outputs remain correct on interfaces not aligned with the global axes.
[Mesh]
[base]
type = GeneratedMeshGenerator
dim = 3
xmax = 1.1
xmin = -0.1
ymax = 1
zmax = 1
nx = 1
ny = 1
nz = 1
[]
[rename_base]
type = RenameBoundaryGenerator
input = base
old_boundary = 'top bottom left right front back'
new_boundary = 'top_base bottom_base left_base right_base front_base back_base'
[]
[base_id]
type = SubdomainIDGenerator
input = rename_base
subdomain_id = 1
[]
[top]
type = GeneratedMeshGenerator
dim = 3
xmax = 1
ymin = 1
ymax = 2
zmax = 1
nx = 1
ny = 1
nz = 1
[]
[rename_top]
type = RenameBoundaryGenerator
input = top
old_boundary = 'top bottom left right front back'
new_boundary = '100 101 102 103 104 105'
[]
[top_id]
type = SubdomainIDGenerator
input = rename_top
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'base_id top_id'
[]
# Rotate the entire mesh 45 degrees about the z-axis
[rotate]
type = TransformGenerator
input = combined
transform = ROTATE
vector_value = '0 0 45'
[]
# Rotated coordinates: (0,2,0) -> (-sqrt(2), sqrt(2), 0)
[top_node]
type = ExtraNodesetGenerator
coord = '-1.4142135623731 1.4142135623731 0'
input = rotate
new_boundary = top_node
[]
# Rotated coordinates: (-0.1,0,0) -> (-0.1*cos45, -0.1*sin45, 0)
[bottom_node]
type = ExtraNodesetGenerator
coord = '-0.070710678118655 -0.070710678118655 0'
input = top_node
new_boundary = bottom_node
[]
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'top_base'
input = bottom_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = '101'
new_block_name = 'primary_lower'
input = secondary
[]
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = TaylorExpansion
block = '1 2'
[]
[]
[]
[]
[AuxVariables]
[damage]
[]
[local_normal_jump]
[]
[cohesive_traction_normal]
[]
[cohesive_traction_tangential_1]
[]
[cohesive_traction_tangential_2]
[]
[cohesive_traction_effective]
[]
[]
[AuxKernels]
[cohesive_damage]
type = CohesiveZoneMortarUserObjectAux
variable = damage
user_object = czm_uo
cohesive_zone_quantity = cohesive_damage
[]
[local_normal_jump]
type = CohesiveZoneMortarUserObjectAux
variable = local_normal_jump
user_object = czm_uo
cohesive_zone_quantity = local_normal_jump
[]
[cohesive_traction_normal]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_normal
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_normal
[]
[cohesive_traction_tangential_1]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_tangential_1
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_tangential_1
[]
[cohesive_traction_tangential_2]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_tangential_2
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_tangential_2
[]
[cohesive_traction_effective]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_effective
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_effective
[]
[]
[BCs]
# Fix bottom boundary in all directions
[fix_bottom_x]
type = DirichletBC
boundary = bottom_base
variable = disp_x
value = 0
preset = true
[]
[fix_bottom_y]
type = DirichletBC
boundary = bottom_base
variable = disp_y
value = 0
preset = true
[]
[anchor_bottom_z]
type = DirichletBC
boundary = bottom_node
variable = disp_z
value = 0
preset = true
[]
# Top boundary: apply displacement along rotated normal [-1/sqrt(2), 1/sqrt(2), 0]
# Original: disp_y = 0.3*t along [0,1,0]
# Rotated: disp_x = -0.3*t/sqrt(2), disp_y = 0.3*t/sqrt(2)
[top_x]
type = FunctionDirichletBC
boundary = 100
variable = disp_x
function = '-0.21213203435596*t'
preset = true
[]
[top_y]
type = FunctionDirichletBC
boundary = 100
variable = disp_y
function = '0.21213203435596*t'
preset = true
[]
[anchor_top_z]
type = DirichletBC
preset = true
boundary = top_node
variable = disp_z
value = 0
[]
[]
[Materials]
[normal_strength]
type = GenericConstantMaterial
prop_names = 'normal_strength'
prop_values = '50'
[]
[shear_strength]
type = GenericConstantMaterial
prop_names = 'shear_strength'
prop_values = '50'
[]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric_isotropic_E_nu
C_ijkl = '1e8 0.3'
block = '1 2'
[]
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
disp_x = disp_x
disp_y = disp_y
disp_z = disp_z
secondary_variable = disp_x
use_physical_gap = true
correct_edge_dropping = true
set_compressive_traction_to_zero = true
friction_coefficient = 0.0
penalty = 0.0
penalty_friction = 0.0
normal_strength = 'normal_strength'
shear_strength = 'shear_strength'
penalty_stiffness = 1000
power_law_parameter = 2.0
GI_c = 5
GII_c = 5
displacements = 'disp_x disp_y disp_z'
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[z]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_z
component = z
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[c_x]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_z]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_z
component = z
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[]
[Postprocessors]
[jump_n]
type = ElementAverageValue
variable = local_normal_jump
block = secondary_lower
[]
[traction_n]
type = ElementAverageValue
variable = cohesive_traction_normal
block = secondary_lower
[]
[traction_t1]
type = ElementAverageValue
variable = cohesive_traction_tangential_1
block = secondary_lower
[]
[traction_t2]
type = ElementAverageValue
variable = cohesive_traction_tangential_2
block = secondary_lower
[]
[traction_eff]
type = ElementAverageValue
variable = cohesive_traction_effective
block = secondary_lower
[]
[cohesive_damage]
type = ElementAverageValue
variable = damage
block = secondary_lower
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
line_search = none
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
automatic_scaling = true
compute_scaling_once = false
off_diagonals_in_auto_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 100
nl_rel_tol = 1e-12
nl_abs_tol = 1e-11
start_time = 0.0
dt = 0.01
end_time = 1.0
dtmin = 0.01
[]
[Outputs]
exodus = false
csv = true
[]
(modules/contact/test/tests/cohesive_zone_model/traction_curve_modeI_3d.i)
# Purpose: Verify mode-I traction-separation behavior in a 3D single-element geometry.
# Assertion: Normal traction follows the CSV gold curve while both tangential traction components remain negligible.
# Why this exists: Protects the 3D local-frame split of normal versus tangential traction components.
[Mesh]
[base]
type = GeneratedMeshGenerator
dim = 3
xmax = 1.1
xmin = -0.1
ymax = 1
zmax = 1
nx = 1
ny = 1
nz = 1
[]
[rename_base]
type = RenameBoundaryGenerator
input = base
old_boundary = 'top bottom left right front back'
new_boundary = 'top_base bottom_base left_base right_base front_base back_base'
[]
[base_id]
type = SubdomainIDGenerator
input = rename_base
subdomain_id = 1
[]
[top]
type = GeneratedMeshGenerator
dim = 3
xmax = 1
ymin = 1
ymax = 2
zmax = 1
nx = 1
ny = 1
nz = 1
[]
[rename_top]
type = RenameBoundaryGenerator
input = top
old_boundary = 'top bottom left right front back'
new_boundary = '100 101 102 103 104 105'
[]
[top_id]
type = SubdomainIDGenerator
input = rename_top
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'base_id top_id'
[]
[top_node]
type = ExtraNodesetGenerator
coord = '0 2 0'
input = combined
new_boundary = top_node
[]
[top_roller_node]
type = ExtraNodesetGenerator
coord = '1 2 0'
input = top_node
new_boundary = top_roller_node
[]
[bottom_node]
type = ExtraNodesetGenerator
coord = '-0.1 0 0'
input = top_roller_node
new_boundary = bottom_node
[]
[bottom_roller_node]
type = ExtraNodesetGenerator
coord = '1.1 0 0'
input = bottom_node
new_boundary = bottom_roller_node
[]
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'top_base'
input = bottom_roller_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = '101'
new_block_name = 'primary_lower'
input = secondary
[]
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = TaylorExpansion
block = '1 2'
[]
[]
[]
[]
[AuxVariables]
[damage]
[]
[local_normal_jump]
[]
[cohesive_traction_normal]
[]
[cohesive_traction_tangential_1]
[]
[cohesive_traction_tangential_2]
[]
[cohesive_traction_effective]
[]
[]
[AuxKernels]
[cohesive_damage]
type = CohesiveZoneMortarUserObjectAux
variable = damage
user_object = czm_uo
cohesive_zone_quantity = cohesive_damage
[]
[local_normal_jump]
type = CohesiveZoneMortarUserObjectAux
variable = local_normal_jump
user_object = czm_uo
cohesive_zone_quantity = local_normal_jump
[]
[cohesive_traction_normal]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_normal
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_normal
[]
[cohesive_traction_tangential_1]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_tangential_1
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_tangential_1
[]
[cohesive_traction_tangential_2]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_tangential_2
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_tangential_2
[]
[cohesive_traction_effective]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_effective
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_effective
[]
[]
[BCs]
[anchor_bottom_x]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_node
variable = disp_x
[]
[anchor_bottom_z]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_node
variable = disp_z
[]
[roller_bottom_z]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_roller_node
variable = disp_z
[]
[roller_bottom_y]
type = DirichletBC
boundary = bottom_base
variable = disp_y
value = 0
preset = true
[]
[anchor_top_x]
type = DirichletBC
preset = true
boundary = top_node
variable = disp_x
value = 0
[]
[anchor_top_z]
type = DirichletBC
preset = true
boundary = top_node
variable = disp_z
value = 0
[]
[roller_top_z]
type = DirichletBC
preset = true
boundary = top_roller_node
variable = disp_z
value = 0
[]
[top]
type = FunctionDirichletBC
boundary = 100
variable = disp_y
function = '0.3*t'
preset = true
[]
[]
[Materials]
[normal_strength]
type = GenericConstantMaterial
prop_names = 'normal_strength'
prop_values = '50'
[]
[shear_strength]
type = GenericConstantMaterial
prop_names = 'shear_strength'
prop_values = '50'
[]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric_isotropic_E_nu
C_ijkl = '1e8 0.3'
block = '1 2'
[]
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
disp_x = disp_x
disp_y = disp_y
disp_z = disp_z
secondary_variable = disp_x
use_physical_gap = true
correct_edge_dropping = true
set_compressive_traction_to_zero = true
friction_coefficient = 0.0
penalty = 0.0
penalty_friction = 0.0
normal_strength = 'normal_strength'
shear_strength = 'shear_strength'
penalty_stiffness = 1000
power_law_parameter = 2.0
GI_c = 5
GII_c = 5
displacements = 'disp_x disp_y disp_z'
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[z]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_z
component = z
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[c_x]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_z]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_z
component = z
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[]
[Postprocessors]
[jump_n]
type = ElementAverageValue
variable = local_normal_jump
block = secondary_lower
[]
[traction_n]
type = ElementAverageValue
variable = cohesive_traction_normal
block = secondary_lower
[]
[traction_t1]
type = ElementAverageValue
variable = cohesive_traction_tangential_1
block = secondary_lower
[]
[traction_t2]
type = ElementAverageValue
variable = cohesive_traction_tangential_2
block = secondary_lower
[]
[traction_eff]
type = ElementAverageValue
variable = cohesive_traction_effective
block = secondary_lower
[]
[cohesive_damage]
type = ElementAverageValue
variable = damage
block = secondary_lower
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
line_search = none
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
automatic_scaling = true
compute_scaling_once = false
off_diagonals_in_auto_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 100
nl_rel_tol = 1e-12
nl_abs_tol = 1e-11
start_time = 0.0
dt = 0.01
end_time = 1.0
dtmin = 0.01
[]
[Outputs]
exodus = false
csv = true
[]
(modules/contact/examples/mortar_czm/progressive_debonding.i)
# Purpose: Verify heterogeneous interface debonding with spatially varying normal strength.
# Setup: Two 5x5 blocks, mode-I pull to t=0.3. Left half (x<0.5) has N=1e4 (weak),
# right half (x>=0.5) has N=1e6 (strong). K=1e6, GI_c=1e3, GII_c=1e2, viscosity=1e-3.
# Exercises: Asymmetric damage initiation, load redistribution (weak->strong), near-brittle
# snap-through with viscous regularization, and full debonding of both sides.
# Outputs: Exodus with interface fields on secondary_lower, CSV with PointValue at x=0.25/0.75.
[Mesh]
[msh]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymax = 1
nx = 5
ny = 5
boundary_name_prefix = bottom
[]
[msh_id]
type = SubdomainIDGenerator
input = msh
subdomain_id = 1
[]
[msh_two]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = 1
ymax = 2
nx = 5
ny = 5
boundary_name_prefix = top
boundary_id_offset = 10
[]
[msh_two_id]
type = SubdomainIDGenerator
input = msh_two
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'msh_id msh_two_id'
[]
[top_node]
type = ExtraNodesetGenerator
coord = '0 2 0'
input = combined
new_boundary = top_node
[]
[bottom_node]
type = ExtraNodesetGenerator
coord = '0 0 0'
input = top_node
new_boundary = bottom_node
[]
# Build subdomains
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'bottom_top'
input = bottom_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = 'top_bottom'
new_block_name = 'primary_lower'
input = secondary
[]
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
generate_output = 'stress_yy'
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = TaylorExpansion
block = '1 2'
[]
[]
[]
[]
[BCs]
[anchor_bottom_x]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_node
variable = disp_x
[]
[anchor_top_x]
type = DirichletBC
preset = true
boundary = top_node
variable = disp_x
value = 0
[]
[top]
type = FunctionDirichletBC
boundary = top_top
variable = disp_y
function = 'if(t<=0.3,t,if(t<=0.6,0.3-(t-0.3),0.6-t))'
preset = true
[]
[roller_bottom_y]
type = DirichletBC
boundary = bottom_bottom
variable = disp_y
value = 0
preset = true
[]
[]
[AuxVariables]
[damage]
block = secondary_lower
[]
[local_normal_jump]
block = secondary_lower
[]
[cohesive_traction_normal]
block = secondary_lower
[]
[cohesive_traction_tangential_magnitude]
block = secondary_lower
[]
[mode_mixity_ratio]
block = secondary_lower
[]
[]
[AuxKernels]
[cohesive_damage]
type = CohesiveZoneMortarUserObjectAux
variable = damage
user_object = czm_uo
cohesive_zone_quantity = cohesive_damage
[]
[local_normal_jump]
type = CohesiveZoneMortarUserObjectAux
variable = local_normal_jump
user_object = czm_uo
cohesive_zone_quantity = local_normal_jump
[]
[cohesive_traction_normal]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_normal
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_normal
[]
[cohesive_traction_tangential_magnitude]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_tangential_magnitude
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_tangential_magnitude
[]
[mode_mixity_ratio]
type = CohesiveZoneMortarUserObjectAux
variable = mode_mixity_ratio
user_object = czm_uo
cohesive_zone_quantity = mode_mixity_ratio
[]
[]
[Materials]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric9
C_ijkl = '1.684e5 0.176e5 0.176e5 1.684e5 0.176e5 1.684e5 0.754e5 0.754e5 0.754e5'
block = '1 2'
[]
[normal_strength]
type = GenericFunctionMaterial
prop_names = 'N'
prop_values = 'if(x<0.5,1,100)*1e4'
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'NEWTON'
line_search = none
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
automatic_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 30
nl_rel_tol = 1e-11
nl_abs_tol = 1e-11
start_time = 0.0
dt = 0.01
end_time = 0.3
dtmin = 1e-4
[]
[Postprocessors]
# Left side (weak: N = 1e4) - sample at x=0.25 on interface
[damage_left]
type = PointValue
variable = damage
point = '0.25 1.0 0'
[]
[jump_n_left]
type = PointValue
variable = local_normal_jump
point = '0.25 1.0 0'
[]
[traction_n_left]
type = PointValue
variable = cohesive_traction_normal
point = '0.25 1.0 0'
[]
# Right side (strong: N = 1e6) - sample at x=0.75 on interface
[damage_right]
type = PointValue
variable = damage
point = '0.75 1.0 0'
[]
[jump_n_right]
type = PointValue
variable = local_normal_jump
point = '0.75 1.0 0'
[]
[traction_n_right]
type = PointValue
variable = cohesive_traction_normal
point = '0.75 1.0 0'
[]
[]
[Outputs]
exodus = true
csv = true
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
disp_x = disp_x
disp_y = disp_y
friction_coefficient = 0.1 # with 2.0 works
secondary_variable = disp_x
penalty = 0e6
penalty_friction = 1e4
use_physical_gap = true
correct_edge_dropping = true
normal_strength = N
shear_strength = 1e3
viscosity = 1e-3
penalty_stiffness = 1e6
power_law_parameter = 2.2
GI_c = 1e3
GII_c = 1e2
displacements = 'disp_x disp_y'
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
correct_edge_dropping = true
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
correct_edge_dropping = true
[]
[c_x]
type = MortarGenericTraction
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
correct_edge_dropping = true
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
correct_edge_dropping = true
[]
[]
(modules/contact/test/tests/cohesive_zone_model/compression_shear_toggle.i)
# Purpose: Verify the compressive-traction toggle behavior under combined compression and shear loading.
# Assertion: Normal cohesive traction is suppressed when enabled and retained when disabled, matching CSV gold baselines.
# Why this exists: Confirms the set_compressive_traction_to_zero switch produces the expected alternate constitutive behavior.
[Mesh]
[base]
type = GeneratedMeshGenerator
dim = 2
xmax = 1.1
ymax = 1
xmin = -0.1
nx = 1
ny = 1
[]
[rename_base]
type = RenameBoundaryGenerator
input = base
old_boundary = 'top bottom left right'
new_boundary = 'top_base bottom_base left_base right_base'
[]
[base_id]
type = SubdomainIDGenerator
input = rename_base
subdomain_id = 1
[]
[top]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = 1
ymax = 2
nx = 1
ny = 1
[]
[rename_top]
type = RenameBoundaryGenerator
input = top
old_boundary = 'top bottom left right'
new_boundary = '100 101 102 103'
[]
[top_id]
type = SubdomainIDGenerator
input = rename_top
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'base_id top_id'
[]
[top_node]
type = ExtraNodesetGenerator
coord = '0 2 0'
input = combined
new_boundary = top_node
[]
[bottom_node]
type = ExtraNodesetGenerator
coord = '-0.1 0 0'
input = top_node
new_boundary = bottom_node
[]
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'top_base'
input = bottom_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = '101'
new_block_name = 'primary_lower'
input = secondary
[]
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = TaylorExpansion
block = '1 2'
[]
[]
[]
[]
[AuxVariables]
[local_normal_jump]
[]
[local_tangential_jump_1]
[]
[cohesive_traction_normal]
[]
[cohesive_traction_tangential_magnitude]
[]
[]
[AuxKernels]
[local_normal_jump]
type = CohesiveZoneMortarUserObjectAux
variable = local_normal_jump
user_object = czm_uo
cohesive_zone_quantity = local_normal_jump
[]
[local_tangential_jump_1]
type = CohesiveZoneMortarUserObjectAux
variable = local_tangential_jump_1
user_object = czm_uo
cohesive_zone_quantity = local_tangential_jump_1
[]
[cohesive_traction_normal]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_normal
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_normal
[]
[cohesive_traction_tangential_magnitude]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_tangential_magnitude
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_tangential_magnitude
[]
[]
[BCs]
[fix_bottom_x]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_base
variable = disp_x
[]
[top_x]
type = FunctionDirichletBC
boundary = 100
variable = disp_x
function = '0.1*t'
preset = true
[]
[top_y]
type = FunctionDirichletBC
boundary = 100
variable = disp_y
function = '-0.1*t'
preset = true
[]
[roller_bottom_y]
type = DirichletBC
boundary = bottom_base
variable = disp_y
value = 0
preset = true
[]
[]
[Materials]
[normal_strength]
type = GenericConstantMaterial
prop_names = 'normal_strength'
prop_values = '50'
[]
[shear_strength]
type = GenericConstantMaterial
prop_names = 'shear_strength'
prop_values = '50'
[]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric_isotropic_E_nu
C_ijkl = '1e8 0.3'
block = '1 2'
[]
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
disp_x = disp_x
disp_y = disp_y
secondary_variable = disp_x
use_physical_gap = true
correct_edge_dropping = true
set_compressive_traction_to_zero = true
friction_coefficient = 0.0
penalty = 0.0
penalty_friction = 0.0
normal_strength = 'normal_strength'
shear_strength = 'shear_strength'
penalty_stiffness = 1000
power_law_parameter = 2.0
GI_c = 5
GII_c = 5
displacements = 'disp_x disp_y'
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[c_x]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[]
[Postprocessors]
[jump_n]
type = ElementAverageValue
variable = local_normal_jump
block = secondary_lower
[]
[jump_t]
type = ElementAverageValue
variable = local_tangential_jump_1
block = secondary_lower
[]
[traction_n]
type = ElementAverageValue
variable = cohesive_traction_normal
block = secondary_lower
[]
[traction_t]
type = ElementAverageValue
variable = cohesive_traction_tangential_magnitude
block = secondary_lower
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
line_search = none
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
automatic_scaling = true
compute_scaling_once = false
off_diagonals_in_auto_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 100
nl_rel_tol = 1e-12
nl_abs_tol = 1e-11
start_time = 0.0
dt = 0.02
end_time = 0.3
dtmin = 0.02
[]
[Outputs]
exodus = false
csv = true
[]
(modules/combined/test/tests/czm_pressure_sidesets/czm_pressure_sideset.i)
# Two blocks bonded by a mortar cohesive zone are peeled apart from x = 0, so the debonded region
# grows from left to right. Due to the model configuration, each secondary face overlaps one or two
# primary faces.
#
# As each secondary face fully separates, it must be added to czm_pressure_secondary, and every
# primary face it fully overlaps in the undisplaced configuration must be added to
# czm_pressure_primary. The length and first moment (integral of x) of each sideset record which
# faces it holds: a sideset holding the first k faces of width h from x = 0 has length k h and
# first moment (k h)^2 / 2. Once the whole interface has failed the blocks carry no load, so the
# void between the two sidesets must equal the integral of the prescribed opening,
# opening * t * (1 - x / 2), over x in [0, 1], which is 3 / 4 * opening * t.
opening = 0.1
[Mesh]
# The interface lies on y = 0 and the upper block is not displaced in x, so the void integrand
# -y n_y vanishes on the vertical closing segments at both ends of the void, which belong to
# neither sideset.
[lower]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = -0.5
ymax = 0
nx = 4
ny = 1
boundary_name_prefix = lower
[]
[lower_id]
type = SubdomainIDGenerator
input = lower
subdomain_id = 1
[]
[upper]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = 0
ymax = 0.5
nx = 3
ny = 1
boundary_name_prefix = upper
boundary_id_offset = 10
[]
[upper_id]
type = SubdomainIDGenerator
input = upper
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'lower_id upper_id'
[]
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = secondary_lower
sidesets = lower_top
input = combined
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
new_block_name = primary_lower
sidesets = upper_bottom
input = secondary
[]
add_sideset_names = 'czm_pressure_primary czm_pressure_secondary'
add_sideset_ids = '1000 1001'
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Physics/SolidMechanics/QuasiStatic/all]
strain = SMALL
add_variables = true
use_automatic_differentiation = true
block = '1 2'
[]
[AuxVariables]
[damage]
[]
[]
[AuxKernels]
[damage]
type = CohesiveZoneMortarUserObjectAux
variable = damage
user_object = czm_uo
cohesive_zone_quantity = cohesive_damage
[]
[]
[BCs]
[fix_x]
type = DirichletBC
boundary = lower_bottom
variable = disp_x
value = 0
[]
[fix_y]
type = DirichletBC
boundary = lower_bottom
variable = disp_y
value = 0
[]
[fix_upper_x]
type = DirichletBC
boundary = upper_top
variable = disp_x
value = 0
[]
[peel_y]
type = FunctionDirichletBC
boundary = upper_top
variable = disp_y
function = '${opening} * t * (1 - x / 2)'
[]
[]
[Materials]
[elasticity_tensor]
type = ADComputeIsotropicElasticityTensor
youngs_modulus = 1e9
poissons_ratio = 0
block = '1 2'
[]
[stress]
type = ADComputeLinearElasticStress
block = '1 2'
[]
[normal_strength]
type = GenericConstantMaterial
prop_names = normal_strength
prop_values = 1e4
[]
[shear_strength]
type = GenericConstantMaterial
prop_names = shear_strength
prop_values = 1e4
[]
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = upper_bottom
secondary_boundary = lower_top
primary_subdomain = primary_lower
secondary_subdomain = secondary_lower
disp_x = disp_x
disp_y = disp_y
secondary_variable = disp_x
use_physical_gap = true
correct_edge_dropping = true
normal_strength = normal_strength
shear_strength = shear_strength
penalty_stiffness = 1e12
power_law_parameter = 2
# Fully separated (damage = 1) once the opening exceeds 2 GI_c / normal_strength = 4e-2, which
# is much larger than the 1e-5 elastic stretch of the blocks at peak traction, so the response
# does not snap back
GI_c = 200
GII_c = 200
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = upper_bottom
secondary_boundary = lower_top
primary_subdomain = primary_lower
secondary_subdomain = secondary_lower
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = upper_bottom
secondary_boundary = lower_top
primary_subdomain = primary_lower
secondary_subdomain = secondary_lower
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[c_x]
type = MortarGenericTraction
primary_boundary = upper_bottom
secondary_boundary = lower_top
primary_subdomain = primary_lower
secondary_subdomain = secondary_lower
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_y]
type = MortarGenericTraction
primary_boundary = upper_bottom
secondary_boundary = lower_top
primary_subdomain = primary_lower
secondary_subdomain = secondary_lower
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[]
[MeshModifiers]
[pressure_sideset_updater]
type = MortarCZMDebondedBoundaryUpdater
block = secondary_lower
damage_variable = damage
cohesive_zone_userobject = czm_uo
primary_boundary = upper_bottom
secondary_boundary = lower_top
primary_subdomain = primary_lower
secondary_subdomain = secondary_lower
debonded_boundary_primary = czm_pressure_primary
debonded_boundary_secondary = czm_pressure_secondary
execute_on = TIMESTEP_END
# Run before the default group so that the postprocessors see faces promoted this step
execution_order_group = -1
[]
[]
[Postprocessors]
# Undisplaced length and first moment of each sideset
[primary_length]
type = AreaPostprocessor
boundary = czm_pressure_primary
[]
[secondary_length]
type = AreaPostprocessor
boundary = czm_pressure_secondary
[]
[primary_moment]
type = FunctionSideIntegral
boundary = czm_pressure_primary
function = x
[]
[secondary_moment]
type = FunctionSideIntegral
boundary = czm_pressure_secondary
function = x
[]
# Area between the two sidesets in the current configuration
[void_area]
type = InternalVolume
boundary = 'czm_pressure_primary czm_pressure_secondary'
component = 1
[]
[expected_void_area_separated]
type = ParsedPostprocessor
expression = '0.75 * ${opening} * t'
use_t = true
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type -pc_factor_mat_solver_type'
petsc_options_value = 'lu superlu_dist'
nl_rel_tol = 1e-10
nl_abs_tol = 1e-10
# Small enough that each secondary face separates in its own time step (t = 0.45, 0.55, 0.6,
# 0.7), so that every intermediate state of both sidesets is checked
dt = 0.05
end_time = 0.8
[]
[Outputs]
csv = true
[]
(modules/contact/test/tests/cohesive_zone_model/traction_curve_modeI_viscous.i)
# Purpose: Verify viscously regularized mode-I traction-separation output.
# Assertion: The viscously regularized traction-separation response matches the CSV gold curve.
# Why this exists: Protects the viscous regularization option of the bilinear cohesive law.
[Mesh]
[base]
type = GeneratedMeshGenerator
dim = 2
xmax = 1.1
ymax = 1
xmin = -0.1
nx = 1
ny = 1
[]
[rename_base]
type = RenameBoundaryGenerator
input = base
old_boundary = 'top bottom left right'
new_boundary = 'top_base bottom_base left_base right_base'
[]
[base_id]
type = SubdomainIDGenerator
input = rename_base
subdomain_id = 1
[]
[top]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = 1
ymax = 2
nx = 1
ny = 1
[]
[rename_top]
type = RenameBoundaryGenerator
input = top
old_boundary = 'top bottom left right'
new_boundary = '100 101 102 103'
[]
[top_id]
type = SubdomainIDGenerator
input = rename_top
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'base_id top_id'
[]
[top_node]
type = ExtraNodesetGenerator
coord = '0 2 0'
input = combined
new_boundary = top_node
[]
[bottom_node]
type = ExtraNodesetGenerator
coord = '-0.1 0 0'
input = top_node
new_boundary = bottom_node
[]
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'top_base'
input = bottom_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = '101'
new_block_name = 'primary_lower'
input = secondary
[]
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = TaylorExpansion
block = '1 2'
[]
[]
[]
[]
[AuxVariables]
[damage]
[]
[local_normal_jump]
[]
[cohesive_traction_normal]
[]
[cohesive_traction_effective]
[]
[]
[AuxKernels]
[cohesive_damage]
type = CohesiveZoneMortarUserObjectAux
variable = damage
user_object = czm_uo
cohesive_zone_quantity = cohesive_damage
[]
[local_normal_jump]
type = CohesiveZoneMortarUserObjectAux
variable = local_normal_jump
user_object = czm_uo
cohesive_zone_quantity = local_normal_jump
[]
[cohesive_traction_normal]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_normal
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_normal
[]
[cohesive_traction_effective]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_effective
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_effective
[]
[]
[BCs]
[anchor_bottom_x]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_node
variable = disp_x
[]
[anchor_top_x]
type = DirichletBC
preset = true
boundary = top_node
variable = disp_x
value = 0
[]
[top]
type = FunctionDirichletBC
boundary = 100
variable = disp_y
function = '0.3*t'
preset = true
[]
[roller_bottom_y]
type = DirichletBC
boundary = bottom_base
variable = disp_y
value = 0
preset = true
[]
[]
[Materials]
[normal_strength]
type = GenericConstantMaterial
prop_names = 'normal_strength'
prop_values = '50'
[]
[shear_strength]
type = GenericConstantMaterial
prop_names = 'shear_strength'
prop_values = '50'
[]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric_isotropic_E_nu
C_ijkl = '1e8 0.3'
block = '1 2'
[]
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
disp_x = disp_x
disp_y = disp_y
secondary_variable = disp_x
use_physical_gap = true
correct_edge_dropping = true
set_compressive_traction_to_zero = true
friction_coefficient = 0.0
penalty = 0.0
penalty_friction = 0.0
normal_strength = 'normal_strength'
shear_strength = 'shear_strength'
penalty_stiffness = 1000
power_law_parameter = 2.0
GI_c = 5
GII_c = 5
viscosity = 1e-2
displacements = 'disp_x disp_y'
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[c_x]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[]
[Postprocessors]
[jump_n]
type = ElementAverageValue
variable = local_normal_jump
block = secondary_lower
[]
[traction_n]
type = ElementAverageValue
variable = cohesive_traction_normal
block = secondary_lower
[]
[traction_eff]
type = ElementAverageValue
variable = cohesive_traction_effective
block = secondary_lower
[]
[cohesive_damage]
type = ElementAverageValue
variable = damage
block = secondary_lower
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
line_search = none
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
automatic_scaling = true
compute_scaling_once = false
off_diagonals_in_auto_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 100
nl_rel_tol = 1e-12
nl_abs_tol = 1e-11
start_time = 0.0
dt = 0.01
end_time = 1.0
dtmin = 0.01
[]
[Outputs]
exodus = false
csv = true
[]
(modules/contact/test/tests/cohesive_zone_model/traction_curve_modeI.i)
# Purpose: Verify monotonic mode-I traction-separation behavior and exposed normal jump/traction outputs.
# Assertion: Mode-I cohesive traction and jump response match CSV gold curves across default and CLI-varied timestep, mesh, mismatch, and parallel runs.
# Why this exists: Protects the bilinear normal traction law output under temporal and spatial changes.
[Mesh]
[base]
type = GeneratedMeshGenerator
dim = 2
xmax = 1.1
ymax = 1
xmin = -0.1
nx = 1
ny = 1
[]
[rename_base]
type = RenameBoundaryGenerator
input = base
old_boundary = 'top bottom left right'
new_boundary = 'top_base bottom_base left_base right_base'
[]
[base_id]
type = SubdomainIDGenerator
input = rename_base
subdomain_id = 1
[]
[top]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = 1
ymax = 2
nx = 1
ny = 1
[]
[rename_top]
type = RenameBoundaryGenerator
input = top
old_boundary = 'top bottom left right'
new_boundary = '100 101 102 103'
[]
[top_id]
type = SubdomainIDGenerator
input = rename_top
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'base_id top_id'
[]
[top_node]
type = ExtraNodesetGenerator
coord = '0 2 0'
input = combined
new_boundary = top_node
[]
[bottom_node]
type = ExtraNodesetGenerator
coord = '-0.1 0 0'
input = top_node
new_boundary = bottom_node
[]
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'top_base'
input = bottom_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = '101'
new_block_name = 'primary_lower'
input = secondary
[]
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = TaylorExpansion
block = '1 2'
[]
[]
[]
[]
[AuxVariables]
[mode_mixity_ratio]
[]
[damage]
[]
[local_normal_jump]
[]
[local_tangential_jump_1]
[]
[cohesive_traction_normal]
[]
[cohesive_traction_tangential_magnitude]
[]
[cohesive_traction_effective]
[]
[]
[AuxKernels]
[mode_mixity_ratio]
type = CohesiveZoneMortarUserObjectAux
variable = mode_mixity_ratio
user_object = czm_uo
cohesive_zone_quantity = mode_mixity_ratio
[]
[cohesive_damage]
type = CohesiveZoneMortarUserObjectAux
variable = damage
user_object = czm_uo
cohesive_zone_quantity = cohesive_damage
[]
[local_normal_jump]
type = CohesiveZoneMortarUserObjectAux
variable = local_normal_jump
user_object = czm_uo
cohesive_zone_quantity = local_normal_jump
[]
[local_tangential_jump_1]
type = CohesiveZoneMortarUserObjectAux
variable = local_tangential_jump_1
user_object = czm_uo
cohesive_zone_quantity = local_tangential_jump_1
[]
[cohesive_traction_normal]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_normal
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_normal
[]
[cohesive_traction_tangential_magnitude]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_tangential_magnitude
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_tangential_magnitude
[]
[cohesive_traction_effective]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_effective
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_effective
[]
[]
[BCs]
[anchor_bottom_x]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_node
variable = disp_x
[]
[anchor_top_x]
type = DirichletBC
preset = true
boundary = top_node
variable = disp_x
value = 0
[]
[top]
type = FunctionDirichletBC
boundary = 100
variable = disp_y
function = '0.3*t'
preset = true
[]
[roller_bottom_y]
type = DirichletBC
boundary = bottom_base
variable = disp_y
value = 0
preset = true
[]
[]
[Materials]
[normal_strength]
type = GenericConstantMaterial
prop_names = 'normal_strength'
prop_values = '50'
[]
[shear_strength]
type = GenericConstantMaterial
prop_names = 'shear_strength'
prop_values = '50'
[]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric_isotropic_E_nu
C_ijkl = '1e8 0.3'
block = '1 2'
[]
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
disp_x = disp_x
disp_y = disp_y
secondary_variable = disp_x
use_physical_gap = true
correct_edge_dropping = true
set_compressive_traction_to_zero = true
friction_coefficient = 0.0
penalty = 0.0
penalty_friction = 0.0
normal_strength = 'normal_strength'
shear_strength = 'shear_strength'
penalty_stiffness = 1000
power_law_parameter = 2.0
GI_c = 5
GII_c = 5
displacements = 'disp_x disp_y'
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[c_x]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[]
[Postprocessors]
[jump_n]
type = ElementAverageValue
variable = local_normal_jump
block = secondary_lower
[]
[jump_t]
type = ElementAverageValue
variable = local_tangential_jump_1
block = secondary_lower
[]
[traction_n]
type = ElementAverageValue
variable = cohesive_traction_normal
block = secondary_lower
[]
[traction_t]
type = ElementAverageValue
variable = cohesive_traction_tangential_magnitude
block = secondary_lower
[]
[traction_eff]
type = ElementAverageValue
variable = cohesive_traction_effective
block = secondary_lower
[]
[cohesive_damage]
type = ElementAverageValue
variable = damage
block = secondary_lower
[]
[energy_rate]
type = ParsedPostprocessor
expression = 'abs(traction_n) * 0.3'
pp_names = 'traction_n'
[]
[energy_dissipated]
type = TimeIntegratedPostprocessor
value = energy_rate
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
line_search = none
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
automatic_scaling = true
compute_scaling_once = false
off_diagonals_in_auto_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 100
nl_rel_tol = 1e-12
nl_abs_tol = 1e-11
start_time = 0.0
dt = 0.01
end_time = 1.0
dtmin = 0.01
[]
[Outputs]
exodus = false
csv = true
[]
(modules/contact/test/tests/cohesive_zone_model/traction_curve_modeIII_3d.i)
# Purpose: Verify pre-debond out-of-plane mode-III CZM outputs in a 3D single-element geometry.
# Assertion: z-shear produces local_tangential_jump_2, local_tangential_jump_effective, and cohesive_traction_tangential_2 matching the CSV gold curve while normal and first tangential traction remain negligible.
# Why this exists: Protects the second local tangential component output before full debond, where post-debond jump output is unstable.
[Mesh]
[base]
type = GeneratedMeshGenerator
dim = 3
xmax = 1.1
xmin = -0.1
ymax = 1
zmax = 1
nx = 1
ny = 1
nz = 1
[]
[rename_base]
type = RenameBoundaryGenerator
input = base
old_boundary = 'top bottom left right front back'
new_boundary = 'top_base bottom_base left_base right_base front_base back_base'
[]
[base_id]
type = SubdomainIDGenerator
input = rename_base
subdomain_id = 1
[]
[top]
type = GeneratedMeshGenerator
dim = 3
xmax = 1
ymin = 1
ymax = 2
zmax = 1
nx = 1
ny = 1
nz = 1
[]
[rename_top]
type = RenameBoundaryGenerator
input = top
old_boundary = 'top bottom left right front back'
new_boundary = '100 101 102 103 104 105'
[]
[top_id]
type = SubdomainIDGenerator
input = rename_top
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'base_id top_id'
[]
[top_node]
type = ExtraNodesetGenerator
coord = '0 2 0'
input = combined
new_boundary = top_node
[]
[bottom_node]
type = ExtraNodesetGenerator
coord = '-0.1 0 0'
input = top_node
new_boundary = bottom_node
[]
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'top_base'
input = bottom_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = '101'
new_block_name = 'primary_lower'
input = secondary
[]
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = TaylorExpansion
block = '1 2'
[]
[]
[]
[]
[AuxVariables]
[damage]
[]
[local_normal_jump]
[]
[cohesive_traction_normal]
[]
[cohesive_traction_tangential_1]
[]
[cohesive_traction_tangential_2]
[]
[cohesive_traction_effective]
[]
[local_tangential_jump_1]
[]
[local_tangential_jump_2]
[]
[local_tangential_jump_effective]
[]
[]
[AuxKernels]
[cohesive_damage]
type = CohesiveZoneMortarUserObjectAux
variable = damage
user_object = czm_uo
cohesive_zone_quantity = cohesive_damage
[]
[local_normal_jump]
type = CohesiveZoneMortarUserObjectAux
variable = local_normal_jump
user_object = czm_uo
cohesive_zone_quantity = local_normal_jump
[]
[cohesive_traction_normal]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_normal
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_normal
[]
[cohesive_traction_tangential_1]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_tangential_1
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_tangential_1
[]
[cohesive_traction_tangential_2]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_tangential_2
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_tangential_2
[]
[cohesive_traction_effective]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_effective
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_effective
[]
[local_tangential_jump_1]
type = CohesiveZoneMortarUserObjectAux
variable = local_tangential_jump_1
user_object = czm_uo
cohesive_zone_quantity = local_tangential_jump_1
[]
[local_tangential_jump_2]
type = CohesiveZoneMortarUserObjectAux
variable = local_tangential_jump_2
user_object = czm_uo
cohesive_zone_quantity = local_tangential_jump_2
[]
[local_tangential_jump_effective]
type = CohesiveZoneMortarUserObjectAux
variable = local_tangential_jump_effective
user_object = czm_uo
cohesive_zone_quantity = local_tangential_jump_effective
[]
[]
[BCs]
[anchor_bottom_x]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_node
variable = disp_x
[]
[roller_bottom_y]
type = DirichletBC
boundary = bottom_base
variable = disp_y
value = 0
preset = true
[]
[fix_bottom_z]
type = DirichletBC
boundary = bottom_base
variable = disp_z
value = 0
preset = true
[]
[anchor_top_x]
type = DirichletBC
preset = true
boundary = top_node
variable = disp_x
value = 0
[]
[roller_top_y]
type = DirichletBC
preset = true
boundary = 100
variable = disp_y
value = 0
[]
[top_z]
type = FunctionDirichletBC
boundary = 100
variable = disp_z
function = '0.3*t'
preset = true
[]
[]
[Materials]
[normal_strength]
type = GenericConstantMaterial
prop_names = 'normal_strength'
prop_values = '50'
[]
[shear_strength]
type = GenericConstantMaterial
prop_names = 'shear_strength'
prop_values = '50'
[]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric_isotropic_E_nu
C_ijkl = '1e8 0.3'
block = '1 2'
[]
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
disp_x = disp_x
disp_y = disp_y
disp_z = disp_z
secondary_variable = disp_x
use_physical_gap = true
correct_edge_dropping = true
set_compressive_traction_to_zero = true
friction_coefficient = 0.0
penalty = 0.0
penalty_friction = 0.0
normal_strength = 'normal_strength'
shear_strength = 'shear_strength'
penalty_stiffness = 1000
power_law_parameter = 2.0
GI_c = 5
GII_c = 5
displacements = 'disp_x disp_y disp_z'
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[z]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_z
component = z
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[c_x]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_z]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_z
component = z
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[]
[Postprocessors]
[jump_n]
type = ElementAverageValue
variable = local_normal_jump
block = secondary_lower
[]
[traction_n]
type = ElementAverageValue
variable = cohesive_traction_normal
block = secondary_lower
[]
[traction_t1]
type = ElementAverageValue
variable = cohesive_traction_tangential_1
block = secondary_lower
[]
[traction_t2]
type = ElementAverageValue
variable = cohesive_traction_tangential_2
block = secondary_lower
[]
[traction_eff]
type = ElementAverageValue
variable = cohesive_traction_effective
block = secondary_lower
[]
[cohesive_damage]
type = ElementAverageValue
variable = damage
block = secondary_lower
[]
[jump_t]
type = ElementAverageValue
variable = local_tangential_jump_2
block = secondary_lower
[]
[jump_t_eff]
type = ElementAverageValue
variable = local_tangential_jump_effective
block = secondary_lower
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
line_search = none
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
automatic_scaling = true
compute_scaling_once = false
off_diagonals_in_auto_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 100
nl_rel_tol = 1e-12
nl_abs_tol = 1e-11
start_time = 0.0
dt = 0.01
# End just before complete debond (damage = 1 at t = 0.67): post-failure the block
# slides freely and the tangential jump output is no longer meaningful.
end_time = 0.66
dtmin = 0.01
[]
[Outputs]
exodus = false
csv = true
[]
(modules/contact/test/tests/cohesive_zone_model/mortar_czm.i)
[Mesh]
[msh]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymax = 1
nx = 5
ny = 5
boundary_name_prefix = bottom
[]
[msh_id]
type = SubdomainIDGenerator
input = msh
subdomain_id = 1
[]
[msh_two]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = 1
ymax = 2
nx = 5
ny = 5
boundary_name_prefix = top
boundary_id_offset = 10
[]
[msh_two_id]
type = SubdomainIDGenerator
input = msh_two
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'msh_id msh_two_id'
[]
[top_node]
type = ExtraNodesetGenerator
coord = '0 2 0'
input = combined
new_boundary = top_node
[]
[bottom_node]
type = ExtraNodesetGenerator
coord = '0 0 0'
input = top_node
new_boundary = bottom_node
[]
# Build subdomains
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'bottom_top'
input = bottom_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = 'top_bottom'
new_block_name = 'primary_lower'
input = secondary
[]
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
generate_output = 'stress_yy'
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = TaylorExpansion
block = '1 2'
[]
[]
[]
[]
[BCs]
[fix_x]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_node
variable = disp_x
[]
[fix_top]
type = DirichletBC
preset = true
boundary = top_top
variable = disp_x
value = 0
[]
[top]
type = FunctionDirichletBC
boundary = top_top
variable = disp_y
function = 'if(t<=0.3,t,if(t<=0.6,0.3-(t-0.3),0.6-t))'
preset = true
[]
[bottom]
type = DirichletBC
boundary = bottom_bottom
variable = disp_y
value = 0
preset = true
[]
[]
[AuxVariables]
[]
[AuxKernels]
[]
[Materials]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric9
C_ijkl = '1.684e5 0.176e5 0.176e5 1.684e5 0.176e5 1.684e5 0.754e5 0.754e5 0.754e5'
block = '1 2'
[]
[normal_strength]
type = GenericFunctionMaterial
prop_names = 'N'
prop_values = 'if(x<0.5,1,100)*1e4'
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'NEWTON'
line_search = none
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
automatic_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 30
nl_rel_tol = 1e-11
nl_abs_tol = 1e-11
start_time = 0.0
dt = 0.01
end_time = 0.05
dtmin = 0.01
[]
[Outputs]
exodus = true
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
disp_x = disp_x
disp_y = disp_y
friction_coefficient = 0.1 # with 2.0 works
secondary_variable = disp_x
penalty = 0e6
penalty_friction = 1e4
use_physical_gap = true
correct_edge_dropping = true
normal_strength = N
shear_strength = 1e3
viscosity = 1e-3
penalty_stiffness = 1e6
power_law_parameter = 2.2
GI_c = 1e3
GII_c = 1e2
displacements = 'disp_x disp_y'
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
correct_edge_dropping = true
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
correct_edge_dropping = true
[]
[c_x]
type = MortarGenericTraction
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
correct_edge_dropping = true
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
correct_edge_dropping = true
[]
[]
(modules/contact/test/tests/cohesive_zone_model/mortar_czm_analysis.i)
[Mesh]
[msh]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymax = 1
nx = 5
ny = 5
boundary_name_prefix = bottom
[]
[msh_id]
type = SubdomainIDGenerator
input = msh
subdomain_id = 1
[]
[msh_two]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = 1
ymax = 2
nx = 5
ny = 5
boundary_name_prefix = top
boundary_id_offset = 10
[]
[msh_two_id]
type = SubdomainIDGenerator
input = msh_two
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'msh_id msh_two_id'
[]
[top_node]
type = ExtraNodesetGenerator
coord = '0 2 0'
input = combined
new_boundary = top_node
[]
[bottom_node]
type = ExtraNodesetGenerator
coord = '0 0 0'
input = top_node
new_boundary = bottom_node
[]
# Build subdomains
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'bottom_top'
input = bottom_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = 'top_bottom'
new_block_name = 'primary_lower'
input = secondary
[]
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
generate_output = 'stress_yy'
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = TaylorExpansion
block = '1 2'
[]
[]
[]
[]
[BCs]
[fix_x]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_node
variable = disp_x
[]
[lateral_top]
type = FunctionDirichletBC
boundary = top_top
variable = disp_x
function = 'if(t<=0.3,t,if(t<=0.6,0.3-(t-0.3),0.6-t))'
preset = true
[]
[top]
type = FunctionDirichletBC
boundary = top_top
variable = disp_y
function = 'if(t<=0.3,t,if(t<=0.6,0.3-(t-0.3),0.6-t))'
preset = true
[]
[bottom]
type = DirichletBC
boundary = bottom_bottom
variable = disp_y
value = 0
preset = true
[]
[]
[AuxVariables]
[mode_mixity_ratio]
[]
[damage]
[]
[local_normal_jump]
[]
[local_tangential_jump]
[]
[]
[AuxKernels]
[mode_mixity_ratio]
type = CohesiveZoneMortarUserObjectAux
variable = mode_mixity_ratio
user_object = czm_uo
cohesive_zone_quantity = mode_mixity_ratio
boundary = 'bottom_top'
[]
[cohesive_damage]
type = CohesiveZoneMortarUserObjectAux
variable = damage
user_object = czm_uo
cohesive_zone_quantity = cohesive_damage
boundary = 'bottom_top'
[]
[local_normal_jump]
type = CohesiveZoneMortarUserObjectAux
variable = local_normal_jump
user_object = czm_uo
cohesive_zone_quantity = local_normal_jump
boundary = 'bottom_top'
[]
[local_tangential_jump]
type = CohesiveZoneMortarUserObjectAux
variable = local_tangential_jump
user_object = czm_uo
cohesive_zone_quantity = local_tangential_jump
boundary = 'bottom_top'
[]
[]
[Materials]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric9
C_ijkl = '1.684e5 0.176e5 0.176e5 1.684e5 0.176e5 1.684e5 0.754e5 0.754e5 0.754e5'
block = '1 2'
[]
[normal_strength]
type = GenericFunctionMaterial
prop_names = 'N'
prop_values = 'if(x<0.5,1,100)*1e4'
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'NEWTON'
line_search = none
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
automatic_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 30
nl_rel_tol = 1e-10
nl_abs_tol = 1e-10
start_time = 0.0
dt = 0.0005
end_time = 0.01
dtmin = 0.0001
[]
[Outputs]
exodus = true
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
disp_x = disp_x
disp_y = disp_y
friction_coefficient = 0.1 # with 2.0 works
secondary_variable = disp_x
penalty = 0e6
penalty_friction = 1e4
use_physical_gap = true
correct_edge_dropping = true
normal_strength = N
shear_strength = 1e3
viscosity = 1e-3
penalty_stiffness = 1e6
mixed_mode_criterion = POWER_LAW
power_law_parameter = 2.2
GI_c = 1e3
GII_c = 1e2
displacements = 'disp_x disp_y'
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
correct_edge_dropping = true
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
correct_edge_dropping = true
[]
[c_x]
type = MortarGenericTraction
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
correct_edge_dropping = true
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 'top_bottom'
secondary_boundary = 'bottom_top'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
correct_edge_dropping = true
[]
[]
(modules/contact/test/tests/cohesive_zone_model/robustness_pure_shear_limit.i)
# Purpose: Verify numerical stability when normal opening is near zero under dominant shear loading.
# Assertion: Dominant-shear output stays finite when the mode-mixity ratio is large, matching the CSV gold baseline.
# Why this exists: Guards the mode-mixity computation against the near-zero normal-opening limit.
[Mesh]
[base]
type = GeneratedMeshGenerator
dim = 2
xmax = 1.1
ymax = 1
xmin = -0.1
nx = 1
ny = 1
[]
[rename_base]
type = RenameBoundaryGenerator
input = base
old_boundary = 'top bottom left right'
new_boundary = 'top_base bottom_base left_base right_base'
[]
[base_id]
type = SubdomainIDGenerator
input = rename_base
subdomain_id = 1
[]
[top]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = 1
ymax = 2
nx = 1
ny = 1
[]
[rename_top]
type = RenameBoundaryGenerator
input = top
old_boundary = 'top bottom left right'
new_boundary = '100 101 102 103'
[]
[top_id]
type = SubdomainIDGenerator
input = rename_top
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'base_id top_id'
[]
[top_node]
type = ExtraNodesetGenerator
coord = '0 2 0'
input = combined
new_boundary = top_node
[]
[bottom_node]
type = ExtraNodesetGenerator
coord = '-0.1 0 0'
input = top_node
new_boundary = bottom_node
[]
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'top_base'
input = bottom_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = '101'
new_block_name = 'primary_lower'
input = secondary
[]
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = TaylorExpansion
block = '1 2'
[]
[]
[]
[]
[AuxVariables]
[damage]
[]
[local_normal_jump]
[]
[local_tangential_jump_1]
[]
[cohesive_traction_normal]
[]
[cohesive_traction_tangential_magnitude]
[]
[cohesive_traction_effective]
[]
[]
[AuxKernels]
[cohesive_damage]
type = CohesiveZoneMortarUserObjectAux
variable = damage
user_object = czm_uo
cohesive_zone_quantity = cohesive_damage
[]
[local_normal_jump]
type = CohesiveZoneMortarUserObjectAux
variable = local_normal_jump
user_object = czm_uo
cohesive_zone_quantity = local_normal_jump
[]
[local_tangential_jump_1]
type = CohesiveZoneMortarUserObjectAux
variable = local_tangential_jump_1
user_object = czm_uo
cohesive_zone_quantity = local_tangential_jump_1
[]
[cohesive_traction_normal]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_normal
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_normal
[]
[cohesive_traction_tangential_magnitude]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_tangential_magnitude
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_tangential_magnitude
[]
[cohesive_traction_effective]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_effective
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_effective
[]
[]
[BCs]
[fix_bottom_x]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_base
variable = disp_x
[]
# Near-zero normal opening to exercise the epsilon_tolerance path in computeModeMixity.
[top_y]
type = FunctionDirichletBC
boundary = 100
variable = disp_y
function = '1e-6*t'
preset = true
[]
# Dominant shear loading drives damage via tangential jump.
[top_x]
type = FunctionDirichletBC
boundary = 100
variable = disp_x
function = '0.3*t'
preset = true
[]
[roller_bottom_y]
type = DirichletBC
boundary = bottom_base
variable = disp_y
value = 0
preset = true
[]
[]
[Materials]
[normal_strength]
type = GenericConstantMaterial
prop_names = 'normal_strength'
prop_values = '50'
[]
[shear_strength]
type = GenericConstantMaterial
prop_names = 'shear_strength'
prop_values = '50'
[]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric_isotropic_E_nu
C_ijkl = '1e8 0.3'
block = '1 2'
[]
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
disp_x = disp_x
disp_y = disp_y
secondary_variable = disp_x
use_physical_gap = true
correct_edge_dropping = true
set_compressive_traction_to_zero = true
friction_coefficient = 0.0
penalty = 0.0
penalty_friction = 0.0
normal_strength = 'normal_strength'
shear_strength = 'shear_strength'
penalty_stiffness = 1000
power_law_parameter = 2.0
GI_c = 5
GII_c = 5
displacements = 'disp_x disp_y'
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[c_x]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[]
[Postprocessors]
[jump_n]
type = ElementAverageValue
variable = local_normal_jump
block = secondary_lower
[]
[jump_t]
type = ElementAverageValue
variable = local_tangential_jump_1
block = secondary_lower
[]
[traction_n]
type = ElementAverageValue
variable = cohesive_traction_normal
block = secondary_lower
[]
[traction_t]
type = ElementAverageValue
variable = cohesive_traction_tangential_magnitude
block = secondary_lower
[]
[traction_eff]
type = ElementAverageValue
variable = cohesive_traction_effective
block = secondary_lower
[]
[cohesive_damage]
type = ElementAverageValue
variable = damage
block = secondary_lower
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
line_search = none
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
automatic_scaling = true
compute_scaling_once = false
off_diagonals_in_auto_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 100
nl_rel_tol = 1e-12
nl_abs_tol = 1e-11
start_time = 0.0
dt = 0.01
end_time = 1.0
dtmin = 0.01
[]
[Outputs]
exodus = false
csv = true
[]
(modules/contact/test/tests/cohesive_zone_model/bilinear_mixed_compare.i)
[Mesh]
[base]
type = GeneratedMeshGenerator
dim = 2
xmax = 1.0
ymax = 1
xmin = -0.0
nx = 1
ny = 1
[]
[rename_base]
type = RenameBoundaryGenerator
input = base
old_boundary = 'top bottom left right'
new_boundary = 'top_base bottom_base left_base right_base'
[]
[base_id]
type = SubdomainIDGenerator
input = rename_base
subdomain_id = 1
[]
[top]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = 1
ymax = 2
nx = 1
ny = 1
[]
[rename_top]
type = RenameBoundaryGenerator
input = top
old_boundary = 'top bottom left right'
new_boundary = '100 101 102 103'
[]
[top_id]
type = SubdomainIDGenerator
input = rename_top
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'base_id top_id'
[]
[top_node]
type = ExtraNodesetGenerator
coord = '0 2 0'
input = combined
new_boundary = top_node
[]
[bottom_node]
type = ExtraNodesetGenerator
coord = '-0.0 0 0'
input = top_node
new_boundary = bottom_node
[]
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'top_base'
input = bottom_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = '101'
new_block_name = 'primary_lower'
input = secondary
[]
patch_update_strategy = auto
patch_size = 20
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
generate_output = 'stress_yy vonmises_stress stress_xy'
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = TaylorExpansion
generate_output = 'vonmises_stress stress_yy stress_xy'
block = '1 2'
[]
[]
[]
[]
[BCs]
[fix_x]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_node
variable = disp_x
[]
[fix_top]
type = DirichletBC
preset = true
boundary = 100
variable = disp_x
value = 0
[]
[top]
type = FunctionDirichletBC
boundary = 100
variable = disp_y
function = 'if(t<=0.3,t,if(t<=0.6,0.3-(t-0.3),0.6-t))'
preset = true
[]
[bottom]
type = DirichletBC
boundary = bottom_base
variable = disp_y
value = 0
preset = true
[]
[]
[Materials]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric9
C_ijkl = '1.684e5 0.176e5 0.176e5 1.684e5 0.176e5 1.684e5 0.754e5 0.754e5 0.754e5'
block = '1 2'
[]
[]
[Postprocessors]
[stress_yy]
type = ElementExtremeValue
variable = stress_yy
value_type = max
block = '1 2'
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
line_search = 'none'
# petsc_options = '-pc_svd_monitor -ksp_monitor_singular_values'
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'svd superlu_dist'
automatic_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 150
nl_rel_tol = 1e-14
nl_abs_tol = 1e-12
start_time = 0.0
dt = 0.01
end_time = 0.85
dtmin = 0.01
[]
[Outputs]
exodus = true
csv = true
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
correct_edge_dropping = true
disp_x = disp_x
disp_y = disp_y
friction_coefficient = 0.0 # with 2.0 works
secondary_variable = disp_x
penalty = 0e6
penalty_friction = 0e4
use_physical_gap = true
# bilinear parameters
normal_strength = 1e4
shear_strength = 1e3
penalty_stiffness = 1e6
power_law_parameter = 2.2
viscosity = 1.0e-3
GI_c = 1e3
GII_c = 1e2
displacements = 'disp_x disp_y'
[]
[]
[Constraints]
[c_x]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[]
(modules/contact/test/tests/cohesive_zone_model/traction_curve_mixed.i)
# Purpose: Verify mixed-mode traction-separation response for POWER_LAW and CLI-overridden BK criteria.
# Assertion: Mixed-mode jump-traction curves for both criteria match their CSV gold baselines within tolerance.
# Why this exists: Confirms criterion-specific mixed-mode damage behavior remains stable across constitutive options.
[Mesh]
[base]
type = GeneratedMeshGenerator
dim = 2
xmax = 1.1
ymax = 1
xmin = -0.1
nx = 1
ny = 1
[]
[rename_base]
type = RenameBoundaryGenerator
input = base
old_boundary = 'top bottom left right'
new_boundary = 'top_base bottom_base left_base right_base'
[]
[base_id]
type = SubdomainIDGenerator
input = rename_base
subdomain_id = 1
[]
[top]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = 1
ymax = 2
nx = 1
ny = 1
[]
[rename_top]
type = RenameBoundaryGenerator
input = top
old_boundary = 'top bottom left right'
new_boundary = '100 101 102 103'
[]
[top_id]
type = SubdomainIDGenerator
input = rename_top
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'base_id top_id'
[]
[top_node]
type = ExtraNodesetGenerator
coord = '0 2 0'
input = combined
new_boundary = top_node
[]
[bottom_node]
type = ExtraNodesetGenerator
coord = '-0.1 0 0'
input = top_node
new_boundary = bottom_node
[]
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'top_base'
input = bottom_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = '101'
new_block_name = 'primary_lower'
input = secondary
[]
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = TaylorExpansion
block = '1 2'
[]
[]
[]
[]
[AuxVariables]
[mode_mixity_ratio]
[]
[damage]
[]
[local_normal_jump]
[]
[local_tangential_jump_1]
[]
[cohesive_traction_normal]
[]
[cohesive_traction_tangential_magnitude]
[]
[cohesive_traction_effective]
[]
[]
[AuxKernels]
[mode_mixity_ratio]
type = CohesiveZoneMortarUserObjectAux
variable = mode_mixity_ratio
user_object = czm_uo
cohesive_zone_quantity = mode_mixity_ratio
[]
[cohesive_damage]
type = CohesiveZoneMortarUserObjectAux
variable = damage
user_object = czm_uo
cohesive_zone_quantity = cohesive_damage
[]
[local_normal_jump]
type = CohesiveZoneMortarUserObjectAux
variable = local_normal_jump
user_object = czm_uo
cohesive_zone_quantity = local_normal_jump
[]
[local_tangential_jump_1]
type = CohesiveZoneMortarUserObjectAux
variable = local_tangential_jump_1
user_object = czm_uo
cohesive_zone_quantity = local_tangential_jump_1
[]
[cohesive_traction_normal]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_normal
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_normal
[]
[cohesive_traction_tangential_magnitude]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_tangential_magnitude
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_tangential_magnitude
[]
[cohesive_traction_effective]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_effective
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_effective
[]
[]
[BCs]
[fix_bottom_x]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_base
variable = disp_x
[]
[top_x]
type = FunctionDirichletBC
boundary = 100
variable = disp_x
function = '0.2*t'
preset = true
[]
[top_y]
type = FunctionDirichletBC
boundary = 100
variable = disp_y
function = '0.2*t'
preset = true
[]
[roller_bottom_y]
type = DirichletBC
boundary = bottom_base
variable = disp_y
value = 0
preset = true
[]
[]
[Materials]
[normal_strength]
type = GenericConstantMaterial
prop_names = 'normal_strength'
prop_values = '50'
[]
[shear_strength]
type = GenericConstantMaterial
prop_names = 'shear_strength'
prop_values = '50'
[]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric_isotropic_E_nu
C_ijkl = '1e8 0.3'
block = '1 2'
[]
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
disp_x = disp_x
disp_y = disp_y
secondary_variable = disp_x
use_physical_gap = true
correct_edge_dropping = true
set_compressive_traction_to_zero = true
friction_coefficient = 0.0
penalty = 0.0
penalty_friction = 0.0
normal_strength = 'normal_strength'
shear_strength = 'shear_strength'
penalty_stiffness = 1000
power_law_parameter = 2.0
mixed_mode_criterion = POWER_LAW
GI_c = 5
GII_c = 2
displacements = 'disp_x disp_y'
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[c_x]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[]
[Postprocessors]
[mode_mixity]
type = ElementAverageValue
variable = mode_mixity_ratio
block = secondary_lower
[]
[jump_n]
type = ElementAverageValue
variable = local_normal_jump
block = secondary_lower
[]
[jump_t]
type = ElementAverageValue
variable = local_tangential_jump_1
block = secondary_lower
[]
[traction_n]
type = ElementAverageValue
variable = cohesive_traction_normal
block = secondary_lower
[]
[traction_t]
type = ElementAverageValue
variable = cohesive_traction_tangential_magnitude
block = secondary_lower
[]
[traction_eff]
type = ElementAverageValue
variable = cohesive_traction_effective
block = secondary_lower
[]
[cohesive_damage]
type = ElementAverageValue
variable = damage
block = secondary_lower
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
line_search = none
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
automatic_scaling = true
compute_scaling_once = false
off_diagonals_in_auto_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 100
nl_rel_tol = 1e-12
nl_abs_tol = 1e-11
start_time = 0.0
dt = 0.01
end_time = 1.0
dtmin = 0.01
[]
[Outputs]
exodus = false
csv = true
[]
(modules/contact/test/tests/cohesive_zone_model/traction_cycle_modeI.i)
# Purpose: Verify irreversible damage behavior in a mode-I load-unload-reload traction cycle.
# Assertion: Reload response follows damaged stiffness and does not heal, matching the CSV gold cycle curve.
# Why this exists: Confirms damage monotonicity and non-healing assumptions in cyclic cohesive loading.
[Mesh]
[base]
type = GeneratedMeshGenerator
dim = 2
xmax = 1.1
ymax = 1
xmin = -0.1
nx = 1
ny = 1
[]
[rename_base]
type = RenameBoundaryGenerator
input = base
old_boundary = 'top bottom left right'
new_boundary = 'top_base bottom_base left_base right_base'
[]
[base_id]
type = SubdomainIDGenerator
input = rename_base
subdomain_id = 1
[]
[top]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = 1
ymax = 2
nx = 1
ny = 1
[]
[rename_top]
type = RenameBoundaryGenerator
input = top
old_boundary = 'top bottom left right'
new_boundary = '100 101 102 103'
[]
[top_id]
type = SubdomainIDGenerator
input = rename_top
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'base_id top_id'
[]
[top_node]
type = ExtraNodesetGenerator
coord = '0 2 0'
input = combined
new_boundary = top_node
[]
[bottom_node]
type = ExtraNodesetGenerator
coord = '-0.1 0 0'
input = top_node
new_boundary = bottom_node
[]
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'top_base'
input = bottom_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = '101'
new_block_name = 'primary_lower'
input = secondary
[]
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = TaylorExpansion
block = '1 2'
[]
[]
[]
[]
[AuxVariables]
[damage]
[]
[local_normal_jump]
[]
[cohesive_traction_normal]
[]
[]
[AuxKernels]
[cohesive_damage]
type = CohesiveZoneMortarUserObjectAux
variable = damage
user_object = czm_uo
cohesive_zone_quantity = cohesive_damage
[]
[local_normal_jump]
type = CohesiveZoneMortarUserObjectAux
variable = local_normal_jump
user_object = czm_uo
cohesive_zone_quantity = local_normal_jump
[]
[cohesive_traction_normal]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_normal
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_normal
[]
[]
[BCs]
[anchor_bottom_x]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_node
variable = disp_x
[]
[anchor_top_x]
type = DirichletBC
preset = true
boundary = top_node
variable = disp_x
value = 0
[]
[top]
type = FunctionDirichletBC
boundary = 100
variable = disp_y
function = 'if(t<=0.3,0.4*t,if(t<=0.6,0.12-0.4*(t-0.3),0.4*(t-0.6)))'
preset = true
[]
[roller_bottom_y]
type = DirichletBC
boundary = bottom_base
variable = disp_y
value = 0
preset = true
[]
[]
[Materials]
[normal_strength]
type = GenericConstantMaterial
prop_names = 'normal_strength'
prop_values = '50'
[]
[shear_strength]
type = GenericConstantMaterial
prop_names = 'shear_strength'
prop_values = '50'
[]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric_isotropic_E_nu
C_ijkl = '1e8 0.3'
block = '1 2'
[]
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
disp_x = disp_x
disp_y = disp_y
secondary_variable = disp_x
use_physical_gap = true
correct_edge_dropping = true
set_compressive_traction_to_zero = true
friction_coefficient = 0.0
penalty = 0.0
penalty_friction = 0.0
normal_strength = 'normal_strength'
shear_strength = 'shear_strength'
penalty_stiffness = 1000
power_law_parameter = 2.0
GI_c = 5
GII_c = 5
displacements = 'disp_x disp_y'
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[c_x]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[]
[Postprocessors]
[jump_n]
type = ElementAverageValue
variable = local_normal_jump
block = secondary_lower
[]
[traction_n]
type = ElementAverageValue
variable = cohesive_traction_normal
block = secondary_lower
[]
[cohesive_damage]
type = ElementAverageValue
variable = damage
block = secondary_lower
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
line_search = none
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
automatic_scaling = true
compute_scaling_once = false
off_diagonals_in_auto_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 100
nl_rel_tol = 1e-12
nl_abs_tol = 1e-11
start_time = 0.0
dt = 0.02
end_time = 0.9
dtmin = 0.02
[]
[Outputs]
exodus = false
csv = true
[]
(modules/contact/test/tests/cohesive_zone_model/robustness_brittle_interface.i)
# Purpose: Verify stable behavior when the softening zone is nearly zero (brittle limit).
# Assertion: Traction output stays finite and damage bounded when delta_f is only marginally larger than delta_0, matching the CSV gold baseline.
# Why this exists: Guards solver and constitutive robustness in the near-brittle limit of the bilinear law.
# Parameters: K=1000, N=50 -> delta_0 = 0.05.
# GI_c = 1.3 -> delta_f = 2*1.3/50 = 0.052, giving delta_f/delta_0 = 1.04 (near-brittle).
# Viscosity = 1e-3 is used to regularize the abrupt damage jump and aid convergence.
[Mesh]
[base]
type = GeneratedMeshGenerator
dim = 2
xmax = 1.1
ymax = 1
xmin = -0.1
nx = 1
ny = 1
[]
[rename_base]
type = RenameBoundaryGenerator
input = base
old_boundary = 'top bottom left right'
new_boundary = 'top_base bottom_base left_base right_base'
[]
[base_id]
type = SubdomainIDGenerator
input = rename_base
subdomain_id = 1
[]
[top]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = 1
ymax = 2
nx = 1
ny = 1
[]
[rename_top]
type = RenameBoundaryGenerator
input = top
old_boundary = 'top bottom left right'
new_boundary = '100 101 102 103'
[]
[top_id]
type = SubdomainIDGenerator
input = rename_top
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'base_id top_id'
[]
[top_node]
type = ExtraNodesetGenerator
coord = '0 2 0'
input = combined
new_boundary = top_node
[]
[bottom_node]
type = ExtraNodesetGenerator
coord = '-0.1 0 0'
input = top_node
new_boundary = bottom_node
[]
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'top_base'
input = bottom_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = '101'
new_block_name = 'primary_lower'
input = secondary
[]
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = TaylorExpansion
block = '1 2'
[]
[]
[]
[]
[AuxVariables]
[damage]
[]
[local_normal_jump]
[]
[cohesive_traction_normal]
[]
[cohesive_traction_effective]
[]
[]
[AuxKernels]
[cohesive_damage]
type = CohesiveZoneMortarUserObjectAux
variable = damage
user_object = czm_uo
cohesive_zone_quantity = cohesive_damage
[]
[local_normal_jump]
type = CohesiveZoneMortarUserObjectAux
variable = local_normal_jump
user_object = czm_uo
cohesive_zone_quantity = local_normal_jump
[]
[cohesive_traction_normal]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_normal
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_normal
[]
[cohesive_traction_effective]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_effective
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_effective
[]
[]
[BCs]
[anchor_bottom_x]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_node
variable = disp_x
[]
[anchor_top_x]
type = DirichletBC
preset = true
boundary = top_node
variable = disp_x
value = 0
[]
[top]
type = FunctionDirichletBC
boundary = 100
variable = disp_y
function = '0.3*t'
preset = true
[]
[roller_bottom_y]
type = DirichletBC
boundary = bottom_base
variable = disp_y
value = 0
preset = true
[]
[]
[Materials]
[normal_strength]
type = GenericConstantMaterial
prop_names = 'normal_strength'
prop_values = '50'
[]
[shear_strength]
type = GenericConstantMaterial
prop_names = 'shear_strength'
prop_values = '50'
[]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric_isotropic_E_nu
C_ijkl = '1e8 0.3'
block = '1 2'
[]
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
disp_x = disp_x
disp_y = disp_y
secondary_variable = disp_x
use_physical_gap = true
correct_edge_dropping = true
set_compressive_traction_to_zero = true
friction_coefficient = 0.0
penalty = 0.0
penalty_friction = 0.0
normal_strength = 'normal_strength'
shear_strength = 'shear_strength'
penalty_stiffness = 1000
power_law_parameter = 2.0
GI_c = 1.3
GII_c = 1.3
viscosity = 1e-3
displacements = 'disp_x disp_y'
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[c_x]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[]
[Postprocessors]
[jump_n]
type = ElementAverageValue
variable = local_normal_jump
block = secondary_lower
[]
[traction_n]
type = ElementAverageValue
variable = cohesive_traction_normal
block = secondary_lower
[]
[traction_eff]
type = ElementAverageValue
variable = cohesive_traction_effective
block = secondary_lower
[]
[cohesive_damage]
type = ElementAverageValue
variable = damage
block = secondary_lower
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
line_search = none
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
automatic_scaling = true
compute_scaling_once = false
off_diagonals_in_auto_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 100
nl_rel_tol = 1e-12
nl_abs_tol = 1e-11
start_time = 0.0
dt = 0.01
end_time = 1.0
dtmin = 0.01
[]
[Outputs]
exodus = false
csv = true
[]
(modules/contact/test/tests/cohesive_zone_model/traction_curve_modeII.i)
# Purpose: Verify mode-II dominant traction-separation behavior with tangential jump and shear traction outputs.
# Assertion: Tangential jump-traction response matches the CSV gold reference within tight tolerance.
# Why this exists: Protects shear-dominant cohesive constitutive behavior against regressions.
[Mesh]
[base]
type = GeneratedMeshGenerator
dim = 2
xmax = 1.1
ymax = 1
xmin = -0.1
nx = 1
ny = 1
[]
[rename_base]
type = RenameBoundaryGenerator
input = base
old_boundary = 'top bottom left right'
new_boundary = 'top_base bottom_base left_base right_base'
[]
[base_id]
type = SubdomainIDGenerator
input = rename_base
subdomain_id = 1
[]
[top]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = 1
ymax = 2
nx = 1
ny = 1
[]
[rename_top]
type = RenameBoundaryGenerator
input = top
old_boundary = 'top bottom left right'
new_boundary = '100 101 102 103'
[]
[top_id]
type = SubdomainIDGenerator
input = rename_top
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'base_id top_id'
[]
[top_left_node]
type = ExtraNodesetGenerator
coord = '0 2 0'
input = combined
new_boundary = top_left_node
[]
[top_right_node]
type = ExtraNodesetGenerator
coord = '1 2 0'
input = top_left_node
new_boundary = top_right_node
[]
[bottom_node]
type = ExtraNodesetGenerator
coord = '-0.1 0 0'
input = top_right_node
new_boundary = bottom_node
[]
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'top_base'
input = bottom_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = '101'
new_block_name = 'primary_lower'
input = secondary
[]
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = TaylorExpansion
block = '1 2'
[]
[]
[]
[]
[AuxVariables]
[mode_mixity_ratio]
[]
[damage]
[]
[local_normal_jump]
[]
[local_tangential_jump_1]
[]
[cohesive_traction_normal]
[]
[cohesive_traction_tangential_magnitude]
[]
[cohesive_traction_effective]
[]
[]
[AuxKernels]
[mode_mixity_ratio]
type = CohesiveZoneMortarUserObjectAux
variable = mode_mixity_ratio
user_object = czm_uo
cohesive_zone_quantity = mode_mixity_ratio
[]
[cohesive_damage]
type = CohesiveZoneMortarUserObjectAux
variable = damage
user_object = czm_uo
cohesive_zone_quantity = cohesive_damage
[]
[local_normal_jump]
type = CohesiveZoneMortarUserObjectAux
variable = local_normal_jump
user_object = czm_uo
cohesive_zone_quantity = local_normal_jump
[]
[local_tangential_jump_1]
type = CohesiveZoneMortarUserObjectAux
variable = local_tangential_jump_1
user_object = czm_uo
cohesive_zone_quantity = local_tangential_jump_1
[]
[cohesive_traction_normal]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_normal
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_normal
[]
[cohesive_traction_tangential_magnitude]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_tangential_magnitude
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_tangential_magnitude
[]
[cohesive_traction_effective]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_effective
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_effective
[]
[]
[BCs]
[fix_bottom_x]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_base
variable = disp_x
[]
[top]
type = FunctionDirichletBC
boundary = 100
variable = disp_x
function = '0.3*t'
preset = true
[]
[top_y_rollers]
type = DirichletBC
boundary = 'top_left_node top_right_node'
variable = disp_y
value = 0
preset = true
[]
[roller_bottom_y]
type = DirichletBC
boundary = bottom_base
variable = disp_y
value = 0
preset = true
[]
[]
[Materials]
[normal_strength]
type = GenericConstantMaterial
prop_names = 'normal_strength'
prop_values = '50'
[]
[shear_strength]
type = GenericConstantMaterial
prop_names = 'shear_strength'
prop_values = '50'
[]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric_isotropic_E_nu
C_ijkl = '1e8 0.3'
block = '1 2'
[]
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
disp_x = disp_x
disp_y = disp_y
secondary_variable = disp_x
use_physical_gap = true
correct_edge_dropping = true
set_compressive_traction_to_zero = true
friction_coefficient = 0.0
penalty = 0.0
penalty_friction = 0.0
normal_strength = 'normal_strength'
shear_strength = 'shear_strength'
penalty_stiffness = 1000
power_law_parameter = 2.0
GI_c = 5
GII_c = 5
displacements = 'disp_x disp_y'
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[c_x]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[]
[Postprocessors]
[mode_mixity]
type = ElementAverageValue
variable = mode_mixity_ratio
block = secondary_lower
[]
[jump_n]
type = ElementAverageValue
variable = local_normal_jump
block = secondary_lower
[]
[jump_t]
type = ElementAverageValue
variable = local_tangential_jump_1
block = secondary_lower
[]
[traction_n]
type = ElementAverageValue
variable = cohesive_traction_normal
block = secondary_lower
[]
[traction_t]
type = ElementAverageValue
variable = cohesive_traction_tangential_magnitude
block = secondary_lower
[]
[traction_eff]
type = ElementAverageValue
variable = cohesive_traction_effective
block = secondary_lower
[]
[cohesive_damage]
type = ElementAverageValue
variable = damage
block = secondary_lower
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
line_search = none
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
automatic_scaling = true
compute_scaling_once = false
off_diagonals_in_auto_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 100
nl_rel_tol = 1e-12
nl_abs_tol = 1e-11
start_time = 0.0
dt = 0.01
end_time = 1.0
dtmin = 0.01
[]
[Outputs]
exodus = false
csv = true
[]
(modules/contact/examples/mortar_czm/debond_recontact.i)
# Purpose: Verify that a fully debonded interface (omega=1) re-establishes compressive contact on closure.
# Verifies: damage stays at 1 after full decohesion; compressive traction is nonzero on recontact via penalty contact.
[Mesh]
[base]
type = GeneratedMeshGenerator
dim = 2
xmax = 1.1
ymax = 1
xmin = -0.1
nx = 1
ny = 1
[]
[rename_base]
type = RenameBoundaryGenerator
input = base
old_boundary = 'top bottom left right'
new_boundary = 'top_base bottom_base left_base right_base'
[]
[base_id]
type = SubdomainIDGenerator
input = rename_base
subdomain_id = 1
[]
[top]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = 1
ymax = 2
nx = 1
ny = 1
[]
[rename_top]
type = RenameBoundaryGenerator
input = top
old_boundary = 'top bottom left right'
new_boundary = '100 101 102 103'
[]
[top_id]
type = SubdomainIDGenerator
input = rename_top
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'base_id top_id'
[]
[top_node]
type = ExtraNodesetGenerator
coord = '0 2 0'
input = combined
new_boundary = top_node
[]
[bottom_node]
type = ExtraNodesetGenerator
coord = '-0.1 0 0'
input = top_node
new_boundary = bottom_node
[]
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'top_base'
input = bottom_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = '101'
new_block_name = 'primary_lower'
input = secondary
[]
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = TaylorExpansion
block = '1 2'
[]
[]
[]
[]
[AuxVariables]
[damage]
[]
[local_normal_jump]
[]
[cohesive_traction_normal]
[]
[]
[AuxKernels]
[cohesive_damage]
type = CohesiveZoneMortarUserObjectAux
variable = damage
user_object = czm_uo
cohesive_zone_quantity = cohesive_damage
[]
[local_normal_jump]
type = CohesiveZoneMortarUserObjectAux
variable = local_normal_jump
user_object = czm_uo
cohesive_zone_quantity = local_normal_jump
[]
[cohesive_traction_normal]
type = CohesiveZoneMortarUserObjectAux
variable = cohesive_traction_normal
user_object = czm_uo
cohesive_zone_quantity = cohesive_traction_normal
[]
[]
[BCs]
[anchor_bottom_x]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_node
variable = disp_x
[]
[anchor_top_x]
type = DirichletBC
preset = true
boundary = top_node
variable = disp_x
value = 0
[]
# Open past delta_f (=0.2mm at K=1000, N=50, GI_c=5) then reverse into compression.
# delta_f = 2*GI_c/N = 0.2; loading reaches 0.24 (full debond) then reverses to -0.096.
[top]
type = FunctionDirichletBC
boundary = 100
variable = disp_y
function = 'if(t<=0.8, 0.3*t, 0.24 - 0.48*(t-0.8))'
preset = true
[]
[roller_bottom_y]
type = DirichletBC
boundary = bottom_base
variable = disp_y
value = 0
preset = true
[]
[]
[Materials]
[normal_strength]
type = GenericConstantMaterial
prop_names = 'normal_strength'
prop_values = '50'
[]
[shear_strength]
type = GenericConstantMaterial
prop_names = 'shear_strength'
prop_values = '50'
[]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric_isotropic_E_nu
C_ijkl = '1e8 0.3'
block = '1 2'
[]
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
disp_x = disp_x
disp_y = disp_y
secondary_variable = disp_x
use_physical_gap = true
correct_edge_dropping = true
set_compressive_traction_to_zero = true
friction_coefficient = 0.0
penalty = 1e4
penalty_friction = 0.0
normal_strength = 'normal_strength'
shear_strength = 'shear_strength'
penalty_stiffness = 1000
power_law_parameter = 2.0
GI_c = 5
GII_c = 5
displacements = 'disp_x disp_y'
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[c_x]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[]
[Postprocessors]
[jump_n]
type = ElementAverageValue
variable = local_normal_jump
block = secondary_lower
[]
[traction_n]
type = ElementAverageValue
variable = cohesive_traction_normal
block = secondary_lower
[]
[cohesive_damage]
type = ElementAverageValue
variable = damage
block = secondary_lower
[]
# Total interface traction (CZM + penalty contact) from bulk stress equilibrium.
# After full debond (omega=1), CZM traction_n = 0, but this captures penalty contact compression.
[total_traction_y]
type = ADSidesetReaction
stress_tensor = stress
boundary = top_base
direction = '0 1 0'
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
line_search = none
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
automatic_scaling = true
compute_scaling_once = false
off_diagonals_in_auto_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 100
nl_rel_tol = 1e-12
nl_abs_tol = 1e-11
start_time = 0.0
dt = 0.02
end_time = 1.36
dtmin = 0.02
[]
[Outputs]
exodus = false
csv = true
[]
(modules/contact/test/tests/cohesive_zone_model/cohesive_friction.i)
[Mesh]
[base]
type = GeneratedMeshGenerator
dim = 2
xmax = 1.1
ymax = 1
xmin = -0.1
nx = 1
ny = 1
[]
[rename_base]
type = RenameBoundaryGenerator
input = base
old_boundary = 'top bottom left right'
new_boundary = 'top_base bottom_base left_base right_base'
[]
[base_id]
type = SubdomainIDGenerator
input = rename_base
subdomain_id = 1
[]
[top]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = 1
ymax = 2
nx = 1
ny = 1
[]
[rename_top]
type = RenameBoundaryGenerator
input = top
old_boundary = 'top bottom left right'
new_boundary = '100 101 102 103'
[]
[top_id]
type = SubdomainIDGenerator
input = rename_top
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'base_id top_id'
[]
[top_node]
type = ExtraNodesetGenerator
coord = '0 2 0'
input = combined
new_boundary = top_node
[]
[bottom_node]
type = ExtraNodesetGenerator
coord = '-0.1 0 0'
input = top_node
new_boundary = bottom_node
[]
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'top_base'
input = bottom_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = '101'
new_block_name = 'primary_lower'
input = secondary
[]
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[AuxVariables]
[penalty_normal_pressure]
order = FIRST
family = LAGRANGE
[]
[mode_mixity_ratio]
[]
[damage]
[]
[local_normal_jump]
[]
[local_tangential_jump]
[]
[]
[AuxKernels]
[penalty_normal_pressure_auxk]
type = MortarUserObjectAux
variable = penalty_normal_pressure
user_object = czm_uo
contact_quantity = normal_pressure
[]
[mode_mixity_ratio]
type = CohesiveZoneMortarUserObjectAux
variable = mode_mixity_ratio
user_object = czm_uo
cohesive_zone_quantity = mode_mixity_ratio
[]
[cohesive_damage]
type = CohesiveZoneMortarUserObjectAux
variable = damage
user_object = czm_uo
cohesive_zone_quantity = cohesive_damage
[]
[local_normal_jump]
type = CohesiveZoneMortarUserObjectAux
variable = local_normal_jump
user_object = czm_uo
cohesive_zone_quantity = local_normal_jump
[]
[local_tangential_jump]
type = CohesiveZoneMortarUserObjectAux
variable = local_tangential_jump
user_object = czm_uo
cohesive_zone_quantity = local_tangential_jump
[]
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
generate_output = 'stress_yy'
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = EigenSolution
generate_output = 'vonmises_stress'
block = '1 2'
[]
[]
[]
[]
[BCs]
[fix_x]
type = DirichletBC
preset = false
value = 0.0
boundary = bottom_node
variable = disp_x
[]
[fix_top]
type = FunctionDirichletBC
preset = false
boundary = 100
variable = disp_x
function = 0.2*t
[]
[top]
type = FunctionDirichletBC
boundary = 100
variable = disp_y
function = -0.01
preset = false
[]
[bottom]
type = DirichletBC
boundary = bottom_base
variable = disp_y
value = 0
preset = false
[]
[]
[Materials]
[normal_strength]
type = GenericConstantMaterial
prop_names = 'normal_strength'
prop_values = '1e3'
[]
[shear_strength]
type = GenericConstantMaterial
prop_names = 'shear_strength'
prop_values = '100'
[]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric_isotropic_E_nu
C_ijkl = '1e6 0.3'
block = '1'
[]
[elasticity_tensor2]
type = ADComputeElasticityTensor
fill_method = symmetric_isotropic_E_nu
C_ijkl = '1e5 0.3'
block = '2'
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
automatic_scaling = true
compute_scaling_once = false
off_diagonals_in_auto_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 20
nl_rel_tol = 1e-12
nl_abs_tol = 1e-9
start_time = 0.0
dt = 0.1
end_time = 1.6
[]
[Outputs]
exodus = true
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
disp_x = disp_x
disp_y = disp_y
secondary_variable = disp_x
use_physical_gap = true
set_compressive_traction_to_zero = true
correct_edge_dropping = true
penalty_friction = 5e6
penalty = 1e5
friction_coefficient = 0.1
# bilinear model parameters
normal_strength = 'normal_strength'
shear_strength = 'shear_strength'
penalty_stiffness = 1e6
power_law_parameter = 2
GI_c = 123
GII_c = 54
displacements = 'disp_x disp_y'
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[t_x]
type = TangentialMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_velocities_uo = czm_uo
[]
[t_y]
type = TangentialMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_velocities_uo = czm_uo
[]
[c_x]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[]
(modules/contact/test/tests/cohesive_zone_model/bilinear_mixed_mortar_only_czm.i)
[Mesh]
[base]
type = GeneratedMeshGenerator
dim = 2
xmax = 1.1
ymax = 1
xmin = -0.1
nx = 1
ny = 1
[]
[rename_base]
type = RenameBoundaryGenerator
input = base
old_boundary = 'top bottom left right'
new_boundary = 'top_base bottom_base left_base right_base'
[]
[base_id]
type = SubdomainIDGenerator
input = rename_base
subdomain_id = 1
[]
[top]
type = GeneratedMeshGenerator
dim = 2
xmax = 1
ymin = 1
ymax = 2
nx = 1
ny = 1
[]
[rename_top]
type = RenameBoundaryGenerator
input = top
old_boundary = 'top bottom left right'
new_boundary = '100 101 102 103'
[]
[top_id]
type = SubdomainIDGenerator
input = rename_top
subdomain_id = 2
[]
[combined]
type = MeshCollectionGenerator
inputs = 'base_id top_id'
[]
[top_node]
type = ExtraNodesetGenerator
coord = '0 2 0'
input = combined
new_boundary = top_node
[]
[bottom_node]
type = ExtraNodesetGenerator
coord = '-0.1 0 0'
input = top_node
new_boundary = bottom_node
[]
[secondary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10001
new_block_name = 'secondary_lower'
sidesets = 'top_base'
input = bottom_node
[]
[primary]
type = LowerDBlockFromSidesetGenerator
new_block_id = 10000
sidesets = '101'
new_block_name = 'primary_lower'
input = secondary
[]
patch_update_strategy = auto
patch_size = 20
allow_renumbering = false
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Physics]
[SolidMechanics]
[QuasiStatic]
generate_output = 'stress_yy'
[all]
strain = FINITE
add_variables = true
use_automatic_differentiation = true
decomposition_method = TaylorExpansion
generate_output = 'vonmises_stress'
block = '1 2'
[]
[]
[]
[]
[BCs]
[fix_x]
type = DirichletBC
preset = true
value = 0.0
boundary = bottom_node
variable = disp_x
[]
[fix_top]
type = DirichletBC
preset = true
boundary = 100
variable = disp_x
value = 0
[]
[top]
type = FunctionDirichletBC
boundary = 100
variable = disp_y
function = 'if(t<=0.3,t,if(t<=0.6,0.3-(t-0.3),0.6-t))'
preset = true
[]
[bottom]
type = DirichletBC
boundary = bottom_base
variable = disp_y
value = 0
preset = true
[]
[]
[Materials]
[normal_strength]
type = GenericConstantMaterial
prop_names = 'normal_strength'
prop_values = '1e3'
[]
[shear_strength]
type = GenericConstantMaterial
prop_names = 'shear_strength'
prop_values = '7.5e2'
[]
[stress]
type = ADComputeFiniteStrainElasticStress
block = '1 2'
[]
[elasticity_tensor]
type = ADComputeElasticityTensor
fill_method = symmetric9
C_ijkl = '1.684e5 0.176e5 0.176e5 1.684e5 0.176e5 1.684e5 0.754e5 0.754e5 0.754e5'
block = '1 2'
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
line_search = none
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
automatic_scaling = true
compute_scaling_once = false
off_diagonals_in_auto_scaling = true
l_max_its = 2
l_tol = 1e-14
nl_max_its = 150
nl_rel_tol = 1e-14
nl_abs_tol = 1e-12
start_time = 0.0
dt = 0.1
end_time = 1.0
dtmin = 0.1
[]
[Outputs]
exodus = true
[]
[UserObjects]
[czm_uo]
type = BilinearMixedModeCohesiveZoneModel
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
correct_edge_dropping = true
disp_x = disp_x
disp_y = disp_y
friction_coefficient = 0.0 # with 2.0 works
secondary_variable = disp_x
penalty = 0e6
penalty_friction = 0e4
use_physical_gap = true
set_compressive_traction_to_zero = true
# bilinear model parameters
normal_strength = 'normal_strength'
shear_strength = 'shear_strength'
penalty_stiffness = 200
power_law_parameter = 0.1
GI_c = 123
GII_c = 54
displacements = 'disp_x disp_y'
[]
[]
[Constraints]
[x]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[y]
type = NormalMortarMechanicalContact
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
weighted_gap_uo = czm_uo
[]
[c_x]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_x
component = x
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[c_y]
type = MortarGenericTraction
primary_boundary = 101
secondary_boundary = 'top_base'
primary_subdomain = 10000
secondary_subdomain = 10001
secondary_variable = disp_y
component = y
use_displaced_mesh = true
compute_lm_residuals = false
cohesive_zone_uo = czm_uo
[]
[]