Maple Questions and Posts

These are Posts and Questions associated with the product, Maple

Dear experts

how can I numerically plot the following integral and have output as csv file.

in this relation, there is a list of omega1 and omega2 for each k1 and k2. for example,

k1 = [1,2,3,4,5]

omega11=[1,2,3,4,5], omega12= [1,2,3,4,5], omgega21= [1,2,3,4,5],omega22= [1,2,3,4,5],

all other coefficients would be calculated based on k values and corresponding omegas

thanks in advance

Given a pde (or a set of pdes) of multiple funcitons, is there a way to look for solutions when one of the functions is kept arbritary.

For example if I have a set of pdes with f1,f2,f3. Is there a way to see if there is a functional form for f2 and f3 such that the equation is satisfied for any f1?

It appears to me that "simplify/siderels" with two arguments is simply some special "simplify" procedure that just makes use of assumptions, but the results of experiments seem to tell an opposite story. 

simplify(sqrt(x**2), [x = 0]);
                               0

simplify(sqrt(x**2), assume = [x = 0]);
                               0

simplify(ln(exp(x)), [x = 0]);
                               0

simplify(ln(exp(x)), assume = [x = 0]);
                               x

`assuming`(simplify(ln(exp(x))), [x = 0]);
                               x

simplify((1 - cos(x)**2 + sin(x)*cos(x))/(sin(x)*cos(x) + cos(x)**2), [sin(x)**2 + cos(x)**2 = 1]);
                            2                
                  1 - cos(x)  + sin(x) cos(x)
                  ---------------------------
                   cos(x) (cos(x) + sin(x))  

simplify((1 - cos(x)**2 + sin(x)*cos(x))/(sin(x)*cos(x) + cos(x)**2), assume = [sin(x)**2 + cos(x)**2 = 1]);
                             tan(x)

How to explain these? 

(While using Maple 2015 this question concerns any other Maple versions)

I hesitated on the title to write and my first idea was to write "How to modify a built-in functions without making a mess?".
I finally changed my mind in order not to orient the answers in a wrong way.

So this question is about the construction of multi-variate distributions and concerns only the Statistics package.
Here are some of the attributes of a univariate random variable that Maple recognizes, and it is quite normal to expect that the construction of a multi-variate random variable (MVRV for short) distribution should get, at least, some of them.

X := RandomVariable(Normal(a, b)):
map(a -> printf("%a\n", a), [exports(attributes(X)[3])]):
Conditions
ParentName
Parameters
CharacteristicFunction
CDF
CGF
HodgesLehmann
Mean
Median
MGF
Mode
PDF
RousseeuwCrouxSn
StandardDeviation
Support
Variance
CDFNumeric
QuantileNumeric
RandomSample
RandomSampleSetup
RandomVariate
MaximumLikelihoodEstimate

If the distribution is continuous the PDF is fundamental in the sense it enables constructs all the other statistics (=attributes) of a MVRV.
But it is nice to use integrated functions, such as Mean, Support, PDF, and so on, to get the expressions or values of these statistics instead of computing them from this PDF.
Let's that I prefer doing this

MyNormal := proc(m, v)
  description "Reparameterized Normal randomvariable, m=mean, v=variance":
  Distribution(
    PDF = (t -> exp(-1/2*(x-m)^2/v)/sqrt(2*Pi*v))
    , Conditions = [Sigma > 0]
    , Mean =m
  )
end proc:

X := RandomVariable(MyNormal(mu, Sigma)):
Mean(X)
                              m

than doing this

MyNormal := proc(m, v)
  description "Reparameterized Normal randomvariable, m=mean, v=variance":
  Distribution(
    PDF = (t -> exp(-1/2*(x-m)^2/v)/sqrt(2*Pi*v))
    , Conditions = [Sigma > 0]
  )
end proc:

X := RandomVariable(MyNormal(mu, Sigma)):
Mean(X);  # of course undefined
mean := int(PDF(X, x), x=-infinity..+infinity) assuming Sigma > 0
                           undefined
                           mean := 1

So, while all the statistics can be recover from the CDF (provided it exists), it's nicer to define these statistics within the Distribution structure (as in the first construction above).

