## 22182 Reputation

15 years, 230 days

## Social Networks and Content at Maplesoft.com

"My friend George Mallory… once did an inexplicable climb on Snowdon. He had left his pipe on a ledge, half-way down one of the Liwedd precipices, and scrambled back by a short cut to retrieve it, then up again by the same route. No one saw what route he took, but when they came to examine it the next day for official record, they found an overhang nearly all the way. By a rule of the Climbers' Club, climbs are never named in honour of their inventors, but only describe natural features. An exception was made here. The climb was recorded as follows: 'Mallory's Pipe, a variation on route 2; see adjoining map. This climb is totally impossible. It has been performed once, in failing light, by Mr G. H. L. Mallory.'." -- "Goodbye to All That", Robert Graves

## numeric error (eg. roundoff)...

Performing exact integration of an expression containing floating-point coefficients is a bad idea in general.

In your example the result contains an imaginary floating-point component due to numeric error (eg. roundoff error). That decreases with high enough working precision for the computation.

The real component also has numeric error, but of course you notice that less in the print-out.

The details vary a little according to the chosen method passed to exact int. (The choice method=contour can result in a purely real result for this example, but suffers from additional numeric error for the real component.)

Using indefinite integration or definite integration with unknown names for the end-points -- possibly along with conversion of float coefficients to exact rationals) one can obtain an expression whose imaginary component can be "simplified away". But the result from that may suffer from even more numeric error when the numeric end-points are finally unitlized.

A generally better choice is to perform purely numeric integration.

```f := -.5241663257+.6125225520*x-0.7468538250e-1*(2*x-3)^2
+0.1319812500e-1*(2*x-3)^3-sin(log(x)):

int(f^2, x = 1 .. 2, numeric);

-8
2.301629098 10  ```

Using Maple 17.02,  continuous_function_ac.mw

## method=ftocms...

```restart;

kernelopts(version);

Maple 17.02, X86 64 LINUX, Sep 5 2013, Build ID 872941

combine(int( sin(log(x+1)),x=-1..1,method=ftocms));

-cos(ln(2)) + sin(ln(2))```

compute_integral_ac.mw

You can also apply FTOC yourself, using the indefinite integral and taking limits (using either limit or MultiSeries:-limit in your Maple 17, it seems).

## possible idea...

Did you try constructing the list entries as function calls, when using seq?

```restart;
alias(a1 = a1(r), a2 = a2(r), a3 = a3(r)):

array3 := [seq(cat(a, i)(r), i = 1 .. 3)];

diff(array3, r);```

A somewhat dubious alternative might be,

```restart;
alias(a1 = a1(r), a2 = a2(r), a3 = a3(r)):

array4 := [seq(parse(cat(a, i)), i = 1 .. 3)];

diff(array4, r);```

I find it difficult to be more precise in suggestions because I don't see an explanation of why you've chosen to use alias.

## alternative...

It is not necessary to store all the (possibly many and large) procedures generated by instantiating at different numeric values for parameter delta, in order to get decent speed-up by avoiding setting the parameter delta inefficiently.

That requires later clean-up (forget) for the memory to be cleared, and for some examples is an extra potential allocation issue. Such memory concerns and the use of advanced printf functionality are unnecessary for accomplishing this plotting task reasonably efficiently. The OP has made no request to store such procedures externally, and that's a pretty rare task.

Instead of the above, either name (here t or delta) can be put on the first axis in the 3D plot, while still getting the speedup of setting values efficiently for parameter delta. The parametric calling sequence of the plot3d command allows for parameter delta to be used more efficiently as the first ("outer") range, while still flexibly allowing either t,delta or delta,t being ascribed to axis,axis respectively.

This below is done in Maple 2017.3, since the OP was using that version a while back.

```restart;
phi:=0:lambda:=0.1:N:=5:M:=sqrt(N*(N+1))*exp(I*phi):omegap:=10:
var:={n(t),u(t),z(t)}:
dsys:={diff(n(t),t)=-2*(n(t)-N)+(u(t)-M)*exp(-2*I*omegap*t/lambda)
+((z(t)-conjugate(M))*exp(2*I*omegap*t/lambda)),
diff(u(t),t)=-2*(1-I*delta)*u(t)
+2*(n(t)-N)*exp(2*I*omegap*t/lambda)+2*M,
diff((z(t),t))=-2*(1+I*delta)*z(t)
+2*(n(t)-N)*exp(-2*I*omegap*t/lambda)
+2*conjugate(M)}:

res1:=dsolve(dsys union {n(0)=0,u(0)=0,z(0)=0},
numeric, parameters=[delta], output=listprocedure):

Nt:=eval(n(t),res1):

Ntd:=proc(t,delta)
if not [args]::list(numeric) then 'procname'(args);
else
if rhs~(Nt(parameters))<>[delta] then
Nt(parameters=[delta]);
end if;
Nt(t);
end if;
end proc:

CodeTools:-Usage(plot3d([t, delta, Ntd(t,delta)],
delta=-10..10, t=0..1, axes=boxed));
memory used=8.44GiB, alloc change=36.00MiB, cpu time=64.89s,
real time=64.91s, gc time=5.06s``` ```CodeTools:-Usage(plot(Ntd(.5,delta), delta=-10..10,
memory used=4.66GiB, alloc change=80.00MiB, cpu time=35.69s,
real time=34.11s, gc time=3.65s``` I wrote the above for Maple 2017.2, but in Maple 2019.2 (and later) some improvement to dsolve's runtime parameter handling allows for even simpler and more understandable code. The conditional test to check whether the parameter value has changed is no longer needed above in order to get the speedup from that use of the parametric call to plot3d. So the following is enough, paring it down while keeping it understandable (IMNSHO):

