Coverity 16571
[u/mrichter/AliRoot.git] / STEER / AliVectorSparse.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /* $Id: AliRun.cxx 30859 2009-02-02 16:24:37Z fca $ */
17
18 #include <string.h>
19 #include "AliVectorSparse.h"
20 #include <TString.h>
21
22 /**********************************************************************************************/
23 /* Sparse vector class, used as row of the AliMatrixSparse class                              */
24 /*                                                                                            */ 
25 /* Author: ruben.shahoyan@cern.ch                                                             */
26 /*                                                                                            */ 
27 /**********************************************************************************************/
28
29 ClassImp(AliVectorSparse)
30
31 //___________________________________________________________
32 AliVectorSparse::AliVectorSparse()
33   : fNElems(0),fIndex(0),fElems(0) {}
34
35 //___________________________________________________________
36 AliVectorSparse::AliVectorSparse(const AliVectorSparse& src)
37   : TObject(src),fNElems(src.fNElems),fIndex(0),fElems(0)
38 {
39   fIndex = new UShort_t[fNElems];
40   fElems = new Double_t[fNElems];
41   memcpy(fIndex,src.fIndex,fNElems*sizeof(UShort_t));
42   memcpy(fElems,src.fElems,fNElems*sizeof(Double_t));
43 }
44
45 //___________________________________________________________
46 void AliVectorSparse::Clear(Option_t*)
47 {
48   delete[] fIndex; fIndex = 0;
49   delete[] fElems; fElems = 0;
50   fNElems = 0;
51 }
52
53 //___________________________________________________________
54 AliVectorSparse& AliVectorSparse::operator=(const AliVectorSparse& src)
55 {
56   if (&src==this) return *this;
57   Clear();
58   TObject::operator=(src);
59   fNElems = src.fNElems;
60   fIndex = new UShort_t[fNElems];
61   fElems = new Double_t[fNElems];
62   memcpy(fIndex,src.fIndex,fNElems*sizeof(UShort_t));
63   memcpy(fElems,src.fElems,fNElems*sizeof(Double_t));
64   //
65   return *this;
66 }
67
68 //___________________________________________________________
69 Double_t AliVectorSparse::FindIndex(Int_t ind) const
70 {
71   // return an element with given index
72   //printf("V: findindex\n");
73   int first = 0;
74   int last = fNElems-1;
75   while (first<=last) {
76     int mid = (first+last)>>1;
77     if (ind>fIndex[mid]) first = mid+1;
78     else if (ind<fIndex[mid]) last = mid-1;
79     else return fElems[mid];
80   }
81   return 0.0;
82 }
83
84 //___________________________________________________________
85 void AliVectorSparse::SetToZero(Int_t ind)
86 {
87   // set element to 0 if it was already defined
88   int first = 0;
89   int last = fNElems-1;
90   while (first<=last) {
91     int mid = (first+last)>>1;
92     if (ind>fIndex[mid]) first = mid+1;
93     else if (ind<fIndex[mid]) last = mid-1;
94     else {fElems[mid] = 0.; return;}
95   }
96 }
97
98 //___________________________________________________________
99 Double_t& AliVectorSparse::FindIndexAdd(Int_t ind)
100 {
101   // increment an element with given index
102   //printf("V: findindexAdd\n");
103   int first = 0;
104   int last = fNElems-1;
105   while (first<=last) {
106     int mid = (first+last)>>1;
107     if (ind>fIndex[mid]) first = mid+1;
108     else if (ind<fIndex[mid]) last = mid-1;
109     else return fElems[mid];
110   }
111   // need to insert a new element
112   UShort_t *arrI = new UShort_t[fNElems+1];
113   memcpy(arrI,fIndex,first*sizeof(UShort_t));
114   arrI[first] = ind;
115   memcpy(arrI+first+1,fIndex+first,(fNElems-first)*sizeof(UShort_t));
116   delete[] fIndex;
117   fIndex = arrI;
118   //
119   Double_t   *arrE = new Double_t[fNElems+1];
120   memcpy(arrE,fElems,first*sizeof(Double_t));
121   arrE[first] = 0;
122   memcpy(arrE+first+1,fElems+first,(fNElems-first)*sizeof(Double_t));
123   delete[] fElems;
124   fElems = arrE;
125   //
126   fNElems++;
127   return fElems[first];
128   //
129 }
130
131 //__________________________________________________________
132 void AliVectorSparse::ReSize(Int_t sz,Bool_t copy)
133 {
134   if (sz<1) {Clear(); return;}
135     // need to insert a new element
136   UShort_t *arrI = new UShort_t[sz];
137   Double_t *arrE = new Double_t[sz];
138   memset(arrI,0,sz*sizeof(UShort_t));
139   memset(arrE,0,sz*sizeof(Double_t));
140   //
141   if (copy && fIndex) {
142     int cpsz = TMath::Min(fNElems,sz);
143     memcpy(arrI,fIndex,cpsz*sizeof(UShort_t));
144     memcpy(arrE,fElems,cpsz*sizeof(Double_t));
145   }
146   delete[] fIndex;
147   delete[] fElems;
148   fIndex = arrI;
149   fElems = arrE;
150   fNElems = sz;
151   //
152 }
153
154 //__________________________________________________________
155 void AliVectorSparse::SortIndices(Bool_t valuesToo)
156 {
157   // sort indices in increasing order. Used to fix the row after ILUk decomposition
158   for (int i=fNElems;i--;) for (int j=i;j--;) if (fIndex[i]<fIndex[j]) { //swap
159         UShort_t tmpI = fIndex[i]; fIndex[i] = fIndex[j]; fIndex[j]=tmpI;
160         if (valuesToo) {Double_t tmpV = fElems[i];fElems[i]=fElems[j];fElems[j]=tmpV;}
161       }
162 }
163
164 //__________________________________________________________
165 void AliVectorSparse::Print(Option_t* opt)  const
166 {
167   TString sopt = opt; sopt.ToLower();
168   int ndig = sopt.Atoi();
169   if (ndig<=1) ndig = 2;
170   sopt = "%2d:%+.";
171   sopt += ndig;
172   sopt += "e |";
173   printf("|");
174   for (int i=0;i<fNElems;i++) printf(sopt.Data(),fIndex[i],fElems[i]);
175   printf("\n");
176 }
177
178 //___________________________________________________________
179 void AliVectorSparse::Add(Double_t *valc,Int_t *indc,Int_t n)
180 {
181   // add indiced array to row. Indices must be in increasing order
182   int indx;
183   int nadd = 0;
184   //
185   int last = fNElems-1;
186   int mid = 0;
187   for (int i=n;i--;) {
188     // if the element with this index is already defined, just add the value
189     int first = 0;
190     Bool_t toAdd = kTRUE;
191     indx = indc[i];
192     while (first<=last) {
193       mid = (first+last)>>1;
194       if (indx>fIndex[mid]) first = mid+1;
195       else if (indx<fIndex[mid]) last = mid-1;
196       else {
197         fElems[mid] += valc[i];
198         indc[i] = -1;
199         toAdd = kFALSE;
200         last = mid-1;   // profit from the indices being ordered
201         break;
202       }
203     }
204     if (toAdd) nadd++;
205   }
206   //
207   if (nadd<1) return; // nothing to do anymore
208   // 
209   // need to expand the row
210   UShort_t *arrI = new UShort_t[fNElems+nadd];
211   Double_t *arrE = new Double_t[fNElems+nadd];
212   // copy old elems embedding the new ones
213   int inew=0,iold=0;
214   for (int i=0;i<n;i++) {  
215     if ( (indx=indc[i])<0) continue;
216     while (iold<fNElems && fIndex[iold]<indx) {
217       arrI[inew]   = fIndex[iold];
218       arrE[inew++] = fElems[iold++];
219     }
220     arrI[inew] = indx;
221     arrE[inew++] = valc[i];
222   }
223   // copy the rest
224   while (iold<fNElems) {
225     arrI[inew]   = fIndex[iold];
226     arrE[inew++] = fElems[iold++];
227   }
228   //
229   delete[] fIndex;
230   delete[] fElems;
231   fIndex = arrI;
232   fElems = arrE;
233   //
234   fNElems += nadd;
235   //
236 }
237