Now some problems appear when you want to construct the Distribution structure for a MVRV.
The attached file contains the construction of MVRV whose ecah components are mutually independent (to keep the things simple) and both have a Unifom distribution.

MV_Uniform.mw

Here are some observations:

  • Defining a multi-variate PDF goes without problems.
  • Defining the Mean (or many other algebraic or numeric statistics) presents a difficulty related to the type of the arguments the build-in function Mean is aimed to recieve.
    But a workaround, not very elegant, can be found.
  • The case of the Support seems unsolvable: I wasn't able to find any workaround to define the support of a MVRV.
  • I did not consider the Conditions attribute, but I'm not sure that, in the case of, let's say, a bi-gaussian random variable I would be capable to set that the variance is a symmetric positive-definite matrix?

I feel like the main restriction to define such MVRV distributions is the types used in the buid-in functions used in the Distribution structure.

Does anyone have an idea to tackle this problem?

  • Are we doomed to use workarounds like the one I used for defining the MVRV mean?
  • Can we modify the calling sequence of some build-in functions without making a mess and keep them working on build-in distributions?
  • Must we overload the construction of these build-in functions?
    Doing for instance:
    restart:
    with(Statistics):
    local Mean:
    Mean := proc(...) ... end proc

Thanks in advance for any suggestion and help.

 For the command LieAlgebras[RootSpaceDecomposition] I don't understand what the command return, I read the help and see the examples but still not understanding.

 

for example it returns:

RSD := RootSpaceDecomposition(CSA);

RSD := table([[-2, -1] = E31, [2, 1] = E13, [1, 2] = E23, [1, -1] = E12, [-1, 1] = E21, [-1, -2] = E32])

I don't understand what means [-2, -1] even they said that is the root but I know that a root is in h* so it must be only a number not a vector.

eta/2*diff(f, eta) + (diff(f, eta, eta) + lambda*diff(g, eta)/g)*g^lambda=0,

g*diff(g, eta) + sigma*eta*diff(g, eta)=0,here lambda,sigma positive parameters initial conditions f'(0)=0,g'(0)=0,f' goes to 1 at eta goes to infinity ,g' goes to 1 at eta goes to infinity

can you solve this coupled nonlinear differential equations  using homotopy analysis method

I cannot view 3d graphics with my version of Ubuntu 20.04. I've updated all my computer's graphics card drivers and the problem persists. If I run without hardware acceleration, nothing changes; no visualization and no production possible.
Do you have any ideas for solving this problem? Maple uses OpenGL libraries for 3D production and visualization, and these libraries are installed on my computer. Would installing mesa solve the problem, for example?

Thanks in advance.

A very easy example: 

solve({a > 0, ln(a) + ln(1 + a) >= 0}, a);
 = 
                        /    2         \ 
                        |---------- < a| 
                       <  (1/2)         >
                        |5      + 1    | 
                        \              / 

The output is . Clearly,  is also a solution to this inequality, so the solutions have been lost. But there is no warning message. Why?

Dear all

I would like if we can solve a problem of a  field homomorphism between two fields. 

Please, if possible a code that help me to define the homomorphism between finite field. 

homomorphism_finite_field.mw

Thank you 

I found another big problem. 

In 2022, I get Error, (in SolveTools:-Polynomial) too many levels of recursion when using alias(seq(c[k] = _C||k, k = 0..10)); at the top of my code and the solution to the ode has c[2],c[3] etc.. as constants of integration.

This solution was given by Kitonum in this answer

I had this for years in the code (i.e. the alias) code.

In Maple 2022 the following gives the above exception error from odetest. I am using 2022, because in 2023 it just hangs on the same code.  

If I remove the alias code, no error shows. (solution is wrong, but that is different story). 

Also, If I remove the alias code, 2023 no long hangs!  

restart;

#kernelopts('assertlevel'=2):

alias(seq(c[k] = _C||k, k = 0..10));

c[0], c[1], c[2], c[3], c[4], c[5], c[6], c[7], c[8], c[9], c[10]

