Virial stress #137
Virial stress #137
Conversation
…ine 109 for virialStress
… energy function for virial stress calculations
…accordingly so we do not recalculate force, needed to add a const version of sliceForce() to avoid compiling errors, virial stress output is still zero for everything.
|
I've made the suggested changes and it is now working with non-fracture problems using the PMB method. (I may have pushed an altered version of elastic wave problem that I used to test it...please disregard this). Will move on to implementing virial stress with fracture. |
Sounds good. One thing that will help quite a bit is using One other note - we'll need to remove the comments tracking incremental changes before we merge (the fact that the force slice was moved outside the kernel is tracked in the git history and would likely be more confusing than not in the future, for example) |
|
1-s2.0-S0020768322000427-main.pdf Added virial stress to pmb with fracture. Could you please doublecheck the volume summation. it should be normalised by the total volume of the horizon regardless of the number of points within the horizon. I'm thinking we might need to adjust the non-fracture volume too as it does not include vol(i) in the division. |
@pabloseleson can you take a look? |
|
Which equation is being implemented exactly? |
Eq. 29. |
|
@kinanbezem7 : I looked at Eq. (29) in the reference here. There is something awkward in that equation. The RHS of the equation is for the stress tensor at X. However, the LHS of the equation doesn't depend on X. Instead, I looked at an alternative reference here and implemented Eq. (24) of that reference. In this reference, there is a dependence on R on both the RHS and LHS, so it makes more sense to me. BTW: this stress tensor is in general nonsymmetric (like the classical Piola stress tensor, see here). So, it seems we would need to compute the full tensor, not only the 6 entries. I only included the 6 entries for now. |
|
@streeve : Please take a look at the current implementation. |
There was a problem hiding this comment.
Couple small things to clean up, but this is close to ready. @kinanbezem7 let me know if it makes sense how to fix the failing test or remaining comments or we can look at it in person
| @@ -397,7 +507,7 @@ class Force<MemorySpace, ForceModel<LinearPMB, Elastic, ModelParams...>> | |||
| template <class PosType, class WType, class ParticleType, | |||
| class ParallelType> | |||
| double computeEnergyFull( WType& W, const PosType& x, const PosType& u, | |||
There was a problem hiding this comment.
@kinanbezem7 Last changes needed are linear PMB + LPS stresses that match PMB
There was a problem hiding this comment.
I'm trying to run the crack_branching problem just to check that stress is working. but i don't see it in the SILO output files. Does it not fall under OtherFields here in particles.hpp?
#ifdef Cabana_ENABLE_SILO
Cabana::Grid::Experimental::SiloParticleOutput::
writePartialRangeTimeStep(
"particles", local_grid->globalGrid(), output_step, output_time,
0, n_local, getPosition( use_reference ), sliceForce(),
sliceDisplacement(), sliceVelocity(),
std::forward<OtherFields>( other )... );
There was a problem hiding this comment.
You've reminded me of a Silo-specific issue with matrix fields. As I recall they just don't show up in paraview. If you open an issue on your hdf5 problems I can try to help there.
Added virial stress to particles.hpp and the HDF5 output for the pmb model.
Added logic for calculating virial stress in compute full energy in force_pmb.hpp.
Compiles and runs kalthoff winkler, however the virial stress in the output is zero.
I think the forces may not be being calculated properly in the energy logic. I'm also unsure of how to handle bond breakage here, but my entire approach may be wrong so this might be totally unnecessary