tadd outline for registering new contacts - granular - granular dynamics simulation
 (HTM) git clone git://src.adamsgaard.dk/granular
 (DIR) Log
 (DIR) Files
 (DIR) Refs
 (DIR) README
 (DIR) LICENSE
       ---
 (DIR) commit a83d57fbe2b6759ec904f4047e0460c83cd89dcf
 (DIR) parent 04d5e7eefa6dad607385b133b32056e06b8d3fd4
 (HTM) Author: Anders Damsgaard <anders@adamsgaard.dk>
       Date:   Thu,  8 Apr 2021 22:11:56 +0200
       
       add outline for registering new contacts
       
       Diffstat:
         M contact.h                           |       5 +++--
         M grain.c                             |      29 +++++++++++++++++++++++++++++
         M grain.h                             |       3 +++
         M simulation.c                        |      26 ++++++++++++--------------
       
       4 files changed, 47 insertions(+), 16 deletions(-)
       ---
 (DIR) diff --git a/contact.h b/contact.h
       t@@ -6,9 +6,10 @@
        struct contact {
                int active;
                size_t i, j;
       -        double centerdist[3];
       -        double overlap;
                double age;
       +        double overlap;
       +        double centerdist[3];
       +        double tandisp[3];
                double stress[3];
        };
        
 (DIR) diff --git a/grain.c b/grain.c
       t@@ -314,3 +314,32 @@ grain_temporal_integration(struct grain *g, double dt)
        {
                grain_temporal_integration_two_term_taylor(g, dt);
        }
       +
       +void
       +grain_register_contact(struct grain *g, size_t i, size_t j,
       +                       double centerdist[3], double overlap)
       +{
       +        int ic, d;
       +
       +        /* first pass: check if contact is already registered */
       +        for (ic = 0; ic <= MAXCONTACTS; ic++)
       +                if (g->contacts[ic].j == j) {
       +                        g->contacts[ic].overlap = overlap;
       +                        for (d = 0; d < 3; d++)
       +                                g->contacts[ic].centerdist[d] = centerdist[d];
       +                        return;
       +                }
       +
       +        /* second pass: register as new contact if overlapping */
       +        if (overlap < 0.0)
       +                for (ic = 0; ic <= MAXCONTACTS; ic++)
       +                        if (!g->contacts[ic].active) {
       +                                contact_new(&g->contacts[ic], i, j);
       +                                g->ncontacts += 1;
       +                                g->contacts[ic].overlap = overlap;
       +                                for (d = 0; d < 3; d++) {
       +                                        g->contacts[ic].centerdist[d] = centerdist[d];
       +                                        g->contacts[ic].tandisp[d] = 0.0;
       +                                }
       +                        }
       +}
 (DIR) diff --git a/grain.h b/grain.h
       t@@ -54,4 +54,7 @@ double grain_rotational_kinetic_energy(const struct grain *g);
        double grain_kinetic_energy(const struct grain *g);
        void grain_temporal_integration(struct grain *g, double dt);
        
       +void grain_register_contact(struct grain *g, size_t i, size_t j,
       +                            double centerdist[3], double overlap);
       +
        #endif
 (DIR) diff --git a/simulation.c b/simulation.c
       t@@ -106,24 +106,22 @@ static void
        sim_check_add_contact(struct simulation *sim, size_t i, size_t j)
        {
        
       -        double overlap;
       -        int ic;
       +        double overlap, centerdist[3], centerdistinv[3];
       +        int d;
        
       -        for (ic = 0; ic <= MAXCONTACTS; ic++) {
       -                continue; /* TODO */
       -        }
       +        if (i >= j)
       +                return;
        
       -        overlap = sim_grain_pair_overlap(sim, i, j);
       -        if (overlap < 0.0) {
       -                for (ic = 0; ic <= MAXCONTACTS; ic++) {
       +        for (d = 0; d < 3; d++) {
       +                centerdist[d] = sim->grains[i].pos[d] - sim->grains[j].pos[d];
       +                centerdistinv[d] = sim->grains[j].pos[d] - sim->grains[i].pos[d];
       +        }
        
       -                        if (ic == MAXCONTACTS)
       -                                errx(1, "MAXCONTACTS (%d) exceeded (grains %zu, %zu)", 
       -                                         MAXCONTACTS, i, j);
       +        overlap = euclidean_norm(centerdist, 3)
       +                        - (sim->grains[i].radius + sim->grains[j].radius);
        
       -                                /* TODO */
       -                }
       -        }
       +        grain_register_contact(&sim->grains[i], j, i, centerdist, overlap);
       +        grain_register_contact(&sim->grains[j], i, j, centerdistinv, overlap);
        }
        
        static void