Carl Love

Carl Love

28020 Reputation

25 Badges

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

MaplePrimes Activity


These are answers submitted by Carl Love

My, what a trivial procedure, hardly worth naming, let alone writing a help page for. Anyway, here it is:

Trapz:= (X, Y)-> (X[2..]-X[..-2]).(Y[2..]+Y[..-2])/2:
trapz:= proc()
local X,Y;
    if nargs=1 then Y:= args[1]; X:= Vector([upperbound(Y)][1], i-> i)
    else X:= args[1]; Y:= args[2] 
    fi;
    if (local d:= rtable_num_dims(Y)) = 1 then Trapz(X,Y) 
    elif d=2 then Vector[row]([upperbound(Y)][2], i-> Trapz(X, Y[..,i]))
    fi
end proc:

I didn't handle the case where has more than 2 dimensions.

You can enter the system, numerically solve it, and plot the solution in 3D like this:

sys:= { 
    diff(x(t),t) = cos(x(t)), x(0) = 1, 
    diff(y(t),t) = sin(z(t)), y(0) = -1,
    diff(z(t),t) = z(t)*y(t), z(0) = 2
}:
sol:= dsolve(sys, numeric):
plots:-odeplot(sol, [x,y,z](t), t= 0..5);

 

If you want to solve it with Picard iteration, you'll need to code that yourself. Here's the plot produced by the above:

Do this:

f:= (x-1)/(x+2);
s:= convert(f, FormalPowerSeries);
c:= unapply((S-> (op(1,S),[op([2,1],S)]))(indets(s, specfunc(Sum))[]));
L:= limit(abs(c(k+1)/c(k)), k= infinity);
solve(L < 1, x); 
                   (-2,2)

Notes:

1. None of your simplify commands were needed, but neither would it be wrong for you to put them back.
2. Your line sum(k-> 1/k^2, k= 1..1000) is nonsense, and its result, 1000*f, is also nonsense. The fact that this doesn't give an error message should be considered a bug in Maple.
3. The line the turns the general term of the series into the function c is now fully automatic, no need to copy-and-paste.
4. The expression for the ratio test is abs(c(k+1)/c(k)). The extra powers that you used are mathematical errors.
5. The final result is given in interval notation, the open interval from -2 to 2.
 

Here is how you can "automatically" add a plot to an existing plot:

AddPlot:= (p, P)-> plots:-display(P, p, _rest):
f:= (x-1)/(x+2):
p1:= plot(f, x= -3..3, discont, view= [-3..3, -5..5]):
p2:= AddPlot(
    plottools:-circle([0,0], 1, color= blue),
    p1,
    scaling= constrained
);

The following simple code is probably adequate for your needs. It could be made slightly more efficient by pre-sieving for the primes.

MyProperty:= proc(q::prime, P::name) 
local p:= ithprime(q - numtheory:-pi(q));
    if p >= q then P:= p; true else false fi
end proc
:    
MySeq:= proc(n::posint)
local k:= 0, R:= Vector(n), p, q:= 1;
    while k < n do 
        q:= nextprime(q);
        if MyProperty(q, 'p') then
            k:= k+1;
            R[k]:= [p,q];
        fi
    od;
    [seq](R)
end proc
:
S:= MySeq(9);
S := [[2, 2], [13, 11], [17, 13], [29, 17], [31, 19], [43, 23], 
  [67, 29], [71, 31], [97, 37], [107, 41]]

(p_seq, q_seq):= map(op~, [1,2], S)[];
   p_seq, q_seq := [2, 13, 17, 29, 31, 43, 67, 71, 97, 107], 
     [2, 11, 13, 17, 19, 23, 29, 31, 37, 41]


 

In Maple, type checking can only be done at run time. This is mostly because types can be arbitrarily complex, and checking them may involve execution of an arbitrary amount of Maple code. Of course, for the example that you show, it'd be theoretically possible before runtime; but I'm sure that you realize that your example is extremely superficial. 

The code below allows you to set the values of the parameters with sliders and have the plot update automatically. For some values of the parameters, there will be warning messages such as "cannot evaluate the solution further (left/right) of ... probably a singularity." These warnings do not prevent you from proceding with the Explore. If you want, the field arrows could be added.

restart:.
A:= a_1*c^3 - b_1:  B:= 3*a_1*c^2 - 3*b_1/c:  C:= d*c - e:
sys:= {diff(x(t), t$2) = (B*x(t)^2+C*x(t)-e)/A, x(0)=f, D(x)(0)=g}:
P:= [a_1, c, b_1, d, e, f, g]:
sol:= dsolve(sys, numeric, parameters= P):
MyPlot:= proc(P)
    sol(parameters= P);
    plots:-odeplot(
        sol,
        [x(t),diff(x(t),t)],
        t= -3..3,
        view= [-6..6, -6..6],
        thickness= 2, axes= boxed, color= red, scaling= constrained
    )
