1//A simple demo for the MUMPS interface, with the return of the schur complement 2//to run it, You just have to execute the instruction within Scilab 3// exec sparse_example.sce; 4 5 6//*********************** MATRIX INITIALISATION ***********************// 7 n=10; 8 mat=sprand(n,n,.5)+speye(n,n); 9 size_schur=3; 10 11// Right Hand side setting 12 RHS = ones(n,1); 13 14 15//****************** Initialisation of the Scilab MUMPS structure ******************// 16timer(); 17[id]=initmumps(); 18 19//Here Job=-1, the next call will only initialise the C and Fortran structure 20[id]=dmumps(id); 21 22id.RHS=RHS; 23id.VAR_SCHUR = [n-size_schur+1:n]; 24 25//******************** CALL TO MUMPS FOR RESOLUTION ON INTERNAL PROBLEM ************// 26job=6; 27id.JOB=job; 28 29[id]=dmumps(id,mat); 30 31// verification of the solution 32solution=id.SOL; 33norm1=norm(mat(1:n-size_schur,1:n-size_schur)*solution(1:n-size_schur) - ones(n-size_schur,1),'inf'); 34if norm1> 10^(-9) then 35 write(%io(2),'WARNING: solution on internal problem may not be OK'); 36else 37 write(%io(2),'SOLUTION on internal problem ok'); 38end 39 40 41//******************* TRY REDUCED RHS FUNCTIONALITY **************// 42id.JOB=3; 43id.ICNTL(26)=1; 44 45// Forward 46[id]=dmumps(id,mat); 47 48// Solve the problem on the Schur complement 49id.REDRHS=id.SCHUR \ id.REDRHS; 50 51// and reinject it to MUMPS 52id.ICNTL(26)=2; 53[id]=dmumps(id,mat); 54solution=id.SOL; 55norm1=norm(mat*solution-RHS,'inf') 56if norm1> 10^(-9) then 57 write(%io(2),'WARNING: solution on complete problem may not be OK'); 58else 59 write(%io(2),'SOLUTION on complete problem ok'); 60end 61 62 63 64//****************** DESTRUCTION OF THE MUMPS INSTANCE ******************// 65job=-2; 66id.JOB=job; 67[id]=dmumps(id); 68t=timer() 69