slizovskiy

25 Reputation

8 Badges

7 years, 151 days

MaplePrimes Activity


These are replies submitted by slizovskiy

since I am in a hurry to see the effects, i've written a rude step interpolation function myself and continue the computation...

But it really looks as I use golden watch as a hammer. 

since I am in a hurry to see the effects, i've written a rude step interpolation function myself and continue the computation...

But it really looks as I use golden watch as a hammer. 

Hello Acer!

Thank you for your reply!   I attach the example, it takes around 5 minutes to generate the data (first integral is done analytically, but the second is only numerical). 

I found that increasing the data size 4 times, makes ArrayInterpolation work around 16 times slower!

About BSpline: I have NOT found in help that it works for multi-dimensional case.

You see in the example, I use many tricks I learned from You, how to work with numerics...

  Maybe there is a faster way to pass arguments to ArrayInterpolation?  I cannot debug it and see what eats time...  
Best wishes,

Sergey

 

 

unassign('Lambda', 'k'); i1 := `assuming`([simplify(convert(int(k^2*sin(theta)/((q^2+2*q*k*cos(theta))*(1/2)+((k^2+q^2+2*q*k*cos(theta))^2-k^4)/(2*Lambda^2)-I*nu), theta = 0 .. Pi), ln))], [m > 0, k > 0, q > 0])

Lambda^2*k*(-ln((4*q*k+Lambda^2+2*q^2+2*k^2+(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))/(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))+ln(-(-(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2)+4*q*k+Lambda^2+2*q^2+2*k^2)/(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))+ln((Lambda^2+2*q^2+2*k^2-4*q*k+(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))/(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))-ln(-(-(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2)+Lambda^2+2*q^2+2*k^2-4*q*k)/(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2)))/(q*(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))

(1.1)

i11 := 2*Re(i1)

2*Re(Lambda^2*k*(-ln((4*q*k+Lambda^2+2*q^2+2*k^2+(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))/(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))+ln(-(-(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2)+4*q*k+Lambda^2+2*q^2+2*k^2)/(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))+ln((Lambda^2+2*q^2+2*k^2-4*q*k+(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))/(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))-ln(-(-(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2)+Lambda^2+2*q^2+2*k^2-4*q*k)/(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2)))/(q*(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2)))

(1.2)

NULL

 

 

i12 := codegen[makeproc]([codegen[optimize](evalf(i11/(2*Pi)^2))], [k, kF2, q, nu, Lambda])

unassign('kF1', 'q', 'omega', 'vF1'); Chi12unint := unapply(('i12')(k, kF2, q, omega, Lambda), k, kF2, q, omega, Lambda, numeric)

proc (k, kF2, q, omega, Lambda) local unnamed; if type(evalf(k), 'numeric') and type(evalf(kF2), 'numeric') and type(evalf(q), 'numeric') and type(evalf(omega), 'numeric') and type(evalf(Lambda), 'numeric') then i12(k, kF2, q, omega, Lambda) elif procname <> 'unknown' and not member(sprintf("%a", procname), {"%", "%%", "%%%"}) then ('procname')(k, kF2, q, omega, Lambda) else unnamed := pointto((Array(1..1, {(1) = 18446884315406287902}))[1]); ('unnamed')(k, kF2, q, omega, Lambda) end if end proc

(1)

Chi12 := proc (kF2, q, omega, Lambda) options operator, arrow; evalf(Int(Chi12unint(k, kF2, q, omega, Lambda), k = 0 .. kF2, epsilon = 1/100, method = _d01ajc)) end proc

proc (kF2, q, omega, Lambda) options operator, arrow; evalf(Int(Chi12unint(k, kF2, q, omega, Lambda), k = 0 .. kF2, epsilon = 1/100, method = _d01ajc)) end proc

(2)

``

Warning,  computation interrupted

 

``

NULL

 

chiTable := [seq([seq([seq(Chi12(kF2, q, omega, .4), omega = 0.1e-2 .. 1, 0.1e-1)], q = 0.1e-2 .. 1, 0.5e-1)], kF2 = 0.5e-1 .. .5, .1)]

chiTableX := [[seq(kF2, kF2 = 0.5e-1 .. .5, .1)], [seq(q, q = 0.1e-2 .. 1, 0.5e-1)], [seq(omega, omega = 0.1e-2 .. 1, 0.1e-1)]]

