hstop.cpp File Reference

#include "SilbertPeriodic.h"

Classes
class	SilbertHstop

Functions
bool	PointIsAboveCurve (Mdouble h, Mdouble a, int study_num)
void	HstopCurve (int study_num)
int	main (int argc, char *argv[])

Variables
fstream	ReportFile

Function Documentation

HstopCurve()

void HstopCurve

(

int

study_num

)

{
    Mdouble hStart = 4; //height at which bisection algorithm starts
    Mdouble aStart = 21; //angle at which bisection algorithm starts
    
    Mdouble hMax = 60; //height at which bisection algorithm stops
    Mdouble dh = 2; //how much I decrease the height after a point is found
    Mdouble da = .5; // how exact I want the angle to be determined
    Mdouble daMax = 1; //how much (in absolute) I increase the angle if the chute still stops
 
    stringstream name;
    name << "Report" << study_num;
    ReportFile.open(name.str().c_str(),fstream::out);
    ReportFile << "height\tangle\tlambda\tmu\tabovemuBottom\tabove\n";
 
    //we start at (aStart,hStart) and find the largest non-flowing angle for hStart
    //first we search in large steps daMax
    
    //find a flowing angle aMax and stopping angle aMin
    Mdouble aMax = aStart;
    while(!PointIsAboveCurve(hStart,aMax,study_num)) aMax+=daMax;
    Mdouble aMin=aMax-daMax;
    if (aMax==aStart) while(PointIsAboveCurve(hStart,aMin,study_num)) {aMin-=daMax; aMax-=daMax;}
    //bisect until aMin is highest stopping angle
    while (aMax-aMin>da*1.99) {
        Mdouble a=(aMin+aMax)/2;
        if (PointIsAboveCurve(hStart,a,study_num)) aMax=a;
        else aMin=a;
    }
    //now increase height gradually; decrease angle if flow starts
    Mdouble a = aMin;
    for (Mdouble h=hStart+dh; h<hMax; h+=dh) {
        while (PointIsAboveCurve(h,a,study_num)) a-=da;
    }
    
    ReportFile.close();
    return; 
}

References a, plotDoE::name, out(), PointIsAboveCurve(), and ReportFile.

Referenced by main().

main()

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

{
    if (argc>1) HstopCurve(atoi(argv[1]));  
    else exit(-1);
}

References HstopCurve().

PointIsAboveCurve()

bool PointIsAboveCurve	(	Mdouble	h,
		Mdouble	a,
		int	study_num
	)

                                                            {
    SilbertHstop problem;
    problem.setInflowHeight(h) ;
    problem.setChuteAngle(a);
    problem.set_study(study_num);
    stringstream name;
    name << "H" << problem.getInflowHeight() 
         << "A" << problem.getChuteAngleDegrees() 
         << "L" << round(100.*problem.getFixedParticleRadius()*2.)/100.
         << "M" << problem.getSlidingFrictionCoefficient()
         << "B" << problem.getSlidingFrictionCoefficientBottom();
    problem.setName(name.str().c_str());
    problem.solve();
    stringstream com;
    ReportFile   << problem.getInflowHeight() 
         << "\t" << problem.getChuteAngleDegrees() 
         << "\t" << round(100.*problem.getFixedParticleRadius()*2.)/100.
         << "\t" << problem.getSlidingFrictionCoefficient()
         << "\t" << problem.getSlidingFrictionCoefficientBottom()
         << "\t" << problem.pointIsAboveCurve
         << endl;
    return problem.pointIsAboveCurve;
}

References a, plotDoE::name, problem, ReportFile, and Eigen::bfloat16_impl::round().

Referenced by HstopCurve().

Variable Documentation

ReportFile

fstream ReportFile

Referenced by HstopCurve(), and PointIsAboveCurve().