Can Epluspar be used for energy vs cost-optimazation?

[mdakyen]

Hi @hongyuanjia,
I am trying to perform energy vs cost-optimization using Epluspar.
Considering the limitations of energyplus in cost analysis, I decided to have my cost variables assigned in an R data frame. My first objective to minimize energy seems ok because it is dependent on Idf file. However, my second objective cost minimization kept prompting the following error:
Error: Objective function 'global_cost' must have an parameter named 'idf'.
In addition: Warning message:
In formals(fun) : argument is not a function

Please I need suggestions to overcome this error.

[hongyuanjia]

Hi @mdakyen . Sorry for the late response. I saw you close the issue. Have you already solved the problem?

[mdakyen]

@hongyuanjia, its OK Sir. Am sure you have been so busy.
I was able to get the program running with cost as my second objective. But it still takes a lot of time to run, because inside my cost function :
eplusr package collects the simulated energy value from EnergyPlus (for energy cost evaluation, since I can only pass an idf to the function)
openxlsx loads the xlsx cost data I have prepared ( which contains the investment cost, annual cost, maintenance and residual value)
tidyverse arranges values as they are computed.
Sir, I wish your package could be a little flexible, by allowing other variables to be passed to objectives functions. So they can be regarded as Global variables.

[hongyuanjia]

I see. Do you know what is the part that costs most of the time?
It would be great if you could provide some reproducible examples or at least some code snippets to let me have an idea on the direction to improve the function flexibility.

[mdakyen]

@hongyuanjia below is my code.

#load packages
library(eplusr)
library(tidyverse)
library(epluspar)
library(openxlsx)

#-------------------------------------------------------------------------------

# read case file
idf <- read_idf(system.file("extdata/1D 2F stepwisemodel04.idf", package = "eplusr"))
# read base case file paths
path_idf <- system.file("extdata/1D 2F stepwisemodel04.idf", package = "eplusr")
path_epw <- file.path(eplus_config(8.9)$dir, "WeatherData/Nicosia_hour.epw")

#-------------------------------------------------------------------------------

