#! ~/octave/RD/Keen-Minsky/KeenMinsky02.m
# An Octave program to verify model of firms and bank accounts
# including use of the the Phillips curve for non linearities, 
# as published by Prof Keen, his Model, MonetaryMinsky03  
# Octave Program author Richard A Lough. Published under the GPL License
# ver .02 published 20120525
######################################################################## 
# Set up the X matrix, initialise with ones
X= ones(12001,16); # Change this to ones(7202,16) to replicate Keen's model 
Y=zeros(12001,10); # Added for graphic output
for j = 1: 10      # this loop added to Keen's model for development work
Minsky = 12/(2+j); # Keen's model had a fixed value of 1.0
X(1,1)= 100; X(1,2)=12; X(1,3)=70; X(1,4)=5; X(1,5)= 13; X(1,6)= 1; X(1,7)=900;
##############constants and initial conditions ############################
# Initial Parameters  "Interest on loan""Interest on Deposit""Omega" "Sigma"
rl= 0.05 ; rd= 0.01 ; Ow= 0.1 ; sigma = 0.27;
# New parameters Output accelerator, population, Phillips variables
v=3; N=300; g=0.03;
x=0.95; y=0; s=2; m=-0.04;
# values for the function Inv
x1=0.04; y1=0.04; s1=2; m1=0; 
x2=0.03; y2=10; s2=100; m2=3;
x3=0.03; y3=2; s3=-50; m3=0.5;
L = 300; K = 900; II=1; Yr=K/v;
###############################################################################
#default values Alpha = 0.015; Beta = 0.02; PIr=0.0593;
Alpha = 0.015; Beta = 0.02; Delta = 0.01; 
Torw=1/26;Torp = 1;Torb=1; 
X(1,8)= 1; X(1,9)= 1.0; X(1,10)= 1.0; X(1,11)= 1.0; X(1,12)=1;
################ main loop starts here ########################################
###############################################################################
# loop through the desired range of time, say 100 years at tenths of months 
for i= 2:12001
  t=1/120;
  h=i-1; 
Yr=X(h,7)/v; L=Yr/X(h,8); Np= N*X(h,9);II=X(h,11);
  Ph = (y-m)*exp((II-x)*s/(y-m)) +m; # Rate of change of wages
  PIr= (Yr*X(h,10)+rd*X(h,3) -X(h,6)*L -rl*X(h,1))/(v*X(h,10)*Yr); 
  Inv = (y1-m1)*exp((PIr-x1)*s1/(y1-m1)) +m1; g= Inv/v - Delta; 
Taurl = (y2-m2)*exp((PIr-x2)*s2/(y2-m2)) +m2; # Firms time lag
Taulc = (y3-m3)*exp((PIr-x3)*s3/(y3-m3)) +m3; # Bank (Capital or Vault) time lag
#############################################################################################################
#############################################################################################################
###############Start simply and add components to understand the effects ####################################
#############################################################################################################
 X(i,1)= -X(h,1)/Taurl +X(h,2)/Taulc +X(h,10)*Yr*Inv*Minsky ;                         # Firms Loan
 X(i,2)= X(h,1)/Taurl  -X(h,2)/Taulc ;                                                # Bank Vault
 X(i,3)=-X(h,1)*(rl+1/Taurl)+X(h,2)/Taulc +rd*X(h,3)+X(h,4)/Torb ;                    # more on next line 
 X(i,3)=X(i,3)+X(h,5)/Torw -X(h,6)*Yr/X(h,8) +X(h,10)*Yr*Inv*Minsky;                  # Firms Deposit
 X(i,4)= rl*X(h,1) -rd*X(h,3) -X(h,4)/Torb  -rd*X(h,5);                               # Bank Capital
 X(i,5)= +(rd -1/Torw)*X(h,5) +X(h,6)*Yr/X(h,8);                                      # Workers Deposit 
 X(i,6)=(Ph+Ow*(g -(Alpha+Beta)) +(X(h,6)/(X(h,10)*X(h,8)*(1-sigma))-1)/Torp)*X(h,6); # Wages 
 X(i,7)= X(h,7)*g ;                                                                   # Capital
 X(i,8)= Alpha*X(h,8);                                                                # Rate of labour productivity growth
 X(i,9)= Beta*X(h,9);                                                                 # Rate of population growth
 X(i,10)= (X(h,6)/(X(h,8)*(1-sigma)) -X(h,10))/Torp;                                  # change of Price Level (inflation)
 X(i,11)= X(h,11)*(g -(Alpha+Beta));                                                  # change of employment Ld(t)
###########################################################################################
###################Increment the X variables by one step in time ##########################
###########################################################################################
 X(i,1)= X(h,1) +X(i,1)*t;
 X(i,2)= X(h,2) +X(i,2)*t;
 X(i,3)= X(h,3) +X(i,3)*t;
 X(i,4)= X(h,4) +X(i,4)*t;
 X(i,5)= X(h,5) +X(i,5)*t;
 X(i,6)= X(h,6) +X(i,6)*t;
 X(i,7)= X(h,7) +X(i,7)*t; 
 X(i,8)= X(h,8) +X(i,8)*t; 
 X(i,9)= X(h,9) +X(i,9)*t;
 X(i,10)= X(h,10) +X(i,10)*t;
 X(i,11)= X(h,11) +X(i,11)*t;
##########################################################################################
################Some things that need watching ###########################################
##########################################################################################
 X(h,13) = Yr;
 X(h,14) = Ph;
 X(h,15) = g ;
 X(h,16) = (X(h,6)/(X(h,10)*X(h,8)*(1-sigma))-1)/Torp;
 Y(h,j)  = X(h,3);
##########################################################################################
#######################Print out some results ############################################
##########################################################################################
endfor
endfor 
printf ("Firms Loan: ");
printf ("%8.4f",X(h,1));
printf ("\n"); 
printf ("Bank Vault: ");
printf ("%8.4f",X(h,2));
printf ("\n");
printf ("Firms Deposit: ");
printf ("%8.4f",X(h,3));
printf ("\n");
printf ("Bank Capital: "); 
printf ("%8.4f",X(h,4));
printf ("\n");
printf ("Workers Deposit: ");
printf ("%8.4f",X(h,5));
printf ("\n");
printf ("wage : "); 
printf ("%8.4f",X(h,6));
printf ("\n");
printf ("Capital = : "); 
printf ("%8.4f",X(h,7));
printf ("\n"); 
printf ("Alpha = : "); 
printf ("%8.4f",X(h,8));
printf ("\n");
printf ("Beta =: "); 
printf ("%8.4f",X(h,9));
printf ("\n");
printf ("change of Price Level (inflation): ");
printf ("%8.4f",X(h,10));
printf ("\n"); 
printf ("Employment ratio : ");
printf ("%8.4f",X(h,11));
printf ("\n");
#printf ("%8.4f, %8.4f, %8.4f,%8.4f, %8.4f, %8.4f,%8.4f, %8.4f, %8.4f,%8.4f, %8.4f, %8.4f\n", X);
#
plot(Y);
fflush (stdout);
##############################################################################################
# This model follows carefully copied equations etc from Prof Keen's MonetaryMinsky03 program. 
# An outer loop was then added to demonstrate Model sensitivity to the Minsky effect.
# The simulation length is increased here to 100 years for development purposes.
# MonetaryMinsky03's length of simulation was undefined 
# The theoretical background to the model is given here:  http://www.youtube.com/watch?v=0SPqMDMbRlo
# and the equations appear around 1h02m into that video presentation. 
# Keen's original MathCad model is available here : 
# http://www.debtdeflation.com/blogs/wp-content/uploads/qed/QED.zip
# This Keen-Minsky model is published to aid Octave code development, and as as a first step 
# in standardising the process of developing the Model and developing software tools. 
# For reference and validation, the figures output by the program are as follows: 
# Firms Loan:   0.9504
# Bank Vault: 331644.3730
# Firms Deposit: -11.2287
# Bank Capital:   0.1791
# Workers Deposit:   0.0000
# wage :   0.0000
# Capital = : 4236.8343
# Alpha = :   4.4807
# Beta =:   7.3866
# change of Price Level (inflation):   0.0000
# Employment ratio :   0.1422 
########################################################################################