PoAB/de2/elements.py

import numpy as np
from pypdevs.DEVS import AtomicDEVS, CoupledDEVS
from pypdevs.infinity import INFINITY

from dataclasses import dataclass, field
import random
import pandas as pd

from de2.routing import get_graph, get_closest_vertex


TUGS = pd.read_excel("20230405_Tugs.xlsx", dtype={"MMSI": str, "NAME": str, "category": str})


@dataclass
class Footprint:
	usage_time: float = 0.0
	idle_time: float = 0.0
	idle_since: float = 0.0

	sailing_time: float = 0.0
	tugging_time: float = 0.0

	fuel_usage: float = 0.0
	velocity_total: float = 0.0
	tasks: int = 0

	@staticmethod
	def get_fuel_efficiency(velocity):
		"""
		Computes how efficient fuel is being used at a certain velocity.

		Function based on the figure from Pieter (see figures/fuel.jpg) and an
		online curve fitter (https://mycurvefit.com/).
		"""
		return 32.89402 + (1.114886 - 32.89402) / (1 + (velocity / 12.20228) ** 6.449555)

	@staticmethod
	def get_fuel(distance, velocity):
		# TODO: find the right value here
		fuel = 5  # litres/meter for 100% efficient fuel usage
		efficiency = Footprint.get_fuel_efficiency(velocity)
		return distance * fuel * efficiency


@dataclass
class Vessel:
	mmsi: str
	velocity: float = 0.0
	source: tuple = None
	target: tuple = None
	task: str = None
	distance_left: float = 0.0
	total_distance: float = 0.0
	time_until_departure: float = 0.0

	footprint: Footprint = field(default_factory=Footprint)
	name: str = ""
	category: str = "weak"

	def tuple(self):
		return self.mmsi, self.name,\
		       self.source, self.target,\
		       self.task, self.total_distance, self.distance_left, self.velocity


class Pool(AtomicDEVS):
	def __init__(self, name):
		super(Pool, self).__init__(name)

		self.state = {
			"waiting": {},
			"should_exit": [],
			"time": 0.0
		}

		self.req_in = self.addInPort("req_in")
		self.vessel_in = self.addInPort("vessel_in")
		self.vessel_out = self.addOutPort("vessel_out")

		# prefill with all vessels of the simulation
		for mmsi in TUGS["MMSI"]:
			if isinstance(mmsi, str):
				self.state["waiting"][mmsi] = Vessel(mmsi)

	def extTransition(self, inputs):
		self.state["time"] += self.elapsed
		if self.req_in in inputs:
			for request in inputs[self.req_in]:
				if request["mmsi"] in self.state["waiting"]:
					vessel = self.state["waiting"][request["mmsi"]]
				else:
					print("VESSEL %s DOES NOT EXIST IN POOL" % str(request["mmsi"]))
					# print("\t", request)
					vessel = Vessel(request["mmsi"])

				vessel.footprint.idle_time += self.state["time"] - vessel.footprint.idle_since

				# vessel.time_until_departure = (request["end"] - request["start"]) / 1000
				# vessel.velocity = request["distance"] / vessel.time_until_departure
				# vessel.distance_left = request["distance"]
				# vessel.total_distance = request["distance"]
				# vessel.source = (request["source_lon"], request["source_lat"])
				# vessel.target = (request["target_lon"], request["target_lat"])
				# vessel.location = request["location"]

				vessel.velocity = request["velocity"]
				if vessel.velocity > 0: # requests for idling should be ignored!
					vessel.distance_left = request["distance"]
					vessel.total_distance = vessel.distance_left
					# vessel.time_until_departure = vessel.distance_left / vessel.velocity
					vessel.source = request["source_lon"], request["source_lat"]
					vessel.target = request["target_lon"], request["target_lat"]
					vessel.task = request["task"]

					vessel.name = request["name"]
					vessel.category = request["category"]

					self.state["should_exit"].append(vessel)
					if vessel.mmsi in self.state["waiting"]:
						del self.state["waiting"][vessel.mmsi]
		if self.vessel_in in inputs:
			for vessel in inputs[self.vessel_in]:
				if vessel.mmsi in self.state["waiting"]:
					raise ValueError("Cannot create duplicate vessels (loc = %s; mmsi = %s; d = %f)" %
					                 (str(vessel.location), str(vessel.mmsi), vessel.total_distance))
				else:
					vessel.source = vessel.target
					vessel.target = None
					if vessel.velocity > 0:
						time = vessel.total_distance / vessel.velocity
						# vessel.footprint.fuel_usage += Footprint.get_fuel(vessel.total_distance, vessel.velocity)
						vessel.footprint.usage_time += time
						if vessel.task == "sailing":
							vessel.footprint.sailing_time += time
						else:
							vessel.footprint.tugging_time += time
						vessel.footprint.velocity_total += vessel.velocity
						vessel.footprint.tasks += 1

					vessel.footprint.idle_since = self.state["time"]

					self.state["waiting"][vessel.mmsi] = vessel
		return self.state

	def timeAdvance(self):
		if len(self.state["should_exit"]) > 0:
			return 0.0
		return INFINITY

	def outputFnc(self):
		if len(self.state["should_exit"]) > 0:
			return {
				self.vessel_out: [self.state["should_exit"][0]]
			}
		return {}

	def intTransition(self):
		self.state["time"] += self.timeAdvance()
		if len(self.state["should_exit"]) > 0:
			self.state["should_exit"].pop(0)
		return self.state