#Optimization section
# Create a GA optimization Job for wall insulation
    ga <- gaoptim_job(path_idf, path_epw)
    
    #define a measure to change the insulation thickness of the wall
    set_insulation <- function (idf, insul_tick){
      # set certain fields in different classes 
      idf$set(.("polystyrene eps") := list(Thickness = insul_tick)) 
      idf
      
    }
    
    # combine all measures into one
    design_options <- function (idf,insulation_thickness){
      idf <- set_insulation (idf, insulation_thickness)
      idf
    }
    
    # specify design space of parameters 
    ga$apply_measure(design_options,insulation_thickness = float_space(0,0.05))
    
    primary_energy <- function (idf) {
      
    primary_energy_per_floor <- as.double((idf$last_job()$tabular_data( 
                  report_name = "AnnualBuildingUtilityPerformanceSummary",
                  table_name = "Site and Source Energy",
                  row_name = "Total Site Energy",
                  column_name = "Total Energy")
                  $value))
    #converting annually consumed energy into primary energy then diving it by
      # building floor area
     primary_energy_per_floor <- ((primary_energy_per_floor*2.6)/190.55)
     primary_energy_per_floor
    
    }

    global_cost <- function(idf){
      
      #Local definition of consumed energy kWh
      
      hold_consumed_energy <- as.double(idf$last_job()$tabular_data( 
      report_name = "AnnualBuildingUtilityPerformanceSummary",
      table_name = "Site and Source Energy",
      row_name = "Total Site Energy",
      column_name = "Total Energy")
      $value)
      
      #Load cost item from IDF file        
      cost_item <-as.double (idf$object("polystyrene")$Cost_per_Each)
     
      #Insulation maximum parameters, note should be set by user
      max_insluation_cost <- cost_item
      height_insulation <- 1.25
      width_insulation <- 0.6
      thickness_insulation <- 0.05
      #Insulation present parameters
      present_thickness_insulation <- idf$object("polystyrene eps")$Thickness
      
      #Insulation manipulation
      volume_max_insluation = height_insulation*width_insulation*
                              thickness_insulation
     
      volume_present_insluation = height_insulation*width_insulation*
                                  present_thickness_insulation
      
      
      cost_per_increase_in_thickness = ((max_insluation_cost*volume_present_insluation)/
                                          (volume_max_insluation))
      # Read the Discount factors excel worksheet
      discount_factors <- openxlsx::read.xlsx("EPBDcostsheet06.xlsx", 
                                    sheet = 1)
      # Read the energy price over 30 years excel worksheet
      energy_price <- openxlsx::read.xlsx("EPBDcostsheet06.xlsx",
                                sheet = 2)
      
      # Evaluate present values for energy, note inflation is included    
      PV_energy <- dplyr::mutate(energy_price,
                                # add new variable for yearly consumed_electricity kWh
                                consumed_electrcity_kWh = hold_consumed_energy,
                                # add new variable for yearly cost value of energy (YCE)
                                YCE_TL = energy_price_inflation*hold_consumed_energy,
                                # add new variable for discount factors
                                discount_factor=discount_factors$discount_factor,
                                # add new variable for energy present value (EPV)
                                energy_present_value_TL = (YCE_TL*discount_factors$discount_factor))  
       
      #Evaluation of the NPV for energy alone = (sum years 1 to 30 )-(day 0) 
      NPV_energy <- PV_energy$energy_present_value_TL
      hold_sum <- sum(NPV_energy[2:31])
      hold_sub <- PV_energy$energy_present_value_TL[1]
      NPV_energy <- (hold_sum-hold_sub)
      NPV_energy
      
      # Read the measures from excel worksheet
      data_base <- openxlsx::read.xlsx("EPBDcostsheet06.xlsx", 
                                        sheet = 3)
      data_base <- dplyr::mutate(data_base,
                      #Purchase cost
                      purchase_cost_roof= 0.00,purchase_cost_wall = 0.00,purchase_cost_floor = 0.00,
                      #Installation cost 
                      installation_cost_roof = 0.00,installation_cost_wall = 0.00,installation_cost_floor = 0.00, 
                      #Initial investment
                      initial_investment_roof = 0.00,initial_investment_wall = 0.00, initial_investment_floor = 0.00,  
                      #Global cost of measures 
                      global_cost_roof= 0.00,global_cost_wall= 0.00,global_cost_floor= 0.00) 
      
#Evaluation of Global Cost = (cost_per_increase_in_thickness*building area)+
      #                        (installation_cost*building area)

#Purchase cost
data_base$purchase_cost_roof =(cost_per_increase_in_thickness*data_base$roof_area)
data_base$purchase_cost_wall =(cost_per_increase_in_thickness*data_base$wall_area)
data_base$purchase_cost_floor =(cost_per_increase_in_thickness*data_base$floor_area)
#Installation cost 
data_base$installation_cost_roof = (data_base$installation_cost*data_base$roof_area)
data_base$installation_cost_wall = (data_base$installation_cost*data_base$wall_area)
data_base$installation_cost_floor = (data_base$installation_cost*data_base$floor_area)
#Initial investment
data_base$initial_investment_roof = (data_base$purchase_cost_roof) + 
                                    (data_base$installation_cost_roof)
data_base$initial_investment_wall = (data_base$purchase_cost_wall) + 
                                    (data_base$installation_cost_wall)
data_base$initial_investment_floor = (data_base$purchase_cost_floor) + 
                                     (data_base$installation_cost_floor)
#Net energy present value
net_energy_present_vaule <- NPV_energy
               
data_base$global_cost_roof= data_base$initial_investment_roof +((net_energy_present_vaule+ 
          data_base$maintenance_cost+data_base$replacement_cost)-data_base$residual_value)
          data_base$global_cost_roof= data_base$global_cost_roof/190.55
          
data_base$global_cost_wall= data_base$initial_investment_wall +((net_energy_present_vaule+ 
          data_base$maintenance_cost+data_base$replacement_cost)-data_base$residual_value)
          data_base$global_cost_wall = data_base$global_cost_wall/190.55
          
data_base$global_cost_floor= data_base$initial_investment_floor +((net_energy_present_vaule+ 
          data_base$maintenance_cost+data_base$replacement_cost)-data_base$residual_value)
          data_base$global_cost_floor = data_base$global_cost_floor/190.55
        
      
data_base$global_cost_wall

      }
    
   
    # set optimization objectives 
    ga$objective(primary_energy,global_cost, .dir = "min")
    
    # specify how to mix solutions 
    ga$recombinator()
    # specify how to change parts of one solution randomly 
    ga$mutator()
    # specify how to select best solutions 
    ga$selector()
    # specify the conditions when to terminate the computation 
    ga$terminator(max_gen = 50L)
    
    # run optimization 
    ga$run(mu = 4, p_recomb = 0.7, p_mut = 0.2)
    
    # get all population 
    population <- ga$population()
    # get Pareto set 
    pareto <- ga$pareto_set()
    
    # plot Pareto front
    p_pareto <- ggplot() +
      geom_point(aes(primary_energy, global_cost), population, color = "darkgoldenrod", alpha = 0.5) +
      geom_line(aes(primary_energy, global_cost), pareto, color = "darkblue", linetype = 2) +
      geom_point(aes(primary_energy, global_cost), pareto, color = "darkblue", size = 2) +
      scale_x_continuous("Primary energy kWh/yr m?", labels = scales::number_format(scale = 0.001)) +
      scale_y_continuous("Global cost TL/ m?",labels = scales::number_format(big.mark = ","))
    
    p_pareto

