~ubuntu-branches/ubuntu/wily/afnix/wily

« back to all changes in this revision

Viewing changes to src/mod/mth/shl/Rmi.cpp

  • Committer: Bazaar Package Importer
  • Author(s): Anibal Monsalve Salazar
  • Date: 2011-03-16 21:31:18 UTC
  • mfrom: (1.1.2 upstream)
  • Revision ID: james.westby@ubuntu.com-20110316213118-gk4k3ez3e5d2huna
Tags: 2.0.0-1
* QA upload.
* New upstream release
* Debian source format is 3.0 (quilt)
* Fix debhelper-but-no-misc-depends
* Fix ancient-standards-version
* Fix package-contains-linda-override
* debhelper compatibility is 7
* Fix dh-clean-k-is-deprecated

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
// ---------------------------------------------------------------------------
 
2
// - Rmi.cpp                                                                 -
 
3
// - afnix:mth module - real matrix interface implementation                 -
 
4
// ---------------------------------------------------------------------------
 
5
// - This program is free software;  you can redistribute it  and/or  modify -
 
6
// - it provided that this copyright notice is kept intact.                  -
 
7
// -                                                                         -
 
8
// - This program  is  distributed in  the hope  that it will be useful, but -
 
9
// - without  any  warranty;  without  even   the   implied    warranty   of -
 
10
// - merchantability or fitness for a particular purpose.  In no event shall -
 
11
// - the copyright holder be liable for any  direct, indirect, incidental or -
 
12
// - special damages arising in any way out of the use of this software.     -
 
13
// ---------------------------------------------------------------------------
 
14
// - copyright (c) 1999-2011 amaury darsch                                   -
 
15
// ---------------------------------------------------------------------------
 
16
 
 
17
#include "Rmi.hpp"
 
18
#include "Real.hpp"
 
19
#include "Math.hpp"
 
20
#include "Vector.hpp"
 
21
#include "Algebra.hpp"
 
22
#include "Boolean.hpp"
 
23
#include "QuarkZone.hpp"
 
24
#include "Exception.hpp"
 
