input <- readLines("path/to/input")
input <- as.integer(input)
library(ompr)
library(ompr.roi)
library(ROI.plugin.glpk)
library(tidyverse)

n <- length(input)
MIPModel() %>% 
  # a variable for each input number
  add_variable(y[i], type = "binary", i = 1:n) %>% 
  set_objective(0) %>% 
  # we want to select exactly two numbers
  add_constraint(sum_expr(y[i], i = 1:n) == 2) %>%
  # the sum of the two selected numbers should be 2020
  add_constraint(sum_expr(input[i] * y[i], i = 1:n) == 2020) %>% 
  solve_model(with_ROI("glpk", verbose = FALSE)) %>% 
  get_solution(y[i]) %>% 
  filter(value == 1) %>%
  pull(i) -> indexes

sum(input[indexes])
#> [1] 2020
prod(input[indexes])
#> [1] 357504