Introduction
Packages
library(tidyverse)
library(plotly)
library(scales)
library(gganimate)
library(knitr)
library(DT)
library(colorRamps)Input Parameters
S <- 1370
A <- 204
B <- 2.17
K <- 3.86
ai <- 0.62
ab <- 0.25
aW <- 0.75
aB <- 0.25
gamma <- 2.2
delta <- 10/gamma
c <- 10
k <- 0.003265*0.75
T0 <- 20
D <- 0.3 # Death Rate
w0 <- 0.5
b0 <- 0.2
u0 <- 1-w0-b0Fit from Textbook Spreadsheet
zones <- c(85,75,65,55,45,35,25,15,6)
coszones <- cospi(zones/180)
sunWt <- c(0.5,0.531,0.624,0.77,0.892,1.021,1.12,1.189,1.219)
df <- data.frame(zones,sunWt)
f <- function(zones,d,n,k){
d*cos(n*zones)^2+k } #cos or cos^2 ??
Fit <- nls(sunWt ~ f(zones,d,n,k),data= df, start=list(d=1,n=0.015,k=0.3))
ggplot(df, aes(zones,sunWt))+geom_point()+geom_smooth(aes(y=f(zones,0.7768699,0.0164348,0.4617747)))+xlab("Latitude")+ylab("Solar Weight Factor")
Functional Forms
Solar Luminosity
\[ S_1(t)= \frac{S}{100}t \] \[ S_2(t)= S(\sin^2(t+90)) \] \[ S_3(t)= S\Big(1-\frac{1}{\sqrt{2\pi}}e^{-(t-50)^2/2}\Big) \] \[ S_4(t)=S\Big(1-\frac{1}{3}\delta(t-50)\Big) \] \[ S_5(t)=\frac{1}{100}(|t-150|+25) \] \[ S_6(t)=S\Big(1+\frac{1}{10} cos(t)\Big) \]
Albedo
\[ a(T)=\frac{e^{\gamma(T+\delta)}}{e^{\gamma(T+\delta)}+1} (\beta-\alpha) + \alpha \]
EBM
\[ T =\frac{ R_{in}(1-a(T))+KT_m-A }{ B+K } \]
Daisyworld
\[ a(T)=wa_w+ba_b+u\Big(\frac{e^{\gamma(T+\delta)}}{e^{\gamma(T+\delta)}+1} (\beta-\alpha) + \alpha\Big)\] \[ T_w=T+c(a-a_w) \] \[ T_b=T+c(a-a_b) \] \[ F_w=1-k(T_0-T_w)^2 \] \[ F_b=1-k(T_0-T_b)^2 \] \[ w'=w+w(uF_w-D) \] \[ b'=b+b(uF_b-D) \]
plot_albedo <- ggplot(data.frame(x= seq(-15,15, by=0.1),y=alb(seq(-15,15, by=0.1),ai,ab,gamma,delta)))+geom_line(aes(x,y),colour="blue")+xlab("Temperature")+ylab("Albedo")
plot_albedo
Run 0
data01 <- ebm01(150,A,B,K,ai,ab,gamma,delta) #from 300
data02 <-ebm02(150,A,B,K,ai,ab,gamma,delta) #from 300
plot01 <- ggplot(data01,aes(Zones)) +
geom_line(aes(y=T, colour = "Temperature",linetype="Temperature")) +
geom_line(aes(y=a*25, colour = "Albedo",linetype="Albedo"))+
geom_line(aes(y=Ti, colour = "Initial Temperature",linetype="Initial Temperature"))+xlab("Latitude")+ylab("Temperature")+scale_colour_manual(name="Legend",values=c("blue","red","red"))+scale_linetype_manual(name="Legend",values=c("Temperature"=1, "Albedo"=1, "Initial Temperature"=2))
plot02 <- ggplot(data02) +
geom_line(aes(t, Temperature),colour = 'green')+xlab("Time")+ylab("Mean Temperature")
plot01
plot02
Run 1
data11 <- ebm11(100,300,A,B,K,ai,ab,gamma,delta)
plot11 <- ggplot(data11,aes(J,Temp))+geom_line(aes(J, Temp),colour = 'red')+xlab("Time")+ylab("Mean Temperature")+ggtitle("Non Dynamic Mean Temperature (T(t) independent from T(t-1))")
plot11
TEMP1 <- ebm12(100,300,A,B,K,ai,ab,gamma,delta)
plot_ly(z=~TEMP1)%>% add_surface() %>% layout(
title = "With Initialization", scene = list(
xaxis = list(title = "S/100"),
yaxis = list(title = "Latitude"),
zaxis = list(title = "T") )) Run 2
data21 <- ebm21(180,300,A,B,K,ai,ab,gamma,delta )
plot21 <- ggplot(data21,aes(J,Temp2))+geom_line(aes(J, Temp2),colour = 'red')+xlab("Time")+ylab("Mean Temperature")+ggtitle("Sinusoidal Perturbation")
plot21
TEMP2 <- ebm22(180,300,A,B,K,ai,ab,gamma,delta )
plot_ly(z=~TEMP2) %>% add_surface() %>% layout(
title = "Without Initialization",scene = list(
xaxis = list(title = "S/100"),
yaxis = list(title = "Latitude"),
zaxis = list(title = "T") ))Run 3 & 4
data31 <- ebm31(200,300,A,B,K,ai,ab,gamma,delta )
data41 <- ebm41(200,300,A,B,K,ai,ab,gamma,delta )
plot31 <- ggplot(data31,aes(J,Temp))+geom_line(aes(J, Temp),colour = 'red')+xlab("Time")+ylab("Mean Temperature")+ggtitle("Gaussian Perturbation")
plot41 <- ggplot(data41,aes(J,Temp))+geom_line(aes(J, Temp),colour = 'red')+xlab("Time")+ylab("Mean Temperature")+ggtitle("Delta Perturbation")
plot31
plot41
TEMP3 <- ebm32(200,300,A,B,K,ai,ab,gamma,delta )
TEMP4 <- ebm42(200,300,A,B,K,ai,ab,gamma,delta)
plot_ly(z=~TEMP3) %>% add_surface()%>% layout(
title = "Without Initialization", scene = list(
xaxis = list(title = "S/100"),
yaxis = list(title = "Latitude"),
zaxis = list(title = "T") ))plot_ly(z=~TEMP4) %>% add_surface()%>% layout(
title = "Without Initialization",scene = list(
xaxis = list(title = "S/100"),
yaxis = list(title = "Latitude"),
zaxis = list(title = "T") ))Hysteresis Cycles
data23 <- ebm23(180,300,A,B,K,ai,ab,gamma,delta )
plot23 <- ggplot(data23,aes(Sarr,Temp2))+geom_point(aes(Sarr, Temp2),colour = 'yellow')+ylab("Mean Temperature")+xlab("S_2(t)")+ggtitle("Mean temperature vs. Sinusoidally Fluctuating Incoming Radiation")
plot23
data51 <- ebm51(300,60,A,B,K,ai,ab,gamma,delta )
plot51 <- ggplot(data51,aes(Sarr,Temp5,group=1))+geom_point(aes(Sarr, Temp5),colour = 'yellow')+geom_segment(aes(x = Sarr[89], y = Temp5[89], xend = Sarr[90], yend = Temp5[90]),colour = 'yellow')+geom_segment(aes(x = Sarr[258], y = Temp5[258], xend = Sarr[259], yend = Temp5[259]),colour = 'yellow')+ylab("Mean Temperature")+xlab("S_5(t)")+ggtitle("Mean temperature vs. Linearly Fluctuating Incoming Radiation")
plot51
Embedding Daisies in EBM
# WITHOUT DAISIES
data_ND1 <- ebm_ND1(500,w0,b0,A,B,K,ai,ab,aW,aB,gamma,delta) #from 500
data_ND2 <- ebm_ND2(500,w0,b0,A,B,K,ai,ab,aW,aB,gamma,delta)
data_ND4a <- ebm_ND4(0,300,1370/920,A,B,K,ai,ab,gamma,delta) #from 300
data_ND4b <- ebm_ND4(1,300,1370/920,A,B,K,ai,ab,gamma,delta) #from 300
TEMP <- ebm_ND3(500,w0,b0,A,B,K,ai,ab,aW,aB,gamma,delta)
plot_ND1 <- ggplot(data_ND1,aes(Zones))+geom_line(aes(y=w, colour = "% White"))+geom_line(aes(y=b, colour = "% Black"))+geom_line(aes(y=u, colour="% Bare Ground"))+geom_line(aes(y=T/5, colour="Temperature"))+ylab("Temperature & Albedo")+xlab("Latitude")+ggtitle("Without Daisies | Time Steps : 500")+scale_colour_manual(name="Legend",values=c("brown","black","white", "red"))
plot_ND2 <- ggplot(data_ND2,aes(I))+geom_line(aes(y=Barr,colour="% Black"))+ geom_line(aes(y=Warr,colour="% White"))+ geom_line(aes(y=Uarr, colour="% Bare Ground"))+geom_line(aes(y=Tarr/25,colour="Temperature"))+ylab("Temperature & Albedo")+xlab("Time")+ggtitle("Without Daisies | Time Steps : 500")+scale_colour_manual(name="Legend",values=c("brown","black","white", "red"))
plot_ND1 
plot_ND2
plot_ly(z=~TEMP)%>% add_surface() %>% layout( title="Temperature Without Daisies",
scene = list(
xaxis = list(title = "Time"),
yaxis = list(title = "Latitude"),
zaxis = list(title = "T") )) #WITH DAISIES
data_D1 <- ebm_D1(500,w0,b0,c,k,D,A,B,K,ai,ab,aW,aB,gamma,delta)
data_D2 <- ebm_D2(500,w0,b0,c,k,D,A,B,K,ai,ab,aW,aB,gamma,delta)
data_D4a <- ebm_D4(2,300,1370/920,w0,b0,c,k,D,A,B,K,ai,ab,aW,aB,gamma,delta)
data_D4b <- ebm_D4(3,300,1370/920,w0,b0,c,k,D,A,B,K,ai,ab,aW,aB,gamma,delta)
TEMP <- ebm_D3(500,w0,b0,c,k,D,A,B,K,ai,ab,aW,aB,gamma,delta)
BLACK <- ebm_Db2(500,w0,b0,c,k,D,A,B,K,ai,ab,aW,aB,gamma,delta)
WHITE <- ebm_Dw2(500,w0,b0,c,k,D,A,B,K,ai,ab,aW,aB,gamma,delta)
plot_D1 <- ggplot(data_D1 ,aes(Zones))+geom_line(aes(y=w, colour = "% White"))+geom_line(aes(y=b, colour = "% Black"))+geom_line(aes(y=u, colour="% Bare Ground"))+geom_line(aes(y=T/5, colour="Temperature"))+ylab("Temperature & Albedo")+xlab("Latitude")+ggtitle("With Daisies")+scale_colour_manual(name="Legend",values=c("brown","black","white", "red"))
plot_D2 <- ggplot(data_D2,aes(I))+geom_line(aes(y=Barr,colour="% Black"))+ geom_line(aes(y=Warr,colour="% White"))+ geom_line(aes(y=Uarr, colour="% Bare Ground"))+geom_line(aes(y=Tarr/25,colour="Temperature"))+ylab("Temperature & Albedo")+xlab("Time")+ggtitle("With Daisies")+scale_colour_manual(name="Legend",values=c("brown","black","white", "red"))
plot_D4a <- ggplot(data_D4a,aes(Sarr,Tarr))+ylab("Temperature")+xlab("Solar Luminosity")+ggtitle("Solar Luminosity vs. Temperature")+geom_line()
plot_D4b <- ggplot(data_D4b,aes(Time,Temp))+ylab("Temperature")+xlab("Time")+ggtitle("Dynamical Global Mean Temperature")+geom_line(colour="red")
plot_D1
plot_D2
plot_ly(z=~TEMP)%>% add_surface() %>% layout(title="Temperature With Daisies",
scene = list(
xaxis = list(title = "Time"),
yaxis = list(title = "Latitude"),
zaxis = list(title = "T") ))plot_ly(z=~BLACK)%>% add_surface() %>% layout(title="Distribution of Black Daisies",
scene = list(
xaxis = list(title = "Time"),
yaxis = list(title = "Latitude"),
zaxis = list(title = "Black") )) plot_ly(z=~WHITE)%>% add_surface() %>% layout(title="Distribution of White Daisies",
scene = list(
xaxis = list(title = "Time"),
yaxis = list(title = "Latitude"),
zaxis = list(title = "White") ))ebm_Db1(150,1370/920,w0,b0,c,k,D,A,B,K,ai,ab,aW,aB,gamma,delta)%>% layout(title="Distribution of Black Daisies",
xaxis = list(title = "Latitude"),
yaxis = list(title = "Solar Luminosity") )ebm_Dw1(150,1370/920,w0,b0,c,k,D,A,B,K,ai,ab,aW,aB,gamma,delta) %>% layout(title="Distribution of White Daisies",
xaxis = list(title = "Latitude"),
yaxis = list(title = "Solar Luminosity") )WHITE <- ebm_D4(0,300,1370/920,w0,b0,c,k,D,A,B,K,ai,ab,aW,aB,gamma,delta)
BLACK <- ebm_D4(1,300,1370/920,w0,b0,c,k,D,A,B,K,ai,ab,aW,aB,gamma,delta)
plot_ly(z=~BLACK)%>% add_surface() %>% layout(title="Distribution of Black Daisies",
scene = list(
xaxis = list(title = "Time"),
yaxis = list(title = "Solar Luminosity"),
zaxis = list(title = "Black") )) plot_ly(z=~WHITE)%>% add_surface() %>% layout(title="Distribution of White Daisies",
scene = list(
xaxis = list(title = "Time"),
yaxis = list(title = "Solar Luminosity"),
zaxis = list(title = "White") )) data_ND4_D4a <- data.frame(data_ND4a,data_D4a)
plot_ND4_D4a <- ggplot(data_ND4_D4a,aes(Sarr))+ylab("Temperature")+xlab("Solar Luminosity")+ggtitle("Solar Luminosity vs. Temperature")+geom_line(aes(y=Tarr, colour = "With Daisies"))+geom_line(aes(y=Tarr1, colour = "Without Daisies"))+scale_colour_manual(name="Legend",values=c("green","brown"))
plot_ND4_D4a
data_ND4_D4b <- data.frame(data_ND4b,data_D4b)
plot_ND4_D4b <- ggplot(data_ND4_D4b,aes(Time))+ylab("Temperature")+xlab("Time")+ggtitle("Dynamical Global Mean Temperature")+geom_line(aes(y=Temp, colour = "With Daisies"))+geom_line(aes(y=Temp1, colour = "Without Daisies"))+scale_colour_manual(name="Legend",values=c("green","brown"))
plot_ND4_D4b