print of double and polynom fixed.
[scilab.git] / scilab / modules / ast / src / cpp / types / singlepoly.cpp
1 /*
2 *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
3 *  Copyright (C) 2008-2008 - DIGITEO - Antoine ELIAS
4 *
5 *  This file must be used under the terms of the CeCILL.
6 *  This source file is licensed as described in the file COPYING, which
7 *  you should have received as part of this distribution.  The terms
8 *  are also available at
9 *  http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
10 *
11 */
12 #include <sstream>
13 #include <math.h>
14 #include "singlepoly.hxx"
15 #include "double.hxx"
16 #include "tostring_common.hxx"
17 #include "configvariable.hxx"
18
19 extern "C"
20 {
21 #include "log.h"
22 #include "exp.h"
23 #include "elem_common.h"
24 }
25
26 using namespace std;
27
28 namespace types
29 {
30 SinglePoly::SinglePoly()
31 {
32     double* pdblCoefR = NULL;
33     int piDims[2] = {1, 1};
34     create(piDims, 2, &pdblCoefR, NULL);
35     pdblCoefR[0] = 0;
36 }
37
38 SinglePoly::SinglePoly(double** _pdblCoefR, int _iRank)
39 {
40     int piDims[2] = {_iRank + 1, 1};
41     create(piDims, 2, _pdblCoefR, NULL);
42 }
43
44 SinglePoly::SinglePoly(double** _pdblCoefR, double** _pdblCoefI, int _iRank)
45 {
46     int piDims[2] = {_iRank + 1, 1};
47     create(piDims, 2, _pdblCoefR, _pdblCoefI);
48 }
49
50 SinglePoly::~SinglePoly()
51 {
52     deleteAll();
53 }
54
55 void SinglePoly::deleteAll()
56 {
57     delete[] m_pRealData;
58     m_pRealData = NULL;
59     deleteImg();
60 }
61
62 void SinglePoly::deleteImg()
63 {
64     if (m_pImgData != NULL)
65     {
66         delete[] m_pImgData;
67         m_pImgData = NULL;
68     }
69 }
70
71 double SinglePoly::getNullValue()
72 {
73     return 0;
74 }
75
76 SinglePoly* SinglePoly::createEmpty(int _iDims, int* _piDims, bool _bComplex)
77 {
78     double* pdblData = NULL;
79     SinglePoly* pSP = new SinglePoly(&pdblData, _piDims[0] - 1);
80     pSP->setComplex(_bComplex);
81
82     return pSP;
83 }
84
85 double* SinglePoly::allocData(int _iSize)
86 {
87     double* pDbl = NULL;
88     try
89     {
90         if (_iSize < 0)
91         {
92             m_pRealData = NULL;
93             m_pImgData = NULL;
94             char message[bsiz];
95             sprintf(message, _("Can not allocate negative size (%d).\n"),  _iSize);
96             ast::ScilabError se(message);
97             se.SetErrorNumber(999);
98             throw (se);
99         }
100         else
101         {
102             pDbl = new double[_iSize];
103         }
104     }
105     catch (std::bad_alloc &e)
106     {
107         char message[bsiz];
108         sprintf(message, _("Can not allocate %.2f MB memory.\n"),  (double) (_iSize * sizeof(double)) / 1.e6);
109         ast::ScilabError se(message);
110         se.SetErrorNumber(999);
111         throw (se);
112     }
113
114     return pDbl;
115 }
116
117 double SinglePoly::copyValue(double _dblData)
118 {
119     return _dblData;
120 }
121
122 int SinglePoly::getRank()
123 {
124     return m_iSize - 1;
125 }
126
127 bool SinglePoly::setRank(int _iRank, bool bSave)
128 {
129     double *pR = NULL;
130     double *pI = NULL;
131     if (bSave == false)
132     {
133         if (getRank() != _iRank)
134         {
135             int piDims[2] = {_iRank + 1, 1};
136             if (m_pImgData == false)
137             {
138                 deleteAll();
139                 create(piDims, 2, &pR, NULL);
140             }
141             else
142             {
143                 deleteAll();
144                 create(piDims, 2, &pR, &pI);
145             }
146
147             return true;
148         }
149
150         return true;
151     }
152     else
153     {
154         double* pdblOldReal = m_pRealData;
155         double* pdblOldImg  = m_pImgData;
156         int iMinSize = Min(m_iSize, _iRank + 1);
157         int piDims[2] = {_iRank + 1, 1};
158
159         if (m_pImgData == false)
160         {
161             create(piDims, 2, &pR, NULL);
162         }
163         else
164         {
165             create(piDims, 2, &pR, &pI);
166             memcpy(m_pImgData, pdblOldImg, iMinSize * sizeof(double));
167         }
168
169         memcpy(m_pRealData, pdblOldReal, iMinSize * sizeof(double));
170
171         if (pdblOldImg)
172         {
173             delete[] pdblOldImg;
174             pdblOldImg = NULL;
175         }
176
177         delete[] pdblOldReal;
178         pdblOldReal = NULL;
179
180         return true;
181     }
182
183     return false;
184 }
185
186 bool SinglePoly::setZeros()
187 {
188     if (m_pRealData != NULL)
189     {
190         memset(m_pRealData, 0x00, m_iSize * sizeof(double));
191     }
192     else
193     {
194         return false;
195     }
196
197     if (isComplex() == true)
198     {
199         if (m_pImgData != NULL)
200         {
201             memset(m_pImgData, 0x00, m_iSize * sizeof(double));
202         }
203         else
204         {
205             return false;
206         }
207     }
208
209     return true;
210 }
211
212 bool SinglePoly::setCoef(Double* _pdblCoefR)
213 {
214     if (m_pRealData == NULL || _pdblCoefR == NULL)
215     {
216         return false;
217     }
218
219     double *pInR    = _pdblCoefR->getReal();
220     double *pInI    = _pdblCoefR->getImg();
221
222     return setCoef(pInR, pInI);
223 }
224
225 bool SinglePoly::setCoef(double* _pdblCoefR, double* _pdblCoefI)
226 {
227     if (_pdblCoefI != NULL && isComplex() == false)
228     {
229         setComplex(true);
230     }
231
232     if (_pdblCoefR != NULL)
233     {
234         memcpy(m_pRealData, _pdblCoefR, m_iSize * sizeof(double));
235     }
236
237     if (_pdblCoefI != NULL)
238     {
239         memcpy(m_pImgData, _pdblCoefI, m_iSize * sizeof(double));
240     }
241
242     return true;
243 }
244
245 void SinglePoly::whoAmI()
246 {
247     std::cout << "types::SinglePoly";
248 }
249
250 bool SinglePoly::evaluate(double _dblInR, double _dblInI, double *_pdblOutR, double *_pdblOutI)
251 {
252     *_pdblOutR = 0;
253     *_pdblOutI = 0;
254     if (m_iSize == 0)
255     {
256         return true;
257     }
258
259     for (int i = 0 ; i < m_iSize ; i++)
260     {
261         //real part
262         *_pdblOutR += m_pRealData[i] * pow(_dblInR, i);
263         //only if variable is complex
264         if (isComplex())
265         {
266             *_pdblOutR -= m_pImgData[i] * pow(_dblInI, i);
267             //img part
268             *_pdblOutI += m_pRealData[i] * pow(_dblInR, i);
269         }
270         *_pdblOutI += m_pRealData[i] * pow(_dblInI, i);
271     }
272
273     return true;
274 }
275
276 void SinglePoly::updateRank(void)
277 {
278     double dblEps = getRelativeMachinePrecision();
279     int iNewRank = getRank();
280     if (m_pImgData)
281     {
282         for (int i = getRank(); i > 0 ; i--)
283         {
284             if (fabs(m_pRealData[i]) <= dblEps && abs(m_pImgData[i]) <= dblEps)
285             {
286                 iNewRank--;
287             }
288             else
289             {
290                 break;
291             }
292         }
293     }
294     else
295     {
296         for (int i = getRank(); i > 0 ; i--)
297         {
298             if (fabs(m_pRealData[i]) <= dblEps)
299             {
300                 iNewRank--;
301             }
302             else
303             {
304                 break;
305             }
306         }
307     }
308
309     if (iNewRank < getRank())
310     {
311         setRank(iNewRank, true);
312     }
313 }
314
315 bool SinglePoly::toString(std::wostringstream& ostr)
316 {
317     ostr << L"FIXME : implement SinglePoly::toString" << std::endl;
318     return true;
319 }
320
321 void SinglePoly::toStringReal(wstring _szVar, list<wstring>* _pListExp , list<wstring>* _pListCoef)
322 {
323     toStringInternal(m_pRealData, _szVar, _pListExp, _pListCoef);
324 }
325
326 void SinglePoly::toStringImg(wstring _szVar, list<wstring>* _pListExp , list<wstring>* _pListCoef)
327 {
328     if (isComplex() == false)
329     {
330         _pListExp->clear();
331         _pListCoef->clear();
332         return;
333     }
334
335     toStringInternal(m_pImgData, _szVar, _pListExp, _pListCoef);
336 }
337
338 bool SinglePoly::subMatrixToString(wostringstream& ostr, int* _piDims, int _iDims)
339 {
340     return false;
341 }
342
343 void SinglePoly::toStringInternal(double *_pdblVal, wstring _szVar, list<wstring>* _pListExp , list<wstring>* _pListCoef)
344 {
345     int iPrecision = ConfigVariable::getFormatSize();
346     int iLineLen = ConfigVariable::getConsoleWidth();
347
348     wostringstream ostemp;
349     wostringstream ostemp2;
350
351     ostemp << L" ";
352     ostemp2 << L" ";
353
354     int iLen = 0;
355     int iLastFlush = 2;
356     for (int i = 0 ; i < m_iSize ; i++)
357     {
358         if (isRealZero(_pdblVal[i]) == false)
359         {
360             DoubleFormat df;
361             getDoubleFormat(_pdblVal[i], &df);
362
363             if (iLen + df.iWidth + df.iSignLen >= iLineLen - 1)
364             {
365                 iLastFlush = i;
366                 _pListExp->push_back(ostemp2.str());
367                 ostemp2.str(L""); //reset stream
368                 addSpaces(&ostemp2, 11); //take from scilab ... why not ...
369
370                 _pListCoef->push_back(ostemp.str());
371                 ostemp.str(L""); //reset stream
372                 addSpaces(&ostemp, 11); //take from scilab ... why not ...
373             }
374
375             bool bFirst = ostemp.str().size() == 1;
376
377             // In scientific notation case bExp == true, so we have to print point (2.000D+10s)
378             // In other case don't print point (2s)
379             df.bPrintPoint = df.bExp;
380             df.bPrintPlusSign = bFirst == false;
381             df.bPrintOne = i == 0;
382             addDoubleValue(&ostemp, _pdblVal[i], &df);
383
384             if (i == 0)
385             {
386                 iLen = static_cast<int>(ostemp.str().size());
387             }
388             else if (i == 1)
389             {
390                 // add polynom name
391                 ostemp << _szVar;
392                 iLen = static_cast<int>(ostemp.str().size());
393             }
394             else
395             {
396                 // add polynom name and exponent
397                 ostemp << _szVar;
398                 iLen = static_cast<int>(ostemp.str().size());
399                 addSpaces(&ostemp2, iLen - static_cast<int>(ostemp2.str().size()));
400                 ostemp2 << i;
401                 int iSize = static_cast<int>(ostemp2.str().size()) - iLen;
402                 addSpaces(&ostemp, iSize);
403             }
404         }
405     }
406
407     if (iLastFlush != 0)
408     {
409         if (ostemp.str() == L" ")
410         {
411             ostemp << L"  0";
412             addSpaces(&ostemp2, static_cast<int>(ostemp.str().size()));
413         }
414
415         if (ostemp2.str() == L" ")
416         {
417             addSpaces(&ostemp2, static_cast<int>(ostemp.str().size()));
418         }
419
420         _pListExp->push_back(ostemp2.str());
421         _pListCoef->push_back(ostemp.str());
422     }
423
424     return;
425 }
426
427 bool SinglePoly::operator==(const InternalType& it)
428 {
429     if (const_cast<InternalType &>(it).isSinglePoly() == false)
430     {
431         return false;
432     }
433
434     SinglePoly* pP = const_cast<InternalType &>(it).getAs<SinglePoly>();
435
436     if (getRank() != pP->getRank())
437     {
438         return false;
439     }
440
441     double *pdblReal = pP->get();
442     for (int i = 0 ; i < getSize() ; i++)
443     {
444         if (m_pRealData[i] != pdblReal[i])
445         {
446             return false;
447         }
448     }
449
450     //both complex
451     if (isComplex() && pP->isComplex())
452     {
453         double *pdblImg = pP->getImg();
454         for (int i = 0 ; i < m_iSize ; i++)
455         {
456             if (m_pImgData[i] != pdblImg[i])
457             {
458                 return false;
459             }
460         }
461     }
462     //pdbl complex check all img values == 0
463     else if (pP->isComplex())
464     {
465         double *pdblImg = pP->getImg();
466         for (int i = 0 ; i < m_iSize ; i++)
467         {
468             if (pdblImg[i])
469             {
470                 return false;
471             }
472         }
473     }
474     //complex check all img values == 0
475     else if (isComplex())
476     {
477         for (int i = 0 ; i < m_iSize ; i++)
478         {
479             if (m_pImgData[i])
480             {
481                 return false;
482             }
483         }
484     }
485
486     return true;
487 }
488
489 bool SinglePoly::operator!=(const InternalType& it)
490 {
491     return !(*this == it);
492 }
493
494 SinglePoly* SinglePoly::clone()
495 {
496     SinglePoly* pPoly = NULL;
497     double *pR = NULL;
498     if (isComplex())
499     {
500         double *pI = NULL;
501         pPoly = new SinglePoly(&pR, &pI, getRank());
502         pPoly->setCoef(m_pRealData, m_pImgData);
503     }
504     else
505     {
506         pPoly = new SinglePoly(&pR, getRank());
507         pPoly->setCoef(m_pRealData, NULL);
508     }
509     return pPoly;
510 }
511
512 SinglePoly* SinglePoly::conjugate()
513 {
514     SinglePoly* pPoly = NULL;
515     if (isComplex())
516     {
517         double *pR = NULL;
518         double *pI = NULL;
519         SinglePoly * pPoly = new SinglePoly(&pR, &pI, getRank());
520
521         Transposition::conjugate(m_iSize, m_pRealData, pR, m_pImgData, pI);
522
523         return pPoly;
524     }
525     else
526     {
527         return clone();
528     }
529 }
530
531 SinglePoly* operator*(const SinglePoly& _lhs, const SinglePoly& _rhs)
532 {
533     SinglePoly& lhs = const_cast<SinglePoly &>(_lhs);
534     SinglePoly& rhs = const_cast<SinglePoly &>(_rhs);
535     SinglePoly* pOut = NULL;
536
537     bool isComplexL = lhs.isComplex();
538     bool isComplexR = rhs.isComplex();
539     bool isComplexOut = isComplexL || isComplexR;
540
541     int iRankL = lhs.getRank();
542     int iRankR = rhs.getRank();
543     int iRankOut = lhs.getRank() + rhs.getRank();
544
545     double* pdblOutR = NULL;
546     double* pdblOutI = NULL;
547     double* pdblLR = lhs.get();
548     double* pdblLI = lhs.getImg();
549     double* pdblRR = rhs.get();
550     double* pdblRI = rhs.getImg();
551
552     if (isComplexOut)
553     {
554         pOut = new SinglePoly(&pdblOutR, &pdblOutI, iRankOut);
555         memset(pdblOutR, 0x00, sizeof(double) * (iRankOut + 1));
556         memset(pdblOutI, 0x00, sizeof(double) * (iRankOut + 1));
557     }
558     else
559     {
560         pOut = new SinglePoly(&pdblOutR, iRankOut);
561         memset(pdblOutR, 0x00, sizeof(double) * (iRankOut + 1));
562     }
563
564     if (isComplexL)
565     {
566         if (isComplexR)
567         {
568             for (int i = 0 ; i < iRankL + 1; ++i)
569             {
570                 for (int j = 0 ; j < iRankR + 1; ++j)
571                 {
572                     pdblOutR[i + j]  += pdblLR[i] * pdblRR[j] - pdblLI[i] * pdblRI[j];
573                     pdblOutI[i + j]  += pdblLI[i] * pdblRR[j] + pdblLR[i] * pdblRI[j];
574                 }
575             }
576         }
577         else
578         {
579             for (int i = 0 ; i < iRankL + 1 ; ++i)
580             {
581                 for (int j = 0 ; j < iRankR + 1 ; ++j)
582                 {
583                     pdblOutR[i + j]  += pdblLR[i] * pdblRR[j];
584                     pdblOutI[i + j]  += pdblLI[i] * pdblRR[j];
585                 }
586             }
587         }
588     }
589     else
590     {
591         if (isComplexR)
592         {
593             for (int i = 0 ; i < iRankL + 1 ; ++i)
594             {
595                 for (int j = 0 ; j < iRankR + 1 ; ++j)
596                 {
597                     pdblOutR[i + j]  += pdblLR[i] * pdblRR[j];
598                     pdblOutI[i + j]  += pdblLR[i] * pdblRI[j];
599                 }
600             }
601         }
602         else
603         {
604             for (int i = 0 ; i < iRankL + 1 ; ++i)
605             {
606                 for (int j = 0 ; j < iRankR + 1 ; ++j)
607                 {
608                     pdblOutR[i + j]  += pdblLR[i] * pdblRR[j];
609                 }
610             }
611         }
612     }
613
614     return pOut;
615 }
616 }
617
618
619