Das Problem ist nicht so einfach wie es aussieht (aber es
ist ein recht interessantes~).
Die Schwierigkeit besteht hier darin, dass die
Gleichungen rational sind.
Im 2.2 Mathematica Buch steht unter 'Solve' im refguide:
"Solve deals primarily with linear and polynomial systems."
^^^^^^^^^^
Der Grund, dass dieser Satz da steht ist der, dass fuer
polynomiale Systeme ein sehr schoener und maechtiger
Algorithmus, der der Groebner Basen existiert, der es
erlaubt, solche Gleichungssysteme zu loesen.
Nun ist es ja relativ einfach, aus rationelen Gleichungen
polynomiale zu machen, man bringt einfach alles auf einen
Hauptnenner
(alles folgende ist in Version 2.2 gerechnet)
In[1]:=
eqn={C12 + (C13*C23)/(C13 + C23) == (C10*C20)/(C10 + C20),
C13 + (C12*C23)/(C12 + C23) == (C10*C30)/(C10 + C30),
C23 + (C12*C13)/(C12 + C23) == (C20*C30)/(C20 + C30)};
In[2]:= eqn1=Together /@ Apply[Subtract, eqn,{1}]
Out[2]= {(C10 C12 C13 - C10 C13 C20 + C12 C13 C20 + C10
C12 C23 +
> C10 C13 C23 - C10 C20 C23 + C12 C20 C23 + C13 C20
C23) /
> ((C10 + C20) (C13 + C23)),
> (C10 C12 C13 + C10 C12 C23 + C10 C13 C23 - C10 C12
C30 + C12 C13 C30 -
> C10 C23 C30 + C12 C23 C30 + C13 C23 C30) / ((C12
+ C23) (C10 + C30)),
2
> (C12 C13 C20 + C12 C20 C23 + C20 C23 + C12 C13 C30
- C12 C20 C30 +
2
> C12 C23 C30 - C20 C23 C30 + C23 C30) / ((C12 +
C23) (C20 + C30))}
... extrahiert die neuen Zaehler ...
In[3]:= {numerators, denominators}=
Factor//@Transpose[{Numerator[#], Denominator[#]}& /@ eqn1]
Out[3]= {{C10 C12 C13 - C10 C13 C20 + C12 C13 C20 + C10
C12 C23 +
> C10 C13 C23 - C10 C20 C23 + C12 C20 C23 + C13 C20 C23,
> C10 C12 C13 + C10 C12 C23 + C10 C13 C23 - C10 C12
C30 + C12 C13 C30 -
> C10 C23 C30 + C12 C23 C30 + C13 C23 C30,
2
> C12 C13 C20 + C12 C20 C23 + C20 C23 + C12 C13 C30
- C12 C20 C30 +
2
> C12 C23 C30 - C20 C23 C30 + C23 C30},
> {(C10 + C20) (C13 + C23), (C12 + C23) (C10 + C30),
> (C12 + C23) (C20 + C30)}}
und nimmt diese als neue Gleichungen in Solve; welches
nach etwa 1 s (auf einem 200 Mz Pentium Pro) die folgenden
Loesungen liefert:
In[4]:= (sol1=Solve[#==0& /@ numerators,
{C10, C20, C30}])//InputForm
Out[4]//InputForm=
{{C10 -> 0, C20 -> 0, C30 -> 0},
{C10 -> (4*C12^4*C13^4 - 2*C12^3*C13^5 +
16*C12^4*C13^3*C23 +
4*C12^3*C13^4*C23 - 4*C12^2*C13^5*C23 +
24*C12^4*C13^2*C23^2 +
32*C12^3*C13^3*C23^2 - 2*C12^2*C13^4*C23^2 -
2*C12*C13^5*C23^2 +
16*C12^4*C13*C23^3 + 52*C12^3*C13^2*C23^3 +
26*C12^2*C13^3*C23^3 -
4*C12*C13^4*C23^3 + 4*C12^4*C23^4 +
34*C12^3*C13*C23^4 +
46*C12^2*C13^2*C23^4 + 10*C12*C13^3*C23^4 -
2*C13^4*C23^4 +
8*C12^3*C23^5 + 26*C12^2*C13*C23^5 +
20*C12*C13^2*C23^5 +
2*C13^3*C23^5 + 4*C12^2*C23^6 + 8*C12*C13*C23^6 +
4*C13^2*C23^6)/
((2*C12*C13 - C13^2 + 2*C12*C23 + C13*C23 + 2*C23^2)*
(C12^2*C13^3 + 2*C12^2*C13^2*C23 + C12*C13^3*C23 +
3*C12^2*C13*C23^2 + C12*C13^2*C23^2 +
2*C12^2*C23^3 +
4*C12*C13*C23^3 + 2*C12*C23^4 + 2*C13*C23^4)),
C20 -> (2*(C12*C13 + C12*C23 + C13*C23)*(C12*C13 +
C12*C23 + C23^2))/
(C13*(C12*C13 + 3*C12*C23 + 2*C23^2)),
C30 -> (2*(C12*C13 + C12*C23 + C13*C23)*(C12*C13 +
C12*C23 + C23^2))/
(C12*(2*C12*C13 - C13^2 + 2*C12*C23 + C13*C23 +
2*C23^2))}}
Wenn man diese Loesungen jetzt wieder in die
urspruenglichen Gleichungen eqn substituiert, dann erhaelt
man:
In[5]:= Together //@ (eqn /. sol1)
(* some messages deleted *)
C12 C13 + C12 C23 + C13 C23
Out[5]= {{--------------------------- == Indeterminate,
C13 + C23
C12 C13 + C12 C23 + C13 C23
> --------------------------- == Indeterminate,
C12 + C23
2
C12 C13 + C12 C23 + C23
> ------------------------ == Indeterminate}, {True,
True, True}}
C12 + C23
D.h. die erste der Loesungen war keine, weil sie nicht
mit den Nennern konsistent war.
Nun ist der soeben beschriebene Weg des Loesens von
Systemen rationaler Gleichungen aber nicht der, wie man
algorithmisch Systeme rationaler Funktionen loest, d.h.
man bringt nicht alles auf einen Nenner und beschaeftigt
sich nur noch mit den Zaehlern, sondern man bringt die
Bedingung, dass die Nenner niemals verschwinden duerfen
mit in die Rechnung ein.
Ein algorithmischer Zugang ist viel saunberer und erlaubt
alle Loesungen mit der richtigen Multiplizitaet zu
finden. (Obwohl im Fall der Gleichungen eqn der
schnellere.) Insbesonders im Falle, dass das durch die
Gleichungen beschriebene Ideal nicht nulldimensional ist
(d.h. wir haben freie Parameter in der Loesung), ist der
algorithmische Zugang der bessere (das System eqn hat fuer
fast alle Realisierungen der Variablen C12, C13, C23 ein
nulldimensionales Ideal).
Technisch macht man dies (das Einbeziehen der Nenner
folgendermassen) folgendermassen
(Fuer Details siehe Kapitel 3.3 von
D. Cox, J. Little, D. O'Shea. Ideals, Varieties and
Algorithms, Springer-Verlag, New York, 1992):
Man nimmt zu den Zaehlergleichungen die Gleichung
1 - Product[nenner_i, {i, anzahlDerGleichunen}]*y
hinzu (hier ist y eine neue, unabhaengige Variable).
Damit wird aber typischerweise der Berechung der
lexikographischen Groebnerbasis (die braucht man um das
Gleichungssytem auf "Dreiecksform" zu bringen) _viel_
schwieriger (y ist eine neue Variable fuer die
lexikographische Order).
(Ab jetzt sind die Mathematica Rechnungen in 3.0, da die
GroebnerBasis hier viel besser ist):
Sehen wir uns also mal explizit die beiden GroebnerBasen
an. Hier die GroebnerBasis, die durch die Zaehler
definiert wird (sie kann recht schnell berechner werden
ist ist von maessiger Groesse):
In[12]:=
(gb1=GroebnerBasis[{numEqn}, {C10, C20, C30}]); // Timing
Out[12]= {0.12 Second, Null}
In[13]:= ByteCount[gb1]
Out[13]= 11800
In[14]:= Length[gb1]
Out[14]= 9
Wenn wir aber jetzt die obige Bedingung - dass die Nenner
nicht verschwinden - mit einbeziehen erhalten wir:
In[19]:= (gb2=GroebnerBasis[
Append[numerators,
1-(Times@@denominators) y], {y, C10, C20, C30}]);
// Timing
Out[19]= {877.1 Second, Null}
In[20]:= ByteCount[gb2]
Out[20]= 5799604
In[21]:= Length[gb2]
Out[21]= 159
D.h. die Zeit diese GroebnerBasis zu brechnen war fast 4
Groessenordnungen groesser und auch die GroebnerBasis
selbst ist fast 500 mal so gross. Mathematica 2.2 ist
nicht in der Lage diese GroebnerBasis in akzeptabler Zeit
zu berechnen.
Mit der Methode, die Gleichungen auf polynomiale Form zu
bekommen, und die Konsistenz mit den Nenner "bei hand" zu
ueberpruefen, kann man das obige Gleichungssystem auch
nach {C12, C13, C23} problemlos aufloesen (jetzt wieder in
Version 2.2).
In[6]:= (sol2=Solve[#==0& /@ numerators,
{C12, C13, C23}]); //Timing
Out[6]= {16.12 Second, Null}
Und nach etwas Vereinfachen sind die Loesungen auch nicht
mehr zu gross.
In[7]:= (sol2 = Simplify[Together //@ sol2]);
In[8]:= ByteCount[sol2]
Out[8]= 21464
Allerdings haben wir auch hier wieder Loesungen dabei,
die auch die Nenner zum Verschwinden bringen.
In[11]:= Together //@ (eqn/.sol2)
(* again some messages deleted *)
C10 C20
C10 C30
Out[11]= {{Indeterminate == ---------, Indeterminate ==
---------,
C10 + C20
C10 + C30
C20 C30
> Indeterminate == ---------},
C20 + C30
C10 C20 C10 C30
> {Indeterminate == ---------, Indeterminate == ---------,
C10 + C20 C10 + C30
C20 C30
> Indeterminate == ---------}, {True, True, True},
{True, True, True}}
C20 + C30
Michael Trott
Wolfram Research
|