C-XSC - A C++ Class Library for Extended Scientific Computing  2.5.4
ivector.cpp
00001 /*
00002 **  CXSC is a C++ library for eXtended Scientific Computing (V 2.5.4)
00003 **
00004 **  Copyright (C) 1990-2000 Institut fuer Angewandte Mathematik,
00005 **                          Universitaet Karlsruhe, Germany
00006 **            (C) 2000-2014 Wiss. Rechnen/Softwaretechnologie
00007 **                          Universitaet Wuppertal, Germany   
00008 **
00009 **  This library is free software; you can redistribute it and/or
00010 **  modify it under the terms of the GNU Library General Public
00011 **  License as published by the Free Software Foundation; either
00012 **  version 2 of the License, or (at your option) any later version.
00013 **
00014 **  This library is distributed in the hope that it will be useful,
00015 **  but WITHOUT ANY WARRANTY; without even the implied warranty of
00016 **  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00017 **  Library General Public License for more details.
00018 **
00019 **  You should have received a copy of the GNU Library General Public
00020 **  License along with this library; if not, write to the Free
00021 **  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
00022 */
00023 
00024 /* CVS $Id: ivector.cpp,v 1.26 2014/01/30 17:23:45 cxsc Exp $ */
00025 
00026 #define _CXSC_CPP
00027 
00028 #include "ivector.hpp"
00029 #include "vector.inl"
00030 #include "ivector.inl"
00031 
00032 #include "idotk.inl"
00033 
00034 namespace cxsc {
00035 
00036 int in ( const ivector& x, const ivector& y )          // Inner inclusion for two ivectors
00037 {                                          //---------------------------------
00038   int i = Lb(x), n = Ub(x), ok = 1;
00039 
00040   while (ok && i <= n) { ok = in(x[i],y[i]); i++; }
00041   return ok;
00042 }
00043 
00044 int in ( int x, ivector& y )     // Inner inclusion of an integer in a ivector
00045 {                                //-------------------------------------------
00046   int    i = Lb(y), n = Ub(y), ok = 1;
00047   double d = x;
00048 
00049   while (ok && i <= n) { ok = in(d,y[i]); i++; }
00050   return ok;
00051 }
00052 
00053 ivector Blow ( const ivector& x, real eps )                     // Epsilon deflation
00054 {                                                         //------------------
00055   int     i;
00056   ivector h(Lb(x),Ub(x));
00057 
00058   for (i = Lb(x); i <= Ub(x); i++) h[i] = Blow(x[i],eps);
00059   return h;
00060 }
00061 
00062 int Disjoint ( ivector& a, ivector& b )             // Test for disjointedness
00063 {                                                   //------------------------
00064   int al = Lb(a), au = Ub(a), bl = Lb(b);
00065   int disjointed, i, d;
00066 
00067   d = bl - al;
00068 
00069   i = al;
00070   disjointed = 0;
00071   do {
00072     if (Disjoint(a[i],b[i+d]))
00073       disjointed = 1;
00074     else
00075       i++;
00076   } while ( !(disjointed || i > au) );
00077   return disjointed;
00078 }
00079 
00080 int Zero ( ivector& x )                               // Check for zero vector
00081 {                                                     //----------------------
00082   int i, ok;
00083   for (i = Lb(x), ok = 1; i <= Ub(x) && ok; i++) ok = (x[i] == 0.0);
00084   return ok;
00085 }
00086 
00087 rvector mid ( ivector& v )                              // Vector of midpoints
00088 {                                                       //--------------------
00089   int i;
00090   int l = Lb(v), u = Ub(v);
00091   rvector w(l,u);
00092 
00093   for (i = l; i <= u; i++) w[i] = mid(v[i]);
00094   return w;
00095 }
00096 
00097 real MaxRelDiam ( const ivector& v )                    // Maximum relative diameter
00098 {                                                 //--------------------------
00099   real r;
00100   int  i, l=Lb(v), u=Ub(v);
00101 
00102   r = 0.0;
00103   for (i = l; i <= u; i++)
00104     if (RelDiam(v[i]) > r) r = RelDiam(v[i]);
00105   return r;
00106 }
00107 
00108 real MaxRelDiam ( const ivector_slice& v )                    // Maximum relative diameter
00109 {                                                 //--------------------------
00110   real r;
00111   int  i, l=Lb(v), u=Ub(v);
00112 
00113   r = 0.0;
00114   for (i = l; i <= u; i++)
00115     if (RelDiam(v[i]) > r) r = RelDiam(v[i]);
00116   return r;
00117 }
00118 
00119 //----------------------------------------------------------------------------
00120 // Checks if the diameter of the interval vector 'v' is less or equal to 'n'
00121 // ulps. Ulp is an abbreviation for: units in the last place of the mantissa.
00122 //----------------------------------------------------------------------------
00123 int UlpAcc ( ivector& v, int n )
00124 {
00125   int k, upper;
00126 
00127   for (k = Lb(v), upper = Ub(v); (k < upper) && UlpAcc(v[k],n); k++);
00128   return UlpAcc(v[k],n);
00129 }
00130 
00131 // The 'DoubleSize' functions double the number of rows of a matrix
00132 // or double the length of a vector preserving existing components.
00133 //------------------------------------------------------------------
00134 void DoubleSize ( ivector& x )
00135 {
00136   int n = Lb(x);
00137   Resize(x,n,2*Ub(x)-n+1);
00138 }
00139 
00140 
00141 //accurate dot product definitions
00142 
00143 
00144         void accumulate(idotprecision &dp, const ivector & rv1, const ivector &rv2)
00145 #if(CXSC_INDEX_CHECK)
00146         throw(OP_WITH_WRONG_DIM)
00147 #else
00148         throw()
00149 #endif
00150         { 
00151 #if(CXSC_INDEX_CHECK)
00152                 if(VecLen(rv1)!=VecLen(rv2)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(idotprecision &, const ivector &, ivector &)"));
00153 #endif
00154                 addDot(dp,rv1,rv2);
00155         }
00156 
00157 
00158         void accumulate(idotprecision &dp, const ivector_slice & sl, const ivector &rv)
00159 #if(CXSC_INDEX_CHECK)
00160         throw(OP_WITH_WRONG_DIM)
00161 #else
00162         throw()
00163 #endif
00164         { 
00165 #if(CXSC_INDEX_CHECK)
00166                 if(VecLen(sl)!=VecLen(rv)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(idotprecision &, const ivector_slice &, ivector &)"));
00167 #endif
00168                 addDot(dp,sl,rv);
00169         }
00170 
00171 
00172          void accumulate(idotprecision &dp, const ivector &rv, const ivector_slice &sl)
00173 #if(CXSC_INDEX_CHECK)
00174         throw(OP_WITH_WRONG_DIM)
00175 #else
00176         throw()
00177 #endif
00178         { 
00179 #if(CXSC_INDEX_CHECK)
00180                 if(VecLen(rv)!=VecLen(sl)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(idotprecision &, const ivector &, ivector_slice &)"));
00181 #endif
00182                 addDot(dp,rv,sl);
00183         }
00184 
00185 
00186         void accumulate(idotprecision &dp, const ivector_slice & sl1, const ivector_slice &sl2)
00187 #if(CXSC_INDEX_CHECK)
00188         throw(OP_WITH_WRONG_DIM)
00189 #else
00190         throw()
00191 #endif
00192         { 
00193 #if(CXSC_INDEX_CHECK)
00194                 if(VecLen(sl1)!=VecLen(sl2)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(idotprecision &, const ivector_slice &, ivector_slice &)"));
00195 #endif
00196                 addDot(dp,sl1,sl2);
00197         }
00198 
00199 
00200          void accumulate(cidotprecision &dp, const ivector & rv1, const ivector &rv2)
00201 #if(CXSC_INDEX_CHECK)
00202         throw(OP_WITH_WRONG_DIM)
00203 #else
00204         throw()
00205 #endif
00206         { 
00207 #if(CXSC_INDEX_CHECK)
00208                 if(VecLen(rv1)!=VecLen(rv2)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(cidotprecision &, const ivector &, ivector &)"));
00209 #endif
00210                 idotprecision tmp = Re(dp);
00211                 tmp.set_k(dp.get_k());
00212                 addDot(tmp,rv1,rv2);
00213                 SetRe(dp,tmp);
00214         }
00215 
00216 
00217          void accumulate(cidotprecision &dp, const ivector_slice & sl, const ivector &rv)
00218 #if(CXSC_INDEX_CHECK)
00219         throw(OP_WITH_WRONG_DIM)
00220 #else
00221         throw()
00222 #endif
00223         { 
00224 #if(CXSC_INDEX_CHECK)
00225                 if(VecLen(sl)!=VecLen(rv)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(cidotprecision &, const ivector_slice &, ivector &)"));
00226 #endif
00227                 idotprecision tmp = Re(dp);
00228                 tmp.set_k(dp.get_k());
00229                 addDot(tmp,sl,rv);
00230                 SetRe(dp,tmp);
00231         }
00232 
00233 
00234          void accumulate(cidotprecision &dp, const ivector &rv, const ivector_slice &sl)
00235 #if(CXSC_INDEX_CHECK)
00236         throw(OP_WITH_WRONG_DIM)
00237 #else
00238         throw()
00239 #endif
00240         { 
00241 #if(CXSC_INDEX_CHECK)
00242                 if(VecLen(rv)!=VecLen(sl)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(cidotprecision &, const ivector &, ivector_slice &)"));
00243 #endif
00244                 idotprecision tmp = Re(dp);
00245                 tmp.set_k(dp.get_k());
00246                 addDot(tmp,rv,sl);
00247                 SetRe(dp,tmp);
00248         }
00249 
00250 
00251          void accumulate(cidotprecision &dp, const ivector_slice & sl1, const ivector_slice &sl2)
00252 #if(CXSC_INDEX_CHECK)
00253         throw(OP_WITH_WRONG_DIM)
00254 #else
00255         throw()
00256 #endif
00257         { 
00258 #if(CXSC_INDEX_CHECK)
00259                 if(VecLen(sl1)!=VecLen(sl2)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(cidotprecision &, const ivector_slice &, ivector_slice &)"));
00260 #endif
00261                 idotprecision tmp = Re(dp);
00262                 tmp.set_k(dp.get_k());
00263                 addDot(tmp,sl1,sl2);
00264                 SetRe(dp,tmp);
00265         }
00266 
00267 
00268          void accumulate(idotprecision &dp, const rvector & rv1, const ivector &rv2)
00269 #if(CXSC_INDEX_CHECK)
00270         throw(OP_WITH_WRONG_DIM)
00271 #else
00272         throw()
00273 #endif
00274         { 
00275 #if(CXSC_INDEX_CHECK)
00276                 if(VecLen(rv1)!=VecLen(rv2)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(idotprecision &, const rvector &, ivector &)"));
00277 #endif
00278                 addDot(dp,rv1,rv2);
00279         }
00280 
00281 
00282          void accumulate(idotprecision &dp, const ivector & rv1, const rvector &rv2)
00283 #if(CXSC_INDEX_CHECK)
00284         throw(OP_WITH_WRONG_DIM)
00285 #else
00286         throw()
00287 #endif
00288         { 
00289 #if(CXSC_INDEX_CHECK)
00290                 if(VecLen(rv1)!=VecLen(rv2)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(idotprecision &, const ivector &, rvector &)"));
00291 #endif
00292                 addDot(dp,rv1,rv2);
00293         }
00294 
00295 
00296          void accumulate(idotprecision &dp, const rvector_slice & sl, const ivector &rv)
00297 #if(CXSC_INDEX_CHECK)
00298         throw(OP_WITH_WRONG_DIM)
00299 #else
00300         throw()
00301 #endif
00302         { 
00303 #if(CXSC_INDEX_CHECK)
00304                 if(VecLen(sl)!=VecLen(rv)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(idotprecision &, const rvector_slice &, ivector &)"));
00305 #endif
00306                 addDot(dp,sl,rv);
00307         }
00308 
00309 
00310          void accumulate(idotprecision &dp,const ivector_slice &sl,const rvector &rv)
00311 #if(CXSC_INDEX_CHECK)
00312         throw(OP_WITH_WRONG_DIM)
00313 #else
00314         throw()
00315 #endif
00316         { 
00317 #if(CXSC_INDEX_CHECK)
00318                 if(VecLen(sl)!=VecLen(rv)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(idotprecision &, const ivector_slice &, rvector &)"));
00319 #endif
00320                 addDot(dp,sl,rv);
00321         }
00322 
00323 
00324          void accumulate(idotprecision &dp, const rvector &rv, const ivector_slice &sl)
00325 #if(CXSC_INDEX_CHECK)
00326         throw(OP_WITH_WRONG_DIM)
00327 #else
00328         throw()
00329 #endif
00330         { 
00331 #if(CXSC_INDEX_CHECK)
00332                 if(VecLen(rv)!=VecLen(sl)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(idotprecision &, const rvector &, ivector_slice &)"));
00333 #endif
00334                 addDot(dp,rv,sl);
00335         }
00336 
00337 
00338          void accumulate(idotprecision &dp,const ivector &rv,const rvector_slice &sl)
00339 #if(CXSC_INDEX_CHECK)
00340         throw(OP_WITH_WRONG_DIM)
00341 #else
00342         throw()
00343 #endif
00344         { 
00345 #if(CXSC_INDEX_CHECK)
00346                 if(VecLen(rv)!=VecLen(sl)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(idotprecision &, const ivector &, rvector_slice &)"));
00347 #endif
00348                 addDot(dp,rv,sl);
00349         }
00350 
00351 
00352          void accumulate(idotprecision &dp, const ivector_slice & sl1, const rvector_slice &sl2)
00353 #if(CXSC_INDEX_CHECK)
00354         throw(OP_WITH_WRONG_DIM)
00355 #else
00356         throw()
00357 #endif
00358         { 
00359 #if(CXSC_INDEX_CHECK)
00360                 if(VecLen(sl1)!=VecLen(sl2)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(idotprecision &, const ivector_slice &, rvector_slice &)"));
00361 #endif
00362                 addDot(dp,sl1,sl2);
00363         }
00364 
00365 
00366          void accumulate(idotprecision &dp, const rvector_slice & sl1, const ivector_slice &sl2)
00367 #if(CXSC_INDEX_CHECK)
00368         throw(OP_WITH_WRONG_DIM)
00369 #else
00370         throw()
00371 #endif
00372         { 
00373 #if(CXSC_INDEX_CHECK)
00374                 if(VecLen(sl1)!=VecLen(sl2)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(idotprecision &, const rvector_slice &, ivector_slice &)"));
00375 #endif
00376                 addDot(dp,sl1,sl2);
00377         }
00378 
00379 
00380          void accumulate(cidotprecision &dp, const rvector & rv1, const ivector &rv2)
00381 #if(CXSC_INDEX_CHECK)
00382         throw(OP_WITH_WRONG_DIM)
00383 #else
00384         throw()
00385 #endif
00386         { 
00387 #if(CXSC_INDEX_CHECK)
00388                 if(VecLen(rv1)!=VecLen(rv2)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(cidotprecision &, const rvector &, ivector &)"));
00389 #endif
00390                 idotprecision tmp = Re(dp);
00391                 tmp.set_k(dp.get_k());
00392                 addDot(tmp,rv1,rv2);
00393                 SetRe(dp,tmp);
00394         }
00395 
00396 
00397          void accumulate(cidotprecision &dp, const ivector & rv1, const rvector &rv2)
00398 #if(CXSC_INDEX_CHECK)
00399         throw(OP_WITH_WRONG_DIM)
00400 #else
00401         throw()
00402 #endif
00403         { 
00404 #if(CXSC_INDEX_CHECK)
00405                 if(VecLen(rv1)!=VecLen(rv2)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(cidotprecision &, const ivector &, rvector &)"));
00406 #endif
00407                 idotprecision tmp = Re(dp);
00408                 tmp.set_k(dp.get_k());
00409                 addDot(tmp,rv1,rv2);
00410                 SetRe(dp,tmp);
00411         }
00412 
00413 
00414          void accumulate(cidotprecision &dp, const rvector_slice & sl, const ivector &rv)
00415 #if(CXSC_INDEX_CHECK)
00416         throw(OP_WITH_WRONG_DIM)
00417 #else
00418         throw()
00419 #endif
00420         { 
00421 #if(CXSC_INDEX_CHECK)
00422                 if(VecLen(sl)!=VecLen(rv)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(cidotprecision &, const rvector_slice &, ivector &)"));
00423 #endif
00424                 idotprecision tmp = Re(dp);
00425                 tmp.set_k(dp.get_k());
00426                 addDot(tmp,sl,rv);
00427                 SetRe(dp,tmp);
00428         }
00429 
00430 
00431          void accumulate(cidotprecision &dp,const ivector_slice &sl,const rvector &rv)
00432 #if(CXSC_INDEX_CHECK)
00433         throw(OP_WITH_WRONG_DIM)
00434 #else
00435         throw()
00436 #endif
00437         { 
00438 #if(CXSC_INDEX_CHECK)
00439                 if(VecLen(sl)!=VecLen(rv)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(cidotprecision &, const ivector_slice &, rvector &)"));
00440 #endif
00441                 idotprecision tmp = Re(dp);
00442                 tmp.set_k(dp.get_k());
00443                 addDot(tmp,sl,rv);
00444                 SetRe(dp,tmp);
00445         }
00446 
00447 
00448          void accumulate(cidotprecision &dp, const rvector &rv, const ivector_slice &sl)
00449 #if(CXSC_INDEX_CHECK)
00450         throw(OP_WITH_WRONG_DIM)
00451 #else
00452         throw()
00453 #endif
00454         { 
00455 #if(CXSC_INDEX_CHECK)
00456                 if(VecLen(rv)!=VecLen(sl)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(cidotprecision &, const rvector &, ivector_slice &)"));
00457 #endif
00458                 idotprecision tmp = Re(dp);
00459                 tmp.set_k(dp.get_k());
00460                 addDot(tmp,rv,sl);
00461                 SetRe(dp,tmp);
00462         }
00463 
00464 
00465          void accumulate(cidotprecision &dp,const ivector &rv,const rvector_slice &sl)
00466 #if(CXSC_INDEX_CHECK)
00467         throw(OP_WITH_WRONG_DIM)
00468 #else
00469         throw()
00470 #endif
00471         { 
00472 #if(CXSC_INDEX_CHECK)
00473                 if(VecLen(rv)!=VecLen(sl)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(cidotprecision &, const ivector &, rvector_slice &)"));
00474 #endif
00475                 idotprecision tmp = Re(dp);
00476                 tmp.set_k(dp.get_k());
00477                 addDot(tmp,rv,sl);
00478                 SetRe(dp,tmp);
00479         }
00480 
00481 
00482          void accumulate(cidotprecision &dp, const ivector_slice & sl1, const rvector_slice &sl2)
00483 #if(CXSC_INDEX_CHECK)
00484         throw(OP_WITH_WRONG_DIM)
00485 #else
00486         throw()
00487 #endif
00488         { 
00489 #if(CXSC_INDEX_CHECK)
00490                 if(VecLen(sl1)!=VecLen(sl2)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(cidotprecision &, const ivector_slice &, rvector_slice &)"));
00491 #endif
00492                 idotprecision tmp = Re(dp);
00493                 tmp.set_k(dp.get_k());
00494                 addDot(tmp,sl1,sl2);
00495                 SetRe(dp,tmp);
00496         }
00497 
00498 
00499          void accumulate(cidotprecision &dp, const rvector_slice & sl1, const ivector_slice &sl2)
00500 #if(CXSC_INDEX_CHECK)
00501         throw(OP_WITH_WRONG_DIM)
00502 #else
00503         throw()
00504 #endif
00505         { 
00506 #if(CXSC_INDEX_CHECK)
00507                 if(VecLen(sl1)!=VecLen(sl2)) cxscthrow(OP_WITH_WRONG_DIM("void accumulate(cidotprecision &, const rvector_slice &, ivector_slice &)"));
00508 #endif
00509                 idotprecision tmp = Re(dp);
00510                 tmp.set_k(dp.get_k());
00511                 addDot(tmp,sl1,sl2);
00512                 SetRe(dp,tmp);
00513         }
00514 
00515 
00516         //Summation accumulates
00517         void accumulate(idotprecision &dp, const ivector& v) {
00518                 addSum(Inf(dp),Inf(v));
00519                 addSum(Sup(dp),Sup(v));
00520         }
00521 
00522         void accumulate(cidotprecision &dp, const ivector& v) {
00523                 addSum(InfRe(dp),Inf(v));
00524                 addSum(SupRe(dp),Sup(v));
00525         }
00526 
00527 } // namespace cxsc
00528