In [1]:
#Robert Driving Sonny Nav
function Axb(A::Array{Float64,2},b::Array{Float64,1})
return A*b
end
Out[1]:
In [2]:
A = rand((2,2))
b = rand((2))
Out[2]:
In [3]:
print(Axb(A,b),A*b)
In [4]:
#Sonny driving, Robert Nav
function Axb_2(A::Array{Float64,2},b::Array{Float64,1})
x = zeros(size(b))
for i in 1:length(A[:,1])
for j in 1:length(A[1,:])
x[i] += A[i,j]*b[j]
end
end
return x
end
Out[4]:
In [5]:
print(Axb_2(A,b),A*b)
In [6]:
#Sonny driving, Robert Nav
function Axb_3(A::Array{Float64,2},b::Array{Float64,1})
x = zeros(size(b))
for j in 1:length(A[1,:])
for i in 1:length(A[:,1])
x[i] += A[i,j]*b[j]
end
end
return x
end
Out[6]:
In [7]:
print(Axb_3(A,b),A*b)
In [8]:
function testing_function(f1,f2,f3,dimension)
A = rand((dimension,dimension))
b = rand((dimension))
println(dimension,"x",dimension, " Matrix")
println("Function 1 Function 2 Function 3")
println((@elapsed f1(A,b))," ",(@elapsed f2(A,b))," ",(@elapsed f2(A,b)))
println()
end
Out[8]:
In [9]:
testing_function(Axb,Axb_2,Axb_3,8)
In [10]:
for i in [4,8,32,264,512,1024,1025,10000]
testing_function(Axb,Axb_2,Axb_3,i)
end
Ooopse we got our terminology mixed-up with or rows and columns. -R
Note the thrashing in the 1024 and 1025 cases. -S
In [11]:
#2b)
#Sonny driving
#Pkg.add("Calculus")
In [12]:
include("CLASS_REPO/Astro585/Lab5/HW5_Q2_lookup_table.jl")
Out[12]:
In [12]:
lookup_table = make_table_linear(x->tanh(x),0.,2.,10);
In [12]:
#2c)
#Sonny driving
using PyPlot
x = linspace(0.,2.,1000);
y1 = tanh(x);
y2 = lookup(lookup_table,x);
res = y2 - y1;
plot(x,y1,label="Function");
plot(x,y2,label="Lookup");
legend(loc="lower right");
#figure();
#plot(x,res,'o',label="Residuals")
#legend();
In [12]:
using PyPlot
plot(x,res,label="Residuals")
legend();
In [13]:
#robert driving
println("Function ", (@elapsed tanh(x)))
println("lookup ", (@elapsed lookup(lookup_table,x)))
In [13]:
#I stole this from stack exchange -R
function fact(n)
if n == 1
n
else
n * fact(n-1)
end
end;
In [14]:
function fn(x::Float64):
@assert x <= 20
@assert x >= 1.0
A = fact(iceil(x))
x -= 1
for i in 1:100
x += 1/i
end
return A*tanh(x)/(sin(x) + cos(x))
end
Out[14]:
In [14]:
LU1 = make_table_linear(x->fn(x),1.1,20.,100);
In [15]:
x = linspace(1.1,19.9,10000)
#robert driving
println("Function ", (@elapsed [fn(i) for i in x]))
println("lookup ", (@elapsed lookup(LU1,x)))
The more absurd we made the function, eventally we got the lookup table to beat it, but it shows the idea -R
In [16]:
include("CLASS_REPO/Astro585/Lab5/HW5_Q2_ecc_anom.jl")
Out[16]:
In [16]:
#Sonny Driving
x = linspace(-2*pi,4*pi,1000);
ecc_1 = ecc_anom(x,0.25);
using PyPlot
plot(x,ecc_1);
xlabel("Mean Anomaly");
ylabel("Eccentric Anomaly");
This would perform better than calculating a lot of eccentricity anomalies (there's a while statement in ecc_anom). -S
In [16]:
ecc_table1 = make_table_linear(x->ecc_anom(x,0.25),0.,2*pi,128);
In [16]:
#Sonny Driving
x = linspace(0.,2*pi,1000);
ecc_2 = ecc_anom(x,0.25);
using PyPlot
plot(x,ecc_2);
plot(x,lookup(ecc_table1,x))
ylim(-7,7);
xlabel("Mean Anomaly");
ylabel("Eccentric Anomaly");
The artifacting has to be something from the interpolation process - otherwise, the ecc_anom call would also show it. -S
In [18]:
#Robert Driving
ecc_table2 = make_table_linear(x->ecc_anom(x,0.25),x->decc_anom_dM(x,0.25),0.,2*pi,128)
Out[18]:
In [18]:
#Robert Driving
x = linspace(0.,2*pi,1000);
ecc_2 = ecc_anom(x,0.25);
using PyPlot
plot(x,ecc_2);
plot(x,lookup(ecc_table2,x))
ylim(-2,7);
xlabel("Mean Anomaly");
ylabel("Eccentric Anomaly");
In [21]:
decc = [decc_anom_dM(i,0.25) for i in x ]
Out[21]:
In [22]:
plot(x,decc)
Out[22]:
In [ ]: