{UTF8}
:MACRO: PINK NOISE(mean, std deviation, correlation time, seed)
PINK NOISE = INTEG(updating pink noise,mean+std deviation*RANDOM NORMAL(6,6,0,1,seed\
))
~ mean
~ Contributed by Ed Anderson, MIT/U. Texas  Austin
Description: The pink noise molecule described generates a simple random
series with autocorrelation. This is useful in representing time series,
like rainfall from day to day, in which today's value has some correlation
with what happened yesterday. This particular formulation will also have
properties such as standard deviation and mean that are insensitive both to
the time step and the correlation (smoothing) time. Finally, the output as
a whole and the difference in values between any two days is guaranteed to
be Gaussian (normal) in distribution.
Behavior: Pink noise series will have both a historical and a random
component during each period. The relative "trendtonoise" ratio is
controlled by the length of the correlation time. As the correlation time
approaches zero, the pink noise output will become more independent of its
historical value and more "noisy." On the other hand, as the correlation
time approaches infinity, the pink noise output will approximate a
continuous time random walk or Brownian motion. Displayed above are two
time series with correlation times of 1 and 8 months. While both series
have approximately the same standard deviation, the 1month correlation time
series is less smooth from period to period than the 8month series, which
is characterized by "sustained" swings in a given direction. Note that this
behavior will be independent of the timestep.
The "pink" in pink noise refers to the power spectrum of the output. A time
series in which each period's observation is independent of the past is
characterized by a flat or "white" power spectrum. Smoothing a time series
attenuates the higher or "bluer" frequencies of the power spectrum, leaving
the lower or "redder" frequencies relatively stronger in the output.
Caveats: This assumes the use of Euler integration with a time step of no
more than 1/4 of the correlation time. Very long correlation times should be
avoided also as the multiplication in the scaled white noise will become
progressively less accurate.
Technical Notes: This particular form of pink noise is superior to that of
Britting presented in Richardson and Pugh (1981) because the Gaussian
(Normal) distribution of the output does not depend on the Central Limit
Theorem. (Dynamo did not have a Gaussian random number generator and hence
R&P had to invoke the CLM to get a normal distribution.) Rather, this
molecule's normal output is a result of the observations being a sum of
Gaussian draws. Hence, the series over short intervals should better
approximate normality than the macro in R&P.
MEAN: This is the desired mean for the pink noise.
STD DEVIATION: This is the desired standard deviation for the pink noise.
CORRELATION TIME: This is the smooth time for the noise, or for the more technically \
minded this is the inverse of the filter's cutoff frequency in radians.
Updated by Tom Fiddaman, 2010, to include a random initial value,
correct units, and use TIME STEP$ keyword

updating pink noise = gap/correlation time
~ mean/correlation time
~ 
gap = scaled white noisePINK NOISE
~ mean
~ 
scaled white noise =mean+white noise*std deviation*SQRT((2time step$/correlation time
)/(time step$/correlation time))
~ mean
~ This adjusts the standard deviation of the white noise to compensate for the time \
step and the
correlation time to produce the appropriate pink noise std \
deviation.

white noise = RANDOM NORMAL(6,6,0,1,seed)
~ dmnl
~ This is an independent, identically distributed random quantity drawn every time \
step. The distribution is gaussian with mean = 0 and variance = 1.
Note that RANDOM NORMAL is truncated +/ 6 standard deviations here.
For Vensim 1.62 syntax, remove the arguments to RANDOM NORMAL.

:END OF MACRO:
mean=
0
~ 1/T
~ 
noise SD=
0.05
~ 1/T [0,1]
~ 
corr time=
1
~ T [0.25,5]
~ 
seed=
0
~ dmnl [0,10000,1]
~ 
slope=
2
~ 1/T [0,2]
~ 
disturbance=
PINK NOISE(mean,noise SD,corr time,seed)
~ 1/T
~ 
a=
1slope*(TimeINITIAL TIME)/(FINAL TIMEINITIAL TIME)
~ 1/T [1,1]
~ 
Flow=
a*SIN(Stock*2*3.14159)+disturbance
~ 1/T
~ 
Initial Stock=
0.5
~ dmnl [0,1]
~ 
Stock= INTEG (
Flow,
Initial Stock)
~ dmnl
~ 
********************************************************
.Control
********************************************************~
Simulation Control Parameters

