Preben Alsholm

13653 Reputation

22 Badges

19 years, 310 days

MaplePrimes Activity


These are replies submitted by Preben Alsholm

@nm I experimented a little. Reducing the degree of the denominator from 3/2 helps to prevent crash:

restart;
u:=Int(1/alpha^(3/2)*sin(1/2*3^(1/2)*ln(alpha))*sin(alpha)*3^(1/2),alpha = 0 .. x);
u1:=Int(1/alpha^(3/2-10^(-6))*sin(1/2*3^(1/2)*ln(alpha))*sin(alpha)*3^(1/2),alpha = 0 .. x);
value(u1); # 1499999/1000000 # No crash

I kept going:
 

restart;
u:=Int(1/alpha^(3/2)*sin(1/2*3^(1/2)*ln(alpha))*sin(alpha)*3^(1/2),alpha = 0 .. x);
u2:=Int(1/alpha^(3/2-10^(-16))*sin(1/2*3^(1/2)*ln(alpha))*sin(alpha)*3^(1/2),alpha = 0 .. x);
value(u2); # 3/2-10^(-16) No crash
u3:=Int(1/alpha^(3/2+10^(-16))*sin(1/2*3^(1/2)*ln(alpha))*sin(alpha)*3^(1/2),alpha = 0 .. x);
value(u3); # 3/2+10^(-16) No crash
value(u); #  3/2 Crash

So the exponent of the denominator being exactly 3/2 seems critical.

@nm I still got the crash after having renamed my maple.ini file.

I notice that C_R uses Windows 10. I use Windows 11, but that could hardly make any difference (?)

I get the crash in Windows and with Physics:-Version() 1715 and also with the previous 1714.
I also get the crash if I rename the Physics Updates folder 2024 to something else, in my case 2024_temp.

If I use the workaround proposed by Thomas Richard I don't get any crash.

@Aung Let us look at your equation after first having replaced N with x as variable of integration and pi with Pi.
Since you say that k is a material constant I suppose we may assume that it is positive:
But we only need k>-2.
 

eq:=int((1 - omega)^(k + 1), omega = 0 .. 1) = C*int((1 - sigma*sin(2*Pi*x))^k, x = 0 .. N);
eq2:=simplify(eq) assuming k>-2;
# What do you want to do with eq2?

 

@acer In this simple case the problem is just one of integration. So I tried this plot:
 

plot('int(Mt,0..t,numeric)',t=-10..10);

@omkardpd The plot produced in my answer, where the plot command picks the points, can also be produced like this, ehere we pick the points:
 

LL:=[seq]([z,temp_proc(z,2)],z=-2..3,.01):
plot(LL);
type(LL,listlist); #true
M:=Matrix(LL);
plot(M); 

Incidentally 3d plots can be produced very easily:
 

plot3d('temp_proc(x,y)',x=-2..3,y=-2..3);
plot3d(temp_proc,-2..3,-2..3); # A procedural version

Finding points first seems rather unneccessary, but could be done:
 

S:=seq(seq([x,y,temp_proc(x,y)],x=-2..3,0.1),y=-2..3,0.1):
S[5];
type([S],listlist);
plots:-matrixplot(Matrix([S]));

@ My problem was that I used your older version of sierp.
With the new one sierp and inv_sierp work fine!
I removed my reply when I realized that I was using the old sierp.

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.

2 3 4 5 6 7 8 Last Page 4 of 229