Carl Love

Carl Love

28035 Reputation

25 Badges

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

MaplePrimes Activity


These are answers submitted by Carl Love

You need to spell out Heaviside and use exp(...) instead of e^.... The command is

inttrans:-laplace(Heaviside(t-1)*exp(3*t-3), t, s);

The Wronskian command is predefined. See ?Wronskian.

The linalg package has long been deprecated. It's functionality has been replaced by the packages VectorCalculus and LinearAlgebra.

The internal reference to h must be stored as z5:-h because h is in z5, so it needs to know where to look up h on re-evaluation of the expression. But you don't need to see the z5:- if you use a `print/h` procedure. Change z5 to this:

z5:= module()
option package;
export
     h:= F-> expand(evalindets(D(F)/F, specfunc(D), d-> op(d)*'h'(op(d))))
;
local
     ModuleLoad:= proc() :-`print/h`:= ()-> ':-h'(args) end proc,
     ModuleUnload:= ()-> unassign(':-`print/h`')
;
     ModuleLoad();
end module:

Of course, you can include any other locals or exports that you want in the module.

The key thing is that you can get any function f to prettyprint differently by defining a global procedure `print/f`. Since the procedure needs to be global, it needs to be assigned in a ModuleLoad procedure.

Symbolic antiderivatives in Maple are allowed to be discontinuous. This especially happens with complex functions. You can see this discontinuity in this plot (the red line is the real part):

F:= I*arctan(sqrt(exp(2*I*t)-1))-I*sqrt(exp(2*I*t)-1):
plot(([Re,Im])~(F), t= 0..Pi, discont);

The int command correctly adjusts for the discontinuity when it's a definite integral.

In the line that defines vNN, you have

vNN= ....

That needs to be

vNN:= ....

In the next line, use eval instead of algsubs:

pressure:= eval(nMom, vNN);

Do you mean something like this?

restart:

interface(rtablesize= 30):

Matrix((30,2), (i,j)-> `if`(j=1, i, 2^i));

Matrix([[1, 2], [2, 4], [3, 8], [4, 16], [5, 32], [6, 64], [7, 128], [8, 256], [9, 512], [10, 1024], [11, 2048], [12, 4096], [13, 8192], [14, 16384], [15, 32768], [16, 65536], [17, 131072], [18, 262144], [19, 524288], [20, 1048576], [21, 2097152], [22, 4194304], [23, 8388608], [24, 16777216], [25, 33554432], [26, 67108864], [27, 134217728], [28, 268435456], [29, 536870912], [30, 1073741824]])

 


Download Powers_of_2.mw

The variable cape very likely has an unassigned variable in it, so it can't be evaluated to a number. If you give the command print(cape), you'll very likely see what that unassigned variable is. If cape is too long to find that variable, the command indets(cape, name) will help you find it.

If you need further help, please use the green up arrow on the toolbar in the MaplePrimes editor to upload your worksheet. It's very difficult to make an accurate diagnosis based on simply reading the code "by eye".

I see what you're trying to do. It's a very common model for iterative numerical algorithms: You want to limit the number of iterations so that you don't have an infinite loop if the method is not converging. You need to use a for-while-do loop.

restart:

PointFixe:= proc(g, d, a, N)
local o, i, u, pm, relerr;
     pm:= 0.222e-15;
     u:= d;
     relerr:= a;
     for i to N while relerr >= a do
          o:= u;
          u:= evalf(eval(g, x= u)); #eval, not subs
          relerr:= abs((u-o)/(abs(u)+pm));
     od;
     #Check if loop ended due to iteration count rather than convergence.
     if i = N+1 then
          error "Relative error is %1 after %2 iterations.", relerr, N
     end if;
     userinfo(1, PointFixe, "Used", i, "iterations.");
     u
end proc:       

 

infolevel[PointFixe]:= 1: #Activates the userinfo statement.

PointFixe(exp(x)-x-5, 2, 1e-8, 100);

Error, (in PointFixe) Relative error is 0.159964237326880e-3 after 100 iterations.


