From 958577cac90a99124cd673fde1926781d966d91f Mon Sep 17 00:00:00 2001 From: Ankit Raj Date: Wed, 21 Jun 2017 10:26:59 +0530 Subject: Updated Scilab2C --- src/c/linearAlgebra/balanc/dbalanca.c | 75 +++++++++++++ src/c/linearAlgebra/includes/balanc.h | 27 +++++ src/c/linearAlgebra/includes/rcond.h | 26 +++++ src/c/linearAlgebra/includes/schur.h | 30 ++++++ src/c/linearAlgebra/interfaces/int_balanc.h | 29 +++++ src/c/linearAlgebra/interfaces/int_rcond.h | 25 +++++ src/c/linearAlgebra/interfaces/int_schur.h | 60 +++++++++++ src/c/linearAlgebra/rcond/drconda.c | 46 ++++++++ src/c/linearAlgebra/schur/dgschura.c | 161 ++++++++++++++++++++++++++++ src/c/linearAlgebra/schur/dschura.c | 126 ++++++++++++++++++++++ 10 files changed, 605 insertions(+) create mode 100644 src/c/linearAlgebra/balanc/dbalanca.c create mode 100644 src/c/linearAlgebra/includes/balanc.h create mode 100644 src/c/linearAlgebra/includes/rcond.h create mode 100644 src/c/linearAlgebra/includes/schur.h create mode 100644 src/c/linearAlgebra/interfaces/int_balanc.h create mode 100644 src/c/linearAlgebra/interfaces/int_rcond.h create mode 100644 src/c/linearAlgebra/interfaces/int_schur.h create mode 100644 src/c/linearAlgebra/rcond/drconda.c create mode 100644 src/c/linearAlgebra/schur/dgschura.c create mode 100644 src/c/linearAlgebra/schur/dschura.c (limited to 'src/c/linearAlgebra') diff --git a/src/c/linearAlgebra/balanc/dbalanca.c b/src/c/linearAlgebra/balanc/dbalanca.c new file mode 100644 index 00000000..a86a1967 --- /dev/null +++ b/src/c/linearAlgebra/balanc/dbalanca.c @@ -0,0 +1,75 @@ +/* Copyright (C) 2017 - IIT Bombay - FOSSEE + + 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 + Author: Siddhesh Wani + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ + +/*Funtion to balance given square matrix to improve its condition number*/ +#include "balanc.h" +#include "lapack.h" +#include "stdlib.h" +#include "string.h" +#include "eye.h" + +void dbalanca(double* in1, int rows, double* in2, double* out1, \ + double* out2, double* out3, double* out4) +{ + char JOB = 'B'; + char RSIDE = 'R'; + char LSIDE = 'L'; + + int ILO, IHI, INFO; + double *buf1, *buf2, *LSCALE, *RSCALE; + double *LWORK; + + + if(in2 == NULL) + { + /*Single input matrix*/ + buf1 = (double*) malloc((double) rows*rows*sizeof(double)); + LSCALE = (double*) malloc((double) rows*sizeof(double)); + + /*copy input to temp buf as lapack function modifies input*/ + memcpy(buf1,in1,rows*rows*sizeof(double)); + + dgebal_(&JOB,&rows,buf1,&rows,&ILO,&IHI,LSCALE,&INFO); + deyea(out2,rows,rows); + dgebak_(&JOB,&RSIDE,&rows,&ILO,&IHI,LSCALE,&rows,out2,&rows,&INFO); + + /*copy computed matrix to output*/ + memcpy(out1,buf1,rows*rows*sizeof(double)); + } + else + { + /*two matrices as inputs*/ + buf1 = (double*) malloc((double) rows*rows*sizeof(double)); + buf2 = (double*) malloc((double) rows*rows*sizeof(double)); + LSCALE = (double*) malloc((double) rows*sizeof(double)); + RSCALE = (double*) malloc((double) rows*sizeof(double)); + LWORK = (double*) malloc((double) 6*rows*sizeof(double)); + + /*copy input to temp buf as lapack function modifies input*/ + memcpy(buf1,in1,rows*rows*sizeof(double)); + memcpy(buf2,in2,rows*rows*sizeof(double)); + + dggbal_(&JOB,&rows,buf1,&rows,buf2,&rows,&ILO,&IHI,LSCALE,RSCALE, \ + LWORK,&INFO); + deyea(out3,rows,rows); + deyea(out4,rows,rows); + + dgebak_(&JOB,&LSIDE,&rows,&ILO,&IHI,LSCALE,&rows,out3,&rows,&INFO); + dgebak_(&JOB,&RSIDE,&rows,&ILO,&IHI,LSCALE,&rows,out4,&rows,&INFO); + + /*copy computed matrix to output*/ + memcpy(out1,buf1,rows*rows*sizeof(double)); + memcpy(out2,buf2,rows*rows*sizeof(double)); + + } + +} \ No newline at end of file diff --git a/src/c/linearAlgebra/includes/balanc.h b/src/c/linearAlgebra/includes/balanc.h new file mode 100644 index 00000000..dcc66b28 --- /dev/null +++ b/src/c/linearAlgebra/includes/balanc.h @@ -0,0 +1,27 @@ + /* Copyright (C) 2017 - IIT Bombay - FOSSEE + + 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 + Author: Siddhesh Wani + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ +#ifndef __BALANC_H__ +#define __BALANC_H__ + + +#ifdef __cplusplus +extern "C" { +#endif + +void dbalanca(double* in1, int rows, double* in2, double* out1, \ + double* out2, double* out3, double* out4); + +#ifdef __cplusplus +} /* extern "C" */ +#endif + +#endif /*__BALANC_H__*/ diff --git a/src/c/linearAlgebra/includes/rcond.h b/src/c/linearAlgebra/includes/rcond.h new file mode 100644 index 00000000..4796f029 --- /dev/null +++ b/src/c/linearAlgebra/includes/rcond.h @@ -0,0 +1,26 @@ + /* Copyright (C) 2017 - IIT Bombay - FOSSEE + + 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 + Author: Siddhesh Wani + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ +#ifndef __RCOND_H__ +#define __RCOND_H__ + + +#ifdef __cplusplus +extern "C" { +#endif + +double drconda(double* in1, int rows); + +#ifdef __cplusplus +} /* extern "C" */ +#endif + +#endif /*__RCOND_H__*/ diff --git a/src/c/linearAlgebra/includes/schur.h b/src/c/linearAlgebra/includes/schur.h new file mode 100644 index 00000000..27b5c2aa --- /dev/null +++ b/src/c/linearAlgebra/includes/schur.h @@ -0,0 +1,30 @@ + /* Copyright (C) 2017 - IIT Bombay - FOSSEE + + 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 + Author: Siddhesh Wani + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ +#ifndef __SCHUR_H__ +#define __SCHUR_H__ + + +#ifdef __cplusplus +extern "C" { +#endif + +double dschura(double* in1, int size, int flag, int nout, double* out1, \ + double* out2); + +double dgschura(double* in1, int size, double* in2, int flag, int nout, \ + double* out1, double* out2, double* out3, double* out4); + +#ifdef __cplusplus +} /* extern "C" */ +#endif + +#endif /*__SCHUR_H__*/ diff --git a/src/c/linearAlgebra/interfaces/int_balanc.h b/src/c/linearAlgebra/interfaces/int_balanc.h new file mode 100644 index 00000000..a16ba8c2 --- /dev/null +++ b/src/c/linearAlgebra/interfaces/int_balanc.h @@ -0,0 +1,29 @@ + /* Copyright (C) 2017 - IIT Bombay - FOSSEE + + 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 + Author: Siddhesh Wani + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ +#ifndef __INT_BALANC_H__ +#define __INT_BALANC_H__ + +#ifdef __cplusplus +extern "C" { +#endif + +#define d2balancd2d2(in1,size1,out1,out2) dbalanca(in1,size1[0],NULL,out1, \ + out2,NULL,NULL) + +#define d2d2balancd2d2d2d2(in1,size1,in2,size2,out1,out2,out3,out4) \ + dbalanca(in1,size1[0],in2,out1,out2,out3,out4) + +#ifdef __cplusplus +} /* extern "C" */ +#endif + +#endif /*__INT_BALANC_H__*/ diff --git a/src/c/linearAlgebra/interfaces/int_rcond.h b/src/c/linearAlgebra/interfaces/int_rcond.h new file mode 100644 index 00000000..6e6a4454 --- /dev/null +++ b/src/c/linearAlgebra/interfaces/int_rcond.h @@ -0,0 +1,25 @@ + /* Copyright (C) 2017 - IIT Bombay - FOSSEE + + 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 + Author: Siddhesh Wani + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ +#ifndef __INT_RCOND_H__ +#define __INT_RCOND_H__ + +#ifdef __cplusplus +extern "C" { +#endif + +#define d2rcondd0(in1,size1) drconda(in1,size1[0]) + +#ifdef __cplusplus +} /* extern "C" */ +#endif + +#endif /*__INT_RCOND_H__*/ diff --git a/src/c/linearAlgebra/interfaces/int_schur.h b/src/c/linearAlgebra/interfaces/int_schur.h new file mode 100644 index 00000000..81324e66 --- /dev/null +++ b/src/c/linearAlgebra/interfaces/int_schur.h @@ -0,0 +1,60 @@ + /* Copyright (C) 2017 - IIT Bombay - FOSSEE + + 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 + Author: Siddhesh Wani + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ +#ifndef __INT_SCHUR_H__ +#define __INT_SCHUR_H__ + +#ifdef __cplusplus +extern "C" { +#endif + +#define d2schurd2(in1,size1,out1) dschura(in1,size1[0],0,1,out1,NULL) + +#define d2schurd2d2(in1,size1,out1,out2) dschura(in1,size1[0],0,2,out1,out2) + +#define d2g2schurd2(in1,size1,in2,size2,out1) (in2[0]=='c')? \ + dschura(in1,size1[0],1,1,out1,NULL): \ + dschura(in1,size1[0],2,1,out1,NULL) + +#define d2g2schurd2d0(in1,size1,in2,size2,out1) (in2[0]=='c')? \ + dschura(in1,size1[0],1,2,out1,NULL): \ + dschura(in1,size1[0],2,2,out1,NULL) + +#define d2g2schurd2d0d2(in1,size1,in2,size2,out1,out2) (in2[0]=='c')? \ + dschura(in1,size1[0],1,3,out1,out2): \ + dschura(in1,size1[0],2,3,out1,out2) + +#define d2d2schurd2d2(in1,size1,in2,size2,out1,out2) dgschura(in1,size1[0], \ + in2,0,2,out1,out2,NULL,NULL) + +#define d2d2schurd2d2d2d2(in1,size1,in2,size2,out1,out2,out3,out4) \ + dgschura(in1,size1[0],in2,0,4,out1,out2,out3,out4) + +#define d2d2g2schurd0(in1,size1,in2,size2,in3,size3) dgschura(in1,size1[0], \ + in2,1,1,NULL,NULL,NULL,NULL) + +#define d2d2g2schurd2d0(in1,size1,in2,size2,in3,size3,out1) \ + dgschura(in1,size1[0],in2,1,2,out1,NULL,NULL,NULL) + +#define d2d2g2schurd2d2d0(in1,size1,in2,size2,in3,size3,out1,out2) \ + dgschura(in1,size1[0],in2,1,3,out1,out2,NULL,NULL) + +#define d2d2g2schurd2d2d2d0(in1,size1,in2,size2,in3,size3,out1,out2,out3) \ + dgschura(in1,size1[0],in2,1,4,out1,out2,out3,NULL) + +#define d2d2g2schurd2d2d2d2d0(in1,size1,in2,size2,in3,size3,out1,out2,out3, \ + out4) dgschura(in1,size1[0],in2,1,5,out1,out2,out3,out4) + +#ifdef __cplusplus +} /* extern "C" */ +#endif + +#endif /*__INT_SCHUR_H__*/ diff --git a/src/c/linearAlgebra/rcond/drconda.c b/src/c/linearAlgebra/rcond/drconda.c new file mode 100644 index 00000000..a203c1e3 --- /dev/null +++ b/src/c/linearAlgebra/rcond/drconda.c @@ -0,0 +1,46 @@ +/* Copyright (C) 2017 - IIT Bombay - FOSSEE + + 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 + Author: Siddhesh Wani + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ + +/*Funtion to find inverse condition number of square matrix*/ + +#include "lapack.h" +#include "stdlib.h" +#include "string.h" + +double drconda(double* in1, int rows) +{ + double *buf, *LDWORK; + int INFO, *IPIV, *LIWORK; + char one = '1'; + double ANORM; + double RCOND = 1; + + buf = (double*) malloc((double) rows*rows*sizeof(double)); + IPIV = (int*) malloc((int) rows*sizeof(int)); + LIWORK = (int*) malloc((int) rows*sizeof(int)); + LDWORK = (double*) malloc((double) 4*rows*sizeof(double)); + + /*Copy input in temp buf, as lapack modifies input*/ + memcpy(buf,in1,rows*rows*sizeof(double)); + + ANORM = dlange_(&one,&rows,&rows,buf,&rows,LDWORK); + + dgetrf_(&rows,&rows,buf,&rows,IPIV,&INFO); + + if(INFO == 0) + { + dgecon_(&one,&rows,buf,&rows,&ANORM,&RCOND,LDWORK,LIWORK,&INFO); + } + + return RCOND; + +} \ No newline at end of file diff --git a/src/c/linearAlgebra/schur/dgschura.c b/src/c/linearAlgebra/schur/dgschura.c new file mode 100644 index 00000000..f17da8c8 --- /dev/null +++ b/src/c/linearAlgebra/schur/dgschura.c @@ -0,0 +1,161 @@ +/* Copyright (C) 2017 - IIT Bombay - FOSSEE + + 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 + Author: Siddhesh Wani + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ + +/*Fucntion to find generalised schur decomposition of given square matrices */ +#include "schur.h" +#include "lapack.h" +#include "stdlib.h" +#include "string.h" +#include + +/*flag --> 0: nothing + --> 1: continuous + --> 2: discrete +*/ + +lapack_logical selctg21( double* in1, double* in2, double* in3); +lapack_logical selctg22( double* in1, double* in2, double* in3); + +double dgschura(double* in1, int size, double* in2, int flag, int nout, \ + double* out1, double* out2, double* out3, double* out4) +{ + char JOBVSL = 'N'; + char JOBVSR = 'N'; + char SORT = 'N'; + int SDIM = 0; + int LDVSL = size, LDVSR = size; + int LWORK = 8*size+16, INFO; + double *ALPHAR, *ALPHAI, *BETA, *VSL, *VSR, *WORK; + int *BWORK; + double ret = 0; + double *buf1, *buf2; /*input is copied to buf, since lapack function directly + modifies the input variable*/ + + /*Used incase of flag > 0*/ + LAPACK_D_SELECT3 selctg = &selctg21; + + if(nout >= 2){ + JOBVSL = 'V'; + JOBVSR = 'V'; + } + if(flag > 0) SORT = 'S'; + + buf1 = (double*) malloc((double) size*size*sizeof(double)); + buf2 = (double*) malloc((double) size*size*sizeof(double)); + ALPHAR = (double*) malloc((double) size*sizeof(double)); + ALPHAI = (double*) malloc((double) size*sizeof(double)); + BETA = (double*) malloc((double) size*sizeof(double)); + VSL = (double*) malloc((double) LDVSL*size*sizeof(double)); + VSR = (double*) malloc((double) LDVSR*size*sizeof(double)); + WORK = (double*) malloc((double) LWORK*sizeof(double)); + BWORK = (int*) malloc((double) size*sizeof(double)); + + + memcpy(buf1,in1,size*size*sizeof(double)); + memcpy(buf2,in2,size*size*sizeof(double)); + + dgges_(&JOBVSL,&JOBVSR,&SORT,selctg,&size,buf1,&size,buf2,&size,&SDIM, \ + ALPHAR,ALPHAI,BETA,VSL,&LDVSL,VSR,&LDVSR,WORK,&LWORK,BWORK,&INFO); + + /*if (INFO != 0) + { + out1 = NULL; + return 0; + }*/ + + if(nout == 1) + { + return(SDIM); + } + else if(nout == 2) + { + if(flag == 0) + { + /*copy in1 to out1 and in2 to out2*/ + memcpy(out1,buf1,size*size*sizeof(double)); + memcpy(out2,buf2,size*size*sizeof(double)); + } + else + { + /*copy VSR to out1 and return SDIM */ + memcpy(out1,VSR,size*size*sizeof(double)); + ret = SDIM; + } + } + else if(nout == 3) + { + /*copy VSL to out1, VSR to out2, return SDIM*/ + memcpy(out1,VSL,size*size*sizeof(double)); + memcpy(out2,VSR,size*size*sizeof(double)); + ret = SDIM; + } + else if(nout == 4) + { + if(flag == 0) + { + /*copy in1 to out1 and in2 to out2*/ + memcpy(out1,buf1,size*size*sizeof(double)); + memcpy(out2,buf2,size*size*sizeof(double)); + /*copy VSL to out3 and VSR to out4*/ + memcpy(out3,VSL,size*size*sizeof(double)); + memcpy(out4,VSR,size*size*sizeof(double)); + } + else + { + /*copy in1 to out1 and in2 to out2*/ + memcpy(out1,buf1,size*size*sizeof(double)); + memcpy(out2,buf2,size*size*sizeof(double)); + /*copy VSR to out3 and return SDIM*/ + memcpy(out3,VSR,size*size*sizeof(double)); + ret = SDIM; + } + } + else /*nout = 5*/ + { + /*copy in1 to out1 and in2 to out2*/ + memcpy(out1,buf1,size*size*sizeof(double)); + memcpy(out2,buf2,size*size*sizeof(double)); + /*copy VSL to out3 and VSR to out4*/ + memcpy(out3,VSL,size*size*sizeof(double)); + memcpy(out4,VSR,size*size*sizeof(double)); + /*return SDIM*/ + ret = SDIM; + } + + free(buf1); + free(buf2); + free(VSL); + free(VSR); + free(ALPHAR); + free(ALPHAI); + free(BETA); + free(WORK); + free(BWORK); + + return ret; +} + +lapack_logical selctg21(double* in1, double* in2, double* in3) +{ + if((sqrt(*in1**in1+*in2**in2)/ *in3) < 1) + return 1; + else + return 0; +} + +lapack_logical selctg22(double* in1, double* in2, double* in3) +{ + if((sqrt(*in1**in1+*in2**in2)/ *in3) < 1) + return 1; + else + return 0; +} \ No newline at end of file diff --git a/src/c/linearAlgebra/schur/dschura.c b/src/c/linearAlgebra/schur/dschura.c new file mode 100644 index 00000000..802caa81 --- /dev/null +++ b/src/c/linearAlgebra/schur/dschura.c @@ -0,0 +1,126 @@ +/* Copyright (C) 2017 - IIT Bombay - FOSSEE + + 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 + Author: Siddhesh Wani + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in + */ + +/*Funtion to find schur decomposition of given square matrix */ +#include "schur.h" +#include "lapack.h" +#include "stdlib.h" +#include "string.h" +#include + +/*flag --> 0: nothing + --> 1: continuous + --> 2: discrete +*/ + +lapack_logical selctg11( double* in1, double* in2); +lapack_logical selctg12( double* in1, double* in2); + +double dschura(double* in1, int size, int flag, int nout, double* out1, \ + double* out2) +{ + char JOBVS = 'N'; + char SORT = 'N'; + int SDIM = 0; + int LVDS = size; + int LWORK = 3*size, INFO; + double *WR, *WI, *VS, *WORK; + int *BWORK; + double ret = 0; + double *buf; /*input is copied to buf, since lapack function direclty + modifies the input variable*/ + + + /*Used incase of flag > 0*/ + LAPACK_D_SELECT2 selctg; + + if(flag == 1 || flag == 0) + selctg = &selctg11; + else if(flag == 2) + selctg = &selctg12; + + if(nout >= 2) JOBVS = 'V'; + if(flag > 0) SORT = 'S'; + + buf = (double*) malloc((double) size*size*sizeof(double)); + WR = (double*) malloc((double) size*sizeof(double)); + WI = (double*) malloc((double) size*sizeof(double)); + VS = (double*) malloc((double) LVDS*size*sizeof(double)); + WORK = (double*) malloc((double) LWORK*sizeof(double)); + BWORK = (int*) malloc((double) size*sizeof(double)); + + + memcpy(buf,in1,size*size*sizeof(double)); + + dgees_(&JOBVS,&SORT,selctg,&size,buf,&size,&SDIM,WR,WI,VS,&LVDS, \ + WORK,&LWORK,BWORK,&INFO); + + /*if (INFO != 0) + { + out1 = NULL; + return 0; + } */ + + if(nout == 1) + { + /*Copy result in in1 to out1*/ + memcpy(out1,buf,size*size*sizeof(double)); + } + else if(nout == 2) + { + if(flag == 0) + { + /*copy in1 to out2 and VS to out1*/ + memcpy(out2,buf,size*size*sizeof(double)); + memcpy(out1,VS,size*size*sizeof(double)); + } + else + { + /*copy VS to out1 and SDIM to out2*/ + memcpy(out1,VS,size*size*sizeof(double)); + ret = SDIM; + } + } + else + { + /*copy VS to out1, SDIM to out2, in1 to out3*/ + memcpy(out1,VS,size*size*sizeof(double)); + memcpy(out2,buf,size*size*sizeof(double)); + ret = SDIM; + } + + free(buf); + free(WI); + free(WR); + free(VS); + free(WORK); + free(BWORK); + + return ret; +} + +lapack_logical selctg11(double* in1, double* in2) +{ + if(*in1 <= 0) + return 1; + else + return 0; +} + +lapack_logical selctg12(double* in1, double* in2) +{ + + if(sqrt(*in1**in1+*in2**in2) < 1) + return 1; + else + return 0; +} \ No newline at end of file -- cgit