// Arup Guha
// 3/6/04
// Solution to Vase problem from 1988 World Finals.
// Uses the disc method to determine the volume of a figure created by
// rotating a curve around the x-axis.

#include <iostream.h>
#include <fstream.h>
#include <stdio.h>
#include <math.h>

double vol(double x, double a, double b, double h, double d);
double invert_vol(double x);
void chartbisect(double area);

const double PI = 3.141592654;
double h, d, a, b;

void main() {

  ifstream fin("vase.in");

  // Read in h, d, a, b
  
  while (fin >> h >> d >> a >> b) {

    printf("h = %.3lf   d = %.3lf   a = %.3lf   b = %.3lf\n\n", h, d, a,b);

    // The area of the circle at the top of the vase. It uses the radius
    // as r(h).
    double area = PI*d*d/4*(1+b);

    // Prints out the chart.
    cout << "Rain\tLine" << endl;
    chartbisect(area);
    cout << endl;
  }
}

// Determines the volume of water needed to fill the vase up to a height
// of x. This is the definite integral of pi*(r(x)^2) evaulated from 0
// to x. The last term in the integral is from subtracting out the value
// of the function when x=0.
double vol(double x) {

  double mult = d*d/4*PI; // Constant to pull out before integrating.

  // This is the definite integral.
  double integral = x - a*h/(2*PI)*cos(2*PI*x/h)+b/(2*h)*x*x+a*h/(2*PI);

  return mult*integral;
}

// Uses the bisection method to determine what x value corresponds to a
// particular volume of water in the vase. NOTE: This bisection method
// only works if the function you are attempting to invert is
// monotonically increasing. ALSO, this function will lead you to an
// infinite loop if you try to pass to it a volume that corresponds to an
// x value that is NOT in between 0 and h.

double invert_vol(double x) {

  // The actual answer must be in between 0 and h, according to the
  // the problem specification.
  double min = 0, max = h;
  double mid=(min+max)/2;
  double volume = vol(mid);

  // Check to see if mid is close to the correct x value we want.
  while (fabs(volume-x) > .0001) {

    // Improve your guess my adjusting either min or max.
    if (volume > x) 
      max = mid;
    else
      min = mid;
    mid = (min+max)/2;
    volume = vol(mid);
  }

  return mid; // Found the x value to return.
}

// Prints out the chart of values.
void chartbisect(double area) {

  double step = 0.1, level;
  double maxvol = vol(h);

  // Calculate the number of lines in the chart.
  int numlines = (int)(10*maxvol/area);

  // Print out each line by passing in the proper parameter to the
  // invert_vol function.
  for (int i=1; i<=numlines; i++) {
    level = invert_vol(step*area); 
    printf("%.1lf\t%.2lf\n",step,level);
    step += 0.1;
  }
}