1 ------------------------------------------------------------------------------
3 -- GNAT RUN-TIME COMPONENTS --
5 -- SYSTEM.GENERIC_COMPLEX_BLAS --
9 -- Copyright (C) 2006-2007, Free Software Foundation, Inc. --
11 -- GNAT is free software; you can redistribute it and/or modify it under --
12 -- terms of the GNU General Public License as published by the Free Soft- --
13 -- ware Foundation; either version 2, or (at your option) any later ver- --
14 -- sion. GNAT is distributed in the hope that it will be useful, but WITH- --
15 -- OUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY --
16 -- or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License --
17 -- for more details. You should have received a copy of the GNU General --
18 -- Public License distributed with GNAT; see file COPYING. If not, write --
19 -- to the Free Software Foundation, 51 Franklin Street, Fifth Floor, --
20 -- Boston, MA 02110-1301, USA. --
22 -- As a special exception, if other files instantiate generics from this --
23 -- unit, or you link this unit with other files to produce an executable, --
24 -- this unit does not by itself cause the resulting executable to be --
25 -- covered by the GNU General Public License. This exception does not --
26 -- however invalidate any other reasons why the executable file might be --
27 -- covered by the GNU Public License. --
29 -- GNAT was originally developed by the GNAT team at New York University. --
30 -- Extensive contributions were provided by Ada Core Technologies Inc. --
32 ------------------------------------------------------------------------------
34 with Ada
.Unchecked_Conversion
; use Ada
;
35 with Interfaces
; use Interfaces
;
36 with Interfaces
.Fortran
; use Interfaces
.Fortran
;
37 with Interfaces
.Fortran
.BLAS
; use Interfaces
.Fortran
.BLAS
;
38 with System
.Generic_Array_Operations
; use System
.Generic_Array_Operations
;
40 package body System
.Generic_Complex_BLAS
is
42 Is_Single
: constant Boolean :=
43 Real
'Machine_Mantissa = Fortran
.Real
'Machine_Mantissa
44 and then Fortran
.Real
(Real
'First) = Fortran
.Real
'First
45 and then Fortran
.Real
(Real
'Last) = Fortran
.Real
'Last;
47 Is_Double
: constant Boolean :=
48 Real
'Machine_Mantissa = Double_Precision
'Machine_Mantissa
50 Double_Precision
(Real
'First) = Double_Precision
'First
52 Double_Precision
(Real
'Last) = Double_Precision
'Last;
54 subtype Complex
is Complex_Types
.Complex
;
58 function To_Double_Precision
(X
: Real
) return Double_Precision
;
59 pragma Inline
(To_Double_Precision
);
61 function To_Double_Complex
(X
: Complex
) return Double_Complex
;
62 pragma Inline
(To_Double_Complex
);
64 function To_Complex
(X
: Double_Complex
) return Complex
;
65 function To_Complex
(X
: Fortran
.Complex
) return Complex
;
66 pragma Inline
(To_Complex
);
68 function To_Fortran
(X
: Complex
) return Fortran
.Complex
;
69 pragma Inline
(To_Fortran
);
73 function To_Double_Complex
is new
74 Vector_Elementwise_Operation
75 (X_Scalar
=> Complex_Types
.Complex
,
76 Result_Scalar
=> Fortran
.Double_Complex
,
77 X_Vector
=> Complex_Vector
,
78 Result_Vector
=> BLAS
.Double_Complex_Vector
,
79 Operation
=> To_Double_Complex
);
81 function To_Complex
is new
82 Vector_Elementwise_Operation
83 (X_Scalar
=> Fortran
.Double_Complex
,
84 Result_Scalar
=> Complex
,
85 X_Vector
=> BLAS
.Double_Complex_Vector
,
86 Result_Vector
=> Complex_Vector
,
87 Operation
=> To_Complex
);
89 function To_Double_Complex
is new
90 Matrix_Elementwise_Operation
92 Result_Scalar
=> Double_Complex
,
93 X_Matrix
=> Complex_Matrix
,
94 Result_Matrix
=> BLAS
.Double_Complex_Matrix
,
95 Operation
=> To_Double_Complex
);
97 function To_Complex
is new
98 Matrix_Elementwise_Operation
99 (X_Scalar
=> Double_Complex
,
100 Result_Scalar
=> Complex
,
101 X_Matrix
=> BLAS
.Double_Complex_Matrix
,
102 Result_Matrix
=> Complex_Matrix
,
103 Operation
=> To_Complex
);
105 function To_Double_Precision
(X
: Real
) return Double_Precision
is
107 return Double_Precision
(X
);
108 end To_Double_Precision
;
110 function To_Double_Complex
(X
: Complex
) return Double_Complex
is
112 return (To_Double_Precision
(X
.Re
), To_Double_Precision
(X
.Im
));
113 end To_Double_Complex
;
115 function To_Complex
(X
: Double_Complex
) return Complex
is
117 return (Real
(X
.Re
), Real
(X
.Im
));
120 function To_Complex
(X
: Fortran
.Complex
) return Complex
is
122 return (Real
(X
.Re
), Real
(X
.Im
));
125 function To_Fortran
(X
: Complex
) return Fortran
.Complex
is
127 return (Fortran
.Real
(X
.Re
), Fortran
.Real
(X
.Im
));
137 Inc_X
: Integer := 1;
139 Inc_Y
: Integer := 1) return Complex
144 type X_Ptr
is access all BLAS
.Complex_Vector
(X
'Range);
145 type Y_Ptr
is access all BLAS
.Complex_Vector
(Y
'Range);
146 function Conv_X
is new Unchecked_Conversion
(Address
, X_Ptr
);
147 function Conv_Y
is new Unchecked_Conversion
(Address
, Y_Ptr
);
149 return To_Complex
(BLAS
.cdotu
(N
, Conv_X
(X
'Address).all, Inc_X
,
150 Conv_Y
(Y
'Address).all, Inc_Y
));
155 type X_Ptr
is access all BLAS
.Double_Complex_Vector
(X
'Range);
156 type Y_Ptr
is access all BLAS
.Double_Complex_Vector
(Y
'Range);
157 function Conv_X
is new Unchecked_Conversion
(Address
, X_Ptr
);
158 function Conv_Y
is new Unchecked_Conversion
(Address
, Y_Ptr
);
160 return To_Complex
(BLAS
.zdotu
(N
, Conv_X
(X
'Address).all, Inc_X
,
161 Conv_Y
(Y
'Address).all, Inc_Y
));
165 return To_Complex
(BLAS
.zdotu
(N
, To_Double_Complex
(X
), Inc_X
,
166 To_Double_Complex
(Y
), Inc_Y
));
175 (Trans_A
: access constant Character;
176 Trans_B
: access constant Character;
180 Alpha
: Complex
:= (1.0, 0.0);
185 Beta
: Complex
:= (0.0, 0.0);
186 C
: in out Complex_Matrix
;
192 subtype A_Type
is BLAS
.Complex_Matrix
(A
'Range (1), A
'Range (2));
193 subtype B_Type
is BLAS
.Complex_Matrix
(B
'Range (1), B
'Range (2));
195 access all BLAS
.Complex_Matrix
(C
'Range (1), C
'Range (2));
197 new Unchecked_Conversion
(Complex_Matrix
, A_Type
);
199 new Unchecked_Conversion
(Complex_Matrix
, B_Type
);
201 new Unchecked_Conversion
(Address
, C_Ptr
);
203 BLAS
.cgemm
(Trans_A
, Trans_B
, M
, N
, K
, To_Fortran
(Alpha
),
204 Conv_A
(A
), Ld_A
, Conv_B
(B
), Ld_B
, To_Fortran
(Beta
),
205 Conv_C
(C
'Address).all, Ld_C
);
211 BLAS
.Double_Complex_Matrix
(A
'Range (1), A
'Range (2));
213 BLAS
.Double_Complex_Matrix
(B
'Range (1), B
'Range (2));
214 type C_Ptr
is access all
215 BLAS
.Double_Complex_Matrix
(C
'Range (1), C
'Range (2));
217 new Unchecked_Conversion
(Complex_Matrix
, A_Type
);
219 new Unchecked_Conversion
(Complex_Matrix
, B_Type
);
220 function Conv_C
is new Unchecked_Conversion
(Address
, C_Ptr
);
222 BLAS
.zgemm
(Trans_A
, Trans_B
, M
, N
, K
, To_Double_Complex
(Alpha
),
223 Conv_A
(A
), Ld_A
, Conv_B
(B
), Ld_B
,
224 To_Double_Complex
(Beta
),
225 Conv_C
(C
'Address).all, Ld_C
);
230 DP_C
: BLAS
.Double_Complex_Matrix
(C
'Range (1), C
'Range (2));
232 if Beta
.Re
/= 0.0 or else Beta
.Im
/= 0.0 then
233 DP_C
:= To_Double_Complex
(C
);
236 BLAS
.zgemm
(Trans_A
, Trans_B
, M
, N
, K
, To_Double_Complex
(Alpha
),
237 To_Double_Complex
(A
), Ld_A
,
238 To_Double_Complex
(B
), Ld_B
, To_Double_Complex
(Beta
),
241 C
:= To_Complex
(DP_C
);
251 (Trans
: access constant Character;
254 Alpha
: Complex
:= (1.0, 0.0);
258 Inc_X
: Integer := 1;
259 Beta
: Complex
:= (0.0, 0.0);
260 Y
: in out Complex_Vector
;
261 Inc_Y
: Integer := 1)
266 subtype A_Type
is BLAS
.Complex_Matrix
(A
'Range (1), A
'Range (2));
267 subtype X_Type
is BLAS
.Complex_Vector
(X
'Range);
268 type Y_Ptr
is access all BLAS
.Complex_Vector
(Y
'Range);
270 new Unchecked_Conversion
(Complex_Matrix
, A_Type
);
272 new Unchecked_Conversion
(Complex_Vector
, X_Type
);
274 new Unchecked_Conversion
(Address
, Y_Ptr
);
276 BLAS
.cgemv
(Trans
, M
, N
, To_Fortran
(Alpha
),
277 Conv_A
(A
), Ld_A
, Conv_X
(X
), Inc_X
, To_Fortran
(Beta
),
278 Conv_Y
(Y
'Address).all, Inc_Y
);
284 BLAS
.Double_Complex_Matrix
(A
'Range (1), A
'Range (2));
286 BLAS
.Double_Complex_Vector
(X
'Range);
287 type Y_Ptr
is access all BLAS
.Double_Complex_Vector
(Y
'Range);
289 new Unchecked_Conversion
(Complex_Matrix
, A_Type
);
291 new Unchecked_Conversion
(Complex_Vector
, X_Type
);
293 new Unchecked_Conversion
(Address
, Y_Ptr
);
295 BLAS
.zgemv
(Trans
, M
, N
, To_Double_Complex
(Alpha
),
296 Conv_A
(A
), Ld_A
, Conv_X
(X
), Inc_X
,
297 To_Double_Complex
(Beta
),
298 Conv_Y
(Y
'Address).all, Inc_Y
);
303 DP_Y
: BLAS
.Double_Complex_Vector
(Y
'Range);
305 if Beta
.Re
/= 0.0 or else Beta
.Im
/= 0.0 then
306 DP_Y
:= To_Double_Complex
(Y
);
309 BLAS
.zgemv
(Trans
, M
, N
, To_Double_Complex
(Alpha
),
310 To_Double_Complex
(A
), Ld_A
,
311 To_Double_Complex
(X
), Inc_X
, To_Double_Complex
(Beta
),
314 Y
:= To_Complex
(DP_Y
);
326 Inc_X
: Integer := 1) return Real
331 subtype X_Type
is BLAS
.Complex_Vector
(X
'Range);
333 new Unchecked_Conversion
(Complex_Vector
, X_Type
);
335 return Real
(BLAS
.scnrm2
(N
, Conv_X
(X
), Inc_X
));
340 subtype X_Type
is BLAS
.Double_Complex_Vector
(X
'Range);
342 new Unchecked_Conversion
(Complex_Vector
, X_Type
);
344 return Real
(BLAS
.dznrm2
(N
, Conv_X
(X
), Inc_X
));
348 return Real
(BLAS
.dznrm2
(N
, To_Double_Complex
(X
), Inc_X
));
352 end System
.Generic_Complex_BLAS
;