ode:=diff(diff(y(x),x),x)+4*diff(y(x),x)+12*y(x) = 80*sin(2*x);
sol:=y(x) = -10*exp(-1/2*2^(1/2)*arctan(sin(2*2^(1/2)*x)/cos(2*2^(1/2)*x)))*(c[3]*cos(2*2^(1/2)*x)+c[2]*sin(2*2^(1/2)*x))*((-1/10*exp(4*I*2^(1/2)*x)+1/10)*c[2]-1/10*I*c[3]*exp(4*I*2^(1/2)*x)-1/10*I*c[3])/((exp(4*I*2^(1/2)*x)-1)*c[2]+I*c[3]*exp(4*I*2^(1/2)*x)+I*c[3])*c[1]-10*I*exp(-1/2*2^(1/2)*arctan(sin(2*2^(1/2)*x)/cos(2*2^(1/2)*x)))*(c[3]*cos(2*2^(1/2)*x)+c[2]*sin(2*2^(1/2)*x))*(cos(2*x)-sin(2*x))*exp(2*I*2^(1/2)*x)*exp(2*I*2^(1/2)*x)^(-1/2*I*2^(1/2))/((exp(4*I*2^(1/2)*x)-1)*c[2]+I*c[3]*exp(4*I*2^(1/2)*x)+I*c[3]);

odetest(sol,ode);

diff(diff(y(x), x), x)+4*(diff(y(x), x))+12*y(x) = 80*sin(2*x)

y(x) = -10*exp(-(1/2)*2^(1/2)*arctan(sin(2*2^(1/2)*x)/cos(2*2^(1/2)*x)))*(c[3]*cos(2*2^(1/2)*x)+c[2]*sin(2*2^(1/2)*x))*((-(1/10)*exp((4*I)*2^(1/2)*x)+1/10)*c[2]-((1/10)*I)*c[3]*exp((4*I)*2^(1/2)*x)-((1/10)*I)*c[3])*c[1]/((exp((4*I)*2^(1/2)*x)-1)*c[2]+I*c[3]*exp((4*I)*2^(1/2)*x)+I*c[3])-(10*I)*exp(-(1/2)*2^(1/2)*arctan(sin(2*2^(1/2)*x)/cos(2*2^(1/2)*x)))*(c[3]*cos(2*2^(1/2)*x)+c[2]*sin(2*2^(1/2)*x))*(cos(2*x)-sin(2*x))*exp((2*I)*2^(1/2)*x)*(exp((2*I)*2^(1/2)*x))^(-((1/2)*I)*2^(1/2))/((exp((4*I)*2^(1/2)*x)-1)*c[2]+I*c[3]*exp((4*I)*2^(1/2)*x)+I*c[3])

Error, (in SolveTools:-Polynomial) too many levels of recursion

 

Download odetest_error_june_15_2023.mw

Why using the alias line above causes this error?  If you remove the alias line, you will see it completes with no error.

But it hangs in 2023. I am no longer using 2023 but went back to 2022 due to too many hangs in 2023. I wonder now if it because of this alias line I had there all the time.  I will remove now and see if this solves some of the hangs I had in 2023.

Anyone can shed some more light on what is going on?

Windows 10.

 Hi,
  Can someone guide me on how to write the code for this program?

sh.docx

I want to plot phase portrait for 4 dimension. help me

4_equation_phase_portrait.mw

thank youvery much

There are some excel sheets that cannot be read w/ Maple 2023 that were read w/ Maple 2022.
I attached an example.ImportFerrites.maple
Does anyone else experienced this issue or has a solution?

Can you change f(eta) to upflow curve and theta(eta) to downflow curve.

In my Problem,Boundary Conditions are

theta(infinity) = 0, (D(f))(infinity) = 1 , (Take, eta =infinity)

Flows will be correct for what value is taken for infinity .

I take  eta = 5. and also tried changing ranges  but could't find it.Please Help to fix the curve.

my code is,

SM.mw

In a recent question/conversation, I had discussed integrating dsolve/numeric-based codes with NLPsolve at 
https://www.mapleprimes.com/questions/236494-Inconsistent-Behavior-With-Dsolvenumeric-And-NLPSolve-sqp

