API tools faq
paste
Advertisement
Jade39

D-N

Jade39
May 1st, 2014
254
0
Never
Add comment
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
Delphi 19.80 KB | None | 0 0
raw download clone embed print report
  1. program NNNNN1;
  2.  
  3. {$APPTYPE CONSOLE}
  4. {$MAXSTACKSIZE 200000000}
  5. {$Q+,I+,R+}
  6.  
  7. uses
  8.   SysUtils;
  9.  
  10. CONST
  11.     eps=0.0001;
  12.     pi=3.1415926535897932384626433832795;
  13.  
  14. TYPE
  15.     vect=Array[1..300] of Real;
  16.     t_vect=Array[1..300] of vect;
  17.  
  18. VAR iter,iter1,i,j,n,m,i1,j1,n1,m1,k,kk,k1,sh1,sh2,ll : Integer;
  19.     norma,norma1,tau,ro,ro_, first,opt,opt1,ss,tt,ro1,ro_1 : Real;
  20.     r,f,y,p,w,tmp,pr1,pr2,pr3,y1,y2,y3 : T_VECT;
  21.     rr,zz,y12,y11,a33p,am,x,yy,pp,yy1,pq,xx,y4,q:vect;
  22. {-------------------------------------------------------------------------------}
  23. FUNCTION scalar(v1,v2: T_VECT;n,m:integer): Real;
  24. Var res : Real;
  25.     i,j : Integer;
  26. Begin
  27.  res:=0;
  28.  For i:=1 To n Do
  29.      For j:=1 To m Do
  30.          res:=res + v1[i][j]*v2[i][j];
  31.  scalar:=res;
  32. End;
  33. {-------------------------------------------------------------------------------}
  34. procedure print(a:t_vect;n,m:integer);
  35. var i,j:integer;
  36. begin
  37.     for i:=1 to n do begin
  38.     for j:=1 to m do begin
  39.       write(a[i,j]:6:4,' ');end;writeln;end;
  40. end;
  41. {-------------------------------------------------------------------------------}
  42. FUNCTION scalar1(v1,v2: vect;n:integer): Real;
  43. Var res : Real;
  44.     i,j : Integer;
  45. Begin
  46.  res:=0;
  47.  For i:=1 To n Do
  48.          res:=res + v1[i]*v2[i];
  49.  scalar1:=res;
  50. End;
  51. {-------------------------------------------------------------------------------}
  52.  function minus(a,b:vect;m:integer):vect;
  53. var i,j,k:integer;
  54. tmp:real;
  55. c:vect;
  56. begin
  57. FillChar(c,sizeof(c),0);
  58. for i:=1 to m do
  59.      c[i]:=a[i]-b[i];
  60.   minus:=c;
  61. end;
  62. {-------------------------------------------------------------------------------}
  63.  function plus(a,b:vect;m:integer):vect;
  64. var i,j,k:integer;
  65. tmp:real;
  66. c:vect;
  67. begin
  68. for i:=1 to m do
  69.      c[i]:=a[i]+b[i];
  70.   plus:=c;
  71. end;
  72. {-------------------------------------------------------------------------------}
  73. function pulti(a:T_VECT;b:vect;m:integer):vect;
  74. var i,j,k:integer;
  75. tmp:real;
  76. c:vect;
  77. begin
  78.  
  79.   for i:=1 to m do begin
  80.      tmp:=0;
  81.      for k:=1 to m do begin
  82.          tmp:=tmp+a[i][k]*b[k];
  83.      end;
  84.      c[i]:=tmp;
  85.    end;
  86.    pulti:=c;
  87. end;
  88. {-------------------------------------------------------------------------------}
  89. function pulti1(a:real;b:vect;m:integer):vect;
  90. var i,j,k:integer;
  91. c:vect;
  92. begin
  93. FillChar(c,sizeof(c),0);
  94.   for i:=1 to m do begin
  95.          c[i]:=a*b[i];
  96.      end;
  97.    pulti1:=c;
  98. end;
  99. {-------------------------------------------------------------------------------}
  100. FUNCTION raz_raz1 (r: T_VECT;n,m:integer): T_VECT;
  101. Var res,z : T_VECT;
  102.     i,j : Integer;
  103. Begin
  104.  FillChar(res,sizeof(res),0);
  105.  
  106.     For i:=1 To n  Do Begin
  107.                RES[i,1]:=2*R[i,1];
  108.                If i<>1 Then
  109.                     RES[i,1]:=RES[i,1]- 0.5*R[i-1,1];
  110.          If i<>n Then
  111.                     RES[i,1]:=RES[i][1] -  0.5*R[i+1,1];
  112.  
  113.                     RES[i,1]:=RES[i][1] - R[i,2];
  114.          End;
  115.         // print(res,n,m);
  116.           //z:=res;
  117. For i:=1 To n  Do Begin
  118.          For j:=2 To m Do Begin
  119.  
  120.                RES[i,j]:=4*r[i,j];
  121.                If i<>1 Then
  122.                     RES[i,j]:=RES[i,j]- r[i-1,j];
  123.          If i<>n Then
  124.                     RES[i,j]:=RES[i][j] -  r[i+1,j];
  125.                If j<>1 Then
  126.                     RES[i,j]:=RES[i][j] - r[i,j-1];
  127.                If j<>m Then
  128.                     RES[i,j]:=RES[i][j] - r[i,j+1];
  129.      End;
  130.  End;
  131.  
  132.  raz_raz1:=RES;
  133. End;
  134. {-------------------------------------------------------------------------------}
  135. FUNCTION GET_TAU1(r: real; p: T_VECT; n,m:integer): Real;
  136. Var s2 : Real;
  137.     tmp : T_VECT;
  138.     i,j : Integer;
  139. Begin
  140.  s2:=0;
  141.  FillChar(tmp,sizeof(tmp),0);
  142.  tmp:=raz_raz1(p,n,m);
  143.  For i:=1 To n  Do
  144.          For j:=1 To m Do
  145.                s2:=s2 + TMP[i,j]*P[i,j];
  146.  GET_TAU1:=r/s2;
  147. End;
  148. {-------------------------------------------------------------------------------}
  149. FUNCTION SSOR_PRE1(f: T_VECT;n,m:integer): T_VECT;
  150. Var y,tmp,r,y2 : T_VECT;
  151.     iter,i,j : Integer;
  152.     norma,prev,s,w : Real;
  153. Begin
  154.  iter:=0;
  155.  norma:=0;
  156.  prev:=0;
  157.  s:=0;
  158.  
  159.  w:=2/(1+sin(pi/(n)));
  160.  //writeln('w= ',w:4:4);
  161.  
  162.  FillChar(y,sizeof(y),0);
  163.  FillChar(r,sizeof(r),0);
  164.  FillChar(tmp,sizeof(tmp),0);
  165.  tmp:=raz_raz1(y,n,m);
  166.  For i:=1 To N Do
  167.          For j:=1 To M Do Begin
  168.                R[i,j]:=F[i,j] - TMP[i,j];
  169.          If norma<abs(R[i,j]) Then
  170.             norma:=abs(R[i,j]);
  171.          End;
  172.  
  173.  For i:=1 To N Do begin
  174.          s:=0;
  175.              Y[i,1]:=(1-w)*Y[i,1];
  176.              If i<>1 Then
  177.                   s:=s + 0.5*Y[i-1,1];
  178.          If i<>n Then
  179.                   s:=s + 0.5*Y[i+1,1];
  180.                   s:=s + Y[i,2];
  181.              s:=s * 0.5*w;
  182.              Y[i,1]:=Y[i][1] + s + 0.5*w*F[i,1];
  183.          //u(i,j) = (1-w)*u(i,j)+0.25*w*(u(i-1,j)+u(i,j-1)+u(i+1,j)+u(i,j+1))+0.25*w*f(i,j);
  184.      End;
  185.  
  186.      For i:=1 To N Do begin
  187.      For j:=2 To M Do Begin
  188.          s:=0;
  189.          //if j<>1 then
  190.              Y[i,j]:=(1-w)*Y[i,j];
  191.              If i<>1 Then
  192.                   s:=s + Y[i-1,j];
  193.          If i<>n Then
  194.                   s:=s + Y[i+1,j];
  195.          If j<>1 Then
  196.                   s:=s + Y[i,j-1];
  197.          If j<>m Then
  198.                   s:=s + Y[i,j+1];
  199.              s:=s * 0.25*w;
  200.         // if j<>1 then
  201.              Y[i,j]:=Y[i][j] + s + 0.25*w*F[i,j];
  202.          //u(i,j) = (1-w)*u(i,j)+0.25*w*(u(i-1,j)+u(i,j-1)+u(i+1,j)+u(i,j+1))+0.25*w*f(i,j);
  203.      End;
  204.      end;
  205.  
  206.  //obratnii hod
  207.  For i:=n DownTo 1 Do
  208.          For j:=n DownTo 2 Do Begin
  209.                s:=0;
  210.          //if j<>1 then
  211.                Y[i,j]:=(1-w)*Y[i,j];
  212.                If i<>1 Then
  213.                     s:=s + Y[i-1,j];
  214.                If i<>n Then
  215.                     s:=s + Y[i+1,j];
  216.                If j<>1 Then
  217.                     s:=s + Y[i,j-1];
  218.                If j<>m Then
  219.                     s:=s + Y[i,j+1];
  220.          s:=s * 0.25*w;
  221.          //if j<>1 then
  222.          Y[i,j]:=Y[i][j] + s + 0.25*w*F[i,j];
  223.                //u(i,j) = (1-w)*u(i,j)+0.25*w*(u(i-1,j)+u(i,j-1)+u(i+1,j)+u(i,j+1))+0.25*w*f(i,j);
  224.          End;
  225.  
  226.         For i:=n DownTo 1 Do begin
  227.                s:=0;
  228.                Y[i,1]:=(1-w)*Y[i,1];
  229.                If i<>1 Then
  230.                     s:=s + 0.5*Y[i-1,1];
  231.                If i<>n Then
  232.                     s:=s + 0.5*Y[i+1,1];
  233.                     s:=s + Y[i,2];
  234.          s:=s * 0.5*w;
  235.          Y[i,1]:=Y[i][1] + s + 0.5*w*F[i,1];
  236.                //u(i,j) = (1-w)*u(i,j)+0.25*w*(u(i-1,j)+u(i,j-1)+u(i+1,j)+u(i,j+1))+0.25*w*f(i,j);
  237.          End;
  238.  
  239.  
  240.     //norma neviazki
  241.  For i:=1 To n Do
  242.         For j:=1 To m Do Begin
  243.             If i<>1 Then
  244.                  norma:=norma + Y[i-1,j];
  245.             If i<>n Then
  246.                  norma:=norma + Y[i+1,j];
  247.             If j<>1 Then
  248.                  norma:=norma + Y[i,j-1];
  249.             If j<>m Then
  250.                  norma:=norma + Y[i,j+1];
  251.             norma:=norma - 4*Y[i,j];
  252.             norma:=abs(F[i,j] + norma);
  253.         End;
  254.  
  255.     SSOR_PRE1:=Y;
  256. End;
  257. {-------------------------------------------------------------------------------}
  258. FUNCTION raz_raz (r: T_VECT;n,m:integer): T_VECT;
  259. Var res : T_VECT;
  260.     i,j : Integer;
  261. Begin
  262.  //FillChar(res,sizeof(res),0);
  263.  For i:=1 To n  Do Begin
  264.          For j:=1 To m Do Begin
  265.                RES[i,j]:=4*R[i,j];
  266.                If i<>1 Then
  267.                     RES[i,j]:=RES[i,j] - R[i-1,j];
  268.          If i<>n Then
  269.                     RES[i,j]:=RES[i][j] -  R[i+1,j];
  270.                If j<>1 Then
  271.                     RES[i,j]:=RES[i][j] - R[i,j-1];
  272.                If j<>m Then
  273.                     RES[i,j]:=RES[i][j] - R[i,j+1];
  274.      End;
  275.  End;
  276.  
  277.  raz_raz:=RES;
  278. End;
  279. {-------------------------------------------------------------------------------}
  280. FUNCTION GET_TAU(r: real; p: T_VECT; n,m:integer): Real;
  281. Var s2 : Real;
  282.     tmp : T_VECT;
  283.     i,j : Integer;
  284. Begin
  285.  s2:=0;
  286.  FillChar(tmp,sizeof(tmp),0);
  287.  tmp:=raz_raz(p,n,m);
  288.  For i:=1 To n  Do
  289.          For j:=1 To m Do
  290.                s2:=s2 + TMP[i,j]*P[i,j];
  291.  GET_TAU:=r/s2;
  292. End;
  293. {-------------------------------------------------------------------------------}
  294. FUNCTION SSOR_PRE(f: T_VECT;n,m:integer): T_VECT;
  295. Var y,tmp,r : T_VECT;
  296.     iter,i,j : Integer;
  297.     norma,prev,s,w : Real;
  298. Begin
  299.  iter:=0;
  300.  norma:=0;
  301.  prev:=0;
  302.  s:=0;
  303.  w:=2/(1+sin(pi/(M+1)));
  304.  //writeln('w= ',w:4:4);
  305.  
  306.  FillChar(y,sizeof(y),0);
  307.  FillChar(r,sizeof(r),0);
  308.  FillChar(tmp,sizeof(tmp),0);
  309.  tmp:=raz_raz(y,n,m);
  310.  For i:=1 To N Do
  311.          For j:=1 To M Do Begin
  312.                R[i,j]:=F[i,j] - TMP[i,j];
  313.          If norma<abs(R[i,j]) Then
  314.             norma:=abs(R[i,j]);
  315.          End;
  316.  For i:=1 To N Do
  317.      For j:=1 To M Do Begin
  318.          s:=0;
  319.              Y[i,j]:=(1-w)*Y[i,j];
  320.              If i<>1 Then
  321.                   s:=s + Y[i-1,j];
  322.          If i<>n Then
  323.                   s:=s + Y[i+1,j];
  324.          If j<>1 Then
  325.                   s:=s + Y[i,j-1];
  326.          If j<>m Then
  327.                   s:=s + Y[i,j+1];
  328.              s:=s * 0.25*w;
  329.              Y[i,j]:=Y[i][j] + s + 0.25*w*F[i,j];
  330.          //u(i,j) = (1-w)*u(i,j)+0.25*w*(u(i-1,j)+u(i,j-1)+u(i+1,j)+u(i,j+1))+0.25*w*f(i,j);
  331.      End;
  332.  //obratnii hod
  333.  For i:=n DownTo 1 Do
  334.          For j:=n DownTo 1 Do Begin
  335.                s:=0;
  336.                Y[i,j]:=(1-w)*Y[i,j];
  337.                If i<>1 Then
  338.                     s:=s + Y[i-1,j];
  339.                If i<>n Then
  340.                     s:=s + Y[i+1,j];
  341.                If j<>1 Then
  342.                     s:=s + Y[i,j-1];
  343.                If j<>m Then
  344.                     s:=s + Y[i,j+1];
  345.          s:=s * 0.25*w;
  346.          Y[i,j]:=Y[i][j] + s + 0.25*w*F[i,j];
  347.                //u(i,j) = (1-w)*u(i,j)+0.25*w*(u(i-1,j)+u(i,j-1)+u(i+1,j)+u(i,j+1))+0.25*w*f(i,j);
  348.          End;
  349.  
  350.     //norma neviazki
  351.  For i:=1 To n Do
  352.         For j:=1 To m Do Begin
  353.             If i<>1 Then
  354.                  norma:=norma + Y[i-1,j];
  355.             If i<>n Then
  356.                  norma:=norma + Y[i+1,j];
  357.             If j<>1 Then
  358.                  norma:=norma + Y[i,j-1];
  359.             If j<>m Then
  360.                  norma:=norma + Y[i,j+1];
  361.             norma:=norma - 4*Y[i,j];
  362.             norma:=abs(F[i,j] + norma);
  363.         End;
  364.  
  365.     SSOR_PRE:=Y;
  366. End;
  367. {-------------------------------------------------------------------------------}
  368. function A33(a,b:integer):T_VECT;
  369. var
  370. i,j:integer;
  371. q:T_VECT;
  372. begin
  373. for i:=1 to a do  begin
  374. for j:=1 to b do  begin
  375.        if i=j then q[i,j]:=4
  376.        else if (i+1=j) or (i=j+1) then q[i,j]:=-1
  377.        else   q[i,j]:=0;
  378.    end;end;
  379.   {  for i:=1 to a do  begin
  380.      for j:=1 to b do  begin
  381.      write(q[i,j]:6:1,' ');
  382.      end; writeln; end; }
  383.  
  384.      A33:=q;
  385. end;
  386. {-------------------------------------------------------------------------------}
  387. procedure Grad0(n,m:integer;r:T_VECT);
  388.  
  389. BEGIN
  390.  
  391.  //FillChar(r,sizeof(r),16);
  392. // FillChar(F,sizeof(F),16);
  393.  //FillChar(y,sizeof(y),0);
  394.  //FillChar(p,sizeof(p),0);
  395.  //FillChar(w,sizeof(w),0);
  396. // FillChar(tmp,sizeof(tmp),0);
  397.  iter:=0;
  398.  //norma:=eps+1;
  399.  tau:=0;
  400.  ro:=0;
  401.  ro_:=0;
  402.  For i:=1 To n Do begin
  403.      For j:=1 To m Do Begin
  404.          //F[i][j]:=16;
  405.          //R[i][j]:=16;
  406.          Y[i][j]:=0;
  407.          P[i][j]:=0;
  408.          W[i][j]:=0;
  409.          tmp[i][j]:=0;
  410.      End;end;
  411.      //f:=r;
  412.  p:=ssor_pre1(r,n,m);
  413.  ro:=scalar(p,r,n,m);
  414.  
  415.  While true Do Begin
  416.      INC(iter);
  417.          norma:=0;
  418.          tau:=get_tau1(ro,p,n,m);
  419.     // writeln('--------------tau---------------');
  420.     // writeln('tau= ',tau);
  421.          For i:=1 To n Do
  422.                For j:=1 To m Do
  423.                      Y[i,j]:=Y[i][j] + tau*P[i,j];
  424.             // writeln('--------------x---------------');
  425.              //print(y,n,m);
  426.      tmp:=raz_raz1(p,n,m);
  427.          For i:=1 To n Do
  428.          For j:=1 To m Do Begin
  429.                      r[i,j]:=r[i][j] - tau*TMP[i,j];
  430.                norma:=norma+abs(r[i,j]);
  431.                     // If norma<abs(r[i,j]) Then
  432.                         //  norma:= abs(r[i,j]);
  433.                End;
  434.          //writeln('--------------r---------------');
  435.          //print(r,n,m);
  436.  
  437.      if iter=1 then first:=norma;
  438.     //writeln('norma_next=',norma:6:3,' norma1=',first:6:3, ' ');
  439.      if abs(norma/first)<eps then break;
  440.  
  441.      w:=ssor_pre1(r,n,m); //pereschet po formyle
  442.      //writeln('w= ',w[1,1]);
  443.  
  444.          ro_:=ro;
  445.          ro:=scalar(w,r,n,m);
  446.  
  447.          For i:=1 To n Do
  448.          For j:=1 To m Do
  449.              P[i,j]:=W[i,j] + (ro/ro_)*P[i,j];
  450.             // writeln('--------------p---------------');
  451.            // print(p,n,m);
  452.  
  453.  End;
  454.  
  455. for i:=1 to n do begin
  456.   for j:=1 to m do
  457.     if j=1 then
  458.       y4[i]:=y[i,j];
  459. end;
  460.  
  461.  
  462. end;
  463. {-------------------------------------------------------------------------------}
  464. procedure Grad(n,m:integer;r:T_VECT);
  465.  
  466. BEGIN
  467.  
  468.  //FillChar(r,sizeof(r),16);
  469. // FillChar(F,sizeof(F),16);
  470.  //FillChar(y,sizeof(y),0);
  471.  //FillChar(p,sizeof(p),0);
  472.  //FillChar(w,sizeof(w),0);
  473. // FillChar(tmp,sizeof(tmp),0);
  474.  iter:=0;
  475.  norma:=eps+1;
  476.  tau:=0;
  477.  ro:=0;
  478.  ro_:=0;
  479.  For i:=1 To n Do begin
  480.      For j:=1 To m Do Begin
  481.          //F[i][j]:=16;
  482.          //R[i][j]:=16;
  483.          Y[i][j]:=0;
  484.          P[i][j]:=0;
  485.          W[i][j]:=0;
  486.          tmp[i][j]:=0;
  487.      End;end;
  488.      //f:=r;
  489.  p:=ssor_pre(r,n,m);
  490.  
  491.  ro:=scalar(p,r,n,m);
  492.  
  493.  While true Do Begin
  494.      INC(iter);
  495.          norma:=0;
  496.          tau:=get_tau(ro,p,n,m);
  497.      //writeln('tau= ',tau);
  498.          For i:=1 To n Do
  499.                For j:=1 To m Do
  500.                      Y[i,j]:=Y[i][j] + tau*P[i,j];
  501.      tmp:=raz_raz(p,n,m);
  502.          For i:=1 To n Do
  503.          For j:=1 To m Do Begin
  504.                      r[i,j]:=r[i][j] - tau*TMP[i,j];
  505.                norma:=norma+abs(r[i,j]);
  506.                     // If norma<abs(r[i,j]) Then
  507.                         //  norma:= abs(r[i,j]);
  508.                End;
  509.  
  510.      if iter=1 then first:=norma;
  511.      if norma/first<eps then break;
  512.      w:=ssor_pre(r,n,m); //pereschet po formyle
  513.      //writeln('w= ',w[1,1]);
  514.          ro_:=ro;
  515.          ro:=scalar(w,r,n,m);
  516.  
  517.          For i:=1 To n Do
  518.          For j:=1 To m Do
  519.              P[i,j]:=W[i,j] + (ro/ro_)*P[i,j];
  520.  End;
  521. { WriteLn('Gradient method with predobusl: ',iter);
  522.  For i:=0 To N-1 Do Begin
  523.      For j:=0 To M-1 Do
  524.          write(Y[i][j]:6:2,'  ');
  525.      WriteLn;
  526.  End;}
  527. end;
  528. {-------------------------------------------------------------------------------}
  529. //òóò íà ëèíèè ñöåïëåíèÿ
  530. procedure nev(n,k,bb,m:integer;z:T_VECT);
  531. var
  532. i,l,i1:integer;
  533. sh:T_VECT;
  534. begin
  535.    Grad(n,k,z);
  536.    sh:=y;
  537.    if(bb=1) then begin
  538.           l:=1;
  539.           for i:=m to n do begin
  540.           zz[l]:=sh[i,n];l:=l+1;//k:=k+1;
  541.           end; end   //end;//end;end;
  542.    else if (bb=2) then begin
  543.           l:=1;
  544.           for i1:=1 to n do begin
  545.             zz[l]:=sh[i1,1];l:=l+1;
  546.           end;end
  547.    else write('ne pravilno vveli bb');
  548. end;
  549. {-------------------------------------------------------------------------------}
  550. //òóò íàõîäèì am
  551. procedure naxod_am(k,n,n1,m,m1:integer;rr,pq:vect);
  552. var
  553. i,j,i1,j1,ll:integer;
  554. begin
  555. //---------------------------------------------------------------
  556. // writeln('-------------1------------');writeln;
  557.               ll:=k;
  558.              for i:=1 to n do begin
  559.  
  560.              for j:=1 to m do begin
  561.               if(i=ll)and(j=m) then begin
  562.                pr1[i,j]:=pq[ll-k+1];
  563.                ll:=ll+1;
  564.               end
  565.              else pr1[i,j]:=0;
  566.              end;end;
  567.  
  568. nev(n,m,1,k,pr1);
  569. y11:=zz;
  570.   //for i1:=1 to n1 do
  571.   //writeln(y11[i1]:2:4,' ');
  572. //writeln;
  573.  //writeln('-------------2------------');writeln;
  574.  
  575.  sh2:=1;
  576.       //for ll:=1 to n1 do begin
  577.         //sh2:=sh2+1;
  578.           for i:=1 to n1 do begin
  579.           for j:=1 to m1 do begin
  580.             if (j=1) then begin
  581.               pr2[i,j]:=pq[sh2];
  582.               sh2:=sh2+1;
  583.             end
  584.             else pr2[i,j]:=0;end;end;
  585. nev(n1,m1,2,k,pr2);
  586. y12:=zz;
  587.   //for i1:=1 to n1 do
  588.   //writeln(y12[i1]:2:4,' ');
  589. //Writeln('------------------A33*p--------------------------');
  590. a33p:=pulti(A33(n1,m1),pq,n1);
  591. //for i:=1 to n1 do
  592.   //Writeln(a33p[i]:6:2,' ');
  593. //Writeln('------------------A33*p-2*p-3*p--------------------------');
  594. am:=minus(a33p,plus(y11,y12,n1),n1);
  595.   //for i:=1 to n1 do
  596.   //Writeln(am[i]:6:2,' '); Writeln;//am=sp
  597. //-----------------------------------------------------------------------
  598.  
  599. end;
  600. {-------------------------------------------------------------------------------}
  601. procedure Grad2(k,n,n1,m,m1:integer;rr:vect);
  602. var
  603. iter:integer;
  604. x:vect;
  605. begin
  606. iter:=0;
  607. norma1:=0;
  608.  FillChar(x,sizeof(x),0);
  609.  kk:=1;
  610.    // ïðàâà ÷àñòü ìàëîãî
  611. for j1:=1 to m1+1 do begin
  612. for i1:=1 to n1 do begin
  613.   if (j1=1) then begin
  614.    pr3[i1,j1]:=rr[kk];
  615.    kk:=kk+1;
  616.   end
  617.   else pr3[i1,j1]:=0;
  618. end;
  619. end;
  620.  
  621. Grad0(n1,m1+1,pr3);
  622.  pq:=y4;
  623.  ro1:=scalar1(pq,rr,n1);
  624.  writeln('-----------------');
  625.  for i:=1 to n1 do
  626.  writeln(pq[i]:6:4,' ');
  627.  writeln('-----------------');
  628.   writeln('ro= ',ro:6:6,' ');
  629.  
  630.  
  631. //while iter>0 do begin
  632. While true Do Begin
  633.  
  634. //iter:=iter-1;
  635. inc(iter);
  636. writeln('---------------------------------------------------- ');
  637. writeln('íîìåð èòåðàöèè:',iter);
  638.  
  639. naxod_am(k,n,n1,m,m1,rr,pq);
  640.  
  641.   Writeln('------------------a_m--------------------------');
  642. for i:=1 to n1 do
  643. writeln(am[i]:6:5,' ');
  644. writeln('---------------------------------------------------- ');
  645.  
  646. tt:=scalar1(am,pq,n1);
  647. writeln('ro=',ro1:6:10,' ');
  648. writeln('scalar1(am,pq,n1)=',tt:6:2,' ');
  649. opt:=ro1/tt;
  650.  
  651.   writeln('lamda opt=',opt:6:10,' ');
  652.   Writeln('------------------x_m+1--------------------------');
  653.  
  654.  yy1:=pulti1(opt,pq,n1);
  655.  
  656.  x:=plus(x,yy1,n1);
  657.  for i:=1 to n1 do
  658.  writeln(x[i]:6:10,' ');writeln;
  659.  
  660.  Writeln('------------------r_m+1--------------------------');
  661.  
  662.  rr:=minus(rr,pulti1(opt,am,n1),n1);
  663.  
  664.    if iter=1 then begin
  665.     For i:=1 To N  Do
  666.        norma1:=norma1+abs(rr[i]);
  667.   end;
  668.  
  669.  for i:=1 to n1 do
  670.  writeln(rr[i]:6:6,' ');writeln;
  671.  
  672.        norma:=0;
  673.      For i:=1 To N  Do
  674.                      norma:=norma+abs(rr[i]);
  675.               ss:=abs(norma/norma1);
  676.           writeln('----------óñëîâèå âûõîà:abs(norma/norma1)<eps--------------');
  677.           writeln(ss:6:10,' < ',eps:6:6);
  678.          If ss<eps Then
  679.                 break;
  680.  
  681.   kk:=1;
  682.    // ïðàâà ÷àñòü ìàëîãî
  683. for j1:=1 to m1+1 do begin
  684. for i1:=1 to n1 do begin
  685.   if (j1=1) then begin
  686.    pr3[i1,j1]:=rr[kk];
  687.    kk:=kk+1;
  688.   end
  689.   else pr3[i1,j1]:=0;
  690. end;
  691. end;
  692.  
  693. Grad0(n1,m1+1,pr3);
  694.  q:=y4;
  695.  ro_1:=ro1;
  696.  ro1:=scalar1(q,rr,n1);
  697.  opt1:=ro1/ro_1;
  698.  pq:=plus(q,pulti1(opt1,pq,n1),n1);
  699.  Writeln('------------------p_m+1--------------------------');
  700.  for i:=1 to n1 do
  701.   writeln(pq[i]:6:6,' ');
  702. //writeln('---------------------------------------------------- ');
  703.  
  704.  
  705.  
  706. end;
  707. xx:=x;
  708. iter1:=iter;
  709. writeln('opa');
  710. end;
  711.  
  712. {-------------------------------------------------------------------------------}
  713. procedure Grad3(n,m,n1,m1:integer;pr1,pr2:T_vect);
  714. begin
  715. Grad(n,m,pr1);
  716. y1:=y;
  717. WriteLn('----------------Ðåøåíèå 1 îáëàñòè------------- ');WriteLn;
  718.  For i:=1 To n Do Begin
  719.      For j:=1 To m Do
  720.          write(y1[i][j]:6:2,'  ');
  721.      WriteLn;
  722.  End;
  723.   WriteLn;
  724. Grad(n1,m1,pr2);
  725. y2:=y;
  726. WriteLn('----------------Ðåøåíèå 2 îáëàñòè------------- ');WriteLn;
  727.  For i:=1 To n1 Do Begin
  728.      For j:=1 To m1 Do
  729.          write(y2[i][j]:6:2,'  ');
  730.      WriteLn;
  731.  End;
  732. end;
  733. {-------------------------------------------------------------------------------}
  734.  
  735. begin
  736. assign(output,'output.txt');
  737. rewrite(output);
  738. n:=7; m:=7; // ðàçìåð áîëüøîãî êâàäðàòà
  739. n1:=3; m1:=3; // ðàçìåð ìåíüøåãî
  740. k:=n-n1+1; k1:=k+n1-1;  // ê - ñ êàêîé ñòðîêè ...
  741. //k:=5;
  742.  
  743. for i:=1 to n do
  744. for j:=1 to m do
  745. pr1[i,j]:=16;   // ïðàâàÿ ÷àñòü áîëüøîãî
  746.  
  747. for i1:=1 to n1 do
  748. for j1:=1 to m1 do
  749. pr2[i1,j1]:=16;  // ïðàâà ÷àñòü ìàëîãî
  750.  
  751. Grad3(n,m,n1,m1,pr1,pr2);
  752.  
  753. writeln('-------------íåâÿçêà íà ëèíèè ñöåïëåíèÿ------------');writeln;
  754.       for i:=1 to n1 do  begin
  755.         rr[i]:=y2[i][1]+y1[k+i-1][n]+16;
  756.         writeln(rr[i]:6:6,' ');
  757.        end;
  758.  
  759. Grad2(k,n,n1,m,m1,rr);
  760.  
  761.  writeln('---------Íàø âåêòîð çà ',iter1,' èòåðàöèé--------------');
  762.  for i:=1 to n1 do
  763.     writeln(xx[i]:6:4,' '); writeln;
  764.  
  765.               ll:=k;
  766.              for i:=1 to n do begin
  767.  
  768.              for j:=1 to n do begin
  769.               if(i=ll)and(j=n) then begin
  770.                pr1[i,j]:=xx[ll-k+1]+16;
  771.                ll:=ll+1;
  772.               end
  773.              else pr1[i,j]:=16;
  774.              end;end;
  775.  sh2:=1;
  776.       //for ll:=1 to n1 do begin
  777.         //sh2:=sh2+1;
  778.           for i:=1 to n1 do begin
  779.           for j:=1 to n1 do begin
  780.             if (j=1) then begin
  781.               pr2[i,j]:=xx[sh2]+16;
  782.               sh2:=sh2+1;
  783.             end
  784.             else pr2[i,j]:=16;end;end;
  785.  
  786. Grad3(n,m,n1,m1,pr1,pr2); writeln;
  787.  
  788. writeln('-------------Ðåøåíèå çàäà÷è------------'); writeln;
  789.     for i:=1 to n do begin
  790.     for j:=1 to n+n1+1 do begin
  791.       if (j<=n) then
  792.       write(y1[i][j]:6:2,' ')
  793.       else if (i>=k)and(i<=k+n1-1) and(j=n+1) then
  794.       write(xx[i-k+1]:6:2,' ')
  795.       else if (i>=k)and(i<=k+n1-1)and (j>n+1)and(j<=n+n1+1) then
  796.       write(y2[i-k+1][j-n-1]:3:2,' ')
  797.       else write(' ')
  798.       end;writeln;end;
  799.  
  800.  
  801.  
  802. END.
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement
Public Pastes
Advertisement
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.  OK, I Understand
Not a member of Pastebin yet?
Sign Up, it unlocks many cool features!