home *** CD-ROM | disk | FTP | other *** search
/ Programmer 7500 / MAX_PROGRAMMERS.iso / PASCAL / NUMPAS.ZIP / NEWTON.PAS < prev    next >
Encoding:
Pascal/Delphi Source File  |  1980-01-01  |  2.2 KB  |  74 lines
  1. program newton(input,output);
  2. (* Program searches for the intersection of two nonlinear functions *)
  3. (* of y1 and y2 by the Newton-Raphson procedure. *)
  4.  
  5. var m,y1,y2,e: real;
  6.     l,max: integer;
  7.  
  8.  
  9. procedure check(x1,x2,d1,d2,e: real; var l: integer; max1: integer); forward;
  10. (* Checks maximum mismatch and continues iteration if maximum *)
  11. (* allowable mismatch is exceeded. *)
  12.  
  13. var m1,m2,f1,f2,j11,j12,j21,j22: real;
  14.  
  15. procedure estimate(var x1,x2: real; max1: integer);
  16. (* Uses Newton-Raphson technique to estimate new solution for x1 and x2 *)
  17. (* from previous estimate. *)
  18.  
  19. var f1,f2,j11,j12,j21,j22,d,d1,d2: real;
  20.    begin
  21.       f1 := x1 + 2.0*x1*x1*x1 + x2 - 2.0;            (* Function # 1 *)
  22.       f2 := x1 + x2 + 4.0*(1.0-exp(-3.0*x2)) - 4.0;  (* Function # 2 *)
  23.       j11 := 1.0 + 6.0*x1*x1;                        (* Jacobian elements *)
  24.       j12 := 1.0;
  25.       j21 := 1.0;
  26.       j22 := 1.0 + 12.0*exp(-3.0*x2);
  27.       d := j11*j22 - j12*j21;
  28.       d1 := (-j22*f1 + j12*f2)/d;
  29.       d2 := (j21*f1 - j11*f2)/d;
  30.       x1 := x1 + d1;            (* New estimate for x1 *)
  31.       x2 := x2 + d2;            (* New estimate for x2 *)
  32.       check(x1,x2,d1,d2,e,l,max1);
  33.    end;
  34.  
  35. procedure check;  (* Comment above *)
  36.  
  37.    begin
  38.       l := l + 1;
  39.       m := 0.0;
  40.       m1 := d1/x1;          (* Mismatch for x1 *)
  41.       m2 := d2/x2;          (* Mismatch for x2 *)
  42.       if abs(m1) > abs(m2) then
  43.          m := m1
  44.       else m := m2;         (* Find maximum mismatch *)
  45.       writeln(l:5,'':5,m:10,'':5,x1:10,'  ',x2:10);
  46.       if (abs(m) > e) and (l < max1) then
  47.          begin
  48.             estimate(x1,x2,max1);
  49.          end;
  50.    end;
  51.  
  52. begin  (* main *)
  53.    writeln('Functions and Jacobians should be in procedure "estimate".');
  54.    writeln;
  55.    write('Input initial estimates for x1 and x2 (with space between): ');
  56.    readln(y1,y2);
  57.    writeln;
  58.    write('Input maximum allowable mismatch: ');
  59.    readln(e);
  60.    writeln;
  61.    write('Input maximum number of iterations: ');
  62.    readln(max);
  63.    writeln;
  64.    writeln('Iter. #   ':10,'Max Mismatch    ':15,'x1  ':12,'x2':10);
  65.    writeln;
  66.    l := 0;
  67.    m := 999.999;
  68.    writeln(l:5,'':5,m:10,'':5,y1:10,'  ',y2:10); (* print initial values *)
  69.    estimate(y1,y2,max);
  70. end.
  71.  
  72.  
  73.  
  74.