@torabi The differential order of the system is seven, and you have one parameter, so you need eight boundary conditions. There are six in bcs, therefore two more must come from extra_bcs.

I took out approxsoln because it's possible to get a solution (in fact, 25 solutions) without it. Note that approxsoln doesn't change the solution at all; it only helps the problem get started. And I guessed that it was just something leftover from a previous problem. Feel free to put it back; it won't make a big difference.

The 25 numbers on the second-to-last line are your values of omega2. You didn't say before that you needed them sorted. I added code to do that. The final line are the indices into res that correspond to the sorted values of omega2. In other words, the last line shows the value of b corresponding to each value of omega2.

warning-3.mw

@Adam Ledger Your request in the original post was

im just curious to know how this one measures up in terms of computational efficiency...if anyone has the time to give it a try and let me know what they think ie faster more logical way about it any feed back is appreciated....

And that's what I did, but you don't seem to appreciate it; rather, you seem to take offense. I showed you how your function could be transcribed into more-normal Maple, still containing essentially a delta function. Then I rewrote it in a much-simpler, more-logical, more-efficient way without a delta function.

I didn't say that the delta function is worthless. I said that it had no practical computational value, meaning that you shouldn't use it in code, and that it had no theoretical computational value, meaning that you shouldn't use it in pseudocode either. The point is to think in terms of control structures, like if statements, which are more logical and less arithmetical.

This is covered on this Wikipedia page.

@Adam Ledger I recommend the book Prime Numbers: A Computational Perspective by Richard Crandall and Carl Pomerance. The algorithms are in pseudocode.

You should narrow your study to one of the following problems:

1. Prime detection: detecting whether an integer is prime
2. Prime counting: counting all the primes less than a certain number
3. Prime listing: listing all the primes less than a certain number

To see the code of a Maple procedure, use showstat:

showstat(isprime);

The command select is in the kernel, so you can't see it. However, it is a triviality which adds no significant time to the command. It is equivalent to

Select:= (B, L::{list,set})->
     `if`(L::list, `[]`, `{}`)(seq(`if`(B(x, _rest), x, NULL), x= L))
:

but it runs faster because it's in the kernel.

To time a Maple procedure and save the time, use CodeTools:-Usage with the output option:

T[n]:= CodeTools:-Usage(P(n), output= realtime);

You can replace realtime with cputime. This'll give you a different timing. Which one is more meaningful is debatable, but it's how they grow with n (asymptotic time complexity) that matters most. That's unlikely to change if you switch between realtime and cputime. Make sure that n is large enough to give you a meaningful timing. Times of 16 ms or less for Maple code are usually meaningless, or at least difficult to interpret correctly.

You should make an empirical study of the asymptotic time complexity of your algorithms. The first step to this is to make a log-log plot of n vs. T[n] using T[n] as defined above. Looking at the ratios of the times of two different algorithms before knowing their time complexity is not very meaningful. Ratios would be more meaningful if you were comparing different implementations of the same algorithm.

The delta function is a convenience for mathematical notation and has little to no practical or even theoretical computational value. You should try to erase it from your mind and instead think in terms of if and select statements and other control structures.

The Legendre prime-counting formula that you gave is mostly of theoretical interest; it has little practical computational value. Read the linked MathWorld page, and follow the links to Meissel's formula and Lehmer's formula. These formulas are more practical.

@aamirkhan Floats are numbers with decimal points.

@Adam Ledger Never mind cleaning up that code or explaining the algorithm: I got it figured out. It was easier to figure out once I saw that the vast majority of your aliases and packages were unused.

Your procedure says, in English, "k is prime if the number of positive divisors of its positive divisors (counted with multiplicity) is 3." Indeed, those three divisors are [1, 1, k], not the [k, 1, 0] that you claimed.

Written in more-normal Maple (yet maintaining your algorithm) your procedure is equivalent to

P:= proc(N::posint)
local k;
uses numtheory;   #For `divisors`
     {seq(k*`if`(nops((op@[op]@divisors)~([divisors(k)[]])) = 3, 1, 0), k= 1..N)} minus {0}
end proc;

