1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194
|
#include "stdafx.h"
#include <iostream>
#include <iomanip>
#include <ctime>
#include <fstream>
#include <vector>
#include <math.h>
#include <stdio.h>
#include <cstdlib>
#include <stdlib.h>
#include "boost/random.hpp"
#include "boost/generator_iterator.hpp"
#include "boost/random/normal_distribution.hpp"
#include <boost/random/uniform_int.hpp>
#include <boost/random/uniform_01.hpp>
using namespace std;
using std::vector;
unsigned histogram[650]={0};
double sums[650]={0};
std::vector<double> frequency_bins[650];
class Gravitar //creates a class called Gravitar
{
public: //public assesors
long double GetRadius();
void SetRadius(long double radius);
long double GetAngle();
void SetAngle(long double angle);
long double GetZed();
void SetZed(long double zed);
long double GetAge();
void SetAge(long double age);
long double GetPeriod();
void SetPeriod(long double period);
private: //private values
long double Radius;
long double Angle;
long double Zed;
long double Age;
long double Period;
};
long double Gravitar::GetRadius() //GetRadius, public assessor function returns value of Radius
{
return Radius; //returns the value Radius
}
long double Gravitar::GetAngle()
{
return Angle; //returns the value of Angle
}
long double Gravitar::GetZed()
{
return Zed;
}
long double Gravitar::GetAge()
{
return Age;
}
long double Gravitar::GetPeriod()
{
return Period;
}
void Gravitar::SetRadius(long double radius)//sets radius to Radius
{
Radius = radius;
}
void Gravitar::SetAngle(long double angle)
{
Angle = angle;
}
void Gravitar::SetZed(long double zed)
{
Zed = zed;
}
void Gravitar::SetAge(long double age)
{
Age = age;
}
void Gravitar::SetPeriod(long double period)
{
Period = period;
}
int main()
{
int m ;
cout << "How many Gravitars?\n\n";
cin >> m;
long double e ;
cout << "\n\nWhat is the ellipticity of the objects? Enter 6, 7, 8, 9, for gravitars of ellipticity of 10^-x.\n\n";
cin >> e;
long double Ellipse = pow (10, -e); //the ellipticity should vary from 10^-6 to 10^-7 to 10^-8 to 10^-9
boost::mt19937 rng( time(0) );
boost::normal_distribution<> nd(5, .69);
boost::variate_generator<boost::mt19937&, boost::normal_distribution<> > var_nor(rng, nd);
boost::minstd_rand intgen;
boost::uniform_01<boost::minstd_rand> gen(intgen);
long double R;
ofstream myfile
//("agedata.txt");
//("positiondata.txt");
("postspindowndata.txt");
//("postspindowndata.dat");
unsigned histogram[650]={0};
for (int R = 0; R < m; R++)
{
long double radius = sqrt(gen());
long double angle = gen();
long double zed = gen();
long double age = gen();
long double period = var_nor()*.001; // normal distribution of periods with a mean of 5 ms and a standard deviation of .69 mseconds
Gravitar RadialComponents;
RadialComponents.SetRadius(radius);
Gravitar AngularComponents;
AngularComponents.SetAngle(angle);
Gravitar ZedComponents;
ZedComponents.SetZed(zed);
Gravitar AgeComponents;
AgeComponents.SetAge(age);
Gravitar PeriodComponents;
PeriodComponents.SetPeriod(period);
long double x = RadialComponents.GetRadius()*cos (AngularComponents.GetAngle()*6.283185307179586476925286766559)*15329.729;
long double y = RadialComponents.GetRadius()*sin (AngularComponents.GetAngle()*6.283185307179586476925286766559)*15329.729;
long double z = ZedComponents.GetZed()*153.29729;
long double a = AgeComponents.GetAge() * 10000000.0 * 31556926.0; // Random number generated between 0 and 1, then multiplied by 10 millions years, multiplied by the number of seconds within a year
long double p = PeriodComponents.GetPeriod(); //p measured in seconds
long double prefrequency = 2.0 / p;
long double xde = x;
long double deltay = (y - 8500);
long double deltay2 = deltay*deltay;
long double yde = sqrt(deltay2);
long double zde = z;
long double xde2 = xde*xde;
long double yde2 = yde*yde;
long double zde2 = zde*zde;
long double dge = sqrt(xde2 + yde2 + zde2);
//WARNING:E is the NOT mathematical constant e.
//It is the const var of ellipticity
long double Ellipse2=Ellipse*Ellipse;
const long double pi=3.1415926535897932384626433832795,
pi4=pi*pi*pi*pi,
pi2=pi*pi,
G=6.673E-11;
long double beta= (32.0 * pi4 * G * 10E+44) / (2.421E+42);
int value = -4;
long double othervalue = -.25;
long double postfrequency = pow( (pow(prefrequency, value) + beta * Ellipse2 * a), othervalue);
long double h = ((4.0 * pi2 * G * 10E+44) / (8.077E+33)) * ((Ellipse * pow( postfrequency, 2)) / dge);
long double hx = h*1000000000;
long double hxr = floor (hx);
double hxrf = hxr/1000000000;
double f = hxrf;
size_t index=size_t((f<1E-10)?0:log10(f)+10);
sums[index]+=f;
frequency_bins[index].push_back(f);
//double f = hxrf;
//size_t index=size_t((f<1E-6)?0:log10(f)+6);
//histogram[index]++;
//std::cout <<"histogram["<<index<<"]="<<histogram[index]<<std::endl;
//printf("Gravitar # %4f %16f %16f %16f %16f %16f %32f %16f %16f %16f %16f \n" , R+1 , x , y , z , a , p , Ellipse , prefrequency , postfrequency, dge, h, hxrf );
}
for (int a=0;a<650;a++){
std::cout <<"Bucket "<<a<<":\n"
"Sum of h values for bucket: "<<sums[a]<<"\n"
"Frequencies in bucket:\n";
for (size_t b=0;b<frequency_bins[a].size();b++)
std::cout <<frequency_bins[a][b]<<" ";
std::cout <<std::endl;
}
myfile.close();
system("PAUSE");
return 0;
}
|