% Code for the plotting of G and H functions % % The following are the computed values of utility and time of the first % departure where there are full cycles. The first set are for the % congested region and the second is for the hypercongested region. % % Value for u0 = 24.376280000001277 @ t_up = 2.460644400008765 % Value for u1 = 24.108100000000000 @ t_up = 2.278937200000062, 2.276374999998793 % Value for u2 = 23.317700000001839 @ t_up = 1.971289999997988 2.690576 s % % Value for u2_1 = 15.681679999993122 @ t_up = 0.919738999994365 % Value for u1_1 = -4.429400000000165 @ t_up = 0.225595599999814 % Value for u0_1 = -66.127459999881751 @ t_up = 0.008035700000128 % % Shorthand Notation for cycle code % % P -> one partial cycle % FP -> one full cycle and one partial cycle % FFP -> two full cycles and one partial cycle % FFFP -> three full cycles and one partial cycle % close all; clear all; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Parameters % Logarithmic Utility Parameters a0 = 15; a1 = .5; b0 = 18; % Other parameters k = 6; % Latest arrival time t = 6 => 11 AM L = 4; % Distance traveled by cars vf = 15; % free flow velocity kj = 10^6; % jam density %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Minimum and Maximum Plot Sizes for time and utility t_min = 0.008; t_max = 2.65; u_min = -100; u_max = 30; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Method 1 for plotting G and H functions (slow, but easier to get higher % accuracy) using contour plots %{ [t,u] = meshgrid(t_min:.001:t_max,u_min:.01:u_max); %[t,u] = meshgrid(-5:.001:0,0:.001:5); T = k - t - exp((1/b0)*(u - a0*log(a1*t))); TT = t+T; TTT = (TT) + k - TT - exp((1/b0)*(u - a0*log(a1*TT))); TTT = real(TTT); TTTT = (TTT) + k - TTT - exp((1/b0)*(u - a0*log(a1*TTT))); TTTT = real(TTTT); TTTTT = (TTTT) + k - TTT - exp((1/b0)*(u - a0*log(a1*TTT))); TTTTT = real(TTTTT); G = ((a0)/(b0))*(1./t).*exp((1/b0)*(u - a0*log(a1*t))).*((a0)/(b0)).*(1./(t+T)).*exp((1/b0)*(u - a0*log(a1*(t+T)))) - 1; G1 = ((a0)/(b0))*(1./t).*exp((1/b0).*(u - a0.*log(a1*t))).*((a0)/(b0)).*(1./(t+T)).*exp((1/b0).*(u - a0*log(a1*(t+T)))).*((a0)/(b0)).*(1./(TTT)).*exp((1/b0).*(u - a0*log(a1*(TTT)))) - 1; G2 = ((a0)/(b0))*(1./t).*exp((1/b0)*(u - a0*log(a1*t))).*((a0)/(b0)).*(1./(t+T)).*exp((1/b0)*(u - a0*log(a1*(t+T)))).*((a0)/(b0)).*(1./(TTT)).*exp((1/b0)*(u - a0*log(a1*(TTT)))).*((a0)/(b0)).*(1./(TTTT)).*exp((1/b0)*(u - a0*log(a1*(TTTT)))) - 1; G3 = ((a0)/(b0))*(1./t).*exp((1/b0)*(u - a0*log(a1*t))).*((a0)/(b0)).*(1./(t+T)).*exp((1/b0)*(u - a0*log(a1*(t+T)))).*((a0)/(b0)).*(1./(TTT)).*exp((1/b0)*(u - a0*log(a1*(TTT)))).*((a0)/(b0)).*(1./(TTTT)).*exp((1/b0)*(u - a0*log(a1*(TTTT)))).*((a0)/(b0)).*(1./(TTTTT)).*exp((1/b0)*(u - a0*log(a1*(TTTTT)))) - 1; G = real(G); G1 = real(G1); G2 = real(G2); G3 = real(G3); %G2 = G2(G2 == real(G2)); %G3 = G3(G3 == real(G3)); %G2 = reshape(G2,[9401,265]); %G3 = reshape(G3,[9401,265]); %{ G = G(G == real(G)); G1 = G1(G1 == real(G1)); G2 = G2(G2 == real(G2)); G3 = G3(G3 == real(G3)); %} %} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Imported values of G that were run a priori, for faster computation G1_times_temp = load('times_G1_9000_ext.mat'); % Working G1_utils_temp = load('utils_G1_9000_ext.mat'); G1_tstar_temp = load('tstars_G1_9000_ext.mat'); G1_times = G1_times_temp.times_G1; G1_utils = G1_utils_temp.utils_G1; G1_stars = G1_tstar_temp.comp_time1; % Computation of Exact solution for one partial cycle case Tprime_und = -1 + ((a0)/(b0))*(1./G1_times).*exp((1/b0)*(G1_utils - a0*log(a1*G1_times))); Tprime_sta = -1 + ((a0)/(b0))*(1./G1_stars).*exp((1/b0)*(G1_utils - a0*log(a1*G1_stars))); CC = -1./(1+(Tprime_und)); NN = (1./(2+Tprime_sta)) + CC.*((1+Tprime_sta)/(2+Tprime_sta)); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Logarithmic Utility Fucntion Equilibrium Values from Algorithm for each % number of cycles % P data p_u_temp_L4 = load('p_u_data.mat'); p_t_temp_L4 = load('p_t_data.mat'); p_T_temp_L4 = load('p_T_data.mat'); p_Tp_temp_L4 = load('p_Tp_data.mat'); p_N_temp_L4 = load('p_N_data.mat'); p_u_L4 = p_u_temp_L4.u_matrix; p_t_L4 = p_t_temp_L4.t_matrix; p_T_L4 = p_T_temp_L4.T_matrix; p_Tp_L4 = p_Tp_temp_L4.Tp_matrix; p_N_L4 = p_N_temp_L4.N_com_matrix; % FP Data fp_u_temp_L4 = load('fp_u_data.mat'); fp_t_temp_L4 = load('fp_t_data.mat'); fp_T_temp_L4 = load('fp_T_data.mat'); fp_Tp_temp_L4 = load('fp_Tp_data.mat'); fp_N_temp_L4 = load('fp_N_data.mat'); fp_u_L4 = fp_u_temp_L4.u_matrix; fp_t_L4 = fp_t_temp_L4.t_matrix; fp_T_L4 = fp_T_temp_L4.T_matrix; fp_Tp_L4 = fp_Tp_temp_L4.Tp_matrix; fp_N_L4 = fp_N_temp_L4.N_com_matrix; % FFP Data fpp_u_temp_L4 = load('fpp_u_data.mat'); fpp_t_temp_L4 = load('fpp_t_data.mat'); %fpp_T_temp_L4 = load('fpp_T_data.mat'); %fpp_Tp_temp_L4 = load('fpp_Tp_data.mat'); fpp_N_temp_L4 = load('fpp_N_data.mat'); fpp_u_L4 = fpp_u_temp_L4.u_matrix; fpp_t_L4 = fpp_t_temp_L4.t_matrix; %fpp_T_L4 = fpp_T_temp_L4.T_matrix; %fpp_Tp_L4 = fpp_Tp_temp_L4.Tp_matrix; fpp_N_L4 = fpp_N_temp_L4.N_com_matrix; % FFFP Data fffp1_u_temp_L4 = load('fffp1_u_data.mat'); fffp1_t_temp_L4 = load('fffp1_t_data.mat'); fffp1_T_temp_L4 = load('fffp1_T_data.mat'); fffp1_Tp_temp_L4 = load('fffp1_Tp_data.mat'); fffp1_N_temp_L4 = load('fffp1_N_data.mat'); fffp1_u_L4 = fffp1_u_temp_L4.u_matrix; fffp1_t_L4 = fffp1_t_temp_L4.t_matrix; fffp1_T_L4 = fffp1_T_temp_L4.T_matrix; fffp1_Tp_L4 = fffp1_Tp_temp_L4.Tp_matrix; fffp1_N_L4 = fffp1_N_temp_L4.N_com_matrix; % t = 1.065421061749116 % u = 17.393399999992759 fffp2_u_temp_L4 = load('fffp2_u_data.mat'); fffp2_t_temp_L4 = load('fffp2_t_data.mat'); fffp2_T_temp_L4 = load('fffp2_T_data.mat'); fffp2_Tp_temp_L4 = load('fffp2_Tp_data.mat'); fffp2_N_temp_L4 = load('fffp2_N_data.mat'); fffp2_u_L4 = fffp2_u_temp_L4.u_matrix; fffp2_t_L4 = fffp2_t_temp_L4.t_matrix; fffp2_T_L4 = fffp2_T_temp_L4.T_matrix; fffp2_Tp_L4 = fffp2_Tp_temp_L4.Tp_matrix; fffp2_N_L4 = fffp2_N_temp_L4.N_com_matrix; % t = 1.018865282924472 % u = 16.880399999992132 fffp3_u_temp_L4 = load('fffp3_u_data.mat'); fffp3_t_temp_L4 = load('fffp3_t_data.mat'); fffp3_T_temp_L4 = load('fffp3_T_data.mat'); fffp3_Tp_temp_L4 = load('fffp3_Tp_data.mat'); fffp3_N_temp_L4 = load('fffp3_N_data.mat'); fffp3_u_L4 = fffp3_u_temp_L4.u_matrix; fffp3_t_L4 = fffp3_t_temp_L4.t_matrix; fffp3_T_L4 = fffp3_T_temp_L4.T_matrix; fffp3_Tp_L4 = fffp3_Tp_temp_L4.Tp_matrix; fffp3_N_L4 = fffp3_N_temp_L4.N_com_matrix; % t = 0.919804957992266 % u = 15.682399999991233 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Backward Bending Segments % FFP data ffpb_u_temp_L4 = load('ffpb_u_data.mat'); ffpb_t_temp_L4 = load('ffpb_t_data.mat'); %ffpb_T_temp_L4 = load('ffpb_T_data.mat'); %ffpb_Tp_temp_L4 = load('ffpb_Tp_data.mat'); ffpb_N_temp_L4 = load('ffpb_N_data.mat'); ffpb_u_L4 = ffpb_u_temp_L4.u_matrix; ffpb_t_L4 = ffpb_t_temp_L4.t_matrix; %ffpb_T_L4 = ffpb_T_temp_L4.T_matrix; %ffpb_Tp_L4 = ffpb_Tp_temp_L4.Tp_matrix; ffpb_N_L4 = ffpb_N_temp_L4.N_com_matrix; % u = 13.685500000001106 % t = 0.801390797988537 ffpb1_u_temp_L4 = load('ffpb1_u_data.mat'); ffpb1_t_temp_L4 = load('ffpb1_t_data.mat'); %ffpb1_T_temp_L4 = load('ffpb1_T_data.mat'); %ffpb1_Tp_temp_L4 = load('ffpb1_Tp_data.mat'); ffpb1_N_temp_L4 = load('ffpb1_N_data.mat'); ffpb1_u_L4 = ffpb1_u_temp_L4.u_matrix; ffpb1_t_L4 = ffpb1_t_temp_L4.t_matrix; %ffpb1_T_L4 = ffpb1_T_temp_L4.T_matrix; %ffpb1_Tp_L4 = ffpb1_Tp_temp_L4.Tp_matrix; ffpb1_N_L4 = ffpb1_N_temp_L4.N_com_matrix; % u = 10.629500000002800 % t = 0.646174220709072 ffpb2_u_temp_L4 = load('ffpb2_u_data.mat'); ffpb2_t_temp_L4 = load('ffpb2_t_data.mat'); %ffpb1_T_temp_L4 = load('ffpb1_T_data.mat'); %ffpb1_Tp_temp_L4 = load('ffpb1_Tp_data.mat'); ffpb2_N_temp_L4 = load('ffpb2_N_data.mat'); ffpb2_u_L4 = ffpb2_u_temp_L4.u_matrix; ffpb2_t_L4 = ffpb2_t_temp_L4.t_matrix; %ffpb1_T_L4 = ffpb1_T_temp_L4.T_matrix; %ffpb1_Tp_L4 = ffpb1_Tp_temp_L4.Tp_matrix; ffpb2_N_L4 = ffpb2_N_temp_L4.N_com_matrix; % u = 4.112500000002958 % t = 0.408178934002925 ffpb3_u_temp_L4 = load('ffpb3_u_data.mat'); ffpb3_t_temp_L4 = load('ffpb3_t_data.mat'); %ffpb1_T_temp_L4 = load('ffpb1_T_data.mat'); %ffpb1_Tp_temp_L4 = load('ffpb1_Tp_data.mat'); ffpb3_N_temp_L4 = load('ffpb3_N_data.mat'); ffpb3_u_L4 = ffpb3_u_temp_L4.u_matrix; ffpb3_t_L4 = ffpb3_t_temp_L4.t_matrix; %ffpb1_T_L4 = ffpb1_T_temp_L4.T_matrix; %ffpb1_Tp_L4 = ffpb1_Tp_temp_L4.Tp_matrix; ffpb3_N_L4 = ffpb3_N_temp_L4.N_com_matrix; % u = 1.715500000003173 % t = 0.345197889160610 ffpb4_u_temp_L4 = load('ffpb4_u_data.mat'); ffpb4_t_temp_L4 = load('ffpb4_t_data.mat'); %ffpb1_T_temp_L4 = load('ffpb1_T_data.mat'); %ffpb1_Tp_temp_L4 = load('ffpb1_Tp_data.mat'); ffpb4_N_temp_L4 = load('ffpb4_N_data.mat'); ffpb4_u_L4 = ffpb4_u_temp_L4.u_matrix; ffpb4_t_L4 = ffpb4_t_temp_L4.t_matrix; %ffpb1_T_L4 = ffpb1_T_temp_L4.T_matrix; %ffpb1_Tp_L4 = ffpb1_Tp_temp_L4.Tp_matrix; ffpb4_N_L4 = ffpb4_N_temp_L4.N_com_matrix; % u = -4.419499999996559 % t = 0.225773063374777 % FP data fpb_u_temp_L4 = load('fpb_u_data.mat'); fpb_t_temp_L4 = load('fpb_t_data.mat'); fpb_T_temp_L4 = load('fpb_T_data.mat'); fpb_Tp_temp_L4 = load('fpb_Tp_data.mat'); fpb_N_temp_L4 = load('fpb_N_data.mat'); fpb_u_L4 = fpb_u_temp_L4.u_matrix; fpb_t_L4 = fpb_t_temp_L4.t_matrix; fpb_T_L4 = fpb_T_temp_L4.T_matrix; fpb_Tp_L4 = fpb_Tp_temp_L4.Tp_matrix; fpb_N_L4 = fpb_N_temp_L4.N_com_matrix; fpb1_u_temp_L4 = load('fpb1_u_data.mat'); fpb1_t_temp_L4 = load('fpb1_t_data.mat'); fpb1_T_temp_L4 = load('fpb1_T_data.mat'); fpb1_Tp_temp_L4 = load('fpb1_Tp_data.mat'); fpb1_N_temp_L4 = load('fpb1_N_data.mat'); fpb1_u_L4 = fpb1_u_temp_L4.u_matrix; fpb1_t_L4 = fpb1_t_temp_L4.t_matrix; fpb1_T_L4 = fpb1_T_temp_L4.T_matrix; fpb1_Tp_L4 = fpb1_Tp_temp_L4.Tp_matrix; fpb1_N_L4 = fpb1_N_temp_L4.N_com_matrix; % u = -22.229600000000726 % t_up = 0.092135919146068 fpb2_u_temp_L4 = load('fpb2_u_data.mat'); fpb2_t_temp_L4 = load('fpb2_t_data.mat'); fpb2_T_temp_L4 = load('fpb2_T_data.mat'); fpb2_Tp_temp_L4 = load('fpb2_Tp_data.mat'); fpb2_N_temp_L4 = load('fpb2_N_data.mat'); fpb2_u_L4 = fpb2_u_temp_L4.u_matrix; fpb2_t_L4 = fpb2_t_temp_L4.t_matrix; fpb2_T_L4 = fpb2_T_temp_L4.T_matrix; fpb2_Tp_L4 = fpb2_Tp_temp_L4.Tp_matrix; fpb2_N_L4 = fpb2_N_temp_L4.N_com_matrix; % u = -28.369600000001686 % t_up = 0.065031652375699 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% al0 = 15; al1 = .5; bl0 = 18; % value of nu nu_m = 25.993549906396904; % Plots for Exact G(t,u) with H(t,u,L=4) plotted figure % Plotting Break points bpH = plot(2.460644400008765,24.376280000001277,'ko','MarkerSize',10), hold on plot(2.278937200000062,24.108100000000000,'ko','MarkerSize',10), hold on plot(1.971289999997988,23.317700000001839,'ko','MarkerSize',10), hold on plot(0.919738999994365,15.681679999993122,'ko','MarkerSize',10), hold on plot(0.225595599999814,-4.429400000000165,'ko','MarkerSize',10), hold on plot(0.008035700000128,-66.127459999881751,'ko','MarkerSize',10), hold on plot(2.460644400008765,24.376280000001277,'ko','MarkerSize',12), hold on plot(2.278937200000062,24.108100000000000,'ko','MarkerSize',12), hold on plot(1.971289999997988,23.317700000001839,'ko','MarkerSize',12), hold on plot(0.919738999994365,15.681679999993122,'ko','MarkerSize',12), hold on plot(0.225595599999814,-4.429400000000165,'ko','MarkerSize',12), hold on plot(0.008035700000128,-66.127459999881751,'ko','MarkerSize',12), hold on % First Segment of H functions H = plot(p_t_L4,p_u_L4,'b','LineWidth',1), hold on plot(fp_t_L4,fp_u_L4,'b','LineWidth',1), hold on plot(fpp_t_L4,fpp_u_L4,'b','LineWidth',1), hold on plot(fffp1_t_L4,fffp1_u_L4,'b','LineWidth',1), hold on plot(fffp2_t_L4,fffp2_u_L4,'b','LineWidth',1), hold on plot(fffp3_t_L4,fffp3_u_L4,'b','LineWidth',1), hold on % Second Segment of H functions plot(ffpb_t_L4,ffpb_u_L4,'b','LineWidth',1), hold on plot(ffpb1_t_L4,ffpb1_u_L4,'b','LineWidth',1), hold on plot(ffpb2_t_L4,ffpb2_u_L4,'b','LineWidth',1), hold on plot(ffpb3_t_L4,ffpb3_u_L4,'b','LineWidth',1), hold on plot(ffpb4_t_L4,ffpb4_u_L4,'b','LineWidth',1), hold on plot(fpb_t_L4,fpb_u_L4,'b','LineWidth',1), hold on plot(fpb1_t_L4,fpb1_u_L4,'b','LineWidth',1), hold on plot(fpb2_t_L4,fpb2_u_L4,'b','LineWidth',1), hold on % Mesh Density (for L_along_Gi is 5000) G0 = fimplicit(@(x,y) ( (al0/bl0)*(al1^(-al0/bl0))*exp(y./bl0).*(x.^((-al0/bl0)-1)) ) - 1, [-6 6 -100 30],'-bs','MeshDensity',1000,'LineWidth',1,'MarkerSize',10) G1 = fimplicit(@(x,y) ( (al0/bl0)*(al1^(-al0/bl0))*exp(y./bl0).*(x.^((-al0/bl0)-1)) ) .* ( (al0/bl0)*(al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*(x.^(-al0/bl0))).^((-al0/bl0)-1)) ) - 1, [-6 6 -100 30],'-md','MeshDensity',200,'LineWidth',1,'MarkerSize',10), hold on % G_2 G2 = fimplicit(@(x,y) ( (al0/bl0)*(al1^(-al0/bl0))*exp(y./bl0).*(x.^((-al0/bl0)-1)) ).*( (al0/bl0)*(al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*(x.^(-al0/bl0))).^((-al0/bl0)-1)) ).*( (al0/bl0)*(al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*(x.^(-al0/bl0))).^(-al0/bl0))).^((-al0/bl0)-1)) ) - 1, [-6 6 -100 30],'-rx','MeshDensity',200,'LineWidth',1,'MarkerSize',10), hold on G3 = fimplicit(@(x,y) ( (al0/bl0)*(al1^(-al0/bl0))*exp(y./bl0).*(x.^((-al0/bl0)-1)) ).*( (al0/bl0)*(al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*(x.^(-al0/bl0))).^((-al0/bl0)-1)) ).*( (al0/bl0)*(al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*(x.^(-al0/bl0))).^(-al0/bl0))).^((-al0/bl0)-1)) ).*( (al0/bl0)*(al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*(x.^(-al0/bl0))).^(-al0/bl0))).^(-al0/bl0))).^((-al0/bl0)-1)) ) - 1, [-6 6 -50 25],'-cs','MeshDensity',200,'LineWidth',1,'MarkerSize',10), hold on G4 = fimplicit(@(x,y) ( (al0/bl0)*(al1^(-al0/bl0))*exp(y./bl0).*(x.^((-al0/bl0)-1)) ).*( (al0/bl0)*(al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*(x.^(-al0/bl0))).^((-al0/bl0)-1)) ).*( (al0/bl0)*(al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*(x.^(-al0/bl0))).^(-al0/bl0))).^((-al0/bl0)-1)) ).*( (al0/bl0)*(al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*(x.^(-al0/bl0))).^(-al0/bl0))).^(-al0/bl0))).^((-al0/bl0)-1)) )* ( (al0/bl0)*(al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).* ... ((k - (al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).*((k - (al1^(-al0/bl0))*exp(y./bl0).* ... (x.^(-al0/bl0))).^(-al0/bl0))).^(-al0/bl0))).^(-al0/bl0))).^(-al0/bl0))).^(-al0/bl0))).^(-al0/bl0))).^((-al0/bl0)-1)) ) - 1, [-6 6 -50 25],'-gv','MeshDensity',800,'LineWidth',1,'MarkerSize',10), hold on % Contour plots for G %contour(t,u,G,[0 0],'k'), hold on %contour(t,u,G1,[0 0],'m'), hold on %contour(t,u,G2,[0 0],'g'), hold on %contour(t,u,G3,[0 0],'c'), hold on %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% lg = legend([bpH,H,G1,G2,G3,G4],'Break Points','$H(\underline{t},\underline{u},\underline{L})$','$G_1(\underline{t},\underline{u})$','$G_2(\underline{t},\underline{u})$','$G_3(\underline{t},\underline{u})$','$G_4(\underline{t},\underline{u})$'); set(lg,'Interpreter','latex'); set(gca,'fontsize',50,'TickLabelInterpreter', 'latex');xlabel('t-underbar') xlabel('time $\underline{t}$','Interpreter','latex') ylabel('Utility $\underline{u}$','Interpreter','latex') title('Functions $G_i(\underline{t},\underline{u})$ and $H(\underline{t},\underline{u},L)$','Interpreter','latex') ylim([-68 26]) xlim([-0.1 3]) grid on % Number of Commuters vs. Utility figure bpN = plot(9.409060516508245e+04,24.376280000001277,'ko','MarkerSize',15), hold on plot(4.015577381961486e+05,24.108100000000000,'ko','MarkerSize',15), hold on plot(1.172702057909408e+06,23.317700000001839,'ko','MarkerSize',15), hold on plot(2.372942129818678e+06,15.681679999993122,'ko','MarkerSize',15), hold on plot(1.877880212005743e+06,-4.429400000000165,'ko','MarkerSize',15), hold on plot(9.961861701928501e+05,-66.127459999881751,'ko','MarkerSize',15), hold on plot(9.409060516508245e+04,24.376280000001277,'ko','MarkerSize',17), hold on plot(4.015577381961486e+05,24.108100000000000,'ko','MarkerSize',17), hold on plot(1.172702057909408e+06,23.317700000001839,'ko','MarkerSize',17), hold on plot(2.372942129818678e+06,15.681679999993122,'ko','MarkerSize',17), hold on plot(1.877880212005743e+06,-4.429400000000165,'ko','MarkerSize',17), hold on plot(9.961861701928501e+05,-66.127459999881751,'ko','MarkerSize',17), hold on % First Segments %plot(N_p_L4,u_p_L4,'r'), hold on N1p = plot(p_N_L4,p_u_L4,'b','LineWidth',3), hold on plot(fp_N_L4,fp_u_L4,'b','LineWidth',3), hold on plot(fpp_N_L4,fpp_u_L4,'b','LineWidth',3), hold on plot(fffp1_N_L4,fffp1_u_L4,'b','LineWidth',3), hold on plot(fffp2_N_L4,fffp2_u_L4,'b','LineWidth',3), hold on plot(fffp3_N_L4,fffp3_u_L4,'b','LineWidth',3), hold on % Second Segments plot(ffpb_N_L4,ffpb_u_L4,'b','LineWidth',3), hold on plot(ffpb1_N_L4,ffpb1_u_L4,'b','LineWidth',3), hold on plot(ffpb2_N_L4,ffpb2_u_L4,'b','LineWidth',3), hold on plot(ffpb3_N_L4,ffpb3_u_L4,'b','LineWidth',3), hold on plot(ffpb4_N_L4,ffpb4_u_L4,'b','LineWidth',3), hold on plot(fpb_N_L4,fpb_u_L4,'b','LineWidth',3), hold on plot(fpb1_N_L4,fpb1_u_L4,'b','LineWidth',3), hold on plot(fpb2_N_L4,fpb2_u_L4,'b','LineWidth',3), hold on plot(NN,G1_utils,'r','LineWidth',3), hold on %plot(N_TEST1_L4,u_TEST1_L4,'r'), hold on %plot(N_TEST2_L4,u_TEST2_L4,'r'), hold on %%%plot(N_TEST3_L4,u_TEST3_L4,'r','LineWidth',3), hold on %plot(N_TEST1a_L4,u_TEST1a_L4,'r'), hold on %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% l1 = legend([bpN,N1p],'Break Points','$N(\underline{u})$'); set(l1,'Interpreter','latex'); set(gca,'fontsize',50,'TickLabelInterpreter', 'latex'); xlabel('$N$ (number of commuters)','Interpreter','latex') ylabel('$\underline{u}$ (utility level)','Interpreter','latex') title('Utility vs Number of Commuters','Interpreter','latex') %ylim([-.5*10^5 10^6]) grid on % Number of Commuters vs. User cost figure bpNc = plot(9.409060516508245e+04,nu_m -24.376280000001277,'ko','MarkerSize',15), hold on plot(4.015577381961486e+05,nu_m -24.108100000000000,'ko','MarkerSize',15), hold on plot(1.172702057909408e+06,nu_m -23.317700000001839,'ko','MarkerSize',15), hold on plot(2.372942129818678e+06,nu_m -15.681679999993122,'ko','MarkerSize',15), hold on plot(1.877880212005743e+06,nu_m --4.429400000000165,'ko','MarkerSize',15), hold on plot(9.961861701928501e+05,nu_m --66.127459999881751,'ko','MarkerSize',15), hold on plot(9.409060516508245e+04,nu_m -24.376280000001277,'ko','MarkerSize',17), hold on plot(4.015577381961486e+05,nu_m -24.108100000000000,'ko','MarkerSize',17), hold on plot(1.172702057909408e+06,nu_m -23.317700000001839,'ko','MarkerSize',17), hold on plot(2.372942129818678e+06,nu_m -15.681679999993122,'ko','MarkerSize',17), hold on plot(1.877880212005743e+06,nu_m --4.429400000000165,'ko','MarkerSize',17), hold on plot(9.961861701928501e+05,nu_m --66.127459999881751,'ko','MarkerSize',17), hold on % First Segments %plot(N_p_L4,u_p_L4,'r'), hold on N1pc = plot(p_N_L4,nu_m -p_u_L4,'b','LineWidth',3), hold on plot(fp_N_L4,nu_m -fp_u_L4,'b','LineWidth',3), hold on plot(fpp_N_L4,nu_m -fpp_u_L4,'b','LineWidth',3), hold on plot(fffp1_N_L4,nu_m -fffp1_u_L4,'b','LineWidth',3), hold on plot(fffp2_N_L4,nu_m -fffp2_u_L4,'b','LineWidth',3), hold on plot(fffp3_N_L4,nu_m -fffp3_u_L4,'b','LineWidth',3), hold on % Second Segments plot(ffpb_N_L4,nu_m -ffpb_u_L4,'b','LineWidth',3), hold on plot(ffpb1_N_L4,nu_m -ffpb1_u_L4,'b','LineWidth',3), hold on plot(ffpb2_N_L4,nu_m -ffpb2_u_L4,'b','LineWidth',3), hold on plot(ffpb3_N_L4,nu_m -ffpb3_u_L4,'b','LineWidth',3), hold on plot(ffpb4_N_L4,nu_m -ffpb4_u_L4,'b','LineWidth',3), hold on plot(fpb_N_L4,nu_m -fpb_u_L4,'b','LineWidth',3), hold on plot(fpb1_N_L4,nu_m -fpb1_u_L4,'b','LineWidth',3), hold on plot(fpb2_N_L4,nu_m -fpb2_u_L4,'b','LineWidth',3), hold on %plot(N_TEST1_L4,u_TEST1_L4,'r'), hold on %plot(N_TEST2_L4,u_TEST2_L4,'r'), hold on %%%plot(N_TEST3_L4,u_TEST3_L4,'r','LineWidth',3), hold on %plot(N_TEST1a_L4,u_TEST1a_L4,'r'), hold on %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% l1 = legend([bpNc,N1pc],'Break Points','$c(N)$'); set(l1,'Interpreter','latex'); set(gca,'fontsize',50,'TickLabelInterpreter', 'latex'); xlabel('$N$ (number of commuters)','Interpreter','latex') ylabel('$c$ (user-cost)','Interpreter','latex') %title('User-Cost vs. Number of Commuters (Log Utility Function)','Interpreter','latex') ylim([0 94]) grid on % N vs time figure % bpN = plot(9.409060516508245e+04,24.376280000001277,'ko','MarkerSize',15), hold on % plot(4.015577381961486e+05,24.108100000000000,'ko','MarkerSize',15), hold on % plot(1.172702057909408e+06,23.317700000001839,'ko','MarkerSize',15), hold on % plot(2.372942129818678e+06,15.681679999993122,'ko','MarkerSize',15), hold on % plot(1.877880212005743e+06,-4.429400000000165,'ko','MarkerSize',15), hold on % plot(9.961861701928501e+05,-66.127459999881751,'ko','MarkerSize',15), hold on % plot(9.409060516508245e+04,24.376280000001277,'ko','MarkerSize',17), hold on % plot(4.015577381961486e+05,24.108100000000000,'ko','MarkerSize',17), hold on % plot(1.172702057909408e+06,23.317700000001839,'ko','MarkerSize',17), hold on % plot(2.372942129818678e+06,15.681679999993122,'ko','MarkerSize',17), hold on % plot(1.877880212005743e+06,-4.429400000000165,'ko','MarkerSize',17), hold on % plot(9.961861701928501e+05,-66.127459999881751,'ko','MarkerSize',17), hold on % First Segments %plot(N_p_L4,u_p_L4,'r'), hold on N1p = plot(p_t_L4,p_N_L4,'b','LineWidth',3), hold on plot(fp_t_L4,fp_N_L4,'b','LineWidth',3), hold on plot(fpp_t_L4,fpp_N_L4,'b','LineWidth',3), hold on plot(fffp1_t_L4,fffp1_N_L4,'b','LineWidth',3), hold on plot(fffp2_t_L4,fffp2_N_L4,'b','LineWidth',3), hold on plot(fffp3_t_L4,fffp3_N_L4,'b','LineWidth',3), hold on % Second Segments plot(ffpb_t_L4,ffpb_N_L4,'b','LineWidth',3), hold on plot(ffpb1_t_L4,ffpb1_N_L4,'b','LineWidth',3), hold on plot(ffpb2_t_L4,ffpb2_N_L4,'b','LineWidth',3), hold on plot(ffpb3_t_L4,ffpb3_N_L4,'b','LineWidth',3), hold on plot(ffpb4_t_L4,ffpb4_N_L4,'b','LineWidth',3), hold on plot(fpb_t_L4,fpb_N_L4,'b','LineWidth',3), hold on plot(fpb1_t_L4,fpb1_N_L4,'b','LineWidth',3), hold on plot(fpb2_t_L4,fpb2_N_L4,'b','LineWidth',3), hold on %plot(N_TEST1_L4,u_TEST1_L4,'r'), hold on %plot(N_TEST2_L4,u_TEST2_L4,'r'), hold on %%%plot(N_TEST3_L4,u_TEST3_L4,'r','LineWidth',3), hold on %plot(N_TEST1a_L4,u_TEST1a_L4,'r'), hold on %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% l1 = legend([bpN,N1p],'Break Points','$N(\underline{u})$'); set(l1,'Interpreter','latex'); set(gca,'fontsize',50,'TickLabelInterpreter', 'latex'); xlabel('$\underline{t}$ time','Interpreter','latex') ylabel('$N$ number of commuters','Interpreter','latex') title('Utility vs Number of Commuters','Interpreter','latex') %ylim([-.5*10^5 10^6]) grid on