Joe Riel

9660 Reputation

23 Badges

20 years, 4 days

MaplePrimes Activity


These are answers submitted by Joe Riel

You presumably don't want B(...), but instead B*(...). With that an another typo fixed, the following works

ode := diff(f(x),x,x,x)+f(x)*diff(f(x),x,x)+B*(1-diff(f(x),x))^2=0;
bcs:= f(0)=0, D(f)(0)=0, D(f)(5)=1;
sys := {bcs, ode};

for b in [-.199,-.18,0,.3] do
    sys0 := eval(sys, B=b);
    dsn := dsolve(sys0, numeric ):
    print(plots:-odeplot(dsn, [x,f(x)], 0..5));
end do:

The first problem is that your points are entered incorrectly, you close the list after the first pair (note the ]]). Change that to a single ].

The next problem is that polynomial interpolation requires distinct independent values; here the x-coordinate 3.43799 is repeated. 

You can sum the terms inside the do loop, say

rtot := 0:
for z in w do
    tot := select(has,z,phi);
    r := 2*tot;
    rtot := rtot + r;
end do:

Simpler is to just do

  2*select(has,w,phi);
select(has, [t], beta);

I don't understand what you want.  You aren't taking alternate terms, that would be

[ [1,1,1,1,1], [], [1,1,1,1,1], [], [1,1], [], ??? ]

One way is to assign the functions to return the cumulative distance. For example,

tom := n -> add(1/k, k=1..n):
jerry := n -> add(2^(2-k), k=1..n):

I used ?add rather than ?sum because it is intended for definite summation, that is, when the limits are known.

Now use ?plots[pointplot] to generate a discrete plot for each snail.

ptom := plots:-pointplot([seq([n,tom(n)], n=1..10)], color=red):
pjer := plots:-pointplot([seq([n,jerry(n)], n=1..10)], color=blue):

Finally, call ?plots[display] to display both plots on the same graph

plots:-display(ptom,pjer);

Connect the sensor between one frame and the prismatic joint.

Look at the ?binomial function. 

Edit: I misinterpreted your question.  Dang, I cannot see the question when editing, so now again am in the dark. Will have to reedit.

After

convert( . , Sum, x);

do

convert(%, binomial);

Unless you explicity assigned mf as an Array and are using parentheses as indices, the assignment you are doing is probably not what you want. Regardless, it isn't clear what you want.  If you just want to create a list of items corresponding to the rhs of your equation, then do

L := [seq(seq(f(Prc)*m(Prc,B), B = 2..200,0.5), Prc = 2..100)];
sort(L);

 

Generally, for digit comparison, it is useful to convert to strings and compare characters.  Here the ?StringTools[CommonPrefix] procedure should be useful.

I believe my hypothesis can be fixed by constraining the applicable functions.  Here's the basic idea

restart;
with(IntegrationTools):

F := 1/deltax*Int(sin(f(x)),x=x1..x1+deltax);
Ftheta1 := Change(F, f(x)=theta, [theta]);
Ftheta2 := eval(Ftheta1, f(x1+deltax)=convert(series(f(x1+deltax),deltax,2),'polynom'));
Ftheta3 := eval(Ftheta2, [f(x1)=2*n*Pi, deltax = Pi/D(f)(x1)]);
Fphi1 := Change(Ftheta3, phi=theta-2*n*Pi, [phi]);
Fphi2 := eval(Fphi1, sin(phi+2*n*Pi) = sin(phi));
Fphi3 := evalindets(Fphi2
                    , specfunc(specfunc(anything,RootOf),D(f))
                    , DF -> eval(series( DF, phi, 2), RootOf(-f(_Z) + 2*n*Pi) = x1)
                   );
                   );
                           Pi
                          /
               D(f)(x1)  |                   sin(phi)
      Fphi3 := --------  |    -------------------------------------- dphi
                  Pi     |                 (2)
                         |               (D   )(f)(x1)            2
                        /     D(f)(x1) + ------------- phi + O(phi )
                          0                D(f)(x1)

This suggests that if f''(x)/f'(x)^2 approaches 0 as x goes to infinity, then the limit goes to 2/Pi.  That explains why ln(x) fails.  However, both x^2 and x*exp(x) meet the critieria.  One can directly compute the limit for f = ln(x):

