Model Description
The authors address the question of how infrastructure design affects SES sustainability in two stages. First, they explore the effects of design variations in shared infrastructure on long-term system behavior in a model system. They examine two types of distribution infrastructure, one with and one without upstream−downstream asymmetry, and different threshold characteristics of infrastructure maintenance. Second, they evaluate how these design variations influence the robustness of system function to an economic shock.
The model assumes that there are N farming households spread across two villages (village 1 and village 2) that manage a shared irrigation infrastructure. Each farmer is endowed with the same amount of available labor (l) each year and the same acreage (a). A farmer appropriates q units of water from the system and allocates labor among three activities: farming, maintaining infrastructure and outside employment at a wage rate w. Governance is represented in the model by the following rules. The expected maintenance labor contribution is proportional to the farmer's acreage (assumed to be same for all farmers in this model). Water allocations are also proportional to acreage, but only among the water rights holders. Only farmers who contributed labor to the infrastructure maintenance before the planting season obtain water rights. Farmers choose between two strategies: group conformist (G) and opportunist (O). Gs follow and enforce rules, and strive to maximize the total welfare of the two villages. Os break the rules and attempt to maximize the individual net income. The model tracks the fraction of Gs in the village i and the resulting performance of the irrigation system.
$\begin{equation}
\frac{dX_i}{dt} = X_i(\pi_i^G - \overline{\pi_i})
\end{equation}$ |
This replicator equation tracks the fractions of the strategies in villages 1 and 2. |
$\begin{equation}
Q = I(L_m)S(t)
\end{equation}$ |
The total irrigation water, Q, depends on the efficiency of the shared irrigation infrastructure and the volume of river discharge. Efficiency of the irrigation infrastructure depends on the sum of maintenance labor provided by all farmers each year. |
$I(L_m) =
\begin{cases}
0 & 0 \leq L_m < \psi - \epsilon \\
\frac{I_{max}}{2\epsilon} (L_m - \psi + \epsilon) & \psi - \epsilon \leq L_m \leq \psi + \epsilon\\
I_{max} & \psi + \epsilon \leq L_m \leq L
\end{cases}$ |
I-max is the maximum infrastructure efficiency. The half-saturation point of labor and half-width of the threshold slope determine the threshold of maintenance. |
$\begin{equation}
l = l_f + l_m + l_e
\end{equation}$ |
A farmer can allocate labor among three activities: farming (l-f), maintaining infrastructure (l-m) and outside employment (l-e) |
$\begin{equation}
\pi_i = pb(l_f)^j(q_i + r_i)^ka_i^(1-j-k) + wl_e
\end{equation}$ |
The equation represents a farmer's total income in village i |
$\begin{equation}
\begin{aligned}
& \underset{L_f,L_m}{\text{max}}
{\Pi = pbL_f^j[I(L_m)S + R]^k A^(1-j-k) + w L_e} \\
&\text{s.t.} L_f + L_m + L_e = L \\
\end{aligned}
\end{equation}$ |
The aggregate income of two villages maximized over the optimal values of total maintenance labor, farming, and employment. L is the total available labor of the farmers in the two villages. The model assumes each farmer is endowed with the same acreage. The uppercase symbols denote aggregate-level quantities. |
$\begin{equation}
\pi_1^G = pb(l_f^G)^j(q_1^G + r_1)^ka_1^(1-j-k) + w(l_e^G) - [\gamma_s(1 - X_1) + \gamma_o(1 - X_2)] \\
\pi_2^G = pb(l_f^G)^j(q_2^G + r_2)^ka_2^(1-j-k) + w(l_e^G) - [\gamma_s(1 - X_2) + \gamma_o(1 - X_1)]
\end{equation}$ |
The payoffs of group conformists in village 1 and 2. The cost of enforcement for a G increases with the frequencies of Os in their own village and the other village. |
$\begin{equation}
\pi_1^O = pb(l_f^O)^j(q_1^O + r_1)^ka_1^(1-j-k) + w(l_e^O) - \delta\bigg(1 - \sigma \frac{Q(L_m)}{Q(L_m*)}\bigg)q_1^O\frac{(X_1 + X_2)}{2}\\
\pi_2^O = pb(l_f^O)^j(q_2^O + r_2)^ka_2^(1-j-k) + w(l_e^O) - \delta\bigg(1 - \sigma \frac{Q(L_m)}{Q(L_m*)}\bigg)q_2^O\frac{(X_1 + X_2)}{2}
\end{equation}$ |
The payoffs of opportunists in the two villages. The probability of opportunists being caught increases with the average frequency of Gs. The penalty varies by situation: it increases with the amount of water stolen, but decreases with water abundance in the system. |
$\begin{equation}
Q(L_m) = I(L_m)S
Q(L_m^*) = I(L_m^*)S
Abdundance = \frac{Q(L_m)}{Q(L_m^*)}
\end{equation}$ |
Q(L-m*) is the optimal amount of water desired by the two villages. Q(L-m) is the total amount of water actually available in the system. |
$\begin{equation}
q_1^O \leq \frac{Q}{N_1} \\
q_1^G = min[\frac{Q}{N^G} , (Q - Q_1^O) / N_1^G]
\end{equation}$ |
Due to their privileged access, Os in village 1 are less constrained by a higher bound on the amount of water they can steal. Gs in village 1 rely on their upstream position to extract water to bring their actual amount as close as possible to allocated amount. |
$\begin{equation}
q_2^O \leq \frac{Q - Q_1^O - Q_1^G)}{N_2}\\
q_2^G = \frac{[Q - Q_1^O - Q_1^G - Q_2^O]}{N_2^G}
\end{equation}$ |
The amount of water that Os in village 2 can steal is constrained on the water available after consumption by Gs and Os in village 1. The second equation is the amount of water available for Gs in village 2. *Note that the model assumes that both Gs and Os in village 1 compete for water at the same time and they are symmetric in their capacity for competing for water. For the same reason, Os in village 2 face an upper bound on the water they obtain. |
#delta=2.9 / 1.0 R=30 / 0 par j=0.3,k=0.4 par p=1,b=1 par N1=50,N2=50 par ind_l=1 par ind_a1=1 par ind_a2=1 par Imax=1,psix=0.2, epsilonx=0.125, S=100 par w=0.2 par R=0 par gamma_s=0.05, gamma_o=0.1 par delta=1.4 par sigma=0.9 A1=(N1*ind_a1) A2=(N2*ind_a2) L=(N1+N2)*ind_l ind_r1=R*(ind_a1/(A1+A2)) ind_r2=R*(ind_a2/(A1+A2)) psi=L*psix epsilon=L*epsilonx #################################################### #functions infra(x)=if(x<=(psi-epsilon)) then (0) else (if(x>=(psi+epsilon)) then (Imax) else (((Imax)/(2*epsilon))*(x-psi+epsilon))) f_income(lfx,qx,rx,ax)=p*b*((lfx)^(j))*((qx+rx)^(k))*((ax)^(1-j-k)) up_bound(x,up)=if(x>=up) then (up) else (x) low_bound(x,low)=if(x<=low) then (low) else (x) bound(x,low,up)=if(x<=low) then (low) else (if(x>=up) then (up) else (x)) step(x,ref)=if(x>ref) then (1) else (0) #################################################### # first block A=A1+A2 Qmax=Imax*S CC=p*b*A^(1-j-k) QR=(Qmax)/(2*epsilon) #################################################### # Region 1 functions Reg1Dec(Lf1tx)=if(Lf1tx<L) then (1) else (2) Reg1Ymax(Dc,Y1max_ax,Y1max_bx)=if(Dc==1) then (Y1max_ax) else (Y1max_bx) Reg1Lf(Dc,Lf1a_starx,Lf1b_starx)=if(Dc==1) then (Lf1a_starx) else (Lf1b_starx) # Region 1 Lm1t=0 Lf1t=(CC*j*R^k/w)^(1/(1-j)) Lf1a_star=Lf1t Y1max_a=CC*(Lf1t^j)*(R)^k+w*(L-Lf1t) Lf1b_star=L Y1max_b=CC*(Lf1b_star)^j*R^k D1=Reg1Dec(Lf1t) Lm_star_1=0 Lf_star_1=Reg1Lf(D1,Lf1a_star,Lf1b_star) Y1max=Reg1Ymax(D1,Y1max_a,Y1max_b) #################################################### # Region 3 functions Reg3Dec(Lf3tx,Lm3tx,x)=if((Lf3tx+Lm3tx)<x) then (1) else (2) Reg3Ymax(Dc,Y3max_ax,Y3max_bx)=if(Dc==1) then (Y3max_ax) else (Y3max_bx) Reg3Lf(Dc,Lf3a_starx,Lf3b_starx)=if(Dc==1) then (Lf3a_starx) else (Lf3b_starx) Reg3Lm(Dc,Lm3a_starx,Lm3b_starx)=if(Dc==1) then (Lm3a_starx) else (Lm3b_starx) # Region 3 Lm3t=psi+epsilon Lf3t=(CC*j*(Qmax+R)^k/w)^(1/(1-j)) Lf3a_star=Lf3t Lm3a_star=Lm3t Y3max_a=CC*Lf3a_star^j*(Qmax+R)^k+w*(L-Lf3a_star- Lm3a_star) Lf3b_star=L-psi-epsilon Lm3b_star=Lm3t Y3max_b=CC*(Lf3b_star)^j*(Qmax+R)^k temp=L D3=Reg3Dec(Lf3t,Lm3t,temp) Lm_star_3=Reg3Lm(D3,Lm3a_star,Lm3b_star) Lf_star_3=Reg3Lf(D3,Lf3a_star,Lf3b_star) Y3max=Reg3Ymax(D3,Y3max_a,Y3max_b) #################################################### # Region 2 functions Reg2Dec(L2sumx,Lm2tx,Lm2Dx,Lm2f_starx,Lf2f_starx)=if(L2sumx<=L) then (DecC(Lm2tx)) else (DecF(Lm2Dx,Lm2f_starx,Lf2f_starx)) DecB(Lm2tx)=if(Lm2tx<(psi-epsilon)) then (2) else (1) DecC(Lm2tx)=if(Lm2tx>(psi+epsilon)) then (3) else (DecB(Lm2tx)) DecE(Lm2Dx)=if(Lm2Dx<(psi-epsilon)) then (5) else (4) DecF(Lm2Dx,Lm2f_starx,Lf2f_starx)=if(Lm2Dx>(psi+epsilon)) then (DecG(Lm2f_starx,Lf2f_starx)) else (DecE(Lm2Dx)) DecG(Lm2f_starx,Lf2f_starx)=if((Lm2f_starx+Lf2f_starx)>=L) then (7) else (6) Reg2Ymax(Dc,Y2max_ax,Y2max_bx,Y2max_cx,Y2max_dx,Y2max_ex,Y2max_fx,Y2max_gx)=if(Dc==1) then (Y2max_ax) else (if(Dc==2) then (Y2max_bx) else (if(Dc==3) then (Y2max_cx) else (if(Dc==4) then (Y2max_dx) else (if(Dc==5) then (Y2max_ex) else (if(Dc==6) then (Y2max_fx) else (Y2max_gx)))))) Reg2Lf(Dc,Lf2a_starx,Lf2b_starx,Lf2c_starx,Lf2d_starx,Lf2e_starx,Lf2f_starx,Lf2g_starx)=if(Dc==1) then (Lf2a_starx) else (if(Dc==2) then (Lf2b_starx) else (if(Dc==3) then (Lf2c_starx) else (if(Dc==4) then (Lf2d_starx) else (if(Dc==5) then (Lf2e_starx) else (if(Dc==6) then (Lf2f_starx) else (Lf2g_starx)))))) Reg2Lm(Dc,Lm2a_starx,Lm2b_starx,Lm2c_starx,Lm2d_starx,Lm2e_starx,Lm2f_starx,Lm2g_starx)=if(Dc==1) then (Lm2a_starx) else (if(Dc==2) then (Lm2b_starx) else (if(Dc==3) then (Lm2c_starx) else (if(Dc==4) then (Lm2d_starx) else (if(Dc==5) then (Lm2e_starx) else (if(Dc==6) then (Lm2f_starx) else (Lm2g_starx)))))) # basics Lf2t=((j*CC/w)*(k*QR/j)^k)^(1/(1-j-k)) Lm2t=psi-epsilon-(R/QR)+(1/QR)*( (j*CC/w)*(QR*k/j)^(1-j) )^(1/(1-j-k)) L2sum=Lf2t+Lm2t # Region 2 Lf2a_star=Lf2t Lm2a_star=Lm2t Y2max_a=(w/j)*(1-j)*Lf2a_star+w*(L-Lm2a_star) Lf2b_star=(CC*j*R^k/w)^(1/(1-j)) Lm2b_star=psi-epsilon Y2max_b=(w/j)*(1-j)*(Lf2b_star)+w*(L-Lm2b_star) Lf2c_star=(CC*j*(Qmax+R)^k/w)^(1/(1-j)) Lm2c_star=psi+epsilon Y2max_c=(w/j)*(1-j)*Lf2c_star + w*(L-Lm2c_star) Lm2D=(k*L+j*(psi-epsilon)-j*R/QR)/(j+k) Lf2D=L-Lm2D Lm2d_star=Lm2D Lf2d_star=Lf2D Qd=R+QR*(Lm2d_star-psi+epsilon) Y2max_d=CC*(Qd)^k*(Lf2d_star)^j Lm2e_star=psi-epsilon Lf2e_star=L-Lm2e_star Y2max_e=CC*(L-Lm2e_star)^j*R^k Lm2f_star=psi+epsilon Lf2f_star=(CC*j*(Qmax+R)^k/w)^(1/(1-j)) Y2max_f=CC*(Lf2f_star)^j*(Qmax+R)^k+w*(L-Lm2f_star-Lf2f_star) Lm2g_star=psi+epsilon Lf2g_star=L-Lm2g_star Y2max_g=CC*(L-psi-epsilon)^j*(Qmax+R)^k D2=Reg2Dec(L2sum,Lm2t,Lm2D,Lm2f_star,Lf2f_star) Lm_star_2=Reg2Lm(D2,Lm2a_star,Lm2b_star,Lm2c_star,Lm2d_star,Lm2e_star,Lm2f_star,Lm2g_star) Lf_star_2=Reg2Lf(D2,Lf2a_star,Lf2b_star,Lf2c_star,Lf2d_star,Lf2e_star,Lf2f_star,Lf2g_star) Y2max=Reg2Ymax(D2,Y2max_a,Y2max_b,Y2max_c,Y2max_d,Y2max_e,Y2max_f,Y2max_g) #################################################### # Final decision functions FinalDec(Y1maxx,Y2maxx,Y3maxx)=if(Y1maxx>Y2maxx) then (if(Y1maxx>Y3maxx) then (1) else (3)) else (2) FinalYmax(Dc,Y1maxx,Y2maxx,Y3maxx)=if(Dc==1) then (Y1maxx) else (if(Dc==2) then (Y2maxx) else (Y3maxx)) FinalLm(Dc,Lm_star_1x,Lm_star_2x,Lm_star_3x)=if(Dc==1) then (Lm_star_1x) else (if(Dc==2) then (Lm_star_2x) else (Lm_star_3x)) FinalLf(Dc,Lf_star_1x,Lf_star_2x,Lf_star_3x)=if(Dc==1) then (Lf_star_1x) else (if(Dc==2) then (Lf_star_2x) else (Lf_star_3x)) # Final decision DF=FinalDec(Y1max,Y2max,Y3max) Ymax=FinalYmax(DF,Y1max,Y2max,Y3max) Lm_star=FinalLm(DF,Lm_star_1,Lm_star_2,Lm_star_3) Lf_star=FinalLm(DF,Lf_star_1,Lf_star_2,Lf_star_3) Lo_star=L-Lm_star-Lf_star aux Y_maxx=Ymax aux Lm_starr=Lm_star aux Lf_tarr=Lf_star aux Lo_starr=Lo_star ############################################### #======maintenance labor lm_g1=min(Lm_star*(ind_a1/(A1+A2)), ind_l) lm_g2=min(Lm_star*(ind_a2/(A1+A2)), ind_l) lm_o1=0 lm_o2=0 #======Actual total maintenance labor Lm=lm_g1*X1*N1+lm_g2*X2*N2 #======Actual water production Q0=infra(Lm)*S Qopt=infra(Lm_star)*S Abund=Q0/Qopt #======water allocation of upstream users temp1=max( ((p*b*(ind_a1)^(1-j-k)*(j/w)^j*(k)^(1-j))/(delta*(1-sigma*Abund)*0.5*(X1+X2))^(1-j))^(1/(1-j-k))-ind_r1, 0 ) q_o1=min(temp1, Q0/N1) Q1_steal=q_o1*(1-X1)*N1 q_g1=min(Q0*(1/(X1*N1+X2*N2)), (Q0-Q1_steal)*(1/(X1*N1))) Q1_rule=q_g1*X1*N1 #======water left to downstream Q1=Q1_rule+Q1_steal Q2=max(Q0-Q1,0) #======water allocation of downstream users temp2=max(((p*b*(ind_a1)^(1-j-k)*(j/w)^j*(k)^(1-j))/(delta*(1-sigma*Abund)*0.5*(X1+X2))^(1-j))^(1/(1-j-k))-ind_r2, 0) q_o2=min(temp2, Q2/N2) Q2_steal=q_o2*(1-X2)*N2 q_g2=max(Q2-Q2_steal, 0)*(1/(X2*N2)) Q2_rule=q_g2*X2*N2 #======labor allocation of group conformist lf_g1=min((p*b*j*((q_g1+ind_r1)^k)*(ind_a1^(1-j-k))/w)^(1/(1-j)), ind_l-lm_g1) le_g1=max(ind_l-lf_g1-lm_g1, 0) lf_g2=min((p*b*j*((q_g2+ind_r2)^k)*(ind_a2^(1-j-k))/w)^(1/(1-j)), ind_l-lm_g2) le_g2=max(ind_l-lf_g2-lm_g2, 0) #======labor allocation of opportunists lf_o1=min((p*b*j*((q_o1+ind_r1)^k)*(ind_a1^(1-j-k))/w)^(1/(1-j)), ind_l-lm_o1) le_o1=max(ind_l-lf_o1-lm_o1, 0) lf_o2=min((p*b*j*((q_o2+ind_r2)^k)*(ind_a2^(1-j-k))/w)^(1/(1-j)), ind_l-lm_o2) le_o2=max(ind_l-lf_o2-lm_o2, 0) #======aggregate LF and LO Lf=(lf_g1*X1+lf_o1*(1-X1))*N1+(lf_g2*X2+lf_o2*(1-X2))*N2 Le=(le_g1*X1+le_o1*(1-X1))*N1+(le_g2*X2+le_o2*(1-X2))*N2 #======upstream user income f_g1=p*b*((lf_g1)^j)*((q_g1+ind_r1)^k)*(ind_a1^(1-j-k)) e_g1=w*le_g1 f_o1=p*b*((lf_o1)^j)*((q_o1+ind_r1)^k)*(ind_a1^(1-j-k)) e_o1=w*le_o1 #======upstream downstream income f_g2=p*b*((lf_g2)^j)*((q_g2+ind_r2)^k)*(ind_a2^(1-j-k)) e_g2=w*le_g2 f_o2=p*b*((lf_o2)^j)*((q_o2+ind_r2)^k)*(ind_a2^(1-j-k)) e_o2=w*le_o2 #======net utility pi_g1=max(f_g1+e_g1-(gamma_s*(1-X1)+gamma_o*(1-X2)), 0) pi_o1=max(f_o1+e_o1-delta*q_o1*(1-sigma*Abund)*0.5*(X1+X2), 0) pi_g2=max(f_g2+e_g2-(gamma_s*(1-X2)+gamma_o*(1-X1)), 0) pi_o2=max(f_o2+e_o2-delta*q_o2*(1-sigma*Abund)*0.5*(X1+X2), 0) #======right hand sides pi_1=pi_g1*X1+pi_o1*(1-X1) pi_2=pi_g2*X2+pi_o2*(1-X2) dX1/dt=X1*(pi_g1-pi_1) dX2/dt=X2*(pi_g2-pi_2) #======aux aux Q0x=Q0 aux v1_avg=(f_g1+e_g1)*X1+(f_o1+e_o1)*(1-X1) aux v2_avg=(f_g2+e_g2)*X2+(f_o2+e_o2)*(1-X2) aux v1_avgx=(f_g1)*X1+(f_o1)*(1-X1) aux v2_avgx=(f_g2)*X2+(f_o2)*(1-X2) aux q_g1x=q_g1 aux q_o1x=q_o1 aux q_g2x=q_g2 aux q_o2x=q_o2 aux q_o1_ax=temp1 aux q_o1_bx=q_o1 aux q_o2_ax=temp2 aux q_o2_bx=q_o2 aux Lfx=Lf aux Lmx=Lm aux Lex=Le aux Lx=Lf+Lm+Le aux lf_g1x=lf_g1 aux lm_g1x=lm_g1 aux le_g1x=le_g1 aux lf_o1x=lf_o1 aux lm_o1x=lm_o1 aux le_o1x=le_o1 aux lf_g2x=lf_g2 aux lm_g2x=lm_g2 aux le_g2x=le_g2 aux lf_o2x=lf_o2 aux lm_o2x=lm_o2 aux le_o2x=le_o2 aux f_g1x=f_g1 aux f_o1x=f_o1 aux f_g2x=f_g2 aux f_o2x=f_o2 aux e_g1x=e_g1 aux e_o1x=e_o1 aux e_g2x=e_g2 aux e_o2x=e_o2 aux pi_g1x=pi_g1 aux pi_o1x=pi_o1 aux pi_g2x=pi_g2 aux pi_o2x=pi_o2 aux pi_1x=pi_1 aux pi_2x=pi_2 yield=(1/p)*(f_g1*((X1)*N1)+f_o1*((1-X1)*N1)+f_g2*((X2)*N2)+f_o2*((1-X2)*N2)) ni=((f_g1+e_g1)*((X1)*N1)+(f_o1+e_o1)*((1-X1)*N1)+(f_g2+e_g2)*((X2)*N2)+(f_o2+e_o2)*((1-X2)*N2)) opt_yield=p*b*(L-Lm_star)^j*(Imax*S)^k*(A1+A2)^(1-j-k) aux yieldx=yield aux yieldxx=yield/opt_yield aux nix=ni s_f_1=(1/p)*(f_g1*((X1)*N1)+f_o1*((1-X1)*N1))/yield s_f_2=(1/p)*(f_g2*((X2)*N2)+f_o2*((1-X2)*N2))/yield gini_yield=abs(s_f_1-s_f_2)/(2*1^2*(s_f_1+s_f_2)/2) aux gini_yieldx=gini_yield s_ni_1=((f_g1+e_g1)*((X1)*N1)+(f_o1+e_o1)*((1-X1)*N1))/ni s_ni_2=((f_g2+e_g2)*((X2)*N2)+(f_o2+e_o2)*((1-X2)*N2))/ni gini_ni=abs(s_ni_1-s_ni_2)/(2*1^2*(s_ni_1+s_ni_2)/2) aux gini_nix=gini_ni aux Q0xx=Q0/100 #=============initial data init X1=0.5, X2=0.5 @ bounds=10000 @ xp=X1, yp=X2 @ xlo=0, xhi=1, ylo=0, yhi=1 @ total=1000 @ dt=1 @ nmesh=90 @ XNC=1 @ YNC=7 @ MAXSTORE=1000000 done #10=purple, 9=blue, 8=light blue, 7=green 6=light green, 5=gold, 4=orange, 3=darker orange, 2=red, 1=dark red, 0=black
Vallury S, Arizona State University.
Effect of infrastructure design on commons dilemmas in social−ecological system dynamics. Proceedings of the National Academy of Sciences. 112(43):13212.
. 2015.