This is a very simple simulation that is not meant to reproduce an exact fit to CMB anisotropy but to get some feel for whether CMB anisotropy occurs only with gravitational redshift as the only source of anisotropy which it does not try to do exactly. A static Einstein universe is an excellent model now that we know that the original reasons for abandoning it were not quite right (Eddington’s instability result). The six parameter Lambda-CDM model fits the anisotropy well which is the model to overthrow. So is it possible to do this? I believe that a static Einstein universe is an infinitely superior model (as it gets rid of the cosmological constant problem since quantum field theory on static Einstein universe is actually close to the measured value of around and the quantum gravity problem is not as difficult here since for example the Wheeler-deWitt equation can be solved explicitly). So here is a little simulation that shows that there is some interesting anisotropy in the simplest type of model for gravitational redshift effect. Distribute a uniform set of points over a round three-sphere and assign them masses by an exponential distribution. Choose a point on the three-sphere as a reference pont and draw a two-sphere around it and histogram the geodesics meeting the 2-sphere. Assume the background radiation has 2.7 K uniform background so that the anisotropy is due only to the gravitational redshift of photons that scatter on the fixed objects (scattered uniformly in the three-sphere). Smoothen the three-sphere histogram and extract the spherical harmonics coefficients for them and plot the power spectrum. We ask: do we see an anisotropy that gives any sense of the anisotropy. The answer to this question is yes (in a relatively quick implementation). This model can be refined to fit the actual CMB anisotropy I am confident. If I don’t succeed, someone else should do this because the Big Bang cosmology is insanely ambitious in claiming anything about the ‘origins of the universe’. A static Einstein universe would be infinitely better science with what we actually can know about the universe without any ability to do any experiments. Of course I am not doing any experiments either but I am doing real work: find a tight fit to the CMB anisotropy in a SIMPLE model. This cannot be that hard.

The code is in octave: The first example of anisotropy for the model is this little graph:

artificial-simple-anisotropy-static-universe-gravitational-redshiftonly

pkg load statistics

pkg load gsl

pkg load image

n=10000

q0=200

X = mvnrnd( [0,0,0,0],diag([1,1,1,1]),n);

for i=1:n,

X(i,:)=X(i,:)/norm(X(i,:),2);

endfor

# stereographic projection + normalize

function z=pwrrnd(k,n)

u=unifrnd(0,1,n)

z=(1-u)^(1/(1-k))*0.001

masses = exprnd(1,n)

function d=direction(x1,x2)

d=x2-x1

d=d/norm(d,2)

function d=distance(x1,x2)

d=acos(dot(x1,x2))

function f=gforce( x1,x2, m1,m2)

f=m1*m2*direction(x1,x2)/distance(x1,x2)^2

#f=0

timestep=0.1

for a=1:n,

force_on_x = 0

x1 = X(a,:)

m1 = masses(a)

for b=1:n,

x2 = X(b,:)

d=distance(x1,x2)

m2=masses(b)

if d>0 and d<pi,

force_on_x = force_on_x + gforce(x1,x2,m1,m2)

endif

endfor

x = x+0.5*force_on_x*timestep^2/m1

x = x/norm(x,2)

endfor

z = zeros(n,3);

for i=1:n,

y=X(i,:);

y1=y(1);

y2=y(2);

y3=y(3);

y4=y(4);

zp=[ y1/(1-y4), y2/(1-y4),y3/(1-y4) ];

z(i,:) = zp/norm(zp,2);

endfor

function t=fitangle(x)

t=x;

while t>2*pi,

t=t-2*pi;

endwhile

while t<0,

t=t+2*pi;

endwhile

endfunction

function b=bin_s2point( z )

theta = atan( z(2)/z(1));

theta=fitangle(theta);

phi = atan( sqrt(z(1)^2+z(2)^2)/z(3));

phi = fitangle(phi);

thetan = floor(theta*100/(2*pi));

phin = floor(phi*100/(2*pi));

b=[thetan,phin];

endfunction

# create a spherical grid in longitudes/latitudes

g = zeros(q0,q0);

refpoint = [1,0,0,0]

for i = 1:n,

cds = bin_s2point(z(i,:));

cds(1)=mod(cds(1),100);

cds(2)=mod(cds(2),100);

sdir = sign(direction(refpoint,X(i,:));

dist = distance(refpoint,X(i,:));

dg = zeros(q0,q0);

dg(cds(1)+1,cds(2)+1) = 1;

dg = imsmooth( dg, “Gaussian”, ceil(m1),ceil(m1));

g = g+dg;

endfor

#g = (g-min(min(g)));

#g=g/sum(sum(g));

h=imsmooth(g,”Gaussian”,10,10);

g=h;

# value at a spherical harmonic

function v=sphCoeff( l,m, g)

v = 0;

q0=100;

for j=1:q0,

for k=1:q0,

w1 = cos(2*pi*j/100);

w2 = exp(1i*2*pi*k/100);

fnval = gsl_sf_legendre_sphPlm( l, m, w1 )*w2;

v = v + fnval*g(j,k);

endfor

endfor

endfunction

function v=sumSphCoeff( l, g)

v=0;

for m=0:l

v = v+sphCoeff(l,m,g);

endfor

v=abs(v)^2

endfunction

valsByL = zeros(1,100);

for k=1:100,

valsByL(1,k)=sumSphCoeff(20*k,g);

endfor

plot( 1:100, valsByL)