@acer was able to help resolve the issue by calling NLPSolve with higher optimality tolerance.

I am starting a new question/post to show the need for evalhf as the previous post takes too long to load.
Code is attached for the same.
 

restart:

currentdir();

"C:\Users\Venkat\OneDrive - The University of Texas at Austin\Documents\trial\Learnningsqpoptim"

Test code written by Dr. Venkat Subramanian at UT Austin, 05/31/2023 and has been updated multiple times. This code uses CVP approach (piecwise constant) to perform optimal control. NLPSolve combination with RK2 approach to integrate ODEs (constant step size) works well for this problem. But dsolve numeric based codes work (with optimality tolerance fix from acer), but cannot handle large values of nvar (optimization variables).

Digits:=15;

15

 

eqodes:=[diff(ca(t),t)=-(u+u^2/2)*1.0*ca(t),diff(cb(t),t)=1.0*u*ca(t)-0.1*cb(t)];

[diff(ca(t), t) = -1.0*(u+(1/2)*u^2)*ca(t), diff(cb(t), t) = 1.0*u*ca(t)-.1*cb(t)]

soln:=dsolve({op(eqodes),ca(0)=alpha,cb(0)=beta},type=numeric,'parameters'=[alpha,beta,u],savebinary=true):

 

 

Note that Vector or Array form can be used for RK2h to implement the procedure for any number of variables/ODEs. But the challenge will be when implicit methods are used to integrate ODEs/DAEs and running them in evalhf form (this can be done with Gauss elimination type linear solver, but this will be limited to small number of ODEs/DAEs say 200 or so).

RK2h:=proc(NN,u,dt,y0::Vector(datatype=float[8]))#this procedure can be made efficient in vector form
local j,c1mid,c2mid,c10,c20,c1,c2;
option hfloat;
c10:=y0[1];c20:=y0[2];
for j from 1 to NN do
  c1mid:=c10+dt/NN*(-(u+u^2/2)*c10):
  c2mid:=c20+dt/NN*(u*c10)-0.1*dt/NN*c20:
  c1:=c10/2+c1mid/2+dt/NN/2*(-(u+u^2/2)*c1mid):
  c2:=c20/2+c2mid/2+dt/NN/2*(u*c1mid)-0.1*dt/NN/2*c2mid:
  c10:=c1:c20:=c2:
  od:
y0[1]:=c10:y0[2]:=c20:
end proc:

 

soln('parameters'=[1,0,0.1]);soln(0.1);

[alpha = 1., beta = 0., u = .1]

[t = .1, ca(t) = HFloat(0.9895549324188543), cb(t) = HFloat(0.009898024276129616)]

 

 

ssdsolve:=proc(x)
interface(warnlevel=0):

#if  type(x,Vector)#if type is not needed for this problem, might be needed for other problems
#then


local z1,n1,i,c10,c20,dt,u;
global soln,nvar;
dt:=evalf(1.0/nvar):
c10:=1.0:c20:=0.0:
for i from 1 to nvar do
u:=x[i]:
soln('parameters'=[c10,c20,u]):
z1:=soln(dt):
c10:=subs(z1,ca(t)):c20:=subs(z1,cb(t)):
od:
-c20;
 #else 'procname'(args):

#end if:

end proc:

 

 

ssRK2:=proc(x)#based on RK2
#interface(warnlevel=0):
option hfloat;
#if  type(x,Vector)
#then

local z1,n1,i,c10,c20,dt,u,NN,y0;
global nvar,RK2h;
y0:=Array(1..2,[1.0,0.0],datatype=float[8]):
#y0[1]:=1.0:y0[2]:=0.0:
dt:=evalf(1.0/nvar):
NN:=256*2/nvar;#NN is hardcode based on the assumption that nvar will be a multiple of 2 <=256


for i from 1 to nvar do
u:=x[i]:
evalhf(RK2h(NN,u,dt,y0)):
od:
-y0[2];
 #else 'procname'(args):

#end if:

end proc:

nvar:=2;

2

ic0:=Vector(nvar,[seq(0.1,i=1..nvar)],datatype=float):

