дело такое. написал код для расчета интеграла перекрывания орбиталей
Код: Выбрать все
#include "stdafx.h"
#include <cmath>
#include <windows.h>
#include <iostream>
using namespace std;
double f(double x, double y, double z,int m)
{ //Подынтегральная функция
double a = 0.5292;
double phi = 2*3.14/(3*m);
double dy = 1.42;
double dx = 0.71;
return (z*z*cos(phi/2)*cos(phi/2)-z*dy*sin(phi)+(dy*dy-y*y)*sin(phi/2)*sin(phi/2))*exp(-1/(2*a)*(sqrt((x+dx)*(x+dx)+(y+dy)*(y+dy)+z*z)+sqrt((x-dx)*(x-dx)+(y-dy)*(y-dy)+z*z)));
}
double rectangle_integrate(double a, double b, int n, double (*f)(double,double,double,int) )
{
long double result, h;
int i,j,k,m;
cin >>m;
h = (b-a)/n; //Шаг сетки
result = 0.0;
for(i=1; i <= n; i++)
{
for (j=1; j<=n; j++)
{
for (k=1; k<=n; k++)
{
result += f( a + h * (i - 0.5), a + h * (j - 0.5),a + h * (k - 0.5), m ); //Вычисляем в средней точке и добавляем в сумму
}
}
}
result *= h*h*h/(32*3.14*pow(a,5));
return result;
}
int _tmain(int argc, _TCHAR* argv[])
{
double integral;
integral=rectangle_integrate(-1,1,100,f);
printf("The value of the integral is: %lf \n", integral);
system("pause");
return 0;
}
проблема вот в чем.
научник сказал что пределы подбирать самостоятельно такие, чтобы экспонента с отрицательным показателем убивалась.
но вот непонятно. характерные размеры это расстояние между соседними атомами графита. то есть 1,42 ангстрем. в рассчетах также используется первый боровский радиус.
я вот пытаюсь считать все в ангстремах, а результат получается отрицательным...это же не хорошо. в коде вроде нигде не напортачил
посоветуйте что нибудь..
спасибо)
