Ticket #10983: trac_10983.patch

File trac_10983.patch, 31.8 KB (added by casamayou, 6 years ago)
  • new file sage/tests/french_book/calculus.py

    # HG changeset patch
    # User Alexandre Casamayou <alexandre.casamayou@orange.fr>
    # Date 1314277883 -7200
    # Node ID 8f9ee6055e9eb5c109fa1d4d6327b14fdab2fef7
    # Parent  2a2abbcad325ccca9399981ceddf5897eb467e64
    #10983 new doctest for french book about Sage
    <alexandre.casamayou@orange.fr> HG: branch 'default' HG: added sage/tests/french_book/calculus.py HG: added sage/tests/french_book/graphique_1.py HG:
    added sage/tests/french_book/graphique_2.py HG: added sage/tests/french_book/graphique_3.py HG: added sage/tests/french_book/sol_calculus.py HG:
    added sage/tests/french_book/sol_graphiques_1.py HG: added sage/tests/french_book/sol_graphiques_2.py
    
    diff -r 2a2abbcad325 -r 8f9ee6055e9e sage/tests/french_book/calculus.py
    - +  
     1r"""
     2sage: var ('a,x') ; y = cos(x+a)*(x+1) ; y
     3(a, x)
     4(x + 1)*cos(a + x)
     5sage: var ('a,x') ; y = cos(x+a)*(x+1) ; y
     6(a, x)
     7(x + 1)*cos(a + x)
     8sage: y.subs(a=-x) ; y.subs(x=pi/2,a=pi/3) ; y.subs(x=0.5,a=2.3)
     9x + 1
     10-1/4*(pi + 2)*sqrt(3)
     11-1.41333351100299
     12sage: y(a=-x) ; y(x=pi/2,a=pi/3) ; y(x=0.5,a=2.3)
     13x + 1
     14-1/4*(pi + 2)*sqrt(3)
     15-1.41333351100299
     16sage: var('x,y,z') ; q = x*y+y*z+z*x
     17(x, y, z)
     18sage: bool (q(x=y,y=z,z=x)==q), bool(q(z=y)(y=x) == 3*x^2)
     19(True, True)
     20sage: var('y z'); f = x^3+y^2+z; f.subs_expr(x^3==y^2, z==1)
     21(y, z)
     222*y^2 + 1
     23sage: f(x)=(2*x+1)^3 ; f(-3)
     24  -125
     25sage: f(x).expand()
     268*x^3 + 12*x^2 + 6*x + 1
     27sage: var('y'); u = sin(x) + x*cos(y)
     28y
     29sage: v = u.function(x, y); v
     30(x, y) |--> x*cos(y) + sin(x)
     31sage: w(x, y) = u; w
     32(x, y) |--> x*cos(y) + sin(x)
     33sage: (x^x/x).simplify()
     34x^(x - 1)
     35sage: f = (e^x-1)/(1+e^(x/2)); f.simplify_exp()
     36e^(1/2*x) - 1
     37sage: f = cos(x)^6 + sin(x)^6 + 3 * sin(x)^2 * cos(x)^2
     38sage: f.simplify_trig()
     391
     40sage: f = cos(x)^6; f.reduce_trig()
     4115/32*cos(2*x) + 3/16*cos(4*x) + 1/32*cos(6*x) + 5/16
     42sage: f = sin(5 * x); f.expand_trig()
     43sin(x)^5 - 10*sin(x)^3*cos(x)^2 + 5*sin(x)*cos(x)^4
     44sage: n = var('n'); f = factorial(n+1)/factorial(n)
     45sage: f.simplify_factorial()
     46n + 1
     47sage: f = sqrt(x^2); f.simplify_radical()
     48abs(x)
     49sage: f = log(x*y); f.simplify_radical()
     50log(x) + log(y)
     51sage: assume(x > 0); bool(sqrt(x^2) == x)
     52True
     53sage: forget(x > 0); bool(sqrt(x^2) == x)
     54False
     55sage: var('n'); assume(n, 'integer'); sin(n*pi).simplify()
     56n
     570
     58sage: z, phi = var('z, phi')
     59sage: solve(z**2 -2 / cos(phi) * z + 5 / cos(phi) ** 2 - 4, z)
     60[z == -(2*sqrt(cos(phi)^2 - 1) - 1)/cos(phi), z == (2*sqrt(cos(phi)^2 - 1) + 1)/cos(phi)]
     61sage: var('y'); solve(y^6==y, y)
     62y
     63[y == e^(2/5*I*pi), y == e^(4/5*I*pi), y == e^(-4/5*I*pi), y == e^(-2/5*I*pi), y == 1, y == 0]
     64sage: solve(x^2-1, x, solution_dict=True)
     65[{x: -1}, {x: 1}]
     66sage: solve([x+y == 3, 2*x+2*y == 6],x,y)
     67[[x == -r1 + 3, y == r1]]
     68sage: solve([cos(x)*sin(x) == 1/2, x+y == 0], x, y)
     69[[x == 1/4*pi + pi*z30, y == -1/4*pi - pi*z30]]
     70sage: solve(x^2+x-1>0,x)
     71[[x < -1/2*sqrt(5) - 1/2], [x > 1/2*sqrt(5) - 1/2]]
     72sage: x, y, z = var('x, y, z')
     73sage: solve([x^2 * y * z == 18, x * y^3 * z == 24,\
     74             x * y * z^4 == 6], x, y, z)
     75[[x == 3, y == 2, z == 1], [x == (1.33721506733 - 2.68548987407*I), y == (-1.70043427146 + 1.05286432575*I), z == (0.932472229404 - 0.361241666187*I)], [x == (1.33721506733 + 2.68548987407*I), y == (-1.70043427146 - 1.05286432575*I), z == (0.932472229404 + 0.361241666187*I)], [x == (-2.55065140719 - 1.57929648863*I), y == (-0.547325980144 + 1.92365128635*I), z == (-0.982973099684 - 0.183749517817*I)], [x == (-2.55065140719 + 1.57929648863*I), y == (-0.547325980144 - 1.92365128635*I), z == (-0.982973099684 + 0.183749517817*I)], [x == (0.27680507839 - 2.98720252889*I), y == (1.47801783444 - 1.34739128729*I), z == (-0.85021713573 - 0.526432162877*I)], [x == (0.27680507839 + 2.98720252889*I), y == (1.47801783444 + 1.34739128729*I), z == (-0.85021713573 + 0.526432162877*I)], [x == (-0.820988970216 + 2.88547692952*I), y == (-1.20526927276 - 1.59603445456*I), z == (0.0922683594633 - 0.995734176295*I)], [x == (-0.820988970216 - 2.88547692952*I), y == (-1.20526927276 + 1.59603445456*I), z == (0.0922683594633 + 0.995734176295*I)], [x == (-1.80790390914 - 2.39405168184*I), y == (0.891476711553 - 1.79032658271*I), z == (0.739008917221 - 0.673695643647*I)], [x == (-1.80790390914 + 2.39405168184*I), y == (0.891476711553 + 1.79032658271*I), z == (0.739008917221 + 0.673695643647*I)], [x == (2.21702675166 + 2.02108693094*I), y == (1.86494445881 + 0.722483332374*I), z == (-0.273662990072 - 0.961825643173*I)], [x == (2.21702675166 - 2.02108693094*I), y == (1.86494445881 - 0.722483332374*I), z == (-0.273662990072 + 0.961825643173*I)], [x == (2.79741668821 - 1.08372499856*I), y == (-1.96594619937 + 0.367499035633*I), z == (-0.602634636379 - 0.79801722728*I)], [x == (2.79741668821 + 1.08372499856*I), y == (-1.96594619937 - 0.367499035633*I), z == (-0.602634636379 + 0.79801722728*I)], [x == (-2.94891929905 + 0.55124855345*I), y == (0.184536718927 + 1.99146835259*I), z == (0.445738355777 - 0.895163291355*I)], [x == (-2.94891929905 - 0.55124855345*I), y == (0.184536718927 - 1.99146835259*I), z == (0.445738355777 + 0.895163291355*I)]]
     76sage: expr = sin(x) + sin(2 * x) + sin(3 * x)
     77sage: solve(expr, x)
     78[sin(3*x) == -sin(2*x) - sin(x)]
     79sage: find_root(expr, 0.1, pi)
     802.0943951023931957
     81sage: f = expr.simplify_trig(); f
     822*(2*cos(x)^2 + cos(x))*sin(x)
     83sage: solve(f, x)
     84[x == 0, x == 2/3*pi, x == 1/2*pi]
     85sage: (x^3+2*x+1).roots(x)
     86[(-1/2*(I*sqrt(3) + 1)*(1/18*sqrt(3)*sqrt(59) - 1/2)^(1/3) - 1/3*(I*sqrt(3) - 1)/(1/18*sqrt(3)*sqrt(59) - 1/2)^(1/3), 1), (-1/2*(-I*sqrt(3) + 1)*(1/18*sqrt(3)*sqrt(59) - 1/2)^(1/3) - 1/3*(-I*sqrt(3) - 1)/(1/18*sqrt(3)*sqrt(59) - 1/2)^(1/3), 1), ((1/18*sqrt(3)*sqrt(59) - 1/2)^(1/3) - 2/3/(1/18*sqrt(3)*sqrt(59) - 1/2)^(1/3), 1)]
     87sage: (x^3+2*x+1).roots(x, ring=RR)
     88 [(-0.453397651516404, 1)]
     89sage: (x^3+2*x+1).roots(x, ring=CC)
     90 [(-0.453397651516404, 1), (0.226698825758202 - 1.46771150871022*I, 1), (0.226698825758202 + 1.46771150871022*I, 1)]
     91sage: var('k n'); sum(k, k, 1, n).factor()
     92 (k, n)
     93 1/2*(n + 1)*n
     94sage: var('n y'); sum(binomial(n,k) * x^k * y^(n-k), k, 0, n)
     95 (n, y)
     96 (x + y)^n
     97sage: var('k n'); sum(binomial(n,k), k, 0, n),\
     98                  sum(k * binomial(n, k), k, 0, n),\
     99                  sum((-1)^k*binomial(n,k), k, 0, n)
     100(k, n)
     101(2^n, n*2^(n - 1), 0)
     102sage: var('a q'); sum(a*q^k, k, 0, n)
     103(a, q)
     104(a*q^(n + 1) - a)/(q - 1)
     105sage: assume(abs(q) < 1)
     106sage: sum(a*q^k, k, 0, infinity)
     107 -a/(q - 1)
     108sage: limit((x**(1/3) -2) / ((x + 19)**(1/3) - 3), x = 8)
     1099/4
     110sage: f(x) = (cos(pi/4-x)-tan(x))/(1-sin(pi/4 + x))
     111sage: limit(f(x), x = pi/4, dir='minus')
     112+Infinity
     113sage: limit(f(x), x = pi/4, dir='plus')
     114-Infinity
     115sage: u(n) = n^100 / 100.^n
     116sage: u(1);u(2);u(3);u(4);u(5);u(6);u(7);u(8);u(9);u(10)
     1170.0100000000000000
     1181.26765060022823e26
     1195.15377520732011e41
     1201.60693804425899e52
     1217.88860905221012e59
     1226.53318623500071e65
     1233.23447650962476e70
     1242.03703597633449e74
     1252.65613988875875e77
     1261.00000000000000e80
     127sage: plot(u(x), x, 1, 40)
     128sage: v(x) = diff(u(x), x); sol = solve(v(x) == 0, x); sol
     129[x == 0, x == (268850/12381)]
     130sage: sol[1].rhs().n().floor()
     13121
     132sage: limit(u(n), n=infinity)
     1330
     134sage: n0 = find_root(u(n) - 1e-8 == 0, 22, 1000); n0
     135105.07496210187252
     136sage: taylor((1+arctan(x))**(1/x), x, 0, 3)
     1371/16*x^3*e + 1/8*x^2*e - 1/2*x*e + e
     138sage: (ln(2* sin(x))).series(x==pi/6, 3)
     139(sqrt(3))*(-1/6*pi + x) + (-2)*(-1/6*pi + x)^2 + Order(-1/216*(pi - 6*x)^3)
     140sage: (ln(2* sin(x))).series(x==pi/6, 3).truncate()
     141-1/18*(pi - 6*x)^2 - 1/6*(pi - 6*x)*sqrt(3)
     142sage: taylor((x**3+x)**(1/3) - (x**3-x)**(1/3), x, infinity, 2)
     1432/3/x
     144sage: tan(4*arctan(1/5)).simplify_trig()
     145120/119
     146sage: tan(pi/4+arctan(1/239)).simplify_trig()
     147120/119
     148sage: f = arctan(x).series(x, 10); f
     1491*x + (-1/3)*x^3 + 1/5*x^5 + (-1/7)*x^7 + 1/9*x^9 + Order(x^10)
     150sage: (16*f.subs(x==1/5) - 4*f.subs(x==1/239)).n(); pi.n()
     1513.14159268240440
     1523.14159265358979
     153sage: var('k'); sum(1/k^2, k, 1, infinity),\
     154                sum(1/k^4, k, 1, infinity),\
     155                sum(1/k^5, k, 1, infinity)
     156k
     157(1/6*pi^2, 1/90*pi^4, zeta(5))
     158sage: s = 2*sqrt(2)/9801*(sum((factorial(4*k))*(1103+26390*k)\
     159     /((factorial(k)) ^ 4 * 396 ^ (4 * k)) for k in (0..11)))
     160sage: (1 / s).n(digits=100); pi.n(digits=100)
     1613.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342121432
     1623.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117068
     163sage: print "%e" % (pi - 1 / s).n(digits=100)
     164-4.364154e-96
     165sage: var('n'); u = sin(pi*(sqrt(4*n^2+1)-2*n))
     166n
     167sage: taylor(u, n, infinity, 3)
     1681/4*pi/n - 1/384*(6*pi + pi^3)/n^3
     169sage: diff(sin(x^2), x)
     1702*x*cos(x^2)
     171sage: function('f',x); function('g',x); diff(f(g(x)), x)
     172f(x)
     173g(x)
     174D[0](f)(g(x))*D[0](g)(x)
     175sage: diff(ln(f(x)), x)
     176D[0](f)(x)/f(x)
     177sage: f(x,y) = x*y + sin(x^2) + e^(-x); derivative(f, x)
     178(x, y) |--> 2*x*cos(x^2) + y - e^(-x)
     179sage: derivative(f, y)
     180(x, y) |--> x
     181sage: var('x y'); f = ln(x**2+y**2) / 2
     182(x, y)
     183sage: delta = diff(f,x,2)+diff(f,y,2)
     184sage: delta.simplify_full()
     1850
     186sage: sin(x).integral(x, 0, pi/2)
     1871
     188sage: integrate(1/(1+x^2), x)
     189arctan(x)
     190sage: integrate(1/(1+x^2), x, -infinity, infinity)
     191pi
     192sage: integrate(exp(-x**2), x, 0, infinity)
     1931/2*sqrt(pi)
     194sage: u = var('u'); f = x * cos(u) / (u^2 + x^2)
     195sage: assume(x>0); f.integrate(u, 0, infinity)
     1961/2*pi*e^(-x)
     197sage: forget(); assume(x<0); f.integrate(u, 0, infinity)
     198-1/2*pi*e^x
     199sage: integral_numerical(sin(x)/x, 0, 1)
     200(0.946083070367182..., 1.050363207929708...e-14)
     201sage: g = integrate(exp(-x**2), x, 0, infinity); g, g.n()
     202(1/2*sqrt(pi), 0.886226925452...)
     203sage: approx = integral_numerical(exp(-x**2), 0, infinity)
     204sage: approx[0]
     2050.88622692545275...
     206sage: A = matrix(QQ,[[1,2],[3,4]]); A
     207[1 2]
     208[3 4]
     209sage: A = matrix(QQ, [[2, 4, 3],[-4,-6,-3],[3,3,1]])
     210sage: A.characteristic_polynomial()
     211x^3 + 3*x^2 - 4
     212sage: A.eigenvalues()
     213[1, -2, -2]
     214sage: A.minimal_polynomial().factor()
     215(x - 1) * (x + 2)^2
     216sage: A.eigenvectors_right()
     217[(1, [
     218(1, -1, 1)
     219], 1), (-2, [
     220(1, -1, 0)
     221], 2)]
     222sage: A.jordan_form(transformation=True)
     223(
     224[ 1| 0  0]
     225[--+-----]  [ 1  1  1]
     226[ 0|-2  1]  [-1 -1  0]
     227[ 0| 0 -2], [ 1  0 -1]
     228)
     229sage: A = matrix(QQ, [[1,-1/2],[-1/2,-1]])
     230sage: A.minimal_polynomial()
     231x^2 - 5/4
     232sage: R = QQ[sqrt(5)]
     233sage: A = A.change_ring(R)
     234sage: A.jordan_form(transformation=True, subdivide=False)
     235(
     236[ 1/2*sqrt5          0]  [         1          1]
     237[         0 -1/2*sqrt5], [-sqrt5 + 2  sqrt5 + 2]
     238)
     239sage: K.<sqrt2> = NumberField(x^2 - 2)
     240sage: L.<sqrt3> = K.extension(x^2 - 3)
     241sage: A = matrix(L, [[2, sqrt2*sqrt3, sqrt2], [sqrt2*sqrt3, 3, sqrt3], [sqrt2, sqrt3, 1]])
     242sage: A.jordan_form(transformation=True)
     243(
     244[6|0|0]
     245[-+-+-]
     246[0|0|0]  [              1               1               0]
     247[-+-+-]  [1/2*sqrt2*sqrt3               0               1]
     248[0|0|0], [      1/2*sqrt2          -sqrt2          -sqrt3]
     249)
     250
     251"""
  • new file sage/tests/french_book/graphique_1.py

    diff -r 2a2abbcad325 -r 8f9ee6055e9e sage/tests/french_book/graphique_1.py
    - +  
     1r"""
     2sage: g = plot(x * sin(1/x), x, -2, 2, plot_points=500)
     3sage: def p(x, n): return(taylor(sin(x), x, 0, n))
     4sage: xmax = 15 ; n = 15
     5sage: g = plot(sin(x), x, -xmax, xmax)
     6sage: for d in range(n):
     7...     g += plot(p(x, 2 * d + 1), x, -xmax, xmax,\
     8...       color=(1.7*d/(2*n), 1.5*d/(2*n), 1-3*d/(4*n)))
     9sage: # g.show(ymin=-2, ymax=2)
     10sage: a = animate([[sin(x), taylor(sin(x), x, 0, 2*k+1)]\
     11...          for k in range(0, 14)], xmin=-14, xmax=14,\
     12...          ymin=-3, ymax=3, figsize=[8, 4])
     13sage: # a.show()
     14sage: f2(x) = 1; f1(x) = -1
     15sage: f = Piecewise([[(-pi,0),f1],[(0,pi),f2]])
     16sage: S = f.fourier_series_partial_sum(20,pi); S
     174/3*sin(3*x)/pi + 4/5*sin(5*x)/pi + 4/7*sin(7*x)/pi + 4/9*sin(9*x)/pi + 4/11*sin(11*x)/pi + 4/13*sin(13*x)/pi + 4/15*sin(15*x)/pi + 4/17*sin(17*x)/pi + 4/19*sin(19*x)/pi + 4*sin(x)/pi
     18sage: g = plot(S, x, -8, 8, color='blue')
     19sage: scie(x) = x - 2 * pi * floor((x + pi) / (2 * pi))
     20sage: g += plot(scie(x) / abs(scie(x)), x, -8, 8, color='red')
     21sage: t = var('t')
     22sage: x = cos(t) + cos(7*t)/2 + sin(17*t)/3
     23sage: y = sin(t) + sin(7*t)/2 + cos(17*t)/3
     24sage: g = parametric_plot( (x, y), (t, 0, 2*pi))
     25sage: # g.show(aspect_ratio=1)
     26sage: t = var('t'); e, n = 2, 20/19
     27sage: g1 = polar_plot(1+e*cos(n*t),(t,0,n*38*pi),plot_points=5000)
     28sage: e, n = 1/3, 20/19
     29sage: g2 = polar_plot(1+e*cos(n*t),(t,0,n*38*pi),plot_points=5000)
     30sage: # g1.show(aspect_ratio=1); g2.show(aspect_ratio=1)
     31sage: z = var('z'); g1 = complex_plot(abs(cos(z^4)) - 1,\
     32...             (-3, 3), (-3, 3), plot_points=400)
     33sage: f = lambda x, y : (abs(cos((x + I * y) ** 4)) - 1)
     34sage: g2 = implicit_plot(f, (-3, 3), (-3, 3), plot_points=400)
     35sage: # g1.show(aspect_ratio=1); g2.show(aspect_ratio=1)
     36sage: f(z) = z^5 + z - 1 + 1/z
     37sage: g = complex_plot(f,\
     38...      (-3, 3), (-3, 3))
     39
     40"""
  • new file sage/tests/french_book/graphique_2.py

    diff -r 2a2abbcad325 -r 8f9ee6055e9e sage/tests/french_book/graphique_2.py
    - +  
     1r"""
     2sage: g = bar_chart([randrange(15) for i in range(20)], color='red')
     3sage: g = bar_chart([x^2 for x in range(1,20)], width=0.2)
     4sage: liste = [10 + floor(10*sin(i)) for i in range(100)]
     5sage: g = bar_chart(liste)
     6sage: g = finance.TimeSeries(liste).plot_histogram(bins=20)
     7sage: from random import *
     8sage: n, l, x, y = 10000, 1, 0, 0; p = [[0, 0]]
     9sage: for k in range(n):
     10...    theta = (2 * pi * random()).n(digits=5)
     11...    x, y = x + l * cos(theta), y + l * sin(theta)
     12...    p.append([x, y])
     13sage: g1 = line([p[n], [0, 0]], color='red', thickness=2)
     14sage: g1 += line(p, thickness=.4); # g1.show(aspect_ratio=1)
     15sage: length = 200; n = var('n')
     16sage: u(n) = n * sqrt(2)
     17sage: z(n) = exp(2 * I * pi * u(n))
     18sage: vertices = [CC(0, 0)]
     19sage: for n in range(1, length):
     20...    vertices.append(vertices[n - 1] + CC(z(n)))
     21sage: g = line(vertices); # g.show(aspect_ratio=1)
     22sage: x = var('x'); y = function('y',x)
     23sage: DE = x*diff(y, x) == 2*y + x^3
     24sage: desolve(DE, [y,x])
     25(c + x)*x^2
     26sage: sol = []
     27sage: for i in srange(-2, 2, 0.2):
     28...    sol.append(desolve(DE, [y, x], ics=[1, i]))
     29...    sol.append(desolve(DE, [y, x], ics=[-1, i]))
     30sage: g = plot(sol, x, -2, 2)
     31sage: y = var('y')
     32sage: g += plot_vector_field((x, 2*y+x^3), (x, -2, 2), (y, -1, 1))
     33sage: # g.show(ymin=-1, ymax=1)
     34sage: # debut variante
     35sage: x = var('x'); y = function('y',x)
     36sage: DE = x*diff(y, x) == 2*y + x^3
     37sage: g = Graphics()
     38sage: for i in srange(-1, 1, 0.1):
     39...       g += line(desolve_rk4(DE, y, ics=[1, i],\
     40...                         step=0.05, end_points=[0,2]))
     41...       g += line(desolve_rk4(DE, y, ics=[-1, i],\
     42...                         step=0.05, end_points=[-2,0]))
     43sage: y = var('y')
     44sage: g = plot_vector_field((x, 2*y + x^3), (x,-2,2), (y,-1,1))
     45sage: # g.show(ymin=-1, ymax=1)
     46sage: # fin variante
     47
     48"""
  • new file sage/tests/french_book/graphique_3.py

    diff -r 2a2abbcad325 -r 8f9ee6055e9e sage/tests/french_book/graphique_3.py
    - +  
     1r"""
     2sage: import scipy; from scipy import integrate
     3sage: f = lambda y, t: - cos(y * t)
     4sage: t = srange(0, 5, 0.1); p = Graphics()
     5sage: for k in srange(0, 10, 0.15):
     6...      y = integrate.odeint(f, k, t)
     7...      p += line(zip(t, flatten(y)))
     8sage: t = srange(0, -5, -0.1); q = Graphics()
     9sage: for k in srange(0, 10, 0.15):
     10...      y = integrate.odeint(f, k, t)
     11...      q += line(zip(t, flatten(y)))
     12sage: y = var('y')
     13sage: v = plot_vector_field((1, -cos(x * y)), (x,-5,5), (y,-2,11))
     14sage: g = p + q + v; # g.show()
     15sage: import scipy; from scipy import integrate
     16sage: a, b, c, d = 1., 0.1, 1.5, 0.75
     17sage: def dX_dt(X, t=0):
     18...    return [ a*X[0] -   b*X[0]*X[1] ,
     19...            -c*X[1] + d*b*X[0]*X[1] ]
     20sage: t = srange(0, 15, .01)                     # echelle de temps
     21sage: X0 = [10, 5]  # conditions initiales : 10 lapins et 5 renards
     22sage: X = integrate.odeint(dX_dt, X0, t)     # resolution numerique
     23sage: lapins, renards =  X.T         # raccourcis de  X.transpose()
     24sage: p = line(zip(t, lapins), color='red') # trace du nb de lapins
     25sage: p += text("Lapins",(12,37), fontsize=10, color='red')
     26sage: p += line(zip(t, renards), color='blue')# idem pr les renards
     27sage: p += text("Renards",(12,7), fontsize=10, color='blue')
     28sage: p.axes_labels(["temps", "population"]); # p.show(gridlines=True)
     29sage: ### Deuxieme graphique :
     30sage: n = 11;  L = srange(6, 18, 12 / n); R=srange(3, 9, 6 / n)
     31sage: CI = zip(L, R)                # liste des conditions initiales
     32sage: def g(x,y):
     33...       v = vector(dX_dt([x, y]))  # pour un trace plus lisible,
     34...       return v/v.norm()          # on norme le champ de vecteurs
     35sage: x, y = var('x, y')
     36sage: q = plot_vector_field(g(x, y), (x, 0, 60), (y, 0, 36))
     37sage: for j in range(n):
     38...       X = integrate.odeint(dX_dt, CI[j], t)        # resolution
     39...       q += line(X, color=hue(.8-float(j)/(1.8*n))) #  graphique
     40sage: q.axes_labels(["lapins","renards"]); # q.show()
     41sage: x, y, t = var('x, y, t')
     42sage: alpha(t) = 1; beta(t) = t/2; gamma(t) = t + t**3/8
     43sage: env = solve([alpha(t)*x+beta(t)*y==gamma(t),\
     44...         diff(alpha(t),t)*x+diff(beta(t),t)*y==diff(gamma(t),t)],\
     45...         [x,y])
     46sage: f = lambda x:x^2 / 4
     47sage: p = plot(f, -8, 8, rgbcolor=(0.2,0.2,0.4)) # trace la parabole
     48sage: for u in srange(0, 8, 0.1): # trace des normales a la parabole
     49...       p += line([[u, f(u)], [-8*u, f(u) + 18]], thickness=.3)
     50...       p += line([[-u, f(u)], [8*u, f(u) + 18]], thickness=.3)
     51sage: p += parametric_plot((env[0][0].rhs(),env[0][1].rhs()),\
     52...       (t, -8, 8),color='red')              # trace la developpee
     53sage: # p.show(xmin=-8, xmax=8, ymin=-1, ymax=12, aspect_ratio=1)
     54sage: t = var('t'); p = 2
     55sage: x(t) = t; y(t) = t^2 / (2 * p)
     56sage: f(t) = [x(t), y(t)]
     57sage: df(t) = [x(t).diff(t), y(t).diff(t)]
     58sage: d2f(t) = [x(t).diff(t, 2), y(t).diff(t, 2)]
     59sage: T(t) = [df(t)[0] / df(t).norm(), df[1](t) / df(t).norm()]
     60sage: N(t) = [-df(t)[1] / df(t).norm(), df[0](t) / df(t).norm()]
     61sage: R(t) =  (df(t).norm())^3 / \
     62...        (df(t)[0] * d2f(t)[1] -df(t)[1] * d2f(t)[0])
     63sage: Omega(t) = [f(t)[0] + R(t)*N(t)[0], f(t)[1] + R(t)*N(t)[1]]
     64sage: g = parametric_plot(f(t), (t, -8, 8), color='green', thickness=2)
     65sage: for u in srange(.4, 4, .2):
     66...    g += line([f(t = u), Omega(t = u)], color='red', alpha = .5)
     67...    g += circle(Omega(t = u), R(t = u), color='blue')
     68sage: # g.show(aspect_ratio=1, xmin=-12, xmax=7, ymin=-3, ymax=12)
     69sage: u, v = var('u, v')
     70sage: h = lambda u,v: u^2 + 2*v^2
     71sage: f = plot3d(h, (u,-1,1), (v,-1,1), aspect_ratio=[1,1,1])
     72sage: f(x, y) = x^2 * y / (x^4 + y^2)
     73sage: t, theta = var('t theta')
     74sage: limit(f(t * cos(theta), t * sin(theta)) / t, t=0)
     75cos(theta)^2/sin(theta)
     76sage: solve(f(x,y) == 1/2, y)
     77[y == x^2]
     78sage: a = var('a'); h = f(x, a*x^2).simplify_rational(); h
     79a/(a^2 + 1)
     80sage: g = plot(h, a, -4, 4)
     81sage: p = plot3d(f(x, y), (x,-2,2), (y,-2,2), plot_points=[150,150])
     82sage: for i in range(1,4):
     83...    p += plot3d(-0.5 + i / 4, (x, -2, 2), (y, -2, 2),\
     84...                 color=hue(i / 10), opacity=.1)
     85sage: x, y, z = var('x, y, z'); a = 1
     86sage: h = lambda x, y, z:(a^2 + x^2 + y^2)^2 - 4*a^2*x^2-z^4
     87sage: f = implicit_plot3d(h, (x, -3, 3), (y, -3, 3), (z, -2, 2),\
     88...                    plot_points=100, adaptative=True)
     89sage: g1 = line3d([(-10*cos(t)-2*cos(5*t)+15*sin(2*t),\
     90...     -15*cos(2*t)+10*sin(t)-2*sin(5*t),\
     91...      10*cos(3*t)) for t in srange(0,6.4,.1)],radius=.5)
     92
     93"""
  • new file sage/tests/french_book/sol_calculus.py

    diff -r 2a2abbcad325 -r 8f9ee6055e9e sage/tests/french_book/sol_calculus.py
    - +  
     1r"""
     2sage: var('n k'); p = 4; s = [n + 1]
     3(n, k)
     4sage: for k in (1..p):
     5...     s = s + [factor((((n + 1)^(k + 1) - sum(binomial(k + 1, j)*s[j]
     6...              for j in (0..k - 1))) / (k + 1)))]
     7sage: s
     8[n + 1, 1/2*(n + 1)*n, 1/6*(n + 1)*(2*n + 1)*n, 1/4*(n + 1)^2*n^2, 1/30*(n + 1)*(2*n + 1)*(3*n^2 + 3*n - 1)*n]
     9sage: var('x h a'); f = function('f', x)
     10(x, h, a)
     11sage: g(x) = taylor(f, x, a, 3)
     12sage: phi(h) = (g(a+3*h)-3*g(a+2*h)+3*g(a+h)-g(a))/h**3
     13sage: phi(h).expand()
     14D[0, 0, 0](f)(a)
     15sage: n = 7; var('x h a'); f = function('f', x)
     16(x, h, a)
     17sage: g(x) = taylor(f, x, a, n)
     18sage: phi(h) = sum(binomial(n,k)*(-1)^(n-k)*g(a+k*h) for k in (0..n))/h**n
     19sage: phi(h).expand()
     20D[0, 0, 0, 0, 0, 0, 0](f)(a)
     21sage: theta = 12 * arctan(1/38) + 20 * arctan(1/57) + 7 * arctan(1/239) + 24 * arctan(1/268)
     22sage: x = tan(theta)
     23sage: y = x.trig_expand()
     24sage: y.trig_simplify()
     251
     26sage: M = 12*(1/38)+20*(1/57)+ 7*(1/239)+24*(1/268)
     27sage: M
     2837735/48039
     29sage: x = var('x')
     30sage: f(x) = taylor(arctan(x), x, 0, 21)
     31sage: approx = 4 * (12 * f(1/38) + 20 * f(1/57) + 7 * f(1/239) + 24 * f(1/268))
     32sage: approx.n(digits = 50); pi.n(digits = 50)
     333.1415926535897932384626433832795028851616168852864
     343.1415926535897932384626433832795028841971693993751
     35sage: approx.n(digits = 50) - pi.n(digits = 50)
     369.6444748591132486785420917537404705292978817080880e-37
     37sage: n = var('n'); phi = lambda x: n*pi +pi/2 - arctan(1/x); x = pi * n
     38sage: for i in range(4): x = taylor(phi(x), n, oo, 2 * i); x
     391/2*pi + pi*n
     401/2*pi + pi*n - 1/(pi*n) + 1/2/(pi*n^2)
     411/2*pi + pi*n - 1/(pi*n) + 1/2/(pi*n^2) - 1/12*(3*pi^2 + 8)/(pi^3*n^3) + 1/8*(pi^2 + 8)/(pi^3*n^4)
     421/2*pi + pi*n - 1/(pi*n) + 1/2/(pi*n^2) - 1/12*(3*pi^2 + 8)/(pi^3*n^3) + 1/8*(pi^2 + 8)/(pi^3*n^4) - 1/240*(15*pi^4 + 240*pi^2 + 208)/(pi^5*n^5) + 1/96*(3*pi^4 + 80*pi^2 + 208)/(pi^5*n^6)
     43sage: f(x, y) = x * y * (x**2 - y**2) / (x**2 + y**2)
     44sage: D1f(x, y) = diff(f(x,y), x)
     45sage: limit((D1f(0,h) - 0) / h, h=0)
     46-1
     47sage: D2f(x, y) = diff(f(x,y), y)
     48sage: limit((D2f(h,0) - 0) / h, h=0)
     491
     50sage: g = plot3d(f(x, y), (x, -3, 3), (y, -3, 3))
     51sage: n, t = var('n, t')
     52sage: v(n) = (4/(8*n+1)-2/(8*n+4)-1/(8*n+5)-1/(8*n+6))*1/16^n
     53sage: assume(8*n+1>0)
     54sage: u(n) = integrate((4*sqrt(2)-8*t^3-4*sqrt(2)*t^4-8*t^5)\
     55...                     * t^(8*n), t, 0, 1/sqrt(2))
     56sage: (u(n)-v(n)).simplify_full()
     570
     58sage: J = integrate((4*sqrt(2)-8*t^3-4*sqrt(2)*t^4-8*t^5)\
     59...                 / (1-t^8), t, 0, 1/sqrt(2))
     60sage: J.simplify_full()
     61pi + 2*log(sqrt(2) - 1) + 2*log(sqrt(2) + 1)
     62sage: ln(exp(J).simplify_log())
     63pi
     64sage: l = sum(v(n) for n in (0..40)); l.n(digits=60); pi.n(digits=60)
     653.14159265358979323846264338327950288419716939937510581474759
     663.14159265358979323846264338327950288419716939937510582097494
     67sage: print "%e" % (l-pi).n(digits=60)
     68-6.227358e-54
     69sage: var('X'); ps = lambda f,g : integral(f * g, X, -pi, pi)
     70X
     71sage: n = 5; Q = sin(X)
     72sage: var('a a0 a1 a2 a3 a4 a5'); a= [a0, a1, a2, a3, a4, a5]
     73(a, a0, a1, a2, a3, a4, a5)
     74sage: P = sum(a[k] * X^k for k in (0..n))
     75sage: equ = [ps(P - Q, X^k) for k in (0..n)]
     76sage: sol = solve(equ, a)
     77sage: P = sum(sol[0][k].rhs() * X^k for k in (0..n))
     78sage: g = plot(P,X,-4,4,color='red') + plot(Q,X,-4,4,color='blue')
     79sage: var('p e theta1 theta2 theta3')
     80(p, e, theta1, theta2, theta3)
     81sage: r(theta) = p / (1-e * cos(theta))
     82sage: r1 = r(theta1); r2 = r(theta2); r3 = r(theta3)
     83sage: R1 = vector([r1 * cos(theta1), r1 * sin(theta1), 0])
     84sage: R2 = vector([r2 * cos(theta2), r2 * sin(theta2), 0])
     85sage: R3 = vector([r3 * cos(theta3), r3 * sin(theta3), 0])
     86sage: D = R1.cross_product(R2) + R2.cross_product(R3) + R3.cross_product(R1)
     87sage: i = vector([1, 0, 0])
     88sage: S = (r1 - r3) * R2 + (r3 - r2) * R1 +   (r2 - r1) * R3
     89sage: V =  S + e * i.cross_product(D)
     90sage: map(lambda x:x.simplify_full(), V) # rep. : [0, 0, 0]
     91[0, 0, 0]
     92sage: map(lambda x:x.simplify_full(), S.cross_product(D))
     93[(e*p^4*sin(theta1)^2*cos(theta2)^2 - 2*e*p^4*sin(theta1)*sin(theta2)*cos(theta1)*cos(theta2) + e*p^4*sin(theta2)^2*cos(theta1)^2 + (e*p^4*cos(theta1)^2 - 2*e*p^4*cos(theta1)*cos(theta2) + e*p^4*cos(theta2)^2)*sin(theta3)^2 + (e*p^4*sin(theta1)^2 - 2*e*p^4*sin(theta1)*sin(theta2) + e*p^4*sin(theta2)^2)*cos(theta3)^2 - 2*(e*p^4*sin(theta1)^2*cos(theta2) + e*p^4*sin(theta2)^2*cos(theta1) - (e*p^4*sin(theta1)*cos(theta1) + e*p^4*sin(theta1)*cos(theta2))*sin(theta2))*cos(theta3) + 2*(e*p^4*sin(theta1)*cos(theta1)*cos(theta2) - e*p^4*sin(theta1)*cos(theta2)^2 - (e*p^4*cos(theta1)^2 - e*p^4*cos(theta1)*cos(theta2))*sin(theta2) - (e*p^4*sin(theta1)*cos(theta1) - e*p^4*sin(theta1)*cos(theta2) - (e*p^4*cos(theta1) - e*p^4*cos(theta2))*sin(theta2))*cos(theta3))*sin(theta3))/(e^2*cos(theta1)^2 + (e^4*cos(theta1)^2 - 2*e^3*cos(theta1) + e^2)*cos(theta2)^2 + (e^4*cos(theta1)^2 - 2*e^3*cos(theta1) + (e^6*cos(theta1)^2 - 2*e^5*cos(theta1) + e^4)*cos(theta2)^2 - 2*(e^5*cos(theta1)^2 - 2*e^4*cos(theta1) + e^3)*cos(theta2) + e^2)*cos(theta3)^2 - 2*(e^3*cos(theta1)^2 - 2*e^2*cos(theta1) + e)*cos(theta2) - 2*(e^3*cos(theta1)^2 + (e^5*cos(theta1)^2 - 2*e^4*cos(theta1) + e^3)*cos(theta2)^2 - 2*e^2*cos(theta1) - 2*(e^4*cos(theta1)^2 - 2*e^3*cos(theta1) + e^2)*cos(theta2) + e)*cos(theta3) - 2*e*cos(theta1) + 1), 0, 0]
     94sage: N = r3 * R1.cross_product(R2) + r1 * R2.cross_product(R3) + r2 * R3.cross_product(R1)
     95sage: W =  p * S + e * i.cross_product(N)
     96sage: print map(lambda x:x.simplify_full(), W)  # rep. : [0, 0, 0]
     97[0, 0, 0]
     98sage: R1=vector([0,1.,0]);R2=vector([2.,2.,0]);R3=vector([3.5,0,0])
     99sage: r1 = R1.norm(); r2 = R2.norm(); r3 = R3.norm()
     100sage: D = R1.cross_product(R2) + R2.cross_product(R3) + R3.cross_product(R1)
     101sage: S = (r1 - r3) * R2 + (r3 - r2) * R1 + (r2 - r1) * R3
     102sage: V =  S + e * i.cross_product(D)
     103sage: N = r3 * R1.cross_product(R2) + r1 * R2.cross_product(R3) \
     104...  + r2 * R3.cross_product(R1)
     105sage: W =  p * S + e * i.cross_product(N)
     106sage: e = S.norm() / D.norm()
     107sage: p = N.norm() / D.norm()
     108sage: a = p/(1-e^2)
     109sage: c = a * e
     110sage: b = sqrt(a^2 - c^2)
     111sage: X = S.cross_product(D)
     112sage: i = X / X.norm()
     113sage: phi = atan2(i[1],i[0]) * 180 / pi.n()
     114sage: print "%.3f %.3f %.3f %.3f %.3f %.3f" % (a, b, c, e, p, phi)
     1152.360 1.326 1.952 0.827 0.746 17.917
     116sage: A = matrix(QQ, [[2, -3, 2, -12, 33],
     117...                   [ 6, 1, 26, -16, 69],
     118...                   [10, -29, -18, -53, 32],
     119...                   [2, 0, 8, -18, 84]])
     120sage: A.right_kernel()
     121Vector space of degree 5 and dimension 2 over Rational Field
     122Basis matrix:
     123[     1      0  -7/34   5/17   1/17]
     124[     0      1  -3/34 -10/17  -2/17]
     125sage: H = A.echelon_form()
     126sage: A.column_space()
     127Vector space of degree 4 and dimension 3 over Rational Field
     128Basis matrix:
     129[       1        0        0 1139/350]
     130[       0        1        0    -9/50]
     131[       0        0        1   -12/35]
     132sage: S.<x,y,z,t>=QQ[]
     133sage: C = matrix(S, 4,1,[x,y,z,t])
     134sage: B = block_matrix([A,C], ncols=2)
     135sage: C = B.echelon_form()
     136sage: C[3,5]*350
     137-1139*x + 63*y + 120*z + 350*t
     138sage: K = A.kernel(); K
     139Vector space of degree 4 and dimension 1 over Rational Field
     140Basis matrix:
     141[        1  -63/1139 -120/1139 -350/1139]
     142sage: matrix(K.0).right_kernel()
     143Vector space of degree 4 and dimension 3 over Rational Field
     144Basis matrix:
     145[       1        0        0 1139/350]
     146[       0        1        0    -9/50]
     147[       0        0        1   -12/35]
     148sage: A = matrix(QQ, [[-2, 1, 1], [8, 1, -5], [4, 3, -3]])
     149sage: C = matrix(QQ, [[1, 2, -1], [2, -1, -1], [-5, 0, 3]])
     150sage: B = C.solve_left(A); B
     151[ 0 -1  0]
     152[ 2  3  0]
     153[ 2  1  0]
     154sage: C.left_kernel()
     155Vector space of degree 3 and dimension 1 over Rational Field
     156Basis matrix:
     157[1 2 1]
     158sage: var('x y z'); v = matrix([[1, 2, 1]])
     159(x, y, z)
     160sage: B = B+(x*v).stack(y*v).stack(z*v); B
     161[      x 2*x - 1       x]
     162[  y + 2 2*y + 3       y]
     163[  z + 2 2*z + 1       z]
     164sage: A == B*C
     165True
     166
     167"""
  • new file sage/tests/french_book/sol_graphiques_1.py

    diff -r 2a2abbcad325 -r 8f9ee6055e9e sage/tests/french_book/sol_graphiques_1.py
    - +  
     1r"""
     2sage: t = var('t'); liste = [a + cos(t) for a in srange(0, 10, 0.1)]
     3sage: g = polar_plot(liste, (t, 0, 2 * pi)); g.show(aspect_ratio=1)
     4sage: f = lambda x: abs(x**2 - 1/4)
     5sage: def liste_pts(u0, n):
     6...     u = u0; liste = [[u0,0]]
     7...     for k in range(n):
     8...         v, u = u, f(u)
     9...         liste.extend([[v,u], [u,u]])
     10...     return(liste)
     11sage: g = line(liste_pts(1.1, 8), rgbcolor = (.9,0,0))
     12sage: g += line(liste_pts(-.4, 8), rgbcolor = (.01,0,0))
     13sage: g += line(liste_pts(1.3, 3), rgbcolor = (.5,0,0))
     14sage: g += plot(f, -1, 3, rgbcolor = 'blue')
     15sage: g += plot(x, -1, 3, rgbcolor = 'green')
     16sage: # g.show(aspect_ratio=1, ymin = -.2, ymax = 3)
     17sage: x = var('x'); y = function('y',x); DE = x^2 * diff(y, x) -  y == 0
     18sage: desolve(DE, [y,x])
     19c*e^(-1/x)
     20sage: g = plot([c*e^(-1/x) for c in srange(-8, 8, 0.4)], (x, -3, 3))
     21sage: y = var('y'); g += plot_vector_field((x^2, y), (x,-3,3), (y,-5,5))
     22sage: # g.show(ymin = -5, ymax = 5)
     23sage: from sage.calculus.desolvers import desolve_system_rk4
     24sage: f = lambda x, y:[a*x-b*x*y,-c*y+d*b*x*y]
     25sage: x,y,t = var('x y t')
     26sage: a, b, c, d = 1., 0.1, 1.5, 0.75
     27sage: P = desolve_system_rk4(f(x,y),[x,y],
     28...        ics=[0,10,5],ivar=t,end_points=15)
     29sage: Ql = [ [i,j] for i,j,k in P]
     30sage: p = line(Ql,color='red')
     31sage: p += text("Lapins",(12,37),fontsize=10,color='red')
     32sage: Qr = [ [i,k] for i,j,k in P]
     33sage: p += line(Qr,color='blue')
     34sage: p += text("Renards",(12,7),fontsize=10,color='blue')
     35sage: p.axes_labels(["temps","population"])
     36sage: # p.show(gridlines=True)
     37sage: ### Deuxieme graphique
     38sage: n = 10;  L = srange(6, 18, 12 / n); R = srange(3, 9, 6 / n)
     39sage: def g(x,y):
     40...     v = vector(f(x, y))
     41...     return v/v.norm()
     42...
     43sage: q = plot_vector_field(g(x, y), (x, 0, 60), (y, 0, 36))
     44sage: for j in range(n):
     45...     P = desolve_system_rk4(f(x,y),[x,y],
     46...            ics=[0,L[j],R[j]],ivar=t,end_points=15)
     47...     Q = [ [j,k] for i,j,k in P]
     48...     q += line(Q, color=hue(.8-j/(2*n)))
     49...
     50sage: q.axes_labels(["nombre de lapins","nombre de renards"])
     51sage: # q.show()
     52sage: import scipy; from scipy import integrate
     53sage: def dX_dt(X, t=0):
     54...     return [X[1], 0.5*X[1] - X[0] - X[1]**3]
     55...
     56sage: t = srange(0, 40, 0.01);  x0 = srange(-2, 2, 0.1)
     57sage: y0 = 2.5
     58sage: CI = [[i, y0] for i in x0] + [[i, -y0] for i in x0]
     59sage: def g(x,y):
     60...     v = vector(dX_dt([x, y]))
     61...     return v/v.norm()
     62...
     63sage: x, y = var('x, y')
     64sage: q = plot_vector_field(g(x, y), (x, -3, 3), (y, -y0, y0))
     65sage: for j in xrange(len(CI)):
     66...     X = integrate.odeint(dX_dt, CI[j], t)
     67...     q += line(X, color=(1.7*j/(4*n),1.5*j/(4*n),1-3*j/(8*n)))
     68...
     69sage: X = integrate.odeint(dX_dt, [0.01,0], t)
     70sage: q += line(X, color = 'red')
     71sage: # q.show()
     72
     73"""
  • new file sage/tests/french_book/sol_graphiques_2.py

    diff -r 2a2abbcad325 -r 8f9ee6055e9e sage/tests/french_book/sol_graphiques_2.py
    - +  
     1r"""
     2sage: import scipy
     3sage: from scipy import integrate
     4sage: # Intervalle de temps :
     5sage: t = srange(0, 40, 0.01)
     6sage: # Conditions initiales en cartesiennes :
     7sage: n = 35
     8sage: CI_cart = [[4, .2 * i] for i in range(n)]
     9sage: # Conditions initiales en polaires :
     10sage: CI = map(lambda x:[sqrt(x[0]**2+x[1]**2),\
     11...        pi - arctan(x[1] / x[0])], CI_cart)
     12sage: # Choix du parametre alpha :
     13sage: alpha = [0.1, 0.5, 1, 1.25]
     14sage: for i in range(len(alpha)):
     15...     # Definition du syteme differentiel :
     16...     def dX_dt(X, t=0):
     17...         return [cos(X[1]) * (1 - 1 / X[0]^2),\
     18...         -sin(X[1]) * (1 / X[0] + 1 / X[0]^3) \
     19...         + 2 * alpha[i] / X[0]^2]
     20...     # trace du disque :
     21...     q = circle((0, 0), 1, fill=True, rgbcolor='purple')
     22...     for j in range(n):
     23...         # resolution du syteme autonome :
     24...         X = integrate.odeint(dX_dt, CI[j], t)
     25...         # passage en cartesiennes :
     26...         Y = [[u[0] * cos(u[1]), u[0] * sin(u[1])] for u in X]
     27...         # trace stocke dans la variable q
     28...         q += line(Y, xmin = -4, xmax = 4, color='blue')
     29...     # Trace final :
     30...     # q.show(aspect_ratio = 1, axes = False)
     31
     32"""