(An astute reader may notice that the [op] is superfluous. I only included it because it corresponds to the outer convert(..., list) in the OP's code.)

Of course that English translation is just a bloated way of saying "k is prime if the number of its positive divisors is 2." A simple Maple translation of that is

P:= N-> select(k-> nops(numtheory:-divisors(k))=2, {$1..N});

As far as I know, Maple's numeric PDE solver doesn't support the use of imaginary coefficients. Often, the pdsolve command itself will complete without error, but any practical attempt to use the solution module will result in error. I think that has happened in this case. Certainly, the direct cause of the problem is the imaginary coefficients.

@Adam Ledger Yes, I'd appreciate it if you made some modifications for readability and re-uploaded it. For my taste, the readability would be greatly improved if there was no use of alias or with whatsoever.

Regarding Statistics:-Count: The number of elements of a set or a list can be found with the top-level commands nops or numelems. These have the benefit of not requiring set-to-list conversion.

Regarding efficiency: You can use CodeTools:-Usage to make some measurements. There is no prepackaged prime-subset command that I'm aware of. But it's trivial to build one from the pre-existing commands. The most trivial I can come up with is

AllPrimes:= (n::posint)-> select(isprime, {$1..n}):

This is also surprisingly efficient. Comparing it to yours:

CodeTools:-Usage(`ℙ`(2^14)):   # `ℙ` is your procedure's name in ASCII.
memory used=0.57GiB, alloc change=0 bytes, cpu time=8.14s, real time=8.18s, gc time=390.62ms

CodeTools:-Usage(AllPrimes(2^14)):
memory used=6.74MiB, alloc change=0 bytes, cpu time=63.00ms, real time=52.00ms, gc time=0ns

(s = seconds, ms = milliseconds)

I have a highly polished Sieve of Eratosthenes procedure that's a little bit faster than AllPrimes. But it's on a different computer, and it'll take me a while to get it.

[A while passes.]

EratosthenesSieve:= proc(N::posint)
# Returns list of primes <= N
   local
      # Only record data for the odds >= 3
      M:= iquo(N-1,2)
     ,prime:= rtable(1..M, fill= true)
     ,p, P, k
   ;
   for p to trunc(evalf(sqrt(N)/2)) do
      if prime[p] then
         P:= 2*p+1;
         for k from (P+1)*p by P to M do  prime[k]:= false  end do
      end if
   end do;
   [`if`(N=1,[][],2), seq(`if`(prime[p], 2*p+1, [][]), p= 1..M)]
end proc:
     
CodeTools:-Usage(EratosthenesSieve(2^14)):
memory used=0.55MiB, alloc change=0 bytes, cpu time=16.00ms, real time=8.00ms, gc time=0ns

The Sieve procedure is written in "pure, native" Maple. I may be able to make it faster with evalhf or compilation.

I combined your three posts into one. 

It's an interesting algorithm that I've never seen before. Some very quick and preliminary tests show that the prime-generation part is correct. The counting part is off by one. 

Your code is exhausting to read due to the overuse of alias (especially) and also overuse of with. I think that there are some packages and aliases that you declare but don't use. Please consider rewriting this with the aliases and package names expanded.

List-to-set conversion of L can be done by {L[]}, and set-to-list conversion of S can be done by [S[]]. I find these easier to read than your single-letter operators. (I'm sure that some others would disagree with that.)

@Kitonum Thanks. While you were voting, I added a significant amount to the Reply, so you'll probably want to reread it.

@Kitonum 

piecewise uses assumptions when evaluating conditions with parameters; `if` doesn't. So, if n isn't manifestly an integer, then the n::odd in the `if` evaluates to false, regardless of assumptions. There is a distinction between types (which are checked with type) and properties (which are checked with is). Many types have the same name as an equivalent property; odd, even, and integer are examples of this. The `::` operator can be used to do a type check or a property check. Which is used depends on the context. piecewise is an is context; `if` is a type context.

The expression piecewise(cond, a, b) is almost equivalent to

(()-> `if`(is(cond)::truefalse, `if`(is(cond), a, b), 'procname'()))()

In other words, is forces a decision to true, false, or FAIL; piecewise just remains unevaluated rather than going to FAIL when the condition can't be decided under the current assumptions.

I say almost equivalent because the call to is from piecewise is moderated through a call to the undocumented procedure PiecewiseTools:-Is, and I don't know the exact details of that. In particular, there are many common circumstances where is(cond) evaluates false, but the piecewise remains unevaluated rather than simplifying to b. The code of PiecewiseTools:-Is is fairly easy to read, easier than the code of is.

@tomleslie 

Yes, I agree that the length check for the purposes of showing the excessive-length message should take into account whether the output is already elided. But that's not the way that it's implemented.

Also, the excessive-length message should include "This is not an error. The result has been computed."

@tomleslie 

Give the command

interface(verboseproc= 3);

Then the procedure definitions returned by dsolve won't be elided. Now compare the output of the following two dsolve commands, and you'll see why the length increases:

sol1:= dsolve(dsys3, numeric, initmesh= 8);
sol2:= dsolve(dsys3, numeric, initmesh= 16);

The OP hasn't provided a reason for such an unusually large value of initmesh, and I doubt that there is any good reason. Inceasing maxmesh to a large value often makes sense. I'd only increase initmesh if the problem seems to have trouble "getting started" (as indicated by an error message), and then I'd only increase it gradually.

@roman_pearce 

This is all-the-more reason why one should be able to set kernelopts(numcpus) before a call to Threads. The current situation---where one can only effectively set it immediately after a restart---is unacceptable.

If I simply change my test list to a million Maple floats rather than 32 million integers, then using any number of threads greater than one leads to a substantial increase in real time! It's about 5 times higher when using two threads and nearly 2 times higher when using eight threads. To run this test, simply change the command the creates the input data to

L:= evalf(RandomTools:-Generate(list(integer, 2^20))):

In another test, I tried to defeat the caching by adding elements in a random order. I changed the input to

L:= evalf(RandomTools:-Generate(list(integer, 2^20))):
J:= RandomTools:-Generate(list(integer(range= 1..2^20), 2^20)):
save L, J, "ThreadsTestData.m":

And I changed the test execution (both places!) to

Threads:-Add(L[k], k= J)

@roman_pearce Thank you for the detailed and very useful information, and thank you for running my test on other CPUs. I think that you've mostly explained the phenomena that I observed. So I should probably ignore CPU times for parallel processes on my CPU. Still, the real-time reduction going from four to eight threads on my CPU is very small.

Something else to consider is the sharing of the caches between threads on the same core.