bl := Vector(nvar,[seq(0.,i=1..nvar)],datatype=float):bu := Vector(nvar,[seq(5.,i=1..nvar)],datatype=float):

infolevel[Optimization]:=15;

15

CodeTools:-Usage(Optimization:-NLPSolve(nvar,evalhf(ssRK2),[],initialpoint=ic0,[bl,bu],optimalitytolerance=1e-6)):

NLPSolve: calling NLP solver

NLPSolve: using method=sqp

NLPSolve: number of problem variables 2

NLPSolve: number of nonlinear inequality constraints 0

NLPSolve: number of nonlinear equality constraints 0

NLPSolve: number of general linear constraints 0

NLPSolve: feasibility tolerance set to 0.1053671213e-7

NLPSolve: optimality tolerance set to 0.1e-5

NLPSolve: iteration limit set to 50

NLPSolve: infinite bound set to 0.10e21

NLPSolve: trying evalhf mode

NLPSolve: trying evalf mode

attemptsolution: number of major iterations taken 11

memory used=0.96MiB, alloc change=0 bytes, cpu time=16.00ms, real time=121.00ms, gc time=0ns

 

Calling the procedures based on RK2 with NLPSolve uses evalf for numerical gradient (fdiff). Providing a procedure for gradient solves this issue.

gradRK2 := proc (x,g)
option hfloat;
local base, i, xnew, del;
global nvar,ssRK2;
xnew:=Array(1..nvar,datatype=float[8]):
base := ssRK2(x);
for i to nvar do
xnew[i] := x[i]; end do;
for i to nvar do
del := max(1e-5,.1e-6*x[i]);
xnew[i] := xnew[i]+del;
g[i] := (ssRK2(xnew)-base)/del;
xnew[i] := xnew[i]-del;
end do;
end proc:

CodeTools:-Usage(Optimization:-NLPSolve(nvar,evalhf(ssRK2),[],initialpoint=ic0,[bl,bu],objectivegradient=gradRK2,optimalitytolerance=1e-6)):

NLPSolve: calling NLP solver

NLPSolve: using method=sqp

NLPSolve: number of problem variables 2

NLPSolve: number of nonlinear inequality constraints 0

NLPSolve: number of nonlinear equality constraints 0

NLPSolve: number of general linear constraints 0

NLPSolve: feasibility tolerance set to 0.1053671213e-7

NLPSolve: optimality tolerance set to 0.1e-5

NLPSolve: iteration limit set to 50

NLPSolve: infinite bound set to 0.10e21

NLPSolve: trying evalhf mode

attemptsolution: number of major iterations taken 11

memory used=184.77KiB, alloc change=0 bytes, cpu time=31.00ms, real time=110.00ms, gc time=0ns

Significant saving in memory used is seen. CPU time is also less, which is more apparent at larger valeus of nvar.

 

dsolvenumeric based codes work for optimization, but evalf is invoked possibly for both the objective and gradient

CodeTools:-Usage(Optimization:-NLPSolve(nvar,(ssdsolve),[],initialpoint=ic0,[bl,bu],optimalitytolerance=1e-6)):

NLPSolve: calling NLP solver

NLPSolve: using method=sqp

NLPSolve: number of problem variables 2

NLPSolve: number of nonlinear inequality constraints 0

NLPSolve: number of nonlinear equality constraints 0

NLPSolve: number of general linear constraints 0

NLPSolve: feasibility tolerance set to 0.1053671213e-7

NLPSolve: optimality tolerance set to 0.1e-5

NLPSolve: iteration limit set to 50

NLPSolve: infinite bound set to 0.10e21

NLPSolve: trying evalhf mode

NLPSolve: trying evalf mode

attemptsolution: number of major iterations taken 11

memory used=14.61MiB, alloc change=37.00MiB, cpu time=94.00ms, real time=213.00ms, gc time=46.88ms

 

Providing gradient procedure doesn't help with evalhf computaiton for dsolve/numeric based procedures.

