summaryrefslogtreecommitdiff
path: root/src/c/linearAlgebra
diff options
context:
space:
mode:
Diffstat (limited to 'src/c/linearAlgebra')
-rw-r--r--src/c/linearAlgebra/includes/schur.h30
-rw-r--r--src/c/linearAlgebra/interfaces/int_schur.h60
-rw-r--r--src/c/linearAlgebra/schur/dgschura.c152
-rw-r--r--src/c/linearAlgebra/schur/dschura.c111
4 files changed, 353 insertions, 0 deletions
diff --git a/src/c/linearAlgebra/includes/schur.h b/src/c/linearAlgebra/includes/schur.h
new file mode 100644
index 0000000..27b5c2a
--- /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_schur.h b/src/c/linearAlgebra/interfaces/int_schur.h
new file mode 100644
index 0000000..81324e6
--- /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/schur/dgschura.c b/src/c/linearAlgebra/schur/dgschura.c
new file mode 100644
index 0000000..56d4845
--- /dev/null
+++ b/src/c/linearAlgebra/schur/dgschura.c
@@ -0,0 +1,152 @@
+/* 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"
+
+
+/*flag --> 0: nothing
+ --> 1: continuous
+ --> 2: discrete
+*/
+
+lapack_logical selctg2( 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 = &selctg2;
+
+ 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 selctg2(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 0000000..c31ddca
--- /dev/null
+++ b/src/c/linearAlgebra/schur/dschura.c
@@ -0,0 +1,111 @@
+/* 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 schur decomposition of given square matrix */
+#include "schur.h"
+#include "lapack.h"
+#include "stdlib.h"
+#include "string.h"
+
+
+/*flag --> 0: nothing
+ --> 1: continuous
+ --> 2: discrete
+*/
+
+lapack_logical selctg1( 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 = &selctg1;
+
+ 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 selctg1(double* in1, double* in2)
+{
+
+ if(sqrt(*in1**in1+*in2**in2) < 1)
+ return 1;
+ else
+ return 0;
+} \ No newline at end of file
generated by cgit (git 2.34.1) at 2024-12-20 22:04:23 +0000