//Circular Pipe Hydraulics
//11/27/2007
//BS
function theta(depth,dia)
{	//half of angle at the center
	var ans=0;
	var r=dia/2;
	if(depth<r)
	{ans=Math.acos(1-depth/r);}
	else
	{ans=Math.PI-Math.acos(depth/r-1);}
	return ans;	
}

function wetted_perimeter(depth,dia)
{
	//wetted perimeter
	var ans=0;
	ans=2*theta(depth,dia)*dia/2;
	return ans;
}


function flow_area(depth,dia)
{
	var ans=0;
	var alpha=theta(depth,dia);
	ans=0.5*((dia*dia/4 *(2*alpha-Math.sin(2*alpha))));
	return ans;
}

function hydraulic_radius(depth,dia)
{
	//Hydraulic Radius
	var ans=0;
	ans=flow_area(depth,dia)/wetted_perimeter(depth,dia);
	return ans;
}

function flow_velocity(depth,dia,slope,n)
{
	//velocity in cfs unit
	//depth and dia in ft
	var ans= 1.486/n*Math.pow(hydraulic_radius(depth,dia),(2/3))*Math.sqrt(slope);
	return ans;
}

function flow_rate(depth,dia,slope,n)
{
	//flow_rate in cfs
	//depth and dia in ft
	var ans=flow_velocity(depth,dia,slope,n)*flow_area(depth,dia);
	return ans;
}

function top_width(depth,dia)
{
	//width of the top free surface
	var ans=0;
	var y=depth;
	var d=dia;
	var r=dia/2;
	if(y<r)
	{var b=r-y;}
	else
	{var b=y-r;}
	ans=2*Math.sqrt(r*r-b*b);
	return ans;
}

function Fr(depth,dia,slope,n)
{ //Froude Number

	var ans=0;
	var g=32.1740;
	ans=flow_velocity(depth,dia,slope,n)/Math.sqrt(g*flow_area(depth,dia)/top_width(depth,dia));
	return ans;
}

function Fr_test(depth_test,dia,knownQ)
{
	//Not exactly for froude calculation
	//Froude1 = (knownQ ^ 2 * topWid(diaFt, yTest) / flowArea(diaFt, yTest) ^ 3 / 32.2) ^ (0.5)	//from VBA
	return Math.pow((knownQ*knownQ*top_width(depth_test,dia)/Math.pow(flow_area(depth_test,dia),3)/32.2),0.5);
}
function normal_depth(q,dia,slope,n)
{	//NOT PERFECT YET
	var max_iterations=1000;
	var counter=0;
	var msg="Message from normal_depth() function :"+'\n';

	var y_max=dia;
	var y_min=0;
	var tolerance=0.00001;
	var error_q=0;
	var y=(y_max+y_min)/2;

	do
	{	counter=counter+1;
		if (counter>max_iterations)
		{
			msg=msg+"No solution found after "+counter+" iterations!";
			alert(msg);
			 return;
		}

		error_q=q-flow_rate(y,dia,slope,n);
		if(error_q>0)
		{	y_min=y;}
		else
		{	y_max=y;}
		y=(y_max+y_min)/2;

		//var msg=msg+(y*12).toFixed(6)+", "+flow_rate(y,dia,slope,n)+", "+error_q+", "+Fr(y,dia,slope,n)+'\n';
		

	}
	while(Math.abs(error_q)>tolerance)
	//document.frm2.textarea1.value=msg;
	return y;
}

	function myCriticalDepth(D,q)
	{
		//RETURNS CRITICAL DEPTH IN FEET
		//D=DIAMETER IN FEET
		//q=FLOWRATE IN CFS
		var counter=0; var xmin=0; var xmax=2*Math.PI; var max_iterations=1000;
		var msg='results : <hr />';

		do
		{	counter++;
			if(counter>max_iterations)
			{
				var errmsg='Critical Depth\nNo solution found\nafter'+counter+' iterations';
				alert(errmsg); return;
			}

			x=(xmax+xmin)/2;		//angle substended at the center
			var A=D*D/8*(x-Math.sin(x));	//flow area
			var y=D/2*(1-Math.cos(x/2));	//flow depth
			var T=Math.sqrt(D*y-y*y)*2;	//top free surface width
			var fx=A*A*A/T-q*q/32.2; 	//fx=0 for critical depth to occur
			if(fx<0){xmin=x;}
			else {xmax=x;}
			msg=msg+'<br/>'+counter+' ======== '+x+' ======== '+fx;
		}while(Math.abs(fx)>0.000001)
		//document.write(msg);
		return D*(1-Math.cos(x/2))/2; 		//critical depth		
	}
/*
function paintWaterLevel(pipedia,depth,scale,canvasname,cx,cy)
{	//NOT USED; IT HAS BEEN MOVED TO PIPE_DRAW.JS 1:18 AM 6/17/2008
	//Paints water flow section in a canvas
	var canvas = document.getElementById(canvasname);
	 var ctx = canvas.getContext("2d");

	 ctx.fillStyle = "rgb(200,0,0)";
	// ctx.fillRect (10, 10, 55, 50);

	 ctx.fillStyle = "rgba(0, 0, 200, 0.5)";
	//One filled arc
	ctx.beginPath();
	var D=pipedia*scale; //pipe dia
	var y=depth*scale; //depth
	var A=Math.asin(2*(D/2-y)/D);	//start angle
	var B=Math.PI-A;		//end angle

	ctx.arc(cx,cy,D,A,B,0);
	ctx.fill();

	var g=canvas.getContext("2d");
	g.beginPath();
	g.fillStyle = "rgba(200,0, 0, 0.5)";
	g.arc(cx,cy,D,0,2*Math.PI,0);

	//draw vertical line
	g.moveTo(cx+D,cy+D);	
	g.lineTo(cx+D+50,cy+D);	//horz tick mark
	g.lineTo(cx+D+50,(cy+D-2*y));	//vertical bar
	g.lineTo(cx+D+50-10,(cy+D-2*y));	//horz tick mark
	g.stroke();

}
*/