計算機は囁く

Mathematica は当然

In[5]:= Reduce[Exists[x, 0 <= x <= 1 && y == x^2 - 2 a x], Reals]

Out[5]= (a <= 0 && 0 <= y <= 1 - 2 a) || (0 < a <= 1/2 && -a^2 <= y <= 1 - 2 a) ||
    (1/2 < a <= 1 && -a^2 <= y <= 0) || (a > 1 && 1 - 2 a <= y <= 0)

と答えます.

一方,qepmax+nnsolvexx は

(%i1) qe([[E,x]],"%and"(0<=x,x<=1,y=x^2-2*a*x));

(%o1) (((a > 0) %and (y+1 > 0) %and (y+2*a-1 <= 0))
  %or ((y >= 0) %and (y+2*a-1 <= 0))
  %or ((y <= 0) %and (y+2*a-1 >= 0))) %and (y+a^2 >= 0)

(%i2) nnsolvexx(%);

(%o2) ((a-1 < 0) %and (a > 0) %and (y+2*a-1 <= 0) %and (y+a^2 >= 0))
  %or ((y >= 0) %and (y+2*a-1 <= 0))
  %or ((y <= 0) %and (y+2*a-1 >= 0))

と囁きます.

set_goal

set_goal([`&0 < b pow 2 - &4 * a * c`],`?x. a * x pow 2 + b * x + c = &0`);;
e(ASM_CASES_TAC `a = &0`);;
e(ASM_CASES_TAC `b = &0`);;
e( X_GET_TAC 0 );;
e( ASM_SIMP_TAC [] );;
e( ARITH_TAC );;
e( EXISTS_TAC `-- c / b` );;
e( ASM_REAL_FIELD_TAC );;
e( EXISTS_TAC `(sqrt (b pow 2 - &4 * a * c) - b) / (&2 * a)` );;
e( POP_ASSUM MP_TAC );;
e( FIRST_X_ASSUM (MP_TAC o (MATCH_MP REAL_LT_IMP_LE) ) );;
e( SIMP_TAC [POW_2;GSYM SQRT_POW2] );;
e( REAL_FIELD_TAC );;
top_thm();;
# set_goal([`&0 < b pow 2 - &4 * a * c`],`?x. a * x pow 2 + b * x + c = &0`);;
val it : goalstack = 1 subgoal (1 total)

  0 [`&0 < b pow 2 - &4 * a * c`]

`?x. a * x pow 2 + b * x + c = &0`

# e(ASM_CASES_TAC `a = &0`);;
val it : goalstack = 2 subgoals (2 total)

  0 [`&0 < b pow 2 - &4 * a * c`]
  1 [`~(a = &0)`]

`?x. a * x pow 2 + b * x + c = &0`

  0 [`&0 < b pow 2 - &4 * a * c`]
  1 [`a = &0`]

`?x. a * x pow 2 + b * x + c = &0`

# e(ASM_CASES_TAC `b = &0`);;
val it : goalstack = 2 subgoals (3 total)

  0 [`&0 < b pow 2 - &4 * a * c`]
  1 [`a = &0`]
  2 [`~(b = &0)`]

`?x. a * x pow 2 + b * x + c = &0`

  0 [`&0 < b pow 2 - &4 * a * c`]
  1 [`a = &0`]
  2 [`b = &0`]

`?x. a * x pow 2 + b * x + c = &0`

# e( X_GET_TAC 0 );;
val it : goalstack = 1 subgoal (3 total)

  0 [`a = &0`]
  1 [`b = &0`]

`&0 < b pow 2 - &4 * a * c ==> (?x. a * x pow 2 + b * x + c = &0)`

# e( ASM_SIMP_TAC []  );;
val it : goalstack = 1 subgoal (3 total)

  0 [`a = &0`]
  1 [`b = &0`]

`&0 < &0 pow 2 - &4 * &0 * c ==> (?x. &0 * x pow 2 + &0 * x + c = &0)`

# e( ARITH_TAC );;
val it : goalstack = 1 subgoal (2 total)

  0 [`&0 < b pow 2 - &4 * a * c`]
  1 [`a = &0`]
  2 [`~(b = &0)`]

`?x. a * x pow 2 + b * x + c = &0`

# e( EXISTS_TAC `-- c / b` );;
val it : goalstack = 1 subgoal (2 total)

  0 [`&0 < b pow 2 - &4 * a * c`]
  1 [`a = &0`]
  2 [`~(b = &0)`]