PointFixe(exp(x)-x-5, 2, 1e-8, 1000);

PointFixe: Used 210 iterations.

-2.45716106995070

 

 

Download PointFixe.mw

Note that the given answer is not correct. The algorithm converged because the step size got too small. So, you'll need to add a check for that. You should have abs(u) being close to 0 to return without error.

You need more square brackets in the plot commands to indicate a parametric plot:

InertiaAxis1 := t->plot([x, V1(t)[2]*x/V1(t)[1], x = -2 .. 2]);
InertiaAxis2 := t->plot([x, V2(t)[2]*x/V2(t)[1], x = -2 .. 2]);

There's a command specifically for this equation. See ?LinearAlgebra,LyapunovSolve.

Here's a one-line procedure that tests whether your lists are valid:

valid:= (L::list(algebraic), V::list(name))->

     not `or`(seq(type(L[k], freeof({V[]} minus {V[k]})), k= 1..nops(V)))
:

See ?type,freeof.

You mentioned monomials. This Question is about polynomials, not monomials.

The coefficients of all the extraneous terms are relatively small. The relative error is ~ 1e-6 (just eyeballing it) over a small range of sigma2. However, it's true that an extra-degree term will ultimately cause unlimited error.  I'm working on some precise error computations, presented in a worksheet below.

The method used for the original determinant is unifloat. (So note that your guess that it was Gaussian elimination was incorrect.) The method used for Axel's method is fracfree. The details of each method that I summarize below can be viewed by setting infolevel[all]:= 5 before doing the computation.

Using the unifloat method, it first determines that the maximum possible degree of the determinant is 9. This estimate is very reasonable based on the maximum degree in every row and column. Both rowwise and columnwise estimates are 9 (as can be seen by inspection). Then it evaluates the matrix's floating-point determinant at 10 points and computes the interpolating polynomial. This method is nearly guaranteed to give you a dense degee-9 polynomial when done with floating-point arithmetic.

The fracfree method proceeds by relatively straightforward fraction-free Gaussian elimination. "Fraction-free" in this context means that no rational functions are used; it doesn't mean that integer coefficients are used.

Here's a worksheet showing the relative errors of the two methods.

restart:

macro(LA= LinearAlgebra):

infolevel[all]:= 5:

A:= <
     180*sigma2^2, -17980.4676072929*sigma2, 60*sigma2^2,
          -1792.74593023400*sigma2, -597.004097753022*sigma2, -60*sigma2;
     -17980.4676072929*sigma2, 60*sigma2^2, -17980.4676072929*sigma2,
          -597.004097753022*sigma2, -597.581976744666*sigma2, 5993.48920243096;
     60*sigma2^2, -17980.4676072929*sigma2, 180*sigma2^2,
          -597.581976744666*sigma2, -1791.01229325907*sigma2, -60*sigma2;
     -1792.74593023400*sigma2, -597.004097753022*sigma2, -597.581976744666*sigma2,
          -60*sigma2, 5993.48920243096, 597.581976744666;
     -597.004097753022*sigma2, -597.581976744666*sigma2, -1791.01229325907*sigma2,
          5993.48920243096, -60*sigma2, 597.004097753022;
     -60*sigma2, 5993.48920243096, -60*sigma2,
          597.581976744666, 597.004097753022, 60
>;

