3
$\begingroup$

Typically, one can restrict symbolic quantities using Assumptions.

Suppose, for instance, that I want to solve $a=x^2$ for $x$, with $a>0$, restricting my answer to the positive root.

r1 = Reduce[a == x^2, x]
Simplify[r1, Assumptions -> x > 0 && a > 0]

$x=-\sqrt{a} \,\,||\, x=\sqrt{a}$

$\sqrt{a}=x$

However, it appears that Assumptions does not work on symbolic quantities when they are expressed as Root objects. For instance, r2 contains four symbolic quantities:

expr=H3O^4 + H3O^3*(KaHA + KaHB) + 
H3O^2*(KaHA*KaHB - (cHA*KaHA + cHB*KaHB) - Kw) - 
H3O*((cHA + cHB)*KaHA*KaHB + Kw*(KaHA + KaHB)) - KaHA*KaHB*Kw;
r2= Reduce[expr==0, H3O]

enter image description here

Now suppose I want to determine if any of these gives a positive value for the variable H3O when all the parameters are postive:

sr2=Simplify[r2, Assumptions->H3O > 0 && KaHA > 0 && KaHB > 0 && cHA > 0 && cHB > 0 && Kw > 0]

The output of sr2 is identical to that of r2. Yet a single numerical test shows that the restrictions given by Assumptions are ignored. Specifically, Simplify/Assumptions returns all four Root objects for H3O, even though only one of them numerically evaluates as positive:

N[sr2 /. {KaHA -> 10^-(375/100), KaHB -> 10^-(4756/1000), cHA -> 5/1000, cHB -> 5/1000, Kw -> 10^-14}, 10]

H3O == -0.001074531007 || H3O == -0.00003186037898 || H3O == -9.999999686*10^-13 || H3O == 0.0009110246413

I tried to apply Simplify/Assumptions to another expression that also contains Root objects, and encountered the same behavior.

So does Assumptions not work for Root objects generally and, if so, why?

It's possible this is happening because Root objects need to be "unpacked" (with ToRadicals) before constraint-testing can be done by Assumptions, and Simplify doesn't unpack Root objects.

One can restrict the answer by instead applying the restrictions within Reduce:

 r3= Reduce[expr==0 && H3O > 0 && KaHA > 0 && KaHB > 0 && cHA > 0 && cHB > 0 && Kw > 0, H3O];
 N[r3 /. {KaHA -> 10^-(375/100), KaHB -> 10^-(4756/1000), cHA -> 5/1000, cHB -> 5/1000, Kw -> 10^-14}, 10]

H3O == 0.0009110246413

As to why I don't just use the latter (r3) as a workaround, r2 runs in a fraction of a second while r3 takes ~10 hours. So I was hoping that r2, followed by Simplify/Assumptions on Root objects, might take less time.

But if unpacking is needed, the time will be long either way. For instance, I could unpack r2 using ToRadicals, and then apply Simplify/Assumptions to that, but then I get back to the same problem: excessive computation time.

Hence it seems I need to either pay on the front end (in Reduce) or on the back end (ToRadicals and then Simplify/Assumptions) because it's the constraint-testing that is time-consuming, and it has to be done one place or the other.

[N.B.: I am using MMA 12.0.0.0 for MacOS.]

$\endgroup$
5
  • $\begingroup$ Block[{KaHA = 10^-(375/100), KaHB = 10^-(4756/1000), cHA = 5/1000, cHB = 5/1000, Kw = 10^-14}, Reduce[expr == 0 && H3O > 0, H3O]] // N returns one positive result nearly instantly... $\endgroup$
    – ciao
    Commented Mar 15, 2020 at 20:26
  • 1
    $\begingroup$ If you know something about the expected solutions, you might be able to use the ordering of algebraic Root objects here. Real roots come first, sorted by value, followed by complex, sorted by real part. So, it you expect three negative roots and one positive, Root[...,4] is what you want. $\endgroup$
    – John Doty
    Commented Mar 22, 2020 at 19:40
  • $\begingroup$ @JohnDoty That's nice to know. I do expect only one positive root. But: The last root isn't generally positive. It is generally positive when all the parameters are positive. For example, if I make the first param. negative, the first 2 roots are neg., and the last 2 are non-real . So, logically, how can I trust that the last root will be pos. under my parameter constraints, when I haven't told MMA what the parameter constraints are? Also, I'd still like MMA to confirm that, generically, the last root is positive, and that only that root is positive, as a check on my expression. $\endgroup$
    – theorist
    Commented Mar 22, 2020 at 19:54
  • $\begingroup$ @ciao. Sure, if I substitute in actual values for all the parameters, I can rapidly solve for the four possibe values of the variable H3O, and ask MMA to pick the positive one. But what I instead want is to more rapidly determine the general symbolic expression for H3O that satisfies my restrictions on the parameters and the variable. I.e., using my labeling, what I want is a quicker way to obtain r3. $\endgroup$
    – theorist
    Commented Mar 22, 2020 at 20:06
  • $\begingroup$ How do you imagine a simplification of Root[-KaHA KaHB Kw + (-cHA KaHA KaHB - cHB KaHA KaHB - KaHA Kw - KaHB Kw) #1 + (-cHA KaHA - cHB KaHB + KaHA KaHB - Kw) #1^2 + (KaHA + KaHB) #1^3 + #1^4 &, 1]? $\endgroup$
    – user64494
    Commented Feb 17, 2021 at 12:12

1 Answer 1

1
$\begingroup$
r1 = Reduce[a == x^2 && x > 0, x]

(* Re[a] > 0 && Im[a] == 0 && x == Sqrt[Re[a]] *)

restrict the parameters that are intended to be restricted.

From the documentation page for Simplify in the section Properties & Relations: Assuming[x > 0, Simplify[Sqrt[x^2]]] (* x *)

Simplify[Sqrt[x^2]]

(* Sqrt[x^2] *)

Simplify[Sqrt[x^2], x > 0]

(* x *)

From the documentation page for Simplify in the section Properties & Relations: Use Solve or Reduce to find solutions of systems of multivariate equations: shows example that do not contain parameters.

This gives a first impression that Reduce or Solve are not intended to give rich parameter solutions. They allow for many multivariations.

This is an experimental equation of three constituents water and A and B. During protolysis there is HA and HB reaching equilibrium. The system is in water at a certain condition {T,p,V} and there is a certain mass or amount of each of theses. The equilibrium constant Kw is a function of {T,p,V}. There are tables filling pages for this. Kw is read-off and @ciao is right use the numerical value, magnitude of it and reduce complexity for Mathematica. Then fastest is

Block[{KaHA = 10^-(375/100), KaHB = 10^-(4756/1000), cHA = 5/1000, 
   cHB = 5/1000, Kw = 10^-14}, 
  NSolve[expr == 0 && H3O > 0, H3O]] // AbsoluteTiming

(* {0.003565, {{H3O -> 0.000911025}}} *)

Block[{KaHA = 10^-(375/100), KaHB = 10^-(4756/1000), cHA = 5/1000, 
    cHB = 5/1000, Kw = 10^-14}, Roots[expr == 0, H3O]] // 
  N // AbsoluteTiming

(* {0.091388, H3O == -0.00107453 + 3.7183910^-21 I || H3O == -0.0000318604 - 5.5478510^-19 I || H3O == -1.10^-12 + 5.6608310^-19 I || H3O == 0.000911025 - 1.50162*10^-20 I} *)

gives the most exact results. Only one is positive and valid and Reals after Chop is applied.

$\endgroup$

Not the answer you're looking for? Browse other questions tagged or ask your own question.