diff --git a/src/main/scala/symsim/examples/concrete/pumping/Pump.scala b/src/main/scala/symsim/examples/concrete/pumping/Pump.scala index c9a2b072..dffb4c69 100644 --- a/src/main/scala/symsim/examples/concrete/pumping/Pump.scala +++ b/src/main/scala/symsim/examples/concrete/pumping/Pump.scala @@ -1,13 +1,8 @@ package symsim package examples.concrete.pumping -/** A case study of control of a water extraction at Marbjerg (operated by HOFOR). - * - * Source: Andreas Holck Høeg-Petersen. Reinforcement Learning for Controlling - * Groundwater Extraction. MSc thesis at the IT University of Copenhagen. - * August 2021. (esp. Chapters 2 and 3) - */ import symsim.concrete.Randomized +import scala.math /** Constants as fit in the model of AHHP (reference above). @@ -15,8 +10,8 @@ import symsim.concrete.Randomized * They are called a, b, c, and d parameters in his report. */ val c1: Double = 1.0 / 12 -val c2: Double = 0.15 / (c1*80*12*24*365) -val c3: Double = 0.1 / (12*24*365) +val c2: Double = 0.15 / (c1 * 80 * 12 * 24 * 365) +val c3: Double = 0.1 / (12 * 24 * 365) val c4: Double = 0.05 /** Size of the history window for the head-level */ @@ -33,6 +28,14 @@ case class PumpState ( w: Double, // Water level in the ground deposits (unobservable) phm: List[Double]): // Past head means (a sliding window of history) + require (f >= FlowMin, s"Flow is too low: $f") + require (f <= FlowMax, s"Flow is too high: $f") + require (h <= HeadMax, s"Head is too high: $h") + require (hm >= HeadMin, s"Head mean is too low: $hm") + require (hm <= HeadMax, s"Head mean is too high: $hm") + require (w >= WaterMin, s"Water level is too low: $w") + require (w <= WaterMax, s"Water level is too high: $w") + override def toString: String = s"[flow=$f, head=$h, head mean=$hm, tank level=$tl, time=$t, " + s"water=$w, past head means=$phm]" @@ -42,9 +45,9 @@ end PumpState /** Discretize cut points for flow, head, and tank variables./ * - * We observe to the first point on the list that is lower than the actual - * value of the variable. So we 'round down'. Values below the last element - * in the cut-point list are rounded up (to the last value). + * We discretize to the first point on the list that is lower than the actual + * value of the variable. So we 'round down'. Values below the last element in + * the cut-point list are rounded up (to the last value). * * For now, we have preconditions about low values in discretization so the * getOrElse (rounding up) should never happen. @@ -57,24 +60,24 @@ def closest (value: Double) (cutPoints: List[Double]): Double = .find { value >= _ } .getOrElse (cutPoints.last) -val HEAD_MIN: Double = 7.0 -val HEAD_MAX: Double = 10.84 -val FLOW_MIN: Double = 0.0 -val FLOW_MAX: Double = 120.0 -val TANK_MIN: Double = 0.0 -val TANK_MAX: Double = 2000.0 -val WATER_MIN: Double = 0.0 -val WATER_MAX: Double = 9.11 +val HeadMin = 7.0 +val HeadMax = 10.84 +val FlowMin = 0.0 +val FlowMax = 120.0 +val TankMin = 0.0 +val TankMax = 2000.0 +val WaterMin = 0.0 +val WaterMax = 9.11 -val flowCutPoints = List (FLOW_MAX, 115.0, 110.0, 105.0, 100.0, 95.0, - 90.0, 85.0, 80.0, 75.0, 70.0, 65.0, 60.0, 55.0, 50.0, FLOW_MIN) +val FlowCutPoints = List (FlowMax, 115.0, 110.0, 105.0, 100.0, 95.0, + 90.0, 85.0, 80.0, 75.0, 70.0, 65.0, 60.0, 55.0, 50.0, FlowMin) -val headCutPoints = List (HEAD_MAX, 10.68, 10.52, 10.36, 10.2, 10.04, +val HeadCutPoints = List (HeadMax, 10.68, 10.52, 10.36, 10.2, 10.04, 9.88, 9.72, 9.56, 9.4, 9.24, 9.08, 8.92, 8.76, 8.6, 8.44, 8.28, - 8.12, 7.96, 7.8, 7.64, 7.48, 7.32, 7.16, HEAD_MIN) + 8.12, 7.96, 7.8, 7.64, 7.48, 7.32, 7.16, HeadMin) -val tankCutPoints = List (TANK_MAX, 1800.0, 1600.0, 1400.0, - 1200.0, 1000.0, 800.0, 600.0, 400.0, 200.0, TANK_MIN) +val TankCutPoints = List (TankMax, 1800.0, 1600.0, 1400.0, + 1200.0, 1000.0, 800.0, 600.0, 400.0, 200.0, TankMin) /** The (discrete) state, observable by the pump controller. */ case class ObservablePumpState ( @@ -83,6 +86,11 @@ case class ObservablePumpState ( hm: Double, // Discretized, derived head mean tl: Double): // Discretized observable tank level + require (FlowCutPoints.contains (f), s"Wrong discrete flow: $f") + require (HeadCutPoints.contains (h), s"Wrong discrete head level: $h") + require (HeadCutPoints.contains (hm), s"Wrong discrete mean head level: $hm") + require (TankCutPoints.contains (tl), s"Wrong discrete tank level: $tl") + override def toString: String = s"(flow=$f, head=$h, head mean=$hm, tank level=$tl)" @@ -93,6 +101,12 @@ type PumpAction = Double type PumpReward = Double +/** A case study of control of a water extraction at Marbjerg (operated by HOFOR). + * + * Source: Andreas Holck Høeg-Petersen. Reinforcement Learning for Controlling + * Groundwater Extraction. MSc thesis at the IT University of Copenhagen. + * August 2021. (esp. Chapters 2 and 3) + */ object Pump extends Agent[PumpState, ObservablePumpState, PumpAction, PumpReward, Randomized], Episodic: @@ -103,22 +117,20 @@ object Pump extends /** If isFinal is true in a state, then the episode is complete */ def isFinal (s: PumpState): Boolean = - s.t >= 4000 || s.h < HEAD_MIN || s.tl > TANK_MAX || s.tl < TANK_MIN + s.t >= 4000 || s.h < HeadMin || s.tl > TankMax || s.tl < TankMin def observe (s: PumpState): ObservablePumpState = - require (s.tl >= TANK_MIN, s"s.tl = ${s.tl} >= $TANK_MIN") - - val df = closest (s.f) (flowCutPoints) - val dh = closest (s.h) (headCutPoints) - val dhm = closest (s.hm) (headCutPoints) - val dtl = closest (s.tl) (tankCutPoints) - + require (s.tl >= TankMin, s"s.tl = ${s.tl} >= $TankMin") + val df = closest (s.f) (FlowCutPoints) + val dh = closest (s.h) (HeadCutPoints) + val dhm = closest (s.hm) (HeadCutPoints) + val dtl = closest (s.tl) (TankCutPoints) ObservablePumpState (df, dh, dhm, dtl) def reward (source: PumpState) (target: PumpState) (a: PumpAction): Double = - if target.h < HEAD_MIN || target.tl < TANK_MIN || target.tl > TANK_MAX + if target.h < HeadMin || target.tl < TankMin || target.tl > TankMax then -9999 else val flowReward = if target.f != source.f then -0.5 else 0 @@ -126,9 +138,9 @@ object Pump extends flowReward - headCost - val amp: Double = 0.41852857594808646 - val freq: Double = 0.0000597030105413 - val phase: Double = -6347.109214682171 + val Amp = 0.41852857594808646 + val Freq = 0.0000597030105413 + val Phase = -6347.109214682171 override def step (s: PumpState) (a: PumpAction) @@ -141,10 +153,10 @@ object Pump extends cd <- getDemand (s.t%24 + 1) d = cd + nd tl1 = s.tl + c1 * (f1 - d) - h1 = s.h + c4 * (s.w + (c1*f1 / Math.PI)) + h1 = s.h + c4 * (s.w + (c1*f1 / math.Pi)) nw <- Randomized.gaussian (0.0, 1.0) w1 = s.w - c2*(c1*f1) + c3 + - (amp*Math.sin (2*Math.PI*(s.t + phase) / freq)) + nw + (Amp * math.sin (2 * math.Pi*(s.t + Phase) / Freq)) + nw hm1 = (1.0 / k)*s.phm.sum phm1 = (s.hm:: s.phm).slice (0, k) s1 = PumpState (f1, h1, hm1, tl1, s.t + 1, w1, phm1) @@ -196,10 +208,10 @@ object PumpInstances given enumState: BoundedEnumerable[ObservablePumpState] = val ss: List[ObservablePumpState] = for - f <- flowCutPoints - h <- headCutPoints - hm <- headCutPoints - tl <- tankCutPoints + f <- FlowCutPoints + h <- HeadCutPoints + hm <- HeadCutPoints + tl <- TankCutPoints yield ObservablePumpState (f, h, hm, tl) BoundedEnumerableFromList (ss*) @@ -213,13 +225,13 @@ object PumpInstances concrete.Randomized.canTestInRandomized lazy val genPumpState: Gen[PumpState] = for - f <- Gen.choose[Double] (FLOW_MIN, FLOW_MAX) - h <- Gen.choose[Double] (HEAD_MIN, HEAD_MAX) - hm <- Gen.choose[Double] (HEAD_MIN, HEAD_MAX) - tl <- Gen.choose[Double] (TANK_MIN, TANK_MAX) + f <- Gen.choose[Double] (FlowMin, FlowMax) + h <- Gen.choose[Double] (HeadMin, HeadMax) + hm <- Gen.choose[Double] (HeadMin, HeadMax) + tl <- Gen.choose[Double] (TankMin, TankMax) t <- Gen.choose[Int] (0, 24) - w <- Gen.choose[Double] (WATER_MIN, WATER_MAX) - phm <- Gen.listOfN (5, Gen.choose (HEAD_MIN, HEAD_MAX)) + w <- Gen.choose[Double] (WaterMin, WaterMax) + phm <- Gen.listOfN (5, Gen.choose (HeadMin, HeadMax)) yield PumpState (f, h, hm, tl, t, w, phm) given arbitraryState: Arbitrary[PumpState] = Arbitrary (genPumpState) diff --git a/src/test/scala/symsim/examples/concrete/pumping/PumpIsAgentSpec.scala b/src/test/scala/symsim/examples/concrete/pumping/PumpIsAgentSpec.scala index ed68f84c..494904be 100644 --- a/src/test/scala/symsim/examples/concrete/pumping/PumpIsAgentSpec.scala +++ b/src/test/scala/symsim/examples/concrete/pumping/PumpIsAgentSpec.scala @@ -1,8 +1,8 @@ package symsim package examples.concrete.pumping -import laws.AgentLaws -import laws.EpisodicLaws +import symsim.laws.AgentLaws +import symsim.laws.EpisodicLaws class PumpIsAgentSpec extends SymSimSpec: diff --git a/src/test/scala/symsim/examples/concrete/pumping/PumpSpec.scala b/src/test/scala/symsim/examples/concrete/pumping/PumpSpec.scala index 92ad6baf..491cbe06 100644 --- a/src/test/scala/symsim/examples/concrete/pumping/PumpSpec.scala +++ b/src/test/scala/symsim/examples/concrete/pumping/PumpSpec.scala @@ -2,8 +2,8 @@ package symsim package examples.concrete.pumping import symsim.concrete.Randomized.given -import CanTestIn.given -import CanTestIn.* +import symsim.CanTestIn.given +import symsim.CanTestIn.* import org.scalacheck.{Arbitrary, Gen, Prop} import org.scalacheck.Arbitrary.* @@ -19,13 +19,13 @@ class PumpSpec extends org.scalacheck.Properties ("Pumping"): // Generators of test data - val flow = Gen.choose[Double] (FLOW_MIN, FLOW_MAX) - val head = Gen.choose[Double] (HEAD_MIN, HEAD_MAX) - val head_mean = Gen.choose[Double] (HEAD_MIN, HEAD_MAX) - val tank = Gen.choose[Double] (TANK_MIN, TANK_MAX) + val flow = Gen.choose[Double] (FlowMin, FlowMax) + val head = Gen.choose[Double] (HeadMin, HeadMax) + val head_mean = Gen.choose[Double] (HeadMin, HeadMax) + val tank = Gen.choose[Double] (TankMin, TankMax) val time = Gen.choose[Int] (0, 24) - val water = Gen.choose[Double] (WATER_MIN, WATER_MAX) - val past_head_mean = Gen.listOfN (5, Gen.choose (HEAD_MIN, HEAD_MAX)) + val water = Gen.choose[Double] (WaterMin, WaterMax) + val past_head_mean = Gen.listOfN (5, Gen.choose (HeadMin, HeadMax)) val actions = Gen.oneOf (Pump.instances.enumAction.membersAscending) val obsStates = Gen.oneOf (Pump.instances.allObservableStates) @@ -36,7 +36,7 @@ class PumpSpec forAll (flow, head, head_mean, tank, time, water, past_head_mean, actions) { (f, h, hm, tl, t, w, phm, a) => for (s1, r) <- Pump.step (PumpState (f, h, hm, tl, t, w, phm)) (a) - yield s1.tl >= TANK_MIN + yield s1.tl >= TankMin } } @@ -48,9 +48,9 @@ class PumpSpec } } - property ("The time can be more than TANK_MAX (comparing is not correct)") = { + property ("The time can be more than TankMax (comparing is not correct)") = { exists (time) { t => - t > TANK_MAX + t > TankMax } } @@ -65,7 +65,7 @@ class PumpSpec (f, h, hm, tl, t, w, phm) => exists(actions) { a => for (s1, r) <- Pump.step(PumpState(f, h, hm, tl, t, w, phm))(a) - yield (s1.tl <= TANK_MAX) + yield (s1.tl <= TankMax) } } }