5164
|
1 /* ========================================================================== */ |
|
2 /* === UMF_scale ============================================================ */ |
|
3 /* ========================================================================== */ |
|
4 |
|
5 /* -------------------------------------------------------------------------- */ |
|
6 /* UMFPACK Version 4.4, Copyright (c) 2005 by Timothy A. Davis. CISE Dept, */ |
|
7 /* Univ. of Florida. All Rights Reserved. See ../Doc/License for License. */ |
|
8 /* web: http://www.cise.ufl.edu/research/sparse/umfpack */ |
|
9 /* -------------------------------------------------------------------------- */ |
|
10 |
|
11 /* Divide a vector of stride 1 by the pivot value. */ |
|
12 |
|
13 #include "umf_internal.h" |
|
14 |
|
15 GLOBAL void UMF_scale |
|
16 ( |
|
17 Int n, |
|
18 Entry pivot, |
|
19 Entry X [ ] |
|
20 ) |
|
21 { |
|
22 Entry x ; |
|
23 double s ; |
|
24 Int i ; |
|
25 |
|
26 /* ---------------------------------------------------------------------- */ |
|
27 /* compute the approximate absolute value of the pivot, and select method */ |
|
28 /* ---------------------------------------------------------------------- */ |
|
29 |
|
30 APPROX_ABS (s, pivot) ; |
|
31 |
|
32 if (s < RECIPROCAL_TOLERANCE || IS_NAN (pivot)) |
|
33 { |
|
34 /* ------------------------------------------------------------------ */ |
|
35 /* tiny, or zero, pivot case */ |
|
36 /* ------------------------------------------------------------------ */ |
|
37 |
|
38 /* The pivot is tiny, or NaN. Do not divide zero by the pivot value, |
|
39 * and do not multiply by 1/pivot, either. */ |
|
40 |
|
41 for (i = 0 ; i < n ; i++) |
|
42 { |
|
43 /* X [i] /= pivot ; */ |
|
44 x = X [i] ; |
|
45 |
|
46 #ifndef NO_DIVIDE_BY_ZERO |
|
47 if (IS_NONZERO (x)) |
|
48 { |
|
49 DIV (X [i], x, pivot) ; |
|
50 } |
|
51 #else |
|
52 /* Do not divide by zero */ |
|
53 if (IS_NONZERO (x) && IS_NONZERO (pivot)) |
|
54 { |
|
55 DIV (X [i], x, pivot) ; |
|
56 } |
|
57 #endif |
|
58 |
|
59 } |
|
60 |
|
61 } |
|
62 else |
|
63 { |
|
64 |
|
65 /* ------------------------------------------------------------------ */ |
|
66 /* normal case. select the x/pivot or x * (1/pivot) method */ |
|
67 /* ------------------------------------------------------------------ */ |
|
68 |
|
69 /* The pivot is not tiny, and is not NaN. Don't bother to check for |
|
70 * zeros in the pivot column, X. */ |
|
71 |
|
72 #if !defined (NRECIPROCAL) && !(defined (__GNUC__) && defined (COMPLEX)) |
|
73 |
|
74 /* -------------------------------------------------------------- */ |
|
75 /* multiply x by (1/pivot) */ |
|
76 /* -------------------------------------------------------------- */ |
|
77 |
|
78 /* Slightly less accurate, but faster. It allows the use of |
|
79 * the level-1 BLAS dscal or zscal routine. This not used when |
|
80 * UMFPACK is used in MATLAB (either as a built-in routine, or as |
|
81 * a mexFunction). |
|
82 * |
|
83 * Using gcc version 3.2 can cause the following code to fail for |
|
84 * some complex matrices (not all), with or without the BLAS. This |
|
85 * was found in Red Hat Linux 7.3 on a Dell Latitude C840 with a |
|
86 * Pentium 4M. Thus, this code is not used when gcc is used, for |
|
87 * the complex case. |
|
88 * |
|
89 * It works just fine with Intel's icc compiler, version 7.0. |
|
90 */ |
|
91 |
|
92 /* pivot = 1 / pivot */ |
|
93 RECIPROCAL (pivot) ; |
|
94 |
|
95 #if defined (USE_NO_BLAS) |
|
96 for (i = 0 ; i < n ; i++) |
|
97 { |
|
98 /* X [i] *= pivot ; */ |
|
99 x = X [i] ; |
|
100 MULT (X [i], x, pivot) ; |
|
101 } |
|
102 #else |
|
103 BLAS_SCAL (n, pivot, X) ; |
|
104 #endif |
|
105 |
|
106 #else |
|
107 |
|
108 /* -------------------------------------------------------------- */ |
|
109 /* divide x by the pivot */ |
|
110 /* -------------------------------------------------------------- */ |
|
111 |
|
112 /* This is slightly more accurate, particularly if the pivot column |
|
113 * consists of only IEEE subnormals. Always do this if UMFPACK is |
|
114 * being compiled as a built-in routine or mexFunction in MATLAB, |
|
115 * or if gcc is being used with complex matrices. */ |
|
116 |
|
117 for (i = 0 ; i < n ; i++) |
|
118 { |
|
119 /* X [i] /= pivot ; */ |
|
120 x = X [i] ; |
|
121 DIV (X [i], x, pivot) ; |
|
122 } |
|
123 |
|
124 #endif |
|
125 |
|
126 } |
|
127 } |