Thursday, April 15, 2010

Solving a Nonlinear System of Four Equations from the Literature

Solving a Nonlinear System of Four Equations from the Literature

Jsun Yui Wong

"Solving systems of nonlinear equations is perhaps the most difficult problem in all of numerical computation," Rice [1, 1993, p. 355].

Originally intended for nonlinear integer programming, the following computer program tries to solve the following nonlinear system, which is based on a system in Shoup [2, pp. 40-41]. These equations are used in line 1141 through line 1150 below. As indicated in line 1157 through line 1160, tolerances of .0001s are used.

5*X(1)+3*X(2)+X(3)+X(4)=15
X(1)*X(2)+X(2)*X(3)+X(3)*X(4)=17
X(1)^2+X(2)^2+X(3)^2-X(4)^2=9
X(1)*X(3)+X(2)*X(4)+X(1)^3=8

0 DEFDBL A-Z
3 DEFINT I,J,K
4 DIM X(42),A(42),L(33),K(33)
5 FOR JJJJ=-32000 TO 32000
14 RANDOMIZE JJJJ
16 M=-1D+17
91 FOR KK=1 TO 4
94 A(KK)=-5+RND*10
95 NEXT KK
126 IMAR=10+FIX(RND*5000)
128 FOR I=1 TO IMAR
129 FOR K=1 TO 4
131 X(K)=A(K)
132 NEXT K
181 J=1+FIX(RND*4)
183 R=(1-RND*2)*A(J)
191 IF RND<.14 THEN X(J)=-5+RND*10 ELSE IF RND<.17 THEN X(J)=A(J)+R ELSE IF RND<.2 THEN X(J)=A(J)+RND*R ELSE IF RND<.25 THEN X(J)=A(J)+RND^2*R ELSE IF RND<.33 THEN X(J)=A(J)+RND^3*R ELSE IF RND<.5 THEN X(J)=A(J)+RND^4*R ELSE X(J)=INT(-5+RND*10)
1141 X(4)=-5*X(1)-3*X(2)-X(3)+15
1145 L1=X(1)*X(2)+X(2)*X(3)+X(3)*X(4)-17
1149 L2=X(1)^2+X(2)^2+X(3)^2-X(4)^2-9
1150 L3=X(1)*X(3)+X(2)*X(4)+X(1)^3-8
1152 L1K=L1
1153 L2K=L2
1154 L3K=L3
1157 IF ABS(L1)>.0001 THEN L1=-ABS(L1)*9999999! ELSE L1=0
1159 IF ABS(L2)>.0001 THEN L2=-ABS(L2)*9999999! ELSE L2=0
1160 IF ABS(L3)>.0001 THEN L3=-ABS(L3)*9999999! ELSE L3=0
1230 P1NEWMAY=L1+L2+L3
1448 P=P1NEWMAY
1451 IF P<=M THEN 1670
1657 FOR KEW=1 TO 4
1658 A(KEW)=X(KEW)
1659 NEXT KEW
1661 M=P
1663 LL1=L1K:LL2=L2K:LL3=L3K
1670 NEXT I
1890 IF M>-5 THEN 1912 ELSE 1999
1912 PRINT A(1),A(2),A(3),A(4),M,JJJJ,LL1,LL2,LL3
1999 NEXT JJJJ

The BASIC computer program above was run with Microsoft's GW BASIC 3.11 interpreter, which is not a compiler. Its complete output through JJJJ=-21946 is presented below. What immediately follows is a manual copy from the computer screen.

1 1 4 3 0
-30258 0 0 0

1 1 4 3 0
-28543 0 0 0

1 1 4 3 0
-21946 0 0 0

Interpreted in accordance with line 1912, the output above was obtained in thirty minutes on a personal computer with an Intel 2.66 GHz. chip and the IBM basica/D interpreter.

References

[1] J. Rice. Numerical methods, software, and analysis, second ed. Academic Press, 1993.

[2] T. E. Shoup. A Practical guide to computer methods for engineers. Prentice-Hall, 1979.

[3] W. H. Press, S. A. Teukolsky, W. T. Vetterling, B. P. Flannery. Numerical recipes: the art of scientific computing, third ed. Cambridge University Press, 2007.

[4] Microsoft Corp. BASIC, second edition (May 1982), Version 1.10. Boca Raton, Florida: IBM Corp., Personal Computer, P. O. Box 1328-C, Boca Raton, Florida 33432, 1981.