R/e2e_run.R

In StrathE2E2: End-to-End Marine Food Web Model

#
# e2e_run.R
#
#' Run the StrathE2E model with a given configuration.
#'
#' Perform a single deterministic run of the StrathE2E model for a configuration
#' defined by an R-list object compiled by a prior call of the e2e_read() function.
#'
#' The function solves a network of Ordinary Differential Equations using the lsoda function
#' in the R deSolve package. The equations represent a shelf-sea food web and its connections to the physical and chemical environment,
#' and to a set of fishing fleets.
#'
#' The outputs from the run are time series (daily intervals) of the masses of each of the state
#' variables, the fluxes between all state variable, and the fluxes in and out of the model including
#' fishery landings. In addition, a range of derived quantities are generated for the final year of the
#' run, including annual averages, maxima and minima of state variables, annual fluxes between state variables, annual landings and
#' discards of each guild of taxa by each fleet of fishing gears, and a set of
#' network indices generated by the R NetIndices package.
#'
#' @param model R-list object defining the model configuration and compiled by the e2e_read() function.
#' @param nyears Number of years (integer) to run the model (default=20).
#' @param quiet Logical. If FALSE then see status outputs during the model run (default=TRUE).
#' @param csv.output Logical. If TRUE then enable writing of csv output files (default=FALSE).
#'
#' @return Model outputs as an R-list object and optionally csv files.
#'
#' @seealso \code{\link{e2e_ls}}, \code{\link{e2e_read}}, \code{\link{e2e_copy}}, \code{\link{e2e_plot_ts}}
#'
#' @importFrom deSolve ode
#'
#' @useDynLib StrathE2E2
#'
#' @export
#'
#' @examples
#' # Load the 2003-2013 version of the North Sea model supplied with the package:
#'     model <- e2e_read("North_Sea", "2003-2013")
#'
#' # Run the model for 2 years and generate the results object
#'     results <- e2e_run(model,nyears=2)
#' # Set csv.output=TRUE to save results to csv files to the directory
#' # specified in the e2e_read() function call.
#'
#' # Time series plot of state varaiables over the full length of the run
#'     e2e_plot_ts(model, results,selection="ECO")
#'     dev.new()
#'     e2e_plot_ts(model, results,selection="CATCH")
#
# ---------------------------------------------------------------------
# |                                                                   |
# | Authors: Mike Heath, Ian Thurlbeck                                |
# | Department of Mathematics and Statistics                          |
# | University of Strathclyde, Glasgow                                |
# |                                                                   |
# | Date of this version: May 2020                                    |
# |                                                                   |
# ---------------------------------------------------------------------

e2e_run <- function(model, nyears=20, quiet=TRUE, csv.output=FALSE) {

   oo <- options()
   on.exit(options(oo))

	pkg.env$quiet <- quiet			# quietens down read/write notes
	pkg.env$csv.output <- csv.output	# controls writing of CSV files

	# Run the Model:
	build <- build_model(model, nyears)
	output <- StrathE2E.run(build)

	# main processed output:
	aggregates			<- aggregate_model_output(model, output)
	total.annual.catch		<- extract_timeseries_annual_landings(model, build, output)
	annual.catch.by.gear		<- disaggregate_landings_discards_by_gear(elt(build, "fleet.output"), total.annual.catch)
	catch.land.disc			<- extract_simulated_catch_land_disc_by_gear_for_given_year(model, annual.catch.by.gear)

	# additional processed output:
	monthly.averages		<- monthly_averages_of_final_year(model, build, output, aggregates)
	annual.results.wholedomain	<- derive_annual_results_wholedomain(model, build, output, aggregates)
	annual.results.offshore		<- derive_annual_results_offshore(model, build, output, aggregates)
	annual.results.inshore		<- derive_annual_results_inshore(model, build, output, aggregates)
	flow.matrices			<- assemble_flow_matrix_from_model_annual_output(model, build, output, aggregates)

	model.path			<- elt(model, "setup", "model.path")
	annual.target.data		<- read_annual_target_data(model.path)
	monthly.target.data		<- read_monthly_target_data(model.path)

	model.target.results		<- derive_model_target_results(model, build, output, aggregates, annual.target.data)
	fit.to.target.data		<- calculate_error_function(model, model.target.results)

	results <- list(
		build			= list(
			model.parameters		= elt(build, "model.parameters"),
			run				= elt(build, "run"),
			drivers				= elt(build, "drivers"),
			forcings			= elt(build, "forcings")
		),
		output			= output,
		aggregates		= aggregates,
		fleet.output		= elt(build, "fleet.output"),
		total.annual.catch	= total.annual.catch,
		annual.catch.by.gear	= annual.catch.by.gear,
		final.year.outputs	= list(
			# catch/landings/discards:
			inshore_catchmat		= elt(catch.land.disc, "inshore_catchmat"),
			inshore_discmat			= elt(catch.land.disc, "inshore_discmat"),
			inshore_landmat			= elt(catch.land.disc, "inshore_landmat"),
			offshore_catchmat		= elt(catch.land.disc, "offshore_catchmat"),
			offshore_landmat		= elt(catch.land.disc, "offshore_landmat"),
			offshore_discmat		= elt(catch.land.disc, "offshore_discmat"),

			monthly.averages		= monthly.averages,

			mass_results_inshore		= elt(annual.results.inshore, "mass_results"),
			maxmass_results_inshore		= elt(annual.results.inshore, "maxmass_results"),
			minmass_results_inshore		= elt(annual.results.inshore, "minmass_results"),
			annual_flux_results_inshore	= elt(annual.results.inshore, "annual_flux_results"),

			mass_results_offshore		= elt(annual.results.offshore, "mass_results"),
			maxmass_results_offshore	= elt(annual.results.offshore, "maxmass_results"),
			minmass_results_offshore	= elt(annual.results.offshore, "minmass_results"),
			annual_flux_results_offshore	= elt(annual.results.offshore, "annual_flux_results"),

			mass_results_wholedomain	= elt(annual.results.wholedomain, "mass_results"),
			maxmass_results_wholedomain	= elt(annual.results.wholedomain, "maxmass_results"),
			minmass_results_wholedomain	= elt(annual.results.wholedomain, "minmass_results"),
			annual_flux_results_wholedomain	= elt(annual.results.wholedomain, "annual_flux_results"),

			# Network fluxes:
			flow_matrix_all_fluxes		= elt(flow.matrices, "flow_matrix_all_fluxes"),
			flow_matrix_excl_spawn_recruit	= elt(flow.matrices, "flow_matrix_excl_spawn_recruit"),
			NetworkIndexResults		= elt(flow.matrices, "NetworkIndexResults"),

			annual.target.data		= annual.target.data,
			monthly.target.data		= monthly.target.data,

			# error function/likelihoods:
			annual_obj			= elt(fit.to.target.data, "annual_obj"),
			partial_chi			= elt(fit.to.target.data, "partial_chi"),
			opt_results			= elt(fit.to.target.data, "opt_results")
		)
	)

	results
}