graddsolve := proc (x,g)
local base, i, xnew, del;
global nvar,ssdsolve;
#xnew:=Vector(nvar,datatype=float[8]):
xnew:=Array(1..nvar,datatype=float[8]):
base := ssdsolve(x);
for i to nvar do
xnew[i] := x[i]; end do;
for i to nvar do
del := max(1e-5,.1e-6*x[i]);
xnew[i] := xnew[i]+del;
g[i] := (ssdsolve(xnew)-base)/del;
xnew[i] := xnew[i]-del;
end do;
end proc:

CodeTools:-Usage(Optimization:-NLPSolve(nvar,(ssdsolve),[],initialpoint=ic0,[bl,bu],objectivegradient=graddsolve,optimalitytolerance=1e-6)):

NLPSolve: calling NLP solver

NLPSolve: using method=sqp

NLPSolve: number of problem variables 2

NLPSolve: number of nonlinear inequality constraints 0

NLPSolve: number of nonlinear equality constraints 0

NLPSolve: number of general linear constraints 0

NLPSolve: feasibility tolerance set to 0.1053671213e-7

NLPSolve: optimality tolerance set to 0.1e-5

NLPSolve: iteration limit set to 50

NLPSolve: infinite bound set to 0.10e21

NLPSolve: trying evalhf mode

NLPSolve: trying evalf mode

attemptsolution: number of major iterations taken 11

memory used=3.82MiB, alloc change=0 bytes, cpu time=31.00ms, real time=129.00ms, gc time=0ns

 

Calling both RK2 and dsolve based procedures again to check the values and computation meterics

s1RK2:=CodeTools:-Usage(Optimization:-NLPSolve(nvar,evalhf(ssRK2),[],initialpoint=ic0,[bl,bu],objectivegradient=gradRK2,optimalitytolerance=1e-6)):

NLPSolve: calling NLP solver

NLPSolve: using method=sqp

NLPSolve: number of problem variables 2

NLPSolve: number of nonlinear inequality constraints 0

NLPSolve: number of nonlinear equality constraints 0

NLPSolve: number of general linear constraints 0

NLPSolve: feasibility tolerance set to 0.1053671213e-7

NLPSolve: optimality tolerance set to 0.1e-5

NLPSolve: iteration limit set to 50

NLPSolve: infinite bound set to 0.10e21

NLPSolve: trying evalhf mode

attemptsolution: number of major iterations taken 11

memory used=185.09KiB, alloc change=0 bytes, cpu time=16.00ms, real time=95.00ms, gc time=0ns

s1dsolve:=CodeTools:-Usage(Optimization:-NLPSolve(nvar,ssdsolve,[],initialpoint=ic0,[bl,bu],objectivegradient=graddsolve,optimalitytolerance=1e-6)):

NLPSolve: calling NLP solver

NLPSolve: using method=sqp

NLPSolve: number of problem variables 2

NLPSolve: number of nonlinear inequality constraints 0

NLPSolve: number of nonlinear equality constraints 0

NLPSolve: number of general linear constraints 0

NLPSolve: feasibility tolerance set to 0.1053671213e-7

NLPSolve: optimality tolerance set to 0.1e-5

NLPSolve: iteration limit set to 50

NLPSolve: infinite bound set to 0.10e21

NLPSolve: trying evalhf mode

NLPSolve: trying evalf mode

attemptsolution: number of major iterations taken 11

memory used=3.82MiB, alloc change=0 bytes, cpu time=31.00ms, real time=122.00ms, gc time=0ns

s1RK2[1];s1dsolve[1];

-.523900304316377463

-.523901163953022553

 

Next, a for loop is written to optimize for increasing values of nvar. One can see that evalhf is important for larger values of nvar. While dsolve/numeric is a superior code, not being able to use it in evalhf format is a significant weakness and should be addressed. Note that dsolve numeric evalutes procedures that are evaluated in evalhf or compiled form, so hopefully this is an easy fix.

infolevel[Optimization]:=0:

for j from 1 to 9 do
nvar:=2^(j-1):
ic0:=Vector(nvar,[seq(0.1,i=1..nvar)],datatype=float):
bl := Vector(nvar,[seq(0.,i=1..nvar)],datatype=float):bu := Vector(nvar,[seq(5.,i=1..nvar)],datatype=float):
soptRK[j]:=CodeTools:-Usage(Optimization:-NLPSolve(nvar,evalhf(ssRK2),[],initialpoint=ic0,[bl,bu],objectivegradient=gradRK2,optimalitytolerance=1e-6)):
print(2^(j-1),soptRK[j][1]);

