Context Navigation

source: branches/HeuristicLab.Problems.GaussianProcessTuning/ILNumerics.2.14.4735.573/Functions/builtin/atan2.cs @ 11316

Visit:

Last change on this file since 11316 was 9102, checked in by gkronber, 12 years ago
#1967: ILNumerics source for experimentation
File size: 87.0 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
40	using System;
41	using ILNumerics;
42	using ILNumerics.Misc;
43	using ILNumerics.Storage;
44	using ILNumerics.Native;
45	using ILNumerics.Exceptions;
46
47
48	namespace ILNumerics {
49
50	public partial class ILMath {
51
52
53
54
55	#region HYCALPER AUTO GENERATED CODE
56
57	/// <summary>Arcus tangens of elements</summary>
58	/// <param name="A">Input array A</param>
59	/// <param name="B">Input array B</param>
60	/// <returns>Array with elementwise arcus tangens of both inputs</returns>
61	/// <remarks><para>On empty input an empty array will be returned.</para>
62	/// <para>A and/or B may be scalar. The scalar value will be applied on all elements of the
63	/// other array.</para>
64	/// <para>If A or B is a colum vector and the other parameter is an array with a matching colum length, the vector is used to operate on all columns of the array.
65	/// Similar, if one parameter is a row vector, it is used to operate along the rows of the other array if its number of columns matches the vector length. This feature
66	/// can be used to replace the (costly) repmat function for most binary operators.</para>
67	/// <para>For all other cases the dimensions of A and B must match.</para></remarks>
68	/// <exception cref="ILNumerics.Exceptions.ILArgumentException">If the size of both arrays does not match any parameter rule.</exception>
69	public unsafe static ILRetArray<float> atan2(ILInArray<float> A, ILInArray<float> B) {
70	using (ILScope.Enter(A,B)) {
71	int outLen;
72	BinOpItMode mode;
73	float [] retArr;
74	float [] arrA = A.GetArrayForRead();
75	float[] arrB = B.GetArrayForRead();
76	ILSize outDims;
77	#region determine operation mode
78	if (A.IsScalar) {
79	outDims = B.Size;
80	if (B.IsScalar) {
81
82	return new ILRetArray<float> (new float [1]{ (float)Math.Atan2 (A.GetValue(0) , B.GetValue(0))}, A.Size);
83	} else if (B.IsEmpty) {
84	return ILRetArray<float>.empty(outDims);
85	} else {
86	outLen = outDims.NumberOfElements;
87	if (!B.TryGetStorage4InplaceOp(out retArr)) {
88	retArr = ILMemoryPool.Pool.New< float > (outLen);
89	mode = BinOpItMode.SAN;
90	} else {
91	mode = BinOpItMode.SAI;
92	}
93	}
94	} else {
95	outDims = A.Size;
96	if (B.IsScalar) {
97	if (A.IsEmpty) {
98	return ILRetArray<float>.empty(A.Size);
99	}
100	outLen = A.S.NumberOfElements;
101	if (!A.TryGetStorage4InplaceOp(out retArr)) {
102	retArr = ILMemoryPool.Pool.New< float > (outLen);
103	mode = BinOpItMode.ASN;
104	} else {
105	mode = BinOpItMode.ASI;
106	}
107	} else {
108	// array + array
109	if (!A.Size.IsSameSize(B.Size)) {
110	return atan2Ex(A,B);
111	}
112	outLen = A.S.NumberOfElements;
113	if (A.TryGetStorage4InplaceOp(out retArr))
114	mode = BinOpItMode.AAIA;
115	else if (B.TryGetStorage4InplaceOp(out retArr))
116	mode = BinOpItMode.AAIB;
117	else {
118	retArr = ILMemoryPool.Pool.New< float > (outLen);
119	mode = BinOpItMode.AAN;
120	}
121	}
122	}
123	#endregion
124	ILDenseStorage<float> retStorage = new ILDenseStorage<float>(retArr, outDims);
125	int i = 0, workerCount = 1;
126	Action<object> worker = data => {
127	Tuple<int, int, IntPtr, IntPtr, IntPtr, BinOpItMode> range
128	= (Tuple<int, int, IntPtr, IntPtr, IntPtr, BinOpItMode>)data;
129	float* cp = (float*)range.Item5 + range.Item1;
130	float scalar;
131	int j = range.Item2;
132	#region loops
133	switch (mode) {
134	case BinOpItMode.AAIA:
135	float* bp = ((float*)range.Item4 + range.Item1);
136	while (j > 20) {
137	cp[0] = (float)Math.Atan2 (cp[0] , bp[0]);
138	cp[1] = (float)Math.Atan2 (cp[1] , bp[1]);
139	cp[2] = (float)Math.Atan2 (cp[2] , bp[2]);
140	cp[3] = (float)Math.Atan2 (cp[3] , bp[3]);
141	cp[4] = (float)Math.Atan2 (cp[4] , bp[4]);
142	cp[5] = (float)Math.Atan2 (cp[5] , bp[5]);
143	cp[6] = (float)Math.Atan2 (cp[6] , bp[6]);
144	cp[7] = (float)Math.Atan2 (cp[7] , bp[7]);
145	cp[8] = (float)Math.Atan2 (cp[8] , bp[8]);
146	cp[9] = (float)Math.Atan2 (cp[9] , bp[9]);
147	cp[10] = (float)Math.Atan2 (cp[10] , bp[10]);
148	cp[11] = (float)Math.Atan2 (cp[11] , bp[11]);
149	cp[12] = (float)Math.Atan2 (cp[12] , bp[12]);
150	cp[13] = (float)Math.Atan2 (cp[13] , bp[13]);
151	cp[14] = (float)Math.Atan2 (cp[14] , bp[14]);
152	cp[15] = (float)Math.Atan2 (cp[15] , bp[15]);
153	cp[16] = (float)Math.Atan2 (cp[16] , bp[16]);
154	cp[17] = (float)Math.Atan2 (cp[17] , bp[17]);
155	cp[18] = (float)Math.Atan2 (cp[18] , bp[18]);
156	cp[19] = (float)Math.Atan2 (cp[19] , bp[19]);
157	cp[20] = (float)Math.Atan2 (cp[20] , bp[20]);
158	cp += 21; bp += 21; j -= 21;
159	}
160	while (j --> 0) {
161
162	cp = (float)Math.Atan2 (cp , *bp);
163	cp++; bp++;
164	}
165	break;
166	case BinOpItMode.AAIB:
167	float* ap = ((float*)range.Item3 + range.Item1);
168	while (j > 20) {
169
170	cp[0] = (float)Math.Atan2 (ap[0] , cp[0]);
171	cp[1] = (float)Math.Atan2 (ap[1] , cp[1]);
172	cp[2] = (float)Math.Atan2 (ap[2] , cp[2]);
173	cp[3] = (float)Math.Atan2 (ap[3] , cp[3]);
174	cp[4] = (float)Math.Atan2 (ap[4] , cp[4]);
175	cp[5] = (float)Math.Atan2 (ap[5] , cp[5]);
176	cp[6] = (float)Math.Atan2 (ap[6] , cp[6]);
177	cp[7] = (float)Math.Atan2 (ap[7] , cp[7]);
178	cp[8] = (float)Math.Atan2 (ap[8] , cp[8]);
179	cp[9] = (float)Math.Atan2 (ap[9] , cp[9]);
180	cp[10] = (float)Math.Atan2 (ap[10] , cp[10]);
181	cp[11] = (float)Math.Atan2 (ap[11] , cp[11]);
182	cp[12] = (float)Math.Atan2 (ap[12] , cp[12]);
183	cp[13] = (float)Math.Atan2 (ap[13] , cp[13]);
184	cp[14] = (float)Math.Atan2 (ap[14] , cp[14]);
185	cp[15] = (float)Math.Atan2 (ap[15] , cp[15]);
186	cp[16] = (float)Math.Atan2 (ap[16] , cp[16]);
187	cp[17] = (float)Math.Atan2 (ap[17] , cp[17]);
188	cp[18] = (float)Math.Atan2 (ap[18] , cp[18]);
189	cp[19] = (float)Math.Atan2 (ap[19] , cp[19]);
190	cp[20] = (float)Math.Atan2 (ap[20] , cp[20]);
191	ap += 21; cp += 21; j -= 21;
192	}
193	while (j --> 0) {
194
195	cp = (float)Math.Atan2 (ap , *cp);
196	ap++; cp++;
197	}
198	break;
199	case BinOpItMode.AAN:
200	ap = ((float*)range.Item3 + range.Item1);
201	bp = ((float*)range.Item4 + range.Item1);
202	while (j > 20) {
203
204	cp[0] = (float)Math.Atan2 (ap[0] , bp[0]);
205	cp[1] = (float)Math.Atan2 (ap[1] , bp[1]);
206	cp[2] = (float)Math.Atan2 (ap[2] , bp[2]);
207	cp[3] = (float)Math.Atan2 (ap[3] , bp[3]);
208	cp[4] = (float)Math.Atan2 (ap[4] , bp[4]);
209	cp[5] = (float)Math.Atan2 (ap[5] , bp[5]);
210	cp[6] = (float)Math.Atan2 (ap[6] , bp[6]);
211	cp[7] = (float)Math.Atan2 (ap[7] , bp[7]);
212	cp[8] = (float)Math.Atan2 (ap[8] , bp[8]);
213	cp[9] = (float)Math.Atan2 (ap[9] , bp[9]);
214	cp[10] = (float)Math.Atan2 (ap[10] , bp[10]);
215	cp[11] = (float)Math.Atan2 (ap[11] , bp[11]);
216	cp[12] = (float)Math.Atan2 (ap[12] , bp[12]);
217	cp[13] = (float)Math.Atan2 (ap[13] , bp[13]);
218	cp[14] = (float)Math.Atan2 (ap[14] , bp[14]);
219	cp[15] = (float)Math.Atan2 (ap[15] , bp[15]);
220	cp[16] = (float)Math.Atan2 (ap[16] , bp[16]);
221	cp[17] = (float)Math.Atan2 (ap[17] , bp[17]);
222	cp[18] = (float)Math.Atan2 (ap[18] , bp[18]);
223	cp[19] = (float)Math.Atan2 (ap[19] , bp[19]);
224	cp[20] = (float)Math.Atan2 (ap[20] , bp[20]);
225	ap+=21; bp+=21; cp+=21; j-=21;
226	}
227	while (j --> 0) {
228
229	cp = (float)Math.Atan2 (ap , *bp);
230	ap++; bp++; cp++;
231	}
232	break;
233	case BinOpItMode.ASI:
234	scalar = ((float)range.Item4);
235	while (j > 20) {
236
237	cp[0] = (float)Math.Atan2 (cp[0] , scalar);
238	cp[1] = (float)Math.Atan2 (cp[1] , scalar);
239	cp[2] = (float)Math.Atan2 (cp[2] , scalar);
240	cp[3] = (float)Math.Atan2 (cp[3] , scalar);
241	cp[4] = (float)Math.Atan2 (cp[4] , scalar);
242	cp[5] = (float)Math.Atan2 (cp[5] , scalar);
243	cp[6] = (float)Math.Atan2 (cp[6] , scalar);
244	cp[7] = (float)Math.Atan2 (cp[7] , scalar);
245	cp[8] = (float)Math.Atan2 (cp[8] , scalar);
246	cp[9] = (float)Math.Atan2 (cp[9] , scalar);
247	cp[10] = (float)Math.Atan2 (cp[10] , scalar);
248	cp[11] = (float)Math.Atan2 (cp[11] , scalar);
249	cp[12] = (float)Math.Atan2 (cp[12] , scalar);
250	cp[13] = (float)Math.Atan2 (cp[13] , scalar);
251	cp[14] = (float)Math.Atan2 (cp[14] , scalar);
252	cp[15] = (float)Math.Atan2 (cp[15] , scalar);
253	cp[16] = (float)Math.Atan2 (cp[16] , scalar);
254	cp[17] = (float)Math.Atan2 (cp[17] , scalar);
255	cp[18] = (float)Math.Atan2 (cp[18] , scalar);
256	cp[19] = (float)Math.Atan2 (cp[19] , scalar);
257	cp[20] = (float)Math.Atan2 (cp[20] , scalar);
258	cp += 21; j -= 21;
259	}
260	while (j --> 0) {
261
262	cp = (float)Math.Atan2 (cp , scalar);
263	cp++;
264	}
265	break;
266	case BinOpItMode.ASN:
267	ap = ((float*)range.Item3 + range.Item1);
268	scalar = ((float)range.Item4);
269	while (j > 20) {
270
271	cp[0] = (float)Math.Atan2 (ap[0] , scalar);
272	cp[1] = (float)Math.Atan2 (ap[1] , scalar);
273	cp[2] = (float)Math.Atan2 (ap[2] , scalar);
274	cp[3] = (float)Math.Atan2 (ap[3] , scalar);
275	cp[4] = (float)Math.Atan2 (ap[4] , scalar);
276	cp[5] = (float)Math.Atan2 (ap[5] , scalar);
277	cp[6] = (float)Math.Atan2 (ap[6] , scalar);
278	cp[7] = (float)Math.Atan2 (ap[7] , scalar);
279	cp[8] = (float)Math.Atan2 (ap[8] , scalar);
280	cp[9] = (float)Math.Atan2 (ap[9] , scalar);
281	cp[10] = (float)Math.Atan2 (ap[10] , scalar);
282	cp[11] = (float)Math.Atan2 (ap[11] , scalar);
283	cp[12] = (float)Math.Atan2 (ap[12] , scalar);
284	cp[13] = (float)Math.Atan2 (ap[13] , scalar);
285	cp[14] = (float)Math.Atan2 (ap[14] , scalar);
286	cp[15] = (float)Math.Atan2 (ap[15] , scalar);
287	cp[16] = (float)Math.Atan2 (ap[16] , scalar);
288	cp[17] = (float)Math.Atan2 (ap[17] , scalar);
289	cp[18] = (float)Math.Atan2 (ap[18] , scalar);
290	cp[19] = (float)Math.Atan2 (ap[19] , scalar);
291	cp[20] = (float)Math.Atan2 (ap[20] , scalar);
292	ap+=21; cp+=21; j -= 21;
293	}
294	while (j --> 0) {
295
296	cp = (float)Math.Atan2 (ap , scalar);
297	ap++; cp++;
298	}
299	break;
300	case BinOpItMode.SAI:
301	scalar = ((float)range.Item3);
302	while (j > 20) {
303
304	cp[0] = (float)Math.Atan2 (scalar , cp[0]);
305	cp[1] = (float)Math.Atan2 (scalar , cp[1]);
306	cp[2] = (float)Math.Atan2 (scalar , cp[2]);
307	cp[3] = (float)Math.Atan2 (scalar , cp[3]);
308	cp[4] = (float)Math.Atan2 (scalar , cp[4]);
309	cp[5] = (float)Math.Atan2 (scalar , cp[5]);
310	cp[6] = (float)Math.Atan2 (scalar , cp[6]);
311	cp[7] = (float)Math.Atan2 (scalar , cp[7]);
312	cp[8] = (float)Math.Atan2 (scalar , cp[8]);
313	cp[9] = (float)Math.Atan2 (scalar , cp[9]);
314	cp[10] = (float)Math.Atan2 (scalar , cp[10]);
315	cp[11] = (float)Math.Atan2 (scalar , cp[11]);
316	cp[12] = (float)Math.Atan2 (scalar , cp[12]);
317	cp[13] = (float)Math.Atan2 (scalar , cp[13]);
318	cp[14] = (float)Math.Atan2 (scalar , cp[14]);
319	cp[15] = (float)Math.Atan2 (scalar , cp[15]);
320	cp[16] = (float)Math.Atan2 (scalar , cp[16]);
321	cp[17] = (float)Math.Atan2 (scalar , cp[17]);
322	cp[18] = (float)Math.Atan2 (scalar , cp[18]);
323	cp[19] = (float)Math.Atan2 (scalar , cp[19]);
324	cp[20] = (float)Math.Atan2 (scalar , cp[20]);
325	cp += 21; j -= 21;
326	}
327	while (j --> 0) {
328
329	cp = (float)Math.Atan2 (scalar , cp);
330	cp++;
331	}
332	break;
333	case BinOpItMode.SAN:
334	scalar = ((float)range.Item3);
335	bp = ((float*)range.Item4 + range.Item1);
336	while (j > 20) {
337
338	cp[0] = (float)Math.Atan2 (scalar , bp[0]);
339	cp[1] = (float)Math.Atan2 (scalar , bp[1]);
340	cp[2] = (float)Math.Atan2 (scalar , bp[2]);
341	cp[3] = (float)Math.Atan2 (scalar , bp[3]);
342	cp[4] = (float)Math.Atan2 (scalar , bp[4]);
343	cp[5] = (float)Math.Atan2 (scalar , bp[5]);
344	cp[6] = (float)Math.Atan2 (scalar , bp[6]);
345	cp[7] = (float)Math.Atan2 (scalar , bp[7]);
346	cp[8] = (float)Math.Atan2 (scalar , bp[8]);
347	cp[9] = (float)Math.Atan2 (scalar , bp[9]);
348	cp[10] = (float)Math.Atan2 (scalar , bp[10]);
349	cp[11] = (float)Math.Atan2 (scalar , bp[11]);
350	cp[12] = (float)Math.Atan2 (scalar , bp[12]);
351	cp[13] = (float)Math.Atan2 (scalar , bp[13]);
352	cp[14] = (float)Math.Atan2 (scalar , bp[14]);
353	cp[15] = (float)Math.Atan2 (scalar , bp[15]);
354	cp[16] = (float)Math.Atan2 (scalar , bp[16]);
355	cp[17] = (float)Math.Atan2 (scalar , bp[17]);
356	cp[18] = (float)Math.Atan2 (scalar , bp[18]);
357	cp[19] = (float)Math.Atan2 (scalar , bp[19]);
358	cp[20] = (float)Math.Atan2 (scalar , bp[20]);
359	bp+=21; cp+=21; j -= 21;
360	}
361	while (j --> 0) {
362
363	cp = (float)Math.Atan2 (scalar , bp);
364	bp++; cp++;
365	}
366	break;
367	default:
368	break;
369	}
370	#endregion
371	System.Threading.Interlocked.Decrement(ref workerCount);
372	//retStorage.PendingEvents.Signal();
373	};
374
375	#region do the work
376	int workItemCount = Settings.s_maxNumberThreads, workItemLength;
377	if (Settings.s_maxNumberThreads > 1 && outLen / 2 > Settings.s_minParallelElement1Count) {
378	if (outLen / workItemCount > Settings.s_minParallelElement1Count) {
379	workItemLength = outLen / workItemCount;
380	//workItemLength = (int)((double)outLen / workItemCount * 1.05);
381	} else {
382	workItemLength = outLen / 2;
383	workItemCount = 2;
384	}
385	} else {
386	workItemLength = outLen;
387	workItemCount = 1;
388	}
389
390	fixed ( float* arrAP = arrA)
391	fixed ( float* arrBP = arrB)
392	fixed ( float* retArrP = retArr) {
393
394	for (; i < workItemCount - 1; i++) {
395	Tuple<int, int, IntPtr, IntPtr, IntPtr, BinOpItMode> range
396	= new Tuple<int, int, IntPtr, IntPtr, IntPtr, BinOpItMode>
397	(i * workItemLength, workItemLength, (IntPtr)arrAP, (IntPtr)arrBP, (IntPtr)retArrP, mode);
398	System.Threading.Interlocked.Increment(ref workerCount);
399	ILThreadPool.QueueUserWorkItem(i, worker, range);
400	}
401	// the last (or may the only) chunk is done right here
402	//System.Threading.Interlocked.Increment(ref retStorage.PendingTasks);
403	worker(new Tuple<int, int, IntPtr, IntPtr, IntPtr, BinOpItMode>
404	(i * workItemLength, outLen - i * workItemLength, (IntPtr)arrAP, (IntPtr)arrBP, (IntPtr)retArrP, mode));
405
406	ILThreadPool.Wait4Workers(ref workerCount);
407	}
408
409	#endregion
410	return new ILRetArray< float>(retStorage);
411	}
412	}
413
414	private static unsafe ILRetArray<float> atan2Ex(ILInArray<float> A, ILInArray<float> B) {
415	//using (ILScope.Enter(A, B)) { we cannot start a new scope here, since this would prevent A and B to be used implace if applicable
416
417	#region parameter checking
418	if (isnull(A) \|\| isnull(B))
419	return empty<float>(ILSize.Empty00);
420	if (A.IsEmpty) {
421	return empty<float>(B.S);
422	} else if (B.IsEmpty) {
423	return empty<float>(A.S);
424	}
425	//if (A.IsScalar \|\| B.IsScalar \|\| A.D.IsSameSize(B.D))
426	// return add(A,B);
427	int dim = -1;
428	for (int l = 0; l < Math.Max(A.S.NumberOfDimensions, B.S.NumberOfDimensions); l++) {
429	if (A.S[l] != B.S[l]) {
430	if (dim >= 0 \|\| (A.S[l] != 1 && B.S[l] != 1)) {
431	throw new ILArgumentException("A and B must have the same size except for one simgleton dimension in A or B");
432	}
433	dim = l;
434	}
435	}
436	if (dim > 1)
437	throw new ILArgumentException("singleton dimension expansion currently is only supported for colum- and row vectors");
438	dim = -(dim - 1); // 0 -> 1, 1 -> 0
439	#endregion
440
441	#region parameter preparation
442
443	float[] retArr;
444
445	float[] arrA = A.GetArrayForRead();
446
447	float[] arrB = B.GetArrayForRead();
448	ILSize outDims;
449	BinOptItExMode mode;
450	int workItemMultiplierLenA;
451	int workItemMultiplierLenB;
452	if (A.IsVector) {
453	outDims = B.S;
454	if (!B.TryGetStorage4InplaceOp(out retArr)) {
455	retArr = ILMemoryPool.Pool.New<float>(outDims.NumberOfElements);
456	mode = BinOptItExMode.VAN;
457	} else {
458	mode = BinOptItExMode.VAI;
459	}
460	workItemMultiplierLenB = outDims[0];
461	workItemMultiplierLenA = dim; // 0 for column, 1 for row vector
462	} else if (B.IsVector) {
463	outDims = A.S;
464	if (!A.TryGetStorage4InplaceOp(out retArr)) {
465	retArr = ILMemoryPool.Pool.New<float>(outDims.NumberOfElements);
466	mode = BinOptItExMode.AVN;
467	} else {
468	mode = BinOptItExMode.AVI;
469	}
470	workItemMultiplierLenB = dim; // 0 for column, 1 for row vector
471	workItemMultiplierLenA = outDims[0];
472	} else {
473	throw new ILArgumentException("A and B must have the same size except for one singleton dimension in either A or B");
474	}
475	int itLen = outDims[0]; // (dim == 0) ? outDims.SequentialIndexDistance(1) : outDims.SequentialIndexDistance(0);
476	#endregion
477
478	#region worker loops definition
479	ILDenseStorage<float> retStorage = new ILDenseStorage<float>(retArr, outDims);
480	int workerCount = 1;
481	Action<object> worker = data => {
482	// expects: iStart, iLen, ap, bp, cp
483	Tuple<int, IntPtr, IntPtr, IntPtr> range =
484	(Tuple<int, IntPtr, IntPtr, IntPtr>)data;
485
486	float* ap;
487
488	float* bp;
489
490	float* cp;
491	switch (mode) {
492	case BinOptItExMode.VAN:
493	if (dim == 0) {
494	bp = (float*)range.Item3;
495	cp = (float*)range.Item4;
496	for (int s = 0; s < range.Item1; s++) {
497	ap = (float*)range.Item2;
498	int l = itLen;
499	while (l > 10) {
500	cp[0] = (float)Math.Atan2 (ap[0] , bp[0]);
501	cp[1] = (float)Math.Atan2 (ap[1] , bp[1]);
502	cp[2] = (float)Math.Atan2 (ap[2] , bp[2]);
503	cp[3] = (float)Math.Atan2 (ap[3] , bp[3]);
504	cp[4] = (float)Math.Atan2 (ap[4] , bp[4]);
505	cp[5] = (float)Math.Atan2 (ap[5] , bp[5]);
506	cp[6] = (float)Math.Atan2 (ap[6] , bp[6]);
507	cp[7] = (float)Math.Atan2 (ap[7] , bp[7]);
508	cp[8] = (float)Math.Atan2 (ap[8] , bp[8]);
509	cp[9] = (float)Math.Atan2 (ap[9] , bp[9]);
510	cp[10] = (float)Math.Atan2 (ap[10] , bp[10]);
511	ap += 11;
512	bp += 11;
513	cp += 11;
514	l -= 11;
515	}
516	while (l-- > 0) {
517	cp++ = (float)Math.Atan2 (ap++ , *bp++);
518	}
519	}
520	} else {
521	// dim == 1
522	ap = (float*)range.Item2;
523	bp = (float*)range.Item3;
524	cp = (float*)range.Item4;
525	for (int s = 0; s < range.Item1; s++) {
526	float val = *ap++;
527	int l = itLen;
528	while (l > 10) {
529	cp[0] = (float)Math.Atan2 (val , bp[0]);
530	cp[1] = (float)Math.Atan2 (val , bp[1]);
531	cp[2] = (float)Math.Atan2 (val , bp[2]);
532	cp[3] = (float)Math.Atan2 (val , bp[3]);
533	cp[4] = (float)Math.Atan2 (val , bp[4]);
534	cp[5] = (float)Math.Atan2 (val , bp[5]);
535	cp[6] = (float)Math.Atan2 (val , bp[6]);
536	cp[7] = (float)Math.Atan2 (val , bp[7]);
537	cp[8] = (float)Math.Atan2 (val , bp[8]);
538	cp[9] = (float)Math.Atan2 (val , bp[9]);
539	cp[10] = (float)Math.Atan2 (val , bp[10]);
540	bp += 11;
541	cp += 11;
542	l -= 11;
543	}
544	while (l-- > 0) {
545	cp++ = (float)Math.Atan2 (val , bp++);
546	}
547	}
548	}
549	break;
550	case BinOptItExMode.VAI:
551	if (dim == 0) {
552	cp = (float*)range.Item4;
553	for (int s = 0; s < range.Item1; s++) {
554	ap = (float*)range.Item2;
555	int l = itLen;
556	while (l > 10) {
557	cp[0] = (float)Math.Atan2 (ap[0] , cp[0]);
558	cp[1] = (float)Math.Atan2 (ap[1] , cp[1]);
559	cp[2] = (float)Math.Atan2 (ap[2] , cp[2]);
560	cp[3] = (float)Math.Atan2 (ap[3] , cp[3]);
561	cp[4] = (float)Math.Atan2 (ap[4] , cp[4]);
562	cp[5] = (float)Math.Atan2 (ap[5] , cp[5]);
563	cp[6] = (float)Math.Atan2 (ap[6] , cp[6]);
564	cp[7] = (float)Math.Atan2 (ap[7] , cp[7]);
565	cp[8] = (float)Math.Atan2 (ap[8] , cp[8]);
566	cp[9] = (float)Math.Atan2 (ap[9] , cp[9]);
567	cp[10] = (float)Math.Atan2 (ap[10] , cp[10]);
568	ap += 11;
569	cp += 11;
570	l -= 11;
571	}
572	while (l-- > 0) {
573	cp = (float)Math.Atan2 (ap++ , *cp);
574	cp++;
575	}
576	}
577	} else {
578	// dim == 1
579	cp = (float*)range.Item4;
580	ap = (float*)range.Item2;
581	for (int s = 0; s < range.Item1; s++) {
582
583	float val = *ap++;
584	int l = itLen;
585	while (l > 10) {
586	cp[0] = (float)Math.Atan2 (val , cp[0]);
587	cp[1] = (float)Math.Atan2 (val , cp[1]);
588	cp[2] = (float)Math.Atan2 (val , cp[2]);
589	cp[3] = (float)Math.Atan2 (val , cp[3]);
590	cp[4] = (float)Math.Atan2 (val , cp[4]);
591	cp[5] = (float)Math.Atan2 (val , cp[5]);
592	cp[6] = (float)Math.Atan2 (val , cp[6]);
593	cp[7] = (float)Math.Atan2 (val , cp[7]);
594	cp[8] = (float)Math.Atan2 (val , cp[8]);
595	cp[9] = (float)Math.Atan2 (val , cp[9]);
596	cp[10] = (float)Math.Atan2 (val , cp[10]);
597	cp += 11;
598	l -= 11;
599	}
600	while (l-- > 0) {
601	cp = (float)Math.Atan2 (val , cp);
602	cp++;
603	}
604	}
605	}
606	break;
607	case BinOptItExMode.AVN:
608	if (dim == 0) {
609	ap = (float*)range.Item2;
610	cp = (float*)range.Item4;
611	for (int s = 0; s < range.Item1; s++) {
612	bp = (float*)range.Item3;
613	int l = itLen;
614	while (l > 10) {
615	cp[0] = (float)Math.Atan2 (ap[0] , bp[0]);
616	cp[1] = (float)Math.Atan2 (ap[1] , bp[1]);
617	cp[2] = (float)Math.Atan2 (ap[2] , bp[2]);
618	cp[3] = (float)Math.Atan2 (ap[3] , bp[3]);
619	cp[4] = (float)Math.Atan2 (ap[4] , bp[4]);
620	cp[5] = (float)Math.Atan2 (ap[5] , bp[5]);
621	cp[6] = (float)Math.Atan2 (ap[6] , bp[6]);
622	cp[7] = (float)Math.Atan2 (ap[7] , bp[7]);
623	cp[8] = (float)Math.Atan2 (ap[8] , bp[8]);
624	cp[9] = (float)Math.Atan2 (ap[9] , bp[9]);
625	cp[10] = (float)Math.Atan2 (ap[10] , bp[10]);
626	ap += 11;
627	bp += 11;
628	cp += 11;
629	l -= 11;
630	}
631	while (l-- > 0) {
632	cp = (float)Math.Atan2 (ap , *bp);
633	ap++;
634	bp++;
635	cp++;
636	}
637	}
638	} else {
639	// dim = 1
640	ap = (float*)range.Item2;
641	bp = (float*)range.Item3;
642	cp = (float*)range.Item4;
643	for (int s = 0; s < range.Item1; s++) {
644	float val = *bp++;
645	int l = itLen;
646	while (l > 10) {
647	cp[0] = (float)Math.Atan2 (ap[0] , val);
648	cp[1] = (float)Math.Atan2 (ap[1] , val);
649	cp[2] = (float)Math.Atan2 (ap[2] , val);
650	cp[3] = (float)Math.Atan2 (ap[3] , val);
651	cp[4] = (float)Math.Atan2 (ap[4] , val);
652	cp[5] = (float)Math.Atan2 (ap[5] , val);
653	cp[6] = (float)Math.Atan2 (ap[6] , val);
654	cp[7] = (float)Math.Atan2 (ap[7] , val);
655	cp[8] = (float)Math.Atan2 (ap[8] , val);
656	cp[9] = (float)Math.Atan2 (ap[9] , val);
657	cp[10] = (float)Math.Atan2 (ap[10] , val);
658	ap += 11;
659	cp += 11;
660	l -= 11;
661	}
662	while (l-- > 0) {
663	cp = (float)Math.Atan2 (ap , val);
664	ap++;
665	cp++;
666	}
667	}
668	}
669	break;
670	case BinOptItExMode.AVI:
671	if (dim == 0) {
672	cp = (float*)range.Item4;
673	for (int s = 0; s < range.Item1; s++) {
674	bp = (float*)range.Item3;
675	int l = itLen;
676	while (l > 10) {
677	cp[0] = (float)Math.Atan2 (cp[0] , bp[0]);
678	cp[1] = (float)Math.Atan2 (cp[1] , bp[1]);
679	cp[2] = (float)Math.Atan2 (cp[2] , bp[2]);
680	cp[3] = (float)Math.Atan2 (cp[3] , bp[3]);
681	cp[4] = (float)Math.Atan2 (cp[4] , bp[4]);
682	cp[5] = (float)Math.Atan2 (cp[5] , bp[5]);
683	cp[6] = (float)Math.Atan2 (cp[6] , bp[6]);
684	cp[7] = (float)Math.Atan2 (cp[7] , bp[7]);
685	cp[8] = (float)Math.Atan2 (cp[8] , bp[8]);
686	cp[9] = (float)Math.Atan2 (cp[9] , bp[9]);
687	cp[10] = (float)Math.Atan2 (cp[10] , bp[10]);
688	bp += 11;
689	cp += 11;
690	l -= 11;
691	}
692	while (l-- > 0) {
693	cp = (float)Math.Atan2 (cp , *bp);
694	bp++;
695	cp++;
696	}
697	}
698	} else {
699	// dim = 1
700	bp = (float*)range.Item3;
701	cp = (float*)range.Item4;
702	for (int s = 0; s < range.Item1; s++) {
703
704	float val = *bp++;
705	int l = itLen;
706	while (l > 10) {
707	cp[0] = (float)Math.Atan2 (cp[0] , val);
708	cp[1] = (float)Math.Atan2 (cp[1] , val);
709	cp[2] = (float)Math.Atan2 (cp[2] , val);
710	cp[3] = (float)Math.Atan2 (cp[3] , val);
711	cp[4] = (float)Math.Atan2 (cp[4] , val);
712	cp[5] = (float)Math.Atan2 (cp[5] , val);
713	cp[6] = (float)Math.Atan2 (cp[6] , val);
714	cp[7] = (float)Math.Atan2 (cp[7] , val);
715	cp[8] = (float)Math.Atan2 (cp[8] , val);
716	cp[9] = (float)Math.Atan2 (cp[9] , val);
717	cp[10] = (float)Math.Atan2 (cp[10] , val);
718	cp += 11;
719	l -= 11;
720	}
721	while (l --> 0) {
722	cp = (float)Math.Atan2 (cp , val);
723	cp++;
724	}
725	}
726	}
727	break;
728	}
729	System.Threading.Interlocked.Decrement(ref workerCount);
730	};
731	#endregion
732
733	#region work distribution
734	int i = 0, workItemCount = Settings.s_maxNumberThreads, workItemLength;
735	if (Settings.s_maxNumberThreads > 1 && outDims.NumberOfElements >= Settings.s_minParallelElement1Count
736	&& outDims[1] > 1) {
737	if (outDims[1] > workItemCount) {
738	workItemLength = outDims[1] / workItemCount;
739	} else {
740	workItemLength = outDims[1] / 2;
741	workItemCount = 2;
742	}
743	} else {
744	workItemLength = outDims[1];
745	workItemCount = 1;
746	}
747
748	fixed ( float* arrAP = arrA)
749	fixed ( float* arrBP = arrB)
750	fixed ( float* retArrP = retArr) {
751
752	for (; i < workItemCount - 1; i++) {
753	Tuple<int, IntPtr, IntPtr, IntPtr> range = new Tuple<int, IntPtr, IntPtr, IntPtr>
754	(workItemLength
755	, (IntPtr)(arrAP + i * workItemMultiplierLenA * workItemLength)
756	, (IntPtr)(arrBP + i * workItemMultiplierLenB * workItemLength)
757	, (IntPtr)(retArrP + i * outDims[0] * workItemLength));
758	System.Threading.Interlocked.Increment(ref workerCount);
759	ILThreadPool.QueueUserWorkItem(i, worker, range);
760	}
761	// the last (or may the only) chunk is done right here
762	//System.Threading.Interlocked.Increment(ref retStorage.PendingTasks);
763	worker(new Tuple<int, IntPtr, IntPtr, IntPtr>
764	(outDims[1] - i * workItemLength
765	, (IntPtr)(arrAP + i * workItemMultiplierLenA * workItemLength)
766	, (IntPtr)(arrBP + i * workItemMultiplierLenB * workItemLength)
767	, (IntPtr)(retArrP + i * outDims[0] * workItemLength)));
768
769	ILThreadPool.Wait4Workers(ref workerCount);
770	}
771	#endregion
772
773	return new ILRetArray<float>(retStorage);
774	//} // no scopes here! it disables implace operations
775	}
776
777	/// <summary>Arcus tangens of elements</summary>
778	/// <param name="A">Input array A</param>
779	/// <param name="B">Input array B</param>
780	/// <returns>Array with elementwise arcus tangens of both inputs</returns>
781	/// <remarks><para>On empty input an empty array will be returned.</para>
782	/// <para>A and/or B may be scalar. The scalar value will be applied on all elements of the
783	/// other array.</para>
784	/// <para>If A or B is a colum vector and the other parameter is an array with a matching colum length, the vector is used to operate on all columns of the array.
785	/// Similar, if one parameter is a row vector, it is used to operate along the rows of the other array if its number of columns matches the vector length. This feature
786	/// can be used to replace the (costly) repmat function for most binary operators.</para>
787	/// <para>For all other cases the dimensions of A and B must match.</para></remarks>
788	/// <exception cref="ILNumerics.Exceptions.ILArgumentException">If the size of both arrays does not match any parameter rule.</exception>
789	public unsafe static ILRetArray<double> atan2(ILInArray<double> A, ILInArray<double> B) {
790	using (ILScope.Enter(A,B)) {
791	int outLen;
792	BinOpItMode mode;
793	double [] retArr;
794	double [] arrA = A.GetArrayForRead();
795	double[] arrB = B.GetArrayForRead();
796	ILSize outDims;
797	#region determine operation mode
798	if (A.IsScalar) {
799	outDims = B.Size;
800	if (B.IsScalar) {
801
802	return new ILRetArray<double> (new double [1]{ Math.Atan2 (A.GetValue(0) , B.GetValue(0))}, A.Size);
803	} else if (B.IsEmpty) {
804	return ILRetArray<double>.empty(outDims);
805	} else {
806	outLen = outDims.NumberOfElements;
807	if (!B.TryGetStorage4InplaceOp(out retArr)) {
808	retArr = ILMemoryPool.Pool.New< double > (outLen);
809	mode = BinOpItMode.SAN;
810	} else {
811	mode = BinOpItMode.SAI;
812	}
813	}
814	} else {
815	outDims = A.Size;
816	if (B.IsScalar) {
817	if (A.IsEmpty) {
818	return ILRetArray<double>.empty(A.Size);
819	}
820	outLen = A.S.NumberOfElements;
821	if (!A.TryGetStorage4InplaceOp(out retArr)) {
822	retArr = ILMemoryPool.Pool.New< double > (outLen);
823	mode = BinOpItMode.ASN;
824	} else {
825	mode = BinOpItMode.ASI;
826	}
827	} else {
828	// array + array
829	if (!A.Size.IsSameSize(B.Size)) {
830	return atan2Ex(A,B);
831	}
832	outLen = A.S.NumberOfElements;
833	if (A.TryGetStorage4InplaceOp(out retArr))
834	mode = BinOpItMode.AAIA;
835	else if (B.TryGetStorage4InplaceOp(out retArr))
836	mode = BinOpItMode.AAIB;
837	else {
838	retArr = ILMemoryPool.Pool.New< double > (outLen);
839	mode = BinOpItMode.AAN;
840	}
841	}
842	}
843	#endregion
844	ILDenseStorage<double> retStorage = new ILDenseStorage<double>(retArr, outDims);
845	int i = 0, workerCount = 1;
846	Action<object> worker = data => {
847	Tuple<int, int, IntPtr, IntPtr, IntPtr, BinOpItMode> range
848	= (Tuple<int, int, IntPtr, IntPtr, IntPtr, BinOpItMode>)data;
849	double* cp = (double*)range.Item5 + range.Item1;
850	double scalar;
851	int j = range.Item2;
852	#region loops
853	switch (mode) {
854	case BinOpItMode.AAIA:
855	double* bp = ((double*)range.Item4 + range.Item1);
856	while (j > 20) {
857	cp[0] = Math.Atan2 (cp[0] , bp[0]);
858	cp[1] = Math.Atan2 (cp[1] , bp[1]);
859	cp[2] = Math.Atan2 (cp[2] , bp[2]);
860	cp[3] = Math.Atan2 (cp[3] , bp[3]);
861	cp[4] = Math.Atan2 (cp[4] , bp[4]);
862	cp[5] = Math.Atan2 (cp[5] , bp[5]);
863	cp[6] = Math.Atan2 (cp[6] , bp[6]);
864	cp[7] = Math.Atan2 (cp[7] , bp[7]);
865	cp[8] = Math.Atan2 (cp[8] , bp[8]);
866	cp[9] = Math.Atan2 (cp[9] , bp[9]);
867	cp[10] = Math.Atan2 (cp[10] , bp[10]);
868	cp[11] = Math.Atan2 (cp[11] , bp[11]);
869	cp[12] = Math.Atan2 (cp[12] , bp[12]);
870	cp[13] = Math.Atan2 (cp[13] , bp[13]);
871	cp[14] = Math.Atan2 (cp[14] , bp[14]);
872	cp[15] = Math.Atan2 (cp[15] , bp[15]);
873	cp[16] = Math.Atan2 (cp[16] , bp[16]);
874	cp[17] = Math.Atan2 (cp[17] , bp[17]);
875	cp[18] = Math.Atan2 (cp[18] , bp[18]);
876	cp[19] = Math.Atan2 (cp[19] , bp[19]);
877	cp[20] = Math.Atan2 (cp[20] , bp[20]);
878	cp += 21; bp += 21; j -= 21;
879	}
880	while (j --> 0) {
881
882	cp = Math.Atan2 (cp , *bp);
883	cp++; bp++;
884	}
885	break;
886	case BinOpItMode.AAIB:
887	double* ap = ((double*)range.Item3 + range.Item1);
888	while (j > 20) {
889
890	cp[0] = Math.Atan2 (ap[0] , cp[0]);
891	cp[1] = Math.Atan2 (ap[1] , cp[1]);
892	cp[2] = Math.Atan2 (ap[2] , cp[2]);
893	cp[3] = Math.Atan2 (ap[3] , cp[3]);
894	cp[4] = Math.Atan2 (ap[4] , cp[4]);
895	cp[5] = Math.Atan2 (ap[5] , cp[5]);
896	cp[6] = Math.Atan2 (ap[6] , cp[6]);
897	cp[7] = Math.Atan2 (ap[7] , cp[7]);
898	cp[8] = Math.Atan2 (ap[8] , cp[8]);
899	cp[9] = Math.Atan2 (ap[9] , cp[9]);
900	cp[10] = Math.Atan2 (ap[10] , cp[10]);
901	cp[11] = Math.Atan2 (ap[11] , cp[11]);
902	cp[12] = Math.Atan2 (ap[12] , cp[12]);
903	cp[13] = Math.Atan2 (ap[13] , cp[13]);
904	cp[14] = Math.Atan2 (ap[14] , cp[14]);
905	cp[15] = Math.Atan2 (ap[15] , cp[15]);
906	cp[16] = Math.Atan2 (ap[16] , cp[16]);
907	cp[17] = Math.Atan2 (ap[17] , cp[17]);
908	cp[18] = Math.Atan2 (ap[18] , cp[18]);
909	cp[19] = Math.Atan2 (ap[19] , cp[19]);
910	cp[20] = Math.Atan2 (ap[20] , cp[20]);
911	ap += 21; cp += 21; j -= 21;
912	}
913	while (j --> 0) {
914
915	cp = Math.Atan2 (ap , *cp);
916	ap++; cp++;
917	}
918	break;
919	case BinOpItMode.AAN:
920	ap = ((double*)range.Item3 + range.Item1);
921	bp = ((double*)range.Item4 + range.Item1);
922	while (j > 20) {
923
924	cp[0] = Math.Atan2 (ap[0] , bp[0]);
925	cp[1] = Math.Atan2 (ap[1] , bp[1]);
926	cp[2] = Math.Atan2 (ap[2] , bp[2]);
927	cp[3] = Math.Atan2 (ap[3] , bp[3]);
928	cp[4] = Math.Atan2 (ap[4] , bp[4]);
929	cp[5] = Math.Atan2 (ap[5] , bp[5]);
930	cp[6] = Math.Atan2 (ap[6] , bp[6]);
931	cp[7] = Math.Atan2 (ap[7] , bp[7]);
932	cp[8] = Math.Atan2 (ap[8] , bp[8]);
933	cp[9] = Math.Atan2 (ap[9] , bp[9]);
934	cp[10] = Math.Atan2 (ap[10] , bp[10]);
935	cp[11] = Math.Atan2 (ap[11] , bp[11]);
936	cp[12] = Math.Atan2 (ap[12] , bp[12]);
937	cp[13] = Math.Atan2 (ap[13] , bp[13]);
938	cp[14] = Math.Atan2 (ap[14] , bp[14]);
939	cp[15] = Math.Atan2 (ap[15] , bp[15]);
940	cp[16] = Math.Atan2 (ap[16] , bp[16]);
941	cp[17] = Math.Atan2 (ap[17] , bp[17]);
942	cp[18] = Math.Atan2 (ap[18] , bp[18]);
943	cp[19] = Math.Atan2 (ap[19] , bp[19]);
944	cp[20] = Math.Atan2 (ap[20] , bp[20]);
945	ap+=21; bp+=21; cp+=21; j-=21;
946	}
947	while (j --> 0) {
948
949	cp = Math.Atan2 (ap , *bp);
950	ap++; bp++; cp++;
951	}
952	break;
953	case BinOpItMode.ASI:
954	scalar = ((double)range.Item4);
955	while (j > 20) {
956
957	cp[0] = Math.Atan2 (cp[0] , scalar);
958	cp[1] = Math.Atan2 (cp[1] , scalar);
959	cp[2] = Math.Atan2 (cp[2] , scalar);
960	cp[3] = Math.Atan2 (cp[3] , scalar);
961	cp[4] = Math.Atan2 (cp[4] , scalar);
962	cp[5] = Math.Atan2 (cp[5] , scalar);
963	cp[6] = Math.Atan2 (cp[6] , scalar);
964	cp[7] = Math.Atan2 (cp[7] , scalar);
965	cp[8] = Math.Atan2 (cp[8] , scalar);
966	cp[9] = Math.Atan2 (cp[9] , scalar);
967	cp[10] = Math.Atan2 (cp[10] , scalar);
968	cp[11] = Math.Atan2 (cp[11] , scalar);
969	cp[12] = Math.Atan2 (cp[12] , scalar);
970	cp[13] = Math.Atan2 (cp[13] , scalar);
971	cp[14] = Math.Atan2 (cp[14] , scalar);
972	cp[15] = Math.Atan2 (cp[15] , scalar);
973	cp[16] = Math.Atan2 (cp[16] , scalar);
974	cp[17] = Math.Atan2 (cp[17] , scalar);
975	cp[18] = Math.Atan2 (cp[18] , scalar);
976	cp[19] = Math.Atan2 (cp[19] , scalar);
977	cp[20] = Math.Atan2 (cp[20] , scalar);
978	cp += 21; j -= 21;
979	}
980	while (j --> 0) {
981
982	cp = Math.Atan2 (cp , scalar);
983	cp++;
984	}
985	break;
986	case BinOpItMode.ASN:
987	ap = ((double*)range.Item3 + range.Item1);
988	scalar = ((double)range.Item4);
989	while (j > 20) {
990
991	cp[0] = Math.Atan2 (ap[0] , scalar);
992	cp[1] = Math.Atan2 (ap[1] , scalar);
993	cp[2] = Math.Atan2 (ap[2] , scalar);
994	cp[3] = Math.Atan2 (ap[3] , scalar);
995	cp[4] = Math.Atan2 (ap[4] , scalar);
996	cp[5] = Math.Atan2 (ap[5] , scalar);
997	cp[6] = Math.Atan2 (ap[6] , scalar);
998	cp[7] = Math.Atan2 (ap[7] , scalar);
999	cp[8] = Math.Atan2 (ap[8] , scalar);
1000	cp[9] = Math.Atan2 (ap[9] , scalar);
1001	cp[10] = Math.Atan2 (ap[10] , scalar);
1002	cp[11] = Math.Atan2 (ap[11] , scalar);
1003	cp[12] = Math.Atan2 (ap[12] , scalar);
1004	cp[13] = Math.Atan2 (ap[13] , scalar);
1005	cp[14] = Math.Atan2 (ap[14] , scalar);
1006	cp[15] = Math.Atan2 (ap[15] , scalar);
1007	cp[16] = Math.Atan2 (ap[16] , scalar);
1008	cp[17] = Math.Atan2 (ap[17] , scalar);
1009	cp[18] = Math.Atan2 (ap[18] , scalar);
1010	cp[19] = Math.Atan2 (ap[19] , scalar);
1011	cp[20] = Math.Atan2 (ap[20] , scalar);
1012	ap+=21; cp+=21; j -= 21;
1013	}
1014	while (j --> 0) {
1015
1016	cp = Math.Atan2 (ap , scalar);
1017	ap++; cp++;
1018	}
1019	break;
1020	case BinOpItMode.SAI:
1021	scalar = ((double)range.Item3);
1022	while (j > 20) {
1023
1024	cp[0] = Math.Atan2 (scalar , cp[0]);
1025	cp[1] = Math.Atan2 (scalar , cp[1]);
1026	cp[2] = Math.Atan2 (scalar , cp[2]);
1027	cp[3] = Math.Atan2 (scalar , cp[3]);
1028	cp[4] = Math.Atan2 (scalar , cp[4]);
1029	cp[5] = Math.Atan2 (scalar , cp[5]);
1030	cp[6] = Math.Atan2 (scalar , cp[6]);
1031	cp[7] = Math.Atan2 (scalar , cp[7]);
1032	cp[8] = Math.Atan2 (scalar , cp[8]);
1033	cp[9] = Math.Atan2 (scalar , cp[9]);
1034	cp[10] = Math.Atan2 (scalar , cp[10]);
1035	cp[11] = Math.Atan2 (scalar , cp[11]);
1036	cp[12] = Math.Atan2 (scalar , cp[12]);
1037	cp[13] = Math.Atan2 (scalar , cp[13]);
1038	cp[14] = Math.Atan2 (scalar , cp[14]);
1039	cp[15] = Math.Atan2 (scalar , cp[15]);
1040	cp[16] = Math.Atan2 (scalar , cp[16]);
1041	cp[17] = Math.Atan2 (scalar , cp[17]);
1042	cp[18] = Math.Atan2 (scalar , cp[18]);
1043	cp[19] = Math.Atan2 (scalar , cp[19]);
1044	cp[20] = Math.Atan2 (scalar , cp[20]);
1045	cp += 21; j -= 21;
1046	}
1047	while (j --> 0) {
1048
1049	cp = Math.Atan2 (scalar , cp);
1050	cp++;
1051	}
1052	break;
1053	case BinOpItMode.SAN:
1054	scalar = ((double)range.Item3);
1055	bp = ((double*)range.Item4 + range.Item1);
1056	while (j > 20) {
1057
1058	cp[0] = Math.Atan2 (scalar , bp[0]);
1059	cp[1] = Math.Atan2 (scalar , bp[1]);
1060	cp[2] = Math.Atan2 (scalar , bp[2]);
1061	cp[3] = Math.Atan2 (scalar , bp[3]);
1062	cp[4] = Math.Atan2 (scalar , bp[4]);
1063	cp[5] = Math.Atan2 (scalar , bp[5]);
1064	cp[6] = Math.Atan2 (scalar , bp[6]);
1065	cp[7] = Math.Atan2 (scalar , bp[7]);
1066	cp[8] = Math.Atan2 (scalar , bp[8]);
1067	cp[9] = Math.Atan2 (scalar , bp[9]);
1068	cp[10] = Math.Atan2 (scalar , bp[10]);
1069	cp[11] = Math.Atan2 (scalar , bp[11]);
1070	cp[12] = Math.Atan2 (scalar , bp[12]);
1071	cp[13] = Math.Atan2 (scalar , bp[13]);
1072	cp[14] = Math.Atan2 (scalar , bp[14]);
1073	cp[15] = Math.Atan2 (scalar , bp[15]);
1074	cp[16] = Math.Atan2 (scalar , bp[16]);
1075	cp[17] = Math.Atan2 (scalar , bp[17]);
1076	cp[18] = Math.Atan2 (scalar , bp[18]);
1077	cp[19] = Math.Atan2 (scalar , bp[19]);
1078	cp[20] = Math.Atan2 (scalar , bp[20]);
1079	bp+=21; cp+=21; j -= 21;
1080	}
1081	while (j --> 0) {
1082
1083	cp = Math.Atan2 (scalar , bp);
1084	bp++; cp++;
1085	}
1086	break;
1087	default:
1088	break;
1089	}
1090	#endregion
1091	System.Threading.Interlocked.Decrement(ref workerCount);
1092	//retStorage.PendingEvents.Signal();
1093	};
1094
1095	#region do the work
1096	int workItemCount = Settings.s_maxNumberThreads, workItemLength;
1097	if (Settings.s_maxNumberThreads > 1 && outLen / 2 > Settings.s_minParallelElement1Count) {
1098	if (outLen / workItemCount > Settings.s_minParallelElement1Count) {
1099	workItemLength = outLen / workItemCount;
1100	//workItemLength = (int)((double)outLen / workItemCount * 1.05);
1101	} else {
1102	workItemLength = outLen / 2;
1103	workItemCount = 2;
1104	}
1105	} else {
1106	workItemLength = outLen;
1107	workItemCount = 1;
1108	}
1109
1110	fixed ( double* arrAP = arrA)
1111	fixed ( double* arrBP = arrB)
1112	fixed ( double* retArrP = retArr) {
1113
1114	for (; i < workItemCount - 1; i++) {
1115	Tuple<int, int, IntPtr, IntPtr, IntPtr, BinOpItMode> range
1116	= new Tuple<int, int, IntPtr, IntPtr, IntPtr, BinOpItMode>
1117	(i * workItemLength, workItemLength, (IntPtr)arrAP, (IntPtr)arrBP, (IntPtr)retArrP, mode);
1118	System.Threading.Interlocked.Increment(ref workerCount);
1119	ILThreadPool.QueueUserWorkItem(i, worker, range);
1120	}
1121	// the last (or may the only) chunk is done right here
1122	//System.Threading.Interlocked.Increment(ref retStorage.PendingTasks);
1123	worker(new Tuple<int, int, IntPtr, IntPtr, IntPtr, BinOpItMode>
1124	(i * workItemLength, outLen - i * workItemLength, (IntPtr)arrAP, (IntPtr)arrBP, (IntPtr)retArrP, mode));
1125
1126	ILThreadPool.Wait4Workers(ref workerCount);
1127	}
1128
1129	#endregion
1130	return new ILRetArray< double>(retStorage);
1131	}
1132	}
1133
1134	private static unsafe ILRetArray<double> atan2Ex(ILInArray<double> A, ILInArray<double> B) {
1135	//using (ILScope.Enter(A, B)) { we cannot start a new scope here, since this would prevent A and B to be used implace if applicable
1136
1137	#region parameter checking
1138	if (isnull(A) \|\| isnull(B))
1139	return empty<double>(ILSize.Empty00);
1140	if (A.IsEmpty) {
1141	return empty<double>(B.S);
1142	} else if (B.IsEmpty) {
1143	return empty<double>(A.S);
1144	}
1145	//if (A.IsScalar \|\| B.IsScalar \|\| A.D.IsSameSize(B.D))
1146	// return add(A,B);
1147	int dim = -1;
1148	for (int l = 0; l < Math.Max(A.S.NumberOfDimensions, B.S.NumberOfDimensions); l++) {
1149	if (A.S[l] != B.S[l]) {
1150	if (dim >= 0 \|\| (A.S[l] != 1 && B.S[l] != 1)) {
1151	throw new ILArgumentException("A and B must have the same size except for one simgleton dimension in A or B");
1152	}
1153	dim = l;
1154	}
1155	}
1156	if (dim > 1)
1157	throw new ILArgumentException("singleton dimension expansion currently is only supported for colum- and row vectors");
1158	dim = -(dim - 1); // 0 -> 1, 1 -> 0
1159	#endregion
1160
1161	#region parameter preparation
1162
1163	double[] retArr;
1164
1165	double[] arrA = A.GetArrayForRead();
1166
1167	double[] arrB = B.GetArrayForRead();
1168	ILSize outDims;
1169	BinOptItExMode mode;
1170	int workItemMultiplierLenA;
1171	int workItemMultiplierLenB;
1172	if (A.IsVector) {
1173	outDims = B.S;
1174	if (!B.TryGetStorage4InplaceOp(out retArr)) {
1175	retArr = ILMemoryPool.Pool.New<double>(outDims.NumberOfElements);
1176	mode = BinOptItExMode.VAN;
1177	} else {
1178	mode = BinOptItExMode.VAI;
1179	}
1180	workItemMultiplierLenB = outDims[0];
1181	workItemMultiplierLenA = dim; // 0 for column, 1 for row vector
1182	} else if (B.IsVector) {
1183	outDims = A.S;
1184	if (!A.TryGetStorage4InplaceOp(out retArr)) {
1185	retArr = ILMemoryPool.Pool.New<double>(outDims.NumberOfElements);
1186	mode = BinOptItExMode.AVN;
1187	} else {
1188	mode = BinOptItExMode.AVI;
1189	}
1190	workItemMultiplierLenB = dim; // 0 for column, 1 for row vector
1191	workItemMultiplierLenA = outDims[0];
1192	} else {
1193	throw new ILArgumentException("A and B must have the same size except for one singleton dimension in either A or B");
1194	}
1195	int itLen = outDims[0]; // (dim == 0) ? outDims.SequentialIndexDistance(1) : outDims.SequentialIndexDistance(0);
1196	#endregion
1197
1198	#region worker loops definition
1199	ILDenseStorage<double> retStorage = new ILDenseStorage<double>(retArr, outDims);
1200	int workerCount = 1;
1201	Action<object> worker = data => {
1202	// expects: iStart, iLen, ap, bp, cp
1203	Tuple<int, IntPtr, IntPtr, IntPtr> range =
1204	(Tuple<int, IntPtr, IntPtr, IntPtr>)data;
1205
1206	double* ap;
1207
1208	double* bp;
1209
1210	double* cp;
1211	switch (mode) {
1212	case BinOptItExMode.VAN:
1213	if (dim == 0) {
1214	bp = (double*)range.Item3;
1215	cp = (double*)range.Item4;
1216	for (int s = 0; s < range.Item1; s++) {
1217	ap = (double*)range.Item2;
1218	int l = itLen;
1219	while (l > 10) {
1220	cp[0] = Math.Atan2 (ap[0] , bp[0]);
1221	cp[1] = Math.Atan2 (ap[1] , bp[1]);
1222	cp[2] = Math.Atan2 (ap[2] , bp[2]);
1223	cp[3] = Math.Atan2 (ap[3] , bp[3]);
1224	cp[4] = Math.Atan2 (ap[4] , bp[4]);
1225	cp[5] = Math.Atan2 (ap[5] , bp[5]);
1226	cp[6] = Math.Atan2 (ap[6] , bp[6]);
1227	cp[7] = Math.Atan2 (ap[7] , bp[7]);
1228	cp[8] = Math.Atan2 (ap[8] , bp[8]);
1229	cp[9] = Math.Atan2 (ap[9] , bp[9]);
1230	cp[10] = Math.Atan2 (ap[10] , bp[10]);
1231	ap += 11;
1232	bp += 11;
1233	cp += 11;
1234	l -= 11;
1235	}
1236	while (l-- > 0) {
1237	cp++ = Math.Atan2 (ap++ , *bp++);
1238	}
1239	}
1240	} else {
1241	// dim == 1
1242	ap = (double*)range.Item2;
1243	bp = (double*)range.Item3;
1244	cp = (double*)range.Item4;
1245	for (int s = 0; s < range.Item1; s++) {
1246	double val = *ap++;
1247	int l = itLen;
1248	while (l > 10) {
1249	cp[0] = Math.Atan2 (val , bp[0]);
1250	cp[1] = Math.Atan2 (val , bp[1]);
1251	cp[2] = Math.Atan2 (val , bp[2]);
1252	cp[3] = Math.Atan2 (val , bp[3]);
1253	cp[4] = Math.Atan2 (val , bp[4]);
1254	cp[5] = Math.Atan2 (val , bp[5]);
1255	cp[6] = Math.Atan2 (val , bp[6]);
1256	cp[7] = Math.Atan2 (val , bp[7]);
1257	cp[8] = Math.Atan2 (val , bp[8]);
1258	cp[9] = Math.Atan2 (val , bp[9]);
1259	cp[10] = Math.Atan2 (val , bp[10]);
1260	bp += 11;
1261	cp += 11;
1262	l -= 11;
1263	}
1264	while (l-- > 0) {
1265	cp++ = Math.Atan2 (val , bp++);
1266	}
1267	}
1268	}
1269	break;
1270	case BinOptItExMode.VAI:
1271	if (dim == 0) {
1272	cp = (double*)range.Item4;
1273	for (int s = 0; s < range.Item1; s++) {
1274	ap = (double*)range.Item2;
1275	int l = itLen;
1276	while (l > 10) {
1277	cp[0] = Math.Atan2 (ap[0] , cp[0]);
1278	cp[1] = Math.Atan2 (ap[1] , cp[1]);
1279	cp[2] = Math.Atan2 (ap[2] , cp[2]);
1280	cp[3] = Math.Atan2 (ap[3] , cp[3]);
1281	cp[4] = Math.Atan2 (ap[4] , cp[4]);
1282	cp[5] = Math.Atan2 (ap[5] , cp[5]);
1283	cp[6] = Math.Atan2 (ap[6] , cp[6]);
1284	cp[7] = Math.Atan2 (ap[7] , cp[7]);
1285	cp[8] = Math.Atan2 (ap[8] , cp[8]);
1286	cp[9] = Math.Atan2 (ap[9] , cp[9]);
1287	cp[10] = Math.Atan2 (ap[10] , cp[10]);
1288	ap += 11;
1289	cp += 11;
1290	l -= 11;
1291	}
1292	while (l-- > 0) {
1293	cp = Math.Atan2 (ap++ , *cp);
1294	cp++;
1295	}
1296	}
1297	} else {
1298	// dim == 1
1299	cp = (double*)range.Item4;
1300	ap = (double*)range.Item2;
1301	for (int s = 0; s < range.Item1; s++) {
1302
1303	double val = *ap++;
1304	int l = itLen;
1305	while (l > 10) {
1306	cp[0] = Math.Atan2 (val , cp[0]);
1307	cp[1] = Math.Atan2 (val , cp[1]);
1308	cp[2] = Math.Atan2 (val , cp[2]);
1309	cp[3] = Math.Atan2 (val , cp[3]);
1310	cp[4] = Math.Atan2 (val , cp[4]);
1311	cp[5] = Math.Atan2 (val , cp[5]);
1312	cp[6] = Math.Atan2 (val , cp[6]);
1313	cp[7] = Math.Atan2 (val , cp[7]);
1314	cp[8] = Math.Atan2 (val , cp[8]);
1315	cp[9] = Math.Atan2 (val , cp[9]);
1316	cp[10] = Math.Atan2 (val , cp[10]);
1317	cp += 11;
1318	l -= 11;
1319	}
1320	while (l-- > 0) {
1321	cp = Math.Atan2 (val , cp);
1322	cp++;
1323	}
1324	}
1325	}
1326	break;
1327	case BinOptItExMode.AVN:
1328	if (dim == 0) {
1329	ap = (double*)range.Item2;
1330	cp = (double*)range.Item4;
1331	for (int s = 0; s < range.Item1; s++) {
1332	bp = (double*)range.Item3;
1333	int l = itLen;
1334	while (l > 10) {
1335	cp[0] = Math.Atan2 (ap[0] , bp[0]);
1336	cp[1] = Math.Atan2 (ap[1] , bp[1]);
1337	cp[2] = Math.Atan2 (ap[2] , bp[2]);
1338	cp[3] = Math.Atan2 (ap[3] , bp[3]);
1339	cp[4] = Math.Atan2 (ap[4] , bp[4]);
1340	cp[5] = Math.Atan2 (ap[5] , bp[5]);
1341	cp[6] = Math.Atan2 (ap[6] , bp[6]);
1342	cp[7] = Math.Atan2 (ap[7] , bp[7]);
1343	cp[8] = Math.Atan2 (ap[8] , bp[8]);
1344	cp[9] = Math.Atan2 (ap[9] , bp[9]);
1345	cp[10] = Math.Atan2 (ap[10] , bp[10]);
1346	ap += 11;
1347	bp += 11;
1348	cp += 11;
1349	l -= 11;
1350	}
1351	while (l-- > 0) {
1352	cp = Math.Atan2 (ap , *bp);
1353	ap++;
1354	bp++;
1355	cp++;
1356	}
1357	}
1358	} else {
1359	// dim = 1
1360	ap = (double*)range.Item2;
1361	bp = (double*)range.Item3;
1362	cp = (double*)range.Item4;
1363	for (int s = 0; s < range.Item1; s++) {
1364	double val = *bp++;
1365	int l = itLen;
1366	while (l > 10) {
1367	cp[0] = Math.Atan2 (ap[0] , val);
1368	cp[1] = Math.Atan2 (ap[1] , val);
1369	cp[2] = Math.Atan2 (ap[2] , val);
1370	cp[3] = Math.Atan2 (ap[3] , val);
1371	cp[4] = Math.Atan2 (ap[4] , val);
1372	cp[5] = Math.Atan2 (ap[5] , val);
1373	cp[6] = Math.Atan2 (ap[6] , val);
1374	cp[7] = Math.Atan2 (ap[7] , val);
1375	cp[8] = Math.Atan2 (ap[8] , val);
1376	cp[9] = Math.Atan2 (ap[9] , val);
1377	cp[10] = Math.Atan2 (ap[10] , val);
1378	ap += 11;
1379	cp += 11;
1380	l -= 11;
1381	}
1382	while (l-- > 0) {
1383	cp = Math.Atan2 (ap , val);
1384	ap++;
1385	cp++;
1386	}
1387	}
1388	}
1389	break;
1390	case BinOptItExMode.AVI:
1391	if (dim == 0) {
1392	cp = (double*)range.Item4;
1393	for (int s = 0; s < range.Item1; s++) {
1394	bp = (double*)range.Item3;
1395	int l = itLen;
1396	while (l > 10) {
1397	cp[0] = Math.Atan2 (cp[0] , bp[0]);
1398	cp[1] = Math.Atan2 (cp[1] , bp[1]);
1399	cp[2] = Math.Atan2 (cp[2] , bp[2]);
1400	cp[3] = Math.Atan2 (cp[3] , bp[3]);
1401	cp[4] = Math.Atan2 (cp[4] , bp[4]);
1402	cp[5] = Math.Atan2 (cp[5] , bp[5]);
1403	cp[6] = Math.Atan2 (cp[6] , bp[6]);
1404	cp[7] = Math.Atan2 (cp[7] , bp[7]);
1405	cp[8] = Math.Atan2 (cp[8] , bp[8]);
1406	cp[9] = Math.Atan2 (cp[9] , bp[9]);
1407	cp[10] = Math.Atan2 (cp[10] , bp[10]);
1408	bp += 11;
1409	cp += 11;
1410	l -= 11;
1411	}
1412	while (l-- > 0) {
1413	cp = Math.Atan2 (cp , *bp);
1414	bp++;
1415	cp++;
1416	}
1417	}
1418	} else {
1419	// dim = 1
1420	bp = (double*)range.Item3;
1421	cp = (double*)range.Item4;
1422	for (int s = 0; s < range.Item1; s++) {
1423
1424	double val = *bp++;
1425	int l = itLen;
1426	while (l > 10) {
1427	cp[0] = Math.Atan2 (cp[0] , val);
1428	cp[1] = Math.Atan2 (cp[1] , val);
1429	cp[2] = Math.Atan2 (cp[2] , val);
1430	cp[3] = Math.Atan2 (cp[3] , val);
1431	cp[4] = Math.Atan2 (cp[4] , val);
1432	cp[5] = Math.Atan2 (cp[5] , val);
1433	cp[6] = Math.Atan2 (cp[6] , val);
1434	cp[7] = Math.Atan2 (cp[7] , val);
1435	cp[8] = Math.Atan2 (cp[8] , val);
1436	cp[9] = Math.Atan2 (cp[9] , val);
1437	cp[10] = Math.Atan2 (cp[10] , val);
1438	cp += 11;
1439	l -= 11;
1440	}
1441	while (l --> 0) {
1442	cp = Math.Atan2 (cp , val);
1443	cp++;
1444	}
1445	}
1446	}
1447	break;
1448	}
1449	System.Threading.Interlocked.Decrement(ref workerCount);
1450	};
1451	#endregion
1452
1453	#region work distribution
1454	int i = 0, workItemCount = Settings.s_maxNumberThreads, workItemLength;
1455	if (Settings.s_maxNumberThreads > 1 && outDims.NumberOfElements >= Settings.s_minParallelElement1Count
1456	&& outDims[1] > 1) {
1457	if (outDims[1] > workItemCount) {
1458	workItemLength = outDims[1] / workItemCount;
1459	} else {
1460	workItemLength = outDims[1] / 2;
1461	workItemCount = 2;
1462	}
1463	} else {
1464	workItemLength = outDims[1];
1465	workItemCount = 1;
1466	}
1467
1468	fixed ( double* arrAP = arrA)
1469	fixed ( double* arrBP = arrB)
1470	fixed ( double* retArrP = retArr) {
1471
1472	for (; i < workItemCount - 1; i++) {
1473	Tuple<int, IntPtr, IntPtr, IntPtr> range = new Tuple<int, IntPtr, IntPtr, IntPtr>
1474	(workItemLength
1475	, (IntPtr)(arrAP + i * workItemMultiplierLenA * workItemLength)
1476	, (IntPtr)(arrBP + i * workItemMultiplierLenB * workItemLength)
1477	, (IntPtr)(retArrP + i * outDims[0] * workItemLength));
1478	System.Threading.Interlocked.Increment(ref workerCount);
1479	ILThreadPool.QueueUserWorkItem(i, worker, range);
1480	}
1481	// the last (or may the only) chunk is done right here
1482	//System.Threading.Interlocked.Increment(ref retStorage.PendingTasks);
1483	worker(new Tuple<int, IntPtr, IntPtr, IntPtr>
1484	(outDims[1] - i * workItemLength
1485	, (IntPtr)(arrAP + i * workItemMultiplierLenA * workItemLength)
1486	, (IntPtr)(arrBP + i * workItemMultiplierLenB * workItemLength)
1487	, (IntPtr)(retArrP + i * outDims[0] * workItemLength)));
1488
1489	ILThreadPool.Wait4Workers(ref workerCount);
1490	}
1491	#endregion
1492
1493	return new ILRetArray<double>(retStorage);
1494	//} // no scopes here! it disables implace operations
1495	}
1496
1497
1498	#endregion HYCALPER AUTO GENERATED CODE
1499
1500	}
1501	}

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Update cookies preferences