alex/CIFS_Storage_Demand_Model
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Added variants

Browse Source
This commit is contained in:
Alex Gebben Work 2026-02-05 16:57:54 -07:00
parent aa6c892ced
commit 71f662d9ed
5 changed files with 272 additions and 14 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

17
Analysis.r
View File

@ -1,19 +1,16 @@
library(tidyverse) library(tidyverse)
Rev_No_Shipping <- readRDS("Data/Cleaned_Data/Model_Estimates.Rds") Rev_No_Shipping <- readRDS("Data/Results/Model_Estimates.Rds")
Rev_No_Shipping_2018 <- readRDS("Data/Cleaned_Data/Model_Estimates_2018.Rds") Rev_No_Shipping_2018 <- readRDS("Data/Results/Model_Estimates_2018.Rds")
Rev_Shipping <- readRDS("Data/Cleaned_Data/Model_Estimates_With_Shipping_Costs.Rds") #Rev_Shipping <- readRDS("Data/Results/Model_Estimates_With_Shipping_Costs.Rds")
Rev_Shipping_2018 <- readRDS("Data/Cleaned_Data/Model_Estimates_With_Shipping_Costs_2018.Rds") #Rev_Shipping_2018 <- readRDS("Data/Results/Model_Estimates_With_Shipping_Costs_2018.Rds")
TEST <- Rev_No_Shipping TEST <- Rev_No_Shipping
#TEST <- Rev_No_Shipping_2018
TEST <- TEST %>% mutate(Current_Profit=Profit/((1+Discount)^(Year-2026))) %>% ungroup TEST <- TEST %>% mutate(Current_Profit=Profit/((1+Discount)^(Year-2026))) %>% ungroup
TEST <- TEST %>% group_by(Year,Location,Discount) %>% filter(Current_Profit==max(Current_Profit)) %>% ungroup TEST <- TEST %>% group_by(Year,Location,Discount) %>% filter(Current_Profit==max(Current_Profit)) %>% ungroup
TEST %>% select(Profit,Current_Profit,Year,Discount) %>% tail
TEST %>% group_by(Location,Discount) %>% filter(Current_Profit==max(Current_Profit)) %>% print(n=100) %>%
TEMP <- TEST %>% group_by(Location,Discount) %>% select(Year,Capacity,Profit,Current_Profit) %>% mutate(Delay_Profit=(lead(Current_Profit)-Current_Profit)/10^6)
TEMP %>% filter(Year>=2025) %>% summarize(Optimal_Year=sum(ifelse(Current_Profit==max(Current_Profit),Year,0)),Max_Profit=max(Current_Profit),Current_Profit=sum(ifelse(Year==2026,Current_Profit,0)),Yearly_Loss=(Max_Profit-Current_Profit)/(2026-Optimal_Year)/10^6,Next_Years_Loss=(Max_Profit-sum(ifelse(Year==(Optimal_Year+1),Current_Profit,0)))/10^6)
TEST %>% group_by(Location,Discount,
ggplot(TEST %>% filter(Discount==0.03) ,aes(x=Year,y=Current_Profit/10^9,color=Capacity))+geom_line()+facet_wrap(Discount~Location,ncol=1)+scale_x_continuous(breaks=seq(1960,2083,by=5)) ggplot(TEST %>% filter(Discount==0.03) ,aes(x=Year,y=Current_Profit/10^9,color=Capacity))+geom_line()+facet_wrap(Discount~Location,ncol=1)+scale_x_continuous(breaks=seq(1960,2083,by=5))
@ -23,5 +20,5 @@ ggplot(TEST %>% filter(Discount %in% c(0.03,0.05,0.07) ,Year>1970) ,aes(x=Year,y
ggplot(TEST %>% filter(Discount %in% c(0.05) ,Year>1970) ,aes(x=Year,y=Profit/10^9,color=Location))+geom_line()+scale_x_continuous(breaks=seq(1960,2083,by=5)) ggplot(TEST %>% filter(Discount %in% c(0.05) ,Year>1970) ,aes(x=Year,y=Profit/10^9,color=Location))+geom_line()+scale_x_continuous(breaks=seq(1960,2083,by=5))
ggplot(TEST %>% filter(Discount %in% c(0.03,0.05,0.07),Year>1970,Year<2050) ,aes(x=Year,y=Marginal/10^6))+ geom_area(aes(y=ifelse(Marginal>0,Marginal/10^6,0)),fill="lightgreen",alpha=0.6)+ geom_area(aes(y=ifelse(Marginal<0,Marginal/10^6,0)),fill="firebrick2",alpha=0.5) +facet_wrap(Discount~Location,ncol=2,scales = "free_y")+scale_x_continuous(breaks=seq(1960,2083,by=5))+scale_y_continuous(breaks=seq(-300,50,by=50))+geom_hline(yintercept=0)+theme_bw()+ylab("Profit Change (Million Dollars)") ggplot(TEST %>% filter(Discount %in% c(0.03,0.05,0.07),Year>1970,Year<2050) ,aes(x=Year,y=Marginal/10^6))+ geom_area(aes(y=ifelse(Marginal>0,Marginal/10^6,0)),fill="lightgreen",alpha=0.6)+ geom_area(aes(y=ifelse(Marginal<0,Marginal/10^6,0)),fill="firebrick2",alpha=0.5) +facet_wrap(Discount~Location,ncol=2,scales = "free_y")+scale_x_continuous(breaks=seq(1960,2083,by=5))+scale_y_continuous(breaks=seq(-300,300,by=50))+geom_hline(yintercept=0)+theme_bw()+ylab("Profit Change (Million Dollars)")

2
Run_All.sh
View File

@ -1,2 +1,4 @@
Rscript ./Scripts/1_Raw_Cost_Data_Clean.r Rscript ./Scripts/1_Raw_Cost_Data_Clean.r
Rscript ./Scripts/2_Compiled_Results_Data.r Rscript ./Scripts/2_Compiled_Results_Data.r
Rscript ./Scripts/3_Half_Cost_Results.r
Rscript ./Scripts/4_Old_Reactor_Cost_Added.r

8
Scripts/2_Compiled_Results_Data.r
View File

@ -4,9 +4,10 @@ library(parallel)
NCORES <- detectCores()-1 NCORES <- detectCores()-1
library(lpSolve) #For solving discrete value maximization for the power plants library(lpSolve) #For solving discrete value maximization for the power plants
####Manual inputs ####Manual inputs
#DISCOUNT_RATE_LIST <- seq(0,1,by=0.0025)
#Range of discount rates to calculate in the model. Each facility will have each rate calculated, so more values slows the results but allows for more discount rates to be reported in the findings. #Range of discount rates to calculate in the model. Each facility will have each rate calculated, so more values slows the results but allows for more discount rates to be reported in the findings.
DISCOUNT_RATE_LIST <- c(0.03,0.0325,0.035,0.0375,0.04,0.045,0.0475,0.05,0.07,0.1) #DISCOUNT_RATE_LIST <- c(0.03,0.05,0.07)
#DISCOUNT_RATE_LIST <- c(0.03,0.0325,0.035,0.0375,0.04,0.045,0.0475,0.05,0.07,0.1)
DISCOUNT_RATE_LIST <- seq(0.01,0.15,by=0.0025)
#The cost per ton of shipping uranium, used to see what can be shipped on day one of the project. #The cost per ton of shipping uranium, used to see what can be shipped on day one of the project.
SHIPPING_COST_PER_TON <- 1.2874*26000 #Inflation adjusted from New Mexico Report SHIPPING_COST_PER_TON <- 1.2874*26000 #Inflation adjusted from New Mexico Report
#The savings per year of having a CIFS at a served reactor (cost to house the CIFS). #The savings per year of having a CIFS at a served reactor (cost to house the CIFS).
@ -106,8 +107,7 @@ TOTAL_VALUE_METRICS <- MULTI_DISCOUNT_RATE_NPV(DISCOUNT_RATE_LIST,TOTAL ,DOLLARS
TOTAL_VALUE_METRICS_ORIG <- MULTI_DISCOUNT_RATE_NPV(DISCOUNT_RATE_LIST,TOTAL_ORIG ,DOLLARS_SAVED_PER_YEAR=CV) %>% left_join(read_csv("Data/Raw_Data/Curie_Spent_Fuel_Site_Totals.csv")) %>% select(Year,Facility,Discount,Total_Tons,Revenue) TOTAL_VALUE_METRICS_ORIG <- MULTI_DISCOUNT_RATE_NPV(DISCOUNT_RATE_LIST,TOTAL_ORIG ,DOLLARS_SAVED_PER_YEAR=CV) %>% left_join(read_csv("Data/Raw_Data/Curie_Spent_Fuel_Site_Totals.csv")) %>% select(Year,Facility,Discount,Total_Tons,Revenue)
####Find the results for each facility size ####Find the results for each facility size
#Main Results #Main Results
CREATE_FULL_RESULTS <- function(REVENUE_RESULTS,COST_RESULTS){return(REVENUE_RESULTS %>% left_join(COST_RESULTS) %>% mutate(Profit=Revenue-Cost) %>% group_by(Discount,Capacity,Location) %>% mutate(Time_Benefit=(1-Discount)*(lead(Profit)-Profit),Op_Cost=Discount*Profit,Marginal=Time_Benefit-Op_Cost) %>% unique)}
CREATE_FULL_RESULTS <- function(REVENUE_RESULTS,COST_RESULTS){return(REVENUE_RESULTS %>% left_join(COST_COST_RESULTS) %>% mutate(Profit=Revenue-Cost) %>% group_by(Discount,Capacity,Location) %>% mutate(Time_Benefit=(1-Discount)*(lead(Profit)-Profit),Op_Cost=Discount*Profit,Marginal=Time_Benefit-Op_Cost) %>% unique)}
Rev_No_Shipping <- do.call(rbind,mclapply(CAPACITY_LIST,function(CAPACITY){do.call(rbind,lapply(DISCOUNT_RATE_LIST,function(x){YEARLY_RESULTS(TOTAL_VALUE_METRICS,x,CAPACITY,0)}))},mc.cores=min(length(CAPACITY_LIST),NCORES))) Rev_No_Shipping <- do.call(rbind,mclapply(CAPACITY_LIST,function(CAPACITY){do.call(rbind,lapply(DISCOUNT_RATE_LIST,function(x){YEARLY_RESULTS(TOTAL_VALUE_METRICS,x,CAPACITY,0)}))},mc.cores=min(length(CAPACITY_LIST),NCORES)))
saveRDS(Rev_No_Shipping,paste0(INTERMEDIATE_DIR ,"Revenue_Estimates.Rds")) saveRDS(Rev_No_Shipping,paste0(INTERMEDIATE_DIR ,"Revenue_Estimates.Rds"))
Rev_No_Shipping <- CREATE_FULL_RESULTS(Rev_No_Shipping,COST_DATA) Rev_No_Shipping <- CREATE_FULL_RESULTS(Rev_No_Shipping,COST_DATA)

125
Scripts/3_Half_Cost_Results.r Normal file
View File

@ -0,0 +1,125 @@
#A script which attempts to pull in all data, and create a data frame with the maximum revenue values for each facility, year and discount rate. The output can then be used to make figures and graphs
library(tidyverse)
library(parallel)
NCORES <- detectCores()-1
library(lpSolve) #For solving discrete value maximization for the power plants
####Manual inputs
#Range of discount rates to calculate in the model. Each facility will have each rate calculated, so more values slows the results but allows for more discount rates to be reported in the findings.
#DISCOUNT_RATE_LIST <- c(0.03,0.05,0.07)
#DISCOUNT_RATE_LIST <- c(0.03,0.0325,0.035,0.0375,0.04,0.045,0.0475,0.05,0.07,0.1)
DISCOUNT_RATE_LIST <- seq(0.01,0.15,by=0.0025)
#The cost per ton of shipping uranium, used to see what can be shipped on day one of the project.
SHIPPING_COST_PER_TON <- 1.2874*26000 #Inflation adjusted from New Mexico Report
#The savings per year of having a CIFS at a served reactor (cost to house the CIFS).
CV <- 1.2874*(6984013)/2 #Data from New Mexico Report, Converted from 2019 to Dec 2025
#Locations to save results
RES_DIR <- "./Data/Results/"
#Directory where the individual CIFS project plan cost data can be found.
#This has data has already been combined with the original files (low,high) and inflation adjusted.
CIFS_INDIVIDUAL_COST_DATA_DIR <- 'Data/Cleaned_Data/'
#Create any need save locations
dir.create(RES_DIR,recursive=TRUE,showWarnings=FALSE)
###################################Cost results
CIFS <- rbind(readRDS(paste0(CIFS_INDIVIDUAL_COST_DATA_DIR,"Texas_CIFS_Costs.Rds")),readRDS(paste0(CIFS_INDIVIDUAL_COST_DATA_DIR,"New_Mexico_CIFS_Costs.Rds")))
#Adjust for inflation
#CIFS[,-1:-5] <- 1.2874*CIFS[,-1:-5]
COSTS <- do.call(rbind,lapply(DISCOUNT_RATE_LIST ,function(DISCOUNT){CIFS %>% group_by(Location,Capacity,Cost_Assumption) %>% mutate(NPC=Total/(1+DISCOUNT)^Year) %>% summarize(Discount=DISCOUNT,Costs=sum(NPC))})) %>% ungroup
TEMP <- COSTS%>% group_by(Location,Cost_Assumption,Discount) %>% summarize(ST_COST=min(Costs),ST_CAPACITY=min(Capacity),M_COST=(max(Costs)-min(Costs))/(max(Capacity)-min(Capacity))) %>% ungroup
CAPACITY_INCREMENT <- 5000
#rm(COST_DATA)
for(i in 1:nrow(TEMP)){
LOC <- as.character(TEMP[i,"Location"])
COST_LEVEL <- as.character(TEMP[i,"Cost_Assumption"])
DISCOUNT <- as.numeric(TEMP[i,"Discount"])
ST_CAP <- as.numeric(TEMP[i,"ST_CAPACITY"])
ST_COST <- as.numeric(TEMP[i,"ST_COST"])
COST_SLOPE <- as.numeric(TEMP[i,]$M_COST)
CAPACITY <- seq(ST_CAP,160000,by=5000)
COST <- ST_COST+COST_SLOPE*CAPACITY_INCREMENT*(0:(length(CAPACITY)-1))
C_RES <- cbind(CAPACITY,COST) %>% as_tibble %>% mutate(Location=LOC,Cost_Assumption=COST_LEVEL,Discount=DISCOUNT) %>% select(Location,Cost_Assumption,Discount,Capacity=CAPACITY,Cost=COST)
if(!exists("COST_DATA")){COST_DATA <- C_RES}else{COST_DATA<- rbind(COST_DATA,C_RES)}
}
COST_DATA <- COST_DATA %>% filter(Cost_Assumption=='High') %>% select(-Cost_Assumption) %>% unique
##All unique capacity levels that the revenues need to be calculated for
CAPACITY_LIST <- COST_DATA %>% pull(Capacity) %>% unique
###
TOTAL <- read_csv("Data/Raw_Data/Curie_Spent_Fuel_Site_Totals.csv") %>% mutate(OP_YEAR=year(Op_Date_Min),CLOSE_YEAR=year(Close_Date_Max))%>% select(Facility,Total_Assemblies,Total_Tons,OP_YEAR,CLOSE_YEAR)
FACILITY_LIST <- TOTAL %>% pull(Facility)
#https://www.nrc.gov/reactors/operating/licensing/renewal/subsequent-license-renewal
SUBMITTED <-rbind(c(FACILITY_LIST[str_detect(FACILITY_LIST,"Duane*" )],2025),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Nine Mile*" )],2026),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Ginna*" )],2026),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Cooper*" )],2026),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Farley*" )],2027),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Prairie*" )],2027),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Brunswick*" )],2027),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Cook" )],2027),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Hope" )],2027),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Salem" )],2027),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Perry" )],2027),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Millstone" )],2028),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Palisades" )],2028),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Beaver" )],2028),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Callaway" )],2029),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Three Mile Island" )],2029),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Davis-Besse" )],2029),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Wolf" )],2030),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Lucie" )],2021),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Robinson" )],2025),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Hatch" )],2025))
SUBMITTED <- SUBMITTED %>% as_tibble
colnames(SUBMITTED ) <- c("Facility","App_Date","Status")
SUBMITTED <- SUBMITTED%>% mutate(Status="Applied",App_Date=as.numeric(App_Date)) %>% select(Facility,Status,App_Date)
#Issued
RENEWED <- rbind(c(FACILITY_LIST[str_detect(FACILITY_LIST,"Turkey" )],"Granted",2018,2033),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Peach" )],"Granted",2019,2034),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Surry" )],"Granted",2020,2033),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"North" )],"Granted",2021,2040),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Monticello" )],"Granted",2024,2030),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Oconee" )],"Granted",2025,2034),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Summer" )],"Granted",2025,2042),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Beach" )],"Granted",2025,2033),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Browns" )],"Granted",2025,2036),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Dresden" )],"Granted",2025,2031))
RENEWED <- RENEWED %>% as_tibble
colnames(RENEWED) <- c("Facility","Status","App_Date","Op_Date")
RENEWED <- RENEWED %>% mutate(Op_Date=as.numeric(Op_Date),App_Date=as.numeric(App_Date))
AVG_LENGTH <- RENEWED %>% mutate(DIFF=Op_Date-App_Date) %>% pull(DIFF) %>% mean %>% round
SUBMITTED <- SUBMITTED %>% mutate(Op_Date=App_Date+AVG_LENGTH)
UPDATE <- rbind(RENEWED,SUBMITTED )
TOTAL_ORIG <- TOTAL
TOTAL <- TOTAL %>% left_join(UPDATE) %>% mutate(CLOSE_YEAR=ifelse(Op_Date>CLOSE_YEAR & !is.na(Status),Op_Date,CLOSE_YEAR)) %>% select(-Status,-App_Date,-Op_Date)
source("./Scripts/Functions/NPV_Functions.r")
#Calculate and individual facilities rate, for all discount rates, and all years
TOTAL_VALUE_METRICS <- MULTI_DISCOUNT_RATE_NPV(DISCOUNT_RATE_LIST,TOTAL ,DOLLARS_SAVED_PER_YEAR=CV)%>% left_join(read_csv("Data/Raw_Data/Curie_Spent_Fuel_Site_Totals.csv")) %>% select(Year,Facility,Discount,Total_Tons,Revenue)
TOTAL_VALUE_METRICS_ORIG <- MULTI_DISCOUNT_RATE_NPV(DISCOUNT_RATE_LIST,TOTAL_ORIG ,DOLLARS_SAVED_PER_YEAR=CV) %>% left_join(read_csv("Data/Raw_Data/Curie_Spent_Fuel_Site_Totals.csv")) %>% select(Year,Facility,Discount,Total_Tons,Revenue)
####Find the results for each facility size
#Main Results
CREATE_FULL_RESULTS <- function(REVENUE_RESULTS,COST_RESULTS){return(REVENUE_RESULTS %>% left_join(COST_RESULTS) %>% mutate(Profit=Revenue-Cost) %>% group_by(Discount,Capacity,Location) %>% mutate(Time_Benefit=(1-Discount)*(lead(Profit)-Profit),Op_Cost=Discount*Profit,Marginal=Time_Benefit-Op_Cost) %>% unique)}
Rev_No_Shipping <- do.call(rbind,mclapply(CAPACITY_LIST,function(CAPACITY){do.call(rbind,lapply(DISCOUNT_RATE_LIST,function(x){YEARLY_RESULTS(TOTAL_VALUE_METRICS,x,CAPACITY,0)}))},mc.cores=min(length(CAPACITY_LIST),NCORES)))
Rev_No_Shipping <- CREATE_FULL_RESULTS(Rev_No_Shipping,COST_DATA)
saveRDS(Rev_No_Shipping,paste0(RES_DIR,"Model_Estimates_Half_Cost_Storage.Rds"))
#Same as main results but the optimization is based on information known in 2018, using unadjusted CURIE map data
Rev_No_Shipping_2018 <- do.call(rbind,mclapply(CAPACITY_LIST,function(CAPACITY){do.call(rbind,lapply(DISCOUNT_RATE_LIST,function(x){YEARLY_RESULTS(TOTAL_VALUE_METRICS_ORIG,x,CAPACITY,0)}))},mc.cores=min(length(CAPACITY_LIST),NCORES)))
Rev_No_Shipping_2018 <- CREATE_FULL_RESULTS(Rev_No_Shipping_2018,COST_DATA)
saveRDS(Rev_No_Shipping_2018,paste0(RES_DIR,"Model_Estimates_2018_Half_Cost_Storage.Rds"))
#Same as main results but only the reactors with a value ready to ship are included. That is to say other reactors may be willing to pay for future rights to hold the spent fuel but only those reactors with current values ready to ship are included, which set the bounds of the starting CIFS size. This is used for demonstration.
Rev_Shipping <- do.call(rbind,mclapply(CAPACITY_LIST,function(CAPACITY){do.call(rbind,lapply(DISCOUNT_RATE_LIST,function(x){YEARLY_RESULTS(TOTAL_VALUE_METRICS,x,CAPACITY,SHIPPING_COST_PER_TON )}))},mc.cores=min(length(CAPACITY_LIST),NCORES)))
Rev_Shipping <- CREATE_FULL_RESULTS(Rev_Shipping,COST_DATA)
saveRDS(Rev_Shipping,paste0(RES_DIR,"Model_Estimates_With_Shipping_Costs_Half_Cost_Storage.Rds"))
#Same as Rev_Shipping, but only information available in 2018 from the CURIE map is used (no updated operation dates).
Rev_Shipping_2018 <- do.call(rbind,mclapply(CAPACITY_LIST,function(CAPACITY){do.call(rbind,lapply(DISCOUNT_RATE_LIST,function(x){YEARLY_RESULTS(TOTAL_VALUE_METRICS_ORIG,x,CAPACITY,SHIPPING_COST_PER_TON )}))},mc.cores=min(length(CAPACITY_LIST),NCORES)))
Rev_Shipping_2018 <- CREATE_FULL_RESULTS(Rev_Shipping_2018,COST_DATA)
saveRDS(Rev_Shipping_2018,paste0(RES_DIR,"Model_Estimates_With_Shipping_Costs_2018_Half_Cost_Storage.Rds"))

