Carl Love

Carl Love

28085 Reputation

25 Badges

13 years, 98 days
Himself
Wayland, Massachusetts, United States
My name was formerly Carl Devore.

MaplePrimes Activity


These are answers submitted by Carl Love

Just change your command from x=0 to h=0:

taylor(y(x+h), h=0);

To have matrix and vector computations done in single precision, create the Matrices and Vectors with the option datatype= float[4].

Example:

V:= Vector(1000, n-> 1/n^3, datatype= float[4]);
V.V;

                     1.01734306201200031

Compare with the result obtained using datatype= float[8], which is double precision.

                     1.01734306198444679

Note that the results begin to differ at the 9th decimal place.

You're right. I don't think that there's any way to get it to accept infinity. Here's a workaround: Use gamma instead of infinity. When the tutor is done, give the command:

limit(subs(gamma= x, %), x= infinity);

For example,

count:= 0:
for i from 1 to 9 do
     for j from 1 to i do
          count:= count+1
     end do
end do;
count;
                                              45

Of course, this example is trivial, and there are much more efficient ways to perform this particular computation.

interface(showassumed= 0):
assume(nu::real, m>0, k>0, T>0):
Int(nu*4*Pi*(m/2/Pi/k/T)^(3/2)*nu^2*exp(-m*nu^2/2/k/T), nu= 0..infinity);

value(%);

Whenever a symbol is redefined by a package, the original value of the symbol can be accessed by prefixing the symbol with :-. So in this case you'd use :-``.

with(Physics):
``(3)*``(2);
Error, (in GetDifferentiationVariables) differentiation variables for dAlembertian, d_ and D_ are not defined; either use Setup(differentiationvariables =  ...) to set them, or indicate the differentiation variables as a second argument to dAlembertian

:-``(3)*:-``(2);
                            (3) (2)

There is a policy stated in the 6th paragraph at ?name : 

Names that start with an underscore (with the exception of _Env) are used as global variable names by Maple and are effectively reserved for use by library code. Do not use a symbol beginning with an underscore. The empty name `` must never be assigned a value.

It seems that Physics violates this.

 

For the constraints to DirectSearch, you could generate all pairs of the variables like this:

map(`<>`@op, combinat:-choose({a,b,c,d}, 2));

You can use VolumeOfRevolution and exchange the roles of x and y:

Student:-Calculus1:-VolumeOfRevolution(
      2*t^3-t^4, t= 0..2,
      axis= vertical, distancefromaxis = -2,
      output= plot, labels= [y,x,z], caption= ``
);

I believe that this is the same surface as Preben plotted.

Is that what you were expecting?

You need to remove the definition of POWER.

restart:
ee:=u^2+v^2-w^2:
applyrule(a::algebraic^b::nonunit(algebraic) = POWER(a,b), ee);

            POWER(u, 2) + POWER(v, 2) - POWER(w, 2)

Note also the presence of the modifier nonunit. This is necessary to prevent an infinite loop (where it treats everything as itself to the power 1).

Another way:

subsindets(ee, `^`, x-> POWER(op(x)));


(**)

f:= x-> ((x^2-x-3)/(x^2+1))*(x^2+x+1);

proc (x) options operator, arrow; (x^2-x-3)*(x^2+x+1)/(x^2+1) end proc

(**)

D(f)(x);

(2*x-1)*(x^2+x+1)/(x^2+1)-2*(x^2-x-3)*(x^2+x+1)*x/(x^2+1)^2+(x^2-x-3)*(2*x+1)/(x^2+1)

(**)

simplify(%);

2*(x^5+2*x^3+2*x^2-2)/(x^2+1)^2

(**)

 


Download diff.mw

You simply need to use the two-argument eval command.


(**)