end proc
:
Explore(MyPlot(P), parameters= (P=~ -5..5));

 

subsindets(B, string, parse)

As Tom admits, explicitly entering the expression to be integrated is "ugly", and I find Axel's modification to be only a slight improvement. And that ugliness is made worse by the fact that the OP actually has a large number of such sets of integrals to do that I removed from the posted worksheet because the syntax issue is identical even though the integrands are not.

I like to build the expression to be integrated directly into the ODE system, like this:

restart;
eq1:=diff(f(y), y$4)+Uhs*diff(E(y),y$3)-(diff(f(y), y$2))+(diff(theta(y), y$1))= 0:
eq2:=diff(theta(y), y$2)+(diff(f(y), y$2)+1)^2+1+diff(theta(y),y$2) = 0:
E:=y->zeta*(cosh(k/2*(h1+h2-2*y)))/(cosh(k/2*(h1-h2))):
bcs:=f(h1) = -(1/2)*(Q-1-d), f(h2) = (1/2)*(Q-1-d), (D(f))(h1) = -1, (D(f))(h2) = -1,theta(h1) = 0, theta(h2) = 1:
epsilon1:=0.1:d:=1:omega:=Pi/6:h1:=-(1+epsilon1*sin(2*Pi*x)):h2:=d+epsilon2*sin(2*Pi*x+omega):F:= Q-1-d:epsilon2:=0.5:x:=1:
alpha:=Pi/6:
de:=eq1,eq2,bcs:
d1 := subs(Uhs =-2, zeta=3,k=1,{de}):
param:= {Uhs =-2, zeta=3,k=1}:
P1:= eval(diff(f(y), y$3)+Uhs*diff(E(y),y$2)-(diff(f(y), y$1)+1)+(theta(y))+sin(alpha),param):
ec:=0.5:
for Q from -3 to 3 by ec do 
    F2[Q]:= dsolve(
        d1 union {diff(J(y),y)=P1, J(h1)=0}, 
        numeric, maxmesh= 25500, continuation= lambda1
    ) 
end do:
for Q from -3 to 3 by ec do 
    P3[Q]:= eval(J(y), F2[Q](1)) - eval(J(y), F2[Q](0))  
end do;
                   P3[-3] := 2.30876985489011
                  P3[-2.5] := 1.36946865887356
                 P3[-2.0] := 0.460748411800256
                 P3[-1.5] := -0.417186862674401
                 P3[-1.0] := -1.26412891537822
                 P3[-0.5] := -2.07986512872413
                  P3[0.] := -2.86417840284371
                  P3[0.5] := -3.61684700714295
                  P3[1.0] := -4.33764445144679
                  P3[1.5] := -5.02633933429334
                  P3[2.0] := -5.68269517762163
                  P3[2.5] := -6.30647027769476
                  P3[3.0] := -6.89741755441151

 

You can place the cursor wherever you want (such as after the x in dx), press backspace, and change it to y. Likewise, place the cursor before the "d" to enter/edit the integrand, then start a new definite integral there. There's no need to use parentheses, but I suppose that you could if you really wanted to. In the hyperlink posted below, I have entered the expression which in standard Maple Notation would be

Int(Int(f(x,y), x= c(y)..d(y)), y= a..b)

https://learn.maplesoft.com/index.html#/?d=DQJKHPFIKPOTBIGPKTGHCFFPPTKSDPBSDGLKFNFUFIGUKNBQAFDLJOOKNHDUERMQKKCSJKPGASGFEIPQLNCUIMMRHKBFMHHGANER

I think that the problem is simply that you've made a mental arithmetic error in determing your expected result, because x^(2*k)/x^k is x^k, not x^2.

An unrelated tip: You shouldn't assume that the bound variable of the Sum is literally (even if that's always true). Instead, extract it via op([2,1], s) or, equivalently, lhs(op(2,s)).

The code that acer wrote for you that uses operators such as #mi(...) and #mo(...was intended to display output via Maple's GUI, not via Latex. I vaguely recall him explicitly saying that, but upon rereading, I can't find it, so maybe he didn't say it. Anyway, I'm certain that that was his intention. Those operators are part of the language HTML (I believe), and not the language Latex.

Now I wouldn't be at all surprised if Latex had some facility for interpretting HTML. If it does, it's likely possible to write a brief Maple procedure to edit Maple's Latex output in such a way that that facility was turned on.

From my personal observation of the expansion of text-formatting languages, there was a rapid divergence about 25 years ago between those primarily intended for paper documents (such as Latex and PDF) and those primarily intended for on-screen documents (such as HTML and XML). Note that I said "primarily": Almost all of them can be used on either medium, but for almost all of them there's one medium that's distinctly much easier. I'd be quite interested in reading here anyone else's opinion on this divergence.

Your main error is attempting to define a function (or procedure) using this syntax:

f(...):= ...;  #bad way