class Sailer(AtomicDEVS):
	def __init__(self, name, v_tugging, v_sailing, v_max):
		super(Sailer, self).__init__(name)

		self.v_tugging = v_tugging
		self.v_sailing = v_sailing
		self.v_max = v_max

		self.state = {
			"vessels": [],
			"time": 0.0
		}

		self.vessel_in = self.addInPort("vessel_in")
		self.update = self.addInPort("update")
		self.vessel_out = self.addOutPort("vessel_out")

	def extTransition(self, inputs):
		self.state["time"] += self.elapsed
		self.update_vessels(self.elapsed)
		if self.vessel_in in inputs:
			for vessel in inputs[self.vessel_in]:
				if vessel.task == "tugging" and self.v_tugging is not None:
					vessel.velocity = -1
					while vessel.velocity <= 0 or vessel.velocity > self.v_max:
						vessel.velocity = random.normalvariate(*self.v_tugging)
				elif vessel.task == "sailing" and self.v_sailing is not None:
					vessel.velocity = -1
					while vessel.velocity <= 0 or vessel.velocity > self.v_max:
						vessel.velocity = random.normalvariate(*self.v_sailing)
				vessel.time_until_departure = vessel.distance_left / vessel.velocity
				self.state["vessels"].append(vessel)
		self.state["vessels"].sort(key=lambda v: v.time_until_departure)
		return self.state

	def update_vessels(self, elapsed):
		for vessel in self.state["vessels"]:
			x = vessel.velocity * elapsed
			vessel.distance_left = max(0.0, vessel.distance_left - x)
			vessel.time_until_departure = round(max(0.0, vessel.time_until_departure - elapsed), 6)

	def timeAdvance(self):
		if len(self.state["vessels"]) > 0:
			v = self.state["vessels"][0]
			return v.time_until_departure
		return INFINITY

	def outputFnc(self):
		if len(self.state["vessels"]) > 0:
			return {
				self.vessel_out: [self.state["vessels"][0]]
			}
		return {}

	def intTransition(self):
		elapsed = self.timeAdvance()
		self.state["time"] += elapsed
		self.state["vessels"].pop(0)
		self.update_vessels(elapsed)
		if len(self.state["vessels"]) > 0:
			self.state["vessels"].sort(key=lambda v: v.time_until_departure)
		return self.state


class Scheduler(AtomicDEVS):
	def __init__(self, name, ivef):
		super(Scheduler, self).__init__(name)

		self.ivef = pd.read_csv(ivef, usecols=["mmsi", "start", "ETA", "source", "target", "task"], dtype={"mmsi": str})
		# self.ivef = self.ivef[self.ivef["mmsi"] == "205501790"]
		self.ivef.sort_values(by=["start"])
		self.starting_time = self.ivef.iloc[0]["start"]
		# self.ivef["start"] -= self.starting_time
		# self.ivef["end"] -= self.starting_time

		self.state = {
			"index": 0,
			"time": 0.0
		}

		self.reqs = self.addOutPort("reqs")

	def timeAdvance(self):
		if self.state["index"] < len(self.ivef):
			start = self.ivef.iloc[self.state["index"]]["start"] - self.starting_time
			return round(start / 1000 - self.state["time"], 6)
		return INFINITY

	def intTransition(self):
		self.state["time"] += self.timeAdvance()
		self.state["index"] += 1
		return self.state

	def outputFnc(self):
		if self.state["index"] < len(self.ivef):
			return {
				self.reqs: [self.ivef.iloc[self.state["index"]]]
			}
		return {}


class RoutePlanner(AtomicDEVS):
	def __init__(self, name):
		super(RoutePlanner, self).__init__(name)

		self.graph = get_graph()

		self.state = {
			"request": []
		}

		self.req_in = self.addInPort("req_in")
		self.req_out = self.addOutPort("req_out")

	def timeAdvance(self):
		if len(self.state["request"]) == 0:
			return INFINITY
		return 0.0

	def extTransition(self, inputs):
		if self.req_in in inputs:
			for request in inputs[self.req_in]:
				tug = TUGS[TUGS["MMSI"] == str(request["mmsi"])].iloc[0]
				request["name"] = tug["NAME"]
				request["category"] = tug["category"]

				request["source_lon"], request["source_lat"] = eval(request["source"])
				request["target_lon"], request["target_lat"] = eval(request["target"])

				src_vertex, _ = get_closest_vertex(self.graph, request["source_lon"], request["source_lat"])
				tgt_vertex, _ = get_closest_vertex(self.graph, request["target_lon"], request["target_lat"])

				request["distance"] = self.graph.distances([src_vertex.index], [tgt_vertex.index], weights="distance", mode="all")[0][0]

				# TODO: check this value with PoAB
				# TODO: make this a distribution instead?
				# request["velocity"] = 6.17  # about 12 knots
				request["velocity"] = request["distance"] / ((request["ETA"] - request["start"]) / 1000)

				# UNUSED!
				# request["end"] = request["start"] + request["distance"] / request["velocity"]

				self.state["request"].append(request)
		return self.state

	def outputFnc(self):
		if len(self.state["request"]) == 0:
			return {}
		return { self.req_out: [self.state["request"][0]] }

	def intTransition(self):
		if len(self.state["request"]) > 0:
			self.state["request"].pop(0)
		return self.state


