Preben Alsholm

13733 Reputation

22 Badges

20 years, 265 days

MaplePrimes Activity


These are replies submitted by Preben Alsholm

We need to know what MuligIndgange, M1, and RemoveList are.

I think, however, that we need to see the whole worksheet since you mention that the code works outside some procedure but not inside.
So use the fat green uparrow in the panel above to upload the worksheet.
This arrow is only visible while editing your question or reply.

@janhardo As far as I see there is basically only the way I gave when using evalindets.
Notice that you cannot simplify this to cos :
 

f:=t->cos(t);
simplify(eval(f));

But you can do:
 

indets(f(t));
evalindets(f(t),function,s->op(0,s));

I think that in the very old days, certainly several versions before Maple 12, f would just be cos.

@janhardo That can be done:
 

DEtools:-DEplot(sys,[x(t),y(t)],t=0..20,[[x(0) = 1, y(0) = 0]],x=-1..1,y=-1..1,stepsize=0.1);
DEtools:-DEplot(sys,[x(t),y(t)],t=0..20,[[x(0) = 1, y(0) = 0]],x=-1..1,y=-1..1,stepsize=0.1,animatecurves=true);
DEtools:-DEplot(sys2,[x(t),y(t)],t=0..200,[[x(0) = 1, y(0) = 0]],x=-20..1,y=-6..6,stepsize=0.1);
DEtools:-DEplot(sys2,[x(t),y(t)],t=0..200,[[x(0) = 1, y(0) = 0]],x=-20..1,y=-6..6,stepsize=0.1,animatecurves=true);

Using sys and sys2 from my answer above.

@janhardo Here is another road to Rome for the simple example above:
 

restart;
eq1 := diff(x(t), t) = y(t);
eq2 := diff(y(t), t) = -x(t);
beginwaarden := x(0) = 1, y(0) = 0;
sol:=dsolve({eq1, eq2, beginwaarden}, {x(t), y(t)});

L:=subs(sol,[x(t),y(t)]); # List with value of x(t) first
XY:=unapply~(L,t);  # You can use any name on the left.
XY(tau);
plot(XY(t),t=0..2*Pi);
plot(XY,0..2*Pi);
## Now the alternative road:
L;
indets(L,function); # Just warming up.
XY2:=evalindets(L,function,s->op(0,s)); # The alternative
plot(XY2(t),t=0..2*Pi);
plot(XY2,0..2*Pi);

 

@janhardo Is this what you want to accomplish:
 

restart;
eq1 := diff(x(t), t) = y(t);
eq2 := diff(y(t), t) = -x(t);
beginwaarden := x(0) = 1, y(0) = 0;
sol:=dsolve({eq1, eq2, beginwaarden}, {x(t), y(t)});

L:=subs(sol,[x(t),y(t)]); # List with value of x(t) first
XY:=unapply~(L,t);  # You can use any name on the left.
XY(tau);
plot(XY(t),t=0..2*Pi);
plot(XY,0..2*Pi);

 

@janhardo Do this:

sol:=dsolve({eq1, eq2, beginwaarden}, {x(t), y(t)});
             
subs(sol,[x(t),y(t)]); # List with value of x(t) first

if you use rhs elementwise as in rhs~(sol) then the result is a set of the right hand sides of the elements in sol.
That leaves it to set ordering to determine the order, not you.

@nm I suspect that you already know this.
 

restart;
integrand:=(2*x^n+1)/(x^(n+1)+x);
res2021:=convert(eval(integrand,n=2021),parfrac): # Takes a while
length(res2021); #19646
simplify(res2021-eval(integrand,n=2021)); # 0

In fact it also works for n=2023..2030.
for n = 2030 it takes between 10 and 15 minutes.

.

@Thomas Richard I tried replacing the number 2022 by n:
 

restart;
integrand:=(2*x^2022+1)/(x^2023+x);
integrand_n:=(2*x^n+1)/(x^(n+1)+x);
infolevel[int]:=3:
int(integrand_n,x);
eval(%,n=2022);

That revealed that a Risch-Norman integrator was successful.
Continuing with:
 

int(integrand_n,x,method=_RETURNVERBOSE);

gave the information that norman, meijerg,and risch were successful.

@Kitonum Your workaround works in 2024.0 too and is very fast.

@ I  have Mathematica 13.2, but I use it extremely little.
Thus I cannot step in right away.

@Axel Vogt I noticed that you have an earlier build.
The results you get, however, are the same in RC4: `Maple 2024.0, X86 64 WINDOWS, Mar 01 2024, Build ID 1794891`

@Rouben Rostamian  That is very nice.

@ Why is the sequence B in inv_sierp being created?  It isn't used for anything.
B is a sequence of zeros and ones.

Another problem. K is a sequence of 0 and 2 only. So never odd.
Thus a part of the code can be commented out (as well as the B's):

for j to k do 
    km := K[k - j + 1]; # either 0 or 2
    t := t*2^d + cd - km; 
    (*if d = 2 then 
      if evalb(type(km, odd)) then t := 1 - t end if; # Never happens
    end if; 
    if evalf(t) < 0 then t := t + 1 end if; *)  # Never
  end do;

You get the same result as before after commenting this out.
### Note:

I tried using the original inv_sierp code with nothing commented out.
The only change to the procedure was that I made it return t ,[K], [B]
 

return t ,[K],[B]
end proc:
d := 2;
k := 20;
t := 0.6164;
P := sierp(d, k, t);

res:= inv_sierp(d, k, P);

res was:
res := 360798224641/2,
[2, 2, 2, 0, 2, 2, 2, 2, 2, 2, 0, 0, 2, 2, 2, 2, 0, 0, 0, 2],
[1, 0, 1, 1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 1, 0, 1, 0, 0, 0, 0, 1, 0, 0, 1, 1]

@ Yes, sierp now works beautifully.
I tried with N = 4, 5, 6, 7, and 8. It is very nice. For N=8 I changed the thickness to 1.

I'm now convinced that somehow this algorithm has to be iterated.
Try this:
 

restart;
with(Fractals:-LSystem);
with(LSystemExamples);
PlotExample(SierpinskiTriangle,7);
for i from 1 to 10 do PlotExample(SierpinskiTriangle,i) end do;

Notice that PlotExample(SierpinskiTriangle,1);  looks like this:

4 5 6 7 8 9 10 Last Page 6 of 230