SingleParticleSegregation.h

Go to the documentation of this file.

// This file is part of the MercuryDPM project (https://www.mercurydpm.org).
// Copyright (c), The MercuryDPM Developers Team. All rights reserved.
// License: BSD 3-Clause License; see the LICENSE file in the root directory.
 
#ifndef MERCURYDPM_SINGLEPARTICLESEGREGATION_H
#define MERCURYDPM_SINGLEPARTICLESEGREGATION_H
 
#include <Boundaries/CubeInsertionBoundary.h>
#include "Chute.h"
#include "Species/LinearViscoelasticSlidingFrictionSpecies.h"
#include "Boundaries/PeriodicBoundary.h"
#include "CG/TimeAveragedCG.h"
using constants::pi;
 
/*
 * Simulate a single intruder particle in a chute flow
 *
 * The simulation is non-dimensionalised such that diameter and mass of the flowing particle and gravity equals unity.
 *
 * To create some default statistics, use MercuryCG:
 * ../MercuryCG/MercuryCG SingleParticleSegregation -coordinates XZ -n 100 -timeaverage -tmin 500 -z 0 22.5 -w 0.5
 */
 
class SingleParticleSegregation : public Chute
{
public:
    // Set simulation parameters
    SingleParticleSegregation() {
 
        // Chute properties
        setGravity(Vec3D(0,0,1));   // Apply unit gravity
        setChuteAngle(30.);                     // Inclination angle of the chute (degrees)
        setChuteLength(15.);                    // Length of chute in x-direction
        setChuteWidth(6.);                      // Width of chute in y-direction
        setInflowHeight(20.);                   // Height of inflow region
        setInflowParticleRadius(0.5);           // Radius of flowing particles
        setFixedParticleRadius(0.425);          // Radius of fixed bottom particles
        setRoughBottomType(MULTILAYER);         // type of bottom used (change to MULTILAYER for a good demo)
 
        // Material properties
        LinearViscoelasticSlidingFrictionSpecies species;
        species.setDensity(6. / constants::pi);             // Normalised density so that bulk particle has mass 1
        double collisionTime = 0.005;                       // Collision time scale and coefficient of restitution
        species.setCollisionTimeAndRestitutionCoefficient(collisionTime, 0.1, 1);
        species.setSlidingFrictionCoefficient(0.5); // Sliding friction coefficient
        species.setSlidingDissipation(2./7.*species.getDissipation()); //Tangential dissipation
        species.setSlidingStiffness(2./7.*species.getStiffness()); //Tangential stiffness
        speciesHandler.copyAndAddObject(species);
 
        // intruder properties
        intruderRadius = 0.85;                  // Diameter of the large intruder particle
        intruderHeight = 0.5*getInflowHeight(); // Initial height of intruder above chute bottom
        springStiffness = 20;                   // Spring stiffness (used in contact model)
 
        // time stepping properties
        setTimeStep(0.05*collisionTime); // Time step (should be 1/50 of collision time for publication!)
        setTimeMax(10);                // Maximum simulation time
        setSaveCount((unsigned)(0.2/getTimeStep())); // Store one snapshot every 0.2 seconds
 
        // set name
        setName("SingleParticleSegregation");
    }
 
    void setupInitialConditions() override
    {
        // remove existing putput files
        removeOldFiles();
 
        // Make chute periodic in y-direction (spanwise)
        makeChutePeriodic();
 
        // Add periodic boundary in x-direction (streamwise)
        PeriodicBoundary b0;
        b0.set(Vec3D(1.0, 0.0, 0.0), getXMin(), getXMax());
        boundaryHandler.copyAndAddObject(b0);
 
        Chute::setupInitialConditions();
        addIntruderParticle();
        addFlowParticlesCompactly();         // Fill the chute with particles at inflow
        setParticlesWriteVTK(true);          // Enable VTK output for visualisation
        boundaryHandler.removeLastObject();  // Remove last boundary (cleanup, optional)
    }
 
    void addIntruderParticle() {
        SphericalParticle p;
        p.setSpecies(speciesHandler.getObject(0));
        p.setRadius(intruderRadius);
        p.setPosition(Vec3D(0.5*getChuteLength(), 0.5*getChuteWidth(), intruderHeight));
        intruderParticle = particleHandler.copyAndAddObject(p);
    }
 
    // add spring force to the intruder particle to force it to remain near a given z-position.
    void computeAdditionalForces() override {
        intruderParticle->addForce({0,0,springStiffness*(intruderHeight-intruderParticle->getPosition().Z)});
    }
 
    // move all particles to keep intruder centred
    void actionsAfterTimeStep() override {
        Vec3D shift = Vec3D(0.5*getChuteLength()-intruderParticle->getPosition().X, 0.5*getChuteWidth()-intruderParticle->getPosition().Y,0);
        for (BaseParticle* p : particleHandler) {
            p->move(shift);
        }
        // record time and shift, so original positions can be reconstructed
        static std::ofstream shiftRecorder (getName()+".shift");
        static Vec3D totalShift {0,0,0};
        totalShift += shift;
        shiftRecorder << getTime() << " " << totalShift << std::endl;
    }
 
    void printTime() const override {
        logger(INFO, "t=%, N=%, ene=%", getTime(), particleHandler.getSize(), getKineticEnergy()/getElasticEnergy());
    }
 
private:
 
    // simulation parameters
    double intruderRadius;      // Diameter of the large intruder particle
    double intruderHeight;      // Initial height of intruder above chute bottom
    double springStiffness;     // Spring stiffness (used in contact model)
 
    // pointers
    SphericalParticle* intruderParticle;
};
 
 
#endif //MERCURYDPM_SINGLEPARTICLESEGREGATION_H