1// =============================================================================
2// Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
3// Copyright (C) 2012 - SE - Sylvestre Ledru
4//
5//  This file is distributed under the same license as the Scilab package.
6// =============================================================================
7
8// <-- CLI SHELL MODE -->
9// <-- LINUX ONLY -->
10nx    = 10;
11
12nev   = 3;
13ncv   = 6;
14bmat  = 'I';
15which = 'LM';
16
17// Local Arrays
18
19iparam  = zeros(11, 1);
20ipntr   = zeros(14, 1);
21_select = zeros(ncv, 1);
22d       = zeros(nev + 1, 1) + 0 * %i;
23z       = zeros(nx, nev) + 0* %i;
24resid   = zeros(nx, 1) + 0 * %i;
25v       = zeros(nx, ncv) + 0 * %i;
26workd   = zeros(3 * nx, 1) + 0 * %i;
27workev  = zeros(2 * ncv, 1) + 0 * %i;
28rwork   = zeros(ncv, 1);
29workl   = zeros(3 * ncv * ncv + 5 *ncv, 1) + 0 * %i;
30
31// Build the complex test matrix
32A            = diag(10 * ones(nx,1) + %i * ones(nx,1));
33A(1:$-1,2:$) = A(1:$-1,2:$) + diag(6 * ones(nx - 1,1));
34A(2:$,1:$-1) = A(2:$,1:$-1) + diag(-6 * ones(nx - 1,1));
35
36tol    = 0;
37ido    = 0;
38
39ishfts = 1;
40maxitr = 300;
41mode1  = 1;
42
43iparam(1) = ishfts;
44iparam(3) = maxitr;
45iparam(7) = mode1;
46
47sigma = complex(0);
48info_znaupd = 0;
49// M A I N   L O O P (Reverse communication)
50while(ido <> 99)
51  // Repeatedly call the routine ZNAUPD and take actions indicated by parameter IDO until
52  // either convergence is indicated or maxitr has been exceeded.
53
54  [ido, resid, v, iparam, ipntr, workd, workl, rwork, info_znaupd] = znaupd(ido, bmat, nx, which, nev, tol, resid, ncv, v, iparam, ipntr, workd, workl, rwork, info_znaupd);
55
56  if(info_znaupd < 0)
57    printf('\nError with znaupd, info = %d\n', info_znaupd);
58    printf('Check the documentation of znaupd\n\n');
59  end
60
61  if(ido == -1 | ido == 1)
62    // Perform matrix vector multiplication
63    workd(ipntr(2):ipntr(2) + nx - 1) = A * workd(ipntr(1):ipntr(1) + nx - 1);
64  end
65end
66
67// Post-Process using ZNEUPD.
68
69rvec    = 1;
70howmany = 'A';
71info_zneupd = 0;
72
73[d, z, resid, iparam, ipntr, workd, workl, rwork, info_zneupd] = zneupd(rvec, howmany, _select, d, z, sigma, workev, bmat, nx, which, nev, tol, resid, ncv, v, ...
74                                                                    iparam, ipntr, workd, workl, rwork, info_zneupd);
75
76assert_checkalmostequal(A * z, z * diag(d(1:3)), sqrt(%eps), 1.e-10);
77