Free cookie consent management tool by TermsFeed Policy Generator

source: branches/HeuristicLab.Problems.GaussianProcessTuning/ILNumerics.2.14.4735.573/Functions/statistics/select.cs @ 9102

Last change on this file since 9102 was 9102, checked in by gkronber, 12 years ago

#1967: ILNumerics source for experimentation

File size: 64.7 KB
Line 
1///
2///    This file is part of ILNumerics Community Edition.
3///
4///    ILNumerics Community Edition - high performance computing for applications.
5///    Copyright (C) 2006 - 2012 Haymo Kutschbach, http://ilnumerics.net
6///
7///    ILNumerics Community Edition is free software: you can redistribute it and/or modify
8///    it under the terms of the GNU General Public License version 3 as published by
9///    the Free Software Foundation.
10///
11///    ILNumerics Community Edition is distributed in the hope that it will be useful,
12///    but WITHOUT ANY WARRANTY; without even the implied warranty of
13///    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14///    GNU General Public License for more details.
15///
16///    You should have received a copy of the GNU General Public License
17///    along with ILNumerics Community Edition. See the file License.txt in the root
18///    of your distribution package. If not, see <http://www.gnu.org/licenses/>.
19///
20///    In addition this software uses the following components and/or licenses:
21///
22///    =================================================================================
23///    The Open Toolkit Library License
24///   
25///    Copyright (c) 2006 - 2009 the Open Toolkit library.
26///   
27///    Permission is hereby granted, free of charge, to any person obtaining a copy
28///    of this software and associated documentation files (the "Software"), to deal
29///    in the Software without restriction, including without limitation the rights to
30///    use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of
31///    the Software, and to permit persons to whom the Software is furnished to do
32///    so, subject to the following conditions:
33///
34///    The above copyright notice and this permission notice shall be included in all
35///    copies or substantial portions of the Software.
36///
37///    =================================================================================
38///   
39
40using System;
41using System.Collections.Generic;
42using System.Text;
43using System.Threading;
44using ILNumerics.Storage;
45using ILNumerics.Misc;
46using ILNumerics.Exceptions;
47 
48
49
50namespace ILNumerics  {
51    public partial class ILMath {
52
53        /// <summary>
54        /// Select the k-th smallest element from an array along a specific dimension
55        /// </summary>
56        /// <param name="A">Input array</param>
57        /// <param name="k">The element to find. If k is smaller 1 or larger than the number of elements in list, the smallest/largest value will be returned.</param>
58        /// <param name="dim">[Optional] Dimension to operate along. If omitted operates along the first non singleton dimension (i.e. != 1).</param>
59        /// <returns><para>Array having the specified dimension reduced to the length 1 with the value of the k-the smallest element along that dimension.</para>
60        /// <para>Exception: If the selected dimension is of size 0 it will remain 0 (an empty set).</para></returns>
61        public static ILRetArray<double> select(ILInArray<double> A, int k, int dim = -1) {
62            using (ILScope.Enter(A)) {
63                if (dim < 0)
64                    dim = A.Size.WorkingDimension();
65
66                if (dim >= A.Size.NumberOfDimensions)
67                    return A.C;
68               
69                if (A.IsScalar) {
70                    return new ILRetArray<double>(new double[] { A.GetValue(0) }, 1, 1);
71                }
72               
73                if (A.S[dim] == 1) return A.C;
74
75                int[] newDims = A.S.ToIntArray();
76                if (A.IsEmpty)
77                {
78                    // If the array is empty, check whether it is empty along the chosen dimension
79                   
80                    if (A.S[dim] > 0)
81                        // no there are potential elements in that dimension, hence we reduce it to 1
82                        newDims[dim] = 1;
83                    else
84                        // yes, empty along chosen dimension so the result is empty, i.e. 0 elements in that dimension
85                        newDims[dim] = 0;
86                    return ILRetArray<double>.empty(new ILSize(newDims));
87                }
88
89                newDims[dim] = 1;
90
91                // Check selected element and replace by useful value
92                if (k < 1)
93                    k = 1;
94                if (k > A.S[dim])
95                    k = A.S[dim];
96
97                ILSize retDimension = new ILSize(newDims);
98
99                 double[] retArr = ILMemoryPool.Pool.New< double>(retDimension.NumberOfElements);
100
101                int inc = A.Size.SequentialIndexDistance(dim);
102                int dimLen = A.Size[dim];
103                int maxRuns = retDimension.NumberOfElements;
104                int modHelp = A.Size.NumberOfElements - 1;
105                int modOut = retDimension.NumberOfElements - 1;
106                int incOut = retDimension.SequentialIndexDistance(dim);
107                int numelA = A.S.NumberOfElements;
108                if (maxRuns == 1) {
109                    int dummy;
110                    retArr[0] = quickselect_worker(A.C.GetArrayForWrite(), 0, A.S[dim] - 1, k, out dummy);
111                } else {
112                    #region may run parallel
113                   
114                    double[] aArray = A.GetArrayForRead();
115                    int i = 0, workItemCount = Settings.s_maxNumberThreads, workItemLength, workerCount = 1;
116                    if (Settings.s_maxNumberThreads > 1 && maxRuns > 1
117                        && numelA / 2 >= Settings.s_minParallelElement1Count) {
118
119                        if (maxRuns >= Settings.s_maxNumberThreads
120                            && numelA / Settings.s_maxNumberThreads > Settings.s_minParallelElement1Count) {
121                            workItemLength = maxRuns / workItemCount;
122                        } else {
123                            workItemLength = maxRuns / 2;
124                            workItemCount = 2;
125                        }
126                       
127                    } else {
128                        workItemLength = maxRuns;
129                        workItemCount = 1;
130                    }
131                    Action<object> action = (data) => {
132                        Tuple<int, int> range = (Tuple<int, int>)data;
133                        int from = range.Item1, to = range.Item2;
134                        for (int c = from; c < to; c++) {
135                            int pos = (int)(((long)dimLen * c * inc) % modHelp);
136                            int posOut = (c * incOut);
137                            if (posOut > modOut)
138                                posOut = ((posOut - 1) % modOut) + 1;
139                           
140                            double[] tmp = ILMemoryPool.Pool.New< double>(dimLen);
141                            int locPos = 0;
142                            int end = pos + dimLen * inc;
143                            for (int j = pos; j < end; j += inc)
144                                tmp[locPos++] = aArray[j];
145                            int dummy;
146                            retArr[posOut] = quickselect_worker(tmp, 0, dimLen - 1, k, out dummy);
147                            ILMemoryPool.Pool.Free(tmp);
148                        }
149                        System.Threading.Interlocked.Decrement(ref workerCount);
150                    };
151                    for (; i < workItemCount - 1; i++) {
152                        Interlocked.Increment(ref workerCount);
153
154                        ILThreadPool.QueueUserWorkItem(i,action, Tuple.Create(i * workItemLength, (i + 1) * workItemLength));
155                    }
156                    action(Tuple.Create(i * workItemLength, maxRuns));
157                    ILThreadPool.Wait4Workers(ref workerCount);
158                    #endregion
159                }
160                return new ILRetArray<double>(retArr, newDims);
161            }
162        }
163
164#region HYCALPER AUTO GENERATED CODE
165
166        /// <summary>
167        /// Select the k-th smallest element from an array along a specific dimension
168        /// </summary>
169        /// <param name="A">Input array</param>
170        /// <param name="k">The element to find. If k is smaller 1 or larger than the number of elements in list, the smallest/largest value will be returned.</param>
171        /// <param name="dim">[Optional] Dimension to operate along. If omitted operates along the first non singleton dimension (i.e. != 1).</param>
172        /// <returns><para>Array having the specified dimension reduced to the length 1 with the value of the k-the smallest element along that dimension.</para>
173        /// <para>Exception: If the selected dimension is of size 0 it will remain 0 (an empty set).</para></returns>
174        public static ILRetArray<Int64> select(ILInArray<Int64> A, int k, int dim = -1) {
175            using (ILScope.Enter(A)) {
176                if (dim < 0)
177                    dim = A.Size.WorkingDimension();
178
179                if (dim >= A.Size.NumberOfDimensions)
180                    return A.C;
181               
182                if (A.IsScalar) {
183                    return new ILRetArray<Int64>(new Int64[] { A.GetValue(0) }, 1, 1);
184                }
185               
186                if (A.S[dim] == 1) return A.C;
187
188                int[] newDims = A.S.ToIntArray();
189                if (A.IsEmpty)
190                {
191                    // If the array is empty, check whether it is empty along the chosen dimension
192                   
193                    if (A.S[dim] > 0)
194                        // no there are potential elements in that dimension, hence we reduce it to 1
195                        newDims[dim] = 1;
196                    else
197                        // yes, empty along chosen dimension so the result is empty, i.e. 0 elements in that dimension
198                        newDims[dim] = 0;
199                    return ILRetArray<Int64>.empty(new ILSize(newDims));
200                }
201
202                newDims[dim] = 1;
203
204                // Check selected element and replace by useful value
205                if (k < 1)
206                    k = 1;
207                if (k > A.S[dim])
208                    k = A.S[dim];
209
210                ILSize retDimension = new ILSize(newDims);
211
212                Int64[] retArr = ILMemoryPool.Pool.New< Int64>(retDimension.NumberOfElements);
213
214                int inc = A.Size.SequentialIndexDistance(dim);
215                int dimLen = A.Size[dim];
216                int maxRuns = retDimension.NumberOfElements;
217                int modHelp = A.Size.NumberOfElements - 1;
218                int modOut = retDimension.NumberOfElements - 1;
219                int incOut = retDimension.SequentialIndexDistance(dim);
220                int numelA = A.S.NumberOfElements;
221                if (maxRuns == 1) {
222                    int dummy;
223                    retArr[0] = quickselect_worker(A.C.GetArrayForWrite(), 0, A.S[dim] - 1, k, out dummy);
224                } else {
225                    #region may run parallel
226                   
227                    Int64[] aArray = A.GetArrayForRead();
228                    int i = 0, workItemCount = Settings.s_maxNumberThreads, workItemLength, workerCount = 1;
229                    if (Settings.s_maxNumberThreads > 1 && maxRuns > 1
230                        && numelA / 2 >= Settings.s_minParallelElement1Count) {
231
232                        if (maxRuns >= Settings.s_maxNumberThreads
233                            && numelA / Settings.s_maxNumberThreads > Settings.s_minParallelElement1Count) {
234                            workItemLength = maxRuns / workItemCount;
235                        } else {
236                            workItemLength = maxRuns / 2;
237                            workItemCount = 2;
238                        }
239                       
240                    } else {
241                        workItemLength = maxRuns;
242                        workItemCount = 1;
243                    }
244                    Action<object> action = (data) => {
245                        Tuple<int, int> range = (Tuple<int, int>)data;
246                        int from = range.Item1, to = range.Item2;
247                        for (int c = from; c < to; c++) {
248                            int pos = (int)(((long)dimLen * c * inc) % modHelp);
249                            int posOut = (c * incOut);
250                            if (posOut > modOut)
251                                posOut = ((posOut - 1) % modOut) + 1;
252                           
253                            Int64[] tmp = ILMemoryPool.Pool.New< Int64>(dimLen);
254                            int locPos = 0;
255                            int end = pos + dimLen * inc;
256                            for (int j = pos; j < end; j += inc)
257                                tmp[locPos++] = aArray[j];
258                            int dummy;
259                            retArr[posOut] = quickselect_worker(tmp, 0, dimLen - 1, k, out dummy);
260                            ILMemoryPool.Pool.Free(tmp);
261                        }
262                        System.Threading.Interlocked.Decrement(ref workerCount);
263                    };
264                    for (; i < workItemCount - 1; i++) {
265                        Interlocked.Increment(ref workerCount);
266
267                        ILThreadPool.QueueUserWorkItem(i,action, Tuple.Create(i * workItemLength, (i + 1) * workItemLength));
268                    }
269                    action(Tuple.Create(i * workItemLength, maxRuns));
270                    ILThreadPool.Wait4Workers(ref workerCount);
271                    #endregion
272                }
273                return new ILRetArray<Int64>(retArr, newDims);
274            }
275        }
276        /// <summary>
277        /// Select the k-th smallest element from an array along a specific dimension
278        /// </summary>
279        /// <param name="A">Input array</param>
280        /// <param name="k">The element to find. If k is smaller 1 or larger than the number of elements in list, the smallest/largest value will be returned.</param>
281        /// <param name="dim">[Optional] Dimension to operate along. If omitted operates along the first non singleton dimension (i.e. != 1).</param>
282        /// <returns><para>Array having the specified dimension reduced to the length 1 with the value of the k-the smallest element along that dimension.</para>
283        /// <para>Exception: If the selected dimension is of size 0 it will remain 0 (an empty set).</para></returns>
284        public static ILRetArray<Int32> select(ILInArray<Int32> A, int k, int dim = -1) {
285            using (ILScope.Enter(A)) {
286                if (dim < 0)
287                    dim = A.Size.WorkingDimension();
288
289                if (dim >= A.Size.NumberOfDimensions)
290                    return A.C;
291               
292                if (A.IsScalar) {
293                    return new ILRetArray<Int32>(new Int32[] { A.GetValue(0) }, 1, 1);
294                }
295               
296                if (A.S[dim] == 1) return A.C;
297
298                int[] newDims = A.S.ToIntArray();
299                if (A.IsEmpty)
300                {
301                    // If the array is empty, check whether it is empty along the chosen dimension
302                   
303                    if (A.S[dim] > 0)
304                        // no there are potential elements in that dimension, hence we reduce it to 1
305                        newDims[dim] = 1;
306                    else
307                        // yes, empty along chosen dimension so the result is empty, i.e. 0 elements in that dimension
308                        newDims[dim] = 0;
309                    return ILRetArray<Int32>.empty(new ILSize(newDims));
310                }
311
312                newDims[dim] = 1;
313
314                // Check selected element and replace by useful value
315                if (k < 1)
316                    k = 1;
317                if (k > A.S[dim])
318                    k = A.S[dim];
319
320                ILSize retDimension = new ILSize(newDims);
321
322                Int32[] retArr = ILMemoryPool.Pool.New< Int32>(retDimension.NumberOfElements);
323
324                int inc = A.Size.SequentialIndexDistance(dim);
325                int dimLen = A.Size[dim];
326                int maxRuns = retDimension.NumberOfElements;
327                int modHelp = A.Size.NumberOfElements - 1;
328                int modOut = retDimension.NumberOfElements - 1;
329                int incOut = retDimension.SequentialIndexDistance(dim);
330                int numelA = A.S.NumberOfElements;
331                if (maxRuns == 1) {
332                    int dummy;
333                    retArr[0] = quickselect_worker(A.C.GetArrayForWrite(), 0, A.S[dim] - 1, k, out dummy);
334                } else {
335                    #region may run parallel
336                   
337                    Int32[] aArray = A.GetArrayForRead();
338                    int i = 0, workItemCount = Settings.s_maxNumberThreads, workItemLength, workerCount = 1;
339                    if (Settings.s_maxNumberThreads > 1 && maxRuns > 1
340                        && numelA / 2 >= Settings.s_minParallelElement1Count) {
341
342                        if (maxRuns >= Settings.s_maxNumberThreads
343                            && numelA / Settings.s_maxNumberThreads > Settings.s_minParallelElement1Count) {
344                            workItemLength = maxRuns / workItemCount;
345                        } else {
346                            workItemLength = maxRuns / 2;
347                            workItemCount = 2;
348                        }
349                       
350                    } else {
351                        workItemLength = maxRuns;
352                        workItemCount = 1;
353                    }
354                    Action<object> action = (data) => {
355                        Tuple<int, int> range = (Tuple<int, int>)data;
356                        int from = range.Item1, to = range.Item2;
357                        for (int c = from; c < to; c++) {
358                            int pos = (int)(((long)dimLen * c * inc) % modHelp);
359                            int posOut = (c * incOut);
360                            if (posOut > modOut)
361                                posOut = ((posOut - 1) % modOut) + 1;
362                           
363                            Int32[] tmp = ILMemoryPool.Pool.New< Int32>(dimLen);
364                            int locPos = 0;
365                            int end = pos + dimLen * inc;
366                            for (int j = pos; j < end; j += inc)
367                                tmp[locPos++] = aArray[j];
368                            int dummy;
369                            retArr[posOut] = quickselect_worker(tmp, 0, dimLen - 1, k, out dummy);
370                            ILMemoryPool.Pool.Free(tmp);
371                        }
372                        System.Threading.Interlocked.Decrement(ref workerCount);
373                    };
374                    for (; i < workItemCount - 1; i++) {
375                        Interlocked.Increment(ref workerCount);
376
377                        ILThreadPool.QueueUserWorkItem(i,action, Tuple.Create(i * workItemLength, (i + 1) * workItemLength));
378                    }
379                    action(Tuple.Create(i * workItemLength, maxRuns));
380                    ILThreadPool.Wait4Workers(ref workerCount);
381                    #endregion
382                }
383                return new ILRetArray<Int32>(retArr, newDims);
384            }
385        }
386        /// <summary>
387        /// Select the k-th smallest element from an array along a specific dimension
388        /// </summary>
389        /// <param name="A">Input array</param>
390        /// <param name="k">The element to find. If k is smaller 1 or larger than the number of elements in list, the smallest/largest value will be returned.</param>
391        /// <param name="dim">[Optional] Dimension to operate along. If omitted operates along the first non singleton dimension (i.e. != 1).</param>
392        /// <returns><para>Array having the specified dimension reduced to the length 1 with the value of the k-the smallest element along that dimension.</para>
393        /// <para>Exception: If the selected dimension is of size 0 it will remain 0 (an empty set).</para></returns>
394        public static ILRetArray<byte> select(ILInArray<byte> A, int k, int dim = -1) {
395            using (ILScope.Enter(A)) {
396                if (dim < 0)
397                    dim = A.Size.WorkingDimension();
398
399                if (dim >= A.Size.NumberOfDimensions)
400                    return A.C;
401               
402                if (A.IsScalar) {
403                    return new ILRetArray<byte>(new byte[] { A.GetValue(0) }, 1, 1);
404                }
405               
406                if (A.S[dim] == 1) return A.C;
407
408                int[] newDims = A.S.ToIntArray();
409                if (A.IsEmpty)
410                {
411                    // If the array is empty, check whether it is empty along the chosen dimension
412                   
413                    if (A.S[dim] > 0)
414                        // no there are potential elements in that dimension, hence we reduce it to 1
415                        newDims[dim] = 1;
416                    else
417                        // yes, empty along chosen dimension so the result is empty, i.e. 0 elements in that dimension
418                        newDims[dim] = 0;
419                    return ILRetArray<byte>.empty(new ILSize(newDims));
420                }
421
422                newDims[dim] = 1;
423
424                // Check selected element and replace by useful value
425                if (k < 1)
426                    k = 1;
427                if (k > A.S[dim])
428                    k = A.S[dim];
429
430                ILSize retDimension = new ILSize(newDims);
431
432                byte[] retArr = ILMemoryPool.Pool.New< byte>(retDimension.NumberOfElements);
433
434                int inc = A.Size.SequentialIndexDistance(dim);
435                int dimLen = A.Size[dim];
436                int maxRuns = retDimension.NumberOfElements;
437                int modHelp = A.Size.NumberOfElements - 1;
438                int modOut = retDimension.NumberOfElements - 1;
439                int incOut = retDimension.SequentialIndexDistance(dim);
440                int numelA = A.S.NumberOfElements;
441                if (maxRuns == 1) {
442                    int dummy;
443                    retArr[0] = quickselect_worker(A.C.GetArrayForWrite(), 0, A.S[dim] - 1, k, out dummy);
444                } else {
445                    #region may run parallel
446                   
447                    byte[] aArray = A.GetArrayForRead();
448                    int i = 0, workItemCount = Settings.s_maxNumberThreads, workItemLength, workerCount = 1;
449                    if (Settings.s_maxNumberThreads > 1 && maxRuns > 1
450                        && numelA / 2 >= Settings.s_minParallelElement1Count) {
451
452                        if (maxRuns >= Settings.s_maxNumberThreads
453                            && numelA / Settings.s_maxNumberThreads > Settings.s_minParallelElement1Count) {
454                            workItemLength = maxRuns / workItemCount;
455                        } else {
456                            workItemLength = maxRuns / 2;
457                            workItemCount = 2;
458                        }
459                       
460                    } else {
461                        workItemLength = maxRuns;
462                        workItemCount = 1;
463                    }
464                    Action<object> action = (data) => {
465                        Tuple<int, int> range = (Tuple<int, int>)data;
466                        int from = range.Item1, to = range.Item2;
467                        for (int c = from; c < to; c++) {
468                            int pos = (int)(((long)dimLen * c * inc) % modHelp);
469                            int posOut = (c * incOut);
470                            if (posOut > modOut)
471                                posOut = ((posOut - 1) % modOut) + 1;
472                           
473                            byte[] tmp = ILMemoryPool.Pool.New< byte>(dimLen);
474                            int locPos = 0;
475                            int end = pos + dimLen * inc;
476                            for (int j = pos; j < end; j += inc)
477                                tmp[locPos++] = aArray[j];
478                            int dummy;
479                            retArr[posOut] = quickselect_worker(tmp, 0, dimLen - 1, k, out dummy);
480                            ILMemoryPool.Pool.Free(tmp);
481                        }
482                        System.Threading.Interlocked.Decrement(ref workerCount);
483                    };
484                    for (; i < workItemCount - 1; i++) {
485                        Interlocked.Increment(ref workerCount);
486
487                        ILThreadPool.QueueUserWorkItem(i,action, Tuple.Create(i * workItemLength, (i + 1) * workItemLength));
488                    }
489                    action(Tuple.Create(i * workItemLength, maxRuns));
490                    ILThreadPool.Wait4Workers(ref workerCount);
491                    #endregion
492                }
493                return new ILRetArray<byte>(retArr, newDims);
494            }
495        }
496        /// <summary>
497        /// Select the k-th smallest element from an array along a specific dimension
498        /// </summary>
499        /// <param name="A">Input array</param>
500        /// <param name="k">The element to find. If k is smaller 1 or larger than the number of elements in list, the smallest/largest value will be returned.</param>
501        /// <param name="dim">[Optional] Dimension to operate along. If omitted operates along the first non singleton dimension (i.e. != 1).</param>
502        /// <returns><para>Array having the specified dimension reduced to the length 1 with the value of the k-the smallest element along that dimension.</para>
503        /// <para>Exception: If the selected dimension is of size 0 it will remain 0 (an empty set).</para></returns>
504        public static ILRetArray<fcomplex> select(ILInArray<fcomplex> A, int k, int dim = -1) {
505            using (ILScope.Enter(A)) {
506                if (dim < 0)
507                    dim = A.Size.WorkingDimension();
508
509                if (dim >= A.Size.NumberOfDimensions)
510                    return A.C;
511               
512                if (A.IsScalar) {
513                    return new ILRetArray<fcomplex>(new fcomplex[] { A.GetValue(0) }, 1, 1);
514                }
515               
516                if (A.S[dim] == 1) return A.C;
517
518                int[] newDims = A.S.ToIntArray();
519                if (A.IsEmpty)
520                {
521                    // If the array is empty, check whether it is empty along the chosen dimension
522                   
523                    if (A.S[dim] > 0)
524                        // no there are potential elements in that dimension, hence we reduce it to 1
525                        newDims[dim] = 1;
526                    else
527                        // yes, empty along chosen dimension so the result is empty, i.e. 0 elements in that dimension
528                        newDims[dim] = 0;
529                    return ILRetArray<fcomplex>.empty(new ILSize(newDims));
530                }
531
532                newDims[dim] = 1;
533
534                // Check selected element and replace by useful value
535                if (k < 1)
536                    k = 1;
537                if (k > A.S[dim])
538                    k = A.S[dim];
539
540                ILSize retDimension = new ILSize(newDims);
541
542                fcomplex[] retArr = ILMemoryPool.Pool.New< fcomplex>(retDimension.NumberOfElements);
543
544                int inc = A.Size.SequentialIndexDistance(dim);
545                int dimLen = A.Size[dim];
546                int maxRuns = retDimension.NumberOfElements;
547                int modHelp = A.Size.NumberOfElements - 1;
548                int modOut = retDimension.NumberOfElements - 1;
549                int incOut = retDimension.SequentialIndexDistance(dim);
550                int numelA = A.S.NumberOfElements;
551                if (maxRuns == 1) {
552                    int dummy;
553                    retArr[0] = quickselect_worker(A.C.GetArrayForWrite(), 0, A.S[dim] - 1, k, out dummy);
554                } else {
555                    #region may run parallel
556                   
557                    fcomplex[] aArray = A.GetArrayForRead();
558                    int i = 0, workItemCount = Settings.s_maxNumberThreads, workItemLength, workerCount = 1;
559                    if (Settings.s_maxNumberThreads > 1 && maxRuns > 1
560                        && numelA / 2 >= Settings.s_minParallelElement1Count) {
561
562                        if (maxRuns >= Settings.s_maxNumberThreads
563                            && numelA / Settings.s_maxNumberThreads > Settings.s_minParallelElement1Count) {
564                            workItemLength = maxRuns / workItemCount;
565                        } else {
566                            workItemLength = maxRuns / 2;
567                            workItemCount = 2;
568                        }
569                       
570                    } else {
571                        workItemLength = maxRuns;
572                        workItemCount = 1;
573                    }
574                    Action<object> action = (data) => {
575                        Tuple<int, int> range = (Tuple<int, int>)data;
576                        int from = range.Item1, to = range.Item2;
577                        for (int c = from; c < to; c++) {
578                            int pos = (int)(((long)dimLen * c * inc) % modHelp);
579                            int posOut = (c * incOut);
580                            if (posOut > modOut)
581                                posOut = ((posOut - 1) % modOut) + 1;
582                           
583                            fcomplex[] tmp = ILMemoryPool.Pool.New< fcomplex>(dimLen);
584                            int locPos = 0;
585                            int end = pos + dimLen * inc;
586                            for (int j = pos; j < end; j += inc)
587                                tmp[locPos++] = aArray[j];
588                            int dummy;
589                            retArr[posOut] = quickselect_worker(tmp, 0, dimLen - 1, k, out dummy);
590                            ILMemoryPool.Pool.Free(tmp);
591                        }
592                        System.Threading.Interlocked.Decrement(ref workerCount);
593                    };
594                    for (; i < workItemCount - 1; i++) {
595                        Interlocked.Increment(ref workerCount);
596
597                        ILThreadPool.QueueUserWorkItem(i,action, Tuple.Create(i * workItemLength, (i + 1) * workItemLength));
598                    }
599                    action(Tuple.Create(i * workItemLength, maxRuns));
600                    ILThreadPool.Wait4Workers(ref workerCount);
601                    #endregion
602                }
603                return new ILRetArray<fcomplex>(retArr, newDims);
604            }
605        }
606        /// <summary>
607        /// Select the k-th smallest element from an array along a specific dimension
608        /// </summary>
609        /// <param name="A">Input array</param>
610        /// <param name="k">The element to find. If k is smaller 1 or larger than the number of elements in list, the smallest/largest value will be returned.</param>
611        /// <param name="dim">[Optional] Dimension to operate along. If omitted operates along the first non singleton dimension (i.e. != 1).</param>
612        /// <returns><para>Array having the specified dimension reduced to the length 1 with the value of the k-the smallest element along that dimension.</para>
613        /// <para>Exception: If the selected dimension is of size 0 it will remain 0 (an empty set).</para></returns>
614        public static ILRetArray<float> select(ILInArray<float> A, int k, int dim = -1) {
615            using (ILScope.Enter(A)) {
616                if (dim < 0)
617                    dim = A.Size.WorkingDimension();
618
619                if (dim >= A.Size.NumberOfDimensions)
620                    return A.C;
621               
622                if (A.IsScalar) {
623                    return new ILRetArray<float>(new float[] { A.GetValue(0) }, 1, 1);
624                }
625               
626                if (A.S[dim] == 1) return A.C;
627
628                int[] newDims = A.S.ToIntArray();
629                if (A.IsEmpty)
630                {
631                    // If the array is empty, check whether it is empty along the chosen dimension
632                   
633                    if (A.S[dim] > 0)
634                        // no there are potential elements in that dimension, hence we reduce it to 1
635                        newDims[dim] = 1;
636                    else
637                        // yes, empty along chosen dimension so the result is empty, i.e. 0 elements in that dimension
638                        newDims[dim] = 0;
639                    return ILRetArray<float>.empty(new ILSize(newDims));
640                }
641
642                newDims[dim] = 1;
643
644                // Check selected element and replace by useful value
645                if (k < 1)
646                    k = 1;
647                if (k > A.S[dim])
648                    k = A.S[dim];
649
650                ILSize retDimension = new ILSize(newDims);
651
652                float[] retArr = ILMemoryPool.Pool.New< float>(retDimension.NumberOfElements);
653
654                int inc = A.Size.SequentialIndexDistance(dim);
655                int dimLen = A.Size[dim];
656                int maxRuns = retDimension.NumberOfElements;
657                int modHelp = A.Size.NumberOfElements - 1;
658                int modOut = retDimension.NumberOfElements - 1;
659                int incOut = retDimension.SequentialIndexDistance(dim);
660                int numelA = A.S.NumberOfElements;
661                if (maxRuns == 1) {
662                    int dummy;
663                    retArr[0] = quickselect_worker(A.C.GetArrayForWrite(), 0, A.S[dim] - 1, k, out dummy);
664                } else {
665                    #region may run parallel
666                   
667                    float[] aArray = A.GetArrayForRead();
668                    int i = 0, workItemCount = Settings.s_maxNumberThreads, workItemLength, workerCount = 1;
669                    if (Settings.s_maxNumberThreads > 1 && maxRuns > 1
670                        && numelA / 2 >= Settings.s_minParallelElement1Count) {
671
672                        if (maxRuns >= Settings.s_maxNumberThreads
673                            && numelA / Settings.s_maxNumberThreads > Settings.s_minParallelElement1Count) {
674                            workItemLength = maxRuns / workItemCount;
675                        } else {
676                            workItemLength = maxRuns / 2;
677                            workItemCount = 2;
678                        }
679                       
680                    } else {
681                        workItemLength = maxRuns;
682                        workItemCount = 1;
683                    }
684                    Action<object> action = (data) => {
685                        Tuple<int, int> range = (Tuple<int, int>)data;
686                        int from = range.Item1, to = range.Item2;
687                        for (int c = from; c < to; c++) {
688                            int pos = (int)(((long)dimLen * c * inc) % modHelp);
689                            int posOut = (c * incOut);
690                            if (posOut > modOut)
691                                posOut = ((posOut - 1) % modOut) + 1;
692                           
693                            float[] tmp = ILMemoryPool.Pool.New< float>(dimLen);
694                            int locPos = 0;
695                            int end = pos + dimLen * inc;
696                            for (int j = pos; j < end; j += inc)
697                                tmp[locPos++] = aArray[j];
698                            int dummy;
699                            retArr[posOut] = quickselect_worker(tmp, 0, dimLen - 1, k, out dummy);
700                            ILMemoryPool.Pool.Free(tmp);
701                        }
702                        System.Threading.Interlocked.Decrement(ref workerCount);
703                    };
704                    for (; i < workItemCount - 1; i++) {
705                        Interlocked.Increment(ref workerCount);
706
707                        ILThreadPool.QueueUserWorkItem(i,action, Tuple.Create(i * workItemLength, (i + 1) * workItemLength));
708                    }
709                    action(Tuple.Create(i * workItemLength, maxRuns));
710                    ILThreadPool.Wait4Workers(ref workerCount);
711                    #endregion
712                }
713                return new ILRetArray<float>(retArr, newDims);
714            }
715        }
716        /// <summary>
717        /// Select the k-th smallest element from an array along a specific dimension
718        /// </summary>
719        /// <param name="A">Input array</param>
720        /// <param name="k">The element to find. If k is smaller 1 or larger than the number of elements in list, the smallest/largest value will be returned.</param>
721        /// <param name="dim">[Optional] Dimension to operate along. If omitted operates along the first non singleton dimension (i.e. != 1).</param>
722        /// <returns><para>Array having the specified dimension reduced to the length 1 with the value of the k-the smallest element along that dimension.</para>
723        /// <para>Exception: If the selected dimension is of size 0 it will remain 0 (an empty set).</para></returns>
724        public static ILRetArray<complex> select(ILInArray<complex> A, int k, int dim = -1) {
725            using (ILScope.Enter(A)) {
726                if (dim < 0)
727                    dim = A.Size.WorkingDimension();
728
729                if (dim >= A.Size.NumberOfDimensions)
730                    return A.C;
731               
732                if (A.IsScalar) {
733                    return new ILRetArray<complex>(new complex[] { A.GetValue(0) }, 1, 1);
734                }
735               
736                if (A.S[dim] == 1) return A.C;
737
738                int[] newDims = A.S.ToIntArray();
739                if (A.IsEmpty)
740                {
741                    // If the array is empty, check whether it is empty along the chosen dimension
742                   
743                    if (A.S[dim] > 0)
744                        // no there are potential elements in that dimension, hence we reduce it to 1
745                        newDims[dim] = 1;
746                    else
747                        // yes, empty along chosen dimension so the result is empty, i.e. 0 elements in that dimension
748                        newDims[dim] = 0;
749                    return ILRetArray<complex>.empty(new ILSize(newDims));
750                }
751
752                newDims[dim] = 1;
753
754                // Check selected element and replace by useful value
755                if (k < 1)
756                    k = 1;
757                if (k > A.S[dim])
758                    k = A.S[dim];
759
760                ILSize retDimension = new ILSize(newDims);
761
762                complex[] retArr = ILMemoryPool.Pool.New< complex>(retDimension.NumberOfElements);
763
764                int inc = A.Size.SequentialIndexDistance(dim);
765                int dimLen = A.Size[dim];
766                int maxRuns = retDimension.NumberOfElements;
767                int modHelp = A.Size.NumberOfElements - 1;
768                int modOut = retDimension.NumberOfElements - 1;
769                int incOut = retDimension.SequentialIndexDistance(dim);
770                int numelA = A.S.NumberOfElements;
771                if (maxRuns == 1) {
772                    int dummy;
773                    retArr[0] = quickselect_worker(A.C.GetArrayForWrite(), 0, A.S[dim] - 1, k, out dummy);
774                } else {
775                    #region may run parallel
776                   
777                    complex[] aArray = A.GetArrayForRead();
778                    int i = 0, workItemCount = Settings.s_maxNumberThreads, workItemLength, workerCount = 1;
779                    if (Settings.s_maxNumberThreads > 1 && maxRuns > 1
780                        && numelA / 2 >= Settings.s_minParallelElement1Count) {
781
782                        if (maxRuns >= Settings.s_maxNumberThreads
783                            && numelA / Settings.s_maxNumberThreads > Settings.s_minParallelElement1Count) {
784                            workItemLength = maxRuns / workItemCount;
785                        } else {
786                            workItemLength = maxRuns / 2;
787                            workItemCount = 2;
788                        }
789                       
790                    } else {
791                        workItemLength = maxRuns;
792                        workItemCount = 1;
793                    }
794                    Action<object> action = (data) => {
795                        Tuple<int, int> range = (Tuple<int, int>)data;
796                        int from = range.Item1, to = range.Item2;
797                        for (int c = from; c < to; c++) {
798                            int pos = (int)(((long)dimLen * c * inc) % modHelp);
799                            int posOut = (c * incOut);
800                            if (posOut > modOut)
801                                posOut = ((posOut - 1) % modOut) + 1;
802                           
803                            complex[] tmp = ILMemoryPool.Pool.New< complex>(dimLen);
804                            int locPos = 0;
805                            int end = pos + dimLen * inc;
806                            for (int j = pos; j < end; j += inc)
807                                tmp[locPos++] = aArray[j];
808                            int dummy;
809                            retArr[posOut] = quickselect_worker(tmp, 0, dimLen - 1, k, out dummy);
810                            ILMemoryPool.Pool.Free(tmp);
811                        }
812                        System.Threading.Interlocked.Decrement(ref workerCount);
813                    };
814                    for (; i < workItemCount - 1; i++) {
815                        Interlocked.Increment(ref workerCount);
816
817                        ILThreadPool.QueueUserWorkItem(i,action, Tuple.Create(i * workItemLength, (i + 1) * workItemLength));
818                    }
819                    action(Tuple.Create(i * workItemLength, maxRuns));
820                    ILThreadPool.Wait4Workers(ref workerCount);
821                    #endregion
822                }
823                return new ILRetArray<complex>(retArr, newDims);
824            }
825        }
826
827#endregion HYCALPER AUTO GENERATED CODE
828
829        #region private helpers
830
831        private static int partition(double[] list, int left, int right, int pivotIndex)
832        {
833           
834            double pivotValue = list[pivotIndex];
835           
836            double tmp = list[right];
837            list[right] = list[pivotIndex];
838            list[pivotIndex] = tmp; // Move pivot to end
839            int storeIndex = left;
840            for (int i = left; i <= right; i++)
841            {
842                if (list[i] < pivotValue)
843                {
844                    tmp = list[storeIndex];
845                    list[storeIndex] = list[i];
846                    list[i] = tmp;
847                    storeIndex++;
848                }
849            }
850            // Move pivot to its final place
851            tmp = list[right];
852            list[right] = list[storeIndex];
853            list[storeIndex] = tmp;
854
855            return storeIndex;
856        }
857
858#region HYCALPER AUTO GENERATED CODE
859
860        private static int partition(Int64[] list, int left, int right, int pivotIndex)
861        {
862           
863            Int64 pivotValue = list[pivotIndex];
864           
865            Int64 tmp = list[right];
866            list[right] = list[pivotIndex];
867            list[pivotIndex] = tmp; // Move pivot to end
868            int storeIndex = left;
869            for (int i = left; i <= right; i++)
870            {
871                if (list[i] < pivotValue)
872                {
873                    tmp = list[storeIndex];
874                    list[storeIndex] = list[i];
875                    list[i] = tmp;
876                    storeIndex++;
877                }
878            }
879            // Move pivot to its final place
880            tmp = list[right];
881            list[right] = list[storeIndex];
882            list[storeIndex] = tmp;
883
884            return storeIndex;
885        }
886        private static int partition(Int32[] list, int left, int right, int pivotIndex)
887        {
888           
889            Int32 pivotValue = list[pivotIndex];
890           
891            Int32 tmp = list[right];
892            list[right] = list[pivotIndex];
893            list[pivotIndex] = tmp; // Move pivot to end
894            int storeIndex = left;
895            for (int i = left; i <= right; i++)
896            {
897                if (list[i] < pivotValue)
898                {
899                    tmp = list[storeIndex];
900                    list[storeIndex] = list[i];
901                    list[i] = tmp;
902                    storeIndex++;
903                }
904            }
905            // Move pivot to its final place
906            tmp = list[right];
907            list[right] = list[storeIndex];
908            list[storeIndex] = tmp;
909
910            return storeIndex;
911        }
912        private static int partition(byte[] list, int left, int right, int pivotIndex)
913        {
914           
915            byte pivotValue = list[pivotIndex];
916           
917            byte tmp = list[right];
918            list[right] = list[pivotIndex];
919            list[pivotIndex] = tmp; // Move pivot to end
920            int storeIndex = left;
921            for (int i = left; i <= right; i++)
922            {
923                if (list[i] < pivotValue)
924                {
925                    tmp = list[storeIndex];
926                    list[storeIndex] = list[i];
927                    list[i] = tmp;
928                    storeIndex++;
929                }
930            }
931            // Move pivot to its final place
932            tmp = list[right];
933            list[right] = list[storeIndex];
934            list[storeIndex] = tmp;
935
936            return storeIndex;
937        }
938        private static int partition(fcomplex[] list, int left, int right, int pivotIndex)
939        {
940           
941            fcomplex pivotValue = list[pivotIndex];
942           
943            fcomplex tmp = list[right];
944            list[right] = list[pivotIndex];
945            list[pivotIndex] = tmp; // Move pivot to end
946            int storeIndex = left;
947            for (int i = left; i <= right; i++)
948            {
949                if (list[i] < pivotValue)
950                {
951                    tmp = list[storeIndex];
952                    list[storeIndex] = list[i];
953                    list[i] = tmp;
954                    storeIndex++;
955                }
956            }
957            // Move pivot to its final place
958            tmp = list[right];
959            list[right] = list[storeIndex];
960            list[storeIndex] = tmp;
961
962            return storeIndex;
963        }
964        private static int partition(float[] list, int left, int right, int pivotIndex)
965        {
966           
967            float pivotValue = list[pivotIndex];
968           
969            float tmp = list[right];
970            list[right] = list[pivotIndex];
971            list[pivotIndex] = tmp; // Move pivot to end
972            int storeIndex = left;
973            for (int i = left; i <= right; i++)
974            {
975                if (list[i] < pivotValue)
976                {
977                    tmp = list[storeIndex];
978                    list[storeIndex] = list[i];
979                    list[i] = tmp;
980                    storeIndex++;
981                }
982            }
983            // Move pivot to its final place
984            tmp = list[right];
985            list[right] = list[storeIndex];
986            list[storeIndex] = tmp;
987
988            return storeIndex;
989        }
990        private static int partition(complex[] list, int left, int right, int pivotIndex)
991        {
992           
993            complex pivotValue = list[pivotIndex];
994           
995            complex tmp = list[right];
996            list[right] = list[pivotIndex];
997            list[pivotIndex] = tmp; // Move pivot to end
998            int storeIndex = left;
999            for (int i = left; i <= right; i++)
1000            {
1001                if (list[i] < pivotValue)
1002                {
1003                    tmp = list[storeIndex];
1004                    list[storeIndex] = list[i];
1005                    list[i] = tmp;
1006                    storeIndex++;
1007                }
1008            }
1009            // Move pivot to its final place
1010            tmp = list[right];
1011            list[right] = list[storeIndex];
1012            list[storeIndex] = tmp;
1013
1014            return storeIndex;
1015        }
1016
1017#endregion HYCALPER AUTO GENERATED CODE
1018
1019
1020        /// <summary>
1021        /// Quick select algorithm: Find the k-th smallest element in list.
1022        /// Will change the list parameter!
1023        /// </summary>
1024        /// <remarks><para>Elements in the array list will be reordered. Make sure to pass a copy if you intend to use that data later</para></remarks>
1025        /// <param name="list">The list to search in</param>
1026        /// <param name="left">The first index in the list to start the search</param>
1027        /// <param name="right">The last index in the list to end the search</param>
1028        /// <param name="k">The k-th smallest element to find in list[left:right]. If k is smaller than 1 or larger than the number of elements the smallest/largest value will be returned.</param>
1029        /// <param name="position">Returns the index in list where the smallest element was found</param>
1030        /// <returns>The k-th smallest element</returns>
1031        private static double quickselect_worker(double[] list, int left, int right, int k, out int position)
1032        {
1033            if ((left < 0) || (right > list.Length - 1))
1034                throw new Exception("Arguments out of range. Left and right must be within the array limits");
1035            position = -1;
1036            while (true)
1037            {
1038                if (left == right) // If the list contains only one element
1039                {
1040                    position = left;
1041                    return list[left];  // Return that element
1042                }
1043                // select pivotIndex between left and right
1044                int pivotIndex = (left + right) / 2;
1045                position = partition(list, left, right, pivotIndex); // = new pivot index
1046                int pivotDist = position - left + 1;
1047                // The pivot is in its final sorted position,
1048                // so pivotDist reflects its 1-based position if list were sorted
1049                if (pivotDist == k)
1050                    return list[position];
1051                else if (k < pivotDist)
1052                {
1053                    //return quickselect(list, left, pivotNewIndex - 1, k);
1054                    right = position - 1;
1055                }
1056                else
1057                {
1058                    //  return quickselect(list, pivotNewIndex + 1, right, k - pivotDist);
1059                    left = position + 1;
1060                    k = k - pivotDist;
1061                }
1062            }
1063        }
1064
1065#region HYCALPER AUTO GENERATED CODE
1066
1067        /// <summary>
1068        /// Quick select algorithm: Find the k-th smallest element in list.
1069        /// Will change the list parameter!
1070        /// </summary>
1071        /// <remarks><para>Elements in the array list will be reordered. Make sure to pass a copy if you intend to use that data later</para></remarks>
1072        /// <param name="list">The list to search in</param>
1073        /// <param name="left">The first index in the list to start the search</param>
1074        /// <param name="right">The last index in the list to end the search</param>
1075        /// <param name="k">The k-th smallest element to find in list[left:right]. If k is smaller than 1 or larger than the number of elements the smallest/largest value will be returned.</param>
1076        /// <param name="position">Returns the index in list where the smallest element was found</param>
1077        /// <returns>The k-th smallest element</returns>
1078        private static Int64 quickselect_worker(Int64[] list, int left, int right, int k, out int position)
1079        {
1080            if ((left < 0) || (right > list.Length - 1))
1081                throw new Exception("Arguments out of range. Left and right must be within the array limits");
1082            position = -1;
1083            while (true)
1084            {
1085                if (left == right) // If the list contains only one element
1086                {
1087                    position = left;
1088                    return list[left];  // Return that element
1089                }
1090                // select pivotIndex between left and right
1091                int pivotIndex = (left + right) / 2;
1092                position = partition(list, left, right, pivotIndex); // = new pivot index
1093                int pivotDist = position - left + 1;
1094                // The pivot is in its final sorted position,
1095                // so pivotDist reflects its 1-based position if list were sorted
1096                if (pivotDist == k)
1097                    return list[position];
1098                else if (k < pivotDist)
1099                {
1100                    //return quickselect(list, left, pivotNewIndex - 1, k);
1101                    right = position - 1;
1102                }
1103                else
1104                {
1105                    //  return quickselect(list, pivotNewIndex + 1, right, k - pivotDist);
1106                    left = position + 1;
1107                    k = k - pivotDist;
1108                }
1109            }
1110        }
1111        /// <summary>
1112        /// Quick select algorithm: Find the k-th smallest element in list.
1113        /// Will change the list parameter!
1114        /// </summary>
1115        /// <remarks><para>Elements in the array list will be reordered. Make sure to pass a copy if you intend to use that data later</para></remarks>
1116        /// <param name="list">The list to search in</param>
1117        /// <param name="left">The first index in the list to start the search</param>
1118        /// <param name="right">The last index in the list to end the search</param>
1119        /// <param name="k">The k-th smallest element to find in list[left:right]. If k is smaller than 1 or larger than the number of elements the smallest/largest value will be returned.</param>
1120        /// <param name="position">Returns the index in list where the smallest element was found</param>
1121        /// <returns>The k-th smallest element</returns>
1122        private static Int32 quickselect_worker(Int32[] list, int left, int right, int k, out int position)
1123        {
1124            if ((left < 0) || (right > list.Length - 1))
1125                throw new Exception("Arguments out of range. Left and right must be within the array limits");
1126            position = -1;
1127            while (true)
1128            {
1129                if (left == right) // If the list contains only one element
1130                {
1131                    position = left;
1132                    return list[left];  // Return that element
1133                }
1134                // select pivotIndex between left and right
1135                int pivotIndex = (left + right) / 2;
1136                position = partition(list, left, right, pivotIndex); // = new pivot index
1137                int pivotDist = position - left + 1;
1138                // The pivot is in its final sorted position,
1139                // so pivotDist reflects its 1-based position if list were sorted
1140                if (pivotDist == k)
1141                    return list[position];
1142                else if (k < pivotDist)
1143                {
1144                    //return quickselect(list, left, pivotNewIndex - 1, k);
1145                    right = position - 1;
1146                }
1147                else
1148                {
1149                    //  return quickselect(list, pivotNewIndex + 1, right, k - pivotDist);
1150                    left = position + 1;
1151                    k = k - pivotDist;
1152                }
1153            }
1154        }
1155        /// <summary>
1156        /// Quick select algorithm: Find the k-th smallest element in list.
1157        /// Will change the list parameter!
1158        /// </summary>
1159        /// <remarks><para>Elements in the array list will be reordered. Make sure to pass a copy if you intend to use that data later</para></remarks>
1160        /// <param name="list">The list to search in</param>
1161        /// <param name="left">The first index in the list to start the search</param>
1162        /// <param name="right">The last index in the list to end the search</param>
1163        /// <param name="k">The k-th smallest element to find in list[left:right]. If k is smaller than 1 or larger than the number of elements the smallest/largest value will be returned.</param>
1164        /// <param name="position">Returns the index in list where the smallest element was found</param>
1165        /// <returns>The k-th smallest element</returns>
1166        private static byte quickselect_worker(byte[] list, int left, int right, int k, out int position)
1167        {
1168            if ((left < 0) || (right > list.Length - 1))
1169                throw new Exception("Arguments out of range. Left and right must be within the array limits");
1170            position = -1;
1171            while (true)
1172            {
1173                if (left == right) // If the list contains only one element
1174                {
1175                    position = left;
1176                    return list[left];  // Return that element
1177                }
1178                // select pivotIndex between left and right
1179                int pivotIndex = (left + right) / 2;
1180                position = partition(list, left, right, pivotIndex); // = new pivot index
1181                int pivotDist = position - left + 1;
1182                // The pivot is in its final sorted position,
1183                // so pivotDist reflects its 1-based position if list were sorted
1184                if (pivotDist == k)
1185                    return list[position];
1186                else if (k < pivotDist)
1187                {
1188                    //return quickselect(list, left, pivotNewIndex - 1, k);
1189                    right = position - 1;
1190                }
1191                else
1192                {
1193                    //  return quickselect(list, pivotNewIndex + 1, right, k - pivotDist);
1194                    left = position + 1;
1195                    k = k - pivotDist;
1196                }
1197            }
1198        }
1199        /// <summary>
1200        /// Quick select algorithm: Find the k-th smallest element in list.
1201        /// Will change the list parameter!
1202        /// </summary>
1203        /// <remarks><para>Elements in the array list will be reordered. Make sure to pass a copy if you intend to use that data later</para></remarks>
1204        /// <param name="list">The list to search in</param>
1205        /// <param name="left">The first index in the list to start the search</param>
1206        /// <param name="right">The last index in the list to end the search</param>
1207        /// <param name="k">The k-th smallest element to find in list[left:right]. If k is smaller than 1 or larger than the number of elements the smallest/largest value will be returned.</param>
1208        /// <param name="position">Returns the index in list where the smallest element was found</param>
1209        /// <returns>The k-th smallest element</returns>
1210        private static fcomplex quickselect_worker(fcomplex[] list, int left, int right, int k, out int position)
1211        {
1212            if ((left < 0) || (right > list.Length - 1))
1213                throw new Exception("Arguments out of range. Left and right must be within the array limits");
1214            position = -1;
1215            while (true)
1216            {
1217                if (left == right) // If the list contains only one element
1218                {
1219                    position = left;
1220                    return list[left];  // Return that element
1221                }
1222                // select pivotIndex between left and right
1223                int pivotIndex = (left + right) / 2;
1224                position = partition(list, left, right, pivotIndex); // = new pivot index
1225                int pivotDist = position - left + 1;
1226                // The pivot is in its final sorted position,
1227                // so pivotDist reflects its 1-based position if list were sorted
1228                if (pivotDist == k)
1229                    return list[position];
1230                else if (k < pivotDist)
1231                {
1232                    //return quickselect(list, left, pivotNewIndex - 1, k);
1233                    right = position - 1;
1234                }
1235                else
1236                {
1237                    //  return quickselect(list, pivotNewIndex + 1, right, k - pivotDist);
1238                    left = position + 1;
1239                    k = k - pivotDist;
1240                }
1241            }
1242        }
1243        /// <summary>
1244        /// Quick select algorithm: Find the k-th smallest element in list.
1245        /// Will change the list parameter!
1246        /// </summary>
1247        /// <remarks><para>Elements in the array list will be reordered. Make sure to pass a copy if you intend to use that data later</para></remarks>
1248        /// <param name="list">The list to search in</param>
1249        /// <param name="left">The first index in the list to start the search</param>
1250        /// <param name="right">The last index in the list to end the search</param>
1251        /// <param name="k">The k-th smallest element to find in list[left:right]. If k is smaller than 1 or larger than the number of elements the smallest/largest value will be returned.</param>
1252        /// <param name="position">Returns the index in list where the smallest element was found</param>
1253        /// <returns>The k-th smallest element</returns>
1254        private static float quickselect_worker(float[] list, int left, int right, int k, out int position)
1255        {
1256            if ((left < 0) || (right > list.Length - 1))
1257                throw new Exception("Arguments out of range. Left and right must be within the array limits");
1258            position = -1;
1259            while (true)
1260            {
1261                if (left == right) // If the list contains only one element
1262                {
1263                    position = left;
1264                    return list[left];  // Return that element
1265                }
1266                // select pivotIndex between left and right
1267                int pivotIndex = (left + right) / 2;
1268                position = partition(list, left, right, pivotIndex); // = new pivot index
1269                int pivotDist = position - left + 1;
1270                // The pivot is in its final sorted position,
1271                // so pivotDist reflects its 1-based position if list were sorted
1272                if (pivotDist == k)
1273                    return list[position];
1274                else if (k < pivotDist)
1275                {
1276                    //return quickselect(list, left, pivotNewIndex - 1, k);
1277                    right = position - 1;
1278                }
1279                else
1280                {
1281                    //  return quickselect(list, pivotNewIndex + 1, right, k - pivotDist);
1282                    left = position + 1;
1283                    k = k - pivotDist;
1284                }
1285            }
1286        }
1287        /// <summary>
1288        /// Quick select algorithm: Find the k-th smallest element in list.
1289        /// Will change the list parameter!
1290        /// </summary>
1291        /// <remarks><para>Elements in the array list will be reordered. Make sure to pass a copy if you intend to use that data later</para></remarks>
1292        /// <param name="list">The list to search in</param>
1293        /// <param name="left">The first index in the list to start the search</param>
1294        /// <param name="right">The last index in the list to end the search</param>
1295        /// <param name="k">The k-th smallest element to find in list[left:right]. If k is smaller than 1 or larger than the number of elements the smallest/largest value will be returned.</param>
1296        /// <param name="position">Returns the index in list where the smallest element was found</param>
1297        /// <returns>The k-th smallest element</returns>
1298        private static complex quickselect_worker(complex[] list, int left, int right, int k, out int position)
1299        {
1300            if ((left < 0) || (right > list.Length - 1))
1301                throw new Exception("Arguments out of range. Left and right must be within the array limits");
1302            position = -1;
1303            while (true)
1304            {
1305                if (left == right) // If the list contains only one element
1306                {
1307                    position = left;
1308                    return list[left];  // Return that element
1309                }
1310                // select pivotIndex between left and right
1311                int pivotIndex = (left + right) / 2;
1312                position = partition(list, left, right, pivotIndex); // = new pivot index
1313                int pivotDist = position - left + 1;
1314                // The pivot is in its final sorted position,
1315                // so pivotDist reflects its 1-based position if list were sorted
1316                if (pivotDist == k)
1317                    return list[position];
1318                else if (k < pivotDist)
1319                {
1320                    //return quickselect(list, left, pivotNewIndex - 1, k);
1321                    right = position - 1;
1322                }
1323                else
1324                {
1325                    //  return quickselect(list, pivotNewIndex + 1, right, k - pivotDist);
1326                    left = position + 1;
1327                    k = k - pivotDist;
1328                }
1329            }
1330        }
1331
1332#endregion HYCALPER AUTO GENERATED CODE
1333       #endregion
1334    }
1335}
Note: See TracBrowser for help on using the repository browser.