with(CurveFitting)

``

ArrayTools[Dimensions](chiTable)

[1 .. 5, 1 .. 20, 1 .. 100]

(3)

ArrayTools[Dimensions](chiTableX[3])

[1 .. 100]

(4)

chiTable1 := Array(chiTable)

time(ArrayInterpolation(chiTableX, chiTable1, Array([[.1, .1, .1]]), verify = false, uniform = true, degree = 2)[1])

0.8e-2

(5)

ChiInterp := proc (kF2, q, omega) ArrayInterpolation(chiTableX, chiTable, Array([[kF2, q, omega]]), verify = false, uniform = true, degree = 2)[1] end proc

proc (kF2, q, omega) CurveFitting:-ArrayInterpolation(chiTableX, chiTable, Array([[kF2, q, omega]]), verify = false, uniform = true, degree = 2)[1] end proc

(6)

``


Download ArrayInterpolationEx.mwArrayInterpolationEx.mw

Hello Acer!

Thank you for your reply!   I attach the example, it takes around 5 minutes to generate the data (first integral is done analytically, but the second is only numerical). 

I found that increasing the data size 4 times, makes ArrayInterpolation work around 16 times slower!

About BSpline: I have NOT found in help that it works for multi-dimensional case.

You see in the example, I use many tricks I learned from You, how to work with numerics...

  Maybe there is a faster way to pass arguments to ArrayInterpolation?  I cannot debug it and see what eats time...  
Best wishes,

Sergey

 

 

unassign('Lambda', 'k'); i1 := `assuming`([simplify(convert(int(k^2*sin(theta)/((q^2+2*q*k*cos(theta))*(1/2)+((k^2+q^2+2*q*k*cos(theta))^2-k^4)/(2*Lambda^2)-I*nu), theta = 0 .. Pi), ln))], [m > 0, k > 0, q > 0])

Lambda^2*k*(-ln((4*q*k+Lambda^2+2*q^2+2*k^2+(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))/(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))+ln(-(-(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2)+4*q*k+Lambda^2+2*q^2+2*k^2)/(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))+ln((Lambda^2+2*q^2+2*k^2-4*q*k+(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))/(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))-ln(-(-(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2)+Lambda^2+2*q^2+2*k^2-4*q*k)/(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2)))/(q*(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))

(1.1)

i11 := 2*Re(i1)

2*Re(Lambda^2*k*(-ln((4*q*k+Lambda^2+2*q^2+2*k^2+(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))/(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))+ln(-(-(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2)+4*q*k+Lambda^2+2*q^2+2*k^2)/(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))+ln((Lambda^2+2*q^2+2*k^2-4*q*k+(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))/(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2))-ln(-(-(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2)+Lambda^2+2*q^2+2*k^2-4*q*k)/(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2)))/(q*(Lambda^4+4*Lambda^2*k^2+4*k^4+(8*I)*nu*Lambda^2)^(1/2)))

(1.2)

NULL

 

 

i12 := codegen[makeproc]([codegen[optimize](evalf(i11/(2*Pi)^2))], [k, kF2, q, nu, Lambda])

unassign('kF1', 'q', 'omega', 'vF1'); Chi12unint := unapply(('i12')(k, kF2, q, omega, Lambda), k, kF2, q, omega, Lambda, numeric)

proc (k, kF2, q, omega, Lambda) local unnamed; if type(evalf(k), 'numeric') and type(evalf(kF2), 'numeric') and type(evalf(q), 'numeric') and type(evalf(omega), 'numeric') and type(evalf(Lambda), 'numeric') then i12(k, kF2, q, omega, Lambda) elif procname <> 'unknown' and not member(sprintf("%a", procname), {"%", "%%", "%%%"}) then ('procname')(k, kF2, q, omega, Lambda) else unnamed := pointto((Array(1..1, {(1) = 18446884315406287902}))[1]); ('unnamed')(k, kF2, q, omega, Lambda) end if end proc

(1)

Chi12 := proc (kF2, q, omega, Lambda) options operator, arrow; evalf(Int(Chi12unint(k, kF2, q, omega, Lambda), k = 0 .. kF2, epsilon = 1/100, method = _d01ajc)) end proc