Unfortunately, the above is not an outright error, and Maple allows you to "get away with" doing this to some extent, and then the new user becomes habituated to it. But it doesn't always work, and, to my mind, should never be used to define a procedure (the above syntax has other valid and perfectly acceptable uses that aren't further discussed in this article).

There are several correct ways to define a procedure:

  1. f:= x-> e #univariate case
    where e is any explicitly specified expression containing x.
  2. f:= (x, y, ...)-> e 
    allows any number (including 0) of parameters (aka variables).
    must be still be an explicitly specified expression of the parameters.
  3. f:= proc(x, y, ...) s1; s2; ...; e end proc
    This allows for arbitrary length and number of instructions. e doesn't quite need to be explicit, but expressions created outside the procedure using the procedure's parameters are not allowed.
  4. f:= unapply(e, [x, y, ...])
    Now e doesn't need to be explicit, and it can be created in a separate place.
  5. f:= subs([x= e1, y= e2, ...], eval(P))
    where is an existing procedure containing undeclared global variables x, y, ..., and e1, e2, ... are any expressions with almost no restrictions (they can even use P's parameters).

And there are even further ways that I'd classify as "meta-programming".

So, here are 4 ways to perform your task, in increasing order of sophistication. Actually, the individual components of each of these are fairly independent, so you can mix-and-match.

Download func_norms.mw

MaplePrimes won't expand my worksheet inline, but you should download it from the link immediately above. Nonetheless, here it is with some of its output reformatted, and the plots removed (simply for brevity; anyway, they're not interesting plots; the interest lies in the code that creates them).

#Method 1: Doesn't create procedures at all.
restart: 
f:= x^2:
n:= int(f^2, x= 0..1):
f:= f/n;
                                          2
                                  f := 5 x 

plot(f, x= 0..1);

#Method 2: Creates procedures with -> and unapply. 
restart:
f:= x-> x^2:
n:= int(f^2, 0..1):
f:= unapply(f(x)/n, x);
            f := x -> 5*x^2

plot(f, 0..1);

#Method 3: Uses a proc procedure to abstract (slightly) the concept of 
#"functional normalization". This procedure optionally displays the plot
#as a side effect with the "normalized" function as the return value.
#
#Attn nitpickers: Yes, I know that this isn't a true norm. I'm just
#working with the OP's formulation. A true norm is covered in Method 4.
restart:
normalize1:= proc(f, ab::range(algebraic), N)
local x, r:= unapply(f(x)/int(N@f, ab), x);
    userinfo(1, procname, NoName, print(plot(r, ab)));
    eval(r)
end proc
:
#Set next line to 1 if you want to see the plot, 0 if you don't:
infolevel[normalize1]:= 1:

f:= x-> x^2:
f:= normalize1(f, 0..1, f-> f^2);
            f := x -> 5*x^2 

#Method 4: Further abstraction of the normalization concept so that the
#whole norm can be specified as a procedure.
restart:
normalize2:= proc(f, N, {Popts::list:= []})
local x, r:= unapply(f(x)/N(f, _rest), x);
    userinfo(1, procname, NoName, print(plot(r, Popts[])));
    eval(r)
end proc
:
infolevel[normalize2]:= 1:
f:= x-> x^2:

#Create a standard p-norm (only standard if p>=1):
N:= (f, ab::range(algebraic), p::algebraic)-> 
    int((f-> f^p)@abs@f, ab, 'continuous')^(1/p)
:
f:= normalize2(f, N, 0..1, 2, Popts= [0..1]);
         f := x -> x^2*5^(1/2) 

infolevel[normalize2]:= 0: #Turn off plotting.
f:= x-> x^2:
f:= normalize2(f, N, 0..1, p); #p is now symbolic.
    f := x-> x^2/(1/(2*p+1))^(1/p)

#Instantiate p (i.e., give it a specific value):
eval(f(x), p= 2);
                                   2  (1/2)
                                  x  5     

#Another way to instantiate: subs into an existing procedure:
#This is not that useful in this particular example. I just wanted to
#briefly illustrate the 5th method of procedure creation that I itemized at the 
#beginning of this article.
f:= subs(p= 2, eval(f));
         f := x -> x^2*5^(1/2)

#But be aware that the "subs" method is quite tricky to master, yet 
#extremely powerful once you do so.


 

See ?readlib. I think that this command has been obsolete for well over 20 years. But older usages, such as some appearing in dsolve, are still in some library code. Calls to readlib do nothing, and they're not needed anyway, so they're harmless.

This just requires a simple wrapper that resets assertlevel twice with the call to FetchAll in between:

MyFetchAll:= proc()
local 
    al:= kernelopts(':-assertlevel'= 1), 
    r:= Database:-SQLite:-FetchAll(args)
;
    kernelopts(':-assertlevel'= al);
    r
end proc
:    
First 79 80 81 82 83 84 85 Last Page 81 of 395