25
 
 
26
namespace afnix {
 
27
 
 
28
  // -------------------------------------------------------------------------
 
29
  // - class section                                                         -
 
30
  // -------------------------------------------------------------------------
 
31
 
 
32
  // create a default matrix
 
33
 
 
34
  Rmi::Rmi (void) {
 
35
    d_rsiz = 0;
 
36
    d_csiz = 0;
 
37
  }
 
38
 
 
39
  // create a square matrix by size
 
40
 
 
41
  Rmi::Rmi (const t_long size) {
 
42
    // check the size
 
43
    if (size < 0) {
 
44
      throw Exception ("size-error", "invalid real matrix size");
 
45
    }
 
46
    d_rsiz = size;
 
47
    d_csiz = size;
 
48
  }
 
49
 
 
50
  // create a matrix by size
 
51
 
 
52
  Rmi::Rmi (const t_long rsiz, const t_long csiz) {
 
53
    // check the size
 
54
    if ((rsiz < 0) || (csiz < 0)) {
 
55
      throw Exception ("size-error", "invalid real matrix size");
 
56
    }
 
57
    d_rsiz = rsiz;
 
58
    d_csiz = csiz;
 
59
  }
 
60
 
 
61
  // clear this matrix
 
62
 
 
63
  void Rmi::clear (void) {
 
64
    wrlock ();
 
65
    try {
 
66
      for (t_long i = 0; i < d_rsiz; i++) {
 
67
        for (t_long j = 0; j < d_csiz; j++) {
 
68
          set (i, j, 0.0);
 
69
        }
 
70
      }
 
71
      unlock ();
 
72
    } catch (...) {
 
73
      unlock ();
 
74
      throw;
 
75
    }
 
76
  }
 
77
 
 
78
  // get the matrix row size
 
79
 
 
80
  t_long Rmi::getrsiz (void) const {
 
81
    rdlock ();
 
82
    try {
 
83
      t_long result = d_rsiz;
 
84
      unlock ();
 
85
      return result;
 
86
    } catch (...) {
 
87
      unlock ();
 
88
      throw;
 
89
    }
 
90
  }
 
91
 
 
92
  // get the matrix column size
 
93
 
 
94
  t_long Rmi::getcsiz (void) const {
 
95
    rdlock ();
 
96
    try {
 
97
      t_long result = d_csiz;
 
98
      unlock ();
 
99
      return result;
 
100
    } catch (...) {
 
101
      unlock ();
 
102
      throw;
 
103
    }
 
104
  }
 
105
 
 
106
  // return true if the matrix is square
 
107
 
 
108
  bool Rmi::issquare (void) const {
 
109
    rdlock ();
 
110
    try {
 
111
      bool result = (d_rsiz == d_csiz);
 
112
      unlock ();
 
113
      return result;
 
114
    } catch (...) {
 
115
      unlock ();
 
116
      throw;
 
117
    }
 
118
  }
 
119
 
 
120
  // compare two matrices by precision
 
121
 
 
122
  bool Rmi::cmp (const Rmi& mx) const {
 
123
    rdlock ();
 
124
    mx.rdlock ();
 
125
    try {
 
126
      // check size first
 
127
      if ((d_rsiz != mx.d_rsiz) || (d_csiz != mx.d_csiz)) {
 
128
        throw Exception ("matrix-error",
 
129
                         "incompatible matrix size with compare");
 
130
      }
 
131
      // initialize result
 
132
      bool result = true;
 
133
      // loop in vector
 
134
      for (t_long i = 0; i < d_rsiz; i++) {
 
135
        for (t_long j = 0; j < d_csiz; j++) {
 
136
          t_real tij = get (i, j);
 
137
          t_real mij = mx.get (i, j);
 
138
          t_real dij = (mij < tij) ? mij - tij : tij - mij;
 
139
          if (dij > Real::d_aeps) {
 
140
            result = false;
 
141
            break;
 
142
          }
 
143
        }
 
144
      }
 
145
      // unlock and return
 
146
      unlock ();
 
147
      mx.unlock ();
 
148
      return result;
 
149
    } catch (...) {
 
150
      unlock ();
 
151
      mx.unlock ();
 
152
      throw;
 
153
    }
 
154
  }
 
155
 
 
156
  // compute the matrix norm
 
157
 
 
158
  t_real Rmi::norm (void) const {
 
159
    rdlock ();
 
160
    try {
 
161
      t_real sum = 0.0;
 
162
      for (t_long i = 0; i < d_rsiz; i++) {
 
163
        for (t_long j = 0; j < d_csiz; j++) {
 
164
          t_real mij = get (i, j);
 
165
          sum += (mij * mij);
 
166
        }
 
167
      }
 
168
      t_real result = Math::sqrt (sum);
 
169
      unlock ();
 
170
      return result;
 
171
    } catch (...) {
 
172
      unlock ();
 
173
      throw;
 
174
    }
 
175
  }
 
176
 
 
177
  // positive multiply a matrix with a vector
 
178
 
 
179
  Rvi& Rmi::pmul (Rvi& r, const Rvi& x) const {
 
180
    rdlock ();
 
181
    try {
 
182
      Algebra::mul (r, *this, x, 1.0);
 
183
      unlock ();
 
184
      return r;
 
185
    } catch (...) {
 
186
      unlock ();
 
187
      throw;
 
188
    }
 
189
  }
 
190
 
 
191
  // negative multiply a matrix with a vector
 
192
 
 
193
  Rvi& Rmi::nmul (Rvi& r, const Rvi& x) const {
 
194
    rdlock ();
 
195
    try {
 
196
      Algebra::mul (r, *this, x, -1.0);
 
197
      unlock ();
 
198
      return r;
 
199
    } catch (...) {
 
200
      unlock ();
 
201
      throw;
 
202
    }
 
203
  }
 
204
 
 
205
  // multiply a matrix with a vector and a scaling factor
 
206
 
 
207
  Rvi& Rmi::smul (Rvi& r, const Rvi& x, const t_real s) const {
 
208
    rdlock ();
 
209
    try {
 
210
      Algebra::mul (r, *this, x, s);
 
211
      unlock ();
 
212
      return r;
 
213
    } catch (...) {
 
214
      unlock ();
 
215
      throw;
 
216
    }
 
217
  }
 
218
 
 
219
  // -------------------------------------------------------------------------
 
220
  // - object section                                                        -
 
221
  // -------------------------------------------------------------------------
 
222
 
 
223
  // the quark zone
 
224
  static const long QUARK_ZONE_LENGTH = 9;
 
225
  static QuarkZone  zone (QUARK_ZONE_LENGTH);
 
226
 
 
227
  // the rmi supported quarks
 
228
  static const long QUARK_QEQ     = zone.intern ("?=");
 
229
  static const long QUARK_SET     = zone.intern ("set");
 
230
  static const long QUARK_GET     = zone.intern ("get");
 
231
  static const long QUARK_NORM    = zone.intern ("norm");
 
232
  static const long QUARK_CLEAR   = zone.intern ("clear");
 
233
  static const long QUARK_TORSO   = zone.intern ("to-row-sparse");
 
234
  static const long QUARK_GETCFM  = zone.intern ("get-cofactor");
 
235
  static const long QUARK_GETRSIZ = zone.intern ("get-row-size");
 
236
  static const long QUARK_GETCSIZ = zone.intern ("get-col-size");
 
237
 
 
238
  // return true if the given quark is defined
 
239
 
 
240
  bool Rmi::isquark (const long quark, const bool hflg) const {
 
241
    rdlock ();
 
242
    if (zone.exists (quark) == true){
 
243
      unlock ();
 
244
      return true;
 
245
    }
 
246
    bool result = hflg ? Object::isquark (quark, hflg) : false;
 
247
    unlock ();
 
248
    return result;
 
249
  }
 
250
 
 
251
  // apply this object with a set of arguments and a quark
 
252
  
 
253
  Object* Rmi::apply (Runnable* robj, Nameset* nset, const long quark,
 
254
                      Vector* argv) {
 
255
    // get the number of arguments
 
256
    long argc = (argv == nilp) ? 0 : argv->length ();
 
257
 
 
258
    // dispatch 0 argument
 
259
    if (argc == 0) {
 
260
      if (quark == QUARK_NORM)    return new Real (norm ());
 
261
      if (quark == QUARK_GETRSIZ) return new Integer (getrsiz ());
 
262
      if (quark == QUARK_GETCSIZ) return new Integer (getcsiz ());
 
263
      if (quark == QUARK_CLEAR) {
 
264
        clear ();
 
265
        return nilp;
 
266
      }
 
267
    }
 
268
    
 
269
    // dispatch 1 argument
 
270
    if (argc == 1) {
 
271
      if (quark == QUARK_QEQ) {
 
272
        Object* obj = argv->get (0);
 
273
        Rmi*    rmo = dynamic_cast <Rmi*> (obj);
 
274
        if (rmo == nilp) {
 
275
          throw Exception ("type-error", "invalid object for compare",
 
276
                           Object::repr (obj));
 
277
        }
 
278
        return new Boolean (cmp (*rmo));
 
279
      }
 
280
      if (quark == QUARK_TORSO) {
 
281
        t_long row = argv->getlong (0);
 
282
        return torso (row);
 
283
      }
 
284
    }
 
285
    // dispatch 2 argument
 
286
    if (argc == 2) {
 
287
      if (quark == QUARK_GET) {
 
288
        t_long lpos = argv->getlong (0);
 
289
        t_long cpos = argv->getlong (1);
 
290
        return new Real (get (lpos, cpos));
 
291
      }
 
292
      if (quark == QUARK_GETCFM) {
 
293
        t_long row = argv->getlong (0);
 
294
        t_long col = argv->getlong (1);
 
295
        return getcfm (row, col);
 
296
      }
 
297
    }
 
298
    // dispatch 3 arguments
 
299
    if (argc == 3) {
 
300
      if (quark == QUARK_SET) {
 
301
        t_long lpos = argv->getlong (0);
 
302
        t_long cpos = argv->getlong (1);
 
303
        t_real  val = argv->getreal (2);
 
304
        set (lpos, cpos, val);
 
305
        return nilp;
 
306
      }
 
307
    }
 
308
    // call the object
 
309
    return Object::apply (robj, nset, quark, argv);
 
310
  }
 
311
}
 
312