Scilab [trunk]: scilab/modules/cacsd/sci_gateway/c/sci

00001 /*-----------------------------------------------------------------------------------*/ 
00002 /* INRIA */
00003 /*-----------------------------------------------------------------------------------*/ 
00004 #include "machine.h"
00005 #include "stack-c.h"
00006 #include <string.h>
00007 #include "sci_contr.h"
00008 /*-----------------------------------------------------------------------------------*/ 
00009 #define MAX(x,y)        (((x)>(y))?(x):(y))
00010 /*-----------------------------------------------------------------------------------*/ 
00011 extern double C2F(dlamch)  __PARAMS((char *CMACH, unsigned long int));
00012 extern int C2F(ab01od)();
00013 /*-----------------------------------------------------------------------------------*/ 
00014 int intab01od(char* fname)
00015 { 
00016         int mA,nA,ptrA, mB,nB, ptrB;
00017         int A,B,U,KSTAIR,V, ptrIWORK,ptrU,ptrTOL;
00018         int ptrKSTAIR,ptrV,ptrDWORK,ptrJUNK,ptrNCONT;
00019         int LDA, LDB, LDU, LDV, LDWORK;
00020         int N, M, mtol, ntol; int un; int INFO, INDCON, NCONT;
00021         char  *JOBU, *JOBV;
00022         double theTOL;
00023         int minlhs=1, minrhs=2, maxlhs=6, maxrhs=3;
00024 
00025         /*     [NCONT,U,KSTAIR,V,A,B]=ab01od(A,B,[TOL])   */
00026 
00027         CheckRhs(minrhs,maxrhs);  CheckLhs(minlhs,maxlhs);
00028         theTOL=(double) C2F(dlamch)("e",1L);
00029         GetRhsVar(1,"d",&mA,&nA,&ptrA);   A=1;        /*     A */
00030         N=mA;
00031         theTOL=0.2*sqrt(2*theTOL)*N;
00032         GetRhsVar(2,"d",&mB,&nB,&ptrB);   B=2;        /*     B */
00033         M=nB;
00034         if (nA != mB || mA != nA )
00035         { Scierror(999,"Invalid A,B matrices \r\n");  return 0; }
00036         if (Rhs == 3) {
00037                 /*    TOL is given:   ab01od(A,B,tol)   */
00038                 GetRhsVar(3,"d",&mtol,&ntol,&ptrTOL);  theTOL=*stk(ptrTOL);    /*     TOL */
00039                 if (theTOL>1.0||theTOL<0.0) {
00040                         Scierror(999,"TOL must be in [0 1]\r\n");  return 0; 
00041                 }
00042         }
00043 
00044         /*     dimensions...    */
00045         LDA=MAX(1,N);  LDB=LDA;  LDU=LDA; LDV=MAX(1,M);
00046         LDWORK = MAX(1, N*M + MAX(N,M) + MAX(N,3*M));
00047 
00048         /*     other parameters of AB01OD   */
00049         JOBU= "N"; if (Lhs >= 2)  JOBU="I";
00050         JOBV= "N"; if (Lhs >= 4)  JOBV="I";
00051 
00052         /*     creating NCONT,U,KSTAIR,V,IWORK,DWORK   */
00053         CreateVar(Rhs+1,"i",(un=1,&un),(un=1,&un),&ptrNCONT);  NCONT=Rhs+1;
00054         CreateVar(Rhs+2,"d",&N,&N,&ptrU);  U=Rhs+2;
00055         CreateVar(Rhs+3,"i",(un=1,&un),&N,&ptrKSTAIR);  KSTAIR=Rhs+3;
00056         CreateVar(Rhs+4,"d",&M,&M,&ptrV);  V=Rhs+4;
00057         CreateVar(Rhs+5,"i",(un=1,&un),&M,&ptrIWORK);
00058         CreateVar(Rhs+6,"d",(un=1,&un),&LDWORK,&ptrDWORK);
00059         C2F(ab01od)( "A", JOBU, JOBV, &N, &M, stk(ptrA), &LDA, 
00060                 stk(ptrB), &LDB, stk(ptrU), &LDU, stk(ptrV), &LDV, 
00061                 istk(ptrNCONT), &INDCON, istk(ptrKSTAIR), &theTOL,  
00062                 istk(ptrIWORK), stk(ptrDWORK), &LDWORK, &INFO );
00063         if (INFO != 0) {
00064                 C2F(errorinfo)("ab01od", &INFO, 6L);
00065                 return 0;
00066         }
00067         if (Lhs >= 3) {
00068                 /*     resizing KSTAIR      */
00069                 CreateVar(Rhs+7,"i",(un=1,&un),&INDCON,&ptrJUNK); 
00070                 KSTAIR=Rhs+7;
00071                 C2F(icopy)(&INDCON,istk(ptrKSTAIR),(un=1,&un),istk(ptrJUNK),(un=1,&un)); }
00072         /*     lhs variables: [NCONT,U,KSTAIR,V,A,B]=ab01od(A,B)   */
00073         LhsVar(1)=NCONT; LhsVar(2)=U; 
00074         LhsVar(3)=KSTAIR; LhsVar(4)=V; 
00075         LhsVar(5)=A; LhsVar(6)=B;
00076         return 0;
00077 }
00078 /*-----------------------------------------------------------------------------------*/
sci_contr.c
