Переведи на язык r:
% ПРОГРАММА РАСЧЕТА ИЗОТЕРМЫ ЭКСТРАКЦИИ 
% ПО МОДЕЛИ ЭКСТРАКЦИОННОГО РАВНОВЕСИЯ А.М.РОЗЕНА
% 9:11 20.04.00
clc;clear;format short e; tic;    j=1;nq=51;ns=10;
%X(1) - U(VI)aq  X(2) - HNO3aq  X(3) - Pu(IV)aq  X(4) - Np(IV)aq
%X(5) - Th(IV)aq  X(6) - U(IV)aq  X(7) - Np(VI)aq  X(8) - Pu(VI)aq
%Y(1) - U(VI)org  Y(2) - HNO3org  Y(3) - Pu(IV)org  Y(4) - Np(IV)org
%Y(5) - Th(IV)org  Y(6) - U(IV)org  Y(7) - Np(VI)org  Y(8) - Pu(VI)org
%n1=55,n2=75:end if
Mu=238.04;  Mpu=239; Mnp=237; Mth=232; n=1; l=1;
X(1)=1e-4/Mu;X(2)=1.5;%X(3)=0.1/Mpu;X(4)=0.001/Mnp;X(5)=1e-9/Mth;X(6)=0.001/Mu;
%   ПОДПРОГРАММА ВЫЧИСЛЕНИЯ РАВНОВЕСНЫХ КОНЦЕНТРАЦИЙ КОМПОНЕНТОВ
while X(1)<1.8;
 S0=300*0.9727/266.32;
A=X(2)+2*X(1);%+4*X(3)+4*X(4)+4*X(5)+4*X(6);
IS=X(2)+3*X(1);%+10*X(3);
Kcu=exp(4.8-5.77*IS^0.5+4.17*IS-0.712*IS^1.5);
Kch1=0.385-0.155*IS+0.024*IS^2;   Kch2=0.0004;
%Kcr=exp(8.745-15.495*IS^0.5+7.117*IS-1.032*IS^1.5);           %? - question
%Kcpu=exp(9.44-13.02*IS^0.5+5.71*IS-0.76*IS^1.5);
%Kcnp=exp(6.33-11.5*IS^0.5+5.61*IS-0.833*IS^1.5);
%lnKcth1=2*log(Kcr)-log(Kcpu);   Kcth1=exp(lnKcth1);
%lnKcth2=3*log(Kcr)-2*log(Kcnp); Kcth2=exp(lnKcth2);
%lnKcth=(lnKcth1+lnKcth2)/2; Kcth=exp(lnKcth);
   fh1=Kch1*X(2)*A;
   fh2=Kch2*X(2)^2*A^2;
   fu=Kcu*A*A*X(1);
   %fpu=Kcpu*A^4*X(3);
  % fnp=Kcnp*A^4*X(4);
  % fth=Kcth*A^4*X(5);
fr=0;
S=2*S0/(1+fh1+fh2+sqrt((1+fh1+fh2)^2+8*S0*(fu+fr)));
   Y(1)=fu*S^2;       %  уран     - VI
   Y(2)=(fh1+fh2)*S;  %  кислота
   %Y(3)=fpu*S^2;      %  плутоний - IV
   %Y(4)=fnp*S^2;      %  нептуний - IV
   %Y(5)=fth*S^2 ;     % торий    - IV
   %Y(6)=fr*S^2;       %  уран     - IV
X(1)=X(1)*1.2; R0(l)=l;             
R1(l)=X(1);
R2(l)=X(2);
R3(l)=Y(1);
R4(l)=Y(2);
%R(5,l)=Y(3);R(6,l)=X(3);
D(1,l)=Kcu;D(2,l)=Kch1;D1(1,l)=log(Kcu);D1(2,l)=log(Kch1);
kt=1:3; %kt=1:7
rt(l,kt)=[l X(1) X(2)];      % X(3) X(4) X(5) X(6)];
%D(1,l)=Kcr;D(2,l)=Kcpu;D(3,l)=Kcnp;D(4,l)=Kcth1;D(5,l)=Kcth2;D(6,l)=Kcth;
%D1(1,l)=log(Kcr);D1(2,l)=log(Kcpu);D1(3,l)=Log(Kcnp);
%D1(4,l)=log(Kcth1);D1(5,l)=log(Kcth2);D1(6,l)=log(Kcth);
S(1,l)=IS;
S(2,l)=IS;
R1(l)=R1(l).*Mu;
R3(l)=R3(l).*Mu;
%R(4,l)=R(4,l).*Mnp;R(6,l)=R(6,l).*Mpu;
    l=l+1;
 end  %end while X(1)>1e-2
% ПОСТРОЕНИЕ РАБОЧЕЙ ЛИНИИ
X=[0,120];X1=[120,120];X2=[43.36,120];X3=[43.36,43.36];X4=[1.6,43.36];
Y=[0,440];Y1=[159,440];Y2=[159,159];Y3=[5.87,159];Y4=[5.87,5.87];
X5=[1.6,1.6];X6=[0.0293,1.6];X7=[0.0293,0.0293];X8=[0,0.0293];
Y5=[0.1073,5.87];Y6=[0.1073,0.1073];Y7=[1.92E-3,0.1073];Y8=[1.92E-3,1.92E-3];

figure(1)
plot(R1,R3,Y,X,Y1,X1,Y2,X2,Y3,X3,Y4,X4,Y5,X5,Y6,X6,Y7,X7,Y8,X8),
grid;title('C(U)в,равн');
xlabel('U(VI)водн');ylabel('U(VI)орг');
title('Изотерма экстракции U(VI)');
gtext('A');gtext('A1');gtext('B');gtext('B1');gtext('C');gtext('D');zoom on;
disp('Нажмите ENTER для продолжения работы');pause