class Clock(AtomicDEVS):
	def __init__(self, name, interval=1):
		super(Clock, self).__init__(name)
		self.interval = interval
		self.state = {
			"time": 0.0
		}
		self.outp = self.addOutPort("outp")

	def timeAdvance(self):
		return self.interval

	def outputFnc(self):
		return {
			self.outp: [True]
		}

	def intTransition(self):
		self.state["time"] += self.timeAdvance()
		return self.state


class Port(CoupledDEVS):
	def __init__(self, name, ivef, v_tugging=None, v_sailing=None, v_max=None):
		super(Port, self).__init__(name)

		self.clock = self.addSubModel(Clock("clock"))
		self.scheduler = self.addSubModel(Scheduler("scheduler", ivef))
		self.planner = self.addSubModel(RoutePlanner("planner"))
		# self.tracer = self.addSubModel(Tracer("tracer", ivef))
		self.pool = self.addSubModel(Pool("pool"))
		self.sailer = self.addSubModel(Sailer("sailer", v_tugging, v_sailing, v_max))

		# self.connectPorts(self.tracer.reqs, self.pool.req_in)
		self.connectPorts(self.scheduler.reqs, self.planner.req_in)
		self.connectPorts(self.planner.req_out, self.pool.req_in)
		self.connectPorts(self.pool.vessel_out, self.sailer.vessel_in)
		self.connectPorts(self.sailer.vessel_out, self.pool.vessel_in)

		self.connectPorts(self.clock.outp, self.sailer.update)

	def select(self, imm_children):
		for child in imm_children:
			if isinstance(child, Sailer):
				return child
		for child in imm_children:
			if isinstance(child, Pool):
				return child
		return imm_children[0]

	def print_statistics(self):
		# Costs // Dummy data for now
		C_Tug_S = {  # euro/hour
			"strong": 150,
			"average": 100,
			"weak": 50
		}
		C_Tug_P = {  # euro/hour
			"strong": 200,
			"average": 133,
			"weak": 67
		}
		C_Tug_F = {  # euro/day
			"strong": 1500,
			"average": 1000,
			"weak": 500
		}
		C_Tug_W = {  # euro/hour
			"strong": 100,
			"average": 67,
			"weak": 33
		}


		S1 = {
			"strong": 0.0,
			"average": 0.0,
			"weak": 0.0
		}
		S2 = {
			"strong": 0.0,
			"average": 0.0,
			"weak": 0.0
		}
		S3 = {
			"strong": 0.0,
			"average": 0.0,
			"weak": 0.0
		}
		S4 = {
			"strong": 0.0,
			"average": 0.0,
			"weak": 0.0
		}
		for vessel in list(self.pool.state["waiting"].values()) + self.sailer.state["vessels"]:
			T_Tug_S = vessel.footprint.sailing_time / 3600
			T_Tug_P = vessel.footprint.tugging_time / 3600
			T_Tug_W = vessel.footprint.idle_time / 3600

			S1[vessel.category] += C_Tug_S[vessel.category] * T_Tug_S
			S2[vessel.category] += C_Tug_P[vessel.category] * T_Tug_P
			S3[vessel.category] += C_Tug_W[vessel.category] * T_Tug_W
			S4[vessel.category] += C_Tug_F[vessel.category] * 31

		for cat in ["strong", "average", "weak"]:
			print("\ncategory:", cat)
			print("OVERALL COST: %f euros" % (S1[cat] + S2[cat] + S3[cat] + S4[cat]))
			print("Travelling: %f euros" % S1[cat])
			print("Tugging: %f euros" % S2[cat])
			print("Waiting: %f euros" % S3[cat])
			print("Fixed: %f euros" % S4[cat])
		# print("AVERAGE USAGES:")
		# S = 0
		# for vessel in list(self.pool.state["waiting"].values()) + self.sailer.state["vessels"]:
		# 	if vessel.footprint.usage_time > 0:
		# 		print("\t%s: %f m/s" % (vessel.mmsi, vessel.footprint.velocity_total / vessel.footprint.tasks))
		# 		S += vessel.footprint.usage_time
		# cost_per_hour = 100
		# print("-------------")
		# print("OVERALL: %f euros" % (S / 3600 * cost_per_hour))