eqns:=[
    S     = (1-alpha)*((1-beta)*A*g*k^(g-1)-beta*(k+b)*A*(g-1)*g*k^(g-2)),
    c_3   = (1-alpha)*A*(1-g)*k^g,
    c_4   = (1-alpha)*(1-beta)*A*g*k^(g-1),
    c_5   = -(1-alpha)*b*(A*g*k^(g-1)-1),
    cv    = c_3*sigma/(sqrt(1+c_4)*sqrt(c_3+c_5)),
    phi_3 = sigma/(sqrt(1+c_4)*sqrt(c_3+c_5))-(1/2)*cv/(c_3+c_5),
    c_3k  = -(1-alpha)*k*A*(g-1)*g*k^(g-2),
    phi_4 = -(1/2)*cv/(1+c_4),
    c_4k  = (1-alpha)*(1-beta)*A*(g-1)*g*k^(g-2),
    phi_5 = -(1/2)*cv/(c_3+c_5),
    c_5k  = -(1-alpha)*b*A*(g-1)*g*k^(g-2),
    c_5b  = (1-alpha)*(1-A*g*k^(g-1)),
    G     = c_3k*phi_3+c_4k*phi_4+c_5k*phi_5,
    B     = -(1+(1-alpha)*(beta*A*g*k^(g-1)-1))/
            (1-(1-alpha)*((1-beta)*A*g*k^(g-1)-beta*(k+b)*A*(g-1)*g*k^(g-2)))
]:

(**)

vals:= [alpha= 0.3, beta= 0.3, g= 0.3, sigma=1, b= 0.85, k= 0.83]:

(**)

eval(eqns, vals);

[S = .269178523926482*A, c_3 = .463361032316227*A, c_4 = .167479891198636*A, c_5 = -.203368439312629*A+.595, cv = c_3/((1+c_4)^(1/2)*(c_3+c_5)^(1/2)), phi_3 = 1/((1+c_4)^(1/2)*(c_3+c_5)^(1/2))-(1/2)*cv/(c_3+c_5), c_3k = .167479891198637*A, phi_4 = -(1/2)*cv/(1+c_4), c_4k = -.141248101010898*A, phi_5 = -(1/2)*cv/(c_3+c_5), c_5k = .171515551227519*A, c_5b = .7-.239256987426622*A, G = c_3k*phi_3+c_4k*phi_4+c_5k*phi_5, B = -(.3+0.717770962279868e-1*A)/(1-.269178523926482*A)]

(**)

 


Download eval.mw

Assuming that you mean p is the conjugate of z, then

restart:
eq:= (2*z-3)*(-3*p+4) - (2*p-3)*(-3*z+4)=0:
solve(eval(eq, [z= x + I*y, p= x - I*y]));

Therefore z can be any real number.

Several pointers in general:

  1. "e" means nothing on input to Maple.The correct name is exp. You had e in one place and exp in another.
  2. Functions should be defined as f:= (x,y)-> ..., not as f(x,y):= ....
  3. If you want numeric integration, use Int, not int. The latter will try symbolic integration first, potentially taking a long time.

Specifically for this problem: I got 6 digits of precision on the integral in mere seconds by passing the optional argument digits= 6 to Int. There seems to be a threshold effect: At 7 digits, it starts to take a very long time. Since you have constants in your problem with only 4 digits, 6 should be enough.

Here's your worksheet with my corrections:

(**)

restart:

(**)

width:= 35:  gap:= 8.89:  B0:= 1.30615:  centralR:= 1.5322550941477154982:

(**)

 

(**)

z:= (x,y)-> sqrt((-0.012029242-0.05060297686 *x+0.2191855734* y)^2+(0.0521043725+0.2191855734* x-0.9493970231* y)^2);

(**)

 

(**)

perpDegradation:= (x,y)-> B0/(1+exp(-6.762031311+17.48081911 *z(x,y)+2.232089510* Pi-0.8171* (7.320276006* z(x,y)+0.9347108501 *Pi-2.936480740)^2+0.2 *(7.320276006 *z(x,y)+0.9347108501* Pi-2.936480740)^3));

(**)

 

(**)

extendedField:= (x,y)-> perpDegradation(x,y)*(Heaviside(-4.331475875*x+7.474536933-y)-Heaviside(-4.331475875*x+6.140092706-y))*Heaviside(-0.2308681911*x-0.05488154218+y);

(**)

 

(**)

By:= (r,theta)-> B0*(Heaviside(r-1.3822551)-Heaviside(r-1.6822551))*Heaviside(.2308681911*r*cos(theta)+0.5488154218e-1-r*sin(theta))/(1+exp(.25028469571727139024+2.388*(11.21653020*theta-2.936480740)+((-1)*.8171)*(11.21653020*theta-2.936480740)^2+.2*(11.21653020*theta-2.936480740)^3));

(**)

 

(**)

