6e6901889305ca3cfab2ce8b4d43a419112da308
[scilab.git] / scilab / modules / ast / src / cpp / operations / types_ldivide.cpp
1 /*
2  *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
3  *  Copyright (C) 2013 - Scilab Enterprises - Cedric Delamarre
4  *
5  * Copyright (C) 2012 - 2016 - Scilab Enterprises
6  *
7  * This file is hereby licensed under the terms of the GNU GPL v2.0,
8  * pursuant to article 5.3.4 of the CeCILL v.2.1.
9  * This file was originally licensed under the terms of the CeCILL v2.1,
10  * and continues to be available under such terms.
11  * For more information, see the COPYING file which you should have received
12  * along with this program.
13  *
14  */
15
16 #include "types_ldivide.hxx"
17 #include "types_divide.hxx"
18 #include "types_finite.hxx"
19
20 extern "C"
21 {
22 #include "matrix_left_division.h"
23 #include "sciprint.h"
24 #include "localization.h"
25 #include "charEncoding.h"
26 }
27
28 using namespace types;
29
30 InternalType *GenericLDivide(InternalType *_pLeftOperand, InternalType *_pRightOperand)
31 {
32     InternalType *pResult       = NULL;
33     GenericType::ScilabType TypeL = _pLeftOperand->getType();
34     GenericType::ScilabType TypeR = _pRightOperand->getType();
35
36     int iResult = 0;
37
38     if (_pLeftOperand->isDouble() && _pLeftOperand->getAs<Double>()->isEmpty())
39     {
40         return Double::Empty();
41     }
42
43     if (_pRightOperand->isDouble() && _pRightOperand->getAs<Double>()->isEmpty())
44     {
45         return Double::Empty();
46     }
47
48     /*
49     ** DOUBLE \ DOUBLE
50     */
51     if (TypeL == GenericType::ScilabDouble && TypeR == GenericType::ScilabDouble)
52     {
53         Double *pL  = _pLeftOperand->getAs<Double>();
54         Double *pR  = _pRightOperand->getAs<Double>();
55
56         iResult = LDivideDoubleByDouble(pL, pR, (Double**)&pResult);
57     }
58
59     /*
60     ** DOUBLE \ SPARSE
61     */
62     if (TypeL == GenericType::ScilabDouble && TypeR == GenericType::ScilabSparse)
63     {
64         Double *pL = _pLeftOperand->getAs<Double>();
65         Sparse *pR = _pRightOperand->getAs<Sparse>();
66
67         iResult = RDivideSparseByDouble(pR, pL, &pResult);
68     }
69
70     //manage errors
71     if (iResult)
72     {
73         switch (iResult)
74         {
75             case 1 :
76                 throw ast::InternalError(_W("Inconsistent row/column dimensions.\n"));
77             case 2 :
78                 throw ast::InternalError(_W("With NaN or Inf a left division by scalar expected.\n"));
79             case 3 :
80                 throw ast::InternalError(_W("Left division by zero...\n"));
81             case 4 :
82                 sciprint(_("Warning : Left division by zero...\n"));
83                 break;
84             default :
85                 sciprint(_("Operator \\ : Error %d not yet managed.\n"), iResult);
86         }
87     }
88
89     /*
90     ** Default case : Return NULL will Call Overloading.
91     */
92     return pResult;
93 }
94
95 int LDivideDoubleByDouble(Double *_pDouble1, Double *_pDouble2, Double **_pDoubleOut)
96 {
97     types::Double* pDblTmp = NULL;
98     int iErr = 0;
99
100     //check finite values of _pDouble1 and _pDouble2
101     if (isDoubleFinite(_pDouble1) == false || isDoubleFinite(_pDouble2) == false)
102     {
103         if (_pDouble1->isScalar() == false)
104         {
105             return 2;
106         }
107     }
108
109     if (_pDouble1->isScalar())
110     {
111         //x \ Y => Y / x
112         return RDivideDoubleByDouble(_pDouble2, _pDouble1, _pDoubleOut);
113     }
114
115     /* if (_pDouble2->isScalar())
116     {
117         //X \ y => X \ (eye() * y)
118         pDblTmp = new types::Double(_pDouble1->getRows(), _pDouble1->getRows(), _pDouble2->isComplex());
119         double dblScalarReal = _pDouble2->get(0);
120         double* pdblReal = pDblTmp->get();
121         int iSize = pDblTmp->getSize();
122         int iRowsP1 = pDblTmp->getRows() + 1;
123
124         memset(pdblReal, 0x00, iSize * sizeof(double));
125         if (_pDouble2->isComplex())
126         {
127             double dblScalarImag = _pDouble2->getImg(0);
128             double* pdblImag = pDblTmp->getImg();
129             memset(pdblImag, 0x00, iSize * sizeof(double));
130             for (int i = 0; i < iSize; i += iRowsP1)
131             {
132                 pdblReal[i] = dblScalarReal;
133                 pdblImag[i] = dblScalarImag;
134             }
135         }
136         else
137         {
138             for (int i = 0; i < iSize; i += iRowsP1)
139             {
140                 pdblReal[i] = dblScalarReal;
141             }
142         }
143
144         _pDouble2 = pDblTmp;
145     }*/
146
147     if (_pDouble1->getDims() > 2 || _pDouble2->getDims() > 2)
148     {
149         //not managed
150         return 0;
151     }
152
153     if (_pDouble1->getRows() != _pDouble2->getRows())
154     {
155         // matrix dimensions do not agree
156         return 1;
157     }
158
159     *_pDoubleOut = new Double(_pDouble1->getCols(), _pDouble2->getCols(), _pDouble1->isComplex() || _pDouble2->isComplex());
160     if ((*_pDoubleOut)->isComplex())
161     {
162         double dblRcond = 0;
163         iErr = iLeftDivisionOfComplexMatrix(
164                    _pDouble1->getReal(), _pDouble1->getImg(), _pDouble1->getRows(), _pDouble1->getCols(),
165                    _pDouble2->getReal(), _pDouble2->getImg(), _pDouble2->getRows(), _pDouble2->getCols(),
166                    (*_pDoubleOut)->getReal(), (*_pDoubleOut)->getImg(), (*_pDoubleOut)->getRows(), (*_pDoubleOut)->getCols(), &dblRcond);
167     }
168     else
169     {
170         double dblRcond = 0;
171         iErr = iLeftDivisionOfRealMatrix(
172                    _pDouble1->getReal(), _pDouble1->getRows(), _pDouble1->getCols(),
173                    _pDouble2->getReal(), _pDouble2->getRows(), _pDouble2->getCols(),
174                    (*_pDoubleOut)->getReal(), (*_pDoubleOut)->getRows(), (*_pDoubleOut)->getCols(), &dblRcond);
175     }
176
177     if (pDblTmp)
178     {
179         delete pDblTmp;
180     }
181
182     return iErr;
183 }
184