Carl Love

Carl Love

28020 Reputation

25 Badges

12 years, 300 days
Himself
Wayland, Massachusetts, United States
My name was formerly Carl Devore.

MaplePrimes Activity


These are answers submitted by Carl Love

Unfortuately, the equals sign (=) cannot be used to test equality of Vectors and Matrices. You need to use LinearAlgebra:-Equal instead. So change (s7 = e) to LinearAlgebra:-Equal(s7, e), and likewise for the rest.

You can use ormap to condense your compound if...elif...end if to a single condition:

if ormap(LinearAlgebra:-Equal, [e,f,g,h], s7) then
     print("Synchronising")
else
     print("Not Synchronising")
end if;

By using `if` and cat this can be reduced to a single line:

print(cat(`if`(ormap(LinearAlgebra:-Equal, [e,f,g,h], s7), "", "Not "),"Synchronising"));

Problem 1 must be a trick question or you transcibed something wrong. There is no square root of 2 modulo 3 (2 is not a quadratic residue modulo 3), so there are also none modulo 3^8. The p-adic Newton-Hensel formula only works if you start with an exact solution modulo the prime base. Roots will confirm that there is no solution. Roots will also give you the solution to problem 2, but clearly the problem was intended to be done with 5-adic lifting (625 = 5^(2^2)).

Roots(x^2-2) mod 3^8;
                               []
Roots(x^3-2) mod 625;
                           [[303, 1]]

Go to the Maple Applications Center and get the paper and worksheet that I wrote on p-adic Newton-Hensel lifting.

Here is the answer to problem 2 via 5-adic lifting:

f:= x-> x^3-2:
Roots(f(x)) mod 5;
                            [[3, 1]]
3 - f(3)/D(f)(3) mod 5^(2^1);
                               3
3 - f(3)/D(f)(3) mod 5^(2^2);
                              303

Is there a chance that you've made an assumption twice or that you've made a second, different assumption on the same variable? Both can have that undesirable effect. Observe:

restart:
assume(m>0); expr:= m^2; assume(m>0);

subs(m=1, expr);

The same thing would happen if I made the second assumption different from the first, as long as it still involved m.

Addressing your first problem, the primitive 4th roots of unity modulo 29: The inert mod command Roots will compute all the roots of a polynomial. We apply it to the polynomial x^p - 1. Then we select the roots whose order modulo 29 is p.

PrimitiveRootsOfUnity:= proc(p::posint, r::And(posint, Not(identical(1))))
local x;
     select(z-> numtheory:-order(z,r)=p, map2(op, 1, Roots(x^p-1) mod r))[]
end proc:

PrimitiveRootsOfUnity(4,29);
                            12, 17

Since there are only two, it is obvious that one is the inverse of the other.
1/12 mod 29;
                               17

It is very easy to convert any CURVES in a plot into POLYGONS. Afterall, a CURVES structure is just a finite collection of points. Addressing your problem, plotting a segment of circle: Any arc of circle determines a segment of circle, right? The following procedure uses the same arguments as plottools:-arc and converts the structure to a polygon.

CircleSegment:= (C::[realcons, realcons], R::positive, AB::range(realcons))->
     plots:-polygonplot(op(plottools:-arc(C, R, AB)), scaling= constrained, _rest):

Example of use (center = [0,0], radius= 1, from 0 to 1 radians):

CircleSegment([0,0], 1, 0..1, color= green);

restart:

Define the general Newton-Raphson operator. This is an operator-valued operator which returns the Newton iteration operator for its input operator.
NR:= f-> (z->z) - f/D(f):

Apply it to our problem:
N:= NR(z-> z^2+1):

Verify that this gives the expected N:
simplify(N(z));

Note that the imaginary unit in Maple is capital I.
T:= z-> (z-I)/(z+I):

Show that T(N(z)) = T(z)^2. (Eq. 1)
T(N(z)) - T(z)^2;

simplify(%);
                               0

Verify that T(N^2(z)) = T(z)^4. Note that here N^2(z) = N(N(z)). In Maple this can also be expressed as (N@@2)(z). You must not enter this as N^2(z). Rather than getting an error message, you'll just get something wrong.
simplify(T((N@@2)(z)) - T(z)^4);
                               0

Q1) What will be the general result?  
Ans: I guess that this means to guess at the answer to Q2 based on the above.

Q2) T(N^k(z)) is what power of T(z)?
Ans: T((N@@k)(z)) = T(N((N@@(k-1))(z))) = T((N@@(k-1))(z))^2, by (Eq. 1) with (N@@(k-1))(z) substituted for z. It follows by induction that T((N@@k)(z)) = T(z)^(2^k). So the answer is 2^k.

Q3) What is T(N^3(z))?
Ans: By the answer to Q2, it is T(z)^8. Verify that:
simplify(T((N@@3)(z)) - T(z)^8);
                               0

The desired upper bound of error should be passed to the procedure in a variable traditionally named epsilon. If no value is passed, then try using the value 10.^(1-Digits) for relative error. This can often be tricky to achieve, so you should also exit your loop after a certain preset maximum number of iterations. The environment variable Digits is preset in Maple.

Here is a simple, yet robust, square root finder that uses all of these principles:

SQRT:= proc(x::positive, {maxiters::posint:= 99}, {epsilon::positive:= 10.^(1-Digits)})
local 
     z:= evalf(x),
     N:= z-> (z + x/z)/2, #Newton-Raphson iteration for sqrt.
     newz,
     err:= epsilon,
     k
;
     for k to maxiters while err >= epsilon do
          newz:= N(z);
          err:= abs((z-newz)/newz); #relative error
          z:= newz
     od;
     if k > maxiters then WARNING("Did not converge.") end if;
     z
end proc;

SQRT(2);
                       
1.41421356237310

Please let me know if that gives you enough information to solve your problem. If not, then you'll need to be more specific about the problem.

If you change 0.6 to 3/5, then you will get the leading logarithmic term that you expect. That being said, I don't know why some logarithmic term is not there when you use 0.6. That might be a bug.

MultiSeries:-series(y, s= 0);

PrimeSum:= proc(X::realcons)
local p:= 2, S:= 1/p, n;
     for n while is(S < X) do
          p:= nextprime(p);
          S:= S+1/p
     end do;
     n
end proc:

PrimeSum(exp(1));
                             10997
ithprime(10997);
                             116423

In terms of efficiency, there is practically no difference between Arrays and Matrices. Arrays can have arbitrary lower bounds and a arbitrary number of dimensions. Matrix lower bounds are always 1, and, of course, they have two dimensions.

There is a profound difference between lists and Vectors. Maple lists are always static, i.e., they cannot be changed. If you try to change a list or just one member of the list, then if Maple allows it at all, then it recreates the entire list from scratch. So if you are going to be changing the elements, use a Vector. If you are not going to change the elements, then it is very slightly more efficient to use lists. Coding lists is substantially more convenient than coding Vectors.

You asked:

How would I specify the allowed types of the elements of that list (for all elements, and for specific ones)?

To specify that all the elements of the list are a specific type, say posint, use ::list(posint). To specify a list whose first element must be a posint, whose second element can be anything, and which must have two elements, use ::[posint, anything]. Note that the word list does not appear in this spec. To specify a list whose second element must be a posint, but whose length is unlimited use 

::satisfies(L-> nops(L) > 1 and L[2]::posint)

Note that any type, no matter how complex, can be specified by satisfies; but it is more efficient to use the built-in types when possible.

And if the element of the list is a list again, how would I specify the type of that's elements? and so forth..

For example, ::list(list(posint)). If you want to specify that all the inner lists have the same length, then use ::listlist(posint). The former concept can be extended to any depth, but there is no listlistlist.

How would I specify the allowed size of the listArray, etc. ?

To specify a list with three elements, use ::[anything, anything, anything]. Note that you can't use to abbreviate this. So an alternative is

TypeTools:-AddType(list3, [anything $ 3]);
P:= proc(L::list3)
....

Arrays (and Vectors and Matrices) are easier to specify a size for: It's simply ::Array(1..3).

How do I specify the types of more than one output, if my proc returns more than one value?

Just specify multiple types separated by commas. For example,

P:= proc(L::list)::list,nonnegint; L,nops(L) end proc;

Note that the return value of procedures and the assigned values of local variables is only checked if you set

kernelopts(assertlevel= 2);

This is useful for debugging during code development. Afterwards, it becomes more efficient when the assertlevel reverts to its default value of 0. The type-checking of arguments does not depend on assertlevel; however, they are only type-checked if they are actually used in the procedure.

Is there also a way to specify that the ineger (or float) argument has to be greater than some given value, e.g. >2, or that it has to lie between two values, or out of a given set?

In general, use the satisfies construct described above. For some limited cases it is better to use a built-in type. For example, to specify an integer greater than 1, use ::And(posint, Not(identical(1))).

 

local D,I,O;
combinat:-randperm([A,B,C,D,E,F,G,H,I,J,K,L,M,N,O,P,Q,R,S]);

If you want a string, then do

convert(combinat:-randperm(19) +~ 64, bytes);

(because "A" is the 65th character in ASCII).

After assigning the above expression to f, do

plot([(numer,denom)(f), z= -1..1]);

This will use the numerator for the horizontal coordinate and the denominator for the vertical. I picked a range for z that seemed to give a good plot, but you might want to change that.

The erf in the solution is correct. The full verification is in the attached worksheet. Also, I solved for y(3).

restart:
ODE:= x^2*(diff(y(x), x, x))+x^3*(diff(y(x), x))+(x^2-2)*y(x) = 0:
Sol:= dsolve(ODE):
Verify solution:
eval(lhs(ODE), Sol);

simplify(%);
                               0
The solution checks.

Solve BVP:
Sol:= dsolve({ODE, y(1)=1, y(2)=2});

eval(Sol, x=3);

 

evalf(%);
                    y(3) = 1.76345278364621

(Whew! I am sick and tired of cut-and-pasting worksheets line by line.)

 

Download erf.mw

MajMaj wrote: pathFILE := "\/MyDirectory\/blabla.txt";

Don't use the backslash to escape the forward slash! That will remove all slashes. That should be

pathFILE := "/MyDirectory/blabla.txt";

As nm says, the forward slash (/) will work with either OS. So, there's no need to use backslash at all.

First 332 333 334 335 336 337 338 Last Page 334 of 395