Byextended:= (r,theta)-> By(r,theta)+extendedField(r*cos(theta),r*sin(theta));

(**)

 

35

8.89

1.30615

1.5322550941477154982

proc (x, y) options operator, arrow; sqrt((-0.12029242e-1-0.5060297686e-1*x+.2191855734*y)^2+(0.521043725e-1+.2191855734*x-.9493970231*y)^2) end proc

proc (x, y) options operator, arrow; B0/(1+exp(-6.762031311+17.48081911*z(x, y)+2.232089510*Pi-.8171*((7.320276006*z(x, y)+.9347108501*Pi-2.936480740)^2)+.2*((7.320276006*z(x, y)+.9347108501*Pi-2.936480740)^3))) end proc

proc (x, y) options operator, arrow; perpDegradation(x, y)*(Heaviside(-4.331475875*x+7.474536933-y)-Heaviside(-4.331475875*x+6.140092706-y))*Heaviside(-.2308681911*x-0.5488154218e-1+y) end proc

proc (r, theta) options operator, arrow; B0*(Heaviside(r-1.3822551)-Heaviside(r-1.6822551))*Heaviside(.2308681911*r*cos(theta)+0.5488154218e-1-r*sin(theta))/(1+exp(.25028469571727139024+2.388*(11.21653020*theta-2.936480740)+(-1)*.8171*((11.21653020*theta-2.936480740)^2)+.2*((11.21653020*theta-2.936480740)^3))) end proc

proc (r, theta) options operator, arrow; By(r, theta)+extendedField(r*cos(theta), r*sin(theta)) end proc

(**)

evalf(Int(Byextended(centralR,theta), theta=0..Pi/6, digits= 6));

.348526

(**)

 


Download 6_digits.mw

 

It's tricky because Maple automatically distributes integer factors over sums. You can do this:

(**)

p:= 18*R*r - 3*s^2 + 9*r^2;

18*R*r+9*r^2-3*s^2

(**)

content(p)*``(primpart(p));

3*``(6*R*r+3*r^2-s^2)

(**)

 

Notice the special quote marks in my command. It's two single backquotes (the key under ESC on my keyboard) together with nothing in between them.

The expand command will put this back to its original form.

Download primpart.mw

 

By starting from the closed-form expression for the Fibonacci numbers, we can generate a one-line algebraic expression for "next Fibonacci".

 

A procedure for the "next Fibonacci" modeled on nextprime. Written by Carl Love 18 September 2013.

(**)

restart:

(**)

Fn:= rsolve({F(n+1) = F(n) + F(n-1), F(1)=1, F(2)=1}, F(n));

-(1/5)*5^(1/2)*(-(1/2)*5^(1/2)+1/2)^n+(1/5)*5^(1/2)*((1/2)*5^(1/2)+1/2)^n

(**)

evalf(%);

-.447213595499958*(-.618033988749900)^n+.447213595499958*1.61803398874990^n

Note that the negative term gets smaller in magnitude rapidly (exponentially) with increasing n. Thus the value of the Fibonacci number can be determined by rounding the positive term to the nearest integer. I manually extract the positive term for this purpose:

(**)

FibOfIndex:= n-> ((1+sqrt(5))/2)^n/sqrt(5);

proc (n) options operator, arrow; (1/2+(1/2)*sqrt(5))^n/sqrt(5) end proc

Note that the above is continuous for any nonegative n. Below is simply the inverse of the above function. It returns the n for a given Fibonacci F. It also is continuous for positive F. That means that we can apply it to the numbers between the Fibonacci numbers!

(**)

IndexOfFib:= unapply(rhs(isolate(F=FibOfIndex(n), n)), F);

proc (F) options operator, arrow; ln(F*5^(1/2))/ln((1/2)*5^(1/2)+1/2) end proc

Now we see that we can make a single algebraic expression for the "next Fibonacci".

(**)

nextFibo:= proc(x)
     option remember;
     local r;
     r:= trunc(

          evalf[3+ilog10(x)](round(FibOfIndex(ceil(IndexOfFib(x)))))

     );
     `if`(r=x, thisproc(x+1), r)
end proc:

Hard code the "difficult" cases.

(**)

(nextFibo(0),nextFibo(1)):= (1,2):

 

Download next_fibonacci.mw

First 347 348 349 350 351 352 353 Last Page 349 of 395