A := Matrix(6, 6, {(1, 1) = 180*sigma2^2, (1, 2) = -17980.4676072929*sigma2, (1, 3) = 60*sigma2^2, (1, 4) = -1792.74593023400*sigma2, (1, 5) = -597.004097753022*sigma2, (1, 6) = -60*sigma2, (2, 1) = -17980.4676072929*sigma2, (2, 2) = 60*sigma2^2, (2, 3) = -17980.4676072929*sigma2, (2, 4) = -597.004097753022*sigma2, (2, 5) = -597.581976744666*sigma2, (2, 6) = 5993.48920243096, (3, 1) = 60*sigma2^2, (3, 2) = -17980.4676072929*sigma2, (3, 3) = 180*sigma2^2, (3, 4) = -597.581976744666*sigma2, (3, 5) = -1791.01229325907*sigma2, (3, 6) = -60*sigma2, (4, 1) = -1792.74593023400*sigma2, (4, 2) = -597.004097753022*sigma2, (4, 3) = -597.581976744666*sigma2, (4, 4) = -60*sigma2, (4, 5) = 5993.48920243096, (4, 6) = 597.581976744666, (5, 1) = -597.004097753022*sigma2, (5, 2) = -597.581976744666*sigma2, (5, 3) = -1791.01229325907*sigma2, (5, 4) = 5993.48920243096, (5, 5) = -60*sigma2, (5, 6) = 597.004097753022, (6, 1) = -60*sigma2, (6, 2) = 5993.48920243096, (6, 3) = -60*sigma2, (6, 4) = 597.581976744666, (6, 5) = 597.004097753022, (6, 6) = 60})

gc():

P[unifloat]:= CodeTools:-Usage(LA:-Determinant(A)):

IndexList: initializing external call for IndexList

IndexList: external IndexList call successfully initialized
ComputeConnect: initializing external call for ComputeConnect
ComputeConnect: external ComputeConnect call successfully initialized
StronglyConnectedBlocks: Decomposed into blocks of sizes  [6]
Determinant/unifloat: dimension is 6 domain is R[sigma2] evaluation/interpolation requires 9 determinants over R
Determinant/unifloat: R[sigma2]: proceeding with determinant 1
Determinant/float: calling external function
Determinant: NAG hw_f07arf
Determinant/unifloat: R[sigma2]: proceeding with determinant 2
Determinant/float: calling external function
Determinant: NAG hw_f07arf
Determinant/unifloat: R[sigma2]: proceeding with determinant 3
Determinant/float: calling external function
Determinant: NAG hw_f07arf
Determinant/unifloat: R[sigma2]: proceeding with determinant 4
Determinant/float: calling external function
Determinant: NAG hw_f07arf
Determinant/unifloat: R[sigma2]: proceeding with determinant 5
Determinant/float: calling external function
Determinant: NAG hw_f07arf
Determinant/unifloat: R[sigma2]: proceeding with determinant 6
Determinant/float: calling external function
Determinant: NAG hw_f07arf
Determinant/unifloat: R[sigma2]: proceeding with determinant 7
Determinant/float: calling external function
Determinant: NAG hw_f07arf
Determinant/unifloat: R[sigma2]: proceeding with determinant 8
Determinant/float: calling external function
Determinant: NAG hw_f07arf
Determinant/unifloat: R[sigma2]: proceeding with determinant 9
Determinant/float: calling external function
Determinant: NAG hw_f07arf
Determinant/unifloat: R[sigma2]: proceeding with determinant 10
Determinant/float: calling external function
Determinant: NAG hw_f07arf
subtype: TESTING: anything <<< nonreal
subtype: AFTER SIMPLIFICATION: anything <<< nonreal
memory used=1.58MiB, alloc change=0 bytes, cpu time=15.00ms, real time=61.00ms, gc time=0ns

gc():

P[fracfree]:= CodeTools:-Usage(LA:-Determinant(convert(A, rational, exact))):

StronglyConnectedBlocks: Decomposed into blocks of sizes  [6]

Determinant/fracfree: gaussian elimination
Determinant/fracfree: elimination at row 1
Determinant/fracfree: pivot selected -60*sigma2
Determinant/fracfree: elimination at row 2
Determinant/fracfree: pivot selected (3/2500000000)*sigma2^2
Determinant/fracfree: elimination at row 3
Determinant/fracfree: pivot selected (2697070141093941/31250000)*sigma2^4
Determinant/fracfree: elimination at row 4
Determinant/fracfree: pivot selected (2689118895350997/312500000000000000000)*sigma2^5+(640658284733379319059892072243256671367577/156250000000000000000000000000000)*sigma2^4
Determinant/fracfree: elimination at row 5
Determinant/fracfree: pivot selected -(4854726253969069605731308228221/15625000000000000000)*sigma2^6-(57732251822958065318061323918368069370263980820815882906453/312500000000000000000000000000000000000000000)*sigma2^5-(85818076159321533072130690706468745958854769498084423957090495217670196273/3125000000000000000000000000000000000000000000000000000000)*sigma2^4
memory used=351.55KiB, alloc change=0 bytes, cpu time=0ns, real time=37.00ms, gc time=0ns