```phi:=0:lambda:=0.1:N:=5:M:=sqrt(N*(N+1))*exp(I*phi):omegap:=10:
dsys:={diff(n(t),t)=-2*(n(t)-N)+(u(t)-M)*exp(-2*I*omegap*t/lambda)
+((z(t)-conjugate(M))*exp(2*I*omegap*t/lambda)),
diff(u(t),t)=-2*(1-I*delta)*u(t)
+2*(n(t)-N)*exp(2*I*omegap*t/lambda)+2*M,
diff((z(t),t))=-2*(1+I*delta)*z(t)
+2*(n(t)-N)*exp(-2*I*omegap*t/lambda)
+2*conjugate(M)}:

res1:=dsolve(dsys union {n(0)=0,u(0)=0,z(0)=0},
numeric, parameters=[delta], output=listprocedure):
Nt:=eval(n(t),res1):
Ntd:=proc(t,delta) Nt(parameters=[delta]); Nt(t); end proc:
plot3d([t, delta, 'Ntd'(t,delta)], delta=-10..10, t=0..1);```

## style=line...

It's not clear whether you want the dashed line segments.

But, instead of shading=none and transparency=1 perhaps you are looking for the effect from style=line and color=black, for the cuboid.

For example,

```fig1:=plots:-textplot3d([[0, 0, 0, "B", 'align' = '{below, left}'],
[0, 0, 1, "A", 'align' = '{above, left}'],
[0, 1, 0, "C", 'align' = '{below, right}'],
[0, 1, 1, "D", 'align' = '{right}'],
[1, 0, 0, "F", 'align' = '{left}'],
[1, 0, 1, "E", 'align' = '{above, left}'],
[1, 1, 0, "G", 'align' = '{below, right}'],
[1, 1, 1, "H", 'align' = '{above, right}']],
'font' = ["times", "roman", 20], 'color' = "Red"):

fig2:=plottools:-cuboid([0, 0, 0], [1, 1, 1],
style = line, color = black):

plots:-display([fig1, fig2], axes = none, size = [350, 350],
orientation = [179, 15, 165]);``` ## is, simplify...

Your two expressions are not always equal.

```restart;

kernelopts(version);

Maple 2021.0, X86 64 LINUX, Mar 5 2021, Build ID 1523359

f__1 := sqrt(4*a^2 + lambda__g^2)*c/(2*lambda__g*a):

f__2 := c*sqrt(1/lambda__g^2 + 1/(4*a^2)):

simplify(f__1-f__2) assuming lambda__g>0, a>0;

0

simplify(f__1-f__2) assuming lambda__g<0, a<0;

0

combine(simplify(f__1-f__2)) assuming real, lambda__g*a>0;

0

is(f__1-f__2=0) assuming real, lambda__g*a>0;

true

# example where not equal
eval([f__1, f__2], [lambda__g=-1, a=1, c=1]);

[  1  (1/2)  1  (1/2)]
[- - 5     , - 5     ]
[  2         2       ]
```

## not easily...

This kind of question has been discussed before, eg. here and here.

Unfortunately the brief answer is no, there is no easy mechanism for that.

But it can be kludged using rescaling and forced tickmarks. One can even write procedures to do that automatically.

```Porig:=plot3d(abs(Zeta(x+y*I)),x=-4..5,y=-10..40,
view=[-4..5,-10..40,0..4],
orientation=[-15,75,0]):
Porig;``` ```xl,xh,yl,yh,ny:=-4,5,-10,40,6:
sc:=(xh-xl)/(yh-yl):
plots:-display(
plottools:-transform((x,y,z)->[x,sc*y,z])(Porig),
axis=[tickmarks=[seq(sc*(yl+(yh-yl)*(i-1)/(ny-1))
=evalf(evalf(
yl+(yh-yl)*(i-1)/(ny-1))),i=1..ny)]],
view=[xl..xh,map(`*`,yl..yh,sc),0..4],
scaling=constrained,
orientation=[-15,75,0]);``` The above rescaling is simpler than the general case, as the x-z aspect ratio was already "nice". The z-range and target aspect ratios can also be programmatically incorporated into such code.

