Free cookie consent management tool by TermsFeed Policy Generator

source: branches/HeuristicLab.Problems.GaussianProcessTuning/ILNumerics.2.14.4735.573/Functions/builtin/conj.cs @ 11194

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

#1967: ILNumerics source for experimentation

File size: 11.2 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 ILNumerics;
44using ILNumerics.Exceptions;
45using ILNumerics.Storage;
46using ILNumerics.Misc;
47
48
49namespace ILNumerics {
50    public partial class ILMath {
51
52        #region convenience functions (real arguments)
53        /// <summary>
54        /// Complex conjugate of A
55        /// </summary>
56        /// <param name="A">Input array</param>
57        /// <returns>The array itself</returns>
58        /// <remarks>This overload is provided for convenience only. It eases the implementation of complex functions, where complex conjugate transposes are needed.</remarks>
59        public static ILRetArray<double> conj (ILInArray<double> A) {
60            using (ILScope.Enter(A))
61            return A.C;
62        }
63        /// <summary>
64        /// Complex conjugate of A
65        /// </summary>
66        /// <param name="A">Input array</param>
67        /// <returns>The array itself</returns>
68        /// <remarks>This overload is provided for convenience only. It eases the implementation of complex functions, where complex conjugate transposes are needed.</remarks>
69        public static ILRetArray<float> conj (ILInArray<float> A) {
70            using (ILScope.Enter(A))
71            return A.C;
72        }
73        #endregion convenience functions (real arguments)
74
75        /// <summary>Complex conjugate of array A</summary>
76        /// <param name="A">Input array</param>
77        /// <returns>Complex conjugate of array A</returns>
78        /// <remarks><para>If the input array is empty, an empty array will be returned.</para>
79        /// <para>The array returned will be a dense array.</para></remarks>
80        public unsafe static ILRetArray<complex> conj(ILInArray<complex> A) {
81            using (ILScope.Enter(A)) {
82                if (A.IsEmpty)
83                    return new ILRetArray<complex>(A.Size);
84                ILSize inDim = A.Size;
85                complex[] arrA = A.GetArrayForRead();
86                complex[] retArr;
87                int outLen = inDim.NumberOfElements;
88                bool inplace = true;
89
90                if (!A.TryGetStorage4InplaceOp(out retArr)) {
91                    retArr = ILMemoryPool.Pool.New<complex>(outLen);
92                    inplace = false;
93                }
94                int i = 0, workItemCount = Settings.s_maxNumberThreads, workItemLength, workerCount = 1;
95                if (Settings.s_maxNumberThreads > 1 && outLen / 2 > Settings.s_minParallelElement1Count) {
96                    if (outLen / workItemCount > Settings.s_minParallelElement1Count) {
97                        workItemLength = outLen / workItemCount;
98                        //workItemLength = (int)((double)outLen / workItemCount * 1.05);
99                    } else {
100                        workItemLength = outLen / 2;
101                        workItemCount = 2;
102                    }
103                } else {
104                    workItemLength = outLen;
105                    workItemCount = 1;
106                }
107                ILDenseStorage<complex> retStorage = new ILDenseStorage<complex>(retArr, inDim);
108
109                Action<object> worker = data => {
110                    Tuple<int, int, IntPtr, IntPtr, bool> range = (Tuple<int, int, IntPtr, IntPtr, bool>)data;
111
112                    complex* cp = ((complex*)range.Item4 + range.Item1);
113
114                    complex* cLast = cp + range.Item2;
115                    if (range.Item5) {
116                        // inplace
117                        while (cp < cLast) {
118                            (*cp).imag = (*cp).imag * -1.0;
119                            cp++;
120                        }
121                    } else {
122                        complex* ap = ((complex*)range.Item3 + range.Item1);
123                        while (cp < cLast) {
124                            (*cp).real = (*ap).real; (*cp).imag = (*ap).imag * -1.0;
125                            ap++;
126                            cp++;
127                        }
128                    }
129                    System.Threading.Interlocked.Decrement(ref workerCount);
130                    //retStorage.PendingEvents.Signal();
131                };
132
133                //retStorage.PendingEvents = new System.Threading.CountdownEvent(workItemCount);
134                fixed (complex* arrAP = arrA)
135                fixed (complex* retArrP = retArr) {
136                    for (; i < workItemCount - 1; i++) {
137                        Tuple<int, int, IntPtr, IntPtr, bool> range
138                            = new Tuple<int, int, IntPtr, IntPtr, bool>
139                                (i * workItemLength, workItemLength, (IntPtr)arrAP, (IntPtr)retArrP, inplace);
140                        System.Threading.Interlocked.Increment(ref workerCount);
141                        ILThreadPool.QueueUserWorkItem(i, worker, range);
142                    }
143                    // the last (or may the only) chunk is done right here
144                    worker(new Tuple<int, int, IntPtr, IntPtr, bool>
145                                (i * workItemLength, outLen - i * workItemLength, (IntPtr)arrAP, (IntPtr)retArrP, inplace));
146
147                    ILThreadPool.Wait4Workers(ref workerCount);
148                }
149                return new ILRetArray<complex>(retStorage);
150            }
151        }
152
153        /// <summary>Complex conjugate of array A</summary>
154        /// <param name="A">Input array</param>
155        /// <returns>Complex conjugate of array A</returns>
156        /// <remarks><para>If the input array is empty, an empty array will be returned.</para>
157        /// <para>The array returned will be a dense array.</para></remarks>
158        public unsafe static  ILRetArray<fcomplex>  conj (ILInArray< fcomplex > A) {
159            using (ILScope.Enter(A)) {
160                if (A.IsEmpty)
161                    return new  ILRetArray<fcomplex>(A.Size);
162                ILSize inDim = A.Size;
163                fcomplex[] arrA = A.GetArrayForRead();
164                fcomplex [] retArr;
165                int outLen = inDim.NumberOfElements;
166                bool inplace = true;
167               
168                if (!A.TryGetStorage4InplaceOp(out retArr)) {
169                    retArr = ILMemoryPool.Pool.New<fcomplex>(outLen);
170                    inplace = false;
171                }
172                int i = 0, workItemCount = Settings.s_maxNumberThreads, workItemLength, workerCount = 1;
173                if (Settings.s_maxNumberThreads > 1 && outLen / 2 > Settings.s_minParallelElement1Count) {
174                    if (outLen / workItemCount > Settings.s_minParallelElement1Count) {
175                        workItemLength = outLen / workItemCount;
176                        //workItemLength = (int)((double)outLen / workItemCount * 1.05);
177                    } else {
178                        workItemLength = outLen / 2;
179                        workItemCount = 2;
180                    }
181                } else {
182                    workItemLength = outLen;
183                    workItemCount = 1;
184                }
185                ILDenseStorage<fcomplex> retStorage = new ILDenseStorage<fcomplex>(retArr, inDim);
186               
187                Action<object> worker = data => {
188                    Tuple<int, int, IntPtr, IntPtr, bool> range = (Tuple<int, int, IntPtr, IntPtr, bool>)data;
189                   
190                    fcomplex* cp = ((fcomplex*)range.Item4 + range.Item1);
191                   
192                    fcomplex* cLast = cp + range.Item2;
193                    if (range.Item5) {
194                        // inplace
195                        while (cp < cLast) {
196                            (*cp).imag = (*cp).imag * -1.0f;
197                            cp++;
198                        }
199                    } else {
200                        fcomplex* ap = ((fcomplex*)range.Item3 + range.Item1);
201                        while (cp < cLast) {
202                            (*cp).real = (*ap).real; (*cp).imag = (*ap).imag * -1.0f;
203                            ap++;
204                            cp++;
205                        }
206                    }
207                    System.Threading.Interlocked.Decrement(ref workerCount);
208                    //retStorage.PendingEvents.Signal();
209                };
210
211                //retStorage.PendingEvents = new System.Threading.CountdownEvent(workItemCount);
212                fixed ( fcomplex* arrAP = arrA)
213                fixed ( fcomplex* retArrP = retArr) {
214                    for (; i < workItemCount - 1; i++) {
215                        Tuple<int, int, IntPtr, IntPtr, bool> range
216                            = new Tuple<int, int, IntPtr, IntPtr, bool>
217                                (i * workItemLength, workItemLength, (IntPtr)arrAP, (IntPtr)retArrP, inplace);
218                        System.Threading.Interlocked.Increment(ref workerCount);
219                        ILThreadPool.QueueUserWorkItem(i,worker, range);
220                    }
221                    // the last (or may the only) chunk is done right here
222                    worker(new Tuple<int, int, IntPtr, IntPtr, bool>
223                                (i * workItemLength, outLen - i * workItemLength, (IntPtr)arrAP, (IntPtr)retArrP, inplace));
224
225                    ILThreadPool.Wait4Workers(ref workerCount);
226                }
227                return new  ILRetArray<fcomplex>(retStorage);
228            }
229        }
230 
231    }
232}
Note: See TracBrowser for help on using the repository browser.