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#include <gsl/gsl_blas_types.h>
30
31#undef __BEGIN_DECLS
32#undef __END_DECLS
33#ifdef __cplusplus
34# define __BEGIN_DECLS extern "C" {
35# define __END_DECLS }
36#else
37# define __BEGIN_DECLS /* empty */
38# define __END_DECLS /* empty */
39#endif
40
41__BEGIN_DECLS
42
43typedef struct
44{
45 size_t size1;
46 size_t size2;
47 size_t tda;
48 long double * data;
49 gsl_block_complex_long_double * block;
50 int owner;
51} gsl_matrix_complex_long_double ;
52
53typedef struct
54{
55 gsl_matrix_complex_long_double matrix;
56} _gsl_matrix_complex_long_double_view;
57
58typedef _gsl_matrix_complex_long_double_view gsl_matrix_complex_long_double_view;
59
60typedef struct
61{
62 gsl_matrix_complex_long_double matrix;
63} _gsl_matrix_complex_long_double_const_view;
64
65typedef const _gsl_matrix_complex_long_double_const_view gsl_matrix_complex_long_double_const_view;
66
67
68/* Allocation */
69
70gsl_matrix_complex_long_double *
71gsl_matrix_complex_long_double_alloc (const size_t n1, const size_t n2);
72
73gsl_matrix_complex_long_double *
74gsl_matrix_complex_long_double_calloc (const size_t n1, const size_t n2);
75
76gsl_matrix_complex_long_double *
77gsl_matrix_complex_long_double_alloc_from_block (gsl_block_complex_long_double * b,
78 const size_t offset,
79 const size_t n1, const size_t n2, const size_t d2);
80
81gsl_matrix_complex_long_double *
82gsl_matrix_complex_long_double_alloc_from_matrix (gsl_matrix_complex_long_double * b,
83 const size_t k1, const size_t k2,
84 const size_t n1, const size_t n2);
85
86gsl_vector_complex_long_double *
87gsl_vector_complex_long_double_alloc_row_from_matrix (gsl_matrix_complex_long_double * m,
88 const size_t i);
89
90gsl_vector_complex_long_double *
91gsl_vector_complex_long_double_alloc_col_from_matrix (gsl_matrix_complex_long_double * m,
92 const size_t j);
93
94void gsl_matrix_complex_long_double_free (gsl_matrix_complex_long_double * m);
95
96/* Views */
97
98_gsl_matrix_complex_long_double_view
99gsl_matrix_complex_long_double_submatrix (gsl_matrix_complex_long_double * m,
100 const size_t i, const size_t j,
101 const size_t n1, const size_t n2);
102
103_gsl_vector_complex_long_double_view
104gsl_matrix_complex_long_double_row (gsl_matrix_complex_long_double * m, const size_t i);
105
106_gsl_vector_complex_long_double_view
107gsl_matrix_complex_long_double_column (gsl_matrix_complex_long_double * m, const size_t j);
108
109_gsl_vector_complex_long_double_view
110gsl_matrix_complex_long_double_diagonal (gsl_matrix_complex_long_double * m);
111
112_gsl_vector_complex_long_double_view
113gsl_matrix_complex_long_double_subdiagonal (gsl_matrix_complex_long_double * m, const size_t k);
114
115_gsl_vector_complex_long_double_view
116gsl_matrix_complex_long_double_superdiagonal (gsl_matrix_complex_long_double * m, const size_t k);
117
118_gsl_vector_complex_long_double_view
119gsl_matrix_complex_long_double_subrow (gsl_matrix_complex_long_double * m,
120 const size_t i, const size_t offset,
121 const size_t n);
122
123_gsl_vector_complex_long_double_view
124gsl_matrix_complex_long_double_subcolumn (gsl_matrix_complex_long_double * m,
125 const size_t j, const size_t offset,
126 const size_t n);
127
128_gsl_matrix_complex_long_double_view
129gsl_matrix_complex_long_double_view_array (long double * base,
130 const size_t n1,
131 const size_t n2);
132
133_gsl_matrix_complex_long_double_view
134gsl_matrix_complex_long_double_view_array_with_tda (long double * base,
135 const size_t n1,
136 const size_t n2,
137 const size_t tda);
138
139_gsl_matrix_complex_long_double_view
140gsl_matrix_complex_long_double_view_vector (gsl_vector_complex_long_double * v,
141 const size_t n1,
142 const size_t n2);
143
144_gsl_matrix_complex_long_double_view
145gsl_matrix_complex_long_double_view_vector_with_tda (gsl_vector_complex_long_double * v,
146 const size_t n1,
147 const size_t n2,
148 const size_t tda);
149
150
151_gsl_matrix_complex_long_double_const_view
152gsl_matrix_complex_long_double_const_submatrix (const gsl_matrix_complex_long_double * m,
153 const size_t i, const size_t j,
154 const size_t n1, const size_t n2);
155
156_gsl_vector_complex_long_double_const_view
157gsl_matrix_complex_long_double_const_row (const gsl_matrix_complex_long_double * m,
158 const size_t i);
159
160_gsl_vector_complex_long_double_const_view
161gsl_matrix_complex_long_double_const_column (const gsl_matrix_complex_long_double * m,
162 const size_t j);
163
164_gsl_vector_complex_long_double_const_view
165gsl_matrix_complex_long_double_const_diagonal (const gsl_matrix_complex_long_double * m);
166
167_gsl_vector_complex_long_double_const_view
168gsl_matrix_complex_long_double_const_subdiagonal (const gsl_matrix_complex_long_double * m,
169 const size_t k);
170
171_gsl_vector_complex_long_double_const_view
172gsl_matrix_complex_long_double_const_superdiagonal (const gsl_matrix_complex_long_double * m,
173 const size_t k);
174
175_gsl_vector_complex_long_double_const_view
176gsl_matrix_complex_long_double_const_subrow (const gsl_matrix_complex_long_double * m,
177 const size_t i, const size_t offset,
178 const size_t n);
179
180_gsl_vector_complex_long_double_const_view
181gsl_matrix_complex_long_double_const_subcolumn (const gsl_matrix_complex_long_double * m,
182 const size_t j, const size_t offset,
183 const size_t n);
184
185_gsl_matrix_complex_long_double_const_view
186gsl_matrix_complex_long_double_const_view_array (const long double * base,
187 const size_t n1,
188 const size_t n2);
189
190_gsl_matrix_complex_long_double_const_view
191gsl_matrix_complex_long_double_const_view_array_with_tda (const long double * base,
192 const size_t n1,
193 const size_t n2,
194 const size_t tda);
195
196_gsl_matrix_complex_long_double_const_view
197gsl_matrix_complex_long_double_const_view_vector (const gsl_vector_complex_long_double * v,
198 const size_t n1,
199 const size_t n2);
200
201_gsl_matrix_complex_long_double_const_view
202gsl_matrix_complex_long_double_const_view_vector_with_tda (const gsl_vector_complex_long_double * v,
203 const size_t n1,
204 const size_t n2,
205 const size_t tda);
206
207/* Operations */
208
209void gsl_matrix_complex_long_double_set_zero (gsl_matrix_complex_long_double * m);
210void gsl_matrix_complex_long_double_set_identity (gsl_matrix_complex_long_double * m);
211void gsl_matrix_complex_long_double_set_all (gsl_matrix_complex_long_double * m, gsl_complex_long_double x);
212
213int gsl_matrix_complex_long_double_fread (FILE * stream, gsl_matrix_complex_long_double * m) ;
214int gsl_matrix_complex_long_double_fwrite (FILE * stream, const gsl_matrix_complex_long_double * m) ;
215int gsl_matrix_complex_long_double_fscanf (FILE * stream, gsl_matrix_complex_long_double * m);
216int gsl_matrix_complex_long_double_fprintf (FILE * stream, const gsl_matrix_complex_long_double * m, const char * format);
217
218int gsl_matrix_complex_long_double_memcpy(gsl_matrix_complex_long_double * dest, const gsl_matrix_complex_long_double * src);
219int gsl_matrix_complex_long_double_swap(gsl_matrix_complex_long_double * m1, gsl_matrix_complex_long_double * m2);
220int gsl_matrix_complex_long_double_tricpy(CBLAS_UPLO_t Uplo, CBLAS_DIAG_t Diag, gsl_matrix_complex_long_double * dest, const gsl_matrix_complex_long_double * src);
221
222int gsl_matrix_complex_long_double_swap_rows(gsl_matrix_complex_long_double * m, const size_t i, const size_t j);
223int gsl_matrix_complex_long_double_swap_columns(gsl_matrix_complex_long_double * m, const size_t i, const size_t j);
224int gsl_matrix_complex_long_double_swap_rowcol(gsl_matrix_complex_long_double * m, const size_t i, const size_t j);
225
226int gsl_matrix_complex_long_double_transpose (gsl_matrix_complex_long_double * m);
227int gsl_matrix_complex_long_double_transpose_memcpy (gsl_matrix_complex_long_double * dest, const gsl_matrix_complex_long_double * src);
228int gsl_matrix_complex_long_double_transpose_tricpy(CBLAS_UPLO_t Uplo_src, CBLAS_DIAG_t Diag, gsl_matrix_complex_long_double * dest, const gsl_matrix_complex_long_double * src);
229
230int gsl_matrix_complex_long_double_conjtrans_memcpy (gsl_matrix_complex_long_double * dest, const gsl_matrix_complex_long_double * src);
231
232int gsl_matrix_complex_long_double_equal (const gsl_matrix_complex_long_double * a, const gsl_matrix_complex_long_double * b);
233
234int gsl_matrix_complex_long_double_isnull (const gsl_matrix_complex_long_double * m);
235int gsl_matrix_complex_long_double_ispos (const gsl_matrix_complex_long_double * m);
236int gsl_matrix_complex_long_double_isneg (const gsl_matrix_complex_long_double * m);
237int gsl_matrix_complex_long_double_isnonneg (const gsl_matrix_complex_long_double * m);
238
239int gsl_matrix_complex_long_double_add (gsl_matrix_complex_long_double * a, const gsl_matrix_complex_long_double * b);
240int gsl_matrix_complex_long_double_sub (gsl_matrix_complex_long_double * a, const gsl_matrix_complex_long_double * b);
241int gsl_matrix_complex_long_double_mul_elements (gsl_matrix_complex_long_double * a, const gsl_matrix_complex_long_double * b);
242int gsl_matrix_complex_long_double_div_elements (gsl_matrix_complex_long_double * a, const gsl_matrix_complex_long_double * b);
243int gsl_matrix_complex_long_double_scale (gsl_matrix_complex_long_double * a, const gsl_complex_long_double x);
244int gsl_matrix_complex_long_double_scale_rows (gsl_matrix_complex_long_double * a, const gsl_vector_complex_long_double * x);
245int gsl_matrix_complex_long_double_scale_columns (gsl_matrix_complex_long_double * a, const gsl_vector_complex_long_double * x);
246int gsl_matrix_complex_long_double_add_constant (gsl_matrix_complex_long_double * a, const gsl_complex_long_double x);
247int gsl_matrix_complex_long_double_add_diagonal (gsl_matrix_complex_long_double * a, const gsl_complex_long_double x);
248
249/***********************************************************************/
250/* The functions below are obsolete */
251/***********************************************************************/
252int gsl_matrix_complex_long_double_get_row(gsl_vector_complex_long_double * v, const gsl_matrix_complex_long_double * m, const size_t i);
253int gsl_matrix_complex_long_double_get_col(gsl_vector_complex_long_double * v, const gsl_matrix_complex_long_double * m, const size_t j);
254int gsl_matrix_complex_long_double_set_row(gsl_matrix_complex_long_double * m, const size_t i, const gsl_vector_complex_long_double * v);
255int gsl_matrix_complex_long_double_set_col(gsl_matrix_complex_long_double * m, const size_t j, const gsl_vector_complex_long_double * v);
256/***********************************************************************/
257
258/* inline functions if you are using GCC */
259
260INLINE_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);
261INLINE_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);
262
263INLINE_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);
264INLINE_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);
265
266#ifdef HAVE_INLINE
267
268INLINE_FUN
269gsl_complex_long_double
270gsl_matrix_complex_long_double_get(const gsl_matrix_complex_long_double * m,
271 const size_t i, const size_t j)
272{
273#if GSL_RANGE_CHECK
274 if (GSL_RANGE_COND(1))
275 {
276 gsl_complex_long_double zero = {{0,0}};
277
278 if (i >= m->size1)
279 {
280 GSL_ERROR_VAL("first index out of range", GSL_EINVAL, zero) ;
281 }
282 else if (j >= m->size2)
283 {
284 GSL_ERROR_VAL("second index out of range", GSL_EINVAL, zero) ;
285 }
286 }
287#endif
288 return *(gsl_complex_long_double *)(m->data + 2*(i * m->tda + j)) ;
289}
290
291INLINE_FUN
292void
293gsl_matrix_complex_long_double_set(gsl_matrix_complex_long_double * m,
294 const size_t i, const size_t j, const gsl_complex_long_double x)
295{
296#if GSL_RANGE_CHECK
297 if (GSL_RANGE_COND(1))
298 {
299 if (i >= m->size1)
300 {
301 GSL_ERROR_VOID("first index out of range", GSL_EINVAL) ;
302 }
303 else if (j >= m->size2)
304 {
305 GSL_ERROR_VOID("second index out of range", GSL_EINVAL) ;
306 }
307 }
308#endif
309 *(gsl_complex_long_double *)(m->data + 2*(i * m->tda + j)) = x ;
310}
311
312INLINE_FUN
313gsl_complex_long_double *
314gsl_matrix_complex_long_double_ptr(gsl_matrix_complex_long_double * m,
315 const size_t i, const size_t j)
316{
317#if GSL_RANGE_CHECK
318 if (GSL_RANGE_COND(1))
319 {
320 if (i >= m->size1)
321 {
322 GSL_ERROR_NULL("first index out of range", GSL_EINVAL) ;
323 }
324 else if (j >= m->size2)
325 {
326 GSL_ERROR_NULL("second index out of range", GSL_EINVAL) ;
327 }
328 }
329#endif
330 return (gsl_complex_long_double *)(m->data + 2*(i * m->tda + j)) ;
331}
332
333INLINE_FUN
334const gsl_complex_long_double *
335gsl_matrix_complex_long_double_const_ptr(const gsl_matrix_complex_long_double * m,
336 const size_t i, const size_t j)
337{
338#if GSL_RANGE_CHECK
339 if (GSL_RANGE_COND(1))
340 {
341 if (i >= m->size1)
342 {
343 GSL_ERROR_NULL("first index out of range", GSL_EINVAL) ;
344 }
345 else if (j >= m->size2)
346 {
347 GSL_ERROR_NULL("second index out of range", GSL_EINVAL) ;
348 }
349 }
350#endif
351 return (const gsl_complex_long_double *)(m->data + 2*(i * m->tda + j)) ;
352}
353
354#endif /* HAVE_INLINE */
355
356__END_DECLS
357
358#endif /* __GSL_MATRIX_COMPLEX_LONG_DOUBLE_H__ */
359