FINAL TIME = 100
~ T
~ The final time for the simulation.

INITIAL TIME = 0
~ T
~ The initial time for the simulation.

SAVEPER =
TIME STEP
~ T [0,?]
~ The frequency with which output is stored.

TIME STEP = 0.125
~ T [0,?]
~ The time step for the simulation.

\\\/// Sketch information  do not modify anything except names
V300 Do not put anything below this section  it will be ignored
*View 1
$192192192,0,Arial100000000025511111196,96,100,0
10,1,Stock,345,185,40,20,3,3,0,0,0,0,0,0
12,2,48,117,187,10,8,0,3,0,0,1,0,0,0
1,3,5,1,4,0,0,22,0,0,0,111,,1(263,187)
1,4,5,2,100,0,0,22,0,0,0,111,,1(168,187)
11,5,48,216,187,6,8,34,3,0,0,1,0,0,0
10,6,Flow,216,204,16,9,40,3,0,0,1,0,0,0
1,7,1,6,1,0,0,0,0,64,0,111,,1(281,248)
10,8,Initial Stock,345,260,36,9,8,3,0,0,1,0,0,0
1,9,8,1,0,0,0,0,0,64,1,111,,1(345,235)
12,10,2495492,945,245,150,150,3,44,0,0,2,0,0,0
Stock,Graph
12,11,2825446,645,245,150,150,3,44,0,0,1,0,0,0
Phase
10,12,a,136,257,6,9,8,3,0,0,1,0,0,0
1,13,12,6,0,0,0,0,0,64,0,111,,1(166,237)
10,14,disturbance,193,295,36,9,8,3,0,0,0,0,0,0
1,15,14,6,0,0,0,0,0,64,0,111,,1(202,256)
10,16,corr time,228,367,31,9,8,3,0,0,1,0,0,0
1,17,16,14,0,0,0,0,0,64,0,111,,1(213,337)
10,18,noise SD,151,366,32,9,8,3,0,0,1,0,0,0
1,19,18,14,0,0,0,0,0,64,0,111,,1(167,336)
10,20,seed,187,411,18,9,8,3,0,0,1,0,0,0
1,21,20,14,0,0,0,0,0,64,0,111,,1(189,359)
10,22,FINAL TIME,65,338,45,9,8,2,0,3,1,0,0,0,128128128,000,10128128128
1,23,22,12,0,0,0,0,0,64,0,111,,1(96,301)
10,24,INITIAL TIME,59,311,48,9,8,2,0,3,1,0,0,0,128128128,000,10128128128
1,25,24,12,0,0,0,0,0,64,0,111,,1(94,285)
10,26,slope,31,279,18,9,8,3,0,0,1,0,0,0
1,27,26,12,0,0,0,0,0,64,0,111,,1(82,268)
10,28,Time,40,250,24,9,8,2,0,3,1,0,0,0,128128128,000,10128128128
1,29,28,12,0,0,0,0,0,64,0,111,,1(90,252)
10,30,mean,268,336,18,9,8,3,0,0,1,0,0,0
1,31,30,14,0,0,0,0,0,64,0,111,,1(236,318)
///\\\
:GRAPH Phase
:TITLE Phase
:XAXIS Stock
:XMIN 0
:XMAX 1
:NOLEGEND 2
:DOTS
:SCALE
:VAR Flow
:YMIN 1
:YMAX 1
:LINEWIDTH 3
:VAR Flow
:DATASET *2
:VAR Flow
:DATASET *3
:L<%^E!@
1:weak.vdf
1:slowing.vdf
9:slowing
22:$,Dollar,Dollars,$s
22:Day,Days
22:Hour,Hours
22:Month,Months
22:Person,People,Persons
22:Unit,Units
22:Week,Weeks
22:Year,Years
18:seed.vsc
20:keyvar.lst
15:0,0,0,1,0,0
19:100,0
27:2,
34:0,
4:Time
5:Stock
35:Date
36:YYYYMMDD
37:2000
38:1
39:1
40:2
41:0
24:0
25:100
26:100