Mercurial > hg > octave-nkf

comparison liboctave/Sparse.cc @ 8752:06b9903a029b

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
fix & clean up complex & sparse sorting
author Jaroslav Hajek <highegg@gmail.com>
date Mon, 16 Feb 2009 10:15:43 +0100
parents 095ae5e0a831
children eb63fbe60fab
comparison
equal deleted inserted replaced
8751:9f7ce4bf7650 8752:06b9903a029b
2072 2072
2073 // Can't use versions of these in Array.cc due to duplication of the 2073 // Can't use versions of these in Array.cc due to duplication of the
2074 // instantiations for Array<double and Sparse<double>, etc 2074 // instantiations for Array<double and Sparse<double>, etc
2075 template <class T> 2075 template <class T>
2076 bool 2076 bool
2077 sparse_ascending_compare (T a, T b) 2077 sparse_ascending_compare (typename ref_param<T>::type a,
2078 typename ref_param<T>::type b)
2078 { 2079 {
2079 return (a < b); 2080 return (a < b);
2080 } 2081 }
2081 2082
2082 template <class T> 2083 template <class T>
2083 bool 2084 bool
2084 sparse_descending_compare (T a, T b) 2085 sparse_descending_compare (typename ref_param<T>::type a,
2086 typename ref_param<T>::type b)
2085 { 2087 {
2086 return (a > b); 2088 return (a > b);
2087 }
2088
2089 template <class T>
2090 bool
2091 sparse_ascending_compare (vec_index<T> *a, vec_index<T> *b)
2092 {
2093 return (a->vec < b->vec);
2094 }
2095
2096 template <class T>
2097 bool
2098 sparse_descending_compare (vec_index<T> *a, vec_index<T> *b)
2099 {
2100 return (a->vec > b->vec);
2101 } 2089 }
2102 2090
2103 template <class T> 2091 template <class T>
2104 Sparse<T> 2092 Sparse<T>
2105 Sparse<T>::sort (octave_idx_type dim, sortmode mode) const 2093 Sparse<T>::sort (octave_idx_type dim, sortmode mode) const
2120 } 2108 }
2121 2109
2122 octave_sort<T> lsort; 2110 octave_sort<T> lsort;
2123 2111
2124 if (mode == ASCENDING) 2112 if (mode == ASCENDING)
2125 lsort.set_compare (sparse_ascending_compare); 2113 lsort.set_compare (sparse_ascending_compare<T>);
2126 else if (mode == DESCENDING) 2114 else if (mode == DESCENDING)
2127 lsort.set_compare (sparse_descending_compare); 2115 lsort.set_compare (sparse_descending_compare<T>);
2128 else 2116 else
2129 abort (); 2117 abort ();
2130 2118
2131 T *v = m.data (); 2119 T *v = m.data ();
2132 octave_idx_type *mcidx = m.cidx (); 2120 octave_idx_type *mcidx = m.cidx ();
2139 2127
2140 octave_idx_type i; 2128 octave_idx_type i;
2141 if (mode == ASCENDING) 2129 if (mode == ASCENDING)
2142 { 2130 {
2143 for (i = 0; i < ns; i++) 2131 for (i = 0; i < ns; i++)
2144 if (sparse_ascending_compare (static_cast<T> (0), v [i])) 2132 if (sparse_ascending_compare<T> (static_cast<T> (0), v [i]))
2145 break; 2133 break;
2146 } 2134 }
2147 else 2135 else
2148 { 2136 {
2149 for (i = 0; i < ns; i++) 2137 for (i = 0; i < ns; i++)
2150 if (sparse_descending_compare (static_cast<T> (0), v [i])) 2138 if (sparse_descending_compare<T> (static_cast<T> (0), v [i]))
2151 break; 2139 break;
2152 } 2140 }
2153 for (octave_idx_type k = 0; k < i; k++) 2141 for (octave_idx_type k = 0; k < i; k++)
2154 mridx [k] = k; 2142 mridx [k] = k;
2155 for (octave_idx_type k = i; k < ns; k++) 2143 for (octave_idx_type k = i; k < ns; k++)
2186 m = m.transpose (); 2174 m = m.transpose ();
2187 nr = m.rows (); 2175 nr = m.rows ();
2188 nc = m.columns (); 2176 nc = m.columns ();
2189 } 2177 }
2190 2178
2191 octave_sort<vec_index<T> *> indexed_sort; 2179 octave_sort<T> indexed_sort;
2192 2180
2193 if (mode == ASCENDING) 2181 if (mode == ASCENDING)
2194 indexed_sort.set_compare (sparse_ascending_compare); 2182 indexed_sort.set_compare (sparse_ascending_compare<T>);
2195 else if (mode == DESCENDING) 2183 else if (mode == DESCENDING)
2196 indexed_sort.set_compare (sparse_descending_compare); 2184 indexed_sort.set_compare (sparse_descending_compare<T>);
2197 else 2185 else
2198 abort (); 2186 abort ();
2199 2187
2200 T *v = m.data (); 2188 T *v = m.data ();
2201 octave_idx_type *mcidx = m.cidx (); 2189 octave_idx_type *mcidx = m.cidx ();
2202 octave_idx_type *mridx = m.ridx (); 2190 octave_idx_type *mridx = m.ridx ();
2203 2191
2204 OCTAVE_LOCAL_BUFFER (vec_index<T> *, vi, nr);
2205 OCTAVE_LOCAL_BUFFER (vec_index<T>, vix, nr);
2206
2207 for (octave_idx_type i = 0; i < nr; i++)
2208 vi[i] = &vix[i];
2209
2210 sidx = Array<octave_idx_type> (dim_vector (nr, nc)); 2192 sidx = Array<octave_idx_type> (dim_vector (nr, nc));
2193 OCTAVE_LOCAL_BUFFER (octave_idx_type, vi, nr);
2211 2194
2212 for (octave_idx_type j = 0; j < nc; j++) 2195 for (octave_idx_type j = 0; j < nc; j++)
2213 { 2196 {
2214 octave_idx_type ns = mcidx [j + 1] - mcidx [j]; 2197 octave_idx_type ns = mcidx [j + 1] - mcidx [j];
2215 octave_idx_type offset = j * nr; 2198 octave_idx_type offset = j * nr;
2220 sidx (offset + k) = k; 2203 sidx (offset + k) = k;
2221 } 2204 }
2222 else 2205 else
2223 { 2206 {
2224 for (octave_idx_type i = 0; i < ns; i++) 2207 for (octave_idx_type i = 0; i < ns; i++)
2225 { 2208 vi[i] = mridx[i];
2226 vi[i]->vec = v[i]; 2209
2227 vi[i]->indx = mridx[i]; 2210 indexed_sort.sort (v, vi, ns);
2228 }
2229
2230 indexed_sort.sort (vi, ns);
2231 2211
2232 octave_idx_type i; 2212 octave_idx_type i;
2233 if (mode == ASCENDING) 2213 if (mode == ASCENDING)
2234 { 2214 {
2235 for (i = 0; i < ns; i++) 2215 for (i = 0; i < ns; i++)
2236 if (sparse_ascending_compare (static_cast<T> (0), 2216 if (sparse_ascending_compare<T> (static_cast<T> (0), v[i]))
2237 vi [i] -> vec))
2238 break; 2217 break;
2239 } 2218 }
2240 else 2219 else
2241 { 2220 {
2242 for (i = 0; i < ns; i++) 2221 for (i = 0; i < ns; i++)
2243 if (sparse_descending_compare (static_cast<T> (0), 2222 if (sparse_descending_compare<T> (static_cast<T> (0), v[i]))
2244 vi [i] -> vec))
2245 break; 2223 break;
2246 } 2224 }
2247 2225
2248 octave_idx_type ii = 0; 2226 octave_idx_type ii = 0;
2249 octave_idx_type jj = i; 2227 octave_idx_type jj = i;
2255 sidx (offset + jj++) = k; 2233 sidx (offset + jj++) = k;
2256 } 2234 }
2257 2235
2258 for (octave_idx_type k = 0; k < i; k++) 2236 for (octave_idx_type k = 0; k < i; k++)
2259 { 2237 {
2260 v [k] = vi [k] -> vec; 2238 sidx (k + offset) = vi [k];
2261 sidx (k + offset) = vi [k] -> indx;
2262 mridx [k] = k; 2239 mridx [k] = k;
2263 } 2240 }
2264 2241
2265 for (octave_idx_type k = i; k < ns; k++) 2242 for (octave_idx_type k = i; k < ns; k++)
2266 { 2243 {
2267 v [k] = vi [k] -> vec; 2244 sidx (k - ns + nr + offset) = vi [k];
2268 sidx (k - ns + nr + offset) = vi [k] -> indx;
2269 mridx [k] = k - ns + nr; 2245 mridx [k] = k - ns + nr;
2270 } 2246 }
2271 2247
2272 v += ns; 2248 v += ns;
2273 mridx += ns; 2249 mridx += ns;