restart;
with(IntegrationTools):
f := ln:
F := 1/deltax*Int(sin(f(x)),x=x1..x1+deltax);
Ftheta := Change(F, f(x)=theta, [theta]);
sol1 := solve(f(x1)=2*n*Pi,[x1]);
sol2 := solve(f(x1+deltax)=(2*n+1)*Pi,[deltax]);
y1 := simplify(eval(eval(Ftheta,sol2[1]),sol1[1]),ln) assuming real;
y2 := Change(y1, theta0=theta-2*n*Pi, [theta0]);
y3 := eval(y2, sin(theta0+2*n*Pi)=sin(theta0));
y4 := value(y3);
y5 := simplify(expand(y4));
                                        1 + exp(Pi)
                             y5 := 1/2 -----------
                                       exp(Pi) - 1
               

Seems like you might be able to handle the digits by using logarithms.  Euclid would then be changed to compute the sum of the logs.  Regardless, this problem quickly becomes intractable.

To determine where the slow parts are, create a procedure and then profile it.  You can use the ?CodeTools[Profiling] package. 

restart;
euclid := proc (s, n)
    if 1 < s and 1 <= n then
        mul(1/(1-ithprime(i)^(-s)), i = 1 .. n);  # changed product to mul; it is faster
    end if
end proc:

Print := proc()
local q,m,s,n,b;
    for q to 5 do       # reduced to 5 to speed up the run
        m := 10^(-q);
        for s from 2 to 10 do
            n := 1;
            b := 1;
            while b >= m do
                b := evalf(abs(Zeta(s)-euclid(s, n)));
                n := n+1
            end do;
            #print(m, s, n-1) # remove output for now
        end do
    end do
end proc:

with(CodeTools:-Profiling):

t := Build('procs'=[euclid, ithprime, Zeta, Print], 'commands'='Print()'):
PrintProfiles(t);
SortBy(t, order=time);

function                 calls      time     time%          words    words%
---------------------------------------------------------------------------
euclid                    2288     8.547     57.90      166746250     85.84
ithprime               1686699     5.386     36.49       13812972      7.11
Print                        1     0.812      5.50       13643772      7.02
Zeta                       129     0.016      0.11          40111      0.02
---------------------------------------------------------------------------
total:                 1689117    14.761    100.00      194243105    100.00

Clearly euclid is the problem.  There is a simple way to significantly improve its performance as used; create a recursive  version that caches its results.  That works because euclid is repeatedly called with increasing values of n, so we can reuse the previous result.

euclid := proc(s,n)
option cache;
    if n = 1 then
        1/(1-1/2^s);
    else
        procname(s,n-1)/(1-1/ithprime(n)^s);
    end if;
end proc:

Reprofiling gives

SortBy(t, order=time);
function                 calls      time     time%          words    words%
---------------------------------------------------------------------------
Print                        1     0.828     88.09       13628634     79.07
ithprime                  5941     0.056      5.96         634754      3.68
euclid                    1963     0.044      4.68        2935657     17.03
Zeta                       129     0.012      1.28          37970      0.22
---------------------------------------------------------------------------
total:                    8034     0.940    100.00       17237015    100.00

I expect both to be 2/Pi. 

F(x) := whatever;

is not equivalent to

F := x -> whatever;

Both are legal, but the first is rarely what you want. It doesn't create a procedure, rather it merely assigns whatever to the function call F(x), which has almost nothing to do with, say, F(3).  As an example, consider

F := x -> x:    # assign a procedure
F(3) := 100: # assign an entry in F's remember table

We can now evaluate

F(z);
            z
F(3);
            100

What has actually happened is that the procedure F is assigned but also has a remember table that bypasses the computation when the argument is precisely 3.  Note that the argument must match exactly for the entry in the remember table to be accessed. For example,

F(3.) := 24:
F(3), F(3.), F(3.0);
              100, 24, 3.0

For your example, what you probably want is

F := unapply(int(x^2, x = 0..x), x);

Clearer might be to use separate variables for the independent variable of the integral and its upper limit:

F := unapply(int(x^2, x=0..y), y);
First 56 57 58 59 60 61 62 Last Page 58 of 114