`a * (--c / b) pow 2 + b * --c / b + c = &0`

# e( ASM_REAL_FIELD_TAC );;
2 basis elements and 0 critical pairs
5 basis elements and 3 critical pairs
5 basis elements and 2 critical pairs
5 basis elements and 1 critical pairs
5 basis elements and 0 critical pairs
val it : goalstack = 1 subgoal (1 total)

  0 [`&0 < b pow 2 - &4 * a * c`]
  1 [`~(a = &0)`]

`?x. a * x pow 2 + b * x + c = &0`

# e( EXISTS_TAC `(sqrt (b pow 2 - &4 * a * c) - b) /  (&2 * a)` );;
val it : goalstack = 1 subgoal (1 total)

  0 [`&0 < b pow 2 - &4 * a * c`]
  1 [`~(a = &0)`]

`a * ((sqrt (b pow 2 - &4 * a * c) - b) / (&2 * a)) pow 2 +
 b * (sqrt (b pow 2 - &4 * a * c) - b) / (&2 * a) +
 c =
 &0`

# e( POP_ASSUM MP_TAC );;
val it : goalstack = 1 subgoal (1 total)

  0 [`&0 < b pow 2 - &4 * a * c`]

`~(a = &0)
 ==> a * ((sqrt (b pow 2 - &4 * a * c) - b) / (&2 * a)) pow 2 +
     b * (sqrt (b pow 2 - &4 * a * c) - b) / (&2 * a) +
     c =
     &0`

# e( FIRST_X_ASSUM (MP_TAC o (MATCH_MP REAL_LT_IMP_LE) )  );;
val it : goalstack = 1 subgoal (1 total)

`&0 <= b pow 2 - &4 * a * c
 ==> ~(a = &0)
 ==> a * ((sqrt (b pow 2 - &4 * a * c) - b) / (&2 * a)) pow 2 +
     b * (sqrt (b pow 2 - &4 * a * c) - b) / (&2 * a) +
     c =
     &0`

# e( SIMP_TAC [POW_2;GSYM SQRT_POW2] );;
val it : goalstack = 1 subgoal (1 total)

`sqrt (b * b - &4 * a * c) * sqrt (b * b - &4 * a * c) = b * b - &4 * a * c
 ==> ~(a = &0)
 ==> a *
     (sqrt (b * b - &4 * a * c) - b) / (&2 * a) *
     (sqrt (b * b - &4 * a * c) - b) / (&2 * a) +
     b * (sqrt (b * b - &4 * a * c) - b) / (&2 * a) +
     c =
     &0`

# e( REAL_FIELD_TAC );;
3 basis elements and 0 critical pairs
4 basis elements and 1 critical pairs
4 basis elements and 0 critical pairs
val it : goalstack = No subgoals

# top_thm () ;;
val it : thm = &0 < b pow 2 - &4 * a * c |- ?x. a * x pow 2 + b * x + c = &0

冪乗と階乗

g(`!n. 3 < n ==> 2 EXP n < FACT n`);;
e(INDUCT_TAC);;
e(ARITH_TAC);;
e(SIMP_TAC[LT;EXP;FACT]);;
e(STRIP_TAC);;
e(POP_ASSUM (SUBST1_TAC o GSYM) THEN ARITH_TAC);;
e(MATCH_MP_TAC LESS_MULT);;
e(ASM_ARITH_TAC);;
prove
 (`!n. 3 < n ==> 2 EXP n < FACT n`,
  INDUCT_TAC THENL
  [ARITH_TAC;
   SIMP_TAC[LT;EXP;FACT] THEN STRIP_TAC THENL
   [POP_ASSUM (SUBST1_TAC o GSYM) THEN ARITH_TAC;
    MATCH_MP_TAC LESS_MULT THEN ASM_ARITH_TAC]]);;
dmtcp_checkpoint (DMTCP + MTCP) 1.2.5
Copyright (C) 2006-2011  Jason Ansel, Michael Rieker, Kapil Arya, and
                                                       Gene Cooperman
This program comes with ABSOLUTELY NO WARRANTY.
This is free software, and you are welcome to redistribute it
under certain conditions; see COPYING file for details.
(Use flag "-q" to hide this message.)

dmtcp_coordinator starting...
    Port: 7779
    Checkpoint Interval: disabled (checkpoint manually instead)
    Exit on last client: 1
Backgrounding...
# g(`!n. 3 < n ==> 2 EXP n < FACT n`);;
val it : goalstack = 1 subgoal (1 total)

