#include <stdio.h>
#include <stdlib.h>
#include <math.h>
/*
Программа проверяет теорему Ньютона о том, что сферически-симметричная оболочка не создает сил тяготения во внутренней полости
Используется метод "Монте-карло" - на каждой итерации выбирается случайная точка внутри слоя оболочки
и вычисляется сила притяжения к ней пробного тела.
Константа gMm для простоты не вычисляется, т.к. на искомое результирующее направление вектора не влияет.
Результирующая вычисляется как приращение к ней проекции на ось х вектора силы притяжения текущей итерации,
т.е. имеющего направление на текущую случайную точку и размер вычисленной силы.
Приглашаю всех заинтересовавшихся к участию!
*/
int main()
{
srand(time(0));
const int R = 6371;
int R1 = R/2; //внешний радиус кольца
int R2 = R/2 - R/10; //внутренний радиус кольца
double Ox = 0; //Ox,Oy - координаты центра оболочки
double Oy = 0; //
double Cx = Ox+R2-R2/10; //Cx,Cy - координаты тестовой точки, отстоящей вправо от центра, на 1/10-ю радиуса полости от внутреннего края.
double Cy = 0;
double F_res=0, F_ca, F, a, Ax, Ay, cosa;
int i, n, k=0;
n = 10000;
for(i = 0 ; i < n; i++)
{
Ax = rand()%((R1+2)*2) - (R1+2); //Ax, Ay - координаты случайной точки
Ay = rand()%((R1+2)*2) - (R1+2);
if ( ((sqrt(Ax*Ax+Ay*Ay))<=R1)&&((sqrt(Ax*Ax+Ay*Ay))>=R2) ) //если точка с данными координатами находится в "толще", то
{
// вычисляем силу тяжести от этой точки к пробной:
a = sqrt( (Cx-Ax)*(Cx-Ax) + (Cy-Ay)*(Cy-Ay) ); //расстояние СА
F_ca = 1/(a*a); // Сила воздействия А на С
//вычисляем угол между текущим результирующим вектором тяготения (лежит на х) и текущим СА
cosa = (Cx*Ax)/( (sqrt(Cx*Cx))*(sqrt(Ax*Ax + Ay*Ay)) ); //
// вычисляем результирующий вектор тяготения к А:
F = cosa/(a*a);
// приращаем вектор от А к общему:
F_res = F_res+F;
}
//system("cls");
cout << F_res*1000;
printf(" ");
}
return 0;
}