diff --git a/.gitignore b/.gitignore index cc169bdf963c348d6dcc099a1c158c074ff831a8..0d1468ba738bc06bc144e8aba2898a0788e31c2e 100644 --- a/.gitignore +++ b/.gitignore @@ -13,7 +13,7 @@ # Ignore everything in "input" folders input/ - +!modules/29_CES_parameters/load/input # Ignore everything in "output" folders output/ diff --git a/config/default.cfg b/config/default.cfg index 9aabbef5273c14822d814e395384cba55d11c3c3..76cbb60d3cd710c0f4f692be01bd20693d9dd738 100755 --- a/config/default.cfg +++ b/config/default.cfg @@ -22,7 +22,7 @@ cfg$title <- "default" cfg$regionmapping <- "config/regionmappingH12.csv" #### Current input data revision (<mainrevision>.<subrevision>) #### -cfg$revision <- 5.945 +cfg$revision <- 5.946 #### Force the model to download new input data #### cfg$force_download <- FALSE diff --git a/modules/30_biomass/magpie_40/bounds.gms b/modules/30_biomass/magpie_40/bounds.gms index 1da85232e1bba1c106c7e0b7a83f08e7a9c0dc18..27579ea13b4efe91e5a5c0b8b14604a1aefb8f52 100644 --- a/modules/30_biomass/magpie_40/bounds.gms +++ b/modules/30_biomass/magpie_40/bounds.gms @@ -22,7 +22,7 @@ $ifthen.edge_esm_transport "%transport%" == "edge_esm" *** Slightly relaxed extraction bounds for biofuels. vm_fuExtr.up(t,regi,"pebios","5")$(t.val ge 2045) = 1.4*p30_datapebio(regi,"pebios","5","maxprod",t); vm_fuExtr.up(t,regi,"pebios","5")$(t.val ge 2055) = p30_datapebio(regi,"pebios","5","maxprod",t); -vm_fuExtr.up(t,regi,"pebioil","5")$(t.val ge 2030) = 1.5*p30_datapebio(regi,"pebioil","5","maxprod",t); +vm_fuExtr.up(t,regi,"pebioil","5")$(t.val ge 2030) = 2*p30_datapebio(regi,"pebioil","5","maxprod",t); vm_fuExtr.up(t,regi,"pebioil","5")$(t.val ge 2050) = p30_datapebio(regi,"pebioil","5","maxprod",t); $endif.edge_esm_transport diff --git a/scripts/output/comparison/EDGEcomparison.R b/scripts/output/comparison/EDGEcomparison.R index deb5cc426e311f2add8972a7a3deeb29819a7650..a93d0376f4982741a6059dfc382c0da71cb52446 100644 --- a/scripts/output/comparison/EDGEcomparison.R +++ b/scripts/output/comparison/EDGEcomparison.R @@ -48,6 +48,7 @@ EJfuelsFrgt_all = NULL emipSource_all = NULL costs_all = NULL pref_FV_all = NULL +demgdpcap_all = NULL scenNames <- getScenNames(outputdirs) EDGEdata_path <- path(outputdirs, paste("EDGE-T/")) @@ -488,6 +489,25 @@ costscompFun = function(newcomp, sharesVS1, EF_shares, pref_FV, capcost4Wall, c } +demgdpcap_Fun = function(demkm, REMIND2ISO_MAPPING) { +GDP_POP = getRMNDGDPcap() +GDP_POP = merge(GDP_POP, REMIND2ISO_MAPPING, by = "iso") +GDP_POP = merge(GDP_POP, REMIND2ISO_MAPPING, by = "iso") + +demcap_gdp = merge(demkm, GDP_POP, by = c("iso", "year")) +demcap_gdp = merge(demcap_gdp, REMIND2ISO_MAPPING, by = "iso") +demcap_gdp = demcap_gdp[,.(dem = sum(demand_F), gdp = sum(weight), pop = sum(POP_val)), by = .(region, sector, year)] +demcap_gdp[, GDP_cap := gdp/pop] +demcap_gdp[, demcap := dem* ## in trillion km + 1e+6/ ## in million km + pop] ## in million km/million people=pkm/person + +return(demcap_gdp) +} + + + + for (outputdir in outputdirs) { ## load mif file name_mif = list.files(path = outputdir, pattern = "REMIND_generic", full.names = F) @@ -562,6 +582,8 @@ for (outputdir in outputdirs) { emipSource = emipSourceFun(miffile) ## calculate costs by component costs = costscompFun(newcomp = newcomp, sharesVS1 = sharesVS1, EF_shares = EF_shares, pref_FV = pref_FV, capcost4Wall = capcost4Wall, capcost4W_BEVFCEV = capcost4W_BEVFCEV, nonf = nonf, totp = totp, REMIND2ISO_MAPPING) + ## per capita demand-gdp per capita + demgdpcap = demgdpcap_Fun(demkm = demandkm, REMIND2ISO_MAPPING) ## add scenario dimension to the results fleet[, scenario := as.character(unique(miffile$scenario))] salescomp[, scenario := unique(miffile$scenario)] @@ -576,6 +598,7 @@ for (outputdir in outputdirs) { emipSource[, scenario := as.character(unique(miffile$scenario))] costs[, scenario := as.character(unique(miffile$scenario))] pref_FV[, scenario := as.character(unique(miffile$scenario))] + demgdpcap[, scenario := as.character(unique(miffile$scenario))] ## rbind scenarios salescomp_all = rbind(salescomp_all, salescomp) fleet_all = rbind(fleet_all, fleet) @@ -590,6 +613,7 @@ for (outputdir in outputdirs) { emipSource_all = rbind(emipSource_all, emipSource) costs_all = rbind(costs_all, costs) pref_FV_all = rbind(pref_FV_all, pref_FV) + demgdpcap_all = rbind(demgdpcap_all, demgdpcap) } ## create string with date and time @@ -614,6 +638,7 @@ saveRDS(EJfuelsFrgt_all, paste0(outdir, "/EJfuelsFrgt_all.RDS")) saveRDS(emipSource_all, paste0(outdir, "/emipSource_all.RDS")) saveRDS(costs_all, paste0(outdir, "/costs_all.RDS")) saveRDS(pref_FV_all, paste0(outdir, "/pref_FV_all.RDS")) +saveRDS(demgdpcap_all, paste0(outdir, "/demgdpcap_all.RDS")) file.copy(file.path("./scripts/output/comparison/notebook_templates", md_template), outdir) rmarkdown::render(path(outdir, md_template), output_format="pdf_document") diff --git a/scripts/output/comparison/notebook_templates/EDGETransportComparison.Rmd b/scripts/output/comparison/notebook_templates/EDGETransportComparison.Rmd index 573697440d1542623bd0555799ccb42062a4cda9..dd12a05f398fe3f64aa1e62aea88bea71d8ffb5f 100644 --- a/scripts/output/comparison/notebook_templates/EDGETransportComparison.Rmd +++ b/scripts/output/comparison/notebook_templates/EDGETransportComparison.Rmd @@ -35,7 +35,7 @@ EJfrgt_all = readRDS("EJfuelsFrgt_all.RDS") emidem_all = readRDS("emidem_all.RDS") costs_all = readRDS("costs_all.RDS") pref_FV_all = readRDS("pref_FV_all.RDS") - +dempkm_cap_all = readRDS("demgdpcap_all.RDS") setConfig(forcecache=T) cols <- c("NG" = "#d11141", @@ -463,7 +463,64 @@ emidem_pf = function(dt){ emidem_pf(emidem_all) ``` +## demand per capita VS gdp per capita +```{r, echo=FALSE, message=FALSE, warning=FALSE, fig.width=14, fig.height=12} +demcapgdpcap_pf = function(dt) { + +dt = dt[year >= 2005] +dt = dt[, year := as.character(year)] +ppass_sm = ggplot()+ + geom_line(data = dt[sector == "trn_pass"], aes(x = GDP_cap, y = demcap, color = region, group = region))+ + geom_point(data = dt[sector == "trn_pass" & year %in% c("2020", "2030", "2050", "2070", "2100")], aes(x = GDP_cap, y = demcap, shape = year, group = region, color = region))+ + theme_minimal()+ + facet_grid(scenario~sector)+ + theme(axis.text.x = element_text(angle = 90, hjust = 1), + axis.text = element_text(size=7), + title = element_text(size=8), + legend.text = element_text(size=8))+ + labs(y = "Per capita demand [km/capita]", x="GDP per capita [$/person]", title = "Passenger short-medium") + +pfreigt_sm = ggplot()+ + geom_line(data = dt[sector == "trn_freight"], aes(x = GDP_cap, y = demcap, color = region, group = region))+ + geom_point(data = dt[sector == "trn_freight" & year %in% c("2020", "2030", "2050", "2070", "2100")], aes(x = GDP_cap, y = demcap, shape = year, group = region, color = region))+ + theme_minimal()+ + facet_grid(scenario~sector)+ + theme(axis.text.x = element_text(angle = 90, hjust = 1), + axis.text = element_text(size=7), + title = element_text(size=8), + legend.text = element_text(size=8))+ + labs(y = "Per capita demand [km/capita]", x="GDP per capita [$/person]", title = "Freight short-medium") + + +ppass_lo = ggplot()+ + geom_line(data = dt[sector == "trn_aviation_intl"], aes(x = GDP_cap, y = demcap, color = region, group = region))+ + geom_point(data = dt[sector == "trn_aviation_intl" & year %in% c("2020", "2030", "2050", "2070", "2100")], aes(x = GDP_cap, y = demcap, shape = year, group = region, color = region))+ + theme_minimal()+ + facet_grid(scenario~sector)+ + theme(axis.text.x = element_text(angle = 90, hjust = 1), + axis.text = element_text(size=7), + title = element_text(size=8), + legend.text = element_text(size=8))+ + labs(y = "Per capita demand [km/capita]", x="GDP per capita [$/person]", title = "International aviation") + +pfreigt_lo = ggplot()+ + geom_line(data = dt[sector == "trn_shipping_intl"], aes(x = GDP_cap, y = demcap, color = region, group = region))+ + geom_point(data = dt[sector == "trn_shipping_intl" & year %in% c("2020", "2030", "2050", "2070", "2100")], aes(x = GDP_cap, y = demcap, shape = year, group = region, color = region))+ + theme_minimal()+ + facet_grid(scenario~sector)+ + theme(axis.text.x = element_text(angle = 90, hjust = 1), + axis.text = element_text(size=7), + title = element_text(size=8), + legend.text = element_text(size=8))+ + labs(y = "Per capita demand [km/capita]", x="GDP per capita [$/person]", title = "International shipping") + +p = list(pfreigt_lo = pfreigt_lo, ppass_lo = ppass_lo, pfreigt_sm = pfreigt_sm, ppass_sm = ppass_sm) + +return(p) +} +demcapgdpcap_pf(dempkm_cap_all) +``` ## Focus on slected region @@ -544,9 +601,9 @@ costspf = function(dt, rp){ dt_inc = dt[type == "inc"] dt_inc = dt_inc[,.(cost = sum(cost)), by = c("region", "year", "type", "technology", "scenario")] dt_tot = rbind(dt_inc[, c("name", "logit_type") := list("Inconvenience cost", "pinc")], dt[type == "real"]) - dt[, technology := factor(technology, levels = legend_ord)] - dt_tot[, technology := factor(technology, levels = legend_ord)] - + dt[, name := factor(name, levels = legend_ord)] + dt_tot[, name := factor(name, levels = legend_ord)] + dt_tot[, alph := ifelse(logit_type=="pinc", 0.6, 1)] plot1 = ggplot()+ geom_bar(data = dt[type == "inc"], aes(x = year, y = cost, group = name, fill = name), position = "stack", stat = "identity")+ theme_minimal()+ @@ -564,7 +621,7 @@ facet_grid(scenario~technology)+ labs(y = "Costs [$/pkm]", x="") plot2 = ggplot()+ - geom_bar(data = dt_tot, aes(x = year, y = cost, group = name, fill = name), position = "stack", stat = "identity")+ + geom_bar(data = dt_tot, aes(x = year, y = cost, group = name, fill = name, alpha = alph), position = "stack", stat = "identity")+ theme_minimal()+ facet_grid(scenario~technology)+ theme(axis.text.x = element_text(angle = 90, size=14, vjust=0.5, hjust=1), @@ -575,9 +632,11 @@ facet_grid(scenario~technology)+ legend.text = element_text(size=14), strip.text = element_text(size=12), strip.background = element_rect(color = "grey"))+ - guides(fill=guide_legend(title="Cost component"))+ scale_fill_manual(values = cols)+ - labs(y = "Costs [$/pkm]", x="") + labs(y = "Costs [$/pkm]", x="")+ + scale_alpha(range=c(0.4,1)) + + guides(alpha=FALSE, linetype=FALSE, + fill=guide_legend(reverse=TRUE, title="Cost component")) plot = list(plot1 = plot1, plot2 = plot2) @@ -782,12 +841,12 @@ ggsave("pTruck.png", pTruck, dpi=500, height = height , width = height * aspect_ ```{r, echo=FALSE, warning=FALSE} prefBusTrucksplotf = function(pref){ -pref = pref[iso=="DEU" & technology %in% c("Electric", "FCEV") & vehicle_type %in% c("Bus_tmp_vehicletype")] -pref = pref[year >= 2020 & year <= 2100 & vehicle_type == "Bus_tmp_vehicletype"] -pref[, vehicle_type := "Buses and Trucks"] +pref = pref[iso=="USA" & technology %in% c("Electric", "FCEV") & vehicle_type %in% c("Bus_tmp_vehicletype", "Truck (0-2.7t)")] +pref = pref[year >= 2020 & year <= 2050] +pref[, vehicle_type := ifelse(vehicle_type == "Bus_tmp_vehicletype", "Large Trucks and Buses", "Small Trucks")] p = ggplot()+ - geom_line(data = pref, aes(x = year, y = value, group = technology, color = technology, linetype = technology))+ + geom_line(data = pref, aes(x = year, y = value, group = interaction(technology, vehicle_type), color = technology, linetype = vehicle_type))+ facet_wrap(~scenario, ncol = 1)+ theme_minimal()+ theme(axis.text.x = element_text(angle = 90, hjust = 1), @@ -796,9 +855,9 @@ p = ggplot()+ legend.text = element_text(size=7), strip.text = element_text(size=7))+ labs(y = "Preference factor [-]", x="")+ - scale_color_manual("Technology", values = cols,labels = c("Electric", "FCEV")) + - scale_linetype_manual("Technology", values = c(1,2), - labels = c("Electric", "FCEV")) + scale_color_manual("Technology", values = cols,labels = c("Electric", "FCEV")) + scale_linetype_manual("Vehicle", values = c(1,2), + labels = c("Large Trucks and Buses", "Small Trucks")) return(p) } p = prefBusTrucksplotf(pref_FV_all) diff --git a/scripts/output/comparison/policyCosts.R b/scripts/output/comparison/policyCosts.R index 6835361d2521b4928dbb3367f783eb27d079de45..7247fdffceee2d38290260305913768cb8308806 100644 --- a/scripts/output/comparison/policyCosts.R +++ b/scripts/output/comparison/policyCosts.R @@ -20,7 +20,8 @@ suppressPackageStartupMessages(library(tidyverse)) -# Function defintions +########################################################################### +# ###### START FUNCTION DEFINITONS ######################################## rm_timestamp <- function(strings, name_timestamp_seperator = "_", timestamp_format = "%Y-%m-%d_%H.%M.%S") { @@ -36,6 +37,7 @@ rm_timestamp <- function(strings, return(my_strings_wo_timeStamp) } + policy_costs_pdf <- function(policy_costs, fileName="PolicyCost.pdf") { @@ -113,13 +115,14 @@ policy_costs_pdf <- function(policy_costs, } + write_new_reporting <- function(mif_path, scen_name, new_polCost_data) { new_mif_path <- paste0(substr(mif_path,1,nchar(mif_path)-4),"_adjustedPolicyCosts.mif") - cat(paste0("A mif file with the name ",crayon::green(paste0("REMIND_generic_",scen_name,"_adjustedPolicyCosts.mif"))," is being created in the ",scen_name," outputfolder.\n")) + cat(paste0("A mif file with the name ",crayon::green(paste0("REMIND_generic_",scen_name,"_adjustedPolicyCosts.mif"))," is being created in the ",scen_name," output folder.\n")) my_data <- magclass::read.report(mif_path) my_variables <- grep("Policy Cost", magclass::getNames(my_data[[1]][[1]]), value = TRUE, invert = T) @@ -134,8 +137,78 @@ write_new_reporting <- function(mif_path, magclass::write.report(my_data, file = new_mif_path, ndigit = 7, skipempty = FALSE) - #magclass::write.report2(my_data, file=new_mif_path, ndigit=7, skipempty = FALSE) + return(new_mif_path) +} + + +report_transfers <- function(pol_mif, ref_mif) { + + # Read in reporting files + pol_run <- magclass::read.report(pol_mif,as.list = F) + ref_run <- magclass::read.report(ref_mif,as.list = F) + + # Get model and scenario names + md <- magclass::getItems(pol_run,3.2) + sc <- magclass::getItems(pol_run,3.1) + + # Tell the user what's going on + cat(paste0("Adding ",crayon::green("transfers")," to ",paste0("REMIND_generic_",sc,"_adjustedPolicyCosts.mif"),".\n")) + + + # Get gdploss + gdploss <- pol_run[,,"Policy Cost|GDP Loss (billion US$2005/yr)"] + # Add rel gdploss (not in percent) + gdploss_rel <- magclass::setNames(pol_run[,,"Policy Cost|GDP Loss|Relative to Reference GDP (percent)"]/100, + "Policy Cost|GDP Loss|Relative to Reference GDP") + # Get gdp + gdp_ref <- ref_run[,,"GDP|MER (billion US$2005/yr)"] + gdp_policy <- pol_run[,,"GDP|MER (billion US$2005/yr)"] + + # Calculate difference to global rel gdploss + delta_gdploss <- gdploss_rel[,,] - gdploss_rel["GLO",,] + # Calculate transfer required to equalize rel gdploss across regions + delta_transfer <- magclass::setNames(delta_gdploss * gdp_ref, + "Policy Cost|Transfers (billion US$2005/yr)") + delta_transfer_rel <- 100*magclass::setNames(delta_transfer/gdp_ref, + "Policy Cost|Transfers|Relative to Reference GDP (percent)") + + + # Calculate new gdp variables + gdp_withtransfers <- magclass::setNames(gdp_policy + delta_transfer, + "GDP|MER|w/ transfers (billion US$2005/yr)") + gdploss_withtransfers <- magclass::setNames(gdp_ref - gdp_withtransfers, + "Policy Cost|GDP Loss|w/ transfers (billion US$2005/yr)") + gdploss_withtransfers_rel <- 100*magclass::setNames(gdploss_withtransfers/gdp_ref, + "Policy Cost|GDP Loss|w/ transfers|Relative to Reference GDP (percent)") + + # Correct sets + magclass::getSets(delta_transfer, fulldim = F)[3] <- "variable" + magclass::getSets(delta_transfer_rel, fulldim = F)[3] <- "variable" + magclass::getSets(gdp_withtransfers, fulldim = F)[3] <- "variable" + magclass::getSets(gdploss_withtransfers, fulldim = F)[3] <- "variable" + magclass::getSets(gdploss_withtransfers_rel, fulldim = F)[3] <- "variable" + + # Bind together + my_transfers <- NULL + my_transfers <- magclass::mbind(my_transfers, delta_transfer) %>% + magclass::mbind(delta_transfer_rel) %>% + magclass::mbind(gdp_withtransfers) %>% + magclass::mbind(gdploss_withtransfers) %>% + magclass::mbind(gdploss_withtransfers_rel) + + pol_run <- magclass::read.report(pol_mif) + pol_run <- magclass::mbind(pol_run[[1]][[1]][,,], my_transfers) %>% + magclass::add_dimension(dim=3.1,add = "model",nm = md) %>% + magclass::add_dimension(dim=3.1,add = "scenario",nm = sc) + + magclass::write.report(pol_run, file = pol_mif, ndigit = 7, skipempty = FALSE) + + return(my_transfers) } +# ###### END FUNCTION DEFINITONS ######################################## +########################################################################### + + # Check for an object called "source_include". If found, that means, this script @@ -144,13 +217,13 @@ write_new_reporting <- function(mif_path, # default values, and made over-writable with the command line. if(!exists("source_include")) { # Set default value - outputdirs <- c("../../../output/default_2020-03-03_09.38.01", - "../../../output/default_slim_2020-03-03_11.37.51", - "../../../output/default_slim_2020-03-03_11.37.51", - "../../../output/default_2020-03-03_09.38.01") + outputdirs <- c("base_noEffChange_2020-03-09_17.16.28/", + "base_allT_lab_1point25_2020-03-27_16.12.35/", + "base_allT_lab_1point25_2020-03-27_16.12.35/", + "base_noEffChange_2020-03-09_17.16.28/") + special_requests <- c("2") # Make over-writtable from command line - lucode::readArgs("outputdirs") - + lucode::readArgs("outputdirs","special_requests") } @@ -192,10 +265,22 @@ while (!happy_with_input) { cat(crayon::green(paste0("\t", pc_pairs ,"\n"))) cat("Is that what you intended?\n") cat(paste0("Type '",crayon::green("y"),"' to continue, '",crayon::blue("r"),"' to reselect output directories, '",crayon::red("n"),"' to abort: ")) + user_input <- get_line() + if(user_input %in% c("y","Y","yes")) { happy_with_input <- TRUE cat(crayon::green("Great!\n")) + + # Get special requests from user + cat(crayon::blue("\nDo you have any special requests?\n")) + cat("1: Skip creation of adjustedPolicyCost reporting\n") + cat("2: Add transfers to adjustedPolicyCost reporting\n") + cat("3: Skip plot creation\n") + cat("4: Plot until 2150 in pdf\n") + cat("Type the number (or numbers seperated by a comma) to choose the special requests, or nothing to continue without any: ") + special_requests <- get_line() %>% str_split(",",simplify = T) %>% as.vector() + } else if(user_input %in% c("r","R","reselect")) { if (exists("choose_folder")) { cat("Remember, the order in which you choose the directories should be:\n") @@ -211,45 +296,71 @@ while (!happy_with_input) { cat(crayon::red("\nStopping execution now.\n\n")) stop("I can't figure this **** out. I give up. ") } + } else { + happy_with_input <- TRUE } } -# Tell the user what is going on -cat(crayon::blue("\nPolicy cost computations:\n")) - # Get Policy costs for every policy-reference pair +cat(crayon::blue("\nComputing Policy costs:\n")) tmp_policy_costs_magpie <- mapply(remind::reportPolicyCosts, pol_gdxs, ref_gdxs, SIMPLIFY = FALSE) +cat(crayon::green("Done!\n")) -tmp_policy_costs <- tmp_policy_costs_magpie %>% - lapply(quitte::as.quitte) %>% - lapply(select, region, period, data, value) - -# Combine results in single tibble, with names like "Pol_w.r.t_Ref" -policy_costs <- rename(tmp_policy_costs[[1]], !!sym(paste0(pol_names[1], "_w.r.t_",ref_names[1])):=value) -if (length(tmp_policy_costs)>1){ - for (i in 2:length(tmp_policy_costs)) { - policy_costs <- tmp_policy_costs[[i]] %>% - rename(!!sym(paste0(pol_names[i], "_w.r.t_",ref_names[i])):=value) %>% - left_join(policy_costs, tmp_policy_costs[[i]], by=c("region", "period", "data")) - } + +# Create "adjustedPolicyCost" reporting file +if (!"1" %in% special_requests) { + cat(crayon::blue("\nCreating new reportings:\n")) + pol_mifs <- paste0(dirname(pol_gdxs), "/REMIND_generic_", pol_names, ".mif") + new_reporting_files <- mapply(write_new_reporting, pol_mifs, pol_names, tmp_policy_costs_magpie) + cat(crayon::green("Done!\n")) } -# and do some pivotting -policy_costs <- policy_costs %>% - pivot_longer(cols = matches(".*w\\.r\\.t.*"), names_to = "Model Output") %>% - pivot_wider(names_from = data) -cat(crayon::green("Done!\n")) -# Create Pdf -cat(crayon::blue("\nPdf creation:\n")) -time_stamp <- format(Sys.time(), "_%Y-%m-%d_%H.%M.%S") -policy_costs_pdf(policy_costs, fileName = paste0("PolicyCost",time_stamp,".pdf")) -cat(crayon::green("Done!\n")) +# Add transfer variables to "adjustedPolicyCost" reporting file +if ("2" %in% special_requests && !"1" %in% special_requests) { + cat(crayon::blue("\nComputing transfers:\n")) + ref_mifs <- paste0(dirname(ref_gdxs), "/REMIND_generic_", ref_names, ".mif") + transfer_info <- mapply(report_transfers, new_reporting_files, ref_mifs, SIMPLIFY = F) + cat(crayon::green("Done!\n")) +} -# Create new reporting file -cat(crayon::blue("\nMif creation:\n")) -pol_mifs <- paste0(dirname(pol_gdxs), "/REMIND_generic_", pol_names, ".mif") -return_check <- mapply(write_new_reporting, pol_mifs, pol_names, tmp_policy_costs_magpie) -cat(crayon::green("Done!\n\n")) +# Create Pdf +if (!"3" %in% special_requests) { + cat(crayon::blue("\nCreating plots:\n")) + + # Add transfers, if they exist + if (exists("transfer_info")) { + tmp_policy_costs_magpie <- mapply(magclass::mbind, tmp_policy_costs_magpie, transfer_info, SIMPLIFY = F) + } + + tmp_policy_costs <- tmp_policy_costs_magpie %>% + lapply(quitte::as.quitte) %>% + lapply(select, region, period, data, value) + + # Combine results in single tibble, with names like "Pol_w.r.t_Ref" + policy_costs <- rename(tmp_policy_costs[[1]], !!sym(paste0(pol_names[1], "_w.r.t_",ref_names[1])):=value) + if (length(tmp_policy_costs)>1){ + for (i in 2:length(tmp_policy_costs)) { + policy_costs <- tmp_policy_costs[[i]] %>% + rename(!!sym(paste0(pol_names[i], "_w.r.t_",ref_names[i])):=value) %>% + left_join(policy_costs, tmp_policy_costs[[i]], by=c("region", "period", "data")) + } + } + # and do some pivotting + policy_costs <- policy_costs %>% + pivot_longer(cols = matches(".*w\\.r\\.t.*"), names_to = "Model Output") %>% + pivot_wider(names_from = data) + + # By default, plots are only created until 2100 + if (!"4" %in% special_requests) { + policy_costs <- policy_costs %>% filter(period<=2100) + } + + time_stamp <- format(Sys.time(), "_%Y-%m-%d_%H.%M.%S") + policy_costs_pdf(policy_costs, fileName = paste0("PolicyCost",time_stamp,".pdf")) + cat(crayon::green("Done!\n")) +} + +cat("\n") \ No newline at end of file diff --git a/scripts/start/prepare_and_run.R b/scripts/start/prepare_and_run.R index 513e12c7c5672b4aecc9de7ba9712cf6d10c8d29..aaaae532cf8218659af6ebd57462063f0527a5b8 100644 --- a/scripts/start/prepare_and_run.R +++ b/scripts/start/prepare_and_run.R @@ -70,6 +70,18 @@ getReportData <- function(path_to_report,inputpath_mag="magpie",inputpath_acc="c map <- rbind(map,data.frame(emimag="Emissions|CH4|Land|Agriculture|+|Rice (Mt CH4/yr)", emirem="ch4rice", factor_mag2rem=1,stringsAsFactors=FALSE)) map <- rbind(map,data.frame(emimag="Emissions|CH4|Land|Agriculture|+|Animal waste management (Mt CH4/yr)", emirem="ch4anmlwst",factor_mag2rem=1,stringsAsFactors=FALSE)) map <- rbind(map,data.frame(emimag="Emissions|CH4|Land|Agriculture|+|Enteric fermentation (Mt CH4/yr)", emirem="ch4animals",factor_mag2rem=1,stringsAsFactors=FALSE)) + } else if("Emissions|N2O-N|Land|Agriculture|+|Animal Waste Management (Mt N2O-N/yr)" %in% getNames(mag)) { + # MAgPIE 4 new + map <- rbind(map,data.frame(emimag="Emissions|CO2|Land (Mt CO2/yr)", emirem="co2luc", factor_mag2rem=1/1000*12/44,stringsAsFactors=FALSE)) + map <- rbind(map,data.frame(emimag="Emissions|N2O-N|Land|Agriculture|+|Animal Waste Management (Mt N2O-N/yr)", emirem="n2oanwstm", factor_mag2rem=28/44,stringsAsFactors=FALSE)) + map <- rbind(map,data.frame(emimag="Emissions|N2O-N|Land|Agriculture|Agricultural Soils|+|Inorganic Fertilizers (Mt N2O-N/yr)", emirem="n2ofertin", factor_mag2rem=28/44,stringsAsFactors=FALSE)) + map <- rbind(map,data.frame(emimag="Emissions|N2O-N|Land|Agriculture|Agricultural Soils|+|Manure applied to Croplands (Mt N2O-N/yr)",emirem="n2oanwstc", factor_mag2rem=28/44,stringsAsFactors=FALSE)) + map <- rbind(map,data.frame(emimag="Emissions|N2O-N|Land|Agriculture|Agricultural Soils|+|Decay of Crop Residues (Mt N2O-N/yr)", emirem="n2ofertcr", factor_mag2rem=28/44,stringsAsFactors=FALSE)) + map <- rbind(map,data.frame(emimag="Emissions|N2O-N|Land|Agriculture|Agricultural Soils|+|Soil Organic Matter Loss (Mt N2O-N/yr)", emirem="n2ofertsom",factor_mag2rem=28/44,stringsAsFactors=FALSE)) + map <- rbind(map,data.frame(emimag="Emissions|N2O-N|Land|Agriculture|Agricultural Soils|+|Pasture (Mt N2O-N/yr)", emirem="n2oanwstp", factor_mag2rem=28/44,stringsAsFactors=FALSE)) + map <- rbind(map,data.frame(emimag="Emissions|CH4|Land|Agriculture|+|Rice (Mt CH4/yr)", emirem="ch4rice", factor_mag2rem=1,stringsAsFactors=FALSE)) + map <- rbind(map,data.frame(emimag="Emissions|CH4|Land|Agriculture|+|Animal waste management (Mt CH4/yr)", emirem="ch4anmlwst",factor_mag2rem=1,stringsAsFactors=FALSE)) + map <- rbind(map,data.frame(emimag="Emissions|CH4|Land|Agriculture|+|Enteric fermentation (Mt CH4/yr)", emirem="ch4animals",factor_mag2rem=1,stringsAsFactors=FALSE)) } else { # MAgPIE 3 map <- rbind(map,data.frame(emimag="Emissions|CO2|Land Use (Mt CO2/yr)", emirem="co2luc", factor_mag2rem=1/1000*12/44,stringsAsFactors=FALSE)) diff --git a/start_bundle_coupled.R b/start_bundle_coupled.R index f29c0633e621264c4986979465fd0fb1f6d526ab..45a13749469da3b4d1dac19d8b1bdd035d958305 100644 --- a/start_bundle_coupled.R +++ b/start_bundle_coupled.R @@ -12,24 +12,26 @@ path_remind <- paste0(getwd(),"/") # provide path to REMIND. Default: the actual path which the script is started from path_magpie <- "/p/projects/piam/runs/coupled-magpie/" +# Paths to the files where scenarios are defined +# path_settings_remind contains the detailed configuration of the REMIND scenarios +# path_settings_coupled defines which runs will be started, coupling infos, and optinal gdx and report inforamtion that overrides path_settings_remind +path_settings_coupled <- paste0(path_remind,"config/scenario_config_coupled_SSPSDP.csv") +path_settings_remind <- paste0(path_remind,"config/scenario_config_SSPSDP.csv") + +# You can put a prefix infront of the names of your runs, this will turn e.g. "SSP2-Base" into "prefix_SSP2-Base". +# This allows storing results of multiple coupled runs (which have the same scenario names) in the same MAgPIE and REMIND output folders. +prefix_runname <- "C_" + # If there are existing runs you would like to take the gdxes (REMIND) or reportings (REMIND or MAgPIE) from provide the path here and the name prefix below. # Note: the sceanrio names of the old runs have to be identical to the runs that are to be started. If they differ please provide the names of the old scenarios in the -# file that you read in below to path_settings_coupled (scenario_config_coupled_xxx.csv). +# file that you specified on path_settings_coupled (scenario_config_coupled_xxx.csv). path_remind_oldruns <- paste0(path_remind,"output/") path_magpie_oldruns <- paste0(path_magpie,"output/") -# The scripts automatically adds a prefix (name of your remind path) to the scenario names. This is useful because it enables -# using the same MAgPIE and REMIND output folders to store results of coupled runs from multiple REMIND revisions (prevents double names) # If you want the script to find gdxs or reports of older runs as starting point for new runs please # provide the prefix of the old run names so the script can find them. prefix_oldruns <- "C_" -# Paths to the files where scenarios are defined -# path_settings_remind contains the detailed configuration of the REMIND scenarios -# path_settings_coupled defines which runs will be started, coupling infos, and optinal gdx and report inforamtion that overrides path_settings_remind -path_settings_coupled <- paste0(path_remind,"config/scenario_config_coupled_SSPSDP.csv") -path_settings_remind <- paste0(path_remind,"config/scenario_config_SSPSDP.csv") - # number of coupling iterations max_iterations <- 5 @@ -95,9 +97,6 @@ if (!identical(common,character(0))) { for(scen in common){ cat(paste0("\n################################\nPreparing run ",scen,"\n")) - prefix_runname <- "C" #strsplit(path_remind,"/")[[1]][length(strsplit(path_remind,"/")[[1]])] - prefix_runname <- paste0(prefix_runname,"_") - runname <- paste0(prefix_runname,scen) # name of the run that is used for the folder names path_report <- NULL # sets the path to the report REMIND is started with in the first loop LU_pricing <- scenarios_coupled[scen, "LU_pricing"] # set the GHG prices to zero up to and including the year specified here