Важная информация
Показано с 1 по 2 из 2

Тема: вычислить параметрический несобственный интеграл

  1. #1 вычислить параметрический несобственный интеграл 
    Новичок
    Регистрация
    04.11.2013
    Сообщений
    3
    Сказал(а) спасибо
    0
    Поблагодарили 0 раз(а) в 0 сообщениях
    Добрый день. Я на форуме dxdy в теме [Ссылки могут видеть только зарегистрированные пользователи. ] спрашивал как со стороны математики вычислить интеграл:
    [math]$$\int_{x_0}^{\infty} (e^{-4x}(32/x+32/x^2+16/x^3+4/x^4)/\sqrt{1-l^2/a^2x^2+2(1+1/x)e^{-2x}}) dx $$[/math]
    Разобрался, преобразовал интеграл, разбил на два, но к сожалению, созданная программа не верно выдаёт ответ. А именно, то что выдаётся в ответе:
    [Ссылки могут видеть только зарегистрированные пользователи. ]

    то что должно быть:
    [Ссылки могут видеть только зарегистрированные пользователи. ]

    Если требуется, то могу подробно описать условие физической задачи и метод преобразования интеграла.

    Вот код программы:
    Код :
    double f (double x, double l);//подкоренное выражение
    double f2 (double t, double c, double l);//преобразованное подкоренное выражение, разложенное в ряд
    double Fun1 (double c, double A, double l);//Метод Симпсона, когда нижний предел  особая точка
    double Fun2 (double A, double B, double l);//Метод Симпсона, когда верхний предел бесконечность
    double F1 (double x, double c, double l);//преобразованная функция
    double F2 (double x, double l);//обычная функция
    double lowLim (double a, double b, double l);//Нахождение нижнего предела интеграла
     
    int main(void)
    {	
    	const short n = 10;
    	double x[n], Ed[n], c, A, B, e0;
    	double a, h = 0.5, x0 = 0.5; 
    	double C = 3 * pow(10.0,8), H = 1.0546 * pow(10.0,-34), e = 1.6 * pow(10.0,-19), m = 9.11 * pow(10.0, -31);
    	//a = pow(H, 2)/ (m * pow(e, 2));
     
    	a=0.529*pow(10.0,-11);
    	//e0 = pow(e,8)/(3 * a * pow(H,3) * pow(C,3));
    	e0 = pow(e,5)/(3 * a * pow(H,3) * pow(C,3));
    	short i;
     
    	// для первой точки
    	x[0]=x0;
    	c = lowLim (0.133758, 0.133760, pow((x[0]), 2)); 
     
    	A = c + 0.05;  
    	B = A + 1; //первое приближение верхнего предела
    	Ed[0] = Fun1 (c, A, pow(x[0], 2)) + Fun2 (A, B, pow(x[0], 2));
     
    	//для последующих точек
    	for (i = 1; i<n; i++)
    	{
    		x[i] = x[i-1] + h; 
    		c = lowLim (0.64, 5.01, pow((x[i]), 2));
     
    		A = c + 0.1;
    		B = A + 1;
    		Ed[i] = Fun1 (c, A, pow(x[0], 2)) + Fun2 (A, B, pow(x[i], 2));
    	}
     
    	ofstream out;
        out.open ("D:\\Задача 151.txt");
    	for (i = 0; i<n; i++)
    		//out << x[i] <<"\t" << Ed[i] <<"\n";
    		out << x[i] <<"\t" << log(Ed[i]/e0) <<"\n";
    	out.close();
    	return 0;
    }
     
    double lowLim (double a, double b, double l)
    {
    	double с, E = 0.00001;
    	do
    	{
    		с = (b + a)/ 2;
    		if (f(a, l) * f(с, l) > 0) a = с;
    		else b = с;
    	}
    	while (fabs(b-a) > E);
    	return (a+b)/2;
    }
     
    double f(double x,double l)
    {
    	return 1-l/pow(x,2)+2*(1+1/x)*exp(-2*x);
    }
     
    double Fun1 (double c, double A, double l)
    {
    	short i, n = 20;
    	double x, I, h;
    	h = (A - c)/ n; //шаг
     
        I = h * (F1(c, c, l) + F1(c + n * h, c, l))/ 3; 
    	for (i = 1; i<n-1; i++)//n-1!
    	{
    	    x = c + i * h; 
    		if (i%2 == 0) I += h * 2 * F1 (x, c, l)/ 3; 
    		if (i%2 == 1) I += h * 4 * F1 (x, c, l)/ 3;
    	}
    	return I;
     
    }
     
    double F1 (double t, double c, double l)
    {
    	return (2*exp((-4)*(pow(t,2)+c))*(16+32*pow((pow(t,2)+c),-1)+32*pow((pow(t,2)+c),-2)+16*pow((pow(t,2)+c),-3)+4*pow((pow(t,2)+c),-4)))/pow(f2(t, c, l), 0.5);
     
     
    }
     
    double f2 (double t, double c, double l)
    {
    	return (2*(l-c*(2*pow(c,2)+2*c+1)*exp(-2*c)))/(pow(c,3))+((2*c*(2*pow(c,3)+2*pow(c,2)+2*c+1)*exp(-2*c)-3*l)*pow(t,2))/(pow(c,4))+(2*(6*l-c*(4*pow(c,4)+4*pow(c,3)+6*pow(c,2)+6*c+3)*exp(-2*c)*pow(t,4)))/(3*pow(c,5))+((exp(-2*c)*(4*pow(c,6)+4*pow(c,5)+8*pow(c,4)+12*pow(c,3)+12*(c,2)-15*exp(2*c)*l+6*c)*pow(t,6))/(3*pow(c,6)));
     
    }
     
    double Fun2 (double A, double B, double l)
    {
     
    	double Ed1, Ed2, /*интегралы-приближения к интегралу правого предела*/ epsilon = 0.0001;
    	short i, n = 20;
    	double x, h;
     
    	Ed1 = 0;
    	do
    	{
    		Ed2 = Ed1;
    		h = (B - A)/ n;
     
    		Ed1 = h * (F2(A, l) + F2(A + n * h, l))/ 3;
    		for (i = 1; i<n; i++)
    		{
    			x = A + i * h; 
    			if (i%2 == 0) Ed1 += h * 2 * F2 (x, l)/ 3;
    			if (i%2 == 1) Ed1 += h * 4 * F2 (x, l)/ 3;
    		}
     
    		n += 10;
    		B += 10*h;
     
    	} while (fabs(Ed2 - Ed1) > epsilon);
     
    	return Ed1;
    }
     
    double F2(double x, double l)
    {
    	return (exp(-4*x)*(32/x+32/(pow(x,2))+16/(pow(x,3))+4/(pow(x,4))))/(pow((1-l*(1/(pow(x,2)))+2*(1+1/x)*exp(-2*x)),0.5));
    }

    Уже около двух недель решаю данное задание, и был-бы очень благодарен если-бы кто то нашел ошибку!
    Ответить с цитированием  
     

  2. #2  
    Профи Аватар для Dimon012
    Регистрация
    09.02.2011
    Адрес
    Владивосток
    Сообщений
    850
    Сказал(а) спасибо
    29
    Поблагодарили 130 раз(а) в 108 сообщениях
    Записей в блоге
    3
    Я не специалист в С, но для поиска ошибок в данном случае могу предложить следующее...
    Сейчас Вы контролируете исходные данные и как я понял некоторые окончательные результаты, которые приведены в таблице выше. Нужно попытаться контролировать более узкие промежуточные области, например у Вас в программе есть цикл, поставьте в нем точку останова и проверьте значения всех переменных в каждом цикле, если программа линейная, нужно поставить несколько точек останова и проверить значения на различных этапах расчета. Я так надеюсь в Вашем случае возможно даже не надо вручную считать расчет т.к. получаемые значения гораздо больше чем ожидаемые, найдите переменную, которая имеет гораздо большее значение чем все остальные и скорее всего в ней и будет содержаться ошибка.

    Удачи!
    Ответить с цитированием  
     

Информация о теме
Пользователи, просматривающие эту тему

Эту тему просматривают: 1 (пользователей: 0 , гостей: 1)

Похожие темы

  1. Нужно вычислить кол-во дней между датами
    от Рафиль в разделе Turbo Pascal
    Ответов: 1
    Последнее сообщение: 19.01.2013, 19:08
  2. Составить массив и вычислить его сумму
    от Кузя21 в разделе Turbo Pascal
    Ответов: 0
    Последнее сообщение: 13.01.2013, 10:41
  3. Вычислить значение функции
    от оля в разделе Turbo Pascal
    Ответов: 9
    Последнее сообщение: 15.12.2011, 23:05
  4. Вычислить значение функции
    от Макс в разделе QBasic
    Ответов: 1
    Последнее сообщение: 06.12.2011, 04:08
  5. Вычислить уравнение
    от Анастасия в разделе QBasic
    Ответов: 1
    Последнее сообщение: 01.12.2011, 22:01
Ваши права
  • Вы не можете создавать новые темы
  • Вы не можете отвечать в темах
  • Вы не можете прикреплять вложения
  • Вы не можете редактировать свои сообщения
  •