Issue with dump files not storing updated values computed for each contact

Submitted by agupta87 on Wed, 07/19/2023 - 22:17

Hey everyone,

I'm a beginner to LIGGGHTS and C++. I'm currently looking to modify LIGGGHTS to compute and store frictional dissipation at the contacts.

To do this, I introduced a new keyword called diss in compute_pair_gran_local.cpp and modified the ComputePairGranLocal::add_pair method accordingly to take in another argument, diss. I created a corresponding flag in this file called dissflag, and if dissflag is 1, then array[ipair][n++] = diss; executes. I then modified calls to this method in pair_gran.h. I never made any changes to pair_gran.cpp.

Because dissipation is computed at each contact, and because I want to store the quantity to be outputted later, I modify contact_interface.h and add a new quantity under the SurfacesIntersectData struct - a mutable double called diss. Then, I modified tangential_model_history.h. Here, I introduced a new history variable using history_offset2 = hsetup->add_history_value("dissVal", "0");. Upon doing this, I then use the following lines of code to attempt to store dissipation:

// under TangentialModel(LAMMPS * lmp, IContactHistorySetup * hsetup,class ContactModelBase *c)
double * const dissValex = &sidata.contact_history[history_offset2];
...
// under inline void surfacesIntersect
double dissPrev;
dissPrev = dissValex[0];
...
const double inelastic1 = shear[0]-(-Ft1/kt);
const double inelastic2 = shear[1]-(-Ft2/kt);
const double inelastic3 = shear[2]-(-Ft3/kt);
const double inelasticmag = sqrt(inelastic1*inelastic1 + inelastic2*inelastic2 + inelastic3*inelastic3);
dissValex[0] = -Ft_friction*inelasticmag;
sidata.diss = dissValex[0];

I print out sidata.diss immediately in tangential_model_history.h, and the values are correct. However, I notice that these values are not transferred appropriately to the dump files, and instead, it just reports "0" as the value, which is obviously incorrect. I notice that the value of sidata.diss is 0 when cpl_add_pair(LCM::SurfacesIntersectData & sidata, LCM::ForceData & i_forces) is called in pair_gran.h. I also notice that there is some copying of SurfacesIntersectData under pair_gran_base.h, so my intuition tells me that I have to somehow modify this file, though I'm unsure exactly how to go about doing this. I notice that a line of interest in pair_gran_base.h is sidata.contact_history = all_contact_hist ? &all_contact_hist[dnum*jj] : NULL;. I'm guessing that this is copying over the history variables that I have previously defined, but I'm unsure how exactly I can access dissVal from that and store in sidata.diss. Perhaps I might be wrong, but I'm very much stuck, so any help would be much appreciated.

Apologies for the long wall of text. I have attached the relevant files, but changed the h and cpp files to txt so they upload.

agupta87 | Fri, 07/21/2023 - 04:15

So it turns out that my suspicions were correct. To ensure that the dissipation value is dumped correctly, all I had to do was to access to all_contact_hist. Because four history variables were added (three for displacement, one for dissipation), I just had to place a statement sidata.diss = all_contact_hist[3] and the dump files report the correct value!