* Copyright (C) 2009 - DIGITEO - Bernard HUGUENEY
* Copyright (C) 2011 - DIGITEO - Cedric DELAMARRE
*
-* This file must be used under the terms of the CeCILL.
-* This source file is licensed as described in the file COPYING, which
-* you should have received as part of this distribution. The terms
-* are also available at
-* http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
+ * Copyright (C) 2012 - 2016 - Scilab Enterprises
+ *
+ * This file is hereby licensed under the terms of the GNU GPL v2.0,
+ * pursuant to article 5.3.4 of the CeCILL v.2.1.
+ * This file was originally licensed under the terms of the CeCILL v2.1,
+ * and continues to be available under such terms.
+ * For more information, see the COPYING file which you should have received
+ * along with this program.
*
*/
/*--------------------------------------------------------------------------*/
#include "function.hxx"
#include "double.hxx"
#include "overload.hxx"
-#include "execvisitor.hxx"
extern "C"
{
types::Function::ReturnValue sci_spec(types::typed_list &in, int _iRetCount, types::typed_list &out)
{
- types::Double* pDblA = NULL;
- types::Double* pDblB = NULL;
double* pDataA = NULL;
double* pDataB = NULL;
bool symmetric = FALSE;
int iRet = 0;
- if ((in.size() != 1) && (in.size() != 2))
+ _iRetCount = std::max(1, _iRetCount);
+
+ if (in.size() != 1 && in.size() != 2)
{
Scierror(77, _("%s: Wrong number of input argument(s): %d to %d expected.\n"), "spec", 1, 2);
return types::Function::Error;
return types::Function::Error;
}
- if ((in[0]->isDouble() == false))
+ if (in[0]->isDouble() == false)
{
- std::wstring wstFuncName = L"%" + in[0]->getShortTypeStr() + L"_spec";
- return Overload::call(wstFuncName, in, _iRetCount, out, new ast::ExecVisitor());
+ std::wstring wstFuncName = L"%" + in[0]->getShortTypeStr() + L"_spec";
+ return Overload::call(wstFuncName, in, _iRetCount, out);
}
- pDblA = in[0]->getAs<types::Double>()->clone()->getAs<types::Double>();
- if (pDblA->getCols() != pDblA->getRows())
+ types::Double* in0 = in[0]->getAs<types::Double>();
+
+ if (in0->getCols() != in0->getRows())
{
Scierror(20, _("%s: Wrong type for input argument #%d: Square matrix expected.\n"), "spec", 1);
return types::Function::Error;
}
- if ((pDblA->getRows() == -1) || (pDblA->getCols() == -1)) // manage eye case
+ if (in0->getRows() == -1 || in0->getCols() == -1) // manage eye case
{
Scierror(271, _("%s: Size varying argument a*eye(), (arg %d) not allowed here.\n"), "spec", 1);
return types::Function::Error;
}
- if ((pDblA->getCols() == 0) || (pDblA->getRows() == 0)) // size null
+ if (in0->getCols() == 0 || in0->getRows() == 0) // size null
{
out.push_back(types::Double::Empty());
for (int i = 1; i < _iRetCount; i++)
return types::Function::OK;
}
+ types::Double* pDblA = in0->clone()->getAs<types::Double>();
+
if (in.size() == 1)
{
types::Double* pDblEigenValues = NULL;
pDataA = (double*)oGetDoubleComplexFromPointer(pDblA->getReal(), pDblA->getImg(), pDblA->getSize());
if (!pDataA)
{
+ pDblA->killMe();
Scierror(999, _("%s: Cannot allocate more memory.\n"), "spec");
return types::Function::Error;
}
int totalSize = pDblA->getSize();
if ((pDblA->isComplex() ? C2F(vfiniteComplex)(&totalSize, (doublecomplex*)pDataA) : C2F(vfinite)(&totalSize, pDataA)) == false)
{
+ if (pDblA->isComplex())
+ {
+ vFreeDoubleComplexFromPointer((doublecomplex*)pDataA);
+ }
+ pDblA->killMe();
Scierror(264, _("%s: Wrong value for input argument %d: Must not contain NaN or Inf.\n"), "spec", 1);
return types::Function::Error;
}
symmetric = isSymmetric(pDblA->getReal(), pDblA->getImg(), pDblA->isComplex(), pDblA->getRows(), pDblA->getCols()) == 1;
- int eigenValuesCols = (_iRetCount == 1) ? 1 : pDblA->getCols();
+ int eigenValuesCols = (_iRetCount <= 1) ? 1 : pDblA->getCols();
if (symmetric)
{
if (iRet < 0)
{
- Scierror(998, _("%s: On entry to ZGEEV parameter number 3 had an illegal value (lapack library problem).\n"), "spec", iRet);
+ vFreeDoubleComplexFromPointer((doublecomplex*)pDataA);
+ pDblA->killMe();
+ Scierror(998, _("%s: On entry to ZGEEV parameter number 3 had an illegal value (lapack library problem).\n"), "spec");
return types::Function::Error;
}
if (iRet > 0)
{
+ vFreeDoubleComplexFromPointer((doublecomplex*)pDataA);
+ pDblA->killMe();
Scierror(24, _("%s: Convergence problem, %d off-diagonal elements of an intermediate tridiagonal form did not converge to zero.\n"), "spec", iRet);
return types::Function::Error;
}
if (_iRetCount == 2)
{
- vGetPointerFromDoubleComplex((doublecomplex*)pDataA, pDblA->getSize() , pDblEigenVectors->getReal(), pDblEigenVectors->getImg());
+ vGetPointerFromDoubleComplex((doublecomplex*)pDataA, pDblA->getSize(), pDblEigenVectors->getReal(), pDblEigenVectors->getImg());
expandToDiagonalOfMatrix(pDblEigenValues->getReal(), pDblA->getCols());
out.push_back(pDblEigenVectors);
}
out.push_back(pDblEigenValues);
+ pDblA->killMe();
+ vFreeDoubleComplexFromPointer((doublecomplex*)pDataA);
}
else // not symmetric
{
doublecomplex* pEigenValues = (doublecomplex*)MALLOC(pDblA->getCols() * sizeof(doublecomplex));
doublecomplex* pEigenVectors = pDblEigenVectors ? (doublecomplex*)MALLOC(sizeof(doublecomplex) * pDblA->getSize()) : NULL;
iRet = iEigen1ComplexM((doublecomplex*)pDataA, pDblA->getCols(), pEigenValues, pEigenVectors);
+ vFreeDoubleComplexFromPointer((doublecomplex*)pDataA);
if (iRet < 0)
{
- Scierror(998, _("%s: On entry to ZHEEV parameter number 3 had an illegal value (lapack library problem).\n"), "spec", iRet);
+ pDblA->killMe();
+ Scierror(998, _("%s: On entry to ZHEEV parameter number 3 had an illegal value (lapack library problem).\n"), "spec");
return types::Function::Error;
}
if (iRet > 0)
{
+ pDblA->killMe();
Scierror(24, _("%s: The QR algorithm failed to compute all the eigenvalues, and no eigenvectors have been computed. Elements and %d+1:N of W contain eigenvalues which have converged.\n"), "spec", iRet);
return types::Function::Error;
}
}
out.push_back(pDblEigenValues);
FREE(pEigenValues);
+ pDblA->killMe();
}
}
else // real
if (symmetric)
{
iRet = iEigen1RealSymmetricM(pDataA, pDblA->getCols(), (_iRetCount == 2), pDblEigenValues->getReal());
-
if (iRet < 0)
{
- Scierror(998, _("%s: On entry to ZGEEV parameter number 3 had an illegal value (lapack library problem).\n"), "spec", iRet);
+ pDblA->killMe();
+ Scierror(998, _("%s: On entry to ZGEEV parameter number 3 had an illegal value (lapack library problem).\n"), "spec");
return types::Function::Error;
}
if (iRet > 0)
{
+ pDblA->killMe();
Scierror(24, _("%s: Convergence problem, %d off-diagonal elements of an intermediate tridiagonal form did not converge to zero.\n"), "spec", iRet);
return types::Function::Error;
}
if (_iRetCount == 2)
{
+ if (pDblEigenVectors)
+ {
+ pDblEigenVectors->killMe();
+ }
expandToDiagonalOfMatrix(pDblEigenValues->getReal(), pDblA->getCols());
out.push_back(pDblA);
}
+ else
+ {
+ pDblA->killMe();
+ }
+
out.push_back(pDblEigenValues);
}
else // not symmetric
if (iRet < 0)
{
- Scierror(998, _("%s: On entry to ZHEEV parameter number 3 had an illegal value (lapack library problem).\n"), "spec", iRet);
+ pDblA->killMe();
+ if (pDblEigenVectors)
+ {
+ pDblEigenVectors->killMe();
+ }
+ Scierror(998, _("%s: On entry to ZHEEV parameter number 3 had an illegal value (lapack library problem).\n"), "spec");
return types::Function::Error;
}
if (iRet > 0)
{
+ pDblA->killMe();
+ if (pDblEigenVectors)
+ {
+ pDblEigenVectors->killMe();
+ }
Scierror(24, _("%s: The QR algorithm failed to compute all the eigenvalues, and no eigenvectors have been computed. Elements and %d+1:N of WR and WI contain eigenvalues which have converged.\n"), "spec", iRet);
return types::Function::Error;
}
expandToDiagonalOfMatrix(pDblEigenValues->getImg(), pDblA->getCols());
out.push_back(pDblEigenVectors);
}
+
out.push_back(pDblEigenValues);
+ pDblA->killMe();
}
}
+
+ return types::Function::OK;
}
if (in.size() == 2)
doublecomplex* pAlpha = NULL;
bool bIsComplex = false;
- if ((in[1]->isDouble() == false))
+ if (in[1]->isDouble() == false)
+ {
+ std::wstring wstFuncName = L"%" + in[1]->getShortTypeStr() + L"_spec";
+ return Overload::call(wstFuncName, in, _iRetCount, out);
+ }
+
+ types::Double* in1 = in[1]->getAs<types::Double>();
+
+ if (in1->getCols() != in1->getRows())
{
- std::wstring wstFuncName = L"%" + in[1]->getShortTypeStr() + L"_spec";
- return Overload::call(wstFuncName, in, _iRetCount, out, new ast::ExecVisitor());
+ Scierror(20, _("%s: Wrong type for input argument #%d: Square matrix expected.\n"), "spec", 2);
+ return types::Function::Error;
}
- pDblB = in[1]->getAs<types::Double>()->clone()->getAs<types::Double>();
- if ((pDblA->getRows() != pDblB->getRows()) && (pDblA->getCols() != pDblB->getCols()))
+
+ if (pDblA->getRows() != in1->getRows() && pDblA->getCols() != in1->getCols())
{
+ pDblA->killMe();
Scierror(999, _("%s: Arguments %d and %d must have equal dimensions.\n"), "spec", 1, 2);
return types::Function::Error;
}
//chekc if A and B are real complex or with imag part at 0
- if (isNoZeroImag(pDblA) == false && isNoZeroImag(pDblB) == false)
+ if (isNoZeroImag(pDblA) == false && isNoZeroImag(in1) == false)
{
//view A and B as real matrix
bIsComplex = false;
}
else
{
- bIsComplex = pDblA->isComplex() || pDblB->isComplex();
+ bIsComplex = pDblA->isComplex() || in1->isComplex();
}
+ types::Double* pDblB = in1->clone()->getAs<types::Double>();
if (bIsComplex)
{
if (pDblA->isComplex() == false)
pDataA = (double*)oGetDoubleComplexFromPointer(pDblA->getReal(), pDblA->getImg(), pDblA->getSize());
pDataB = (double*)oGetDoubleComplexFromPointer(pDblB->getReal(), pDblB->getImg(), pDblB->getSize());
- if (!pDblA || !pDblB)
+ if (!pDataA && !pDataB)
{
Scierror(999, _("%s: Cannot allocate more memory.\n"), "spec");
return types::Function::Error;
}
+
+ if (!pDataA)
+ {
+ vFreeDoubleComplexFromPointer((doublecomplex*)pDataB);
+ Scierror(999, _("%s: Cannot allocate more memory.\n"), "spec");
+ return types::Function::Error;
+ }
+
+ if (!pDataB)
+ {
+ vFreeDoubleComplexFromPointer((doublecomplex*)pDataA);
+ Scierror(999, _("%s: Cannot allocate more memory.\n"), "spec");
+ return types::Function::Error;
+ }
}
else
{
if ((pDblA->isComplex() ? C2F(vfiniteComplex)(&totalSize, (doublecomplex*)pDataA) : C2F(vfinite)(&totalSize, pDataA)) == false)
{
+ pDblA->killMe();
+ pDblB->killMe();
Scierror(264, _("%s: Wrong value for input argument %d: Must not contain NaN or Inf.\n"), "spec", 1);
return types::Function::Error;
}
if ((pDblB->isComplex() ? C2F(vfiniteComplex)(&totalSize, (doublecomplex*)pDataB) : C2F(vfinite)(&totalSize, pDataB)) == false)
{
+ pDblA->killMe();
+ pDblB->killMe();
Scierror(264, _("%s: Wrong value for input argument %d: Must not contain NaN or Inf.\n"), "spec", 2);
return types::Function::Error;
}
if (iRet < 0)
{
- Scierror(998, _("%s: On entry to ZHEEV parameter number 3 had an illegal value (lapack library problem).\n"), "spec", iRet);
+ pDblA->killMe();
+ pDblB->killMe();
+ Scierror(998, _("%s: On entry to ZHEEV parameter number 3 had an illegal value (lapack library problem).\n"), "spec");
return types::Function::Error;
}
}
else
{
- if (iRet == pDblA->getCols() + 1)
+ if (iRet == pDblA->getCols() + 1)
{
Scierror(999, _("%s: Other than QZ iteration failed in DHGEQZ.\n"), "spec");
}
- if (iRet == pDblA->getCols() + 2)
+ if (iRet == pDblA->getCols() + 2)
{
Scierror(999, _("%s: Error return from DTGEVC.\n"), "spec");
}
{
Scierror(24, _("%s: The QR algorithm failed to compute all the eigenvalues, and no eigenvectors have been computed. Elements and %d+1:N of W contain eigenvalues which have converged.\n"), "spec", iRet);
}
+
+ pDblA->killMe();
+ pDblB->killMe();
+ if (pDataA)
+ {
+ vFreeDoubleComplexFromPointer((doublecomplex*)pDataA);
+ }
+
+ if (pDataB)
+ {
+ vFreeDoubleComplexFromPointer((doublecomplex*)pDataB);
+ }
return types::Function::Error;
}
{
vFreeDoubleComplexFromPointer(pR);
}
- if (pDblB->isComplex())
+ if (bIsComplex && pDblB->isComplex())
{
vFreeDoubleComplexFromPointer((doublecomplex*)pDataB);
}
+ pDblB->killMe();
} // if(in.size() == 2)
projects / scilab.git / blobdiff