134
Scripts/4_Old_Reactor_Cost_Added.r Normal file
View File

@ -0,0 +1,134 @@
#A script which attempts to pull in all data, and create a data frame with the maximum revenue values for each facility, year and discount rate. The output can then be used to make figures and graphs
library(tidyverse)
library(parallel)
NCORES <- detectCores()-1
library(lpSolve) #For solving discrete value maximization for the power plants
####Manual inputs
#Range of discount rates to calculate in the model. Each facility will have each rate calculated, so more values slows the results but allows for more discount rates to be reported in the findings.
#DISCOUNT_RATE_LIST <- c(0.03,0.05,0.07)
#DISCOUNT_RATE_LIST <- c(0.03,0.0325,0.035,0.0375,0.04,0.045,0.0475,0.05,0.07,0.1)
DISCOUNT_RATE_LIST <- seq(0.01,0.15,by=0.0025)
#The cost per ton of shipping uranium, used to see what can be shipped on day one of the project.
SHIPPING_COST_PER_TON <- 1.2874*26000 #Inflation adjusted from New Mexico Report
#The savings per year of having a CIFS at a served reactor (cost to house the CIFS).
CV <- 1.2874*(6984013) #Data from New Mexico Report, Converted from 2019 to Dec 2025
#Locations to save results
RES_DIR <- "./Data/Results/"
#Directory where the individual CIFS project plan cost data can be found.
#This has data has already been combined with the original files (low,high) and inflation adjusted.
CIFS_INDIVIDUAL_COST_DATA_DIR <- 'Data/Cleaned_Data/'
#Create any need save locations
dir.create(RES_DIR,recursive=TRUE,showWarnings=FALSE)
###################################Cost results
CIFS <- rbind(readRDS(paste0(CIFS_INDIVIDUAL_COST_DATA_DIR,"Texas_CIFS_Costs.Rds")),readRDS(paste0(CIFS_INDIVIDUAL_COST_DATA_DIR,"New_Mexico_CIFS_Costs.Rds")))
#Adjust for inflation
#CIFS[,-1:-5] <- 1.2874*CIFS[,-1:-5]
COSTS <- do.call(rbind,lapply(DISCOUNT_RATE_LIST ,function(DISCOUNT){CIFS %>% group_by(Location,Capacity,Cost_Assumption) %>% mutate(NPC=Total/(1+DISCOUNT)^Year) %>% summarize(Discount=DISCOUNT,Costs=sum(NPC))})) %>% ungroup
TEMP <- COSTS%>% group_by(Location,Cost_Assumption,Discount) %>% summarize(ST_COST=min(Costs),ST_CAPACITY=min(Capacity),M_COST=(max(Costs)-min(Costs))/(max(Capacity)-min(Capacity))) %>% ungroup
CAPACITY_INCREMENT <- 5000
#rm(COST_DATA)
for(i in 1:nrow(TEMP)){
LOC <- as.character(TEMP[i,"Location"])
COST_LEVEL <- as.character(TEMP[i,"Cost_Assumption"])
DISCOUNT <- as.numeric(TEMP[i,"Discount"])
ST_CAP <- as.numeric(TEMP[i,"ST_CAPACITY"])
ST_COST <- as.numeric(TEMP[i,"ST_COST"])
COST_SLOPE <- as.numeric(TEMP[i,]$M_COST)
CAPACITY <- seq(ST_CAP,160000,by=5000)
COST <- ST_COST+COST_SLOPE*CAPACITY_INCREMENT*(0:(length(CAPACITY)-1))
C_RES <- cbind(CAPACITY,COST) %>% as_tibble %>% mutate(Location=LOC,Cost_Assumption=COST_LEVEL,Discount=DISCOUNT) %>% select(Location,Cost_Assumption,Discount,Capacity=CAPACITY,Cost=COST)
if(!exists("COST_DATA")){COST_DATA <- C_RES}else{COST_DATA<- rbind(COST_DATA,C_RES)}
}
#COST_DATA <- COST_DATA %>% filter(Cost_Assumption=='Average') %>% select(-Cost_Assumption) %>% unique
COST_DATA <- COST_DATA %>% filter(Cost_Assumption=='High') %>% select(-Cost_Assumption) %>% unique
##All unique capacity levels that the revenues need to be calculated for
CAPACITY_LIST <- COST_DATA %>% pull(Capacity) %>% unique
###
TOTAL <- read_csv("Data/Raw_Data/Curie_Spent_Fuel_Site_Totals.csv") %>% mutate(OP_YEAR=year(Op_Date_Min),CLOSE_YEAR=year(Close_Date_Max))%>% select(Facility,Total_Assemblies,Total_Tons,OP_YEAR,CLOSE_YEAR)
FACILITY_LIST <- TOTAL %>% pull(Facility)
#https://www.nrc.gov/reactors/operating/licensing/renewal/subsequent-license-renewal
SUBMITTED <-rbind(c(FACILITY_LIST[str_detect(FACILITY_LIST,"Duane*" )],2025),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Nine Mile*" )],2026),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Ginna*" )],2026),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Cooper*" )],2026),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Farley*" )],2027),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Prairie*" )],2027),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Brunswick*" )],2027),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Cook" )],2027),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Hope" )],2027),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Salem" )],2027),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Perry" )],2027),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Millstone" )],2028),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Palisades" )],2028),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Beaver" )],2028),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Callaway" )],2029),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Three Mile Island" )],2029),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Davis-Besse" )],2029),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Wolf" )],2030),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Lucie" )],2021),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Robinson" )],2025),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Hatch" )],2025))
SUBMITTED <- SUBMITTED %>% as_tibble
colnames(SUBMITTED ) <- c("Facility","App_Date","Status")
SUBMITTED <- SUBMITTED%>% mutate(Status="Applied",App_Date=as.numeric(App_Date)) %>% select(Facility,Status,App_Date)
#Issued
RENEWED <- rbind(c(FACILITY_LIST[str_detect(FACILITY_LIST,"Turkey" )],"Granted",2018,2033),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Peach" )],"Granted",2019,2034),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Surry" )],"Granted",2020,2033),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"North" )],"Granted",2021,2040),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Monticello" )],"Granted",2024,2030),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Oconee" )],"Granted",2025,2034),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Summer" )],"Granted",2025,2042),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Beach" )],"Granted",2025,2033),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Browns" )],"Granted",2025,2036),
c(FACILITY_LIST[str_detect(FACILITY_LIST,"Dresden" )],"Granted",2025,2031))
RENEWED <- RENEWED %>% as_tibble
colnames(RENEWED) <- c("Facility","Status","App_Date","Op_Date")
RENEWED <- RENEWED %>% mutate(Op_Date=as.numeric(Op_Date),App_Date=as.numeric(App_Date))
AVG_LENGTH <- RENEWED %>% mutate(DIFF=Op_Date-App_Date) %>% pull(DIFF) %>% mean %>% round
SUBMITTED <- SUBMITTED %>% mutate(Op_Date=App_Date+AVG_LENGTH)
UPDATE <- rbind(RENEWED,SUBMITTED )
TOTAL_ORIG <- TOTAL
TOTAL <- TOTAL %>% left_join(UPDATE) %>% mutate(CLOSE_YEAR=ifelse(Op_Date>CLOSE_YEAR & !is.na(Status),Op_Date,CLOSE_YEAR)) %>% select(-Status,-App_Date,-Op_Date)
source("./Scripts/Functions/NPV_Functions.r")
#Custom function overide to add in a cost of 25 million for old sites
VALUE_ADD <- function(r,CURRENT_YEAR,CLOSE_YEAR){
Years_Until_Close <- max(CLOSE_YEAR-CURRENT_YEAR+1,0)
VALUES <- (1+r)^-(1:10^4)
CV <- 1.2874*(6984013) #Data from New Mexico Report, Converted from 2019 to Dec 2025
OLD_ADJUST <- (25*10^6)/(CV)
if(Years_Until_Close==0){return(sum(VALUES))}else if(Years_Until_Close>40){return(0)}else if(Years_Until_Close>10){return(sum(VALUES[-1:-Years_Until_Close])-OLD_ADJUST)} else{return(sum(VALUES[-1:-Years_Until_Close]))}
}
#Calculate and individual facilities rate, for all discount rates, and all years
TOTAL_VALUE_METRICS <- MULTI_DISCOUNT_RATE_NPV(DISCOUNT_RATE_LIST,TOTAL ,DOLLARS_SAVED_PER_YEAR=CV)%>% left_join(read_csv("Data/Raw_Data/Curie_Spent_Fuel_Site_Totals.csv")) %>% select(Year,Facility,Discount,Total_Tons,Revenue)
TOTAL_VALUE_METRICS_ORIG <- MULTI_DISCOUNT_RATE_NPV(DISCOUNT_RATE_LIST,TOTAL_ORIG ,DOLLARS_SAVED_PER_YEAR=CV) %>% left_join(read_csv("Data/Raw_Data/Curie_Spent_Fuel_Site_Totals.csv")) %>% select(Year,Facility,Discount,Total_Tons,Revenue)
####Find the results for each facility size
#Main Results
CREATE_FULL_RESULTS <- function(REVENUE_RESULTS,COST_RESULTS){return(REVENUE_RESULTS %>% left_join(COST_RESULTS) %>% mutate(Profit=Revenue-Cost) %>% group_by(Discount,Capacity,Location) %>% mutate(Time_Benefit=(1-Discount)*(lead(Profit)-Profit),Op_Cost=Discount*Profit,Marginal=Time_Benefit-Op_Cost) %>% unique)}
Rev_No_Shipping <- do.call(rbind,mclapply(CAPACITY_LIST,function(CAPACITY){do.call(rbind,lapply(DISCOUNT_RATE_LIST,function(x){YEARLY_RESULTS(TOTAL_VALUE_METRICS,x,CAPACITY,0)}))},mc.cores=min(length(CAPACITY_LIST),NCORES)))
Rev_No_Shipping <- CREATE_FULL_RESULTS(Rev_No_Shipping,COST_DATA)
saveRDS(Rev_No_Shipping,paste0(RES_DIR,"Model_Estimates_Old_Sites_Cost_Addition.Rds"))
#Same as main results but the optimization is based on information known in 2018, using unadjusted CURIE map data
Rev_No_Shipping_2018 <- do.call(rbind,mclapply(CAPACITY_LIST,function(CAPACITY){do.call(rbind,lapply(DISCOUNT_RATE_LIST,function(x){YEARLY_RESULTS(TOTAL_VALUE_METRICS_ORIG,x,CAPACITY,0)}))},mc.cores=min(length(CAPACITY_LIST),NCORES)))
Rev_No_Shipping_2018 <- CREATE_FULL_RESULTS(Rev_No_Shipping_2018,COST_DATA)
saveRDS(Rev_No_Shipping_2018,paste0(RES_DIR,"Model_Estimates_2018_Old_Sites_Cost_Addition.Rds"))
#Same as main results but only the reactors with a value ready to ship are included. That is to say other reactors may be willing to pay for future rights to hold the spent fuel but only those reactors with current values ready to ship are included, which set the bounds of the starting CIFS size. This is used for demonstration.
Rev_Shipping <- do.call(rbind,mclapply(CAPACITY_LIST,function(CAPACITY){do.call(rbind,lapply(DISCOUNT_RATE_LIST,function(x){YEARLY_RESULTS(TOTAL_VALUE_METRICS,x,CAPACITY,SHIPPING_COST_PER_TON )}))},mc.cores=min(length(CAPACITY_LIST),NCORES)))
Rev_Shipping <- CREATE_FULL_RESULTS(Rev_Shipping,COST_DATA)
saveRDS(Rev_Shipping,paste0(RES_DIR,"Model_Estimates_With_Shipping_Costs_Old_Sites_Cost_Addition.Rds"))
#Same as Rev_Shipping, but only information available in 2018 from the CURIE map is used (no updated operation dates).
Rev_Shipping_2018 <- do.call(rbind,mclapply(CAPACITY_LIST,function(CAPACITY){do.call(rbind,lapply(DISCOUNT_RATE_LIST,function(x){YEARLY_RESULTS(TOTAL_VALUE_METRICS_ORIG,x,CAPACITY,SHIPPING_COST_PER_TON )}))},mc.cores=min(length(CAPACITY_LIST),NCORES)))
Rev_Shipping_2018 <- CREATE_FULL_RESULTS(Rev_Shipping_2018,COST_DATA)
saveRDS(Rev_Shipping_2018,paste0(RES_DIR,"Model_Estimates_With_Shipping_Costs_2018_Old_Sites_Cost_Addition.Rds"))