 /*-----------------------------------------------------------------------------------------
 LIBRARY: REACTOR
 COMPONENT: reactorAB
 PARTITION: partition1
 EXPERIMENT: OptimSin
 TEMPLATE: TRANSIENT
 CREATION DATE: 21/12/2017
-----------------------------------------------------------------------------------------*/
USE OPTIM_METHODS


CONST INTEGER numC = 3				-- número de restricciones del problema de optimización
CONST INTEGER numU = 2				-- número de variables de decisión

FUNCTION INTEGER coste_y_restricciones (IN REAL esnopt_x[], IN INTEGER needF, OUT REAL esnopt_F[], IN INTEGER explicit_derivatives, IN INTEGER needG, OUT REAL esnopt_G[])
	DECLS


	BODY

	RESET_VARIABLES ()
	
      q  = esnopt_x[1]
		Fr = esnopt_x[2]
		
		SET_INIT_ACTIVE(TRUE)
	
		TIME = 0
		TSTOP = 10
		CINT = 1
		IMETHOD = IDAS	
		INTEG()		
		-- Introducción de los valores de las funcines a SNOPT
		--	si la funcion es de camino, hay que introducir la cuadratura ( arr_quad )
		-- en caso contrario, si la restricción es de punto final, al usar F_optim se usa el valor de la función.
		esnopt_F[1] = F_optim[1]
		esnopt_F[2] = F_optim[2]
		esnopt_F[3] = F_optim[3]
		esnopt_F[4] = F_optim[4]

		RETURN 0
END FUNCTION


EXPERIMENT OptimSin ON reactorAB.partition1
   DECLS
		REAL dec_var[numU]      	-- valor inicial de las variables de decision, size numU. Necesario para las llamadas a ESnopt
		
	  	REAL xlow[numU]				-- valor inferior de las variables de decision, size numU. Necesario para las llamadas a ESnopt
   	REAL xupp[numU]				-- valor superior de las variables de decision, size numU. Necesario para las llamadas a ESnopt
  
   	REAL Flow[numC + 1]			-- valor inferior de la funcion objetivo y las restricciones, size (numC + función de coste). Necesario para las llamadas a ESnopt
   	REAL Fupp[numC + 1]			-- valor superior de la funcion objetivo y las restricciones, size (numC + función de coste). Necesario para las llamadas a ESnopt

		INTEGER infoESnopt = 0
		
   OBJECTS
	

   INIT
      -- initial values for state variables
      cA = 2.891
      Tr = 51.5
		T = T0

   BOUNDS
      -- Set equations for boundaries: boundVar = f(TIME;...)
      Fr = 50
      q = 2

   BODY
      -- creates an ASCII file with the results in table format
      -- REPORT_TABLE("results.rpt", "*",TRUE)

      -- set the debug level (valid range [0,4])
      DEBUG_LEVEL = 1
      -- select default integration solver. Valid methods are IDAS (_SPARSE), DASSL(_SPARSE), CVODE_BDF(_SPARSE), CVODE_AM, RK4, EULER, AM1, AM2 and AM4
      IMETHOD = IDAS -- default is DASSL, recommended is either IDAS or IDAS_SPARSE  
      -- Settings for different actions. Valid actions are: SEV_DISABLE, SEV_NONE, SEV_WARNING, SEV_ERROR, SEV_KILLPOINT, SEV_FATAL
      eSetErrorAction(ERR_BAD_OPER,SEV_KILLPOINT) -- Detect bad numerical operations (eg division by zero), default is SEV_NONE, recommended is SEV_KILLPOINT
      eSetErrorAction(ERR_NAN_INF, SEV_KILLPOINT) -- Detect NaN or Inf values, default is SEV_NONE, recommended is SEV_KILLPOINT
      eSetConfig(CFG_FORCE_STOP_CINT,TRUE) -- Force or not to stop the solver each CINT. Sometimes you can speed up the simulation if you select FALSE
      -- set tolerances and other important inputs
      REL_ERROR = 1e-06 -- transient solver relative tolerance
      ABS_ERROR = 1e-06 -- transient solver absolute tolerance
      TOLERANCE = 1e-06 -- steady solver relative tolerance
      INIT_INTEG_STEP = -1 -- initial integration step size (-1 means use the solver estimation) 
      MAX_INTEG_STEP = -1 -- maximum integration step size (-1 means use the solver estimation) 
      NSTEPS = 1 -- Only for explicit solvers use CINT/NSTEPS as integration step size
