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