`!n. 3 < n ==> 2 EXP n < FACT n`

# e(INDUCT_TAC);;
val it : goalstack = 2 subgoals (2 total)

  0 [`3 < n ==> 2 EXP n < FACT n`]

`3 < SUC n ==> 2 EXP SUC n < FACT (SUC n)`

`3 < 0 ==> 2 EXP 0 < FACT 0`

# e(ARITH_TAC);;
val it : goalstack = 1 subgoal (1 total)

  0 [`3 < n ==> 2 EXP n < FACT n`]

`3 < SUC n ==> 2 EXP SUC n < FACT (SUC n)`

# e(SIMP_TAC[LT;EXP;FACT]);;
val it : goalstack = 1 subgoal (1 total)

  0 [`3 < n ==> 2 EXP n < FACT n`]

`3 = n \/ 3 < n ==> 2 * 2 EXP n < SUC n * FACT n`

# e(STRIP_TAC);;
val it : goalstack = 2 subgoals (2 total)

  0 [`3 < n ==> 2 EXP n < FACT n`]
  1 [`3 < n`]

`2 * 2 EXP n < SUC n * FACT n`

  0 [`3 < n ==> 2 EXP n < FACT n`]
  1 [`3 = n`]

`2 * 2 EXP n < SUC n * FACT n`

# e(POP_ASSUM (SUBST1_TAC o GSYM) THEN ARITH_TAC);;
val it : goalstack = 1 subgoal (1 total)

  0 [`3 < n ==> 2 EXP n < FACT n`]
  1 [`3 < n`]

`2 * 2 EXP n < SUC n * FACT n`

# e(MATCH_MP_TAC LESS_MULT);;
val it : goalstack = 1 subgoal (1 total)

  0 [`3 < n ==> 2 EXP n < FACT n`]
  1 [`3 < n`]

`2 < SUC n /\ 2 EXP n < FACT n`

# e(ASM_ARITH_TAC);;
val it : goalstack = No subgoals

fourier_elim を懐柔する(その2)

前回の qfq では,例えば

(%i1) f:abs(x-a)=abs(x-b)+c;
(%o1)                     abs(x - a) = abs(x - b) + c
(%i2) qfq(f);
(%o2) ((0 < c) %and (a = b - c) %and (x = b))
 %or ((0 < c) %and (b = a - c) %and (x = a - c))
 %or ((a = b) %and (c = 0) %and (x = b))
 %or ((a = b - c) %and (max(b, b - c) < x))
 %or ((a = b - c) %and (c < 0) %and (x = b - c))
 %or ((a = c + b) %and (c < 0) %and (x = c + b))
 %or ((a = c + b) %and (x < min(b, c + b)))
 %or ((a = 2 x - c - b) %and (x < min(b, c + b)))
 %or ((a = 2 x + c - b) %and (max(b, b - c) < x))

のように出力に max,min が入る場合があります.

原因は

(%i3) listofvars(f);
(%o3)                            [a, x, c, b]

として変数を抽出している点,つまり,上の例の場合,fourier_elim が先ず a について解き,次に x について解き,…と進んでしまう点にあります.

そこで,変数の順序を指定できるように

qfq(f):=block([g],
if op(f)="%and" then (g:args(f),f2q(fourier_elim(g,
(if listp(fev) and length(fev)>0 then fev else listofvars(g))
)))
else f2q(fourier_elim(f,
(if listp(fev) and length(fev)>0 then fev else listofvars(f))
)))$

と改め,例えば

(%i4) fev:[c,b,a];
(%o4)                              [c, b, a]

とすれば

(%i5) qfq(f);
(%o5) ((a = x) %and (b = x) %and (c = 0))
 %or ((a = x) %and (b = x - c) %and (c < 0))
 %or ((a = x) %and (b = x + c) %and (c < 0))
 %or ((a < x) %and (b = x) %and (c = x - a))
 %or ((a < x) %and (b < x) %and (c = b - a))
 %or ((a < x) %and (c = 2 x - b - a) %and (x < b))
 %or ((b = x) %and (c = a - x) %and (x < a))
 %or ((b < x) %and (c = (- 2 x) + b + a) %and (x < a))
 %or ((c = a - b) %and (x < a) %and (x < b))

となり,これなら

(%i6) g:fqe([[E,x]],f);fqe([],(g) %eq (abs(c)<=abs(a-b)));
(%o6) ((c - b + a >= 0) %and (c + b - a <= 0))
                                   %or ((c - b + a <= 0) %and (c + b - a >= 0))