--      REPORT_MODE = IS_EVENT -- by default it reports results at every CINT and event detection. Other valid options are IS_STEP, IS_CINT and IS_MANUAL_REFRESH 

		dec_var[1] = q
		dec_var[2] = Fr
		
		xlow[1] = 0.5
		xupp[1] = 5

		xlow[2] = 5
		xupp[2] = 90

		TIME = 0
      TSTOP = 15
      CINT = 0.1
		
	
		-- inicialización de los límites de las restricciones y la función de coste
		Flow[1] = -1.0e6
		Fupp[1] =  1.0e6			
		
		Flow[2] = Liminfx
		Fupp[2] = 1			
		
		Flow[3] = 0
		Fupp[3] = 0

		Flow[4] = 0
		Fupp[4] = 0
					
		-- Representacion en el monitor de los valores previos a la optimizacion
		WRITE ("\n\n\t\t------------------- limits before optimization \n")
		FOR ( i IN 1,numU )
				WRITE ("\t\t var_dec[%d]:\t %g \t %g \t %g \n", i, xlow[i], dec_var[i], xupp[i])
		END FOR
		FOR ( i IN 1,numC + 1 )
				WRITE ("\t\t F_optim[%d]:\t %g \t %g \t %g \n", i, Flow[i], F_optim[i], Fupp[i])		
		END FOR
		WRITE ("\n\t\t------------------- END of limits  \n\n")

   	--Optimization extern routine call
		setSilentMode(TRUE)
		SET_REPORT_ACTIVE("#MONITOR",FALSE)
		
		esnopt_init (numU, numC)
		esnopt_set_variables_bounds_and_initial_values (xlow, xupp, dec_var)
		esnopt_set_constraints_bounds_and_initial_values (Flow, Fupp, F_optim)
		esnopt_set_cost_function_and_constraints (coste_y_restricciones)
		esnopt_set_explicit_derivatives ( 0 )
		esnopt_set_function_precision ( 1.0e-5 )
		esnopt_set_iterations_limit (200) 
		infoESnopt = esnopt ()
		
		-- Final de la optimización, obtención de los resultados para la simulación.
		setSilentMode(FALSE)
		SET_REPORT_ACTIVE("#MONITOR",TRUE)
		
		esnopt_print_data ()
		esnopt_get_variables_values(dec_var)
		esnopt_free ()	
		
		-- Llamada al integrador
		RESET()
--			SET_INIT_ACTIVE(TRUE)

		-- Modificación de los parámetros con los nuevos valores obtenidos		
		q = dec_var[1]
		Fr = dec_var[2]

      TIME = 0
		CINT = 0.1

  	   INTEG()
		
		-- Representacion en el monitor de los valores posteriores a la optimizacion
		WRITE ("\n\n\t\t------------------- The optimiser succeed, returned the value: %d. \n", infoESnopt)
		WRITE ("\n\n\t\t------------------- constraints and decision variables after optimization \n")
		FOR ( i IN 1,numU )
			WRITE ("\t\t var_dec[%d]:\t %g \t %g \t %g \n", i, xlow[i], dec_var[i], xupp[i])
		END FOR
		WRITE ("\n\n\t\t------------------- cost function and constraints after optimization \n")
		FOR ( i IN 1,numC + 1 )
			WRITE ("\t\t F_optim[%d]:\t %g \t %g \t %g \n", i, Flow[i], F_optim[i], Fupp[i])		
		END FOR
		WRITE ("\n\t\t------------------- end   \n\n")
END EXPERIMENT
