* *--- sample program illustrating the use of DGEXPV * implicit none external dgcoov double precision tic, tac, clock *--- matrix data ... *--- BEWARE: these values must match those in dgmatv.f integer n, nz, nmax, nzmax parameter( nmax = 5500, nzmax = 600000 ) integer ia(nzmax), ja(nzmax) double precision a(nzmax) common /RMAT/ a, ia, ja, nz, n *--- arguments variables ... integer m, mmax, lwsp, liwsp parameter( mmax = 50 ) parameter( lwsp = nmax*(mmax+2)+5*(mmax+2)**2+7, liwsp = mmax+2 ) integer iwsp(liwsp) double precision t, tol, anorm double precision v(nmax), w(nmax), wsp(lwsp) integer i, itrace, iflag, io double precision ZERO, ONE parameter( ZERO=0.0d0, ONE=1.0d0 ) intrinsic ABS, DBLE * *--- Executable statements ... *--- load matrix in COOrdinates ... open( UNIT=7, STATUS='old', IOstat=io, FILE='mycoo.dat' ) if ( io.ne.0 ) stop 'Could not access matrix' read( UNIT=7,FMT=* ) n, nz print*,'order :',n,' number of nonzero :',nz if ( n.gt.nmax ) stop 'Please increase nmax.' if ( nz.gt.nzmax ) stop 'Please increase nzmax.' *--- read data... read( UNIT=7,FMT=* ) (ia(i), ja(i), a(i), i = 1,nz) *--- compute the infinite norm of A ... do i = 1,n wsp(i) = ZERO enddo do i = 1,nz wsp(ia(i)) = wsp(ia(i)) + ABS( a(i) ) enddo anorm = wsp(1) do i = 2,n if ( anorm.lt.DBLE(wsp(i)) ) anorm = wsp(i) enddo write(UNIT=*,FMT='(A,E8.2)') '||A||_inf= ',anorm *--- the operand vector v is set to e_1 ... do i = 1,n v(i) = ZERO enddo v(1) = ONE *--- set other input arguments ... t = 1.0d0 tol = 1.0d-7 m = 30 itrace = 0 *--- compute w = exp(t*A)v with DGEXPV ... tic = clock() call DGEXPV( n, m, t,v,w, tol, anorm, . wsp,lwsp, iwsp,liwsp, dgcoov, itrace, iflag ) tac = clock() *--- print 9001,'----------------------------------------------------' print 9001,'DGEXPV has completed:' print 9001,'----------------------------------------------------' print 9001,'w(1:n) =' do i = 1,n print*,w(i) enddo *--- display some statistics if desired ... print 9001,'final report----------------------------------------' print 9002,'runtime = ',tac-tic print 9002,'||A||_inf = ',anorm print 9003,'nz =',nz print 9003,'n =',n print 9003,'m =',m print 9003,'itrace =',itrace print 9003,'iflag =',iflag print 9003,'ibrkflag =',iwsp(6) print 9003,'mbrkdwn =',iwsp(7) print 9003,'nstep =',iwsp(4) print 9003,'nreject =',iwsp(5) print 9003,'nmult =',iwsp(1) print 9003,'nexph =',iwsp(2) print 9003,'nscale =',iwsp(3) print 9002,'tol = ',tol print 9002,'t = ',t print 9002,'tbrkdwn = ',DBLE( wsp(7) ) print 9002,'step_min = ',DBLE( wsp(1) ) print 9002,'step_max = ',DBLE( wsp(2) ) print 9002,'max_round = ',DBLE( wsp(3) ) print 9002,'sum_round = ',DBLE( wsp(4) ) print 9002,'max_error = ',DBLE( wsp(5) ) print 9002,'sum_error = ',DBLE( wsp(6) ) print 9002,'hump = ',DBLE( wsp(9) ) print 9002,'scale-norm= ',DBLE( wsp(10) ) 9001 format(A) 9002 format(A,E8.2) 9003 format(A,I9) END *----------------------------------------------------------------------| *----------------------------------------------------------------------|