IFPACK  Development
 All Classes Namespaces Files Functions Variables Enumerations Friends Pages
blas_dh.c
1 /*@HEADER
2 // ***********************************************************************
3 //
4 // Ifpack: Object-Oriented Algebraic Preconditioner Package
5 // Copyright (2002) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ***********************************************************************
40 //@HEADER
41 */
42 
43 #include "blas_dh.h"
44 
45 #undef __FUNC__
46 #define __FUNC__ "matvec_euclid_seq"
47 void
48 matvec_euclid_seq (int n, int *rp, int *cval, double *aval, double *x,
49  double *y)
50 {
51  START_FUNC_DH int i, j;
52  int from, to, col;
53  double sum;
54 
55  if (np_dh > 1)
56  SET_V_ERROR ("only for sequential case!\n");
57 
58  {
59  for (i = 0; i < n; ++i)
60  {
61  sum = 0.0;
62  from = rp[i];
63  to = rp[i + 1];
64  for (j = from; j < to; ++j)
65  {
66  col = cval[j];
67  sum += (aval[j] * x[col]);
68  }
69  y[i] = sum;
70  }
71  }
72 END_FUNC_DH}
73 
74 #undef __FUNC__
75 #define __FUNC__ "Axpy"
76 void
77 Axpy (int n, double alpha, double *x, double *y)
78 {
79  START_FUNC_DH int i;
80 
81  for (i = 0; i < n; ++i)
82  {
83  y[i] = alpha * x[i] + y[i];
84  }
85 END_FUNC_DH}
86 
87 
88 #undef __FUNC__
89 #define __FUNC__ "CopyVec"
90 void
91 CopyVec (int n, double *xIN, double *yOUT)
92 {
93  START_FUNC_DH int i;
94 
95  for (i = 0; i < n; ++i)
96  {
97  yOUT[i] = xIN[i];
98  }
99 END_FUNC_DH}
100 
101 
102 #undef __FUNC__
103 #define __FUNC__ "ScaleVec"
104 void
105 ScaleVec (int n, double alpha, double *x)
106 {
107  START_FUNC_DH int i;
108 
109  for (i = 0; i < n; ++i)
110  {
111  x[i] *= alpha;
112  }
113 END_FUNC_DH}
114 
115 #undef __FUNC__
116 #define __FUNC__ "InnerProd"
117 double
118 InnerProd (int n, double *x, double *y)
119 {
120  START_FUNC_DH double result, local_result = 0.0;
121 
122  int i;
123 
124  for (i = 0; i < n; ++i)
125  {
126  local_result += x[i] * y[i];
127  }
128 
129  if (np_dh > 1)
130  {
131  MPI_Allreduce (&local_result, &result, 1, MPI_DOUBLE, MPI_SUM, comm_dh);
132  }
133  else
134  {
135  result = local_result;
136  }
137 
138 END_FUNC_VAL (result)}
139 
140 #undef __FUNC__
141 #define __FUNC__ "Norm2"
142 double
143 Norm2 (int n, double *x)
144 {
145  START_FUNC_DH double result, local_result = 0.0;
146  int i;
147 
148  for (i = 0; i < n; ++i)
149  {
150  local_result += (x[i] * x[i]);
151  }
152 
153  if (np_dh > 1)
154  {
155  MPI_Allreduce (&local_result, &result, 1, MPI_DOUBLE, MPI_SUM, comm_dh);
156  }
157  else
158  {
159  result = local_result;
160  }
161  result = sqrt (result);
162 END_FUNC_VAL (result)}