1 /* matrix/gsl_matrix_complex_long_double.h
2  *
3  * Copyright (C) 1996, 1997, 1998, 1999, 2000, 2007 Gerard Jungman, Brian Gough
4  *
5  * This program is free software; you can redistribute it and/or modify
6  * it under the terms of the GNU General Public License as published by
7  * the Free Software Foundation; either version 3 of the License, or (at
8  * your option) any later version.
9  *
10  * This program is distributed in the hope that it will be useful, but
11  * WITHOUT ANY WARRANTY; without even the implied warranty of
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
13  * General Public License for more details.
14  *
15  * You should have received a copy of the GNU General Public License
16  * along with this program; if not, write to the Free Software
17  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
18  */
19 
20 #ifndef __GSL_MATRIX_COMPLEX_LONG_DOUBLE_H__
21 #define __GSL_MATRIX_COMPLEX_LONG_DOUBLE_H__
22 
23 #include <stdlib.h>
24 #include <gsl/gsl_types.h>
25 #include <gsl/gsl_errno.h>
26 #include <gsl/gsl_complex.h>
27 #include <gsl/gsl_check_range.h>
28 #include <gsl/gsl_vector_complex_long_double.h>
29 
30 #undef __BEGIN_DECLS
31 #undef __END_DECLS
32 #ifdef __cplusplus
33 # define __BEGIN_DECLS extern "C" {
34 # define __END_DECLS }
35 #else
36 # define __BEGIN_DECLS /* empty */
37 # define __END_DECLS /* empty */
38 #endif
39 
40 __BEGIN_DECLS
41 
42 typedef struct
43 {
44   size_t size1;
45   size_t size2;
46   size_t tda;
47   long double * data;
48   gsl_block_complex_long_double * block;
49   int owner;
50 } gsl_matrix_complex_long_double ;
51 
52 typedef struct
53 {
54   gsl_matrix_complex_long_double matrix;
55 } _gsl_matrix_complex_long_double_view;
56 
57 typedef _gsl_matrix_complex_long_double_view gsl_matrix_complex_long_double_view;
58 
59 typedef struct
60 {
61   gsl_matrix_complex_long_double matrix;
62 } _gsl_matrix_complex_long_double_const_view;
63 
64 typedef const _gsl_matrix_complex_long_double_const_view gsl_matrix_complex_long_double_const_view;
65 
66 
67 /* Allocation */
68 
69 gsl_matrix_complex_long_double *
70 gsl_matrix_complex_long_double_alloc (const size_t n1, const size_t n2);
71 
72 gsl_matrix_complex_long_double *
73 gsl_matrix_complex_long_double_calloc (const size_t n1, const size_t n2);
74 
75 gsl_matrix_complex_long_double *
76 gsl_matrix_complex_long_double_alloc_from_block (gsl_block_complex_long_double * b,
77                                            const size_t offset,
78                                            const size_t n1, const size_t n2, const size_t d2);
79 
80 gsl_matrix_complex_long_double *
81 gsl_matrix_complex_long_double_alloc_from_matrix (gsl_matrix_complex_long_double * b,
82                                             const size_t k1, const size_t k2,
83                                             const size_t n1, const size_t n2);
84 
85 gsl_vector_complex_long_double *
86 gsl_vector_complex_long_double_alloc_row_from_matrix (gsl_matrix_complex_long_double * m,
87                                                 const size_t i);
88 
89 gsl_vector_complex_long_double *
90 gsl_vector_complex_long_double_alloc_col_from_matrix (gsl_matrix_complex_long_double * m,
91                                                 const size_t j);
92 
93 void gsl_matrix_complex_long_double_free (gsl_matrix_complex_long_double * m);
94 
95 /* Views */
96 
97 _gsl_matrix_complex_long_double_view
98 gsl_matrix_complex_long_double_submatrix (gsl_matrix_complex_long_double * m,
99                             const size_t i, const size_t j,
100                             const size_t n1, const size_t n2);
101 
102 _gsl_vector_complex_long_double_view
103 gsl_matrix_complex_long_double_row (gsl_matrix_complex_long_double * m, const size_t i);
104 
105 _gsl_vector_complex_long_double_view
106 gsl_matrix_complex_long_double_column (gsl_matrix_complex_long_double * m, const size_t j);
107 
108 _gsl_vector_complex_long_double_view
109 gsl_matrix_complex_long_double_diagonal (gsl_matrix_complex_long_double * m);
110 
111 _gsl_vector_complex_long_double_view
112 gsl_matrix_complex_long_double_subdiagonal (gsl_matrix_complex_long_double * m, const size_t k);
113 
114 _gsl_vector_complex_long_double_view
115 gsl_matrix_complex_long_double_superdiagonal (gsl_matrix_complex_long_double * m, const size_t k);
116 
117 _gsl_vector_complex_long_double_view
118 gsl_matrix_complex_long_double_subrow (gsl_matrix_complex_long_double * m,
119                                  const size_t i, const size_t offset,
120                                  const size_t n);
121 
122 _gsl_vector_complex_long_double_view
123 gsl_matrix_complex_long_double_subcolumn (gsl_matrix_complex_long_double * m,
124                                     const size_t j, const size_t offset,
125                                     const size_t n);
126 
127 _gsl_matrix_complex_long_double_view
128 gsl_matrix_complex_long_double_view_array (long double * base,
129                              const size_t n1,
130                              const size_t n2);
131 
132 _gsl_matrix_complex_long_double_view
133 gsl_matrix_complex_long_double_view_array_with_tda (long double * base,
134                                       const size_t n1,
135                                       const size_t n2,
136                                       const size_t tda);
137 
138 _gsl_matrix_complex_long_double_view
139 gsl_matrix_complex_long_double_view_vector (gsl_vector_complex_long_double * v,
140                               const size_t n1,
141                               const size_t n2);
142 
143 _gsl_matrix_complex_long_double_view
144 gsl_matrix_complex_long_double_view_vector_with_tda (gsl_vector_complex_long_double * v,
145                                        const size_t n1,
146                                        const size_t n2,
147                                        const size_t tda);
148 
149 
150 _gsl_matrix_complex_long_double_const_view
151 gsl_matrix_complex_long_double_const_submatrix (const gsl_matrix_complex_long_double * m,
152                                   const size_t i, const size_t j,
153                                   const size_t n1, const size_t n2);
154 
155 _gsl_vector_complex_long_double_const_view
156 gsl_matrix_complex_long_double_const_row (const gsl_matrix_complex_long_double * m,
157                             const size_t i);
158 
159 _gsl_vector_complex_long_double_const_view
160 gsl_matrix_complex_long_double_const_column (const gsl_matrix_complex_long_double * m,
161                                const size_t j);
162 
163 _gsl_vector_complex_long_double_const_view
164 gsl_matrix_complex_long_double_const_diagonal (const gsl_matrix_complex_long_double * m);
165 
166 _gsl_vector_complex_long_double_const_view
167 gsl_matrix_complex_long_double_const_subdiagonal (const gsl_matrix_complex_long_double * m,
168                                     const size_t k);
169 
170 _gsl_vector_complex_long_double_const_view
171 gsl_matrix_complex_long_double_const_superdiagonal (const gsl_matrix_complex_long_double * m,
172                                       const size_t k);
173 
174 _gsl_vector_complex_long_double_const_view
175 gsl_matrix_complex_long_double_const_subrow (const gsl_matrix_complex_long_double * m,
176                                        const size_t i, const size_t offset,
177                                        const size_t n);
178 
179 _gsl_vector_complex_long_double_const_view
180 gsl_matrix_complex_long_double_const_subcolumn (const gsl_matrix_complex_long_double * m,
181                                           const size_t j, const size_t offset,
182                                           const size_t n);
183 
184 _gsl_matrix_complex_long_double_const_view
185 gsl_matrix_complex_long_double_const_view_array (const long double * base,
186                                    const size_t n1,
187                                    const size_t n2);
188 
189 _gsl_matrix_complex_long_double_const_view
190 gsl_matrix_complex_long_double_const_view_array_with_tda (const long double * base,
191                                             const size_t n1,
192                                             const size_t n2,
193                                             const size_t tda);
194 
195 _gsl_matrix_complex_long_double_const_view
196 gsl_matrix_complex_long_double_const_view_vector (const gsl_vector_complex_long_double * v,
197                                     const size_t n1,
198                                     const size_t n2);
199 
200 _gsl_matrix_complex_long_double_const_view
201 gsl_matrix_complex_long_double_const_view_vector_with_tda (const gsl_vector_complex_long_double * v,
202                                              const size_t n1,
203                                              const size_t n2,
204                                              const size_t tda);
205 
206 /* Operations */
207 
208 void gsl_matrix_complex_long_double_set_zero (gsl_matrix_complex_long_double * m);
209 void gsl_matrix_complex_long_double_set_identity (gsl_matrix_complex_long_double * m);
210 void gsl_matrix_complex_long_double_set_all (gsl_matrix_complex_long_double * m, gsl_complex_long_double x);
211 
212 int gsl_matrix_complex_long_double_fread (FILE * stream, gsl_matrix_complex_long_double * m) ;
213 int gsl_matrix_complex_long_double_fwrite (FILE * stream, const gsl_matrix_complex_long_double * m) ;
214 int gsl_matrix_complex_long_double_fscanf (FILE * stream, gsl_matrix_complex_long_double * m);
215 int gsl_matrix_complex_long_double_fprintf (FILE * stream, const gsl_matrix_complex_long_double * m, const char * format);
216 
217 int gsl_matrix_complex_long_double_memcpy(gsl_matrix_complex_long_double * dest, const gsl_matrix_complex_long_double * src);
218 int gsl_matrix_complex_long_double_swap(gsl_matrix_complex_long_double * m1, gsl_matrix_complex_long_double * m2);
219 
220 int gsl_matrix_complex_long_double_swap_rows(gsl_matrix_complex_long_double * m, const size_t i, const size_t j);
221 int gsl_matrix_complex_long_double_swap_columns(gsl_matrix_complex_long_double * m, const size_t i, const size_t j);
222 int gsl_matrix_complex_long_double_swap_rowcol(gsl_matrix_complex_long_double * m, const size_t i, const size_t j);
223 
224 int gsl_matrix_complex_long_double_transpose (gsl_matrix_complex_long_double * m);
225 int gsl_matrix_complex_long_double_transpose_memcpy (gsl_matrix_complex_long_double * dest, const gsl_matrix_complex_long_double * src);
226 
227 int gsl_matrix_complex_long_double_equal (const gsl_matrix_complex_long_double * a, const gsl_matrix_complex_long_double * b);
228 
229 int gsl_matrix_complex_long_double_isnull (const gsl_matrix_complex_long_double * m);
230 int gsl_matrix_complex_long_double_ispos (const gsl_matrix_complex_long_double * m);
231 int gsl_matrix_complex_long_double_isneg (const gsl_matrix_complex_long_double * m);
232 int gsl_matrix_complex_long_double_isnonneg (const gsl_matrix_complex_long_double * m);
233 
234 int gsl_matrix_complex_long_double_add (gsl_matrix_complex_long_double * a, const gsl_matrix_complex_long_double * b);
235 int gsl_matrix_complex_long_double_sub (gsl_matrix_complex_long_double * a, const gsl_matrix_complex_long_double * b);
236 int gsl_matrix_complex_long_double_mul_elements (gsl_matrix_complex_long_double * a, const gsl_matrix_complex_long_double * b);
237 int gsl_matrix_complex_long_double_div_elements (gsl_matrix_complex_long_double * a, const gsl_matrix_complex_long_double * b);
238 int gsl_matrix_complex_long_double_scale (gsl_matrix_complex_long_double * a, const gsl_complex_long_double x);
239 int gsl_matrix_complex_long_double_add_constant (gsl_matrix_complex_long_double * a, const gsl_complex_long_double x);
240 int gsl_matrix_complex_long_double_add_diagonal (gsl_matrix_complex_long_double * a, const gsl_complex_long_double x);
241 
242 /***********************************************************************/
243 /* The functions below are obsolete                                    */
244 /***********************************************************************/
245 int gsl_matrix_complex_long_double_get_row(gsl_vector_complex_long_double * v, const gsl_matrix_complex_long_double * m, const size_t i);
246 int gsl_matrix_complex_long_double_get_col(gsl_vector_complex_long_double * v, const gsl_matrix_complex_long_double * m, const size_t j);
247 int gsl_matrix_complex_long_double_set_row(gsl_matrix_complex_long_double * m, const size_t i, const gsl_vector_complex_long_double * v);
248 int gsl_matrix_complex_long_double_set_col(gsl_matrix_complex_long_double * m, const size_t j, const gsl_vector_complex_long_double * v);
249 /***********************************************************************/
250 
251 /* inline functions if you are using GCC */
252 
253 INLINE_DECL gsl_complex_long_double gsl_matrix_complex_long_double_get(const gsl_matrix_complex_long_double * m, const size_t i, const size_t j);
254 INLINE_DECL void gsl_matrix_complex_long_double_set(gsl_matrix_complex_long_double * m, const size_t i, const size_t j, const gsl_complex_long_double x);
255 
256 INLINE_DECL gsl_complex_long_double * gsl_matrix_complex_long_double_ptr(gsl_matrix_complex_long_double * m, const size_t i, const size_t j);
257 INLINE_DECL const gsl_complex_long_double * gsl_matrix_complex_long_double_const_ptr(const gsl_matrix_complex_long_double * m, const size_t i, const size_t j);
258 
259 #ifdef HAVE_INLINE
260 
261 INLINE_FUN
262 gsl_complex_long_double
gsl_matrix_complex_long_double_get(const gsl_matrix_complex_long_double * m,const size_t i,const size_t j)263 gsl_matrix_complex_long_double_get(const gsl_matrix_complex_long_double * m,
264                      const size_t i, const size_t j)
265 {
266 #if GSL_RANGE_CHECK
267   if (GSL_RANGE_COND(1))
268     {
269       gsl_complex_long_double zero = {{0,0}};
270 
271       if (i >= m->size1)
272         {
273           GSL_ERROR_VAL("first index out of range", GSL_EINVAL, zero) ;
274         }
275       else if (j >= m->size2)
276         {
277           GSL_ERROR_VAL("second index out of range", GSL_EINVAL, zero) ;
278         }
279     }
280 #endif
281   return *(gsl_complex_long_double *)(m->data + 2*(i * m->tda + j)) ;
282 }
283 
284 INLINE_FUN
285 void
gsl_matrix_complex_long_double_set(gsl_matrix_complex_long_double * m,const size_t i,const size_t j,const gsl_complex_long_double x)286 gsl_matrix_complex_long_double_set(gsl_matrix_complex_long_double * m,
287                      const size_t i, const size_t j, const gsl_complex_long_double x)
288 {
289 #if GSL_RANGE_CHECK
290   if (GSL_RANGE_COND(1))
291     {
292       if (i >= m->size1)
293         {
294           GSL_ERROR_VOID("first index out of range", GSL_EINVAL) ;
295         }
296       else if (j >= m->size2)
297         {
298           GSL_ERROR_VOID("second index out of range", GSL_EINVAL) ;
299         }
300     }
301 #endif
302   *(gsl_complex_long_double *)(m->data + 2*(i * m->tda + j)) = x ;
303 }
304 
305 INLINE_FUN
306 gsl_complex_long_double *
gsl_matrix_complex_long_double_ptr(gsl_matrix_complex_long_double * m,const size_t i,const size_t j)307 gsl_matrix_complex_long_double_ptr(gsl_matrix_complex_long_double * m,
308                              const size_t i, const size_t j)
309 {
310 #if GSL_RANGE_CHECK
311   if (GSL_RANGE_COND(1))
312     {
313       if (i >= m->size1)
314         {
315           GSL_ERROR_NULL("first index out of range", GSL_EINVAL) ;
316         }
317       else if (j >= m->size2)
318         {
319           GSL_ERROR_NULL("second index out of range", GSL_EINVAL) ;
320         }
321     }
322 #endif
323   return (gsl_complex_long_double *)(m->data + 2*(i * m->tda + j)) ;
324 }
325 
326 INLINE_FUN
327 const gsl_complex_long_double *
gsl_matrix_complex_long_double_const_ptr(const gsl_matrix_complex_long_double * m,const size_t i,const size_t j)328 gsl_matrix_complex_long_double_const_ptr(const gsl_matrix_complex_long_double * m,
329                                    const size_t i, const size_t j)
330 {
331 #if GSL_RANGE_CHECK
332   if (GSL_RANGE_COND(1))
333     {
334       if (i >= m->size1)
335         {
336           GSL_ERROR_NULL("first index out of range", GSL_EINVAL) ;
337         }
338       else if (j >= m->size2)
339         {
340           GSL_ERROR_NULL("second index out of range", GSL_EINVAL) ;
341         }
342     }
343 #endif
344   return (const gsl_complex_long_double *)(m->data + 2*(i * m->tda + j)) ;
345 }
346 
347 #endif /* HAVE_INLINE */
348 
349 __END_DECLS
350 
351 #endif /* __GSL_MATRIX_COMPLEX_LONG_DOUBLE_H__ */
352