#define _USE_MATH_DEFINES
#include <cmath>
#include <iostream>
#include <math.h>
#include <fstream>
using namespace std;
const double A = 0.7; //константы, даны по условию задачи
const double B = 0.3;
const double bm_step = 0.001; //шаг
double bm;
double F(double y) //подынтегральная функция
{
return sqrt(fabs(A*exp(-(y*y)) + B*exp(-(y/2)) - bm));
}


//
//
//Вычисление интеграла МЕТОДОМ ТРАПЕЦИЙ
//
//
//



double trap(double a, double b, int N)
{
double h = (b - a)/N;
double tmpA, tmpB;
double s;

tmpA = a;
tmpB = a + h;
s = 0; //обнуляем площадь, для того что бы потом просуммировать

for (int i = 0; i < N; i++) //цикл
{
s = s + (F(tmpA) + F(tmpB))*h/2; //формула для нахождения площади
//ограниченную сверху кривой

tmpA = tmpB;
tmpB = tmpB + h;
}

return s;
}



//Функция для нахождения верхнего предела интегрирования
//
//
//МЕТОД ДИХОТОМИИ
//
//
//

double f (double Ym, double bm) //функция f(y)
{
return (0.7*exp(-(Ym*Ym))+0.3*exp(-(Ym/2))-bm);
}
double Ym(double a, double b, double eps, double bm)
{
double c;
do
{
c=(a+b)/2; //объявляем шаг вычисления
if (f(c,bm)*f(a,bm)<0) b=c; //вызов функции
else a=c; //если f(a)*f(c) > 0
//то к значению a присваиваем значение c
}
while (fabs(a-b)>=eps); //«fabs» - берём модуль разности f(b)-f(a)
//и сравниваем с заданной точностью
return c;
}


int main() //заголовок главной функции
{
ofstream file; //открытие файла, file переменная описанная как ofstream
double V,a,b,eps,epsp; //описание вещественных переменных
int N, N1, m;
bm = 0.05;
cout <<"a="; cin >>a; //задаем левую границу интервала
cout <<"b="; cin >>b; //задаем правую границу интервала
cout <<"eps="; cin >>eps; //задаем точность вычислений
cout <<"epsp="; cin >>epsp;
cout<<"N="; //вывод строки символа N= с помощью cout
cin>>N; //ввод N с помощью cin
cout<<endl;
cout<<"Ym="<<Ym(a,b,eps,0.05)<<endl;
file.open("results.txt"); //открывает файл,
//results –полное имя файла на диске,
//txt- режим работы с открытым файлом


for (m = 0; m < 6; m++) //цикл
{
bm = 0.05; //начальное условие
while (bm <= 0.9) //пока bm меньше
// или равно 0.9 выполняется
{
N1=N+1;
while (fabs(trap(a,b,N)-trap(a,b,N1))>epsp)
{
N1=N1*2;
}


V = (m*M_PI + 3*M_PI/4)/trap(0, Ym(a,b,eps,bm), N1);
cout<<"m="<<m<<" "<<"bm="<<bm<<" "<<"trap="<<trap(0, Ym(a,b,eps,bm),N1)<<" "<< "V="<<V<<endl;

if (ceil(V) <= 40) //ceil -округляет дробное значение V
//в большую сторону
//и если меньше или равно 40, то выводит //соответствующие значения на экран
{
cout.setf(ios::left); //выводит значения m=
//использовав флаг форматирования //выравнивание по левой границе
cout<<"m="; //функции вывода результатов на экран
cout.width(10);
cout<<m<<"bm=";
cout.width(10); //устанавливает ширину поля вывода,
//10 позиций
cout<<bm<<"V=";
cout.width(10);
cout<<V<<endl;
cout<<Ym(a,b,eps,bm)<<"Ym=";
cout.width(10);
cout<<Ym(a,b,eps,bm)<<endl;
file.setf(ios::left); //записывает значения в файл
file<<"m=";
file.width(10);
file<<m<<"bm=";
file.width(10);
file<<bm<<"V=";
file.width(10);
file<<V<<endl;
}

bm += bm_step; //сложение с присваиванием
}
}

return 0;
}