idaholab / moose

Multiphysics Object Oriented Simulation Environment
https://www.mooseframework.org
GNU Lesser General Public License v2.1
1.77k stars 1.05k forks source link

SingularShapeTensorEliminatorUserObjectPD is not called when other UserObjects are active #24239

Open andrisfreimanis opened 1 year ago

andrisfreimanis commented 1 year ago

Bug Description

SingularShapeTensorEliminatorUserObjectPD is not called when other UserObjects, such as NodalDamageIndexPD and NodalNumIntactBondsPD, are active.

Steps to Reproduce

I modified Peridynamics module test "failure_tests/2D_singular_shape_tensor_H1NOSPD.i" by adding either NodalDamageIndexPD, or NodalNumIntactBondsPD UserObject. Enabling "singular_shape_tensor" and either "damage" or "intact_bonds" shows that the problem is with both UOs. I copied the two blocks below and the full input files are attached.

[AuxVariables]
  [critical_stress]
    family = MONOMIAL
    order = CONSTANT
  []
  [damage]
  []
  [intact_bonds_num]
  []
[]

[UserObjects]
  active = 'singular_shape_tensor damage'
  [singular_shape_tensor]
    type = SingularShapeTensorEliminatorUserObjectPD
  []
  [damage]
    type = NodalDamageIndexPD
    variable = damage
  []
  [intact_bonds]
    type = NodalNumIntactBondsPD
    variable = intact_bonds_num
  []
[]

Running both files with "show_actions=true" shows that UOs are added

[DBG][ACT] TASK (         add_user_object) TYPE (             AddUserObjectAction) NAME (singular_shape_tensor) Memory usage 142MB
[DBG][ACT] TASK (         add_user_object) TYPE (             AddUserObjectAction) NAME (    intact_bonds) Memory usage 142MB

and

[DBG][ACT] TASK (         add_user_object) TYPE (             AddUserObjectAction) NAME (singular_shape_tensor) Memory usage 142MB
[DBG][ACT] TASK (         add_user_object) TYPE (             AddUserObjectAction) NAME (          damage) Memory usage 142MB

Enabling " show_execution_order = 'TIMESTEP_END'" shows that "bond_status" and "damage" is executed, but singular_shape_tensor UO is not. The SingularShapeTensorEliminatorUserObjectPD is not run, the shape tensor becomes singular, and an error is thrown. The same (save "damage" changes to "intact_bonds") output is given when "intact_bonds" are enabled.

Currently Executing
  Finished Solving                                                                       [ 41.41 s] [  157 MB]
[DBG] Executing auxiliary kernels on elements on TIMESTEP_END
[DBG] Ordering of AuxKernels on block 1000
[DBG] bond_status
[DBG] Computing elemental user objects on TIMESTEP_END
[DBG] Execution order of objects types on each element then its sides:
[DBG] - element user objects
[DBG] - domain user objects
[DBG] - element user objects contributing to the Jacobian
[DBG] - side user objects
[DBG] - domain user objects executing on sides
[DBG] - side user objects contributing to the Jacobian
[DBG] - internal side user objects
[DBG] - domain user objects executing on internal sides
[DBG] - interface user objects
[DBG] - domain user objects executing at interfaces
[DBG] Ordering of User Objects on block 1000
[DBG] Executing ElementUserObject on TIMESTEP_END
[DBG] Order of execution:
[DBG] damage

[DBG] Only user objects active on local element/sides are executed

*** ERROR ***
The following error occurred in the object "strain", of type "ComputeFiniteStrainNOSPD".

Singular shape tensor is detected in ComputeStrainBaseNOSPD! Use SingularShapeTensorEliminatorUserObjectPD to avoid singular shape tensor!

The expected output of SingularShapeTensorEliminatorUserObjectPD, which h

The expected output with only SingularShapeTensorEliminatorUserObjectPD enabled is:

[DBG] Executing auxiliary kernels on elements on TIMESTEP_END
[DBG] Ordering of AuxKernels on block 1000
[DBG] bond_status
[DBG] Initializing, executing & finalizing general UO 'singular_shape_tensor' on TIMESTEP_END
 Singular shape tensor detected! Elimination in process ... 
     Loop: 1, Singularities: 10
     Loop: 2, Singularities: 0
 Elimination done!
[DBG] Initializing, executing & finalizing general UO 'allGhostElemPD' on TIMESTEP_END

Impact

Not enabling SingularShapeTensorEliminatorUserObjectPD means peridynamic correspondence models will fail with singular shape tensor, when damage is high in some local area. This highly reduces model usability.

NodalDamageIndexPD UO automatically calculates damage at each node. Technically, it could be done as a post-processing step after Moose simulation finishes by extracting node and element data from the results file and calculating it separately. However, it is much more convenient for visualization and saves time if the damage results are available already in the exodus output file.

Not working together with NodalNumIntactBondsPD is annoying, but I doubt critical.

singular_and_dmg_ind_UO.txt singular_and_intact_bonds_UO.txt

lindsayad commented 1 year ago

@bwspenc @dschwen as you guys were kind enough to volunteer for peridynamics ownership, do you have any thoughts on this or who we might dispatch this issue to?

andrisfreimanis commented 1 year ago

I looked into this a bit. The problem here is the execution order of UOs. If at the timestep end the SingularShapeTensorEliminatorUserObjectPD would always be executed first, then it would set the shape tensor to identity before other UOs use the singular one. Instead, PPs and UOs are executed in (I guess) the order they were created, so some other UO forces evaluation of a singular shape tensor, before the SingularShapeTensorEliminatorUserObjectPD can fix it.

Setting execution_order_group = 0 to SingularShapeTensorEliminatorUserObjectPD and execution_order_group = 1 for all other UOs and PPs seems to make it work.

Not sure if this is a bug or working as intended.

dschwen commented 1 year ago

They only get executed in the order they are created if there is no explicit dependency between the UOs. Without looking at the code a getUserObject< SingularShapeTensorEliminatorUserObjectPD >("nameparameter") issued by the UOs that need it to execute first, could possible resolve this dependency automatically. If not, then the execution_order_group flag will be a viable solution, too. We should just be explicit about that in the docs.