tadd acceleration lock setting - granular - granular dynamics simulation
 (HTM) git clone git://src.adamsgaard.dk/granular
 (DIR) Log
 (DIR) Files
 (DIR) Refs
 (DIR) README
 (DIR) LICENSE
       ---
 (DIR) commit 8617d6d7191d79d1fef9852751a40b5d40e9b936
 (DIR) parent b60ee39cf043b8660bf185e351259f9186c4c331
 (HTM) Author: Anders Damsgaard <anders@adamsgaard.dk>
       Date:   Wed, 21 Apr 2021 12:01:53 +0200
       
       add acceleration lock setting
       
       Diffstat:
         M grain.c                             |      49 +++++++++++++++++++++----------
         M grain.h                             |       1 +
         M grains.c                            |       1 +
         M granular.5                          |       6 ++++++
         M simulation.c                        |       2 +-
       
       5 files changed, 42 insertions(+), 17 deletions(-)
       ---
 (DIR) diff --git a/grain.c b/grain.c
       t@@ -21,7 +21,7 @@ grain_new(void)
        void
        grain_defaults(struct grain *g)
        {
       -        size_t i;
       +        int i;
        
                g->radius = 0.5;
                for (i = 0; i < 3; i++) {
       t@@ -37,7 +37,8 @@ grain_defaults(struct grain *g)
                        = g->forceext[i]
                        = g->contact_stress[i]
                        = 0.0;
       -                g->gridpos[i] = 0;
       +                g->gridpos[i]
       +                = g->acc_lock[i] = 0;
                }
                g->fixed = 0;
                g->rotating = 1;
       t@@ -58,16 +59,25 @@ grain_defaults(struct grain *g)
        static void
        print_padded_nd_double(FILE *stream, const double *arr)
        {
       -        size_t i;
       +        int i;
                
                for (i = 0; i < 3; i++)
                        fprintf(stream, "%.*g\t", FLOATPREC, arr[i]);
        }
        
        static void
       -print_padded_nd_int(FILE *stream, const size_t *arr)
       +print_padded_nd_int(FILE *stream, const int *arr)
        {
       -        size_t i;
       +        int i;
       +        
       +        for (i = 0; i < 3; i++)
       +                fprintf(stream, "%d\t", arr[i]);
       +}
       +
       +static void
       +print_padded_nd_uint(FILE *stream, const size_t *arr)
       +{
       +        int i;
                
                for (i = 0; i < 3; i++)
                        fprintf(stream, "%zu\t", arr[i]);
       t@@ -80,6 +90,7 @@ grain_print(FILE *stream, const struct grain *g)
                print_padded_nd_double(stream, g->pos);
                print_padded_nd_double(stream, g->vel);
                print_padded_nd_double(stream, g->acc);
       +        print_padded_nd_int(stream, g->acc_lock);
                print_padded_nd_double(stream, g->force);
                print_padded_nd_double(stream, g->angpos);
                print_padded_nd_double(stream, g->angvel);
       t@@ -97,7 +108,7 @@ grain_print(FILE *stream, const struct grain *g)
                fprintf(stream, "%.*g\t", FLOATPREC, g->tensile_strength);
                fprintf(stream, "%.*g\t", FLOATPREC, g->shear_strength);
                fprintf(stream, "%.*g\t", FLOATPREC, g->fracture_toughness);
       -        print_padded_nd_int(stream, g->gridpos);
       +        print_padded_nd_uint(stream, g->gridpos);
                fprintf(stream, "%zu\t", g->ncontacts);
                print_padded_nd_double(stream, g->contact_stress);
                fprintf(stream, "%.*g\t", FLOATPREC, g->thermal_energy);
       t@@ -117,6 +128,7 @@ grain_read(char *line)
                                 "%lg\t%lg\t%lg\t" /* pos */
                                 "%lg\t%lg\t%lg\t" /* vel */
                                 "%lg\t%lg\t%lg\t" /* acc */
       +                         "%d\t%d\t%d\t"    /* acc_lock */
                                 "%lg\t%lg\t%lg\t" /* force */
                                 "%lg\t%lg\t%lg\t" /* angpos */
                                 "%lg\t%lg\t%lg\t" /* angvel */
       t@@ -125,9 +137,9 @@ grain_read(char *line)
                                 "%lg\t%lg\t%lg\t" /* disp */
                                 "%lg\t%lg\t%lg\t" /* forceext */
                                 "%lg\t"           /* density */
       -                         "%d\t"           /* fixed */
       -                         "%d\t"           /* rotating */
       -                         "%d\t"           /* enabled */
       +                         "%d\t"            /* fixed */
       +                         "%d\t"            /* rotating */
       +                         "%d\t"            /* enabled */
                                 "%lg\t"           /* youngs_modulus */
                                 "%lg\t"           /* poissons_ratio */
                                 "%lg\t"           /* friction_coeff */
       t@@ -138,11 +150,12 @@ grain_read(char *line)
                                 "%zu\t"           /* ncontacts */
                                 "%lg\t%lg\t%lg\t" /* contact_stress */
                                 "%lg\t"           /* thermal_energy */
       -                         "%d\n",          /* color */
       +                         "%d\n",           /* color */
                                 &diameter,
                                 &g->pos[0], &g->pos[1], &g->pos[2], 
                                 &g->vel[0], &g->vel[1], &g->vel[2], 
                                 &g->acc[0], &g->acc[1], &g->acc[2], 
       +                         &g->acc_lock[0], &g->acc_lock[1], &g->acc_lock[2], 
                                 &g->force[0], &g->force[1], &g->force[2], 
                                 &g->angpos[0], &g->angpos[1], &g->angpos[2], 
                                 &g->angvel[0], &g->angvel[1], &g->angvel[2], 
       t@@ -164,7 +177,7 @@ grain_read(char *line)
                                 &g->ncontacts,
                                 &g->contact_stress[0], &g->contact_stress[1], &g->contact_stress[2], 
                                 &g->thermal_energy,
       -                         &g->color) != 50)
       +                         &g->color) != 53)
                        errx(1, "%s: could not read line: %s", __func__, line);
        
                g->radius = 0.5 * diameter;
       t@@ -178,7 +191,7 @@ grain_read(char *line)
        int
        grain_check_values(const struct grain *g)
        {
       -        size_t i;
       +        int i;
                int status = 0;
        
                check_float_positive("grain->radius", g->radius, &status);
       t@@ -187,6 +200,7 @@ grain_check_values(const struct grain *g)
                        check_float("grain->pos", g->pos[i], &status);
                        check_float("grain->vel", g->vel[i], &status);
                        check_float("grain->acc", g->acc[i], &status);
       +                check_int_bool("grain->acc_lock", g->acc_lock[i], &status);
                        check_float("grain->force", g->force[i], &status);
                        check_float("grain->angpos", g->angpos[i], &status);
                        check_float("grain->angvel", g->angvel[i], &status);
       t@@ -251,7 +265,7 @@ grain_moment_of_inertia(const struct grain *g)
        void
        grain_zero_kinematics(struct grain *g)
        {
       -        size_t i;
       +        int i;
        
                for (i = 0; i < 3; i++) {
                        g->vel[i] = 0.0;
       t@@ -286,13 +300,16 @@ grain_kinetic_energy(const struct grain *g)
        static void
        grain_temporal_integration_two_term_taylor(struct grain *g, double dt)
        {
       -        size_t i;
       +        int i;
                double dx, mass = grain_mass(g);
                double moment_of_inertia = grain_moment_of_inertia(g);
        
                for (i = 0; i < 3; i++) {
       -                g->acc[i] = g->force[i] / mass;
       -                g->angacc[i] = g->torque[i] / moment_of_inertia;
       +                if (!g->acc_lock[i])
       +                        g->acc[i] = g->force[i] / mass;
       +
       +                if (g->rotating)
       +                        g->angacc[i] = g->torque[i] / moment_of_inertia;
                }
        
                if (g->fixed)
 (DIR) diff --git a/grain.h b/grain.h
       t@@ -12,6 +12,7 @@ struct grain {
                double pos[3];
                double vel[3];
                double acc[3];
       +        int acc_lock[3];
                double force[3];
                double angpos[3];
                double angvel[3];
 (DIR) diff --git a/grains.c b/grains.c
       t@@ -70,6 +70,7 @@ grains_print_vtk(FILE *stream, const struct grain *grains, size_t ng)
        
                VTK_XML_VECTOR(vel, "Velocity [m/s]", "Float64", VTK_FLOAT_FMT);
                VTK_XML_VECTOR(acc, "Acceleration [m/s^2]", "Float64", VTK_FLOAT_FMT);
       +        VTK_XML_VECTOR(acc_lock, "Acceleration lock [-]", "Int64", "%d ");
                VTK_XML_VECTOR(force, "Force [N]", "Float64", VTK_FLOAT_FMT);
                VTK_XML_VECTOR(angpos, "Angular position [rad]", "Float64", VTK_FLOAT_FMT);
                VTK_XML_VECTOR(angvel, "Angular velocity [rad/s]", "Float64", VTK_FLOAT_FMT);
 (DIR) diff --git a/granular.5 b/granular.5
       t@@ -36,6 +36,12 @@ acceleration, y [m]
        .It
        acceleration, z [m]
        .It
       +acceleration lock, (0: false, 1: true) [-]
       +.It
       +acceleration lock, (0: false, 1: true) [-]
       +.It
       +acceleration lock, (0: false, 1: true) [-]
       +.It
        force, x [m]
        .It
        force, y [m]
 (DIR) diff --git a/simulation.c b/simulation.c
       t@@ -148,7 +148,7 @@ sim_resolve_interactions(struct simulation *sim)
                for (i = 0; i < sim->ng; i++)
                        for (ic = 0; ic < MAXCONTACTS; ic++)
                                if (sim->grains[i].contacts[ic].active)
       -                                continue;
       +                                grain_interact(i);
        }
        
        void