epogrebnyak
commented
6 years ago Honestly, I never understood why you had a second round of distributions, this seems unusual for me. The purpose of Monte Carlo is to seed random input to a rather complex system (usually, where to an analytic solution is possible) and obtain/aggregate/interpret the results.
At the interpretation stage, you might have a question "can resulting empiric distribution be explained by a simpler law?", eg when the result looks like normal distribution you can substitute running the experiments with a proxy distribution, with some degree of certainty. That can help if (a) running the simulations takes a long time (b) their result is not solvable analytically, (c) you need to plug dist to somewhere else to approximate the system inspected.
But this is a step after the MC experiments are run and only if you need reducing your result to a distribution (which in turn can be a good or a bad approximation of the result). Otherwise, you just seed the inputs and get many results, which form a discrete experimental distribution.
The equations you have for the geothermal system are quite simple and you can have an analytic result of system cost distribution even without MC simulation. It is still a good exercise, but the result is very much as expected.
Generally, what is the research question that you try to solve by looking at different systems? There must be a real-world problem that calls for examining with these experiments. The source of randomness can be weather (costs of surviving extreme temperatures), adding pump repair costs, interlink between electricity market and geothermal (running a heat pump at night time when electricity cost is low), among others.
Werderaner
I just had a closer look at your new code and I fully understand it now. This code is kind of similar to the code I initially created, but way shorter and easier.
cop_gw = seed_triangular(DISTRIBUTION_PARAM_COP_GWP) cop_hp = seed_triangular(DISTRIBUTION_PARAM_COP_HP) p = seed_triangular(DISTRIBUTION_PARAM_PRICE)
ec_gw = (ENERGY_HEATING_MWH + ENERGY_COOLING_MWH) / cop_gw
ec_hp = ENERGY_HEATING_MWH / (cop_hp - 1)
total_costs_np = (ec_gw + ec_hp) * p
However (maybe its because I am too tired), I am not sure at the moment if we are simulating a "real" Monte Carlo simulation. Because what we do is just summing the triangular and "kind of" normal distributions. The size of the distributions is defined by N = 10^6.
But if you consider my original code, you maybe remember that I also included some iterations over the summation at the end:
def Total_electricity_costs(): PowercostsHP = np.random.normal(loc=MEAN2, scale=STD2) PowercostsPC = np.random.normal(loc=MEAN3, scale=STD3) PowercostsPH = np.random.normal(loc=MEAN4, scale=STD4) return PowercostsHP +PowercostsPH + PowercostsPC sim5 = np.zeros(N) for i in range(N): sim5[i] = Total_electricity_costs() Plot = plt.hist(sim5, bins=75, density=True)
To illustrate this, a little exaple: If you consider 3 parameters (a,b,c) and each of them represents a set of values. a = [1,2,...10], b = [10,11,...20], c = [100,101,...110]. Each value has a certain exceedens probability.
To get the result of a + b + c, what we do right now is: 1+10+100, 2+11+101, 3+12+102 ... and so on. Am I wrong?
But what we should do is (10^6 times): 1 + 15 + 109, 5 + 20 + 100, 9 + 12 + 105 , so randomly chosen values of a,b,c based on the probability of the values.
Do you understand what I mean? What do you think? OR am I totally wrong? Maybe I do not understand the concept of numpy as a whole.