infolevel[all]:= 0:

P[fracfree]:= evalf(P[fracfree]):

R:= LA:-RandomVector(99, generator= -2.0..2.0):

d[e]:= map(r-> LA:-Determinant(eval(A, sigma2= r)), R):

d[u]:= map(r-> eval(P[unifloat], sigma2= r), R):

d[ff]:= map(r-> eval(P[fracfree], sigma2= r), R):

So, over the 99 evaluation points, the maximum relative error of the unifloat method is

LA:-Norm((d[e] -~ d[u]) /~ d[e], infinity);

HFloat(1.083377744650285e-6)

(That's amazingly close to my "eyeballing" estimate of 1e-6.)

And the maximum relative error of the fracfree method is

LA:-Norm((d[e] -~ d[ff]) /~ d[e], infinity);

HFloat(1.6430641959088563e-14)

Considering these errors and the time and memory used for each determinant computation, the convert(..., rational, exact) method is the clear winner for a 6x6 univariate polynomial matrix, but the default method isn't grossly erroneous or ineffcient either.

 

Download poly_det_error.mw

plots:-implicitplot3d(
     (x-y)^2+(x-z)^2+(z-x)^2=3, x= -2..2, y= -2..2, z= -2..2,
     grid= [20$3], style= surface
);

Here is the solution. If you had Maple 2016, I would've used a DataFrame rather than a Matrix for the final output. That would make it eaier to read with column and row labels.

restart:

macro(ST= StringTools, PD= StringTools:-PatternDictionary, LT= ListTools, St= Statistics):
bid:= PD:-Create(ospd3):
words:= select(
     w-> ST:-AndMap(ST:-IsAlpha, w), #Remove words with non-alpha characters.
     [seq(ST:-UpperCase(PD:-Get(bid,i)), i= 1..PD:-Size(bid) - 1)]
):

nops(words);

80458

L:= sort~((St:-Tally@op~@[op])~(LT:-Classify(nops, ST:-Explode~(words))), key= -rhs):

c:= max(indices(Data, nolist));

39

interface(rtablesize= max(26, c)):

Matrix(
     (26,c),
     (i,j)-> `if`(assigned(L[j]) and nops(L[j]) >= i, convert(parse(lhs(L[j][i])), `local`), ``)
);

Matrix([[I, A, A, A, E, E, E, E, E, E, E, E, E, E, E, I, I, E, E, E, E, E, A, E, I, E, I, A, ``, ``, E, ``, I, ``, ``, O, ``, I, I], [A, I, E, E, A, A, A, A, A, I, I, I, I, I, I, E, T, I, O, O, R, O, N, I, E, I, A, I, ``, ``, O, ``, T, ``, ``, T, ``, T, T], [``, O, O, O, S, S, S, S, I, A, T, A, O, R, T, N, E, O, N, A, O, A, I, R, A, N, E, L, ``, ``, R, ``, A, ``, ``, I, ``, A, A], [``, E, I, S, O, R, I, I, T, O, R, T, T, O, R, T, N, T, T, I, I, L, E, T, O, S, O, E, ``, ``, T, ``, L, ``, ``, M, ``, L, E], [``, M, T, L, R, I, R, R, R, T, A, R, R, T, O, O, R, L, A, N, P, R, T, L, N, L, L, O, ``, ``, F, ``, E, ``, ``, E, ``, E, R], [``, N, N, I, L, O, O, O, O, R, O, O, A, A, N, R, A, R, R, L, S, T, L, C, T, O, R, R, ``, ``, N, ``, O, ``, ``, N, ``, M, S], [``, T, S, R, I, L, N, N, N, N, N, N, N, N, A, S, O, S, I, R, A, P, R, P, R, G, N, T, ``, ``, D, ``, R, ``, ``, A, ``, O, L], [``, S, P, N, N, N, L, L, S, S, S, S, S, S, S, A, D, N, L, T, C, C, S, U, M, X, T, P, ``, ``, I, ``, U, ``, ``, L, ``, R, O], [``, D, D, T, T, D, T, T, L, L, C, C, C, C, C, C, M, A, D, D, T, I, O, S, L, T, M, S, ``, ``, A, ``, D, ``, ``, P, ``, U, P], [``, C, U, D, D, T, D, C, C, C, L, L, P, L, L, L, C, H, P, M, N, N, M, A, U, F, P, M, ``, ``, M, ``, M, ``, ``, X, ``, D, U], [``, R, M, U, C, C, C, D, U, U, P, P, L, P, U, P, P, C, S, S, L, M, C, M, P, D, C, U, ``, ``, L, ``, P, ``, ``, C, ``, P, D], [``, U, R, M, U, U, U, U, D, P, M, M, M, M, P, U, S, P, M, W, H, H, F, V, S, A, F, N, ``, ``, U, ``, B, ``, ``, R, ``, B, M], [``, W, B, P, M, M, G, M, M, M, U, U, U, U, M, G, L, D, C, C, U, Y, U, Y, Y, ``, U, V, ``, ``, Q, ``, N, ``, ``, B, ``, N, B], [``, F, L, B, P, G, M, G, P, D, H, H, H, H, H, D, G, G, V, G, M, U, G, H, D, ``, D, Y, ``, ``, X, ``, S, ``, ``, J, ``, S, N], [``, P, G, C, Y, P, P, P, H, H, D, D, D, D, B, H, U, U, G, P, W, V, P, ``, V, ``, V, ``, ``, ``, C, ``, V, ``, ``, Y, ``, V, V], [``, L, C, H, B, B, B, B, G, G, B, G, G, B, G, M, H, F, U, U, G, X, Y, ``, Q, ``, Y, ``, ``, ``, P, ``, Y, ``, ``, ``, ``, Y, Y], [``, H, H, G, H, H, H, H, B, B, G, B, Y, Y, D, Y, F, M, X, V, V, G, X, ``, F, ``, B, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``], [``, B, F, K, G, Y, F, F, Y, Y, Y, Y, B, G, V, X, V, B, F, Y, ``, B, B, ``, Z, ``, K, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``], [``, Y, Y, F, K, F, Y, Y, F, V, F, F, F, F, Y, F, X, Y, B, H, ``, S, D, ``, H, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``], [``, G, W, Y, F, K, K, K, W, F, V, V, V, V, J, B, W, Q, Y, F, ``, D, H, ``, W, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``], [``, K, K, W, W, W, W, V, K, W, K, K, K, X, F, W, B, X, Q, B, ``, ``, ``, ``, C, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``], [``, V, V, V, V, V, V, W, V, K, W, X, W, W, W, Q, Z, W, ``, K, ``, ``, ``, ``, B, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``], [``, X, X, J, J, Z, Z, Z, X, X, X, W, Q, K, Q, V, J, V, ``, Q, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``], [``, J, J, Z, Z, J, J, X, Z, Q, Q, Q, X, Z, X, K, Y, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``], [``, Z, Z, X, X, X, X, J, Q, J, Z, Z, Z, Q, K, J, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``], [``, Q, Q, Q, Q, Q, Q, Q, J, Z, J, J, J, ``, Z, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``, ``]])

All the numeric data is stored in L:

LTT:= table~(L):

LTT[7]["E"];

14338

So there 14,338 Es in seven-letter words.

``


Download length_by_letter.mw

Another way, if you have trouble remembering the @:

evalf(Int(x-> h(g(x)), x= 1..2);

Using unapply with quotes on its first argument is pointless (and slightly wasteful). You could just as well define g via

g:= a-> fsolve(a*y^2 - sin(y), y=2);

So, I wonder where this post that recommends unapply is so that I can go correct it.

First 208 209 210 211 212 213 214 Last Page 210 of 395