Merge remote-tracking branch 'origin/master' into YaSp
[scilab.git] / scilab / modules / arnoldi / sci_gateway / c / sci_dneupd.c
1 /*
2  * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
3  * Copyright (C) ????-2008 - INRIA
4  *
5  * This file must be used under the terms of the CeCILL.
6  * This source file is licensed as described in the file COPYING, which
7  * you should have received as part of this distribution.  The terms
8  * are also available at
9  * http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
10  *
11  */
12
13 #include <math.h>
14 #include <string.h>
15 #include "api_scilab.h"
16 #include "core_math.h"
17 #include "gw_arnoldi.h"
18 #include "localization.h"
19 #include "Scierror.h"
20 /*--------------------------------------------------------------------------*/
21 extern int C2F(dneupd)(int *rvec, char *howmny, int *select, double *dr,
22                        double *di, double *z, int *ldz, double *sigmar,
23                        double *sigmai, double *workev, char *bmat, int *n,
24                        char *which, int *nev, double *tol, double *resid,
25                        int *ncv, double *v, int *ldv, int *iparam, int *ipntr,
26                        double *workd, double *workl, int *lworkl, int *info);
27 /*--------------------------------------------------------------------------*/
28 int sci_dneupd(char *fname, void *pvApiCtx)
29 {
30     SciErr sciErr;
31
32     int* piAddrpRVEC    = NULL;
33     int* pRVEC          = NULL;
34     int* piAddrpHOWMANY = NULL;
35     char* pHOWMANY      = NULL;
36     int* piAddrpSELECT  = NULL;
37     int* pSELECT        = NULL;
38     int* piAddrpDr      = NULL;
39     double* pDr         = NULL;
40     int* piAddrpDi      = NULL;
41     double* pDi         = NULL;
42     int* piAddrpZ       = NULL;
43     double* pZ          = NULL;
44     int* piAddrpSIGMAr  = NULL;
45     double* pSIGMAr     = NULL;
46     int* piAddrpSIGMAi  = NULL;
47     double* pSIGMAi     = NULL;
48     int* piAddrpWORKev  = NULL;
49     double* pWORKev     = NULL;
50     int* piAddrpBMAT    = NULL;
51     char* pBMAT         = NULL;
52     int* piAddrpN       = NULL;
53     int* pN             = NULL;
54     int* piAddrpWHICH   = NULL;
55     char* pWHICH        = NULL;
56     int* piAddrpNEV     = NULL;
57     int* pNEV           = NULL;
58     int* piAddrpTOL     = NULL;
59     double* pTOL        = NULL;
60     int* piAddrpRESID   = NULL;
61     double* pRESID      = NULL;
62     int* piAddrpNCV     = NULL;
63     int* pNCV           = NULL;
64     int* piAddrpV       = NULL;
65     double* pV          = NULL;
66     int* piAddrpIPARAM  = NULL;
67     int* pIPARAM        = NULL;
68     int* piAddrpIPNTR   = NULL;
69     int* pIPNTR         = NULL;
70     int* piAddrpWORKD   = NULL;
71     double* pWORKD      = NULL;
72     int* piAddrpWORKL   = NULL;
73     double* pWORKL      = NULL;
74     int* piAddrpINFO    = NULL;
75     int* pINFO          = NULL;
76
77     int mRVEC,     nRVEC;
78     int mHOWMANY,  nHOWMANY;
79     int mSELECT,   nSELECT;
80     int Dr,        mDr,       nDr;
81     int Di,        mDi,       nDi;
82     int Z,         mZ,        nZ;
83     int mSIGMAr,   nSIGMAr;
84     int mSIGMAi,   nSIGMAi;
85     int mWORKev,   nWORKev;
86     int mBMAT,     nBMAT;
87     int mN,        nN;
88     int mWHICH,    nWHICH;
89     int mNEV,      nNEV;
90     int mTOL,      nTOL;
91     int RESID,    mRESID,    nRESID;
92     int mNCV,      nNCV;
93     int V,        mV,        nV;
94     int IPARAM,   mIPARAM,   nIPARAM;
95     int IPNTR,    mIPNTR,    nIPNTR;
96     int WORKD,    mWORKD,    nWORKD;
97     int WORKL,    mWORKL,    nWORKL;
98     int INFO,     mINFO,     nINFO;
99
100     int minlhs = 1, minrhs = 22, maxlhs = 10, maxrhs = 22;
101     int LDZ, LDV, LWORKL;
102     int sizeWORKL = 0;
103
104     CheckInputArgument(pvApiCtx, minrhs, maxrhs);
105     CheckOutputArgument(pvApiCtx, minlhs, maxlhs);
106
107     /*                                                  VARIABLE = NUMBER   */
108     sciErr = getVarAddressFromPosition(pvApiCtx, 1, &piAddrpRVEC);
109     if (sciErr.iErr)
110     {
111         printError(&sciErr, 0);
112         return 1;
113     }
114
115     // Retrieve a matrix of double at position 1.
116     sciErr = getMatrixOfDoubleAsInteger(pvApiCtx, piAddrpRVEC, &mRVEC, &nRVEC, &pRVEC);
117     if (sciErr.iErr)
118     {
119         printError(&sciErr, 0);
120         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 1);
121         return 1;
122     }
123
124     sciErr = getVarAddressFromPosition(pvApiCtx, 3, &piAddrpSELECT);
125     if (sciErr.iErr)
126     {
127         printError(&sciErr, 0);
128         return 1;
129     }
130
131     // Retrieve a matrix of double at position 3.
132     sciErr = getMatrixOfDoubleAsInteger(pvApiCtx, piAddrpSELECT, &mSELECT, &nSELECT, &pSELECT);
133     if (sciErr.iErr)
134     {
135         printError(&sciErr, 0);
136         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 3);
137         return 1;
138     }
139
140     sciErr = getVarAddressFromPosition(pvApiCtx, 4, &piAddrpDr);
141     if (sciErr.iErr)
142     {
143         printError(&sciErr, 0);
144         return 1;
145     }
146
147     // Retrieve a matrix of double at position 4.
148     sciErr = getMatrixOfDouble(pvApiCtx, piAddrpDr, &mDr, &nDr, &pDr);
149     if (sciErr.iErr)
150     {
151         printError(&sciErr, 0);
152         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 4);
153         return 1;
154     }
155
156     Dr =  4;
157     sciErr = getVarAddressFromPosition(pvApiCtx, 5, &piAddrpDi);
158     if (sciErr.iErr)
159     {
160         printError(&sciErr, 0);
161         return 1;
162     }
163
164     // Retrieve a matrix of double at position 5.
165     sciErr = getMatrixOfDouble(pvApiCtx, piAddrpDi, &mDi, &nDi, &pDi);
166     if (sciErr.iErr)
167     {
168         printError(&sciErr, 0);
169         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 5);
170         return 1;
171     }
172
173     Di =  5;
174     sciErr = getVarAddressFromPosition(pvApiCtx, 6, &piAddrpZ);
175     if (sciErr.iErr)
176     {
177         printError(&sciErr, 0);
178         return 1;
179     }
180
181     // Retrieve a matrix of double at position 6.
182     sciErr = getMatrixOfDouble(pvApiCtx, piAddrpZ, &mZ, &nZ, &pZ);
183     if (sciErr.iErr)
184     {
185         printError(&sciErr, 0);
186         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 6);
187         return 1;
188     }
189
190     Z =  6;
191     sciErr = getVarAddressFromPosition(pvApiCtx, 7, &piAddrpSIGMAr);
192     if (sciErr.iErr)
193     {
194         printError(&sciErr, 0);
195         return 1;
196     }
197
198     // Retrieve a matrix of double at position 7.
199     sciErr = getMatrixOfDouble(pvApiCtx, piAddrpSIGMAr, &mSIGMAr, &nSIGMAr, &pSIGMAr);
200     if (sciErr.iErr)
201     {
202         printError(&sciErr, 0);
203         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 7);
204         return 1;
205     }
206
207     sciErr = getVarAddressFromPosition(pvApiCtx, 8, &piAddrpSIGMAi);
208     if (sciErr.iErr)
209     {
210         printError(&sciErr, 0);
211         return 1;
212     }
213
214     // Retrieve a matrix of double at position 8.
215     sciErr = getMatrixOfDouble(pvApiCtx, piAddrpSIGMAi, &mSIGMAi, &nSIGMAi, &pSIGMAi);
216     if (sciErr.iErr)
217     {
218         printError(&sciErr, 0);
219         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 8);
220         return 1;
221     }
222
223     sciErr = getVarAddressFromPosition(pvApiCtx, 9, &piAddrpWORKev);
224     if (sciErr.iErr)
225     {
226         printError(&sciErr, 0);
227         return 1;
228     }
229
230     // Retrieve a matrix of double at position 9.
231     sciErr = getMatrixOfDouble(pvApiCtx, piAddrpWORKev, &mWORKev, &nWORKev, &pWORKev);
232     if (sciErr.iErr)
233     {
234         printError(&sciErr, 0);
235         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 9);
236         return 1;
237     }
238
239     sciErr = getVarAddressFromPosition(pvApiCtx, 11, &piAddrpN);
240     if (sciErr.iErr)
241     {
242         printError(&sciErr, 0);
243         return 1;
244     }
245
246     // Retrieve a matrix of double at position 11.
247     sciErr = getMatrixOfDoubleAsInteger(pvApiCtx, piAddrpN, &mN, &nN, &pN);
248     if (sciErr.iErr)
249     {
250         printError(&sciErr, 0);
251         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 11);
252         return 1;
253     }
254
255     sciErr = getVarAddressFromPosition(pvApiCtx, 13, &piAddrpNEV);
256     if (sciErr.iErr)
257     {
258         printError(&sciErr, 0);
259         return 1;
260     }
261
262     // Retrieve a matrix of double at position 13.
263     sciErr = getMatrixOfDoubleAsInteger(pvApiCtx, piAddrpNEV, &mNEV, &nNEV, &pNEV);
264     if (sciErr.iErr)
265     {
266         printError(&sciErr, 0);
267         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 13);
268         return 1;
269     }
270
271     sciErr = getVarAddressFromPosition(pvApiCtx, 14, &piAddrpTOL);
272     if (sciErr.iErr)
273     {
274         printError(&sciErr, 0);
275         return 1;
276     }
277
278     // Retrieve a matrix of double at position 14.
279     sciErr = getMatrixOfDouble(pvApiCtx, piAddrpTOL, &mTOL, &nTOL, &pTOL);
280     if (sciErr.iErr)
281     {
282         printError(&sciErr, 0);
283         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 14);
284         return 1;
285     }
286
287     sciErr = getVarAddressFromPosition(pvApiCtx, 15, &piAddrpRESID);
288     if (sciErr.iErr)
289     {
290         printError(&sciErr, 0);
291         return 1;
292     }
293
294     // Retrieve a matrix of double at position 15.
295     sciErr = getMatrixOfDouble(pvApiCtx, piAddrpRESID, &mRESID, &nRESID, &pRESID);
296     if (sciErr.iErr)
297     {
298         printError(&sciErr, 0);
299         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 15);
300         return 1;
301     }
302
303     RESID = 15;
304     sciErr = getVarAddressFromPosition(pvApiCtx, 16, &piAddrpNCV);
305     if (sciErr.iErr)
306     {
307         printError(&sciErr, 0);
308         return 1;
309     }
310
311     // Retrieve a matrix of double at position 16.
312     sciErr = getMatrixOfDoubleAsInteger(pvApiCtx, piAddrpNCV, &mNCV, &nNCV, &pNCV);
313     if (sciErr.iErr)
314     {
315         printError(&sciErr, 0);
316         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 16);
317         return 1;
318     }
319
320     sciErr = getVarAddressFromPosition(pvApiCtx, 17, &piAddrpV);
321     if (sciErr.iErr)
322     {
323         printError(&sciErr, 0);
324         return 1;
325     }
326
327     // Retrieve a matrix of double at position 17.
328     sciErr = getMatrixOfDouble(pvApiCtx, piAddrpV, &mV, &nV, &pV);
329     if (sciErr.iErr)
330     {
331         printError(&sciErr, 0);
332         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 17);
333         return 1;
334     }
335
336     V = 17;
337     sciErr = getVarAddressFromPosition(pvApiCtx, 18, &piAddrpIPARAM);
338     if (sciErr.iErr)
339     {
340         printError(&sciErr, 0);
341         return 1;
342     }
343
344     // Retrieve a matrix of double at position 18.
345     sciErr = getMatrixOfDoubleAsInteger(pvApiCtx, piAddrpIPARAM, &mIPARAM, &nIPARAM, &pIPARAM);
346     if (sciErr.iErr)
347     {
348         printError(&sciErr, 0);
349         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 18);
350         return 1;
351     }
352
353     IPARAM = 18;
354     sciErr = getVarAddressFromPosition(pvApiCtx, 19, &piAddrpIPNTR);
355     if (sciErr.iErr)
356     {
357         printError(&sciErr, 0);
358         return 1;
359     }
360
361     // Retrieve a matrix of double at position 19.
362     sciErr = getMatrixOfDoubleAsInteger(pvApiCtx, piAddrpIPNTR, &mIPNTR, &nIPNTR, &pIPNTR);
363     if (sciErr.iErr)
364     {
365         printError(&sciErr, 0);
366         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 19);
367         return 1;
368     }
369
370     IPNTR = 19;
371     sciErr = getVarAddressFromPosition(pvApiCtx, 20, &piAddrpWORKD);
372     if (sciErr.iErr)
373     {
374         printError(&sciErr, 0);
375         return 1;
376     }
377
378     // Retrieve a matrix of double at position 20.
379     sciErr = getMatrixOfDouble(pvApiCtx, piAddrpWORKD, &mWORKD, &nWORKD, &pWORKD);
380     if (sciErr.iErr)
381     {
382         printError(&sciErr, 0);
383         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 20);
384         return 1;
385     }
386
387     WORKD = 20;
388     sciErr = getVarAddressFromPosition(pvApiCtx, 21, &piAddrpWORKL);
389     if (sciErr.iErr)
390     {
391         printError(&sciErr, 0);
392         return 1;
393     }
394
395     // Retrieve a matrix of double at position 21.
396     sciErr = getMatrixOfDouble(pvApiCtx, piAddrpWORKL, &mWORKL, &nWORKL, &pWORKL);
397     if (sciErr.iErr)
398     {
399         printError(&sciErr, 0);
400         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 21);
401         return 1;
402     }
403
404     WORKL = 21;
405     sciErr = getVarAddressFromPosition(pvApiCtx, 22, &piAddrpINFO);
406     if (sciErr.iErr)
407     {
408         printError(&sciErr, 0);
409         return 1;
410     }
411
412     // Retrieve a matrix of double at position 22.
413     sciErr = getMatrixOfDoubleAsInteger(pvApiCtx, piAddrpINFO, &mINFO, &nINFO, &pINFO);
414     if (sciErr.iErr)
415     {
416         printError(&sciErr, 0);
417         Scierror(202, _("%s: Wrong type for argument %d: A real expected.\n"), fname, 22);
418         return 1;
419     }
420
421     INFO = 22;
422
423     LWORKL = mWORKL * nWORKL;
424     LDV = Max(1, *(int*)(pN));
425     LDZ = LDV;
426
427     /* Check some sizes */
428     if (mIPARAM*nIPARAM != 11)
429     {
430         Scierror(999, _("%s: Wrong size for input argument %s: An array of size %d expected.\n"), fname, "IPARAM", 11);
431         return 0;
432     }
433
434     if (mIPNTR*nIPNTR != 14)
435     {
436         Scierror(999, _("%s: Wrong size for input argument %s: An array of size %d expected.\n"), fname, "IPNTR", 14);
437         return 0;
438     }
439
440     if (mRESID*nRESID != *(int*)(pN))
441     {
442         Scierror(999, _("%s: Wrong size for input argument %s: An array of size %d expected.\n"), fname, "RESID", *(int*)(pN));
443         return 0;
444     }
445
446     if (mWORKD * nWORKD < 3 * *(int*)(pN))
447     {
448         Scierror(999, _("%s: Wrong size for input argument %s: An array of size %d expected.\n"), fname, "WORKD", 3 * *(int*)(pN));
449         return 0;
450     }
451
452     if (mSELECT*nSELECT != *(int*)(pNCV))
453     {
454         Scierror(999, _("%s: Wrong size for input argument %s: An array of size %d expected.\n"), fname, "SELECT", *(int*)(pNCV));
455         return 0;
456     }
457
458     if (mDr*nDr != (*(int*)(pNEV) + 1))
459     {
460         Scierror(999, _("%s: Wrong size for input argument %s: An array of size %d expected.\n"), fname, "Dr", *(int*)(pNEV) + 1);
461         return 0;
462     }
463
464     if (mDi*nDi != (*(int*)(pNEV) + 1))
465     {
466         Scierror(999, _("%s: Wrong size for input argument %s: An array of size %d expected.\n"), fname, "Di", *(int*)(pNEV) + 1);
467         return 0;
468     }
469
470     if ((mZ != *(int*)(pN)) || (nZ != *(int*)(pNEV) + 1))
471     {
472         Scierror(999, _("%s: Wrong size for input argument %s: A matrix of size %dx%d expected.\n"), fname, "Z", *(int*)(pN), *(int*)(pNEV) + 1);
473         return 0;
474     }
475
476     if (mWORKev*nWORKev != 3 * *(int*)(pNCV))
477     {
478         Scierror(999, _("%s: Wrong size for input argument %s: An array of size %d expected.\n"), fname, "WORKev", 3 * *(int*)(pNCV));
479         return 0;
480     }
481
482     if ((mV != *(int*)(pN)) || (nV != *(int*)(pNCV)))
483     {
484         Scierror(999, _("%s: Wrong size for input argument %s: A matrix of size %dx%d expected.\n"), fname, "V", *(int*)(pN), *(int*)(pNCV));
485         return 0;
486     }
487
488     sizeWORKL = 3 * *(int*)(pNCV) * *(int*)(pNCV) + 6 * *(int*)(pNCV);
489
490     if ((mWORKL * nWORKL < sizeWORKL))
491     {
492         Scierror(999, _("%s: Wrong size for input argument %s: An array of size %d expected.\n"), fname, "WORKL", sizeWORKL);
493         return 0;
494     }
495
496     sciErr = getVarAddressFromPosition(pvApiCtx, 2, &piAddrpHOWMANY);
497     if (sciErr.iErr)
498     {
499         printError(&sciErr, 0);
500         return 1;
501     }
502
503     // Retrieve a matrix of double at position 2.
504     if (getAllocatedSingleString(pvApiCtx, piAddrpHOWMANY, &pHOWMANY))
505     {
506         Scierror(202, _("%s: Wrong type for argument #%d: A string expected.\n"), fname, 2);
507         return 1;
508     }
509
510     sciErr = getVarAddressFromPosition(pvApiCtx, 10, &piAddrpBMAT);
511     if (sciErr.iErr)
512     {
513         freeAllocatedSingleString(pHOWMANY);
514         printError(&sciErr, 0);
515         return 1;
516     }
517
518     // Retrieve a matrix of double at position 10.
519     if (getAllocatedSingleString(pvApiCtx, piAddrpBMAT, &pBMAT))
520     {
521         freeAllocatedSingleString(pHOWMANY);
522         Scierror(202, _("%s: Wrong type for argument #%d: A string expected.\n"), fname, 10);
523         return 1;
524     }
525
526     sciErr = getVarAddressFromPosition(pvApiCtx, 12, &piAddrpWHICH);
527     if (sciErr.iErr)
528     {
529         freeAllocatedSingleString(pHOWMANY);
530         freeAllocatedSingleString(pBMAT);
531         printError(&sciErr, 0);
532         return 1;
533     }
534
535     // Retrieve a matrix of double at position 12.
536     if (getAllocatedSingleString(pvApiCtx, piAddrpWHICH, &pWHICH))
537     {
538         freeAllocatedSingleString(pHOWMANY);
539         freeAllocatedSingleString(pBMAT);
540         Scierror(202, _("%s: Wrong type for argument #%d: A string expected.\n"), fname, 12);
541         return 1;
542     }
543
544     C2F(dneupd)(pRVEC, pHOWMANY, pSELECT,
545                 pDr, pDi, pZ, &LDZ,
546                 pSIGMAr, pSIGMAi, pWORKev,
547                 pBMAT, pN, pWHICH,
548                 pNEV, pTOL, pRESID,
549                 pNCV, pV, &LDV,
550                 pIPARAM, pIPNTR,
551                 pWORKD, pWORKL, &LWORKL,
552                 pINFO);
553
554     freeAllocatedSingleString(pHOWMANY);
555     freeAllocatedSingleString(pBMAT);
556     freeAllocatedSingleString(pWHICH);
557
558     if (pINFO[0] < 0)
559     {
560         C2F(errorinfo)("dneupd", (int*)(pINFO), 6L);
561         return 0;
562     }
563
564     AssignOutputVariable(pvApiCtx, 1)  = Dr;
565     AssignOutputVariable(pvApiCtx, 2)  = Di;
566     AssignOutputVariable(pvApiCtx, 3)  = Z;
567     AssignOutputVariable(pvApiCtx, 4)  = RESID;
568     AssignOutputVariable(pvApiCtx, 5)  = V;
569     AssignOutputVariable(pvApiCtx, 6)  = IPARAM;
570     AssignOutputVariable(pvApiCtx, 7)  = IPNTR;
571     AssignOutputVariable(pvApiCtx, 8)  = WORKD;
572     AssignOutputVariable(pvApiCtx, 9)  = WORKL;
573     AssignOutputVariable(pvApiCtx, 10) = INFO;
574
575     ReturnArguments(pvApiCtx);
576
577     return 0;
578 }
579 /*--------------------------------------------------------------------------*/