%figure(2)
%subplot(2,2,1),plot(X(1),Y(1)),grid;zoom on;title('C(U)в,равн');%axis([2 nh 0 0.3]);
%subplot(2,2,2),plot(Y(1)),grid;zoom on;title('C(U)o,равн');%axis([2 nh 0 0.3]);
%subplot(2,2,3),plot(X(2)),grid;zoom on;title('C(H)в,равн');%axis([2 nh 0 0.3]);
%subplot(2,2,4),plot(Y(2)),grid;zoom on;title('C(H)о,равн');%axis([2 nh 0 0.3]);
%disp('Нажмите ENTER для продолжения работы');pause


   clc
%input"Если надо С(U) в г/л,введите 1,если в моль/л,нажмите ENTER ";q3:
%if q3=1 then:for i=0 to n2

  end;
q1=0;q2=0;
while q1~=4;
while q2<5;
q1=menu('<Меню выдачи результатов>-Вы можете выдать результаты',' На экран ',...
'Запись в файл','Возврат в начало программы',' Завершение работы ');
   if q1==4;q2=5;end
   if q1==1;

q2=menu('<МЕНЮ ВЫБОРА ТАБЛИЦ>-Вы можете выбрать','Концентрации','Константы',...
'Логарифмы констант','Все таблицы','Возврат в меню выбора устройства');
if (q2==1)+(q2==4);
clc; disp('     РАВНОВЕСНЫЕ КОНЦЕНТРАЦИИ ');
disp('   ')
disp(' -----------------------------------------------------------------------------');
disp('     N     U(VI)в         Hв       Pu(IV)в     U(VI)о        Hо        Pu(IV)о');
disp(' -----------------------------------------------------------------------------');
for l=1:ns:nq,
rc=sprintf('%5u  %5.3e  %5.3e ',l, rt(l,2),rt(l,3));disp(rc);
  end
disp('Нажмите ENTER для продолжения работы');pause
end       %end if (q2==1)+(q2==4);
if (q2==2)+(q2==4);pause
clc; disp('        КОНСТАНТЫ РАВНОВЕСИЯ An(IV)  ');
disp('   ')
disp('------------------------------------------------------------------------');
disp('   N     Kcr       Kcpu      Kcnp       Kcth1      Kcth2        Kcth ');
disp('  -----------------------------------------------------------------------');
for l=1:ns:nq,
cr=sprintf('%5u  %5.3e  %5.3e  %5.3e  %5.3e  %5.3e  %5.3e',l, ct(l,2),...
ct(l,3),ct(l,4),ct(l,5),ct(l,6),ct(l,7));disp(cr);
end
disp('Нажмите ENTER для продолжения работы');pause
end       %end if (q2==2)+(q2==4);
if (q2==3)+(q2==4);
clc; disp('  ЛОГАРИФМЫ КОНСТАНТ РАВНОВЕСИЯ');
disp('   ')
disp('------------------------------------------------------------------------');
disp('   N     LnKcu         LnKcnp   LnKcpu     LnKcth1    LnKcth2      LnKcth ');
disp(' -------------------------------------------------------------------------');
for l=1:ns:nq, 
lc=sprintf('%5u  %5.3e  %5.3e  %5.3e  %5.3e  %5.3e  %5.3e',l, lt(l,2),...
lt(l,3),lt(l,4),lt(l,5),lt(l,6),lt(l,7));disp(lc);
end
disp('Нажмите ENTER для продолжения работы');pause
end       %end if (q2==3)+(q2==4);


 %cls:print"   запись в файл "
  %input"Номер начального сечения: ";n1
%  input"Номер конечного сечения: ";n2: if n2=0 then
% n2=nh:end if
%input"Шаг таблицы: ";st:if st=0 then
%  st=fix(nh/5):end if:  cls
%rem пример записи информации в файл
%OPEN"CCC.dta" for output as #1
                    pause    
q0=menu('  <Меню способа выдачи результатов> - Вы можете выбрать',...
'На экран','На печать','В виде графиков','В электронные таблицы',...
'Запись в файл','Возврат в начало программы','Завершение работы');
     if q0==7, q1=4; q2=8; q3=4; end;    
   if q0==1;
q1=menu('  <Меню выбора группы таблиц или переменных> - Вы можете выбрать',...
'макрокомпоненты','Pu,Np,Th',...
'D,K и вспомогательные переменные',...
'Возврат в меню способа выдачи результатов');

      if q1==1;
q2=menu('  <Меню выбора таблиц или переменных> - Вы можете выбрать',...
'Равновесные значения','Рабочие значения','Шестивалентный уран',...
'Азотная кислота','Вспомогательные значения','Коэффициенты распределения',...
'Все таблицы (переменные)','Возврат в меню способа выдачи результатов');
               % if (q2==1)+(q2==13);
  if (q2==1)+(q2==7);
clc;disp('            РАВНОВЕСНЫЕ КОНЦЕНТРАЦИИ ');
disp('  N    U(VI)в     H(водн.)    U(VI)о     H(орг.)    U(IV)в    U(IV)о')           
disp('---------------------------------------------------------------------');
     for l=1:ns:nh;
sce=sprintf('%3u  %5.3e %5.3e %5.3e %5.3e',l,cet(l,2),...  
cet(l,3),cet(l,4),cet(l,5));
disp (sce);
     end;
disp('Нажмите ENTER для продолжения работы');pause
         end ; %end  if (q2==1)+(q2==4);
         %if (q2==2) %~= (q2==13) 
      end     %end  if q1==1;
   end   %end if q0==1;
end;end
end