View Javadoc
1   /*
2    * Licensed to the Apache Software Foundation (ASF) under one or more
3    * contributor license agreements.  See the NOTICE file distributed with
4    * this work for additional information regarding copyright ownership.
5    * The ASF licenses this file to You under the Apache License, Version 2.0
6    * (the "License"); you may not use this file except in compliance with
7    * the License.  You may obtain a copy of the License at
8    *
9    *      https://www.apache.org/licenses/LICENSE-2.0
10   *
11   * Unless required by applicable law or agreed to in writing, software
12   * distributed under the License is distributed on an "AS IS" BASIS,
13   * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14   * See the License for the specific language governing permissions and
15   * limitations under the License.
16   */
17  
18  /*
19   * This is not the original file distributed by the Apache Software Foundation
20   * It has been modified by the Hipparchus project
21   */
22  package org.hipparchus.complex;
23  
24  import java.io.Serializable;
25  
26  import org.hipparchus.exception.LocalizedCoreFormats;
27  import org.hipparchus.exception.MathIllegalArgumentException;
28  import org.hipparchus.exception.MathIllegalStateException;
29  import org.hipparchus.util.FastMath;
30  import org.hipparchus.util.SinCos;
31  
32  /**
33   * A helper class for the computation and caching of the {@code n}-th roots
34   * of unity.
35   */
36  public class RootsOfUnity implements Serializable {
37  
38      /** Serializable version id. */
39      private static final long serialVersionUID = 20120201L;
40  
41      /** Number of roots of unity. */
42      private int omegaCount;
43  
44      /** Real part of the roots. */
45      private double[] omegaReal;
46  
47      /**
48       * Imaginary part of the {@code n}-th roots of unity, for positive values
49       * of {@code n}. In this array, the roots are stored in counter-clockwise
50       * order.
51       */
52      private double[] omegaImaginaryCounterClockwise;
53  
54      /**
55       * Imaginary part of the {@code n}-th roots of unity, for negative values
56       * of {@code n}. In this array, the roots are stored in clockwise order.
57       */
58      private double[] omegaImaginaryClockwise;
59  
60      /**
61       * {@code true} if {@link #computeRoots(int)} was called with a positive
62       * value of its argument {@code n}. In this case, counter-clockwise ordering
63       * of the roots of unity should be used.
64       */
65      private boolean isCounterClockWise;
66  
67      /**
68       * Build an engine for computing the {@code n}-th roots of unity.
69       */
70      public RootsOfUnity() {
71          omegaCount = 0;
72          omegaReal = null;
73          omegaImaginaryCounterClockwise = null;
74          omegaImaginaryClockwise = null;
75          isCounterClockWise = true;
76      }
77  
78      /**
79       * Returns {@code true} if {@link #computeRoots(int)} was called with a
80       * positive value of its argument {@code n}. If {@code true}, then
81       * counter-clockwise ordering of the roots of unity should be used.
82       *
83       * @return {@code true} if the roots of unity are stored in counter-clockwise order
84       * @throws MathIllegalStateException if no roots of unity have been computed yet
85       */
86      public boolean isCounterClockWise()
87              throws MathIllegalStateException {
88          synchronized (this) {
89              if (omegaCount == 0) {
90                  throw new MathIllegalStateException(LocalizedCoreFormats.ROOTS_OF_UNITY_NOT_COMPUTED_YET);
91              }
92              return isCounterClockWise;
93          }
94      }
95  
96      /**
97       * Computes the {@code n}-th roots of unity.
98       * <p>
99       * The roots are stored in {@code omega[]}, such that {@code omega[k] = w ^ k},
100      * where {@code k = 0, ..., n - 1}, {@code w = exp(2 * pi * i / n)} and
101      * {@code i = sqrt(-1)}.
102      * <p>
103      * Note that {@code n} can be positive of negative
104      * <ul>
105      * <li>{@code abs(n)} is always the number of roots of unity.</li>
106      * <li>If {@code n > 0}, then the roots are stored in counter-clockwise order.</li>
107      * <li>If {@code n < 0}, then the roots are stored in clockwise order.</li>
108      * </ul>
109      *
110      * @param n the (signed) number of roots of unity to be computed
111      * @throws MathIllegalArgumentException if {@code n = 0}
112      */
113     public void computeRoots(int n) throws MathIllegalArgumentException {
114 
115         if (n == 0) {
116             throw new MathIllegalArgumentException(LocalizedCoreFormats.CANNOT_COMPUTE_0TH_ROOT_OF_UNITY);
117         }
118 
119         synchronized (this) {
120             isCounterClockWise = n > 0;
121 
122             // avoid repetitive calculations
123             final int absN = FastMath.abs(n);
124 
125             if (absN == omegaCount) {
126                 return;
127             }
128 
129             // calculate everything from scratch
130             final double t  = 2.0 * FastMath.PI / absN;
131             final SinCos sc = FastMath.sinCos(t);
132             omegaReal = new double[absN];
133             omegaImaginaryCounterClockwise = new double[absN];
134             omegaImaginaryClockwise = new double[absN];
135             omegaReal[0] = 1.0;
136             omegaImaginaryCounterClockwise[0] = 0.0;
137             omegaImaginaryClockwise[0] = 0.0;
138             for (int i = 1; i < absN; i++) {
139                 omegaReal[i] = omegaReal[i - 1] * sc.cos() -
140                                 omegaImaginaryCounterClockwise[i - 1] * sc.sin();
141                 omegaImaginaryCounterClockwise[i] = omegaReal[i - 1] * sc.sin() +
142                                 omegaImaginaryCounterClockwise[i - 1] * sc.cos();
143                 omegaImaginaryClockwise[i] = -omegaImaginaryCounterClockwise[i];
144             }
145             omegaCount = absN;
146         }
147     }
148 
149     /**
150      * Get the real part of the {@code k}-th {@code n}-th root of unity.
151      *
152      * @param k index of the {@code n}-th root of unity
153      * @return real part of the {@code k}-th {@code n}-th root of unity
154      * @throws MathIllegalStateException if no roots of unity have been computed yet
155      * @throws MathIllegalArgumentException if {@code k} is out of range
156      */
157     public double getReal(int k)
158             throws MathIllegalArgumentException, MathIllegalStateException {
159 
160         synchronized (this) {
161             if (omegaCount == 0) {
162                 throw new MathIllegalStateException(LocalizedCoreFormats.ROOTS_OF_UNITY_NOT_COMPUTED_YET);
163             }
164             if ((k < 0) || (k >= omegaCount)) {
165                 throw new MathIllegalArgumentException(LocalizedCoreFormats.OUT_OF_RANGE_ROOT_OF_UNITY_INDEX,
166                                                        k, 0, omegaCount - 1);
167             }
168 
169             return omegaReal[k];
170         }
171     }
172 
173     /**
174      * Get the imaginary part of the {@code k}-th {@code n}-th root of unity.
175      *
176      * @param k index of the {@code n}-th root of unity
177      * @return imaginary part of the {@code k}-th {@code n}-th root of unity
178      * @throws MathIllegalStateException if no roots of unity have been computed yet
179      * @throws MathIllegalArgumentException if {@code k} is out of range
180      */
181     public double getImaginary(int k)
182             throws MathIllegalArgumentException, MathIllegalStateException {
183 
184         synchronized (this) {
185             if (omegaCount == 0) {
186                 throw new MathIllegalStateException(LocalizedCoreFormats.ROOTS_OF_UNITY_NOT_COMPUTED_YET);
187             }
188             if ((k < 0) || (k >= omegaCount)) {
189                 throw new MathIllegalArgumentException(LocalizedCoreFormats.OUT_OF_RANGE_ROOT_OF_UNITY_INDEX,
190                                                        k, 0, omegaCount - 1);
191             }
192 
193             return isCounterClockWise ?
194                    omegaImaginaryCounterClockwise[k] :
195                    omegaImaginaryClockwise[k];
196         }
197     }
198 
199     /**
200      * Returns the number of roots of unity currently stored.
201      * <p>
202      * If {@link #computeRoots(int)} was called with {@code n}, then this method
203      * returns {@code abs(n)}. If no roots of unity have been computed yet, this
204      * method returns 0.
205      *
206      * @return the number of roots of unity currently stored
207      */
208     public int getNumberOfRoots() {
209         synchronized (this) {
210             return omegaCount;
211         }
212     }
213 }