(%o7)                                true

と処理できます.

fourier_elim を懐柔する

abs や max,min も気軽に使えるよう,fourier_elim を内蔵した fqe を作りました.

f2q(f):=subst([universalset=true,emptyset=false,"or"="%or","["="%and"],f)$

qfq(f):=block([g],
if op(f)="%and" then (g:args(f),f2q(fourier_elim(g,listofvars(g))))
else f2q(fourier_elim(f,listofvars(f))))$

fqe(L,F):=
qe(L,scanmap(lambda([f],if not(atom(f) or freeof(abs,max,min,dispform(f,all))) then qfq(f) else f),F))$

動かない式を動かす(nnsolve その4)

まず,次をご覧ください.

(%i1) f:((x-1/x)^2+(y^2-4*y+5)^3+(z+1)^2/z^4<=1)%and(x<0);
                                2
                         (z + 1)      2           3        1 2
(%o1)      (x < 0) %and (-------- + (y  - 4 y + 5)  + (x - -)  <= 1)
                             4                             x
                            z

人間が見るとバレバレですが,ストレートな qepmax では

(%i2) qe([],f);
                                  2  6  4       2  5  4       2  4  4
(%o2) (x < 0) %and (z < 0) %and (x  y  z  - 12 x  y  z  + 63 x  y  z
        2  3  4        2  2  4        2    4    4  4        2  4    4    2  2
 - 184 x  y  z  + 315 x  y  z  - 300 x  y z  + x  z  + 122 x  z  + z  + x  z
      2      2
 + 2 x  z + x  = 0)

までです.

一方,nnsolve では

(%i3) nnsolve(f);
(%o3)            (x + 1 = 0) %and (y - 2 = 0) %and (z + 1 = 0)

のように解けます.

より端的な例は

(%i6) qe([],a^2+b^2+c^2+d^2=0);
                              2    2    2    2
(%o6)                        d  + c  + b  + a  = 0
(%i7) nnsolve(a^2+b^2+c^2+d^2=0);
(%o7)           (a = 0) %and (b = 0) %and (c = 0) %and (d = 0)

といったものです.

こうした「動かない式を動かす」ために作ったのが ncond です.

ncond(f):=(f)%and(
substpart("%and",map(lambda([q],qe(q,f)),
fullmap(lambda([x],[E,x]),
delete([],full_listify(powerset(setify(listofvars(f))))))),0))$

これは,所謂,必要条件を付加する関数で

(%i8) ncond(a^2+b^2+c^2+d^2=0);
                                   2    2                          2    2
(%o8) (a = 0) %and (b = 0) %and (b  + a  = 0) %and (c = 0) %and (c  + a  = 0)
        2    2             2    2    2                          2    2
 %and (c  + b  = 0) %and (c  + b  + a  = 0) %and (d = 0) %and (d  + a  = 0)
        2    2             2    2    2             2    2
 %and (d  + b  = 0) %and (d  + b  + a  = 0) %and (d  + c  = 0)
        2    2    2             2    2    2             2    2    2    2
 %and (d  + c  + a  = 0) %and (d  + c  + b  = 0) %and (d  + c  + b  + a  = 0)
(%i9) nns(qe([],%));
(%o9)          (a = 0) %and (b = 0) %and (c = 0) %and (d = 0)
(%i10) ncond(f);
(%o10) (x < 0) %and (x + 1 = 0) %and (y - 2 = 0)
        2  6       2  5       2  4        2  3        2  2        2      4
 %and (x  y  - 12 x  y  + 63 x  y  - 184 x  y  + 315 x  y  - 300 x  y + x
        2
 + 122 x  + 1 = 0) %and (z < 0) %and (z + 1 = 0)
        4  4      2  4    4    2  2      2      2
 %and (x  z  - 2 x  z  + z  + x  z  + 2 x  z + x  = 0)
        6  4       5  4       4  4        3  4        2  4          4        4
 %and (y  z  - 12 y  z  + 63 y  z  - 184 y  z  + 315 y  z  - 300 y z  + 124 z
                                  2
    2                      (z + 1)      2           3        1 2
 + z  + 2 z + 1 = 0) %and (-------- + (y  - 4 y + 5)  + (x - -)  <= 1)
                               4                             x
                              z
(%i11) nns(qe([],%));
(%o11)            (x + 1 = 0) %and (y - 2 = 0) %and (z + 1 = 0)

となります.