[hongyuanjia]

@mdakyen Thanks for the code. It helps me to understand the issue a lot.

Here are some suggestions:

• openxlsx loads the xlsx cost data I have prepared ( which contains the investment cost, annual cost, maintenance and residual value)

  You can put the reading data process outside of your global_cost() function to make them global, and directly use them inside global_cost()

  # Read the Discount factors excel worksheet
  discount_factors <- openxlsx::read.xlsx("EPBDcostsheet06.xlsx", sheet = 1)
  # Read the energy price over 30 years excel worksheet
  energy_price <- openxlsx::read.xlsx("EPBDcostsheet06.xlsx", sheet = 2)
  # Read the measures from excel worksheet
  data_base <- openxlsx::read.xlsx("EPBDcostsheet06.xlsx", sheet = 3)
  
  height_insulation <- 1.25
  width_insulation <- 0.6
  thickness_insulation <- 0.05
  
  global_cost <- function(idf) {
    ......
    # Evaluate present values for energy, note inflation is included    
    PV_energy <- dplyr::mutate(energy_price,
                              # add new variable for yearly consumed_electricity kWh
                              consumed_electrcity_kWh = hold_consumed_energy,
                              # add new variable for yearly cost value of energy (YCE)
                              YCE_TL = energy_price_inflation*hold_consumed_energy,
                              # add new variable for discount factors
                              discount_factor=discount_factors$discount_factor,
                              # add new variable for energy present value (EPV)
                              energy_present_value_TL = (YCE_TL*discount_factors$discount_factor))
    ......
  }

• In global_cost() you can give another argument named param, like global_cost <- function(idf, param) {...}. param gives you values of all the optimization parameter (i.e. those you specified in ga$apply_measure()) for current individual in a list. So in your case, instead of

  present_thickness_insulation <- idf$object("polystyrene eps")$Thickness

  You can use

  present_thickness_insulation <- param$insulation_thickness

[hongyuanjia]

I was able to get the program running with cost as my second objective. But it still takes a lot of time to run, because inside my cost function :

• eplusr package collects the simulated energy value from EnergyPlus (for energy cost evaluation, since I can only pass an idf to the function)
• openxlsx loads the xlsx cost data I have prepared ( which contains the investment cost, annual cost, maintenance and residual value)
• tidyverse arranges values as they are computed.

I believe all those processes will not take too much time. The most time-consuming part is still EnergyPlus simulation itself.

[mdakyen]

@hongyuanjia well explained sir.
Thanks a lot for your recommendations. I will implement your suggestions ASAP. Then keep you posted.
Thanks once more.

[mdakyen]

@hongyuanjia good day.
I tried to implement the recommendations you made, but I was faced with the following error:
Error in dplyr::mutate(energy_price, consumed_electrcity_kWh = hold_consumed_energy, :
object 'energy_price' not found
The error is caused by reading data process outside of the global_cost() function, which is not globally recognized.

[hongyuanjia]

@mdakyen Thanks for the report. I will test on my side to locate the bug.

BTW, I did not come out a good way about how to add an extra parameter in the GAOptim$objective() to let the user specify extra arguments to the objective functions. If you have any suggestions, please let me know.

[mdakyen]

@hongyuanjia Ok will think about possible suggestions towards specifying extra arguments to the objective function. For now, I have been able to minimize two objectives (energy consumption and cost) simultaneously. However, am thinking of adding a third objective ( either carbon emission or thermal comfort) as I develop my model.
I will keep you posted once I get to that phase. Please let me know if your able to locate the bug that limits the definition of a variable as global.