clear % octave -q --no-gui % wear_loop vhe_model='calc01'; % Tsim input data model to be used in this wear prediction loop vhe_model='calc01'; % R600 stiff bogie %vhe_model='calc02'; % Tangent track; No track irr.; With traction iloop= 0; % Start wheel profile # nloop= 3; % End wheel profile # string= sprintf('mkdir -p qout'); [istat, output]= system(string); if (istat ~= 0); error(output); exit; end string= sprintf('mkdir -p kpf/qout'); [istat, output]= system(string); if (istat ~= 0); error(output); exit; end ident= sprintf('%s_%3.3d', vhe_model, iloop); if (iloop==0) string= ['cp -p kpf/w_prof/outer_wheel_000.wheel kpf/w_prof/outer_wheel_', ident, '.wheel']; [istat, output]= system(string); if (istat ~= 0); error(output); exit; end string= ['cp -p kpf/w_prof/inner_wheel_000.wheel kpf/w_prof/inner_wheel_', ident, '.wheel']; [istat, output]= system(string); if (istat ~= 0); error(output); exit; end end outer_wheel_prof= sprintf('outer_wheel_%s', ident); inner_wheel_prof= sprintf('inner_wheel_%s', ident); string= ['[outer_wheel_prof_X, outer_wheel_prof_Y]= read_2col_file ("kpf/w_prof/',outer_wheel_prof,'.wheel");']; % Read the initial profiles try eval(string) catch error("Exit status from command %s\n not equal zero\n%s\n", string, lasterr ()); end string= ['[inner_wheel_prof_X, inner_wheel_prof_Y]= read_2col_file ("kpf/w_prof/',inner_wheel_prof,'.wheel");']; try eval(string) catch error("Exit status from command %s\n not equal zero\n%s\n", string, lasterr ()); end % % Wear step loop % ------------------------ while( iloop kpff/',outer_wheel_prof,'.kpff']; [istat, output]= system(string); if (istat ~= 0); error(output); exit; end % string= ['kpf kpff/',outer_wheel_prof,'.kpff > qout/',outer_wheel_prof,'.out']; %% Run program KPF outer_wheel_prof [istat, output]= system(string); if (istat ~= 0) printf(' ***ERROR*** In wear_loop.m\n'); printf(' Exit status from program KPF not equal 0 (zero)\n'); printf(' Please open and read the error message in file:\n'); printf(' kpf/qout/%s.out\n', outer_wheel_prof); exit end string= ['sed -e "s!WHEEL!w_prof/',inner_wheel_prof,'.wheel!" -e "s!RAIL!\$genkpf/../r_prof/uic60/uic60i40.rail!" kpff/Master_same_profiles.kpff > kpff/',inner_wheel_prof,'.kpff']; [istat, output]= system(string); if (istat ~= 0); error(output); exit; end % string= ['kpf kpff/',inner_wheel_prof,'.kpff > qout/',inner_wheel_prof,'.out']; %% Run program KPF inner_wheel_prof [istat, output]= system(string); if (istat ~= 0) printf(' ***ERROR*** In wear_loop.m\n'); printf(' Exit status from program KPF not equal 0 (zero)\n'); printf(' Please open and read the error message in file:\n'); printf(' kpf/qout/%s.out\n', inner_wheel_prof); exit end cd .. t_kpfr= toc; disp(['t_kpfr= ',num2str(round(t_kpfr)),' s']) disp('End kpfr }}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}') % Adjust the wear depending on contact angle disp('Beg gamma_cos }}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}') string= ['(cd kpf/kpfr; kpf_conv_kpfr2matlab ',outer_wheel_prof,'.kpfr ',outer_wheel_prof,'.m)'] [istat, output]= system(string); if (istat ~= 0) printf(' ***ERROR*** In wear_loop.m\n'); printf(' Exit status from program "kpf_conv_kpfr2matlab" not equal 0 (zero)\n'); printf(' output= %s\n', output); exit end eta= [-.040:.001:.040]; cd kpf/kpfr eval (outer_wheel_prof) cd ../.. cp1_zfn= interp1 (cp1_zfn_(:,1),cp1_zfn_(:,2), eta, "extrap"); cp2_zfn= interp1 (cp2_zfn_(:,1),cp2_zfn_(:,2), eta, "extrap"); cp3_zfn= interp1 (cp3_zfn_(:,1),cp3_zfn_(:,2), eta, "extrap"); cp1_gamfn = interp1 (cp1_gamfn_(:,1),cp1_gamfn_(:,2), eta, "extrap"); cp2_gamfn = interp1 (cp2_gamfn_(:,1),cp2_gamfn_(:,2), eta, "extrap"); cp3_gamfn = interp1 (cp3_gamfn_(:,1),cp3_gamfn_(:,2), eta, "extrap"); cp1_poswfn= interp1 (cp1_poswfn_(:,1),cp1_poswfn_(:,2), eta, "extrap"); cp2_poswfn= interp1 (cp2_poswfn_(:,1),cp2_poswfn_(:,2), eta, "extrap"); cp3_poswfn= interp1 (cp3_poswfn_(:,1),cp3_poswfn_(:,2), eta, "extrap"); % plot(eta,cp1_zfn) % axis([ -.2 .2 -.005 .005 ]) for ii= 1:length(eta) cpa_gamfn(ii) = cp1_gamfn(ii); cpa_poswfn(ii)= cp1_poswfn(ii); if (cp2_zfn(ii) < cp1_zfn(ii)) cpa_gamfn(ii) = cp2_gamfn(ii); cpa_poswfn(ii)= cp2_poswfn(ii); endif if (cp3_zfn(ii) < cp2_zfn(ii)) cpa_gamfn(ii) = cp3_gamfn(ii); cpa_poswfn(ii)= cp3_poswfn(ii); endif endfor % plot(eta,cpa_gamfn,';cpa;', eta,cp1_gamfn,';cp1;', eta,cp2_gamfn,';cp2;', eta,cp3_gamfn,';cp3;') % plot(eta,cp1_gamfn,';cp1;', eta,cp2_gamfn,';cp2;', eta,cp3_gamfn,';cp3;', eta,cpa_gamfn,';cpa;') % plot(eta,cp1_poswfn,';cp1;', eta,cp2_poswfn,';cp2;', eta,cp3_poswfn,';cp3;', eta,cpa_poswfn,';cpa;') % plot(eta,cpa_poswfn) for ii= 2:length(eta) ;% Avoid backsteps if (cpa_poswfn(ii) <= cpa_poswfn(ii-1)) cpa_poswfn(ii)= cpa_poswfn(ii-1)+eps; endif endfor posw_111l_1D= [ -4.475000E-02 -4.425000E-02 -4.375000E-02 -4.325000E-02 -4.275000E-02 -4.225000E-02 -4.175000E-02 -4.125000E-02 -4.075000E-02 -4.025000E-02 -3.975000E-02 -3.925000E-02 -3.875000E-02 -3.825000E-02 -3.775000E-02 -3.725000E-02 -3.675000E-02 -3.625000E-02 -3.575000E-02 -3.525000E-02 -3.475000E-02 -3.425000E-02 -3.375000E-02 -3.325000E-02 -3.275000E-02 -3.225000E-02 -3.175000E-02 -3.125000E-02 -3.075000E-02 -3.025000E-02 -2.975000E-02 -2.925000E-02 -2.875000E-02 -2.825000E-02 -2.775000E-02 -2.725000E-02 -2.675000E-02 -2.625000E-02 -2.575000E-02 -2.525000E-02 -2.475000E-02 -2.425000E-02 -2.375000E-02 -2.325000E-02 -2.275000E-02 -2.225000E-02 -2.175000E-02 -2.125000E-02 -2.075000E-02 -2.025000E-02 -1.975000E-02 -1.925000E-02 -1.875000E-02 -1.825000E-02 -1.775000E-02 -1.725000E-02 -1.675000E-02 -1.625000E-02 -1.575000E-02 -1.525000E-02 -1.475000E-02 -1.425000E-02 -1.375000E-02 -1.325000E-02 -1.275000E-02 -1.225000E-02 -1.175000E-02 -1.125000E-02 -1.075000E-02 -1.025000E-02 -9.750001E-03 -9.250000E-03 -8.749999E-03 -8.250002E-03 -7.750001E-03 -7.250000E-03 -6.749999E-03 -6.250001E-03 -5.750000E-03 -5.249999E-03 -4.749998E-03 -4.250001E-03 -3.750000E-03 -3.249999E-03 -2.749998E-03 -2.250001E-03 -1.750000E-03 -1.249999E-03 -7.500015E-04 -2.500005E-04 2.500005E-04 7.500015E-04 1.249999E-03 1.750000E-03 2.250001E-03 2.750002E-03 3.249999E-03 3.750000E-03 4.250001E-03 4.750002E-03 5.249999E-03 5.750000E-03 6.250001E-03 6.750003E-03 7.250000E-03 7.750001E-03 8.250002E-03 8.749999E-03 9.250000E-03 9.750001E-03 1.025000E-02 1.075000E-02 1.125000E-02 1.175000E-02 1.225000E-02 1.275000E-02 1.325000E-02 1.375000E-02 1.425000E-02 1.475000E-02 1.525000E-02 1.575000E-02 1.625000E-02 1.675000E-02 1.725000E-02 1.775000E-02 1.825000E-02 1.875000E-02 1.925000E-02 1.975000E-02 2.025000E-02 2.075000E-02 2.125000E-02 2.175000E-02 2.225000E-02 2.275001E-02 2.325000E-02 2.375000E-02 2.425000E-02 2.475000E-02 2.525000E-02 2.575000E-02 2.625000E-02 2.675000E-02 2.725000E-02 2.775000E-02 2.825000E-02 2.875000E-02 2.925000E-02 2.975000E-02 3.025001E-02 3.075000E-02 3.125000E-02 3.175000E-02 3.225000E-02 3.275000E-02 3.325000E-02 3.375001E-02 3.425000E-02 3.475000E-02 3.525000E-02 3.575000E-02 3.625000E-02 3.675000E-02 3.725000E-02 3.775001E-02 3.825000E-02 3.875000E-02 3.925000E-02 3.975000E-02 4.025000E-02 4.075000E-02 4.125001E-02 4.175000E-02 4.225000E-02 4.275000E-02 4.325000E-02 4.375000E-02 4.425000E-02 4.475001E-02]; posw_eta= interp1 (cpa_poswfn,eta, posw_111l_1D, "extrap"); % for ii= 1:length(posw_111l_1D) % if (posw_eta(ii) < -.040) posw_eta(ii)= -.040; endif % if (posw_eta(ii) > .040) posw_eta(ii)= .040; endif % endfor gamma= interp1 (eta,cpa_gamfn, posw_eta, "extrap"); ;% Pick up gamma for ii= 1:length(eta) ;% Avoid too big angles in gamma if (gamma(ii) > 1.5) cpa_poswfn(ii)= 1.5; endif endfor gamma_cos= cos(gamma); % plot(posw_eta,gamma) % plot(posw_eta,gamma_cos) % axis([ -40 40 -1 1 ]) disp('End gamma_cos }}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}') % disp('Run tsim {{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{'); tic; # awk '/jack/ {c+=1} {if(c==2){sub("jack","jill",$0)};print}' file # awk '/Wheel-rail geometry functions/ {c+=1} { if (c==2) { # printf("insert file kpf/kpfr/%s.kpfr \\\\n", outer_wheel_prof); # }}' file string= sprintf("awk '/Wheel-rail geometry functions/ {c+=1} {print; if (c==2) {++c;"); strin2= sprintf("printf(\"insert file kpf/kpfr/%s.kpfr\\n\");", outer_wheel_prof); string= strcat(string,strin2); strin2= sprintf("printf(\"insert file kpf/kpfr/%s.kpfr\\n\");", inner_wheel_prof); string= strcat(string,strin2); strin2= sprintf("printf(\"in_substruct %s [ 111 ] \\n\");", outer_wheel_prof); string= strcat(string,strin2); strin2= sprintf("printf(\"in_substruct %s [ 112 ] \\n\");", inner_wheel_prof); string= strcat(string,strin2); strin2= sprintf("printf(\"in_substruct %s [ 121 ] \\n\");", inner_wheel_prof); string= strcat(string,strin2); strin2= sprintf("printf(\"in_substruct %s [ 122 ] \\n\")}}'", outer_wheel_prof); string= strcat(string,strin2); strin2= sprintf(' runf/%s.tsimf > runf/%s.tsimf', vhe_model, ident); string= strcat(string,strin2) # strin2= sprintf("printf(\"runf/%s.tsimf > runf/%s.tsimf", vhe_model, ident); string= strcat(string,strin2); # strin2= sprintf("printf(\"insert file kpf/kpfr/%s.kpfr \\\\n\", inner_wheel_prof)"); string= strcat(string,strin2); # strin2= sprintf(' runf/%s.tsimf > runf/%s.tsimf', vhe_model, ident); string= strcat(string,strin2); # string= sprintf('sed "/Wheel-rail geometry functions/ a\\\n'); # strin2= sprintf('insert file kpf/kpfr/%s.kpfr \\\\n', outer_wheel_prof); string= strcat(string,strin2); # strin2= sprintf('insert file kpf/kpfr/%s.kpfr \\\\n', inner_wheel_prof); string= strcat(string,strin2); # strin2= sprintf('in_substruct %s [ 111 ] \\\\n', outer_wheel_prof); string= strcat(string,strin2); # strin2= sprintf('in_substruct %s [ 112 ] \\\\n', inner_wheel_prof); string= strcat(string,strin2); # strin2= sprintf('in_substruct %s [ 121 ] \\\\n', inner_wheel_prof); string= strcat(string,strin2); # strin2= sprintf('in_substruct %s [ 122 ]"', outer_wheel_prof); string= strcat(string,strin2); # strin2= sprintf(' runf/%s.tsimf > runf/%s.tsimf', vhe_model, ident); string= strcat(string,strin2); [istat, output]= system(string); if (istat ~= 0); error(output); exit; end string= sprintf('tsim -no_addarg -overwrite runf/%s.tsimf > qout/tsim_%s.out', ... ident, ident); [istat, output]= system(string); if (istat ~= 0) printf(' ***ERROR*** In wear_loop.m\n'); printf(' Exit status from program TSIM not equal 0 (zero)\n'); printf(' Please open and read the error message in file:\n'); printf(' qout/tsim_%s.out\n', ident); exit end t_tsim= toc; disp(['t_tsim= ',num2str(round(t_tsim)),' s']) disp('End tsim }}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}') % disp('Run mplot {{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{') tic; [istat, output]= system('del -q mplotr/wear_111l.print'); if (istat ~= 0); error(output); exit; end [istat, output]= system('del -q mplotr/wear_111r.print'); if (istat ~= 0); error(output); exit; end [istat, output]= system('del -q mplotr/wear_112l.print'); if (istat ~= 0); error(output); exit; end [istat, output]= system('del -q mplotr/wear_112r.print'); if (istat ~= 0); error(output); exit; end [istat, output]= system('del -q mplotr/wear_121l.print'); if (istat ~= 0); error(output); exit; end [istat, output]= system('del -q mplotr/wear_121r.print'); if (istat ~= 0); error(output); exit; end [istat, output]= system('del -q mplotr/wear_122l.print'); if (istat ~= 0); error(output); exit; end [istat, output]= system('del -q mplotr/wear_122r.print'); if (istat ~= 0); error(output); exit; end string= sprintf('mplot -no_addarg -no_interactive -no_use_MPfile mplotf/wear.mplotf %s > qout/mplot_%s.out', ident, ident); [istat, output]= system(string); if (istat ~= 0) printf(' ***ERROR*** In wear_loop.m\n'); printf(' Exit status from program MPLOT not equal 0 (zero)\n'); printf(' Please open and read the error message in file:\n'); printf(' qout/mplot_%s.out\n', ident); exit end t_mplot=toc; disp(['t_mplot= ',num2str(round(t_mplot)),' s']) disp('End mplot }}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}') % disp('Update wheel profile {{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{{') load -ascii mplotr/wear_111l.print load -ascii mplotr/wear_111r.print load -ascii mplotr/wear_112l.print load -ascii mplotr/wear_112r.print load -ascii mplotr/wear_121l.print load -ascii mplotr/wear_121r.print load -ascii mplotr/wear_122l.print load -ascii mplotr/wear_122r.print % wear_111l = importdata('mplotr/wear_111l.print', ' ', 12); % wear_111r = importdata('mplotr/wear_111r.print', ' ', 12); % wear_112l = importdata('mplotr/wear_112l.print', ' ', 12); % wear_112r = importdata('mplotr/wear_112r.print', ' ', 12); % wear_121l = importdata('mplotr/wear_121l.print', ' ', 12); % wear_121r = importdata('mplotr/wear_121r.print', ' ', 12); % wear_122l = importdata('mplotr/wear_122l.print', ' ', 12); % wear_122r = importdata('mplotr/wear_122r.print', ' ', 12); % wear_111l = wear_111l.data; % wear_111r = wear_111r.data; % wear_112l = wear_112l.data; % wear_112r = wear_112r.data; % wear_121l = wear_121l.data; % wear_121r = wear_121r.data; % wear_122l = wear_122l.data; % wear_122r = wear_122r.data; wear_outer= (wear_111l(:,2) + wear_111r(:,2) + wear_122l(:,2) + wear_122r(:,2)) /4.*1000.*100.; % Convert wear to mm wear_inner= (wear_112l(:,2) + wear_112r(:,2) + wear_121l(:,2) + wear_121r(:,2)) /4.*1000.*100.; % Multiply with 100 in order to % fake a longer simulation outer_wheel_prof_Y= outer_wheel_prof_Y - wear_outer; % update profile inner_wheel_prof_Y= inner_wheel_prof_Y - wear_inner; ident= sprintf('%s_%3.3d', vhe_model, iloop+1); write_2col_file(['kpf/w_prof/outer_wheel_', ident, '.wheel'], outer_wheel_prof_X, outer_wheel_prof_Y) write_2col_file(['kpf/w_prof/inner_wheel_', ident, '.wheel'], inner_wheel_prof_X, inner_wheel_prof_Y) disp('End wheel profile }}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}') iloop= iloop + 1; end % } =================================================================================== % figure % hold on % plot(outer_wheel_prof_X,outer_wheel_prof_Y) % hold off