od:

memory used=88.40KiB, alloc change=0 bytes, cpu time=0ns, real time=10.00ms, gc time=0ns

1, HFloat(-0.5008626988793192)

memory used=155.66KiB, alloc change=0 bytes, cpu time=16.00ms, real time=30.00ms, gc time=0ns

2, -.523900304316377463

memory used=175.36KiB, alloc change=0 bytes, cpu time=62.00ms, real time=80.00ms, gc time=0ns

4, -.535152497919956782

memory used=243.69KiB, alloc change=0 bytes, cpu time=156.00ms, real time=239.00ms, gc time=0ns

8, -.540546896131004706

memory used=260.09KiB, alloc change=0 bytes, cpu time=141.00ms, real time=260.00ms, gc time=0ns

16, -.542695734426874465

memory used=385.84KiB, alloc change=0 bytes, cpu time=313.00ms, real time=545.00ms, gc time=0ns

32, -.542932877726400531

memory used=0.65MiB, alloc change=0 bytes, cpu time=734.00ms, real time=1.18s, gc time=0ns

64, -.543011976841572652

memory used=1.24MiB, alloc change=0 bytes, cpu time=1.55s, real time=2.40s, gc time=0ns

128, -.543035276649319276

memory used=2.99MiB, alloc change=0 bytes, cpu time=3.45s, real time=5.92s, gc time=0ns

256, -.543041496228812814

for j from 1 to 6 do
nvar:=2^(j-1):
ic0:=Vector(nvar,[seq(0.1,i=1..nvar)],datatype=float):
bl := Vector(nvar,[seq(0.,i=1..nvar)],datatype=float):bu := Vector(nvar,[seq(5.,i=1..nvar)],datatype=float):
soptdsolve[j]:=CodeTools:-Usage(Optimization:-NLPSolve(nvar,evalf(ssdsolve),[],initialpoint=ic0,[bl,bu],objectivegradient=graddsolve,optimalitytolerance=1e-6)):
print(2^(j-1),soptdsolve[j][1]);

od:

memory used=0.66MiB, alloc change=0 bytes, cpu time=16.00ms, real time=11.00ms, gc time=0ns

1, HFloat(-0.5008631947224957)

memory used=3.79MiB, alloc change=0 bytes, cpu time=15.00ms, real time=52.00ms, gc time=0ns

2, -.523901163953022553

memory used=21.28MiB, alloc change=0 bytes, cpu time=78.00ms, real time=236.00ms, gc time=0ns

4, -.535153942647626391

memory used=144.82MiB, alloc change=-4.00MiB, cpu time=469.00ms, real time=1.60s, gc time=46.88ms

8, -.540549407239521607

memory used=347.30MiB, alloc change=16.00MiB, cpu time=1.27s, real time=3.45s, gc time=140.62ms

16, -.542699055038344258

memory used=1.33GiB, alloc change=-8.00MiB, cpu time=7.66s, real time=15.92s, gc time=750.00ms

32, -.542936165630524603

 

SS:=[seq(soptRK[j][1],j=1..9)];

[HFloat(-0.5008626988793192), -.523900304316377463, -.535152497919956782, -.540546896131004706, -.542695734426874465, -.542932877726400531, -.543011976841572652, -.543035276649319276, -.543041496228812814]

 

E1:=[seq(SS[i]-SS[i+1],i=1..nops(SS)-1)];

[HFloat(0.02303760543705824), 0.11252193603580e-1, 0.5394398211048e-2, 0.2148838295869e-2, 0.237143299527e-3, 0.79099115172e-4, 0.23299807746e-4, 0.6219579494e-5]

To get 6 Digits accuracy we need nvar = 256 which may not be attainable with dsolve/numeric approach unless we are able to call it evalhf format.

 

 


 

Download ImportanceofevalhfNLPSolve.mw

First 37 38 39 40 41 42 43 Last Page 39 of 2218