RE: Sovling non-linear eq-sys (beginners question)

*To*: mathgroup at smc.vnet.net*Subject*: [mg71233] RE: [mg71207] Sovling non-linear eq-sys (beginners question)*From*: "David Park" <djmp at earthlink.net>*Date*: Sat, 11 Nov 2006 03:38:19 -0500 (EST)

You can solve many conic section problems with the free package from my web site below. Here is an example of finding a general conic section from 5 points on the curve. Needs["ConicSections`ConicSections`"] The general equation for a conic section is generaleqn = a*x^2 + b*x*y + c*y^2 + d*x + e*y == 1; We need 5 points for a conic in general position. Here is a test set. pt1 = {x -> 2, y -> 5}; pt2 = {x -> -1, y -> 4}; pt3 = {x -> -2, y -> -2}; pt4 = {x -> 4, y -> 5}; pt5 = {x -> 6, y -> -7}; We substitute into the general equation and solve the resulting equations. {generaleqn /. pt1, generaleqn /. pt2, generaleqn /. pt3, generaleqn /. pt4, generaleqn /. pt5}; sol = First@Solve[%] {a -> 1113/38866, b -> 867/38866, c -> 1001/19433, d -> -(11013/38866), e -> -(228/19433)} We then substitute the solution into the general equation to obtain the specific conic equation. I switched to approximate numbers to keep the following expressions short. coniceqn = generaleqn /. sol // N -0.28335820511501053*x + 0.02863685483455977*x^2 - 0.01173261977049349*y + 0.022307415221530385*x*y + 0.051510317501157823*y^2 == 1. Now, I finally use some routines from the ConicSections package. The first routine, ParseConic, will parse the equation and return the a semi-axis size, the eccentricity, and a parametrized curve for the conic. P is an orthogonal matrix; T is a vector translation; R is either an IdentityMatrix or a matrix that reflects in the line x == y. They tell how to go from the conic in standard position to a conic in general position. I don't show all the output. {{avalue, eccentricity}, curve[t_], {P, T, R}} = ParseConic[coniceqn]; The curve parametrization that was calculated is curve[t] = {5.35468 + 7.89892 Cos[t] + 2.10747 Sin[t], -1.04558 - 3.21402 Cos[t] + 5.17941 Sin[t]} ParametricPlot[Evaluate[curve[t]], {t, 0, 2Pi}, AspectRatio -> Automatic, Frame -> True]; We can obtain information about the conic in standard position with the command... information = StandardConic[{avalue, eccentricity}] {conictype -> "Ellipse", conicequation -> 0.013750816629065068*x^2 + 0.031981904691988035*y^2 == 1, coniccurve -> {8.527775420556369*Cos[t], 5.591751174503861*Sin[t]}, coniccurvedomain -> {-Pi, Pi}, coniccenter -> {0, 0}, conicfocus -> {{6.438576894460455, 0}, {-6.438576894460455, 0}}, conicdirectrix -> {x == -11.294880035682084, x == 11.294880035682084}, conicvertex -> {{8.527775420556369, 0}, {-8.527775420556369, 0}}} For the corresponding information about the conic in its actual position we can use... TransformEllipseRules[P, T, R][information] // Simplify {conictype -> "Ellipse", conicequation -> 1.*x^2 + x*(-9.89487870619946 + 0.7789757412398916*y) + 1.7987421383647793*(-4.5214428589005795 + y)* (4.293670631128352 + y) == 0, coniccurve -> {5.354680644977494 + 7.898924351991414*Cos[t] + 2.1074687543949118*Sin[t], -1.045581448350521 - 3.214023600560889*Cos[t] + 5.179406978295339*Sin[t]}, coniccurvedomain -> {-Pi, Pi}, coniccenter -> {5.354680644977494, -1.045581448350521}, conicfocus -> {{11.31846713274418, -3.472209387325002}, {-0.6091058427891927, 1.3810464906239601}}, conicdirectrix -> {6.413969458840636 + 1.*x == 0.4068938322913037*y, 1.*x == 17.97421203377969 + 0.4068938322913037*y}, conicvertex -> {{13.253604996968907, -4.259605048911411}, {-2.54424370701392, 2.168442152210368}}} This should work for any type of conic section in any position. David Park djmp at earthlink.net http://home.earthlink.net/~djmp/ From: 17538 at student.hhs.se [mailto:17538 at student.hhs.se] To: mathgroup at smc.vnet.net I wonder how I can use Mathematica to solve a system of non-linear equations. If I want to find the parameters of an ellipse, given coordinates of some points on it, then I have a system of four non-linear equations. Where should I go from here in order to make Mathematica help me with the substitutions necessary to express each parameter as a function of the coordinates of the points? Here follows an example: The formula of an ellipse (with an axis parallell to the x-axis of the coordinate system) is: (x-ox)^2 / s^2 + (y-oy)^2 / b^2 = 1 where s: half length of small axis. b: half length of big axis. ox: x-coordinate of the center of the ellipse. oy: y-coordinate of the center of the ellipse. x,y: coordinates of a point on perimeter of the ellipse. I solve for each of the 4 parameters s, b, ox, oy and insert unique {x,y}-values in each equation: s = f(x1,y1,b,ox,oy) b = f(x2,y2,s,ox,oy) ox = f(x3,y3,b,s,oy) oy = f(x4,y4,b,s,ox) [Of course, there are two functions for each parameter, one describing the upper half and the other describing the lower half, so to speak.] Now I should be able to substitute my way to a solution like this: Inserting b = f(x2,y2,s,ox,oy) in s = f(x1,y1,b,ox,oy) gives: s = f(x1,y1,f(x2,y2,s,ox,oy),ox,oy) Solving for s gives: s = f(x1,y1,x2,y2,ox,oy) And then continue with inserting oy = f(...) in s' = f(...) above and again solving for s, then the same with ox so that I arrive at s = f(x1,y1,x2,y2,x3,y3,x4,y4). Then the procedure is repeated with each of the other parameters b, ox, oy. The thing is that it gets quite tedious even for a pretty simple problem like this. And what about doing it for an ellipse with a tilted axis (not parallell with the x-axis of the coordinate system)? That would introduce a new parameter, a new equation in the system and significantly longer equations too (maybe it is pointless to work with algebraic solutions in that case anyway). So, what tricks can Mathematica perform in order to do the work for me? I'm grateful for any comments or reference that could help me!