- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
I am trying to use a standard IMSL routine LINDS to invert a positive definite matrix,.
I installed the IMSL library based on these instructions.
https://software.intel.com/en-us/articles/installing-and-using-the-imsl-libraries
Then I included the lines
INCLUDE 'link_fnl_shared.h'
use linds_int
in my program, before setting up the matrix to be inverted
The program compiles without error, but when I debug it it, it exits without inverting the matrix. It does not trigger a break point, but it just stops.
Any help will be greatly appreciated. Thank you very much
Link Copied
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
This would likely happen if you didn't set PATH to cover all dll you need at run time. If your IMSL is the one bought as an option with ifort, that should happen if you open the ifort cmd window. dependencywalker may be informative.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
No, I am using Intel Visual Fortran for Windows. Not sure where to look..
Thanks!
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
When you run a Fortran program in the Visual Studio debugger, you must first set a breakpoint on an executable statement. Otherwise, the program will run to completion and exit. If the program had a console window, that goes away too. My guess is that this is what is happening.
Try this - select Debug > Start without debugging. What happens?
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
I ran the example given in the IMSL documentation for LINDS. I obtained the correct results and the test code ran with no errors.
As Steve has stated, the normal sequence of events is for the program to run in the debugger to completion, i.e., until a STOP or END statement is executed, if no keyboard input is required. If you did run into difficulties when you ran your code, you will need to provide details of the test case, etc.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
When I start without debugging, I get a message
*** TERMINAL ERROR 2 FROM LINDS. The required storage cannot be allocated. The workspace is based on N, where N=1258655008.
The matrix I am inverting is 11 by 11.
Steve Lionel (Intel) wrote:
When you run a Fortran program in the Visual Studio debugger, you must first set a breakpoint on an executable statement. Otherwise, the program will run to completion and exit. If the program had a console window, that goes away too. My guess is that this is what is happening.
Try this - select Debug > Start without debugging. What happens?
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Sounds as if you passed an uninitialized variable. Please show us a small but complete test case.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Here is the program
If I write the whole program in the main program it runs fine. The problem arises when I do it in the subroutine.
Thank you again for your time.
______________
program Test
INCLUDE 'link_fnl_shared.h'
use linds_int
implicit none
integer nf, nt, i, j
parameter (nf=338, nt=11)
real a(nf,nt), apa(nt, nt), apainv(nt, nt), eye(nt, nt)
do i=1,nf
do j=1,nt
a(i,j)=real(i/j)
end do
end do
call getinv(a, nf, nt, apa, apainv)
eye=matmul(apa, apainv)
do i=1, nt
write(*,*) eye(i,:)
pause
end do
end program Test
subroutine getinv(a, nf, nt, apa, apainv)
external linds
integer nf, nt
real a(nf, nt), ap(nt, nf), apa(nt,nt), apainv(nt, nt)
ap=transpose(a)
apa=matmul(ap, a)
call linds (apa, apainv)
end
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
You have some rather basic errors in your code, and it is to your benefit to read Fortran books and understand why your program did not work.
You call a subroutine (LINDS) that requires an explicit interface, but you did not provide one. The "use linds_int" is not useful in a routine from which LINDS is not called. On the other hand, it must be included in any routine that calls LINDS.
Move the USE LINDS_INT statement from the main program to the subroutine GETINV, remove the EXTERNAL statement, and your program will work.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Thank you for your help. It does indeed work.
Re. Fortran books, do you have a particular recommendation?
Finally, I find that the computation of the inverse is accurate when everything is declared as double precision (real*8), rather than real. When I declare as real, the inverse multiplied by the matrix was off the identity matrix (diagonal was 1.007 for example, off diagonals were also non zero from the third decimal place ).
In comparison, the Numerical recipes routines ludcmp and lubksb performed well for real variables and I did not have this inaccuracy. .
Is that something one needs to be careful about in all IMSL routines or just this particular one?
Thank you again for your time and patience in answering my rather elementary question.
Best wishes
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Fortran Books:
Lawrences - Compaq Visual Fortran is not bad - a bit dated
The Metcalf and Reid books are all good - a bit to academic at times.
If you want to learn to program well - try Abelson and Sussman or Winston and Horne - all LISP - but if you LISP well - you are a long way forward.
The old MICROSOFT Fortran manuals are wonderful if you can get them.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
In comparison, the Numerical recipes routines ludcmp and lubksb performed well for real variables and I did not have this inaccuracy. .
-----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
I have been playng with three inverters, one from Harrison in 1973, one from Conte and deBoore which is mid 70's and I just got through getting PARDISO to run with these problems.
All gave the inversions to an accuracy of at least six figures and agreed.
!=============================================================================== ! Copyright 2004-2016 Intel Corporation All Rights Reserved. ! ! The source code, information and material ("Material") contained herein is ! owned by Intel Corporation or its suppliers or licensors, and title to such ! Material remains with Intel Corporation or its suppliers or licensors. The ! Material contains proprietary information of Intel or its suppliers and ! licensors. The Material is protected by worldwide copyright laws and treaty ! provisions. No part of the Material may be used, copied, reproduced, ! modified, published, uploaded, posted, transmitted, distributed or disclosed ! in any way without Intel's prior express written permission. No license under ! any patent, copyright or other intellectual property rights in the Material ! is granted to or conferred upon you, either expressly, by implication, ! inducement, estoppel or otherwise. Any license under such intellectual ! property rights must be express and approved by Intel in writing. ! ! Unless otherwise agreed by Intel in writing, you may not remove or alter this ! notice or any other notice embedded in Materials by Intel or Intel's ! suppliers or licensors in any way. !=============================================================================== !=============================================================================== ! Copyright 2004-2016 Intel Corporation All Rights Reserved. ! ! The source code, information and material ("Material") contained herein is ! owned by Intel Corporation or its suppliers or licensors, and title to such ! Material remains with Intel Corporation or its suppliers or licensors. The ! Material contains proprietary information of Intel or its suppliers and ! licensors. The Material is protected by worldwide copyright laws and treaty ! provisions. No part of the Material may be used, copied, reproduced, ! modified, published, uploaded, posted, transmitted, distributed or disclosed ! in any way without Intel's prior express written permission. No license under ! any patent, copyright or other intellectual property rights in the Material ! is granted to or conferred upon you, either expressly, by implication, ! inducement, estoppel or otherwise. Any license under such intellectual ! property rights must be express and approved by Intel in writing. ! ! Unless otherwise agreed by Intel in writing, you may not remove or alter this ! notice or any other notice embedded in Materials by Intel or Intel's ! suppliers or licensors in any way. !=============================================================================== * Content : MKL PARDISO Fortran example * ******************************************************************************** C---------------------------------------------------------------------- C Example program to show the use of the "PARDISO" routine C for symmetric linear systems C--------------------------------------------------------------------- C This program can be downloaded from the following site: C www.pardiso-project.org C C (C) Olaf Schenk, Department of Computer Science, C University of Basel, Switzerland. C Email: olaf.schenk@unibas.ch C C--------------------------------------------------------------------- MODULE Solver INTEGER, PARAMETER :: dp = selected_real_kind(15, 307) REAL (KIND=dp) :: g = 9.806, pi = 3.14159265D0 contains subroutine linsolve(n, error, nnz, a, ja, ia, b) implicit none include 'mkl_pardiso.fi' INTEGER :: n, nnz REAL (KIND=dp) :: x(n) ! Pardiso needs X(N) even if sol is returned in b REAL (KIND=dp) :: b(n) C.. Internal solver memory pointer for 64-bit architectures C.. INTEGER*8 pt(64) C.. Internal solver memory pointer for 32-bit architectures C.. INTEGER*4 pt(64) C.. This is OK in both cases TYPE(MKL_PARDISO_HANDLE) pt(64) C.. All other variables INTEGER maxfct, mnum, mtype, phase, nrhs, msglvl INTEGER iparm(64) INTEGER ia(nnz) INTEGER ja(nnz) INTEGER i, idum(1) REAL (KIND=dp) :: ddum(1) integer error REAL (KIND=dp) a(nnz) C.. Fill all arrays containing matrix data. DATA nrhs /1/, maxfct /1/, mnum /1/ C.. C.. Setup PARDISO control parameter C.. DO i = 1, 64 iparm(i) = 0 END DO iparm(1) = 1 ! no solver default iparm(2) = 2 ! fill-in reordering from METIS iparm(3) = 1 ! numbers of processors iparm(4) = 0 ! no iterative-direct algorithm iparm(5) = 0 ! no user fill-in reducing permutation iparm(6) = 0 ! =0 solution on the first n components of x iparm(7) = 0 ! not in use iparm(8) = 9 ! numbers of iterative refinement steps iparm(9) = 0 ! not in use iparm(10) = 13 ! perturb the pivot elements with 1E-13 iparm(11) = 1 ! use nonsymmetric permutation and scaling MPS iparm(12) = 0 ! not in use iparm(13) = 1 ! maximum weighted matching algorithm is switched-on (default for non-symmetric) iparm(14) = 0 ! Output: number of perturbed pivots iparm(15) = 0 ! not in use iparm(16) = 0 ! not in use iparm(17) = 0 ! not in use iparm(18) = -1 ! Output: number of nonzeros in the factor LU iparm(19) = -1 ! Output: Mflops for LU factorization iparm(20) = 0 ! Output: Numbers of CG Iterations error = 0 ! initialize error flag msglvl = 1 ! print statistical information mtype = 11 ! real unsymmetric C.. Initialize the internal solver memory pointer. This is only C necessary for the FIRST call of the PARDISO solver. DO i = 1, 64 pt(i)%DUMMY = 0 END DO C.. Reordering and Symbolic Factorization, This step also allocates C all memory that is necessary for the factorization phase = 11 ! only reordering and symbolic factorization CALL pardiso (pt, maxfct, mnum, mtype, phase, n, a, ia, ja, &idum, nrhs, iparm, msglvl, ddum, ddum, error) WRITE(*,*) 'Reordering completed ... ' IF (error .NE. 0) THEN WRITE(*,*) 'The following ERROR was detected: ', error STOP 1 END IF WRITE(*,*) 'Number of nonzeros in factors = ',iparm(18) WRITE(*,*) 'Number of factorization MFLOPS = ',iparm(19) C.. Factorization. phase = 22 ! only factorization CALL pardiso (pt, maxfct, mnum, mtype, phase, n, a, ia, ja, &idum, nrhs, iparm, msglvl, ddum, ddum, error) WRITE(*,*) 'Factorization completed ... ' IF (error .NE. 0) THEN WRITE(*,*) 'The following ERROR was detected: ', error STOP 1 END IF C.. Back substitution and iterative refinement iparm(8) = 2 ! max numbers of iterative refinement steps phase = 33 ! only factorization CALL pardiso (pt, maxfct, mnum, mtype, phase, n, a, ia, ja, &idum, nrhs, iparm, msglvl, b, x, error) WRITE(*,*) 'Solve completed ... ' WRITE(*,*) 'The solution of the system is ' DO i = 1, n WRITE(*,*) ' x(',i,') = ', x(i) END DO C.. Termination and release of memory phase = -1 ! release internal memory CALL pardiso (pt, maxfct, mnum, mtype, phase, n, ddum, idum, &idum, idum, nrhs, iparm, msglvl, ddum, ddum, error) Return End subroutine linsolve end module
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
The Lawrence book is NOT about Fortran - it's how to create Windows API applications in Compaq (and Intel) Visual Fortran.
I mentioned some books in https://software.intel.com/en-us/blogs/2013/12/30/doctor-fortran-in-its-a-modern-fortran-world
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
The driver routine for PARDISO to take a square matrix and put it into the correct form is shown - it is not pretty - merely allows me to test LINSOLVE before I use it for real.
ludcmp and lubksb -- PARDISO is a lot faster -- a.txt is the sample from the INTEL examples.
Lawrence is a good book --
!------------------------------------------------------------------------------------------------------------------------ ! ! Test Problem ! !------------------------------------------------------------------------------------------------------------------------ ! ja(1) = 1 ! ja(2) = 2 ! ja(3) = 4 ! ja(4) = 1 ! ja(5) = 2 ! ja(6) = 3 ! ja(7) = 4 ! ja(8) = 5 ! ja(9) = 1 ! ja(10) = 3 ! ja(11) = 4 ! ja(12) = 2 ! ja(13) = 5 ! ia(1) = 1 ! ia(2) = 4 ! ia(3) = 6 ! ia(4) = 9 ! ia(5) = 12 ! ia(6) = 14 ! a(1) = 1.d0 ! a(2) = -1.d0 ! a(3) = -3.d0 ! a(4) = -2.d0 ! a(5) = 5.d0 ! a(6) = 4.d0 ! a(7) = 6.d0 ! a(8) = 4.d0 ! a(9) = -4.d0 ! a(10)= 2.d0 ! a(11) = -7.d0 ! a(12) = 8.d0 ! a(13) = -5.d0 !DATA a ! 1 / 1.d0,-1.d0, -3.d0, ! 2 -2.d0, 5.d0, ! 3 4.d0, 6.d0, 4.d0, ! 4 -4.d0, 2.d0, 7.d0, ! 5 8.d0, -5.d0/ PROGRAM pardiso_unsym use Solver IMPLICIT NONE include 'mkl_pardiso.fi' REAL (KIND=dp), ALLOCATABLE :: a(:), b(:), c(:,:), atemp(:) integer, allocatable :: ja(:), ia(:), iatemp(:), jatemp(:) INTEGER :: nnodes = 5 integer error, iaN integer :: nnz = 200, istat, i, j , count integer :: size = 0 integer :: size1 = 0 integer :: size2 = 0 integer :: flag = 0 REAL (KIND=dp) :: delta = 0.000000000000001 open(1, file="a.txt", STATUS = 'old') read(1,*)count nnodes = count iaN = count+1 ALLOCATE (a(nnz), ja(nnz), ia(iaN), jatemp(nnz), iatemp(iaN), b(count), c(count,count), atemp(nnz), STAT=istat) IF (istat.NE.0) THEN WRITE (*, *) '*** Could not allocate some arrays in LINSOLVE' STOP END IF b = 1.0D0 a = 0.0d0 c = 0.0d0 iatemp = 0 jatemp = 0 atemp = 0 ia = 0 ja = 0 do 200 i = 1, count flag = 0 do 300 j = 1, count read(1,*)c(i,j) if(abs(c(i,j)) .gt. delta) then write(*,*)c(i,j) size = size + 1 jatemp(size) = j iatemp(size1 + 1) = size+1 atemp(size) = c(i,j) if(flag .eq. 0) then size1 = size1 + 1 iatemp(size1) = size write(*,*)size1, iatemp(size1) flag = 1 endif endif 300 end do ! write(*,400)(c(i,j),j=1,count) 400 format(5(1x,F10.3)) 200 end do do 500 i = 1, count read(1,*)b(i) 500 end do write(*,100)count, size 100 Format(' Matrix Count :: ',i5, ' Size :: ',i5) call LinSolve(nnodes, nnz, error, atemp, jatemp, iatemp, b) end
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
for symmetric linear systems
Steve:
There is a small error in the example for the non-symmetric matrix in the INTEL examples - it says symmetric in the FORTRAN text file.
John
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
John, please mention the example error in the MKL forum.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Sangeeta P. wrote:
I find that the computation of the inverse is accurate when everything is declared as double precision (real*8), rather than real. When I declare as real, the inverse multiplied by the matrix was off the identity matrix (diagonal was 1.007 for example, off diagonals were also non zero from the third decimal place ).
In comparison, the Numerical recipes routines ludcmp and lubksb performed well for real variables and I did not have this inaccuracy. .
Is that something one needs to be careful about in all IMSL routines or just this particular one?
The finding (regarding results with real*8 versus real) is somewhat unexpected. Please present the specific test case in detail. Most of the time, such discrepancies are caused by passing incorrect arguments to routines such as LINDS.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
mecej4 wrote:
Quote:
Sangeeta P. wrote:
I find that the computation of the inverse is accurate when everything is declared as double precision (real*8), rather than real. When I declare as real, the inverse multiplied by the matrix was off the identity matrix (diagonal was 1.007 for example, off diagonals were also non zero from the third decimal place ).
In comparison, the Numerical recipes routines ludcmp and lubksb performed well for real variables and I did not have this inaccuracy. .
Is that something one needs to be careful about in all IMSL routines or just this particular one?
The finding (regarding results with real*8 versus real) is somewhat unexpected. Please present the specific test case in detail. Most of the time, such discrepancies are caused by passing incorrect arguments to routines such as LINDS.
Given the pure mathematics behind the different inversion routines and the long history of their use from the 1950's I would not believe the original statement -- I agree with mecej4, one needs to present a solid working example to make such an assertion and then the error can be located. I doubt it is an error in the supplied rotuines.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
There is *always* some round-off error when inverting a matrix using a finite-resolution algorithm. If it is ill-conditioned then the error will show up under a smaller resolution than when it is not. The observation that round-off is serious under REAL*4 than REAL*8 tells me that the reason might be that the subject matrix is ill conditioned, rather than an algorithm error.
Any good text on numerical methods should discuss the conditioning concept.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Dear all
Thanks for the comments and suggestions, especially on good reference books.
Here is my program using comparing LINDS and the routines LUDCMP& LUBKSB.
They are not terribly different, but the latter is marginally better. Using REAL and REAL*8 creates a substantial difference. To get the latter, please just change the declarations on the main program and the LUDCMP& LUBKSB subroutines.
Thanks again everyone, for engaging with my problem.
Sangeeta
______________________________________________
program Test
INCLUDE 'link_fnl_shared.h'
use linds_int
implicit none
integer nf, nt, i, j
parameter (nf=338, nt=3)
integer indx(nt)
real a(nf,nt), ap(nt, nf), apa(nt, nt), apainv(nt, nt), eye(nt, nt), b(nf,nt), bp(nt, nf), bpb(nt, nt), bpbinv(nt, nt), eyenew(nt,nt), d, perfecteye(nt, nt), norm_linds, norm_nr
do i=1,nf
do j=1,nt
a(i,j)=real(i/j)
b(i,j)=a(i,j)
end do
end do
ap=transpose(a)
apa=matmul(ap, a)
!Constructing inverse with LINDS
call linds (apa, apainv)
eye=matmul(apa, apainv)
!Constructing inverse with LUDCMP & LUBKSB
bp=transpose(b)
bpb=matmul(bp, b)
!Initializing inverse and creating the accurate identity matrix perfecteye
do i=1, nt
do j=1,nt
bpbinv(i,j)=0.d0
perfecteye(i,j)=0.d0
end do
bpbinv(i,i)=1.d0
perfecteye(i,i)=1.d0
end do
call ludcmp(bpb, nt, nt, indx, d)
do j=1, nt
call lubksb(bpb, nt, nt, indx, bpbinv(1, j))
end do
eyenew=matmul(apa, bpbinv)
!Comparing the distance between perfecteye and the identity matrix got by multiplying the matrix and the inverse
norm_linds=sum(abs(eye-perfecteye))
norm_nr=sum(abs(eyenew-perfecteye))
write(*,*) norm_linds, norm_nr
pause
end program Test
SUBROUTINE ludcmp(a,n,np,indx,d)
INTEGER n,np,indx(n),NMAX
REAL d,a(np,np),TINY
PARAMETER (NMAX=500,TINY=1.0e-20)
INTEGER i,imax,j,k
REAL aamax,dum,sum,vv(NMAX)
d=1.
do 12 i=1,n
aamax=0.
do 11 j=1,n
if (abs(a(i,j)).gt.aamax) aamax=abs(a(i,j))
11 continue
if (aamax.eq.0.) pause 'singular matrix in ludcmp'
vv(i)=1./aamax
12 continue
do 19 j=1,n
do 14 i=1,j-1
sum=a(i,j)
do 13 k=1,i-1
sum=sum-a(i,k)*a(k,j)
13 continue
a(i,j)=sum
14 continue
aamax=0.
do 16 i=j,n
sum=a(i,j)
do 15 k=1,j-1
sum=sum-a(i,k)*a(k,j)
15 continue
a(i,j)=sum
dum=vv(i)*abs(sum)
if (dum.ge.aamax) then
imax=i
aamax=dum
endif
16 continue
if (j.ne.imax)then
do 17 k=1,n
dum=a(imax,k)
a(imax,k)=a(j,k)
a(j,k)=dum
17 continue
d=-d
vv(imax)=vv(j)
endif
indx(j)=imax
if(a(j,j).eq.0.)a(j,j)=TINY
if(j.ne.n)then
dum=1./a(j,j)
do 18 i=j+1,n
a(i,j)=a(i,j)*dum
18 continue
endif
19 continue
return
END
SUBROUTINE lubksb(a,n,np,indx,b)
INTEGER n,np,indx(n)
REAL a(np,np),b(n)
INTEGER i,ii,j,ll
REAL sum
ii=0
do 12 i=1,n
ll=indx(i)
sum=b(ll)
b(ll)=b(i)
if (ii.ne.0)then
do 11 j=ii,i-1
sum=sum-a(i,j)*b(j)
11 continue
else if (sum.ne.0.) then
ii=i
endif
b(i)=sum
12 continue
do 14 i=n,1,-1
sum=b(i)
do 13 j=i+1,n
sum=sum-a(i,j)*b(j)
13 continue
b(i)=sum/a(i,i)
14 continue
return
END
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
do i=1,nf
do j=1,nt
a(i,j)=real(i/j)
b(i,j)=a(i,j)
end do
end do
your inversion matrix is so dense and not banded that it is not funny - this is your a matrix run through PARDISO - adjust your problem to something that is realistic and maybe it will work as you want.
1.00000000000000
1 1
2.00000000000000
2 2
1.00000000000000
3.00000000000000
3 4
1.00000000000000
1.00000000000000
4.00000000000000
4 7
2.00000000000000
1.00000000000000
1.00000000000000
5.00000000000000
5 11
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
6.00000000000000
6 16
3.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
7.00000000000000
7 22
3.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
8.00000000000000
8 29
4.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
9.00000000000000
9 37
4.00000000000000
3.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
10.0000000000000
10 46
5.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
11.0000000000000
11 56
5.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
12.0000000000000
12 67
6.00000000000000
4.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
13.0000000000000
13 79
6.00000000000000
4.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
14.0000000000000
14 92
7.00000000000000
4.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
15.0000000000000
15 106
7.00000000000000
5.00000000000000
3.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
16.0000000000000
16 121
8.00000000000000
5.00000000000000
4.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
17.0000000000000
17 137
8.00000000000000
5.00000000000000
4.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
18.0000000000000
18 154
9.00000000000000
6.00000000000000
4.00000000000000
3.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
19.0000000000000
19 172
9.00000000000000
6.00000000000000
4.00000000000000
3.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
20.0000000000000
20 191
10.0000000000000
6.00000000000000
5.00000000000000
4.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
21.0000000000000
21 211
10.0000000000000
7.00000000000000
5.00000000000000
4.00000000000000
3.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
22.0000000000000
22 232
11.0000000000000
7.00000000000000
5.00000000000000
4.00000000000000
3.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
23.0000000000000
23 254
11.0000000000000
7.00000000000000
5.00000000000000
4.00000000000000
3.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
24.0000000000000
24 277
12.0000000000000
8.00000000000000
6.00000000000000
4.00000000000000
4.00000000000000
3.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
25.0000000000000
25 301
12.0000000000000
8.00000000000000
6.00000000000000
5.00000000000000
4.00000000000000
3.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
26.0000000000000
26 326
13.0000000000000
8.00000000000000
6.00000000000000
5.00000000000000
4.00000000000000
3.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
27.0000000000000
27 352
13.0000000000000
9.00000000000000
6.00000000000000
5.00000000000000
4.00000000000000
3.00000000000000
3.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
28.0000000000000
28 379
14.0000000000000
9.00000000000000
7.00000000000000
5.00000000000000
4.00000000000000
4.00000000000000
3.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
29.0000000000000
29 407
14.0000000000000
9.00000000000000
7.00000000000000
5.00000000000000
4.00000000000000
4.00000000000000
3.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
30.0000000000000
30 436
15.0000000000000
10.0000000000000
7.00000000000000
6.00000000000000
5.00000000000000
4.00000000000000
3.00000000000000
3.00000000000000
3.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
2.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
1.00000000000000
Matrix Count :: 30 Size :: 465
=== PARDISO: solving a real nonsymmetric system ===
1-based array indexing is turned ON
PARDISO double precision computation is turned ON
METIS algorithm at reorder step is turned ON
Scaling is turned ON
Matching is turned ON
Summary: ( reordering phase )
================
Times:
======
Time spent in calculations of symmetric matrix portrait (fulladj): 0.000008 s
Time spent in reordering of the initial matrix (reorder) : 0.000149 s
Time spent in symbolic factorization (symbfct) : 0.000755 s
Time spent in data preparations for factorization (parlist) : 0.000001 s
Time spent in allocation of internal data structures (malloc) : 0.001801 s
Time spent in additional calculations : 0.000024 s
Total time spent : 0.002738 s
Statistics:
===========
Parallel Direct Factorization is running on 4 OpenMP
< Linear system Ax = b >
number of equations: 30
number of non-zeros in A: 465
number of non-zeros in A (%): 51.666667
number of right-hand sides: 1
< Factors L and U >
number of columns for each panel: 128
number of independent subgraphs: 0
< Preprocessing with state of the art partitioning metis>
number of supernodes: 1
size of largest supernode: 30
number of non-zeros in L: 900
number of non-zeros in U: 1
number of non-zeros in L+U: 901
Reordering completed ...
Number of nonzeros in factors = 901
Number of factorization MFLOPS = 0
=== PARDISO is running in In-Core mode, because iparam(60)=0 ===
Percentage of computed non-zeros for LL^T factorization
100 %
=== PARDISO: solving a real nonsymmetric system ===
Single-level factorization algorithm is turned ON
Summary: ( factorization phase )
================
Times:
======
Time spent in copying matrix to internal data structure (A to LU): 0.000000 s
Time spent in factorization step (numfct) : 0.001749 s
Time spent in allocation of internal data structures (malloc) : 0.001294 s
Time spent in additional calculations : 0.000004 s
Total time spent : 0.003048 s
Statistics:
===========
Parallel Direct Factorization is running on 4 OpenMP
< Linear system Ax = b >
number of equations: 30
number of non-zeros in A: 465
number of non-zeros in A (%): 51.666667
number of right-hand sides: 1
< Factors L and U >
number of columns for each panel: 128
number of independent subgraphs: 0
< Preprocessing with state of the art partitioning metis>
number of supernodes: 1
size of largest supernode: 30
number of non-zeros in L: 900
number of non-zeros in U: 1
number of non-zeros in L+U: 901
gflop for the numerical factorization: 0.000018
gflop/s for the numerical factorization: 0.010048
Factorization completed ...
=== PARDISO: solving a real nonsymmetric system ===
Summary: ( solution phase )
================
Times:
======
Time spent in direct solver at solve step (solve) : 0.000239 s
Time spent in additional calculations : 0.000730 s
Total time spent : 0.000969 s
Statistics:
===========
Parallel Direct Factorization is running on 4 OpenMP
< Linear system Ax = b >
number of equations: 30
number of non-zeros in A: 465
number of non-zeros in A (%): 51.666667
number of right-hand sides: 1
< Factors L and U >
number of columns for each panel: 128
number of independent subgraphs: 0
< Preprocessing with state of the art partitioning metis>
number of supernodes: 1
size of largest supernode: 30
number of non-zeros in L: 900
number of non-zeros in U: 1
number of non-zeros in L+U: 901
gflop for the numerical factorization: 0.000018
gflop/s for the numerical factorization: 0.010048
Solve completed ...
The solution of the system is
x( 1 ) = 0.000000000000000E+000
x( 2 ) = -21.0000000000000
x( 3 ) = 21.0000000000000
x( 4 ) = 21.0000000000000
x( 5 ) = 1.878240250526247E-032
x( 6 ) = 7.771561172376110E-017
x( 7 ) = -1.436853791652556E-031
x( 8 ) = -21.0000000000000
x( 9 ) = -7.888609052210118E-031
x( 10 ) = 21.0000000000000
x( 11 ) = 8.545993139894292E-031
x( 12 ) = -21.0000000000000
x( 13 ) = 2.253888300631437E-032
x( 14 ) = 21.0000000000000
x( 15 ) = -21.0000000000000
x( 16 ) = 21.0000000000000
x( 17 ) = 0.000000000000000E+000
x( 18 ) = 7.027711745877242E-015
x( 19 ) = 0.000000000000000E+000
x( 20 ) = -21.0000000000000
x( 21 ) = -21.0000000000000
x( 22 ) = 21.0000000000000
x( 23 ) = 0.000000000000000E+000
x( 24 ) = 21.0000000000000
x( 25 ) = 1.972152263052530E-031
x( 26 ) = 21.0000000000000
x( 27 ) = -21.0000000000000
x( 28 ) = -21.0000000000000
x( 29 ) = 0.000000000000000E+000
x( 30 ) = -5.995204332975852E-015
Press any key to continue . . .
- Subscribe to RSS Feed
- Mark Topic as New
- Mark Topic as Read
- Float this Topic for Current User
- Bookmark
- Subscribe
- Printer Friendly Page