Inheritance and Polymorphism in the Source Code of the Pair_Style granular

Dear all,

I am currently attempting to modify the “hertz/material” model within the pair_style command, but I am encountering difficulties in understanding the source code.
Specifically, in the gran_sub_mod_normal.cpp, there is the GranSubModNormalHertzMaterial class, which inherits from the GranSubModNormalHertz class. However, I cannot locate where the GranSubModNormalHertzMaterial class can be applied. Could someone assist me with this?

I am using the August 2, 2023 version of LAMMPS.

Thank you for your help.

These classes are instantiated in the PairGranular::coeff() function.
They are managed by the GranularModel class which contains factories that create the requested instances of the sub-models.

It is a bit difficult to sort out, since it involved C preprocessing and C++ templates. The mechanism, however, is not different from how other LAMMPS styles are mapped to their names.

If you are doing modifications, you best create a new submodel it will then be automatically included during the build process and can be referenced by its name in the granular pair style’s coeff command.

If you haven’t already found it, there’s also a documentation describing how to add new contact models:

https://docs.lammps.org/Modify_gran_sub_mod.html

Dear akohlmey

Thank you for your response. With your guidance and the approach suggested by jtclemm, I successfully created a new submodel and generated the executable file. Now, I have a new question.

When calculating the normal contact force, I need to determine the contact force formula based on the value of the delta variable. This variable is defined in the GranularModel class in the granular_model.h file and represents the particles overlap. Therefore, I need to know the situation of delta from the previous time step and compare it with the current time step. Is there a good way to address this issue?

Well, what have you tried up to now? I would have imagined something like this would work:

// ...
// in relevant loop
delta_old = delta;
delta = calculate_new_delta();
result = compare(delta_old, delta);
// ...

To save information from a previous timestep (e.g. deltaold) for each pairwise contact, you’ll have to create a new ‘history’ entry.

This is mostly done in frictional submodels to track the sliding/rolling/twisting history of the contact. The classes in gran_sub_mod_tangential.cpp are good examples. They first define the size of the history array in the constructor (size_history = 1 for your case) then grab data from gm->history[history_index] where history_index is an automatically generated index that points to the location of data for that specific submodel (in case more than one submodel store history). LAMMPS then handles everything else (repopulating arrays when neighbor lists are built, communicating across processors, saving to restarts, etc.).

Two details to keep in mind:

  1. LAMMPS sometimes calls force calculations outside of the normal timestep loop. The variable gm->history_update indicates when history values should be updated.
  2. By default, LAMMPS assumes the history stores geometric vectors. So two neighboring atoms A and B swap their order in the neighbor list to B and A, values in the history array are negated. This can be overridden by defining a nondefault history transfer factor. See GranSubModTangentialMindlinRescale for an example. It stores 4 values, the first three are a vector (so the factor is -1, the default) while the 4th is a scalar (the factor is +1 so it doesn’t change sign if atoms swap).

Hopefully this helps. The granular pairstyle recently got an overhaul to simplify the process of creating new contact models. While this automated many tasks and cut out a web of if statements which was difficult to navigate, it created a lot of obscure references to the gm class.

Dear jtclemm

Thank you for your patient response. Your answers have been very helpful to me. I will get back to you after understanding and trying the methods you mentioned.