proc (kF2, q, omega, Lambda) options operator, arrow; evalf(Int(Chi12unint(k, kF2, q, omega, Lambda), k = 0 .. kF2, epsilon = 1/100, method = _d01ajc)) end proc

(2)

``

Warning,  computation interrupted

 

``

NULL

 

chiTable := [seq([seq([seq(Chi12(kF2, q, omega, .4), omega = 0.1e-2 .. 1, 0.1e-1)], q = 0.1e-2 .. 1, 0.5e-1)], kF2 = 0.5e-1 .. .5, .1)]

chiTableX := [[seq(kF2, kF2 = 0.5e-1 .. .5, .1)], [seq(q, q = 0.1e-2 .. 1, 0.5e-1)], [seq(omega, omega = 0.1e-2 .. 1, 0.1e-1)]]

with(CurveFitting)

``

ArrayTools[Dimensions](chiTable)

[1 .. 5, 1 .. 20, 1 .. 100]

(3)

ArrayTools[Dimensions](chiTableX[3])

[1 .. 100]

(4)

chiTable1 := Array(chiTable)

time(ArrayInterpolation(chiTableX, chiTable1, Array([[.1, .1, .1]]), verify = false, uniform = true, degree = 2)[1])

0.8e-2

(5)

ChiInterp := proc (kF2, q, omega) ArrayInterpolation(chiTableX, chiTable, Array([[kF2, q, omega]]), verify = false, uniform = true, degree = 2)[1] end proc

proc (kF2, q, omega) CurveFitting:-ArrayInterpolation(chiTableX, chiTable, Array([[kF2, q, omega]]), verify = false, uniform = true, degree = 2)[1] end proc

(6)

``


Download ArrayInterpolationEx.mwArrayInterpolationEx.mw

@karahosen 

pow is known in older versions, but it's unknown in Maple 16.

evalf(pow(1, 1))  returns unevaluated. 

The working compiled procedure seem to involve conversionof double  to M_INT type in some places...  Can it be dangerous???  I mean, is it rounding to integer?

The working compiled procedure seem to involve conversionof double  to M_INT type in some places...  Can it be dangerous???  I mean, is it rounding to integer?

Thank you, Acer.  That was helpful!

Thank you, Acer.  That was helpful!

@PatrickT 

Comments:  virgin installation of Ubuntu 11.10.  Maple 15 works ok, Maple 16 gives a bit stange plots...

 

Formally, only up to Ubuntu 11.10  is supported.  But it's strange anyway.

In the meantime I have tried Maple 16,  it computes integrals better, but still some part of my old worksheet return unevaluated.

I will also try the same under Windows.

Formally, only up to Ubuntu 11.10  is supported.  But it's strange anyway.

In the meantime I have tried Maple 16,  it computes integrals better, but still some part of my old worksheet return unevaluated.

I will also try the same under Windows.

@Axel Vogt 

Thank you for help, Axel!

  Anyhow, it may be better to combine Maple with Matlab if I want good numerical efficiency. 

@Axel Vogt 

Thank you for help, Axel!

  Anyhow, it may be better to combine Maple with Matlab if I want good numerical efficiency. 

 

@Axel Vogt 

Dear Axel, thank you for trying to help.  I again posted an oversimplified problem. Here is a fuller version:

This works
> Compiler:-Compile(proc (U, kF1, kF2, m, k, omega, q) ((-1)*2.*I*omega*m+k^2+q^2+2.*q*k)/sqrt(6.*k^2*q^2+4.*q*k^3+4.*q^3*k+k^4+q^4+4.*omega^2*m^2) end proc, optimize = true);
 
This does not:
> prc := codegen[makeproc]((-(2.*I)*omega*m+k^2+q^2+2.*q*k)/sqrt(6.*k^2*q^2+4.*q*k^3+4.*q^3*k+k^4+q^4+4.*omega^2*m^2), [U, kF1, kF2, m, k, omega, q]);
> Compiler:-Compile(prc);
Error, (in Print) rational numbers and arithmetic are not yet supported
 
You see that any innocent shuffling of expression with sqrt  (like makeproc) automatically converts it to the form that is not accepted by Compile. 
I have reported a similar problem of this type  here  http://www.mapleprimes.com/questions/135011-Compiling-Errors
It seems like one needs some Wizard tricks to make Compile work!
1 2 Page 1 of 2