There is also no very nice way to programmatically remove that extra white space around that last plot. You can resize and manually zoom, but the only fully programmatic way (of which I know) involves using DocumentTools, eg. DocumentTools:-Layout:-Plot and is more effort.

## escaped locals...

Each separate evaluated call to your procedure f will return a different and distinct escaped local whose name appears as x.

When you enter f()+f() then that is automatically simplified to  2*f()  before the call to f is evaluated. Hence the result of 2*x.  It works differently in the handling of  f()-f()  and two distinct calls to f get evaluated.

## something...

I cannot tell whether you want results in terms of (only) conjugates, or whether abs calls are acceptable.

But one idea might be to yank out the determinant, eg.

or,

## StringTools...

```restart;
with(StringTools,PatternDictionary):

bid:=PatternDictionary:-Create(`builtin`):
words:=[seq(PatternDictionary:-Get(bid,i),
i=1..PatternDictionary:-Size(bid)-1)]:

# and now...

StringTools:-Anagrams("trun",words);

"runt", "turn"

StringTools:-Anagrams("crate",words);

"cater", "carte", "caret", "crate", "trace"```

If you are pulling the candidate string (ie. the first argument passed to the StringTools:-Anagrams command) from somewhere else then you may wish to first hit it with StringTools:-LowerCase so that it matches the strings in the wordlist.

You could also use ospd3 (a large Scrabble dictionary) instead of `builtin`. The former contains some words not present in the builtin list, eg. "gherkins". Or you could construct your own, from some imported collection.

```restart;
with(StringTools,PatternDictionary):

bid:=PatternDictionary:-Create(`builtin`):
words:=[seq(PatternDictionary:-Get(bid,i),
i=1..PatternDictionary:-Size(bid)-1)]:
bidospd3:=PatternDictionary:-Create(`ospd3`):
wordsospd3:=[seq(PatternDictionary:-Get(bidospd3,i),
i=1..PatternDictionary:-Size(bidospd3)-1)]:

#...and now, using an example for comparison

cand := StringTools:-LowerCase("SAET");

cand := "saet"

[{StringTools:-Anagrams(cand,words)}[]];

["east", "seat"]

[{StringTools:-Anagrams(cand,wordsospd3)}[]];

["ates", "east", "eats", "etas", "seat"]
```

## package calls, inside a procedure...

The with command will not work properly from inside a procedure, and that is documented.

with(LinearAlgebra);
inside your procedure GaußKronrodQuadraturKurz, you could utulize the uses facility. You can add this line immediately following the local declarations,
uses LinearAlgebra;

Another choice would be to utilize the full-form name of the package's exported commands, eg.
LinearAlgebra:-Eigenvectors(...);

## don't use print...

Simply call it as,

Tabulate(DF):

without the print. And then the terminating full colon should have the desired effect.

(Your call to print will force the explicit printing of the returned value, which shouldn't be very surprising because that is the purpose of that command...)

## reformulation, and numeric error...

It is worthwhile looking at the result from your call  PDF(0.5*X1 + 0.5*X2, t)  , ie. looking at what is actually being passed to the plot command, and comparing with other formulations.

The integral can be computed (even with floating-point coefficients), without as much numeric damage to coefficients and ensuing evaluation at plotting points. Loss of accuracy can be incurred for this particular example during reformulation (eg. expand, normal, simplify ...).

PDFexpandedpwpoly.mw

I am not suggesting that this is a better workaround that using exact rational coefficients for the original mix. I am trying to address the OP's question about why it can go wrong.

## assume...

@Maria2212 That call to assume inside procedure S is the cause of the problem.

Remove the call to assume that you have inside S. Then try changing the call to dsolve inside S to this, instead:

dsolve(L,y(v)) assuming lambda<0;

Also, inside procedure S there is a problematic call to subs which contains this:

indets(%%,name)

What name are you trying to extract, using that indets call? That method of getting at the name is very poor and error-prone. It should be replaced with something better, but you should tell us first what name it is trying to obtain. Is it one of _C1, _C2, lambda, or v? (Judging from your supplied output I suspect that it is a roundabout attempt to get the name lambda. But perhaps that is another programming mistake.) This should be fixed.

You should also get rid of all the use of % and %% inside S, and instead use assignments to additional local variables.

The use of sum and name-concatenation inside procedure Galerkin is poor. Better would be add and indexed names.

You should not set Digits:=2 . That is a very poor way to display results with only 2 digits, since it incurs significant roundoff error.

## notes...

There are a couple of difficulties. One is the radical in the denominator, which can be handled using evala or rationalize. Another is the conjugated term.

example_symbolic_ac.mw

﻿