HertzianSinterForceUnitTest.cpp File Reference

#include "DPMBase.h"
#include <Species/HertzianSinterSpecies.h>
#include <Logger.h>

Classes
class	HertzianSinterForceUnitTest
	This code tests our plastic force model, as published in Luding 2008. More...
class	LongHertzianSinterForceUnitTest

Functions
int	main (int argc UNUSED, char *argv[] UNUSED)

Function Documentation

main()

int main	(	int argc	UNUSED,
		char *argv[]	UNUSED
	)

{   
    HertzianSinterForceUnitTest sf;
    double k1=100.0;
    sf.species->setPlasticParameters(k1, 10.0*k1, k1, 1.0);
    sf.species->setDissipation(0);
    sf.setParticleDimensions(3);
    sf.setSystemDimensions(3);
 
    logger(INFO,"Testing particle particles collision for elastic plastic forces. \n"
        "This will be done for several values of scaled relative velocity chi");
    
    //sf.restartFile.getFstream().precision(20);
    
    //Set up constant data that will be used
    const std::vector<double> chi = {0.34, 0.69, 1.1, 1.37};
    const std::vector<Vec3D> leftFinalVecloity = {
        Vec3D(-0.032721738352012,0.0,0.0),
        Vec3D(-0.0138683231953154,0.0,0.0),
        Vec3D(-0.0204655358555405,0.0,0.0),
        Vec3D(-0.163049415300304,0.0,0.0)};
    const std::vector<Vec3D> leftFinalPosition = {
        Vec3D(0.995546292935715,1.0,1.0),
        Vec3D(1.00695193269955,1.0,1.0),
        Vec3D(1.00840467123501,1.0,1.0),
        Vec3D(0.969386085767181,1.0,1.0)};
 
    //Loop over all test cases
    for (int i=0; i<4; i++)
    {
        logger(INFO, "Running for chi=%", chi[i]);
        sf.set_chi(chi[i]);
        std::stringstream ss("");
        ss << "HertzianSinterForceUnitTest" << sf.get_chi();
        sf.setName(ss.str().c_str());
        sf.solve();
        sf.writeRestartFile();
 
//        //Now check the particles are in the right place for each of the 4 cases
//        auto pIt = sf.particleHandler.begin();
//        if (!(*pIt)->getPosition().isEqualTo(leftFinalPosition[i], 1e-10))
//                logger(FATAL,"Left particle is in the wrong position. It is at % and should be %",(*pIt)->getPosition(),leftFinalPosition[i]);
//        if (!(*pIt)->getVelocity().isEqualTo(leftFinalVecloity[i]  , 1e-10))
//                logger(FATAL,"Left particle has the wrong velocity. It is at % and should be %",(*pIt)->getVelocity(),leftFinalVecloity[i]);
    }
 
    //A longer simulation, with sintering activated
    LongHertzianSinterForceUnitTest lsf;
    lsf.species->setPlasticParameters(k1, 10.0*k1, k1, 1.0);
    lsf.species->setSinterRate(1e-3);
    lsf.species->setDissipation(2e-1);
    lsf.setParticleDimensions(3);
    lsf.setSystemDimensions(3);
    logger(INFO, "Running longer");
    lsf.set_chi(0.47);
    lsf.setName("LongHertzianSinterForceUnitTest");
    lsf.solve();
    lsf.writeRestartFile();
    
    
    
    //(1-((k1+kc)/k2)^1.5) / (1-(k1/k2)^1.5)   
    
    logger(INFO,"Execute 'gnuplot HertzianSinterForceUnitTest.gnu' to view "
                       "output");
    helpers::writeToFile("HertzianSinterForceUnitTest.gnu",
                         "set xlabel 'displacement'\n"
                         "set ylabel 'force'\n"
                         "plot 'HertzianSinterForceUnitTest1.1.fstat' u 7:9, 4./3.*100*x**1.5, 4./3.*1000*(x-0.1)**1.5-4./3.*100*x**1.5, -4./3.*100*x**1.5/(x<0.1)\n"
                         );
    
    logger(INFO,"Execute 'gnuplot LongHertzianSinterForceUnitTest.gnu' to view "
                  "output");
    helpers::writeToFile("LongHertzianSinterForceUnitTest.gnu",
                         "set xlabel 'time'\n"
                         "set ylabel 'displacement'\n"
                         //"plot 'LongHertzianSinterForceUnitTest1.1.fstat' u 7:9 w lp\n"
                         "plot 'LongHertzianSinterForceUnitTest.fstat' u 1:7, sqrt(0.001*(x+120))\n"
                         );
}

References HertzianSinterForceUnitTest::get_chi(), constants::i, INFO, logger, HertzianSinterForceUnitTest::set_chi(), DPMBase::setName(), DPMBase::setParticleDimensions(), DPMBase::setSystemDimensions(), DPMBase::solve(), HertzianSinterForceUnitTest::species, DPMBase::writeRestartFile(), and helpers::writeToFile().